]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Simplify code logic. 16208/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Mon, 23 Oct 2023 23:30:54 +0000 (17:30 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 27 Oct 2023 22:13:28 +0000 (16:13 -0600)
source/fe/fe_system.cc

index da84c7809fd40f9f8316f953e84e4aa53ba3aec4..4df28e839ce1a68b883daac961bc0ec1d5e8da8d 100644 (file)
@@ -905,9 +905,6 @@ FESystem<dim, spacedim>::get_restriction_matrix(
           this->n_dofs_per_cell())
         return this->restriction[refinement_case - 1][child];
 
-      // Check if some of the matrices of the base elements are void.
-      bool do_restriction = true;
-
       // shortcut for accessing local restrictions further down
       std::vector<const FullMatrix<double> *> base_matrices(
         this->n_base_elements());
@@ -916,57 +913,50 @@ FESystem<dim, spacedim>::get_restriction_matrix(
         {
           base_matrices[i] =
             &base_element(i).get_restriction_matrix(child, refinement_case);
-          if (base_matrices[i]->n() != base_element(i).n_dofs_per_cell())
-            do_restriction = false;
-        }
-      Assert(do_restriction,
-             (typename FiniteElement<dim, spacedim>::ExcProjectionVoid()));
 
-      // if we did not encounter void matrices, initialize the matrix sizes
-      if (do_restriction)
-        {
-          FullMatrix<double> restriction(this->n_dofs_per_cell(),
-                                         this->n_dofs_per_cell());
+          Assert(base_matrices[i]->n() == base_element(i).n_dofs_per_cell(),
+                 (typename FiniteElement<dim, spacedim>::ExcProjectionVoid()));
+        }
 
-          // distribute the matrices of the base finite elements to the
-          // matrices of this object. for this, loop over all degrees of
-          // freedom and take the respective entry of the underlying base
-          // element.
-          //
-          // note that we by definition of a base element, they are
-          // independent, i.e. do not couple. only DoFs that belong to the
-          // same instance of a base element may couple
-          for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
-            for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
-              {
-                // first find out to which base element indices i and j
-                // belong, and which instance thereof in case the base element
-                // has a multiplicity greater than one. if they should not
-                // happen to belong to the same instance of a base element,
-                // then they cannot couple, so go on with the next index
-                if (this->system_to_base_table[i].first !=
-                    this->system_to_base_table[j].first)
-                  continue;
-
-                // so get the common base element and the indices therein:
-                const unsigned int base =
-                  this->system_to_base_table[i].first.first;
-
-                const unsigned int base_index_i =
-                                     this->system_to_base_table[i].second,
-                                   base_index_j =
-                                     this->system_to_base_table[j].second;
-
-                // if we are sure that DoFs i and j may couple, then copy
-                // entries of the matrices:
-                restriction(i, j) =
-                  (*base_matrices[base])(base_index_i, base_index_j);
-              }
+      FullMatrix<double> restriction(this->n_dofs_per_cell(),
+                                     this->n_dofs_per_cell());
+
+      // distribute the matrices of the base finite elements to the
+      // matrices of this object. for this, loop over all degrees of
+      // freedom and take the respective entry of the underlying base
+      // element.
+      //
+      // note that we by definition of a base element, they are
+      // independent, i.e. do not couple. only DoFs that belong to the
+      // same instance of a base element may couple
+      for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
+        for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
+          {
+            // first find out to which base element indices i and j
+            // belong, and which instance thereof in case the base element
+            // has a multiplicity greater than one. if they should not
+            // happen to belong to the same instance of a base element,
+            // then they cannot couple, so go on with the next index
+            if (this->system_to_base_table[i].first !=
+                this->system_to_base_table[j].first)
+              continue;
+
+            // so get the common base element and the indices therein:
+            const unsigned int base = this->system_to_base_table[i].first.first;
+
+            const unsigned int base_index_i =
+                                 this->system_to_base_table[i].second,
+                               base_index_j =
+                                 this->system_to_base_table[j].second;
+
+            // if we are sure that DoFs i and j may couple, then copy
+            // entries of the matrices:
+            restriction(i, j) =
+              (*base_matrices[base])(base_index_i, base_index_j);
+          }
 
-          const_cast<FullMatrix<double> &>(
-            this->restriction[refinement_case - 1][child]) =
-            std::move(restriction);
-        }
+      const_cast<FullMatrix<double> &>(
+        this->restriction[refinement_case - 1][child]) = std::move(restriction);
     }
 
   return this->restriction[refinement_case - 1][child];
@@ -997,45 +987,38 @@ FESystem<dim, spacedim>::get_prolongation_matrix(
           this->n_dofs_per_cell())
         return this->prolongation[refinement_case - 1][child];
 
-      bool                                    do_prolongation = true;
       std::vector<const FullMatrix<double> *> base_matrices(
         this->n_base_elements());
       for (unsigned int i = 0; i < this->n_base_elements(); ++i)
         {
           base_matrices[i] =
             &base_element(i).get_prolongation_matrix(child, refinement_case);
-          if (base_matrices[i]->n() != base_element(i).n_dofs_per_cell())
-            do_prolongation = false;
+
+          Assert(base_matrices[i]->n() == base_element(i).n_dofs_per_cell(),
+                 (typename FiniteElement<dim, spacedim>::ExcEmbeddingVoid()));
         }
-      Assert(do_prolongation,
-             (typename FiniteElement<dim, spacedim>::ExcEmbeddingVoid()));
 
-      if (do_prolongation)
-        {
-          FullMatrix<double> prolongate(this->n_dofs_per_cell(),
-                                        this->n_dofs_per_cell());
+      FullMatrix<double> prolongate(this->n_dofs_per_cell(),
+                                    this->n_dofs_per_cell());
 
-          for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
-            for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
-              {
-                if (this->system_to_base_table[i].first !=
-                    this->system_to_base_table[j].first)
-                  continue;
-                const unsigned int base =
-                  this->system_to_base_table[i].first.first;
-
-                const unsigned int base_index_i =
-                                     this->system_to_base_table[i].second,
-                                   base_index_j =
-                                     this->system_to_base_table[j].second;
-                prolongate(i, j) =
-                  (*base_matrices[base])(base_index_i, base_index_j);
-              }
+      for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
+        for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
+          {
+            if (this->system_to_base_table[i].first !=
+                this->system_to_base_table[j].first)
+              continue;
+            const unsigned int base = this->system_to_base_table[i].first.first;
+
+            const unsigned int base_index_i =
+                                 this->system_to_base_table[i].second,
+                               base_index_j =
+                                 this->system_to_base_table[j].second;
+            prolongate(i, j) =
+              (*base_matrices[base])(base_index_i, base_index_j);
+          }
 
-          const_cast<FullMatrix<double> &>(
-            this->prolongation[refinement_case - 1][child]) =
-            std::move(prolongate);
-        }
+      const_cast<FullMatrix<double> &>(
+        this->prolongation[refinement_case - 1][child]) = std::move(prolongate);
     }
 
   return this->prolongation[refinement_case - 1][child];

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.