]> https://gitweb.dealii.org/ - dealii.git/commitdiff
ReferenceCell: consolidate some tables. 14710/head
authorDavid Wells <drwells@email.unc.edu>
Sun, 22 Jan 2023 00:30:12 +0000 (19:30 -0500)
committerDavid Wells <drwells@email.unc.edu>
Sun, 22 Jan 2023 20:21:49 +0000 (15:21 -0500)
include/deal.II/grid/reference_cell.h
source/grid/reference_cell.cc

index 5ff5140e4d54b0220a64624222fddffe391ac50f..bbd9682be668b9b166f88cc81ac08a7026bcb911 100644 (file)
@@ -806,6 +806,37 @@ private:
    */
   constexpr ReferenceCell(const std::uint8_t kind);
 
+  /**
+   * Table containing all vertex permutations for a line.
+   */
+  static constexpr ndarray<unsigned int, 2, 2> line_vertex_permutations = {
+    {{{1, 0}}, {{0, 1}}}};
+
+
+  /**
+   * Table containing all vertex permutations for a triangle.
+   */
+  static constexpr ndarray<unsigned int, 6, 3> triangle_vertex_permutations = {
+    {{{0, 2, 1}},
+     {{0, 1, 2}},
+     {{2, 1, 0}},
+     {{1, 2, 0}},
+     {{1, 0, 2}},
+     {{2, 0, 1}}}};
+
+  /**
+   * Table containing all vertex permutations for a quadrilateral.
+   */
+  static constexpr ndarray<unsigned int, 8, 4>
+    quadrilateral_vertex_permutations = {{{{0, 2, 1, 3}},
+                                          {{0, 1, 2, 3}},
+                                          {{2, 3, 0, 1}},
+                                          {{2, 0, 3, 1}},
+                                          {{3, 1, 2, 0}},
+                                          {{3, 2, 1, 0}},
+                                          {{1, 0, 3, 2}},
+                                          {{1, 3, 0, 2}}}};
+
   /**
    * A kind of constructor -- not quite private because it can be
    * called by anyone, but at least hidden in an internal namespace.
@@ -1865,80 +1896,23 @@ ReferenceCell::standard_to_real_face_vertex(
         break;
       case ReferenceCells::Triangle:
       case ReferenceCells::Quadrilateral:
-        {
-          static constexpr ndarray<unsigned int, 2, 2> table = {
-            {{{1, 0}}, {{0, 1}}}};
-
-          return table[face_orientation][vertex];
-        }
+        return line_vertex_permutations[face_orientation][vertex];
       case ReferenceCells::Tetrahedron:
-        {
-          static constexpr ndarray<unsigned int, 6, 3> table = {{{{0, 2, 1}},
-                                                                 {{0, 1, 2}},
-                                                                 {{2, 1, 0}},
-                                                                 {{1, 2, 0}},
-                                                                 {{1, 0, 2}},
-                                                                 {{2, 0, 1}}}};
-
-          return table[face_orientation][vertex];
-        }
+        return triangle_vertex_permutations[face_orientation][vertex];
       case ReferenceCells::Pyramid:
-        {
-          if (face == 0) // The quadrilateral face
-            {
-              return GeometryInfo<3>::standard_to_real_face_vertex(
-                vertex,
-                Utilities::get_bit(face_orientation, 0),
-                Utilities::get_bit(face_orientation, 2),
-                Utilities::get_bit(face_orientation, 1));
-            }
-          else // One of the triangular faces
-            {
-              static const ndarray<unsigned int, 6, 3> table = {{{{0, 2, 1}},
-                                                                 {{0, 1, 2}},
-                                                                 {{2, 1, 0}},
-                                                                 {{1, 2, 0}},
-                                                                 {{1, 0, 2}},
-                                                                 {{2, 0, 1}}}};
-
-              return table[face_orientation][vertex];
-            }
-        }
+        // face 0 is a quadrilateral
+        if (face == 0)
+          return quadrilateral_vertex_permutations[face_orientation][vertex];
+        else
+          return triangle_vertex_permutations[face_orientation][vertex];
       case ReferenceCells::Wedge:
-        {
-          if (face > 1) // One of the quadrilateral faces
-            {
-              return GeometryInfo<3>::standard_to_real_face_vertex(
-                vertex,
-                Utilities::get_bit(face_orientation, 0),
-                Utilities::get_bit(face_orientation, 2),
-                Utilities::get_bit(face_orientation, 1));
-            }
-          else // One of the triangular faces
-            {
-              static const ndarray<unsigned int, 6, 3> table = {{{{0, 2, 1}},
-                                                                 {{0, 1, 2}},
-                                                                 {{2, 1, 0}},
-                                                                 {{1, 2, 0}},
-                                                                 {{1, 0, 2}},
-                                                                 {{2, 0, 1}}}};
-
-              return table[face_orientation][vertex];
-            }
-        }
+        // faces 0 and 1 are triangles
+        if (face > 1)
+          return quadrilateral_vertex_permutations[face_orientation][vertex];
+        else
+          return triangle_vertex_permutations[face_orientation][vertex];
       case ReferenceCells::Hexahedron:
-        {
-          static constexpr ndarray<unsigned int, 8, 4> table = {
-            {{{0, 2, 1, 3}},
-             {{0, 1, 2, 3}},
-             {{2, 3, 0, 1}},
-             {{2, 0, 3, 1}},
-             {{3, 1, 2, 0}},
-             {{3, 2, 1, 0}},
-             {{1, 0, 3, 2}},
-             {{1, 3, 0, 2}}}};
-          return table[face_orientation][vertex];
-        }
+        return quadrilateral_vertex_permutations[face_orientation][vertex];
       default:
         Assert(false, ExcNotImplemented());
     }
@@ -1958,6 +1932,14 @@ ReferenceCell::standard_to_real_face_line(
   AssertIndexRange(face, n_faces());
   AssertIndexRange(line, face_reference_cell(face).n_lines());
 
+  static constexpr ndarray<unsigned int, 6, 3> triangle_table = {{{{2, 1, 0}},
+                                                                  {{0, 1, 2}},
+                                                                  {{1, 0, 2}},
+                                                                  {{1, 2, 0}},
+                                                                  {{0, 2, 1}},
+                                                                  {{2, 0, 1}}}};
+
+
   switch (this->kind)
     {
       case ReferenceCells::Vertex:
@@ -1967,62 +1949,31 @@ ReferenceCell::standard_to_real_face_line(
         Assert(false, ExcNotImplemented());
         break;
       case ReferenceCells::Tetrahedron:
-        {
-          static constexpr ndarray<unsigned int, 6, 3> table = {{{{2, 1, 0}},
-                                                                 {{0, 1, 2}},
-                                                                 {{1, 0, 2}},
-                                                                 {{1, 2, 0}},
-                                                                 {{0, 2, 1}},
-                                                                 {{2, 0, 1}}}};
-
-          return table[face_orientation][line];
-        }
+        return triangle_table[face_orientation][line];
       case ReferenceCells::Pyramid:
-        {
-          if (face == 0) // The quadrilateral face
-            {
-              return GeometryInfo<3>::standard_to_real_face_line(
-                line,
-                Utilities::get_bit(face_orientation, 0),
-                Utilities::get_bit(face_orientation, 2),
-                Utilities::get_bit(face_orientation, 1));
-            }
-          else // One of the triangular faces
-            {
-              static constexpr ndarray<unsigned int, 6, 3> table = {
-                {{{2, 1, 0}},
-                 {{0, 1, 2}},
-                 {{1, 0, 2}},
-                 {{1, 2, 0}},
-                 {{0, 2, 1}},
-                 {{2, 0, 1}}}};
-
-              return table[face_orientation][line];
-            }
-        }
+        if (face == 0) // The quadrilateral face
+          {
+            return GeometryInfo<3>::standard_to_real_face_line(
+              line,
+              Utilities::get_bit(face_orientation, 0),
+              Utilities::get_bit(face_orientation, 2),
+              Utilities::get_bit(face_orientation, 1));
+          }
+        else // One of the triangular faces
+          {
+            return triangle_table[face_orientation][line];
+          }
       case ReferenceCells::Wedge:
-        {
-          if (face > 1) // One of the quadrilateral faces
-            {
-              return GeometryInfo<3>::standard_to_real_face_line(
-                line,
-                Utilities::get_bit(face_orientation, 0),
-                Utilities::get_bit(face_orientation, 2),
-                Utilities::get_bit(face_orientation, 1));
-            }
-          else // One of the triangular faces
-            {
-              static constexpr ndarray<unsigned int, 6, 3> table = {
-                {{{2, 1, 0}},
-                 {{0, 1, 2}},
-                 {{1, 0, 2}},
-                 {{1, 2, 0}},
-                 {{0, 2, 1}},
-                 {{2, 0, 1}}}};
-
-              return table[face_orientation][line];
-            }
-        }
+        if (face > 1) // One of the quadrilateral faces
+          {
+            return GeometryInfo<3>::standard_to_real_face_line(
+              line,
+              Utilities::get_bit(face_orientation, 0),
+              Utilities::get_bit(face_orientation, 2),
+              Utilities::get_bit(face_orientation, 1));
+          }
+        else // One of the triangular faces
+          return triangle_table[face_orientation][line];
       case ReferenceCells::Hexahedron:
         {
           static constexpr ndarray<unsigned int, 8, 4> table = {
index 3ae63c455624116fd71ec5d6538a21497410475d..0fb41857259ba5e39b39bfd49a8b8f237429682b 100644 (file)
@@ -76,6 +76,13 @@ namespace
 
 } // namespace
 
+constexpr ndarray<unsigned int, 2, 2> ReferenceCell::line_vertex_permutations;
+
+constexpr ndarray<unsigned int, 6, 3>
+  ReferenceCell::triangle_vertex_permutations;
+
+constexpr ndarray<unsigned int, 8, 4>
+  ReferenceCell::quadrilateral_vertex_permutations;
 
 std::string
 ReferenceCell::to_string() 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.