From acc0a241e3914056379514261dafd4a37cda6eff Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 26 Aug 2021 15:27:53 -0600 Subject: [PATCH] Change the return type of hp::FECollection::hp_*_dof_identities. --- include/deal.II/hp/fe_collection.h | 16 +++++++++---- source/dofs/dof_handler_policy.cc | 36 ++++++++++++++++++------------ source/hp/fe_collection.cc | 30 ++++++++++++++++--------- 3 files changed, 53 insertions(+), 29 deletions(-) diff --git a/include/deal.II/hp/fe_collection.h b/include/deal.II/hp/fe_collection.h index a5d4ba6326..428552e0f1 100644 --- a/include/deal.II/hp/fe_collection.h +++ b/include/deal.II/hp/fe_collection.h @@ -311,23 +311,31 @@ namespace hp * the vector is a set of pairs `(fe_index,dof_index)` that identifies * the `fe_index` (an element of the `fes` argument to this function) of * an element and the `dof_index` indicates the how-manyth degree of freedom - * of that element on a vertex participates in this identity. + * of that element on a vertex participates in this identity. Now, + * every `fe_index` can appear only once in these sets (for each identity, + * only one degree of freedom of a finite element can be involved -- + * otherwise we would have identities between different DoFs of the same + * element, which would make the element not unisolvent), and as a + * consequence the function does not actually return a set of + * `(fe_index,dof_index)` pairs for each identity, but instead a `std::map` + * from `fe_index` to `dof_index`, which is conceptually of course + * equivalent to a `std::set` of pairs, but in practice is easier to query. */ - std::vector>> + std::vector> hp_vertex_dof_identities(const std::set &fes) const; /** * Same as hp_vertex_dof_indices(), except that the function treats degrees * of freedom on lines. */ - std::vector>> + std::vector> hp_line_dof_identities(const std::set &fes) const; /** * Same as hp_vertex_dof_indices(), except that the function treats degrees * of freedom on quads. */ - std::vector>> + std::vector> hp_quad_dof_identities(const std::set &fes, const unsigned int face_no = 0) const; diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index 03ee305f16..8ce4f1ad4f 100644 --- a/source/dofs/dof_handler_policy.cc +++ b/source/dofs/dof_handler_policy.cc @@ -158,7 +158,7 @@ namespace internal // exists if (identities.get() == nullptr) { - std::vector>> + std::vector> complete_identities; switch (structdim) @@ -194,29 +194,37 @@ namespace internal // pairs (fe_index,dof_index). Because we put in exactly // two fe indices, we know that each entry of the outer // vector needs to contain a set of exactly two such - // pairs. Check this. + // pairs. Check this. While there, also check that + // the two entries actually reference fe_index_1 and + // fe_index_2: for (const auto &complete_identity : complete_identities) - Assert(complete_identity.size() == 2, ExcInternalError()); + { + Assert(complete_identity.size() == 2, ExcInternalError()); + Assert(complete_identity.find(fe_index_1) != + complete_identity.end(), + ExcInternalError()); + Assert(complete_identity.find(fe_index_2) != + complete_identity.end(), + ExcInternalError()); + } #endif // Next reduce these sets of two pairs by removing the // fe_index parts: We know which indices we have. But we // have to make sure in which order we consider the // pair, by considering whether the fe_index part we are - // throwing away matched fe_index_1 or fe_index_2. + // throwing away matched fe_index_1 or fe_index_2. Fortunately, + // this is easy to do because we can ask the std::map for the + // dof_index that matches a given fe_index: DoFIdentities reduced_identities; for (const auto &complete_identity : complete_identities) { - std::pair reduced_identity( - complete_identity.begin()->second, - (++complete_identity.begin())->second); - if (complete_identity.begin()->first == fe_index_1) - reduced_identities.emplace_back(reduced_identity); - else if (complete_identity.begin()->first == fe_index_2) - reduced_identities.emplace_back(reduced_identity.second, - reduced_identity.first); - else - Assert(false, ExcInternalError()); + const unsigned int dof_index_1 = + complete_identity.at(fe_index_1); + const unsigned int dof_index_2 = + complete_identity.at(fe_index_2); + + reduced_identities.emplace_back(dof_index_1, dof_index_2); } #ifdef DEBUG diff --git a/source/hp/fe_collection.cc b/source/hp/fe_collection.cc index 0089db25a9..f8022758fc 100644 --- a/source/hp/fe_collection.cc +++ b/source/hp/fe_collection.cc @@ -296,7 +296,7 @@ namespace hp * Implement the action of the hp_*_dof_identities() functions * in a generic way. */ - std::vector>> + std::vector> compute_hp_dof_identities( const std::set &fes, const std::function>( @@ -318,17 +318,16 @@ namespace hp const auto reduced_identities = query_identities(fe_index_1, fe_index_2); - std::vector>> - complete_identities; + std::vector> complete_identities; for (const auto &reduced_identity : reduced_identities) { // Each identity returned by query_identities() is a pair of // dof indices. Prefix each with its fe index and put the result // into a vector - std::set> - complete_identity = {{fe_index_1, reduced_identity.first}, - {fe_index_2, reduced_identity.second}}; + std::map complete_identity = { + {fe_index_1, reduced_identity.first}, + {fe_index_2, reduced_identity.second}}; complete_identities.emplace_back(std::move(complete_identity)); } @@ -416,7 +415,7 @@ namespace hp // run through because it modifies the graph and thus invalidates // iterators. But because SG stores all of these edges, we can remove them // all from G after collecting the edges in SG.) - std::vector> identities; + std::vector> identities; while (identities_graph.size() > 0) { Graph sub_graph; // SG @@ -470,7 +469,16 @@ namespace hp // At this point we're sure that we have extracted a complete // sub-graph ("clique"). The DoFs involved are all identical then, and // we will store this identity so we can return it later. - identities.emplace_back(std::move(sub_graph_nodes)); + // + // The sub-graph is given as a set of Node objects, which is just + // a collection of (fe_index,dof_index) pairs. Because each + // fe_index can only appear once there, a better data structure + // is a std::map from fe_index to dof_index, which can conveniently + // be initialized from a range of iterators to pairs: + identities.emplace_back(sub_graph_nodes.begin(), + sub_graph_nodes.end()); + Assert(identities.back().size() == sub_graph_nodes.size(), + ExcInternalError()); } return identities; @@ -480,7 +488,7 @@ namespace hp template - std::vector>> + std::vector> FECollection::hp_vertex_dof_identities( const std::set &fes) const { @@ -494,7 +502,7 @@ namespace hp template - std::vector>> + std::vector> FECollection::hp_line_dof_identities( const std::set &fes) const { @@ -508,7 +516,7 @@ namespace hp template - std::vector>> + std::vector> FECollection::hp_quad_dof_identities( const std::set &fes, const unsigned int face_no) const -- 2.39.5