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