]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move some more functions from internal::ReferenceCell::* classes to ReferenceCell.
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 10 Feb 2021 04:18:29 +0000 (21:18 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 11 Feb 2021 20:10:55 +0000 (13:10 -0700)
include/deal.II/grid/reference_cell.h
include/deal.II/grid/tria_accessor.templates.h
source/fe/fe.cc
source/grid/tria_accessor.cc

index ef703d929ec4ad2de2a2a2892cddc7526e4fef92..53f4e9ccad87ac756a6940176a4ad36fc434b523 100644 (file)
@@ -281,6 +281,42 @@ public:
    * @{
    */
 
+  /**
+   * Return which child cells are adjacent to a certain face of the
+   * mother cell.
+   *
+   * For example, in 2D the layout of a quadrilateral cell is as follows:
+   * @verbatim
+   * .      3
+   * .   2-->--3
+   * .   |     |
+   * . 0 ^     ^ 1
+   * .   |     |
+   * .   0-->--1
+   * .      2
+   * @endverbatim
+   * Vertices and faces are indicated with their numbers, faces also with
+   * their directions.
+   *
+   * Now, when refined, the layout is like this:
+   * @verbatim
+   * *--*--*
+   * | 2|3 |
+   * *--*--*
+   * | 0|1 |
+   * *--*--*
+   * @endverbatim
+   *
+   * Thus, the child cells on face 0 are (ordered in the direction of the
+   * face) 0 and 2, on face 3 they are 2 and 3, etc.
+   *
+   * For three spatial dimensions, the exact order of the children is laid
+   * down in the general documentation of this class.
+   */
+  unsigned int
+  child_cell_on_face(const unsigned int face_n,
+                     const unsigned int subface_n) const;
+
   /**
    * For a given vertex in a cell, return a pair of a face index and a
    * vertex index within this face.
@@ -305,6 +341,38 @@ public:
   std::array<unsigned int, 2>
   standard_line_to_face_and_line_index(const unsigned int line) const;
 
+  /**
+   * Map face line number to cell line number.
+   */
+  unsigned int
+  face_to_cell_lines(const unsigned int  face,
+                     const unsigned int  line,
+                     const unsigned char face_orientation) const;
+
+  /**
+   * Map face vertex number to cell vertex number.
+   */
+  unsigned int
+  face_to_cell_vertices(const unsigned int  face,
+                        const unsigned int  vertex,
+                        const unsigned char face_orientation) const;
+
+  /**
+   * Correct vertex index depending on face orientation.
+   */
+  unsigned int
+  standard_to_real_face_vertex(const unsigned int  vertex,
+                               const unsigned int  face,
+                               const unsigned char face_orientation) const;
+
+  /**
+   * Correct line index depending on face orientation.
+   */
+  unsigned int
+  standard_to_real_face_line(const unsigned int  line,
+                             const unsigned int  face,
+                             const unsigned char face_orientation) const;
+
   /**
    * @}
    */
@@ -793,6 +861,53 @@ ReferenceCell::face_reference_cell(const unsigned int face_no) const
 
 
 
+inline unsigned int
+ReferenceCell::child_cell_on_face(const unsigned int face,
+                                  const unsigned int subface) const
+{
+  AssertIndexRange(face, n_faces());
+
+  if (*this == ReferenceCells::Vertex)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Line)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Triangle)
+    {
+      static constexpr unsigned int subcells[3][2] = {{0, 1}, {1, 2}, {2, 0}};
+
+      return subcells[face][subface];
+    }
+  else if (*this == ReferenceCells::Quadrilateral)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Tetrahedron)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Pyramid)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Wedge)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Hexahedron)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+
+  Assert(false, ExcNotImplemented());
+  return {};
+}
+
+
+
 inline std::array<unsigned int, 2>
 ReferenceCell::standard_vertex_to_face_and_vertex_index(
   const unsigned int vertex) const
@@ -919,6 +1034,323 @@ ReferenceCell::standard_line_to_face_and_line_index(
 
 
 
+inline unsigned int
+ReferenceCell::face_to_cell_lines(const unsigned int  face,
+                                  const unsigned int  line,
+                                  const unsigned char face_orientation) const
+{
+  AssertIndexRange(face, n_faces());
+  AssertIndexRange(line, face_reference_cell(face).n_lines());
+
+  if (*this == ReferenceCells::Vertex)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Line)
+    {
+      return GeometryInfo<1>::face_to_cell_lines(
+        face,
+        line,
+        Utilities::get_bit(face_orientation, 0),
+        Utilities::get_bit(face_orientation, 2),
+        Utilities::get_bit(face_orientation, 1));
+    }
+  else if (*this == ReferenceCells::Triangle)
+    {
+      return face;
+    }
+  else if (*this == ReferenceCells::Quadrilateral)
+    {
+      return GeometryInfo<2>::face_to_cell_lines(
+        face,
+        line,
+        Utilities::get_bit(face_orientation, 0),
+        Utilities::get_bit(face_orientation, 2),
+        Utilities::get_bit(face_orientation, 1));
+    }
+  else if (*this == ReferenceCells::Tetrahedron)
+    {
+      const static unsigned int table[4][3] = {{0, 1, 2},
+                                               {0, 3, 4},
+                                               {2, 5, 3},
+                                               {1, 4, 5}};
+
+      return table[face]
+                  [standard_to_real_face_line(line, face, face_orientation)];
+    }
+  else if (*this == ReferenceCells::Pyramid)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Wedge)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Hexahedron)
+    {
+      return GeometryInfo<3>::face_to_cell_lines(
+        face,
+        line,
+        Utilities::get_bit(face_orientation, 0),
+        Utilities::get_bit(face_orientation, 2),
+        Utilities::get_bit(face_orientation, 1));
+    }
+
+  Assert(false, ExcNotImplemented());
+  return numbers::invalid_unsigned_int;
+}
+
+
+
+inline unsigned int
+ReferenceCell::face_to_cell_vertices(const unsigned int  face,
+                                     const unsigned int  vertex,
+                                     const unsigned char face_orientation) const
+{
+  AssertIndexRange(face, n_faces());
+  AssertIndexRange(vertex, face_reference_cell(face).n_vertices());
+
+  if (*this == ReferenceCells::Vertex)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Line)
+    {
+      return GeometryInfo<1>::face_to_cell_vertices(
+        face,
+        vertex,
+        Utilities::get_bit(face_orientation, 0),
+        Utilities::get_bit(face_orientation, 2),
+        Utilities::get_bit(face_orientation, 1));
+    }
+  else if (*this == ReferenceCells::Triangle)
+    {
+      static const unsigned int table[3][2] = {{0, 1}, {1, 2}, {2, 0}};
+
+      return table[face][face_orientation ? vertex : (1 - vertex)];
+    }
+  else if (*this == ReferenceCells::Quadrilateral)
+    {
+      return GeometryInfo<2>::face_to_cell_vertices(
+        face,
+        vertex,
+        Utilities::get_bit(face_orientation, 0),
+        Utilities::get_bit(face_orientation, 2),
+        Utilities::get_bit(face_orientation, 1));
+    }
+  else if (*this == ReferenceCells::Tetrahedron)
+    {
+      static const unsigned int table[4][3] = {{0, 1, 2},
+                                               {1, 0, 3},
+                                               {0, 2, 3},
+                                               {2, 1, 3}};
+
+      return table[face][standard_to_real_face_vertex(
+        vertex, face, face_orientation)];
+    }
+  else if (*this == ReferenceCells::Pyramid)
+    {
+      constexpr auto            X           = numbers::invalid_unsigned_int;
+      static const unsigned int table[5][4] = {
+        {0, 1, 2, 3}, {0, 2, 4, X}, {3, 1, 4, X}, {1, 0, 4, X}, {2, 3, 4, X}};
+
+      return table[face][standard_to_real_face_vertex(
+        vertex, face, face_orientation)];
+    }
+  else if (*this == ReferenceCells::Wedge)
+    {
+      constexpr auto            X           = numbers::invalid_unsigned_int;
+      static const unsigned int table[6][4] = {
+        {1, 0, 2, X}, {3, 4, 5, X}, {0, 1, 3, 4}, {1, 2, 4, 5}, {2, 0, 5, 3}};
+
+      return table[face][standard_to_real_face_vertex(
+        vertex, face, face_orientation)];
+    }
+  else if (*this == ReferenceCells::Hexahedron)
+    {
+      return GeometryInfo<3>::face_to_cell_vertices(
+        face,
+        vertex,
+        Utilities::get_bit(face_orientation, 0),
+        Utilities::get_bit(face_orientation, 2),
+        Utilities::get_bit(face_orientation, 1));
+    }
+
+  Assert(false, ExcNotImplemented());
+  return numbers::invalid_unsigned_int;
+}
+
+
+
+inline unsigned int
+ReferenceCell::standard_to_real_face_vertex(
+  const unsigned int  vertex,
+  const unsigned int  face,
+  const unsigned char face_orientation) const
+{
+  AssertIndexRange(face, n_faces());
+  AssertIndexRange(vertex, face_reference_cell(face).n_vertices());
+
+  if (*this == ReferenceCells::Vertex)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Line)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Triangle)
+    {
+      static const unsigned int table[2][2] = {{1, 0}, {0, 1}};
+
+      return table[face_orientation][vertex];
+    }
+  else if (*this == ReferenceCells::Quadrilateral)
+    {
+      return GeometryInfo<2>::standard_to_real_line_vertex(vertex,
+                                                           face_orientation);
+    }
+  else if (*this == ReferenceCells::Tetrahedron)
+    {
+      static const unsigned int table[6][3] = {
+        {0, 2, 1}, {0, 1, 2}, {2, 1, 0}, {1, 2, 0}, {1, 0, 2}, {2, 0, 1}};
+
+      return table[face_orientation][vertex];
+    }
+  else if (*this == ReferenceCells::Pyramid)
+    {
+      if (face == 0) // Quad
+        {
+          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 // Tri
+        {
+          static const unsigned int table[6][3] = {
+            {0, 2, 1}, {0, 1, 2}, {2, 1, 0}, {1, 2, 0}, {1, 0, 2}, {2, 0, 1}};
+
+          return table[face_orientation][vertex];
+        }
+    }
+  else if (*this == ReferenceCells::Wedge)
+    {
+      if (face > 1) // QUAD
+        {
+          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 // TRI
+        {
+          static const unsigned int table[6][3] = {
+            {0, 2, 1}, {0, 1, 2}, {2, 1, 0}, {1, 2, 0}, {1, 0, 2}, {2, 0, 1}};
+
+          return table[face_orientation][vertex];
+        }
+    }
+  else if (*this == ReferenceCells::Hexahedron)
+    {
+      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));
+    }
+
+  Assert(false, ExcNotImplemented());
+  return numbers::invalid_unsigned_int;
+}
+
+
+
+inline unsigned int
+ReferenceCell::standard_to_real_face_line(
+  const unsigned int  line,
+  const unsigned int  face,
+  const unsigned char face_orientation) const
+{
+  AssertIndexRange(face, n_faces());
+  AssertIndexRange(line, face_reference_cell(face).n_lines());
+
+  if (*this == ReferenceCells::Vertex)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Line)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Triangle)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Quadrilateral)
+    {
+      Assert(false, ExcNotImplemented());
+    }
+  else if (*this == ReferenceCells::Tetrahedron)
+    {
+      static const unsigned int table[6][3] = {
+        {2, 1, 0}, {0, 1, 2}, {1, 0, 2}, {1, 2, 0}, {0, 2, 1}, {2, 0, 1}};
+
+      return table[face_orientation][line];
+    }
+  else if (*this == ReferenceCells::Pyramid)
+    {
+      if (face == 0) // QUAD
+        {
+          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 // TRI
+        {
+          static const unsigned int table[6][3] = {
+            {2, 1, 0}, {0, 1, 2}, {1, 0, 2}, {1, 2, 0}, {0, 2, 1}, {2, 0, 1}};
+
+          return table[face_orientation][line];
+        }
+    }
+  else if (*this == ReferenceCells::Wedge)
+    {
+      if (face > 1) // QUAD
+        {
+          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 // TRI
+        {
+          static const unsigned int table[6][3] = {
+            {2, 1, 0}, {0, 1, 2}, {1, 0, 2}, {1, 2, 0}, {0, 2, 1}, {2, 0, 1}};
+
+          return table[face_orientation][line];
+        }
+    }
+  else if (*this == ReferenceCells::Hexahedron)
+    {
+      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));
+    }
+
+  Assert(false, ExcNotImplemented());
+  return numbers::invalid_unsigned_int;
+}
+
+
+
 namespace ReferenceCells
 {
   template <int dim>
@@ -1345,7 +1777,6 @@ namespace internal
         return {{0, 0}};
       }
 
-    public:
       /**
        * Correct vertex index depending on face orientation.
        */
@@ -1380,6 +1811,7 @@ namespace internal
         return 0;
       }
 
+    public:
       /**
        * Combine face and line orientation.
        */
@@ -1411,7 +1843,6 @@ namespace internal
         return ReferenceCells::Invalid;
       }
 
-    public:
       /**
        * Map face line number to cell line number.
        */
@@ -1444,6 +1875,7 @@ namespace internal
         return 0;
       }
 
+    public:
       /**
        * Map an ExodusII vertex number to a deal.II vertex number.
        */
@@ -1468,6 +1900,7 @@ namespace internal
         return 0;
       }
 
+    private:
       /**
        * Indices of child cells that are adjacent to a certain face of the
        * mother cell.
index 8d2d26fc2d858362f7663fdb907f74fecf4b199a..b579ba684f6e8c51f50619d1d03bb103ec6a1afc 100644 (file)
@@ -624,7 +624,7 @@ namespace internal
           accessor.reference_cell().standard_line_to_face_and_line_index(i);
         const auto quad_index = pair[0];
         const auto line_index =
-          accessor.reference_cell_info().standard_to_real_face_line(
+          accessor.reference_cell().standard_to_real_face_line(
             pair[1], pair[0], face_orientation_raw(accessor, quad_index));
 
         return accessor.quad(quad_index)->line_index(line_index);
@@ -836,7 +836,7 @@ namespace internal
           accessor.reference_cell().standard_line_to_face_and_line_index(line);
         const auto quad_index = pair[0];
         const auto line_index =
-          accessor.reference_cell_info().standard_to_real_face_line(
+          accessor.reference_cell().standard_to_real_face_line(
             pair[1], pair[0], face_orientation_raw(accessor, quad_index));
 
         return accessor.reference_cell_info().combine_face_and_line_orientation(
@@ -1022,7 +1022,7 @@ namespace internal
             corner);
         const auto line_index = pair[0];
         const auto vertex_index =
-          accessor.reference_cell_info().standard_to_real_face_vertex(
+          accessor.reference_cell().standard_to_real_face_vertex(
             pair[1], pair[0], accessor.line_orientation(line_index));
 
         return accessor.line(line_index)->vertex_index(vertex_index);
@@ -1039,7 +1039,7 @@ namespace internal
             corner);
         const auto face_index = pair[0];
         const auto vertex_index =
-          accessor.reference_cell_info().standard_to_real_face_vertex(
+          accessor.reference_cell().standard_to_real_face_vertex(
             pair[1], pair[0], face_orientation_raw(accessor, face_index));
 
         return accessor.quad(face_index)->vertex_index(vertex_index);
index 9d40d929e33300333ee684e2710c5bdc25494509..610e3afbb9572c479c3fe08be4f75ca5fa9b7994 100644 (file)
@@ -571,9 +571,6 @@ FiniteElement<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
                                                  const bool face_flip,
                                                  const bool face_rotation) const
 {
-  const auto &refence_cell =
-    internal::ReferenceCell::get_cell(this->reference_cell());
-
   AssertIndexRange(face_index, this->n_dofs_per_face(face));
   AssertIndexRange(face, this->reference_cell().n_faces());
 
@@ -612,11 +609,11 @@ 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 (refence_cell.face_to_cell_vertices(face,
-                                                 face_vertex,
-                                                 face_orientation +
-                                                   2 * face_rotation +
-                                                   4 * face_flip) *
+      return (this->reference_cell().face_to_cell_vertices(face,
+                                                           face_vertex,
+                                                           face_orientation +
+                                                             2 * face_rotation +
+                                                             4 * face_flip) *
                 this->n_dofs_per_vertex() +
               dof_index_on_vertex);
     }
@@ -632,11 +629,11 @@ 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() +
-              refence_cell.face_to_cell_lines(face,
-                                              face_line,
-                                              face_orientation +
-                                                2 * face_rotation +
-                                                4 * face_flip) *
+              this->reference_cell().face_to_cell_lines(face,
+                                                        face_line,
+                                                        face_orientation +
+                                                          2 * face_rotation +
+                                                          4 * face_flip) *
                 this->n_dofs_per_line() +
               dof_index_on_line);
     }
index 29663e79991c74cd5c1e5966440f59b9c5244dda..6f4e270ac5a0bcb21aec750e4e60963cbeadbb4c 100644 (file)
@@ -2868,11 +2868,10 @@ CellAccessor<dim, spacedim>::neighbor_child_on_subface(
                   (1 - subface) :
                   subface;
 
-              const auto &info =
-                internal::ReferenceCell::get_cell(ReferenceCells::Triangle);
               const unsigned int neighbor_child_index =
-                info.child_cell_on_face(neighbor_face, neighbor_subface);
-              TriaIterator<CellAccessor<dim, spacedim>> sub_neighbor =
+                ReferenceCells::Triangle.child_cell_on_face(neighbor_face,
+                                                            neighbor_subface);
+              const TriaIterator<CellAccessor<dim, spacedim>> sub_neighbor =
                 neighbor_cell->child(neighbor_child_index);
 
               // neighbor's child is not allowed to be further refined for the

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.