Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-2981-gc9b2d476df 2025-04-01 01:00:00+00:00
\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}} \newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=} \newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]} \newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
fe_collection.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2003 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
17
20
21#include <deque>
22#include <limits>
23#include <set>
24
25
26
28
29namespace hp
30{
31 template <int dim, int spacedim>
33 {
34 set_default_hierarchy();
35 }
36
37
38
39 template <int dim, int spacedim>
46
47
48
49 template <int dim, int spacedim>
51 const std::vector<const FiniteElement<dim, spacedim> *> &fes)
52 : FECollection()
53 {
54 Assert(fes.size() > 0,
55 ExcMessage("Need to pass at least one finite element."));
56
57 for (unsigned int i = 0; i < fes.size(); ++i)
58 push_back(*fes[i]);
59 }
60
61
62
63 template <int dim, int spacedim>
64 void
66 const FiniteElement<dim, spacedim> &new_fe)
67 {
68 // check that the new element has the right number of components. only check
69 // with the first element, since all the other elements have already passed
70 // the test against the first element
71 Assert(this->empty() ||
72 new_fe.n_components() == this->operator[](0).n_components(),
73 ExcMessage("All elements inside a collection need to have the "
74 "same number of vector components!"));
75
77 }
78
79
80
81 template <int dim, int spacedim>
84 {
85 Assert(this->size() > 0, ExcNoFiniteElements());
86
87 // Since we can only add elements to an FECollection, we are safe comparing
88 // the sizes of this object and the MappingCollection. One circumstance that
89 // might lead to their sizes diverging is this:
90 // - An FECollection is created and then this function is called. The shared
91 // map is now initialized.
92 // - A second FECollection is made as a copy of this one. The shared map is
93 // not recreated.
94 // - The second FECollection is then resized by adding a new FE. The shared
95 // map is thus invalid for the second instance.
96 if (!reference_cell_default_linear_mapping ||
97 reference_cell_default_linear_mapping->size() != this->size())
98 {
99 auto &this_nc = const_cast<FECollection<dim, spacedim> &>(*this);
100
102 std::make_shared<MappingCollection<dim, spacedim>>();
103
104 for (const auto &fe : *this)
105 this_nc.reference_cell_default_linear_mapping->push_back(
106 fe.reference_cell()
107 .template get_default_linear_mapping<dim, spacedim>());
108 }
109
110 return *reference_cell_default_linear_mapping;
112
113
114
115 template <int dim, int spacedim>
116 std::set<unsigned int>
118 const std::set<unsigned int> &fes,
119 const unsigned int codim) const
120 {
121 if constexpr (running_in_debug_mode())
122 {
123 // Validate user inputs.
124 Assert(codim <= dim, ExcImpossibleInDim(dim));
125 Assert(this->size() > 0, ExcEmptyObject());
126 for (const auto &fe : fes)
127 AssertIndexRange(fe, this->size());
128 }
129
130 // Check if any element of this FECollection is able to dominate all
131 // elements of @p fes. If one was found, we add it to the set of
132 // dominating elements.
133 std::set<unsigned int> dominating_fes;
134 for (unsigned int current_fe = 0; current_fe < this->size(); ++current_fe)
135 {
136 // Check if current_fe can dominate all elements in @p fes.
139 for (const auto &other_fe : fes)
140 domination =
141 domination &
142 this->operator[](current_fe)
143 .compare_for_domination(this->operator[](other_fe), codim);
144
145 // If current_fe dominates, add it to the set.
148 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
149 dominating_fes.insert(current_fe);
150 }
151 return dominating_fes;
152 }
153
154
155
156 template <int dim, int spacedim>
157 std::set<unsigned int>
159 const std::set<unsigned int> &fes,
160 const unsigned int codim) const
161 {
162 if constexpr (running_in_debug_mode())
163 {
164 // Validate user inputs.
165 Assert(codim <= dim, ExcImpossibleInDim(dim));
166 Assert(this->size() > 0, ExcEmptyObject());
167 for (const auto &fe : fes)
168 AssertIndexRange(fe, this->size());
169 }
170
171 // Check if any element of this FECollection is dominated by all
172 // elements of @p fes. If one was found, we add it to the set of
173 // dominated elements.
174 std::set<unsigned int> dominated_fes;
175 for (unsigned int current_fe = 0; current_fe < this->size(); ++current_fe)
176 {
177 // Check if current_fe is dominated by all other elements in @p fes.
180 for (const auto &other_fe : fes)
181 domination =
182 domination &
183 this->operator[](current_fe)
184 .compare_for_domination(this->operator[](other_fe), codim);
185
186 // If current_fe is dominated, add it to the set.
189 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
190 dominated_fes.insert(current_fe);
191 }
192 return dominated_fes;
193 }
194
195
196
197 template <int dim, int spacedim>
198 unsigned int
200 const std::set<unsigned int> &fes,
201 const unsigned int codim) const
202 {
203 // If the set of elements contains only a single element,
204 // then this very element is considered to be the dominating one.
205 if (fes.size() == 1)
206 return *fes.begin();
207
208 if constexpr (running_in_debug_mode())
209 {
210 // Validate user inputs.
211 Assert(codim <= dim, ExcImpossibleInDim(dim));
212 Assert(this->size() > 0, ExcEmptyObject());
213 for (const auto &fe : fes)
214 AssertIndexRange(fe, this->size());
215 }
216
217 // There may also be others, in which case we'll check if any of these
218 // elements is able to dominate all others. If one was found, we stop
219 // looking further and return the dominating element.
220 for (const auto &current_fe : fes)
221 {
222 // Check if current_fe can dominate all elements in @p fes.
225 for (const auto &other_fe : fes)
226 if (current_fe != other_fe)
227 domination =
228 domination &
229 this->operator[](current_fe)
230 .compare_for_domination(this->operator[](other_fe), codim);
231
232 // If current_fe dominates, return its index.
235 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
236 return current_fe;
237 }
238
239 // If we couldn't find the dominating object, return an invalid one.
241 }
242
243
244
245 template <int dim, int spacedim>
246 unsigned int
248 const std::set<unsigned int> &fes,
249 const unsigned int codim) const
250 {
251 // If the set of elements contains only a single element,
252 // then this very element is considered to be the dominated one.
253 if (fes.size() == 1)
254 return *fes.begin();
255
256 if constexpr (running_in_debug_mode())
257 {
258 // Validate user inputs.
259 Assert(codim <= dim, ExcImpossibleInDim(dim));
260 Assert(this->size() > 0, ExcEmptyObject());
261 for (const auto &fe : fes)
262 AssertIndexRange(fe, this->size());
263 }
264
265 // There may also be others, in which case we'll check if any of these
266 // elements is dominated by all others. If one was found, we stop
267 // looking further and return the dominated element.
268 for (const auto &current_fe : fes)
269 {
270 // Check if current_fe is dominated by all other elements in @p fes.
273 for (const auto &other_fe : fes)
274 if (current_fe != other_fe)
275 domination =
276 domination &
277 this->operator[](current_fe)
278 .compare_for_domination(this->operator[](other_fe), codim);
279
280 // If current_fe is dominated, return its index.
283 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
284 return current_fe;
285 }
286
287 // If we couldn't find the dominated object, return an invalid one.
290
291
292
293 template <int dim, int spacedim>
294 unsigned int
296 const std::set<unsigned int> &fes,
297 const unsigned int codim) const
298 {
299 unsigned int fe_index = find_dominating_fe(fes, codim);
300
301 if (fe_index == numbers::invalid_fe_index)
302 {
303 const std::set<unsigned int> dominating_fes =
304 find_common_fes(fes, codim);
305 fe_index = find_dominated_fe(dominating_fes, codim);
306 }
307
308 return fe_index;
309 }
310
311
312
313 template <int dim, int spacedim>
314 unsigned int
316 const std::set<unsigned int> &fes,
317 const unsigned int codim) const
318 {
319 unsigned int fe_index = find_dominated_fe(fes, codim);
320
321 if (fe_index == numbers::invalid_fe_index)
322 {
323 const std::set<unsigned int> dominated_fes =
324 find_enclosing_fes(fes, codim);
325 fe_index = find_dominating_fe(dominated_fes, codim);
326 }
327
328 return fe_index;
329 }
330
331
332
333 namespace
334 {
339 std::vector<std::map<unsigned int, unsigned int>>
340 compute_hp_dof_identities(
341 const std::set<unsigned int> &fes,
342 const std::function<std::vector<std::pair<unsigned int, unsigned int>>(
343 const unsigned int,
344 const unsigned int)> &query_identities)
346 // Let's deal with the easy cases first. If the set of fe indices is empty
347 // or has only one entry, then there are no identities:
348 if (fes.size() <= 1)
349 return {};
350
351 // If the set has two entries, then the
352 // FiniteElement::hp_*_dof_identities() function directly returns what we
353 // need. We just need to prefix its output with the respective fe indices:
354 if (fes.size() == 2)
355 {
356 const unsigned int fe_index_1 = *fes.begin();
357 const unsigned int fe_index_2 = *(++fes.begin());
358 const auto reduced_identities =
359 query_identities(fe_index_1, fe_index_2);
360
361 std::vector<std::map<unsigned int, unsigned int>> complete_identities;
362
363 for (const auto &reduced_identity : reduced_identities)
364 {
365 // Each identity returned by query_identities() is a pair of
366 // dof indices. Prefix each with its fe index and put the result
367 // into a vector
368 std::map<unsigned int, unsigned int> complete_identity = {
369 {fe_index_1, reduced_identity.first},
370 {fe_index_2, reduced_identity.second}};
371 complete_identities.emplace_back(std::move(complete_identity));
372 }
373
374 return complete_identities;
375 }
376
377 // Now for the general case of three or more elements:
378 //
379 // Consider all degrees of freedom of the identified elements (represented
380 // via (fe_index,dof_index) pairs) as the nodes in a graph. Then draw
381 // edges for all DoFs that are identified based on what the elements
382 // selected in the argument say. Let us first build this graph, where we
383 // only store the edges of the graph, and as a consequence ignore nodes
384 // (DoFs) that simply don't show up at all in any of the identities:
385 using Node = std::pair<unsigned int, unsigned int>;
386 using Edge = std::pair<Node, Node>;
387 using Graph = std::set<Edge>;
388
389 Graph identities_graph;
390 for (const unsigned int fe_index_1 : fes)
391 for (const unsigned int fe_index_2 : fes)
392 if (fe_index_1 != fe_index_2)
393 for (const auto &identity :
394 query_identities(fe_index_1, fe_index_2))
395 identities_graph.emplace(Node(fe_index_1, identity.first),
396 Node(fe_index_2, identity.second));
397
398 if constexpr (running_in_debug_mode())
399 {
400 // Now verify that indeed the graph is symmetric: If one element
401 // declares that certain ones of its DoFs are to be unified with those
402 // of the other, then the other one should agree with this. As a
403 // consequence of this test succeeding, we know that the graph is
404 // actually undirected.
405 for (const auto &edge : identities_graph)
406 Assert(identities_graph.find({edge.second, edge.first}) !=
407 identities_graph.end(),
409 }
411 // The next step is that we ought to verify that if there is an identity
412 // between (fe1,dof1) and (fe2,dof2), as well as with (fe2,dof2) and
413 // (fe3,dof3), then there should also be an identity between (fe1,dof1)
414 // and (fe3,dof3). The same logic should apply to chains of four
415 // identities.
416 //
417 // This means that the graph we have built above must be composed of a
418 // collection of complete sub-graphs (complete = each possible edge in the
419 // sub-graph exists) -- or, using a different term, that the graph
420 // consists of a number of "cliques". Each of these cliques is then one
421 // extended identity between two or more DoFs, and these are the ones that
422 // we will want to return.
423 //
424 // To ascertain that this is true, and to figure out what we want to
425 // return, we need to divide the graph into its sub-graphs and then ensure
426 // that each sub-graph is indeed a clique. This is slightly cumbersome,
427 // but can be done as follows:
428 // - pick one edge 'e' of G
429 // - add e=(n1,n2) to the sub-graph SG
430 // - set N={n1,n2}
431 // - loop over the remainder of the edges 'e' of the graph:
432 // - if 'e' has one or both nodes in N:
433 // - add 'e' to SG and
434 // - add its two nodes to N (they may already be in there)
435 // - remove 'e' from G
436 //
437 // In general, this may not find the entire sub-graph if the edges are
438 // stored in random order. For example, if the graph consisted of the
439 // following edges in this order:
440 // (a,b)
441 // (c,d)
442 // (a,c)
443 // (a,d)
444 // (b,c)
445 // (b,d)
446 // then the graph itself is a clique, but the algorithm outlined above
447 // would skip the edge (c,d) because neither of the nodes are already
448 // in the set N which at that point is still (a,b).
449 //
450 // But, we store the graph with a std::set, which stored edges in sorted
451 // order where the order is the lexicographic order of nodes. This ensures
452 // that we really capture all edges that correspond to a sub-graph (but
453 // we will assert this as well).
454 //
455 // (For programmatic reasons, we skip the removal of 'e' from G in a first
456 // run through because it modifies the graph and thus invalidates
457 // iterators. But because SG stores all of these edges, we can remove them
458 // all from G after collecting the edges in SG.)
459 std::vector<std::map<unsigned int, unsigned int>> identities;
460 while (identities_graph.size() > 0)
461 {
462 Graph sub_graph; // SG
463 std::set<Node> sub_graph_nodes; // N
464
465 sub_graph.emplace(*identities_graph.begin());
466 sub_graph_nodes.emplace(identities_graph.begin()->first);
467 sub_graph_nodes.emplace(identities_graph.begin()->second);
468
469 for (const Edge &e : identities_graph)
470 if ((sub_graph_nodes.find(e.first) != sub_graph_nodes.end()) ||
471 (sub_graph_nodes.find(e.second) != sub_graph_nodes.end()))
472 {
473 sub_graph.insert(e);
474 sub_graph_nodes.insert(e.first);
475 sub_graph_nodes.insert(e.second);
476 }
477
478 // We have now obtained a sub-graph from the overall graph.
479 // Now remove it from the bigger graph
480 for (const Edge &e : sub_graph)
481 identities_graph.erase(e);
483 if constexpr (running_in_debug_mode())
484 {
485 // There are three checks we ought to perform:
486 // - That the sub-graph is undirected, i.e. that every edge
487 // appears
488 // in both directions
489 for (const auto &edge : sub_graph)
490 Assert(sub_graph.find({edge.second, edge.first}) !=
491 sub_graph.end(),
493
494 // - None of the nodes in the sub-graph should have appeared in
495 // any of the other sub-graphs. If they did, then we have a bug
496 // in extracting sub-graphs. This is actually more easily
497 // checked the other way around: none of the nodes of the
498 // sub-graph we just extracted should be in any of the edges of
499 // the *remaining* graph
500 for (const Node &n : sub_graph_nodes)
501 for (const Edge &e : identities_graph)
502 Assert((n != e.first) && (n != e.second), ExcInternalError());
503 // - Second, the sub-graph we just extracted needs to be complete,
504 // i.e.,
505 // be a "clique". We check this by counting how many edges it
506 // has. for 'n' nodes in 'N', we need to have n*(n-1) edges (we
507 // store both directed edges).
508 Assert(sub_graph.size() ==
509 sub_graph_nodes.size() * (sub_graph_nodes.size() - 1),
511 }
512
513 // At this point we're sure that we have extracted a complete
514 // sub-graph ("clique"). The DoFs involved are all identical then, and
515 // we will store this identity so we can return it later.
516 //
517 // The sub-graph is given as a set of Node objects, which is just
518 // a collection of (fe_index,dof_index) pairs. Because each
519 // fe_index can only appear once there, a better data structure
520 // is a std::map from fe_index to dof_index, which can conveniently
521 // be initialized from a range of iterators to pairs:
522 identities.emplace_back(sub_graph_nodes.begin(),
523 sub_graph_nodes.end());
524 Assert(identities.back().size() == sub_graph_nodes.size(),
526 }
527
528 return identities;
529 }
530 } // namespace
531
532
533
534 template <int dim, int spacedim>
535 std::vector<std::map<unsigned int, unsigned int>>
537 const std::set<unsigned int> &fes) const
538 {
539 auto query_vertex_dof_identities = [this](const unsigned int fe_index_1,
540 const unsigned int fe_index_2) {
541 return (*this)[fe_index_1].hp_vertex_dof_identities((*this)[fe_index_2]);
542 };
543 return compute_hp_dof_identities(fes, query_vertex_dof_identities);
544 }
545
546
547
548 template <int dim, int spacedim>
549 std::vector<std::map<unsigned int, unsigned int>>
551 const std::set<unsigned int> &fes) const
552 {
553 auto query_line_dof_identities = [this](const unsigned int fe_index_1,
554 const unsigned int fe_index_2) {
555 return (*this)[fe_index_1].hp_line_dof_identities((*this)[fe_index_2]);
556 };
557 return compute_hp_dof_identities(fes, query_line_dof_identities);
558 }
559
560
561
562 template <int dim, int spacedim>
563 std::vector<std::map<unsigned int, unsigned int>>
565 const std::set<unsigned int> &fes,
566 const unsigned int face_no) const
567 {
568 auto query_quad_dof_identities = [this,
569 face_no](const unsigned int fe_index_1,
570 const unsigned int fe_index_2) {
571 return (*this)[fe_index_1].hp_quad_dof_identities((*this)[fe_index_2],
572 face_no);
573 };
574 return compute_hp_dof_identities(fes, query_quad_dof_identities);
575 }
576
577
579 template <int dim, int spacedim>
580 void
582 const std::function<
583 unsigned int(const typename hp::FECollection<dim, spacedim> &,
584 const unsigned int)> &next,
585 const std::function<
586 unsigned int(const typename hp::FECollection<dim, spacedim> &,
587 const unsigned int)> &prev)
588 {
589 // copy hierarchy functions
590 hierarchy_next = next;
591 hierarchy_prev = prev;
592 }
593
594
595
596 template <int dim, int spacedim>
597 void
599 {
600 // establish hierarchy corresponding to order of indices
601 set_hierarchy(&DefaultHierarchy::next_index,
602 &DefaultHierarchy::previous_index);
603 }
604
605
606
607 template <int dim, int spacedim>
608 std::vector<unsigned int>
610 const unsigned int fe_index) const
611 {
612 AssertIndexRange(fe_index, this->size());
613
614 std::deque<unsigned int> sequence = {fe_index};
615
616 // get predecessors
617 {
618 unsigned int front = sequence.front();
619 unsigned int previous;
620 while ((previous = previous_in_hierarchy(front)) != front)
621 {
622 sequence.push_front(previous);
623 front = previous;
625 Assert(sequence.size() <= this->size(),
627 "The registered hierarchy is not terminated: "
628 "previous_in_hierarchy() does not stop at a final index."));
629 }
630 }
631
632 // get successors
633 {
634 unsigned int back = sequence.back();
635 unsigned int next;
636 while ((next = next_in_hierarchy(back)) != back)
637 {
638 sequence.push_back(next);
639 back = next;
640
641 Assert(sequence.size() <= this->size(),
643 "The registered hierarchy is not terminated: "
644 "next_in_hierarchy() does not stop at a final index."));
645 }
646 }
647
648 return {sequence.begin(), sequence.end()};
649 }
650
651
652
653 template <int dim, int spacedim>
654 unsigned int
656 const unsigned int fe_index) const
657 {
658 AssertIndexRange(fe_index, this->size());
659
660 const unsigned int new_fe_index = hierarchy_next(*this, fe_index);
661 AssertIndexRange(new_fe_index, this->size());
662
663 return new_fe_index;
664 }
665
666
667
668 template <int dim, int spacedim>
669 unsigned int
671 const unsigned int fe_index) const
672 {
673 AssertIndexRange(fe_index, this->size());
675 const unsigned int new_fe_index = hierarchy_prev(*this, fe_index);
676 AssertIndexRange(new_fe_index, this->size());
677
678 return new_fe_index;
679 }
680
681
682
683 template <int dim, int spacedim>
686 const FEValuesExtractors::Scalar &scalar) const
687 {
688 Assert(this->size() > 0,
689 ExcMessage("This collection contains no finite element."));
690
691 // get the mask from the first element of the collection
692 const ComponentMask mask = (*this)[0].component_mask(scalar);
693
694 // but then also verify that the other elements of the collection
695 // would return the same mask
696 for (unsigned int c = 1; c < this->size(); ++c)
697 Assert(mask == (*this)[c].component_mask(scalar), ExcInternalError());
698
699 return mask;
700 }
701
702
703 template <int dim, int spacedim>
706 const FEValuesExtractors::Vector &vector) const
707 {
708 Assert(this->size() > 0,
709 ExcMessage("This collection contains no finite element."));
710
711 // get the mask from the first element of the collection
712 const ComponentMask mask = (*this)[0].component_mask(vector);
713
714 // but then also verify that the other elements of the collection
715 // would return the same mask
716 for (unsigned int c = 1; c < this->size(); ++c)
717 Assert(mask == (*this)[c].component_mask(vector), ExcInternalError());
719 return mask;
720 }
721
722
723 template <int dim, int spacedim>
726 const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
727 {
728 Assert(this->size() > 0,
729 ExcMessage("This collection contains no finite element."));
730
731 // get the mask from the first element of the collection
732 const ComponentMask mask = (*this)[0].component_mask(sym_tensor);
733
734 // but then also verify that the other elements of the collection
735 // would return the same mask
736 for (unsigned int c = 1; c < this->size(); ++c)
737 Assert(mask == (*this)[c].component_mask(sym_tensor), ExcInternalError());
738
739 return mask;
740 }
741
742
743 template <int dim, int spacedim>
746 {
747 Assert(this->size() > 0,
748 ExcMessage("This collection contains no finite element."));
749
750 // get the mask from the first element of the collection
751 const ComponentMask mask = (*this)[0].component_mask(block_mask);
752
753 // but then also verify that the other elements of the collection
754 // would return the same mask
755 for (unsigned int c = 1; c < this->size(); ++c)
756 Assert(mask == (*this)[c].component_mask(block_mask),
757 ExcMessage("Not all elements of this collection agree on what "
758 "the appropriate mask should be."));
759
760 return mask;
761 }
762
763
764 template <int dim, int spacedim>
767 const FEValuesExtractors::Scalar &scalar) const
768 {
769 Assert(this->size() > 0,
770 ExcMessage("This collection contains no finite element."));
771
772 // get the mask from the first element of the collection
773 const BlockMask mask = (*this)[0].block_mask(scalar);
774
775 // but then also verify that the other elements of the collection
776 // would return the same mask
777 for (unsigned int c = 1; c < this->size(); ++c)
778 Assert(mask == (*this)[c].block_mask(scalar),
779 ExcMessage("Not all elements of this collection agree on what "
780 "the appropriate mask should be."));
781
782 return mask;
783 }
784
785
786 template <int dim, int spacedim>
789 const FEValuesExtractors::Vector &vector) const
790 {
791 Assert(this->size() > 0,
792 ExcMessage("This collection contains no finite element."));
793
794 // get the mask from the first element of the collection
795 const BlockMask mask = (*this)[0].block_mask(vector);
796
797 // but then also verify that the other elements of the collection
798 // would return the same mask
799 for (unsigned int c = 1; c < this->size(); ++c)
800 Assert(mask == (*this)[c].block_mask(vector),
801 ExcMessage("Not all elements of this collection agree on what "
802 "the appropriate mask should be."));
803
804 return mask;
805 }
806
807
808 template <int dim, int spacedim>
811 const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
812 {
813 Assert(this->size() > 0,
814 ExcMessage("This collection contains no finite element."));
815
816 // get the mask from the first element of the collection
817 const BlockMask mask = (*this)[0].block_mask(sym_tensor);
818
819 // but then also verify that the other elements of the collection
820 // would return the same mask
821 for (unsigned int c = 1; c < this->size(); ++c)
822 Assert(mask == (*this)[c].block_mask(sym_tensor),
823 ExcMessage("Not all elements of this collection agree on what "
824 "the appropriate mask should be."));
825
826 return mask;
827 }
828
830
831 template <int dim, int spacedim>
834 const ComponentMask &component_mask) const
835 {
836 Assert(this->size() > 0,
837 ExcMessage("This collection contains no finite element."));
838
839 // get the mask from the first element of the collection
840 const BlockMask mask = (*this)[0].block_mask(component_mask);
841
842 // but then also verify that the other elements of the collection
843 // would return the same mask
844 for (unsigned int c = 1; c < this->size(); ++c)
845 Assert(mask == (*this)[c].block_mask(component_mask),
846 ExcMessage("Not all elements of this collection agree on what "
847 "the appropriate mask should be."));
848
849 return mask;
850 }
851
852
853
854 template <int dim, int spacedim>
855 unsigned int
857 {
858 Assert(this->size() > 0, ExcNoFiniteElements());
859
860 const unsigned int nb = this->operator[](0).n_blocks();
861 for (unsigned int i = 1; i < this->size(); ++i)
862 Assert(this->operator[](i).n_blocks() == nb,
863 ExcMessage("Not all finite elements in this collection have "
864 "the same number of components."));
865
866 return nb;
867 }
868} // namespace hp
869
870
871
872// explicit instantiations
873#include "hp/fe_collection.inst"
874
875
unsigned int n_components() const
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const =0
std::vector< std::map< unsigned int, unsigned int > > hp_vertex_dof_identities(const std::set< unsigned int > &fes) const
unsigned int previous_in_hierarchy(const unsigned int fe_index) const
std::vector< unsigned int > get_hierarchy_sequence(const unsigned int fe_index=0) const
unsigned int find_dominating_fe_extended(const std::set< unsigned int > &fes, const unsigned int codim=0) const
const MappingCollection< dim, spacedim > & get_reference_cell_default_linear_mapping() const
std::set< unsigned int > find_common_fes(const std::set< unsigned int > &fes, const unsigned int codim=0) const
void push_back(const FiniteElement< dim, spacedim > &new_fe)
unsigned int next_in_hierarchy(const unsigned int fe_index) const
std::shared_ptr< MappingCollection< dim, spacedim > > reference_cell_default_linear_mapping
unsigned int find_dominating_fe(const std::set< unsigned int > &fes, const unsigned int codim=0) const
std::set< unsigned int > find_enclosing_fes(const std::set< unsigned int > &fes, const unsigned int codim=0) const
unsigned int find_dominated_fe(const std::set< unsigned int > &fes, const unsigned int codim=0) const
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
std::vector< std::map< unsigned int, unsigned int > > hp_line_dof_identities(const std::set< unsigned int > &fes) const
void set_hierarchy(const std::function< unsigned int(const typename hp::FECollection< dim, spacedim > &, const unsigned int)> &next, const std::function< unsigned int(const typename hp::FECollection< dim, spacedim > &, const unsigned int)> &prev)
unsigned int find_dominated_fe_extended(const std::set< unsigned int > &fes, const unsigned int codim=0) const
unsigned int n_blocks() const
std::vector< std::map< unsigned int, unsigned int > > hp_quad_dof_identities(const std::set< unsigned int > &fes, const unsigned int face_no=0) const
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
constexpr bool running_in_debug_mode()
Definition config.h:78
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
Point< 2 > second
Definition grid_out.cc:4633
Point< 2 > first
Definition grid_out.cc:4632
static ::ExceptionBase & ExcEmptyObject()
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
std::size_t size
Definition mpi.cc:745
Definition hp.h:117
constexpr types::fe_index invalid_fe_index
Definition types.h:254
void prev(std::tuple< I1, I2 > &t)