]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Move several free functions into class ReferenceCell::Type.
authorWolfgang Bangerth <bangerth@colostate.edu>
Fri, 15 Jan 2021 16:07:33 +0000 (09:07 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Tue, 19 Jan 2021 23:49:47 +0000 (16:49 -0700)
include/deal.II/grid/connectivity.h
include/deal.II/grid/reference_cell.h
source/base/qprojector.cc
source/fe/fe_nothing.cc
source/fe/mapping_fe.cc
source/fe/mapping_fe_field.cc
source/grid/manifold.cc
source/grid/reference_cell.cc
tests/simplex/unit_tangential_vectors_01.cc

index 07b1fc376affb855d5d7acf818771969f91953c0..d049fa5f0935a39fe72a9b30350eb6d740051aba 100644 (file)
@@ -1164,9 +1164,8 @@ namespace internal
             {
               col_d[offset_i] = counter;
               orientations[offset_i] =
-                ReferenceCell::compute_orientation(ad_entity_types[offset_i],
-                                                   ad_entity_vertices[offset_i],
-                                                   ref_indices);
+                ad_entity_types[offset_i].compute_orientation(
+                  ad_entity_vertices[offset_i], ref_indices);
             }
         }
       ptr_0.push_back(col_0.size());
index b7386fed1b9b49f07d9bf640e160ab31ccd0a5bf..636e671abb808e0f14a5197c2f02aef2baa779e4 100644 (file)
@@ -20,6 +20,7 @@
 
 #include <deal.II/base/array_view.h>
 #include <deal.II/base/geometry_info.h>
+#include <deal.II/base/tensor.h>
 
 
 DEAL_II_NAMESPACE_OPEN
@@ -73,6 +74,66 @@ namespace ReferenceCell
      */
     Type(const CellKinds kind);
 
+    /**
+     * Return the dimension of the reference cell represented by the current
+     * object.
+     */
+    unsigned int
+    get_dimension() const;
+
+    /**
+     * Compute the value of the $i$-th linear shape function at location $\xi$
+     * for the current reference-cell type.
+     */
+    template <int dim>
+    double
+    d_linear_shape_function(const Point<dim> &xi, const unsigned int i) const;
+
+    /**
+     * Compute the gradient of the $i$-th linear shape function at location
+     * $\xi$ for the current reference-cell type.
+     */
+    template <int dim>
+    Tensor<1, dim>
+    d_linear_shape_function_gradient(const Point<dim> & xi,
+                                     const unsigned int i) const;
+
+    /*
+     * Return $i$-th unit tangential vector of a face of the reference cell.
+     * The vectors are arranged such that the
+     * cross product between the two vectors returns the unit normal vector.
+     *
+     * @pre $i$ must be between zero and `dim-1`.
+     */
+    template <int dim>
+    Tensor<1, dim>
+    unit_tangential_vectors(const unsigned int face_no,
+                            const unsigned int i) const;
+
+    /**
+     * Determine the orientation of the current entity described by its
+     * vertices @p var_1 relative to an entity described by @p var_0.
+     */
+    template <typename T, std::size_t N>
+    unsigned char
+    compute_orientation(const std::array<T, N> &vertices_0,
+                        const std::array<T, N> &vertices_1) const;
+
+    /**
+     * Inverse function of compute_orientation().
+     */
+    template <typename T, std::size_t N>
+    std::array<T, N>
+    permute_according_orientation(const std::array<T, N> &vertices,
+                                  const unsigned int      orientation) const;
+
+    /**
+     * Return a text representation of the reference cell represented by the
+     * current object.
+     */
+    std::string
+    to_string() const;
+
     /**
      * Conversion operator to an integer.
      */
@@ -259,13 +320,10 @@ namespace ReferenceCell
 
 
 
-  /**
-   * Return the dimension of the given reference-cell type @p type.
-   */
   inline unsigned int
-  get_dimension(const Type &type)
+  Type::get_dimension() const
   {
-    switch (type)
+    switch (kind)
       {
         case Type::Vertex:
           return 0;
@@ -284,13 +342,15 @@ namespace ReferenceCell
       }
   }
 
+
+
   /**
    * Convert the given reference cell type to a string.
    */
   inline std::string
-  to_string(const Type &type)
+  Type::to_string() const
   {
-    switch (type)
+    switch (kind)
       {
         case Type::Vertex:
           return "Vertex";
@@ -393,21 +453,17 @@ namespace ReferenceCell
   }
 
 
-  /**
-   * Compute the value of the $i$-th linear shape function at location $\xi$ for
-   * a given reference-cell type.
-   */
+
   template <int dim>
   inline double
-  d_linear_shape_function(const Type &       reference_cell,
-                          const Point<dim> & xi,
-                          const unsigned int i)
+  Type::d_linear_shape_function(const Point<dim> & xi,
+                                const unsigned int i) const
   {
-    if (reference_cell == get_hypercube(dim))
+    AssertDimension(dim, get_dimension());
+    if (*this == get_hypercube(dim))
       return GeometryInfo<dim>::d_linear_shape_function(xi, i);
 
-    if (reference_cell ==
-        Type::Tri) // see also Simplex::ScalarPolynomial::compute_value
+    if (*this == Type::Tri) // see also Simplex::ScalarPolynomial::compute_value
       {
         switch (i)
           {
@@ -420,8 +476,7 @@ namespace ReferenceCell
           }
       }
 
-    if (reference_cell ==
-        Type::Tet) // see also Simplex::ScalarPolynomial::compute_value
+    if (*this == Type::Tet) // see also Simplex::ScalarPolynomial::compute_value
       {
         switch (i)
           {
@@ -437,19 +492,18 @@ namespace ReferenceCell
           }
       }
 
-    if (reference_cell ==
+    if (*this ==
         Type::Wedge) // see also Simplex::ScalarWedgePolynomial::compute_value
       {
-        return d_linear_shape_function<2>(Type::Tri,
-                                          Point<2>(xi[std::min(0, dim - 1)],
-                                                   xi[std::min(1, dim - 1)]),
-                                          i % 3) *
-               d_linear_shape_function<1>(Type::Line,
-                                          Point<1>(xi[std::min(2, dim - 1)]),
-                                          i / 3);
+        return Type(Type::Tri).d_linear_shape_function<2>(
+                 Point<2>(xi[std::min(0, dim - 1)], xi[std::min(1, dim - 1)]),
+                 i % 3) *
+               Type(Type::Line)
+                 .d_linear_shape_function<1>(Point<1>(xi[std::min(2, dim - 1)]),
+                                             i / 3);
       }
 
-    if (reference_cell ==
+    if (*this ==
         Type::Pyramid) // see also
                        // Simplex::ScalarPyramidPolynomial::compute_value
       {
@@ -486,21 +540,18 @@ namespace ReferenceCell
     return 0.0;
   }
 
-  /**
-   * Compute the gradient of the $i$-th linear shape function at location $\xi$
-   * for a given reference-cell type.
-   */
+
+
   template <int dim>
   inline Tensor<1, dim>
-  d_linear_shape_function_gradient(const Type &       reference_cell,
-                                   const Point<dim> & xi,
-                                   const unsigned int i)
+  Type::d_linear_shape_function_gradient(const Point<dim> & xi,
+                                         const unsigned int i) const
   {
-    if (reference_cell == get_hypercube(dim))
+    AssertDimension(dim, get_dimension());
+    if (*this == get_hypercube(dim))
       return GeometryInfo<dim>::d_linear_shape_function_gradient(xi, i);
 
-    if (reference_cell ==
-        Type::Tri) // see also Simplex::ScalarPolynomial::compute_grad
+    if (*this == Type::Tri) // see also Simplex::ScalarPolynomial::compute_grad
       {
         switch (i)
           {
@@ -518,26 +569,21 @@ namespace ReferenceCell
     return Point<dim>(+0.0, +0.0, +0.0);
   }
 
-  /**
-   * Return i-th unit tangential vector of a face of the reference cell.
-   * The vectors are arranged such that the
-   * cross product between the two vectors returns the unit normal vector.
-   */
+
   template <int dim>
   inline Tensor<1, dim>
-  unit_tangential_vectors(const Type &       reference_cell,
-                          const unsigned int face_no,
-                          const unsigned int i)
+  Type::unit_tangential_vectors(const unsigned int face_no,
+                                const unsigned int i) const
   {
-    AssertDimension(dim, get_dimension(reference_cell));
+    AssertDimension(dim, get_dimension());
     AssertIndexRange(i, dim - 1);
 
-    if (reference_cell == get_hypercube(dim))
+    if (*this == get_hypercube(dim))
       {
         AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
         return GeometryInfo<dim>::unit_tangential_vectors[face_no][i];
       }
-    else if (reference_cell == Type::Tri)
+    else if (*this == Type::Tri)
       {
         AssertIndexRange(face_no, 3);
         static const std::array<Tensor<1, dim>, 3> table = {
@@ -547,7 +593,7 @@ namespace ReferenceCell
 
         return table[face_no];
       }
-    else if (reference_cell == Type::Tet)
+    else if (*this == Type::Tet)
       {
         AssertIndexRange(face_no, 4);
         static const std::array<std::array<Tensor<1, dim>, 2>, 4> table = {
@@ -563,7 +609,7 @@ namespace ReferenceCell
 
         return table[face_no][i];
       }
-    else if (reference_cell == Type::Wedge)
+    else if (*this == Type::Wedge)
       {
         AssertIndexRange(face_no, 5);
         static const std::array<std::array<Tensor<1, dim>, 2>, 5> table = {
@@ -576,7 +622,7 @@ namespace ReferenceCell
 
         return table[face_no][i];
       }
-    else if (reference_cell == Type::Pyramid)
+    else if (*this == Type::Pyramid)
       {
         AssertIndexRange(face_no, 5);
         static const std::array<std::array<Tensor<1, dim>, 2>, 5> table = {
@@ -605,7 +651,7 @@ namespace ReferenceCell
   inline Tensor<1, dim>
   unit_normal_vectors(const Type &reference_cell, const unsigned int face_no)
   {
-    AssertDimension(dim, get_dimension(reference_cell));
+    AssertDimension(dim, reference_cell.get_dimension());
 
     if (reference_cell == get_hypercube(dim))
       {
@@ -615,7 +661,7 @@ namespace ReferenceCell
     else if (dim == 2)
       {
         const auto tangential =
-          unit_tangential_vectors<dim>(reference_cell, face_no, 0);
+          reference_cell.unit_tangential_vectors<dim>(face_no, 0);
 
         Tensor<1, dim> result;
 
@@ -627,8 +673,8 @@ namespace ReferenceCell
     else if (dim == 3)
       {
         return cross_product_3d(
-          unit_tangential_vectors<dim>(reference_cell, face_no, 0),
-          unit_tangential_vectors<dim>(reference_cell, face_no, 1));
+          reference_cell.unit_tangential_vectors<dim>(face_no, 0),
+          reference_cell.unit_tangential_vectors<dim>(face_no, 1));
       }
 
     Assert(false, ExcNotImplemented());
@@ -1945,19 +1991,14 @@ namespace ReferenceCell
 
 
 
-  /**
-   * Determine the orientation of an entity of @p type described by its
-   * vertices @p var_1 relative to an entity described by @p var_0.
-   */
   template <typename T, std::size_t N>
   inline unsigned char
-  compute_orientation(const ReferenceCell::Type entity_type,
-                      const std::array<T, N> &  vertices_0,
-                      const std::array<T, N> &  vertices_1)
+  Type::compute_orientation(const std::array<T, N> &vertices_0,
+                            const std::array<T, N> &vertices_1) const
   {
     AssertIndexRange(
-      ReferenceCell::internal::Info::get_cell(entity_type).n_vertices(), N + 1);
-    if (entity_type == ReferenceCell::Type::Line)
+      ReferenceCell::internal::Info::get_cell(*this).n_vertices(), N + 1);
+    if (*this == ReferenceCell::Type::Line)
       {
         const std::array<T, 2> i{{vertices_0[0], vertices_0[1]}};
         const std::array<T, 2> j{{vertices_1[0], vertices_1[1]}};
@@ -1970,7 +2011,7 @@ namespace ReferenceCell
         if (i == std::array<T, 2>{{j[1], j[0]}})
           return 0;
       }
-    else if (entity_type == ReferenceCell::Type::Tri)
+    else if (*this == ReferenceCell::Type::Tri)
       {
         const std::array<T, 3> i{{vertices_0[0], vertices_0[1], vertices_0[2]}};
         const std::array<T, 3> j{{vertices_1[0], vertices_1[1], vertices_1[2]}};
@@ -1999,7 +2040,7 @@ namespace ReferenceCell
         if (i == std::array<T, 3>{{j[1], j[0], j[2]}})
           return 4;
       }
-    else if (entity_type == ReferenceCell::Type::Quad)
+    else if (*this == ReferenceCell::Type::Quad)
       {
         const std::array<T, 4> i{
           {vertices_0[0], vertices_0[1], vertices_0[2], vertices_0[3]}};
@@ -2039,25 +2080,22 @@ namespace ReferenceCell
           return 6;
       }
 
-    Assert(
-      false,
-      (internal::NoPermutation<T, N>(entity_type, vertices_0, vertices_1)));
+    Assert(false,
+           (internal::NoPermutation<T, N>(*this, vertices_0, vertices_1)));
 
     return -1;
   }
 
-  /**
-   * Inverse function of compute_orientation().
-   */
+
+
   template <typename T, std::size_t N>
   inline std::array<T, N>
-  permute_according_orientation(const ReferenceCell::Type entity_type,
-                                const std::array<T, N> &  vertices,
-                                const unsigned int        orientation)
+  Type::permute_according_orientation(const std::array<T, N> &vertices,
+                                      const unsigned int      orientation) const
   {
     std::array<T, 4> temp;
 
-    if (entity_type == ReferenceCell::Type::Line)
+    if (*this == ReferenceCell::Type::Line)
       {
         switch (orientation)
           {
@@ -2071,7 +2109,7 @@ namespace ReferenceCell
               Assert(false, ExcNotImplemented());
           }
       }
-    else if (entity_type == ReferenceCell::Type::Tri)
+    else if (*this == ReferenceCell::Type::Tri)
       {
         switch (orientation)
           {
@@ -2097,7 +2135,7 @@ namespace ReferenceCell
               Assert(false, ExcNotImplemented());
           }
       }
-    else if (entity_type == ReferenceCell::Type::Quad)
+    else if (*this == ReferenceCell::Type::Quad)
       {
         switch (orientation)
           {
index cdcef09bcd2540c75298da0d73b0cf8426a03a5e..85a997f514aa9bd836b984fd3796979828ef037e 100644 (file)
@@ -647,10 +647,8 @@ QProjector<3>::project_to_all_faces(
     [](const auto &face, const auto &orientation) -> std::vector<Point<3>> {
     std::array<Point<3>, 3> vertices;
     std::copy_n(face.first.begin(), face.first.size(), vertices.begin());
-    const auto temp =
-      ReferenceCell::permute_according_orientation(ReferenceCell::Type::Tri,
-                                                   vertices,
-                                                   orientation);
+    const auto temp = ReferenceCell::Type(ReferenceCell::Type::Tri)
+                        .permute_according_orientation(vertices, orientation);
     return std::vector<Point<3>>(temp.begin(),
                                  temp.begin() + face.first.size());
   };
@@ -659,10 +657,8 @@ QProjector<3>::project_to_all_faces(
     [](const auto &face, const auto &orientation) -> std::vector<Point<3>> {
     std::array<Point<3>, 4> vertices;
     std::copy_n(face.first.begin(), face.first.size(), vertices.begin());
-    const auto temp =
-      ReferenceCell::permute_according_orientation(ReferenceCell::Type::Quad,
-                                                   vertices,
-                                                   orientation);
+    const auto temp = ReferenceCell::Type(ReferenceCell::Type::Quad)
+                        .permute_according_orientation(vertices, orientation);
     return std::vector<Point<3>>(temp.begin(),
                                  temp.begin() + face.first.size());
   };
index 44009a9e489c0979d77b6d34c8387d3620bc737a..99b7e325be7c679be05e44146b17eb9c07f1bcd3 100644 (file)
@@ -76,8 +76,7 @@ FE_Nothing<dim, spacedim>::get_name() const
 
   std::vector<std::string> name_components;
   if (this->reference_cell_type() != ReferenceCell::get_hypercube(dim))
-    name_components.push_back(
-      ReferenceCell::to_string(this->reference_cell_type()));
+    name_components.push_back(this->reference_cell_type().to_string());
   if (this->n_components() > 1)
     name_components.push_back(std::to_string(this->n_components()));
   if (dominate)
index be1982c505c5785ffa69098fd47f6ab8d4ad0e4f..3b8415b4d865ba5e7fa4160dfcc93a399b1d2768 100644 (file)
@@ -158,17 +158,17 @@ MappingFE<dim, spacedim>::InternalData::initialize_face(
       for (unsigned int i = 0; i < n_faces; ++i)
         {
           unit_tangentials[i].resize(n_original_q_points);
-          std::fill(unit_tangentials[i].begin(),
-                    unit_tangentials[i].end(),
-                    ReferenceCell::unit_tangential_vectors<dim>(
-                      reference_cell_type, i, 0));
+          std::fill(
+            unit_tangentials[i].begin(),
+            unit_tangentials[i].end(),
+            reference_cell_type.template unit_tangential_vectors<dim>(i, 0));
           if (dim > 2)
             {
               unit_tangentials[n_faces + i].resize(n_original_q_points);
               std::fill(unit_tangentials[n_faces + i].begin(),
                         unit_tangentials[n_faces + i].end(),
-                        ReferenceCell::unit_tangential_vectors<dim>(
-                          reference_cell_type, i, 1));
+                        reference_cell_type
+                          .template unit_tangential_vectors<dim>(i, 1));
             }
         }
     }
@@ -871,9 +871,8 @@ MappingFE<dim, spacedim>::MappingFE(const FiniteElement<dim, spacedim> &fe)
   for (unsigned int point = 0; point < n_points; ++point)
     for (unsigned int i = 0; i < n_shape_functions; ++i)
       mapping_support_point_weights(point, i) =
-        ReferenceCell::d_linear_shape_function(reference_cell_type,
-                                               mapping_support_points[point],
-                                               i);
+        reference_cell_type.d_linear_shape_function(
+          mapping_support_points[point], i);
 }
 
 
index d440769a9a95f559a1d925d38fb059ab816cc58e..bb21ab4fa8d4e33429693db180d98720196f1dfc 100644 (file)
@@ -580,16 +580,16 @@ MappingFEField<dim, spacedim, VectorType, void>::compute_face_data(
               data.unit_tangentials[i].resize(n_original_q_points);
               std::fill(data.unit_tangentials[i].begin(),
                         data.unit_tangentials[i].end(),
-                        ReferenceCell::unit_tangential_vectors<dim>(
-                          reference_cell_type, i, 0));
+                        reference_cell_type
+                          .template unit_tangential_vectors<dim>(i, 0));
               if (dim > 2)
                 {
                   data.unit_tangentials[n_faces + i].resize(
                     n_original_q_points);
                   std::fill(data.unit_tangentials[n_faces + i].begin(),
                             data.unit_tangentials[n_faces + i].end(),
-                            ReferenceCell::unit_tangential_vectors<dim>(
-                              reference_cell_type, i, 1));
+                            reference_cell_type
+                              .template unit_tangential_vectors<dim>(i, 1));
                 }
             }
         }
index 2fe79620a4d64a217933ff472cd655385c3fb9f7..4dba26b001103160499274dfeba2fdb639255b8f 100644 (file)
@@ -913,16 +913,17 @@ FlatManifold<dim, spacedim>::normal_vector(
     {
       Point<spacedim> F;
       for (const unsigned int v : face->vertex_indices())
-        F += face->vertex(v) * ReferenceCell::d_linear_shape_function(
-                                 face_reference_cell_type, xi, v);
+        F += face->vertex(v) *
+             face_reference_cell_type.d_linear_shape_function(xi, v);
 
       for (unsigned int i = 0; i < facedim; ++i)
         {
           grad_F[i] = 0;
           for (const unsigned int v : face->vertex_indices())
             grad_F[i] +=
-              face->vertex(v) * ReferenceCell::d_linear_shape_function_gradient(
-                                  face_reference_cell_type, xi, v)[i];
+              face->vertex(v) *
+              face_reference_cell_type.d_linear_shape_function_gradient(xi,
+                                                                        v)[i];
         }
 
       Tensor<1, facedim> J;
index 8ae0afa4f4d0af99b1ed8084962d4f0581e9240a..d7b2e33f1ab849a66f524d500d0e66875f4b4074 100644 (file)
@@ -37,7 +37,7 @@ namespace ReferenceCell
   make_triangulation(const Type &                  reference_cell,
                      Triangulation<dim, spacedim> &tria)
   {
-    AssertDimension(dim, get_dimension(reference_cell));
+    AssertDimension(dim, reference_cell.get_dimension());
 
     if (reference_cell == get_hypercube(dim))
       {
@@ -110,7 +110,7 @@ namespace ReferenceCell
   std::unique_ptr<Mapping<dim, spacedim>>
   get_default_mapping(const Type &reference_cell, const unsigned int degree)
   {
-    AssertDimension(dim, get_dimension(reference_cell));
+    AssertDimension(dim, reference_cell.get_dimension());
 
     if (reference_cell == get_hypercube(dim))
       return std::make_unique<MappingQGeneric<dim, spacedim>>(degree);
@@ -136,7 +136,7 @@ namespace ReferenceCell
   const Mapping<dim, spacedim> &
   get_default_linear_mapping(const Type &reference_cell)
   {
-    AssertDimension(dim, get_dimension(reference_cell));
+    AssertDimension(dim, reference_cell.get_dimension());
 
     if (reference_cell == get_hypercube(dim))
       {
@@ -195,7 +195,7 @@ namespace ReferenceCell
   get_gauss_type_quadrature(const Type &   reference_cell,
                             const unsigned n_points_1D)
   {
-    AssertDimension(dim, get_dimension(reference_cell));
+    AssertDimension(dim, reference_cell.get_dimension());
 
     if (reference_cell == get_hypercube(dim))
       return QGauss<dim>(n_points_1D);
@@ -215,7 +215,7 @@ namespace ReferenceCell
   Quadrature<dim> &
   get_nodal_type_quadrature(const Type &reference_cell)
   {
-    AssertDimension(dim, get_dimension(reference_cell));
+    AssertDimension(dim, reference_cell.get_dimension());
 
     const auto create_quadrature = [](const Type &reference_cell) {
       Triangulation<dim> tria;
index fae80a25baf4121059fa5d86cca9ead67d5cd0e7..99f5a36ba465423688af05d51419ead4b05ca3ea 100644 (file)
@@ -33,9 +33,7 @@ test(const ReferenceCell::Type &reference_cell)
   for (const auto face_no :
        ReferenceCell::internal::Info::get_cell(reference_cell).face_indices())
     {
-      deallog << ReferenceCell::unit_normal_vectors<dim>(reference_cell,
-                                                         face_no)
-              << std::endl;
+      deallog << reference_cell.unit_normal_vectors<dim>(face_no) << std::endl;
 
       for (unsigned int i = 0; i < dim - 1; ++i)
         deallog << ReferenceCell::unit_tangential_vectors<dim>(reference_cell,

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.