]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Update FE::face_to_cell_index() to use combined orientations.
authorDavid Wells <drwells@email.unc.edu>
Thu, 21 Dec 2023 22:54:49 +0000 (17:54 -0500)
committerDavid Wells <drwells@email.unc.edu>
Fri, 22 Dec 2023 13:42:34 +0000 (08:42 -0500)
15 files changed:
doc/news/changes/incompatibilities/20231221DavidWells [new file with mode: 0644]
include/deal.II/fe/fe.h
include/deal.II/fe/fe_q_base.h
include/deal.II/fe/fe_system.h
include/deal.II/numerics/vector_tools_boundary.templates.h
source/dofs/dof_tools_constraints.cc
source/fe/fe.cc
source/fe/fe_q_base.cc
source/fe/fe_system.cc
tests/fe/face_to_cell_q1_2d.cc
tests/fe/face_to_cell_q2_2d.cc
tests/fe/face_to_cell_q2xq2_2d.cc
tests/fe/face_to_cell_q3_2d.cc
tests/fe/face_to_cell_q3xq4_2d.cc
tests/fe/face_to_cell_q4_2d.cc

diff --git a/doc/news/changes/incompatibilities/20231221DavidWells b/doc/news/changes/incompatibilities/20231221DavidWells
new file mode 100644 (file)
index 0000000..4c1e189
--- /dev/null
@@ -0,0 +1,6 @@
+Changed: FiniteElement::adjust_quad_dof_index_for_face_orientation() and
+FiniteElement::face_to_cell_index() now use the standardized
+<tt>combined_orientation</tt> orientation encoding as input arguments rather
+than three booleans.
+<br>
+(David Wells, 2023/12/21)
index a2334ae281057dfe7f22d769ced159fdb3e51081..b9e093aef7bc75d03d114b53b21e73355d661b57 100644 (file)
@@ -1497,15 +1497,8 @@ public:
    * index must be between zero and dofs_per_face.
    * @param face The number of the face this degree of freedom lives on. This
    * number must be between zero and GeometryInfo::faces_per_cell.
-   * @param face_orientation One part of the description of the orientation of
-   * the face. See
-   * @ref GlossFaceOrientation.
-   * @param face_flip One part of the description of the orientation of the
-   * face. See
-   * @ref GlossFaceOrientation.
-   * @param face_rotation One part of the description of the orientation of
-   * the face. See
-   * @ref GlossFaceOrientation.
+   * @param combined_orientation The combined orientation flag containing the
+   * orientation, rotation, and flip of the face. See @ref GlossFaceOrientation.
    * @return The index of this degree of freedom within the set of degrees of
    * freedom on the entire cell. The returned value will be between zero and
    * dofs_per_cell.
@@ -1527,11 +1520,11 @@ public:
    * freedom actually represent.
    */
   virtual unsigned int
-  face_to_cell_index(const unsigned int face_dof_index,
-                     const unsigned int face,
-                     const bool         face_orientation = true,
-                     const bool         face_flip        = false,
-                     const bool         face_rotation    = false) const;
+  face_to_cell_index(
+    const unsigned int  face_dof_index,
+    const unsigned int  face,
+    const unsigned char combined_orientation =
+      ReferenceCell::default_combined_face_orientation()) const;
 
   /**
    * For lines with non-standard line_orientation in 3d, the dofs on lines
index 6eb213aa448dd07846712341cd0b754f9361d5d7..880bb488fc13db0771399192ec639775e59bb752 100644 (file)
@@ -184,25 +184,18 @@ public:
    * index must be between zero and dofs_per_face.
    * @param face The number of the face this degree of freedom lives on. This
    * number must be between zero and GeometryInfo::faces_per_cell.
-   * @param face_orientation One part of the description of the orientation of
-   * the face. See
-   * @ref GlossFaceOrientation.
-   * @param face_flip One part of the description of the orientation of the
-   * face. See
-   * @ref GlossFaceOrientation.
-   * @param face_rotation One part of the description of the orientation of
-   * the face. See
-   * @ref GlossFaceOrientation.
+   * @param combined_orientation The combined orientation flag containing the
+   * orientation, rotation, and flip of the face. See @ref GlossFaceOrientation.
    * @return The index of this degree of freedom within the set of degrees of
    * freedom on the entire cell. The returned value will be between zero and
    * dofs_per_cell.
    */
   virtual unsigned int
-  face_to_cell_index(const unsigned int face_dof_index,
-                     const unsigned int face,
-                     const bool         face_orientation = true,
-                     const bool         face_flip        = false,
-                     const bool         face_rotation = false) const override;
+  face_to_cell_index(
+    const unsigned int  face_dof_index,
+    const unsigned int  face,
+    const unsigned char combined_orientation =
+      ReferenceCell::default_combined_face_orientation()) const override;
 
   /**
    * Return a list of constant modes of the element. For this element, the
index 1271c4a3e60d1c4d4359376f3ab0d280d7eac019..7f055dc7312f4f833610c12388faa9efe4ff4cd0 100644 (file)
@@ -880,25 +880,18 @@ public:
    * index must be between zero and dofs_per_face.
    * @param face The number of the face this degree of freedom lives on. This
    * number must be between zero and GeometryInfo::faces_per_cell.
-   * @param face_orientation One part of the description of the orientation of
-   * the face. See
-   * @ref GlossFaceOrientation.
-   * @param face_flip One part of the description of the orientation of the
-   * face. See
-   * @ref GlossFaceOrientation.
-   * @param face_rotation One part of the description of the orientation of
-   * the face. See
-   * @ref GlossFaceOrientation.
+   * @param combined_orientation The combined orientation flag containing the
+   * orientation, rotation, and flip of the face. See @ref GlossFaceOrientation.
    * @return The index of this degree of freedom within the set of degrees of
    * freedom on the entire cell. The returned value will be between zero and
    * dofs_per_cell.
    */
   virtual unsigned int
-  face_to_cell_index(const unsigned int face_dof_index,
-                     const unsigned int face,
-                     const bool         face_orientation = true,
-                     const bool         face_flip        = false,
-                     const bool         face_rotation = false) const override;
+  face_to_cell_index(
+    const unsigned int  face_dof_index,
+    const unsigned int  face,
+    const unsigned char combined_orientation =
+      ReferenceCell::default_combined_face_orientation()) const override;
 
   /**
    * Implementation of the respective function in the base class.
index 7b8b9d4a42168f5a9cdf2ba2458e3aef1ec7c7a2..c2cfa6bb076d554c6fb81f89bf5269ed2a3dc65e 100644 (file)
@@ -2192,11 +2192,8 @@ namespace VectorTools
             dof_values(i) +=
               tmp * (normals[q_point] *
                      fe_values[vec].value(
-                       fe.face_to_cell_index(i,
-                                             face,
-                                             cell->face_orientation(face),
-                                             cell->face_flip(face),
-                                             cell->face_rotation(face)),
+                       fe.face_to_cell_index(
+                         i, face, cell->combined_face_orientation(face)),
                        q_point));
         }
 
@@ -2213,9 +2210,7 @@ namespace VectorTools
             fe.get_nonzero_components(fe.face_to_cell_index(
               i,
               face,
-              cell->face_orientation(face),
-              cell->face_flip(face),
-              cell->face_rotation(face)))[first_vector_component])
+              cell->combined_face_orientation(face)))[first_vector_component])
           constraints.add_constraint(face_dof_indices[i],
                                      {},
                                      (std::abs(dof_values[i]) > 1e-14 ?
@@ -2303,11 +2298,8 @@ namespace VectorTools
             dof_values_local(i) +=
               tmp * (normals[q_point] *
                      fe_values[vec].value(
-                       fe.face_to_cell_index(i,
-                                             face,
-                                             cell->face_orientation(face),
-                                             cell->face_flip(face),
-                                             cell->face_rotation(face)),
+                       fe.face_to_cell_index(
+                         i, face, cell->combined_face_orientation(face)),
                        q_point));
         }
 
@@ -2322,9 +2314,7 @@ namespace VectorTools
             fe.get_nonzero_components(fe.face_to_cell_index(
               i,
               face,
-              cell->face_orientation(face),
-              cell->face_flip(face),
-              cell->face_rotation(face)))[first_vector_component])
+              cell->combined_face_orientation(face)))[first_vector_component])
           {
             dof_values[face_dof_indices[i]]     = dof_values_local(i);
             projected_dofs[face_dof_indices[i]] = fe.degree;
index 392eeff2f34ad044b852c74fb581e5ee82690f23..518a2ecf6309896566ba62392aedc56c444dc514 100644 (file)
@@ -1869,6 +1869,11 @@ namespace DoFTools
       static const int dim      = FaceIterator::AccessorType::dimension;
       static const int spacedim = FaceIterator::AccessorType::space_dimension;
 
+      const unsigned char combined_face_orientation =
+        ::dealii::internal::combined_face_orientation(face_orientation,
+                                                      face_rotation,
+                                                      face_flip);
+
       const bool use_mg = (level != numbers::invalid_unsigned_int);
 
       // If we don't use multigrid, we should be in the case where face_1 is
@@ -2008,15 +2013,10 @@ namespace DoFTools
       // Build up a cell to face index for face_2:
       for (unsigned int i = 0; i < dofs_per_face; ++i)
         {
-          const unsigned int cell_index =
-            fe.face_to_cell_index(i,
-                                  0, /* It doesn't really matter, just
-                                      * assume we're on the first face...
-                                      */
-                                  true,
-                                  false,
-                                  false // default orientation
-            );
+          const unsigned int cell_index = fe.face_to_cell_index(
+            i,
+            // It doesn't really matter, just assume we're on the first face...
+            0);
           cell_to_rotated_face_index[cell_index] = i;
         }
 
@@ -2100,7 +2100,7 @@ namespace DoFTools
                   // given orientation:
                   const unsigned int j =
                     cell_to_rotated_face_index[fe.face_to_cell_index(
-                      jj, 0, face_orientation, face_flip, face_rotation)];
+                      jj, 0, combined_face_orientation)];
 
                   if (std::abs(transformation(i, jj)) > eps)
                     constraint_entries.emplace_back(dofs_1[j],
@@ -2123,7 +2123,7 @@ namespace DoFTools
           // orientation:
           const unsigned int j =
             cell_to_rotated_face_index[fe.face_to_cell_index(
-              target, 0, face_orientation, face_flip, face_rotation)];
+              target, 0, combined_face_orientation)];
 
           auto dof_left  = dofs_1[j];
           auto dof_right = dofs_2[i];
index 5238f2cae45d816f756975796b14df76808c56d0..e7e54d3f6ad892af75b5920ea4a42b2d30c544c6 100644 (file)
@@ -562,11 +562,10 @@ FiniteElement<dim, spacedim>::block_mask(
 
 template <int dim, int spacedim>
 unsigned int
-FiniteElement<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
-                                                 const unsigned int face,
-                                                 const bool face_orientation,
-                                                 const bool face_flip,
-                                                 const bool face_rotation) const
+FiniteElement<dim, spacedim>::face_to_cell_index(
+  const unsigned int  face_index,
+  const unsigned int  face,
+  const unsigned char combined_orientation) const
 {
   AssertIndexRange(face_index, this->n_dofs_per_face(face));
   AssertIndexRange(face, this->reference_cell().n_faces());
@@ -583,9 +582,9 @@ FiniteElement<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
   // see the function's documentation for an explanation of this
   // assertion -- in essence, derived classes have to implement
   // an overloaded version of this function if we are to use any
-  // other than standard orientation
-  if ((face_orientation != true) || (face_flip != false) ||
-      (face_rotation != false))
+  // other than default (standard) orientation
+  if (combined_orientation !=
+      ReferenceCell::default_combined_face_orientation())
     Assert((this->n_dofs_per_line() <= 1) && (this->n_dofs_per_quad(face) <= 1),
            ExcMessage(
              "The function in this base class can not handle this case. "
@@ -607,11 +606,7 @@ FiniteElement<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
       // then get the number of this vertex on the cell and translate
       // this to a DoF number on the cell
       return (this->reference_cell().face_to_cell_vertices(
-                face,
-                face_vertex,
-                internal::combined_face_orientation(face_orientation,
-                                                    face_rotation,
-                                                    face_flip)) *
+                face, face_vertex, combined_orientation) *
                 this->n_dofs_per_vertex() +
               dof_index_on_vertex);
     }
@@ -627,12 +622,9 @@ FiniteElement<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
       const unsigned int dof_index_on_line = index % this->n_dofs_per_line();
 
       return (this->get_first_line_index() +
-              this->reference_cell().face_to_cell_lines(
-                face,
-                face_line,
-                internal::combined_face_orientation(face_orientation,
-                                                    face_rotation,
-                                                    face_flip)) *
+              this->reference_cell().face_to_cell_lines(face,
+                                                        face_line,
+                                                        combined_orientation) *
                 this->n_dofs_per_line() +
               dof_index_on_line);
     }
index 890cd39dfd952884626249a8f1eef944ebeea469..40b21bade6d97f86262b59fc60eb29ff5680169a 100644 (file)
@@ -1113,11 +1113,10 @@ FE_Q_Base<dim, spacedim>::initialize_quad_dof_index_permutation()
 
 template <int dim, int spacedim>
 unsigned int
-FE_Q_Base<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
-                                             const unsigned int face,
-                                             const bool face_orientation,
-                                             const bool face_flip,
-                                             const bool face_rotation) const
+FE_Q_Base<dim, spacedim>::face_to_cell_index(
+  const unsigned int  face_index,
+  const unsigned int  face,
+  const unsigned char combined_orientation) const
 {
   AssertIndexRange(face_index, this->n_dofs_per_face(face));
   AssertIndexRange(face, GeometryInfo<dim>::faces_per_cell);
@@ -1144,10 +1143,10 @@ FE_Q_Base<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
 
       // then get the number of this vertex on the cell and translate
       // this to a DoF number on the cell
-      return (GeometryInfo<dim>::face_to_cell_vertices(
-                face, face_vertex, face_orientation, face_flip, face_rotation) *
-                this->n_dofs_per_vertex() +
-              dof_index_on_vertex);
+      return this->reference_cell().face_to_cell_vertices(
+               face, face_vertex, combined_orientation) *
+               this->n_dofs_per_vertex() +
+             dof_index_on_vertex;
     }
   else if (face_index < this->get_first_face_quad_index(face))
     // DoF is on a face
@@ -1170,9 +1169,8 @@ FE_Q_Base<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
             break;
 
           case 2:
-            // in 2d, only face_orientation has a meaning. if it is false (i.e.,
-            // the non-default case), then consider dofs in reverse order
-            if (face_orientation == true)
+            if (combined_orientation ==
+                ReferenceCell::default_combined_face_orientation())
               adjusted_dof_index_on_line = dof_index_on_line;
             else
               adjusted_dof_index_on_line =
@@ -1188,8 +1186,8 @@ FE_Q_Base<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
             // that said, the Q2 case is easy enough to implement, as is the
             // case where everything is in standard orientation
             Assert((this->n_dofs_per_line() <= 1) ||
-                     ((face_orientation == true) && (face_flip == false) &&
-                      (face_rotation == false)),
+                     combined_orientation ==
+                       ReferenceCell::default_combined_face_orientation(),
                    ExcNotImplemented());
             adjusted_dof_index_on_line = dof_index_on_line;
             break;
@@ -1199,8 +1197,9 @@ FE_Q_Base<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
         }
 
       return (this->get_first_line_index() +
-              GeometryInfo<dim>::face_to_cell_lines(
-                face, face_line, face_orientation, face_flip, face_rotation) *
+              this->reference_cell().face_to_cell_lines(face,
+                                                        face_line,
+                                                        combined_orientation) *
                 this->n_dofs_per_line() +
               adjusted_dof_index_on_line);
     }
@@ -1217,8 +1216,8 @@ FE_Q_Base<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
       // just have to draw a bunch of pictures. in the meantime,
       // we can implement the Q2 case in which it is simple
       Assert((this->n_dofs_per_quad(face) <= 1) ||
-               ((face_orientation == true) && (face_flip == false) &&
-                (face_rotation == false)),
+               combined_orientation ==
+                 ReferenceCell::default_combined_face_orientation(),
              ExcNotImplemented());
       return (this->get_first_quad_index(face) + index);
     }
index 4df28e839ce1a68b883daac961bc0ec1d5e8da8d..587e4b272a30a495bb21ed65af831253f3c036ec 100644 (file)
@@ -1027,11 +1027,10 @@ FESystem<dim, spacedim>::get_prolongation_matrix(
 
 template <int dim, int spacedim>
 unsigned int
-FESystem<dim, spacedim>::face_to_cell_index(const unsigned int face_dof_index,
-                                            const unsigned int face,
-                                            const bool         face_orientation,
-                                            const bool         face_flip,
-                                            const bool face_rotation) const
+FESystem<dim, spacedim>::face_to_cell_index(
+  const unsigned int  face_dof_index,
+  const unsigned int  face,
+  const unsigned char combined_orientation) const
 {
   // we need to ask the base elements how they want to translate
   // the DoFs within their own numbering. thus, translate to
@@ -1041,11 +1040,7 @@ FESystem<dim, spacedim>::face_to_cell_index(const unsigned int face_dof_index,
 
   const unsigned int base_face_to_cell_index =
     this->base_element(face_base_index.first.first)
-      .face_to_cell_index(face_base_index.second,
-                          face,
-                          face_orientation,
-                          face_flip,
-                          face_rotation);
+      .face_to_cell_index(face_base_index.second, face, combined_orientation);
 
   // it would be nice if we had a base_to_system_index function, but
   // all that exists is a component_to_system_index function. we can't do
index d982f0b45bd6df2fed3818a417e35051c3073f79..603f72da423cb03856c29205f12f2791f5f8e641 100644 (file)
@@ -35,18 +35,13 @@ test()
     {
       deallog << "Face=" << face << std::endl;
 
-      for (int orientation = 0; orientation < 2; ++orientation)
+      for (unsigned char orientation = 0; orientation < 2; ++orientation)
         {
           deallog << "  orientation=" << (orientation == 0 ? "false" : "true")
                   << std::endl
                   << "    ";
           for (unsigned int i = 0; i < dofs_per_face; ++i)
-            deallog << fe.face_to_cell_index(i,
-                                             face,
-                                             (orientation == 0 ? false : true),
-                                             false,
-                                             false)
-                    << " - ";
+            deallog << fe.face_to_cell_index(i, face, orientation) << " - ";
           deallog << std::endl;
         }
     }
index 4ae0330b4e6e97b55634f4a8855c9c65b712ae6c..87bc64166b09ae0d1ff9deb3d900bc8c0a32c1ae 100644 (file)
@@ -35,18 +35,13 @@ test()
     {
       deallog << "Face=" << face << std::endl;
 
-      for (int orientation = 0; orientation < 2; ++orientation)
+      for (unsigned char orientation = 0; orientation < 2; ++orientation)
         {
           deallog << "  orientation=" << (orientation == 0 ? "false" : "true")
                   << std::endl
                   << "    ";
           for (unsigned int i = 0; i < dofs_per_face; ++i)
-            deallog << fe.face_to_cell_index(i,
-                                             face,
-                                             (orientation == 0 ? false : true),
-                                             false,
-                                             false)
-                    << " - ";
+            deallog << fe.face_to_cell_index(i, face, orientation) << " - ";
           deallog << std::endl;
         }
     }
index 2dbb7382288eaa2b0b71ac0795acf22746800ac5..50309ba64189715e00d129512264cae39383a1fc 100644 (file)
@@ -36,18 +36,13 @@ test()
     {
       deallog << "Face=" << face << std::endl;
 
-      for (int orientation = 0; orientation < 2; ++orientation)
+      for (unsigned char orientation = 0; orientation < 2; ++orientation)
         {
           deallog << "  orientation=" << (orientation == 0 ? "false" : "true")
                   << std::endl
                   << "    ";
           for (unsigned int i = 0; i < dofs_per_face; ++i)
-            deallog << fe.face_to_cell_index(i,
-                                             face,
-                                             (orientation == 0 ? false : true),
-                                             false,
-                                             false)
-                    << " - ";
+            deallog << fe.face_to_cell_index(i, face, orientation) << " - ";
           deallog << std::endl;
         }
     }
index 4357d1f31af5a029b7976495d00363584af54a96..ee7721ebba704b0885ef63e25227e2a40f6896e8 100644 (file)
@@ -35,18 +35,13 @@ test()
     {
       deallog << "Face=" << face << std::endl;
 
-      for (int orientation = 0; orientation < 2; ++orientation)
+      for (unsigned char orientation = 0; orientation < 2; ++orientation)
         {
           deallog << "  orientation=" << (orientation == 0 ? "false" : "true")
                   << std::endl
                   << "    ";
           for (unsigned int i = 0; i < dofs_per_face; ++i)
-            deallog << fe.face_to_cell_index(i,
-                                             face,
-                                             (orientation == 0 ? false : true),
-                                             false,
-                                             false)
-                    << " - ";
+            deallog << fe.face_to_cell_index(i, face, orientation) << " - ";
           deallog << std::endl;
         }
     }
index 562934aaec0715959d2b3399254eb8c7317d1c63..874d3394da62ed20ce9d07acd4547ca703fe6e64 100644 (file)
@@ -36,18 +36,13 @@ test()
     {
       deallog << "Face=" << face << std::endl;
 
-      for (int orientation = 0; orientation < 2; ++orientation)
+      for (unsigned char orientation = 0; orientation < 2; ++orientation)
         {
           deallog << "  orientation=" << (orientation == 0 ? "false" : "true")
                   << std::endl
                   << "    ";
           for (unsigned int i = 0; i < dofs_per_face; ++i)
-            deallog << fe.face_to_cell_index(i,
-                                             face,
-                                             (orientation == 0 ? false : true),
-                                             false,
-                                             false)
-                    << " - ";
+            deallog << fe.face_to_cell_index(i, face, orientation) << " - ";
           deallog << std::endl;
         }
     }
index 3bdde4eadaf76498e0f2a6d9a7de3e405b19a805..56efd4e29897a859079f45d61763d85f02af465c 100644 (file)
@@ -35,18 +35,13 @@ test()
     {
       deallog << "Face=" << face << std::endl;
 
-      for (int orientation = 0; orientation < 2; ++orientation)
+      for (unsigned char orientation = 0; orientation < 2; ++orientation)
         {
           deallog << "  orientation=" << (orientation == 0 ? "false" : "true")
                   << std::endl
                   << "    ";
           for (unsigned int i = 0; i < dofs_per_face; ++i)
-            deallog << fe.face_to_cell_index(i,
-                                             face,
-                                             (orientation == 0 ? false : true),
-                                             false,
-                                             false)
-                    << " - ";
+            deallog << fe.face_to_cell_index(i, face, orientation) << " - ";
           deallog << std::endl;
         }
     }

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.