Reference documentation for deal.II version 9.5.0
\(\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// Copyright (C) 2003 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16
18
21
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>
42 : FECollection()
43 {
44 push_back(fe);
45 }
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->size() == 0 ||
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;
111 }
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#ifdef DEBUG
122 // Validate user inputs.
123 Assert(codim <= dim, ExcImpossibleInDim(dim));
124 Assert(this->size() > 0, ExcEmptyObject());
125 for (const auto &fe : fes)
126 AssertIndexRange(fe, this->size());
127#endif
128
129 // Check if any element of this FECollection is able to dominate all
130 // elements of @p fes. If one was found, we add it to the set of
131 // dominating elements.
132 std::set<unsigned int> dominating_fes;
133 for (unsigned int current_fe = 0; current_fe < this->size(); ++current_fe)
134 {
135 // Check if current_fe can dominate all elements in @p fes.
138 for (const auto &other_fe : fes)
139 domination =
140 domination &
141 this->operator[](current_fe)
142 .compare_for_domination(this->operator[](other_fe), codim);
143
144 // If current_fe dominates, add it to the set.
147 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
148 dominating_fes.insert(current_fe);
149 }
150 return dominating_fes;
151 }
152
153
154
155 template <int dim, int spacedim>
156 std::set<unsigned int>
158 const std::set<unsigned int> &fes,
159 const unsigned int codim) const
160 {
161#ifdef DEBUG
162 // Validate user inputs.
163 Assert(codim <= dim, ExcImpossibleInDim(dim));
164 Assert(this->size() > 0, ExcEmptyObject());
165 for (const auto &fe : fes)
166 AssertIndexRange(fe, this->size());
167#endif
168
169 // Check if any element of this FECollection is dominated by all
170 // elements of @p fes. If one was found, we add it to the set of
171 // dominated elements.
172 std::set<unsigned int> dominated_fes;
173 for (unsigned int current_fe = 0; current_fe < this->size(); ++current_fe)
174 {
175 // Check if current_fe is dominated by all other elements in @p fes.
178 for (const auto &other_fe : fes)
179 domination =
180 domination &
181 this->operator[](current_fe)
182 .compare_for_domination(this->operator[](other_fe), codim);
183
184 // If current_fe is dominated, add it to the set.
187 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
188 dominated_fes.insert(current_fe);
190 return dominated_fes;
191 }
192
193
194
195 template <int dim, int spacedim>
196 unsigned int
198 const std::set<unsigned int> &fes,
199 const unsigned int codim) const
200 {
201 // If the set of elements contains only a single element,
202 // then this very element is considered to be the dominating one.
203 if (fes.size() == 1)
204 return *fes.begin();
205
206#ifdef DEBUG
207 // Validate user inputs.
208 Assert(codim <= dim, ExcImpossibleInDim(dim));
209 Assert(this->size() > 0, ExcEmptyObject());
210 for (const auto &fe : fes)
211 AssertIndexRange(fe, this->size());
212#endif
213
214 // There may also be others, in which case we'll check if any of these
215 // elements is able to dominate all others. If one was found, we stop
216 // looking further and return the dominating element.
217 for (const auto &current_fe : fes)
218 {
219 // Check if current_fe can dominate all elements in @p fes.
222 for (const auto &other_fe : fes)
223 if (current_fe != other_fe)
224 domination =
225 domination &
226 this->operator[](current_fe)
227 .compare_for_domination(this->operator[](other_fe), codim);
229 // If current_fe dominates, return its index.
232 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
233 return current_fe;
234 }
235
236 // If we couldn't find the dominating object, return an invalid one.
238 }
239
240
241
242 template <int dim, int spacedim>
243 unsigned int
245 const std::set<unsigned int> &fes,
246 const unsigned int codim) const
247 {
248 // If the set of elements contains only a single element,
249 // then this very element is considered to be the dominated one.
250 if (fes.size() == 1)
251 return *fes.begin();
252
253#ifdef DEBUG
254 // Validate user inputs.
255 Assert(codim <= dim, ExcImpossibleInDim(dim));
256 Assert(this->size() > 0, ExcEmptyObject());
257 for (const auto &fe : fes)
258 AssertIndexRange(fe, this->size());
259#endif
260
261 // There may also be others, in which case we'll check if any of these
262 // elements is dominated by all others. If one was found, we stop
263 // looking further and return the dominated element.
264 for (const auto &current_fe : fes)
265 {
266 // Check if current_fe is dominated by all other elements in @p fes.
269 for (const auto &other_fe : fes)
270 if (current_fe != other_fe)
271 domination =
272 domination &
273 this->operator[](current_fe)
274 .compare_for_domination(this->operator[](other_fe), codim);
275
276 // If current_fe is dominated, return its index.
279 /*covers cases like {Q2,Q3,Q1,Q1} with fes={2,3}*/))
280 return current_fe;
281 }
282
283 // If we couldn't find the dominated object, return an invalid one.
285 }
286
287
289 template <int dim, int spacedim>
290 unsigned int
292 const std::set<unsigned int> &fes,
293 const unsigned int codim) const
294 {
295 unsigned int fe_index = find_dominating_fe(fes, codim);
296
297 if (fe_index == numbers::invalid_fe_index)
298 {
299 const std::set<unsigned int> dominating_fes =
300 find_common_fes(fes, codim);
301 fe_index = find_dominated_fe(dominating_fes, codim);
302 }
303
304 return fe_index;
305 }
306
307
308
309 template <int dim, int spacedim>
310 unsigned int
312 const std::set<unsigned int> &fes,
313 const unsigned int codim) const
314 {
315 unsigned int fe_index = find_dominated_fe(fes, codim);
316
317 if (fe_index == numbers::invalid_fe_index)
318 {
319 const std::set<unsigned int> dominated_fes =
320 find_enclosing_fes(fes, codim);
321 fe_index = find_dominating_fe(dominated_fes, codim);
322 }
323
324 return fe_index;
325 }
326
327
328
329 namespace
330 {
335 std::vector<std::map<unsigned int, unsigned int>>
336 compute_hp_dof_identities(
337 const std::set<unsigned int> &fes,
338 const std::function<std::vector<std::pair<unsigned int, unsigned int>>(
339 const unsigned int,
340 const unsigned int)> & query_identities)
341 {
342 // Let's deal with the easy cases first. If the set of fe indices is empty
343 // or has only one entry, then there are no identities:
344 if (fes.size() <= 1)
345 return {};
346
347 // If the set has two entries, then the
348 // FiniteElement::hp_*_dof_identities() function directly returns what we
349 // need. We just need to prefix its output with the respective fe indices:
350 if (fes.size() == 2)
352 const unsigned int fe_index_1 = *fes.begin();
353 const unsigned int fe_index_2 = *(++fes.begin());
354 const auto reduced_identities =
355 query_identities(fe_index_1, fe_index_2);
356
357 std::vector<std::map<unsigned int, unsigned int>> complete_identities;
359 for (const auto &reduced_identity : reduced_identities)
360 {
361 // Each identity returned by query_identities() is a pair of
362 // dof indices. Prefix each with its fe index and put the result
363 // into a vector
364 std::map<unsigned int, unsigned int> complete_identity = {
365 {fe_index_1, reduced_identity.first},
366 {fe_index_2, reduced_identity.second}};
367 complete_identities.emplace_back(std::move(complete_identity));
368 }
369
370 return complete_identities;
371 }
372
373 // Now for the general case of three or more elements:
374 //
375 // Consider all degrees of freedom of the identified elements (represented
376 // via (fe_index,dof_index) pairs) as the nodes in a graph. Then draw
377 // edges for all DoFs that are identified based on what the elements
378 // selected in the argument say. Let us first build this graph, where we
379 // only store the edges of the graph, and as a consequence ignore nodes
380 // (DoFs) that simply don't show up at all in any of the identities:
381 using Node = std::pair<unsigned int, unsigned int>;
382 using Edge = std::pair<Node, Node>;
383 using Graph = std::set<Edge>;
385 Graph identities_graph;
386 for (const unsigned int fe_index_1 : fes)
387 for (const unsigned int fe_index_2 : fes)
388 if (fe_index_1 != fe_index_2)
389 for (const auto &identity :
390 query_identities(fe_index_1, fe_index_2))
391 identities_graph.emplace(Node(fe_index_1, identity.first),
392 Node(fe_index_2, identity.second));
393
394#ifdef DEBUG
395 // Now verify that indeed the graph is symmetric: If one element
396 // declares that certain ones of its DoFs are to be unified with those
397 // of the other, then the other one should agree with this. As a
398 // consequence of this test succeeding, we know that the graph is actually
399 // undirected.
400 for (const auto &edge : identities_graph)
401 Assert(identities_graph.find({edge.second, edge.first}) !=
402 identities_graph.end(),
404#endif
405
406 // The next step is that we ought to verify that if there is an identity
407 // between (fe1,dof1) and (fe2,dof2), as well as with (fe2,dof2) and
408 // (fe3,dof3), then there should also be an identity between (fe1,dof1)
409 // and (fe3,dof3). The same logic should apply to chains of four
410 // identities.
411 //
412 // This means that the graph we have built above must be composed of a
413 // collection of complete sub-graphs (complete = each possible edge in the
414 // sub-graph exists) -- or, using a different term, that the graph
415 // consists of a number of "cliques". Each of these cliques is then one
416 // extended identity between two or more DoFs, and these are the ones that
417 // we will want to return.
418 //
419 // To ascertain that this is true, and to figure out what we want to
420 // return, we need to divide the graph into its sub-graphs and then ensure
421 // that each sub-graph is indeed a clique. This is slightly cumbersome,
422 // but can be done as follows:
423 // - pick one edge 'e' of G
424 // - add e=(n1,n2) to the sub-graph SG
425 // - set N={n1,n2}
426 // - loop over the remainder of the edges 'e' of the graph:
427 // - if 'e' has one or both nodes in N:
428 // - add 'e' to SG and
429 // - add its two nodes to N (they may already be in there)
430 // - remove 'e' from G
431 //
432 // In general, this may not find the entire sub-graph if the edges are
433 // stored in random order. For example, if the graph consisted of the
434 // following edges in this order:
435 // (a,b)
436 // (c,d)
437 // (a,c)
438 // (a,d)
439 // (b,c)
440 // (b,d)
441 // then the graph itself is a clique, but the algorithm outlined above
442 // would skip the edge (c,d) because neither of the nodes are already
443 // in the set N which at that point is still (a,b).
444 //
445 // But, we store the graph with a std::set, which stored edges in sorted
446 // order where the order is the lexicographic order of nodes. This ensures
447 // that we really capture all edges that correspond to a sub-graph (but
448 // we will assert this as well).
449 //
450 // (For programmatic reasons, we skip the removal of 'e' from G in a first
451 // run through because it modifies the graph and thus invalidates
452 // iterators. But because SG stores all of these edges, we can remove them
453 // all from G after collecting the edges in SG.)
454 std::vector<std::map<unsigned int, unsigned int>> identities;
455 while (identities_graph.size() > 0)
456 {
457 Graph sub_graph; // SG
458 std::set<Node> sub_graph_nodes; // N
459
460 sub_graph.emplace(*identities_graph.begin());
461 sub_graph_nodes.emplace(identities_graph.begin()->first);
462 sub_graph_nodes.emplace(identities_graph.begin()->second);
463
464 for (const Edge &e : identities_graph)
465 if ((sub_graph_nodes.find(e.first) != sub_graph_nodes.end()) ||
466 (sub_graph_nodes.find(e.second) != sub_graph_nodes.end()))
467 {
468 sub_graph.insert(e);
469 sub_graph_nodes.insert(e.first);
470 sub_graph_nodes.insert(e.second);
471 }
472
473 // We have now obtained a sub-graph from the overall graph.
474 // Now remove it from the bigger graph
475 for (const Edge &e : sub_graph)
476 identities_graph.erase(e);
477
478#ifdef DEBUG
479 // There are three checks we ought to perform:
480 // - That the sub-graph is undirected, i.e. that every edge appears
481 // in both directions
482 for (const auto &edge : sub_graph)
483 Assert(sub_graph.find({edge.second, edge.first}) != sub_graph.end(),
485
486 // - None of the nodes in the sub-graph should have appeared in
487 // any of the other sub-graphs. If they did, then we have a bug
488 // in extracting sub-graphs. This is actually more easily checked
489 // the other way around: none of the nodes of the sub-graph we
490 // just extracted should be in any of the edges of the *remaining*
491 // graph
492 for (const Node &n : sub_graph_nodes)
493 for (const Edge &e : identities_graph)
494 Assert((n != e.first) && (n != e.second), ExcInternalError());
495 // - Second, the sub-graph we just extracted needs to be complete,
496 // i.e.,
497 // be a "clique". We check this by counting how many edges it has.
498 // for 'n' nodes in 'N', we need to have n*(n-1) edges (we store
499 // both directed edges).
500 Assert(sub_graph.size() ==
501 sub_graph_nodes.size() * (sub_graph_nodes.size() - 1),
503#endif
504
505 // At this point we're sure that we have extracted a complete
506 // sub-graph ("clique"). The DoFs involved are all identical then, and
507 // we will store this identity so we can return it later.
508 //
509 // The sub-graph is given as a set of Node objects, which is just
510 // a collection of (fe_index,dof_index) pairs. Because each
511 // fe_index can only appear once there, a better data structure
512 // is a std::map from fe_index to dof_index, which can conveniently
513 // be initialized from a range of iterators to pairs:
514 identities.emplace_back(sub_graph_nodes.begin(),
515 sub_graph_nodes.end());
516 Assert(identities.back().size() == sub_graph_nodes.size(),
518 }
519
520 return identities;
521 }
522 } // namespace
523
524
525
526 template <int dim, int spacedim>
527 std::vector<std::map<unsigned int, unsigned int>>
529 const std::set<unsigned int> &fes) const
530 {
531 auto query_vertex_dof_identities = [this](const unsigned int fe_index_1,
532 const unsigned int fe_index_2) {
533 return (*this)[fe_index_1].hp_vertex_dof_identities((*this)[fe_index_2]);
534 };
535 return compute_hp_dof_identities(fes, query_vertex_dof_identities);
536 }
537
538
539
540 template <int dim, int spacedim>
541 std::vector<std::map<unsigned int, unsigned int>>
543 const std::set<unsigned int> &fes) const
544 {
545 auto query_line_dof_identities = [this](const unsigned int fe_index_1,
546 const unsigned int fe_index_2) {
547 return (*this)[fe_index_1].hp_line_dof_identities((*this)[fe_index_2]);
548 };
549 return compute_hp_dof_identities(fes, query_line_dof_identities);
550 }
551
552
553
554 template <int dim, int spacedim>
555 std::vector<std::map<unsigned int, unsigned int>>
557 const std::set<unsigned int> &fes,
558 const unsigned int face_no) const
559 {
560 auto query_quad_dof_identities = [this,
561 face_no](const unsigned int fe_index_1,
562 const unsigned int fe_index_2) {
563 return (*this)[fe_index_1].hp_quad_dof_identities((*this)[fe_index_2],
564 face_no);
565 };
566 return compute_hp_dof_identities(fes, query_quad_dof_identities);
567 }
568
569
570
571 template <int dim, int spacedim>
572 void
574 const std::function<
575 unsigned int(const typename hp::FECollection<dim, spacedim> &,
576 const unsigned int)> &next,
577 const std::function<
578 unsigned int(const typename hp::FECollection<dim, spacedim> &,
579 const unsigned int)> &prev)
580 {
581 // copy hierarchy functions
582 hierarchy_next = next;
583 hierarchy_prev = prev;
584 }
585
586
587
588 template <int dim, int spacedim>
589 void
591 {
592 // establish hierarchy corresponding to order of indices
593 set_hierarchy(&DefaultHierarchy::next_index,
594 &DefaultHierarchy::previous_index);
595 }
596
597
598
599 template <int dim, int spacedim>
600 std::vector<unsigned int>
602 const unsigned int fe_index) const
603 {
604 AssertIndexRange(fe_index, this->size());
605
606 std::deque<unsigned int> sequence = {fe_index};
607
608 // get predecessors
609 {
610 unsigned int front = sequence.front();
611 unsigned int previous;
612 while ((previous = previous_in_hierarchy(front)) != front)
613 {
614 sequence.push_front(previous);
615 front = previous;
616
617 Assert(sequence.size() <= this->size(),
619 "The registered hierarchy is not terminated: "
620 "previous_in_hierarchy() does not stop at a final index."));
621 }
622 }
624 // get successors
625 {
626 unsigned int back = sequence.back();
627 unsigned int next;
628 while ((next = next_in_hierarchy(back)) != back)
629 {
630 sequence.push_back(next);
631 back = next;
632
633 Assert(sequence.size() <= this->size(),
635 "The registered hierarchy is not terminated: "
636 "next_in_hierarchy() does not stop at a final index."));
637 }
638 }
639
640 return {sequence.begin(), sequence.end()};
641 }
642
643
644
645 template <int dim, int spacedim>
646 unsigned int
648 const unsigned int fe_index) const
649 {
650 AssertIndexRange(fe_index, this->size());
651
652 const unsigned int new_fe_index = hierarchy_next(*this, fe_index);
653 AssertIndexRange(new_fe_index, this->size());
655 return new_fe_index;
656 }
657
658
659
660 template <int dim, int spacedim>
661 unsigned int
663 const unsigned int fe_index) const
664 {
665 AssertIndexRange(fe_index, this->size());
666
667 const unsigned int new_fe_index = hierarchy_prev(*this, fe_index);
668 AssertIndexRange(new_fe_index, this->size());
669
670 return new_fe_index;
671 }
672
674
675 template <int dim, int spacedim>
678 const FEValuesExtractors::Scalar &scalar) const
679 {
680 Assert(this->size() > 0,
681 ExcMessage("This collection contains no finite element."));
682
683 // get the mask from the first element of the collection
684 const ComponentMask mask = (*this)[0].component_mask(scalar);
685
686 // but then also verify that the other elements of the collection
687 // would return the same mask
688 for (unsigned int c = 1; c < this->size(); ++c)
689 Assert(mask == (*this)[c].component_mask(scalar), ExcInternalError());
690
691 return mask;
692 }
694
695 template <int dim, int spacedim>
698 const FEValuesExtractors::Vector &vector) const
699 {
700 Assert(this->size() > 0,
701 ExcMessage("This collection contains no finite element."));
702
703 // get the mask from the first element of the collection
704 const ComponentMask mask = (*this)[0].component_mask(vector);
705
706 // but then also verify that the other elements of the collection
707 // would return the same mask
708 for (unsigned int c = 1; c < this->size(); ++c)
709 Assert(mask == (*this)[c].component_mask(vector), ExcInternalError());
710
711 return mask;
712 }
713
714
715 template <int dim, int spacedim>
718 const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
719 {
720 Assert(this->size() > 0,
721 ExcMessage("This collection contains no finite element."));
722
723 // get the mask from the first element of the collection
724 const ComponentMask mask = (*this)[0].component_mask(sym_tensor);
725
726 // but then also verify that the other elements of the collection
727 // would return the same mask
728 for (unsigned int c = 1; c < this->size(); ++c)
729 Assert(mask == (*this)[c].component_mask(sym_tensor), ExcInternalError());
730
731 return mask;
732 }
733
734
735 template <int dim, int spacedim>
738 {
739 Assert(this->size() > 0,
740 ExcMessage("This collection contains no finite element."));
741
742 // get the mask from the first element of the collection
743 const ComponentMask mask = (*this)[0].component_mask(block_mask);
744
745 // but then also verify that the other elements of the collection
746 // would return the same mask
747 for (unsigned int c = 1; c < this->size(); ++c)
748 Assert(mask == (*this)[c].component_mask(block_mask),
749 ExcMessage("Not all elements of this collection agree on what "
750 "the appropriate mask should be."));
751
752 return mask;
753 }
754
755
756 template <int dim, int spacedim>
759 const FEValuesExtractors::Scalar &scalar) const
760 {
761 Assert(this->size() > 0,
762 ExcMessage("This collection contains no finite element."));
763
764 // get the mask from the first element of the collection
765 const BlockMask mask = (*this)[0].block_mask(scalar);
766
767 // but then also verify that the other elements of the collection
768 // would return the same mask
769 for (unsigned int c = 1; c < this->size(); ++c)
770 Assert(mask == (*this)[c].block_mask(scalar),
771 ExcMessage("Not all elements of this collection agree on what "
772 "the appropriate mask should be."));
773
774 return mask;
775 }
776
777
778 template <int dim, int spacedim>
781 const FEValuesExtractors::Vector &vector) const
782 {
783 Assert(this->size() > 0,
784 ExcMessage("This collection contains no finite element."));
785
786 // get the mask from the first element of the collection
787 const BlockMask mask = (*this)[0].block_mask(vector);
788
789 // but then also verify that the other elements of the collection
790 // would return the same mask
791 for (unsigned int c = 1; c < this->size(); ++c)
792 Assert(mask == (*this)[c].block_mask(vector),
793 ExcMessage("Not all elements of this collection agree on what "
794 "the appropriate mask should be."));
795
796 return mask;
798
799
800 template <int dim, int spacedim>
803 const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
804 {
805 Assert(this->size() > 0,
806 ExcMessage("This collection contains no finite element."));
807
808 // get the mask from the first element of the collection
809 const BlockMask mask = (*this)[0].block_mask(sym_tensor);
810
811 // but then also verify that the other elements of the collection
812 // would return the same mask
813 for (unsigned int c = 1; c < this->size(); ++c)
814 Assert(mask == (*this)[c].block_mask(sym_tensor),
815 ExcMessage("Not all elements of this collection agree on what "
816 "the appropriate mask should be."));
817
818 return mask;
819 }
820
821
822
823 template <int dim, int spacedim>
826 const ComponentMask &component_mask) const
827 {
828 Assert(this->size() > 0,
829 ExcMessage("This collection contains no finite element."));
830
831 // get the mask from the first element of the collection
832 const BlockMask mask = (*this)[0].block_mask(component_mask);
833
834 // but then also verify that the other elements of the collection
835 // would return the same mask
836 for (unsigned int c = 1; c < this->size(); ++c)
837 Assert(mask == (*this)[c].block_mask(component_mask),
838 ExcMessage("Not all elements of this collection agree on what "
839 "the appropriate mask should be."));
840
841 return mask;
842 }
843
844
845
846 template <int dim, int spacedim>
847 unsigned int
849 {
850 Assert(this->size() > 0, ExcNoFiniteElements());
851
852 const unsigned int nb = this->operator[](0).n_blocks();
853 for (unsigned int i = 1; i < this->size(); ++i)
854 Assert(this->operator[](i).n_blocks() == nb,
855 ExcMessage("Not all finite elements in this collection have "
856 "the same number of components."));
857
858 return nb;
859 }
860} // namespace hp
861
862
863
864// explicit instantiations
865#include "fe_collection.inst"
866
867
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 2 > second
Definition grid_out.cc:4616
Point< 2 > first
Definition grid_out.cc:4615
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:118
const types::fe_index invalid_fe_index
Definition types.h:228