]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Change the return type of hp::FECollection::hp_*_dof_identities. 12716/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 26 Aug 2021 21:27:53 +0000 (15:27 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 27 Aug 2021 17:14:58 +0000 (11:14 -0600)
include/deal.II/hp/fe_collection.h
source/dofs/dof_handler_policy.cc
source/hp/fe_collection.cc

index a5d4ba63260fb83f9e1c65df41513d49d424ceb6..428552e0f11e2777aa90f49933fb8c2e76bbb738 100644 (file)
@@ -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::set<std::pair<unsigned int, unsigned int>>>
+    std::vector<std::map<unsigned int, unsigned int>>
     hp_vertex_dof_identities(const std::set<unsigned int> &fes) const;
 
     /**
      * Same as hp_vertex_dof_indices(), except that the function treats degrees
      * of freedom on lines.
      */
-    std::vector<std::set<std::pair<unsigned int, unsigned int>>>
+    std::vector<std::map<unsigned int, unsigned int>>
     hp_line_dof_identities(const std::set<unsigned int> &fes) const;
 
     /**
      * Same as hp_vertex_dof_indices(), except that the function treats degrees
      * of freedom on quads.
      */
-    std::vector<std::set<std::pair<unsigned int, unsigned int>>>
+    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;
 
index 03ee305f1604a64b77fc64248badc7c71048f8ef..8ce4f1ad4fc6142dccbd9659a20dc9f60836bf74 100644 (file)
@@ -158,7 +158,7 @@ namespace internal
           // exists
           if (identities.get() == nullptr)
             {
-              std::vector<std::set<std::pair<unsigned int, unsigned int>>>
+              std::vector<std::map<unsigned int, unsigned int>>
                 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<unsigned int, unsigned int> 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
index 0089db25a9ae482826e942c2b7deab9f4b1a6cf6..f8022758fc8db01d6cabd8e1b9c8497431df03e1 100644 (file)
@@ -296,7 +296,7 @@ namespace hp
      * Implement the action of the hp_*_dof_identities() functions
      * in a generic way.
      */
-    std::vector<std::set<std::pair<unsigned int, unsigned int>>>
+    std::vector<std::map<unsigned int, unsigned int>>
     compute_hp_dof_identities(
       const std::set<unsigned int> &fes,
       const std::function<std::vector<std::pair<unsigned int, unsigned int>>(
@@ -318,17 +318,16 @@ namespace hp
           const auto         reduced_identities =
             query_identities(fe_index_1, fe_index_2);
 
-          std::vector<std::set<std::pair<unsigned int, unsigned int>>>
-            complete_identities;
+          std::vector<std::map<unsigned int, unsigned int>> 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<std::pair<unsigned int, unsigned int>>
-                complete_identity = {{fe_index_1, reduced_identity.first},
-                                     {fe_index_2, reduced_identity.second}};
+              std::map<unsigned int, unsigned int> 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<std::set<Node>> identities;
+      std::vector<std::map<unsigned int, unsigned int>> 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 <int dim, int spacedim>
-  std::vector<std::set<std::pair<unsigned int, unsigned int>>>
+  std::vector<std::map<unsigned int, unsigned int>>
   FECollection<dim, spacedim>::hp_vertex_dof_identities(
     const std::set<unsigned int> &fes) const
   {
@@ -494,7 +502,7 @@ namespace hp
 
 
   template <int dim, int spacedim>
-  std::vector<std::set<std::pair<unsigned int, unsigned int>>>
+  std::vector<std::map<unsigned int, unsigned int>>
   FECollection<dim, spacedim>::hp_line_dof_identities(
     const std::set<unsigned int> &fes) const
   {
@@ -508,7 +516,7 @@ namespace hp
 
 
   template <int dim, int spacedim>
-  std::vector<std::set<std::pair<unsigned int, unsigned int>>>
+  std::vector<std::map<unsigned int, unsigned int>>
   FECollection<dim, spacedim>::hp_quad_dof_identities(
     const std::set<unsigned int> &fes,
     const unsigned int            face_no) const

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.