]> https://gitweb.dealii.org/ - dealii.git/commitdiff
All previous issues were addressed 18576/head
authorNatalia Nebulishvili <natalia.nebulishvili@rub.de>
Sun, 22 Jun 2025 15:20:10 +0000 (17:20 +0200)
committerNatalia Nebulishvili <natalia.nebulishvili@rub.de>
Sun, 22 Jun 2025 15:20:10 +0000 (17:20 +0200)
include/deal.II/fe/fe_raviart_thomas.h
include/deal.II/matrix_free/shape_info.templates.h
source/fe/fe_raviart_thomas.cc
source/fe/fe_raviart_thomas_nodal.cc

index 46541cdd61e99a8bb06d4ba4adc2ed99eee7c0cb..153f89d29c428441c0f53d66558ed77c4a627c60 100644 (file)
@@ -148,9 +148,8 @@ public:
    * Compute the lexicographic to hierarchic numbering underlying this class,
    * necessary for the creation of the respective vector polynomial space.
    */
-  std::vector<unsigned int>
-  get_lexicographic_numbering(const unsigned int normal_degree,
-                              const unsigned int tangential_degree) const;
+  static std::vector<unsigned int>
+  get_lexicographic_numbering(const unsigned int degree);
 
   /**
    * This function returns @p true, if the shape function @p shape_index has
@@ -360,14 +359,6 @@ public:
   virtual std::unique_ptr<FiniteElement<dim, dim>>
   clone() const override;
 
-  /**
-   * Compute the lexicographic to hierarchic numbering underlying this class,
-   * necessary for the creation of the respective vector polynomial space.
-   */
-  std::vector<unsigned int>
-  get_lexicographic_numbering(const unsigned int normal_degree,
-                              const unsigned int tangential_degree) const;
-
   virtual void
   get_face_interpolation_matrix(const FiniteElement<dim> &source,
                                 FullMatrix<double>       &matrix,
index 8c852e76ffee351bce22a6f38a2985991699bc3e..8619ee71d7a2243dc9e3896d937828cb4363133b 100644 (file)
@@ -268,8 +268,8 @@ namespace internal
           const unsigned int dofs_per_face_normal = fe_in.n_dofs_per_face();
 
           lexicographic_numbering =
-            PolynomialsRaviartThomas<dim>::get_lexicographic_numbering(
-              fe_in.degree - 1);
+            FE_RaviartThomas<dim>::get_lexicographic_numbering(fe_in.degree -
+                                                               1);
 
           // To get the right shape_values of the RT element
           std::vector<unsigned int> lex_normal, lex_tangent;
index cadd78a4ac67bdd8b746af88239ef41642a08b94..6911821c4bad568c0164946f24f938d5680c3eeb 100644 (file)
@@ -45,8 +45,7 @@ FE_RaviartThomas<dim>::FE_RaviartThomas(const unsigned int deg)
   : FE_PolyTensor<dim>(
       PolynomialsVectorAnisotropic<dim>(deg + 1,
                                         deg,
-                                        get_lexicographic_numbering(deg + 1,
-                                                                    deg)),
+                                        get_lexicographic_numbering(deg)),
       FiniteElementData<dim>(get_dpo_vector(deg),
                              dim,
                              deg + 1,
@@ -700,45 +699,37 @@ FE_RaviartThomas<dim>::has_support_on_face(const unsigned int shape_index,
 
 template <int dim>
 std::vector<unsigned int>
-FE_RaviartThomas<dim>::get_lexicographic_numbering(
-  const unsigned int normal_degree,
-  const unsigned int tangential_degree) const
+FE_RaviartThomas<dim>::get_lexicographic_numbering(const unsigned int degree)
 {
-  const unsigned int n_dofs_face =
-    Utilities::pow(tangential_degree + 1, dim - 1);
+  const unsigned int        n_dofs_face = Utilities::pow(degree + 1, dim - 1);
   std::vector<unsigned int> lexicographic_numbering;
   // component 1
   for (unsigned int j = 0; j < n_dofs_face; ++j)
     {
       lexicographic_numbering.push_back(j);
-      if (normal_degree > 1)
-        for (unsigned int i = n_dofs_face * 2 * dim;
-             i < n_dofs_face * 2 * dim + normal_degree - 1;
-             ++i)
-          lexicographic_numbering.push_back(i + j * (normal_degree - 1));
+      for (unsigned int i = n_dofs_face * 2 * dim;
+           i < n_dofs_face * 2 * dim + degree;
+           ++i)
+        lexicographic_numbering.push_back(i + j * degree);
       lexicographic_numbering.push_back(n_dofs_face + j);
     }
 
   // component 2
-  unsigned int layers = (dim == 3) ? tangential_degree + 1 : 1;
+  unsigned int layers = (dim == 3) ? degree + 1 : 1;
   for (unsigned int k = 0; k < layers; ++k)
     {
-      unsigned int k_add = k * (tangential_degree + 1);
-      for (unsigned int j = n_dofs_face * 2;
-           j < n_dofs_face * 2 + tangential_degree + 1;
+      unsigned int k_add = k * (degree + 1);
+      for (unsigned int j = n_dofs_face * 2; j < n_dofs_face * 2 + degree + 1;
            ++j)
         lexicographic_numbering.push_back(j + k_add);
 
-      if (normal_degree > 1)
-        for (unsigned int i = n_dofs_face * (2 * dim + (normal_degree - 1));
-             i < n_dofs_face * (2 * dim + (normal_degree - 1)) +
-                   (normal_degree - 1) * (tangential_degree + 1);
-             ++i)
-          {
-            lexicographic_numbering.push_back(i + k_add * tangential_degree);
-          }
-      for (unsigned int j = n_dofs_face * 3;
-           j < n_dofs_face * 3 + tangential_degree + 1;
+      for (unsigned int i = n_dofs_face * (2 * dim + degree);
+           i < n_dofs_face * (2 * dim + degree) + degree * (degree + 1);
+           ++i)
+        {
+          lexicographic_numbering.push_back(i + k_add * degree);
+        }
+      for (unsigned int j = n_dofs_face * 3; j < n_dofs_face * 3 + degree + 1;
            ++j)
         lexicographic_numbering.push_back(j + k_add);
     }
@@ -748,12 +739,10 @@ FE_RaviartThomas<dim>::get_lexicographic_numbering(
     {
       for (unsigned int i = 4 * n_dofs_face; i < 5 * n_dofs_face; ++i)
         lexicographic_numbering.push_back(i);
-      if (normal_degree > 1)
-        for (unsigned int i =
-               6 * n_dofs_face + n_dofs_face * 2 * (normal_degree - 1);
-             i < 6 * n_dofs_face + n_dofs_face * 3 * (normal_degree - 1);
-             ++i)
-          lexicographic_numbering.push_back(i);
+      for (unsigned int i = 6 * n_dofs_face + n_dofs_face * 2 * degree;
+           i < 6 * n_dofs_face + n_dofs_face * 3 * degree;
+           ++i)
+        lexicographic_numbering.push_back(i);
       for (unsigned int i = 5 * n_dofs_face; i < 6 * n_dofs_face; ++i)
         lexicographic_numbering.push_back(i);
     }
index d25986923628c8671e22b72ccdccd3979f217ba0..28f3f4134214ed0ebb05845d10b7933ed6a34b2a 100644 (file)
@@ -67,10 +67,10 @@ namespace
 template <int dim>
 FE_RaviartThomasNodal<dim>::FE_RaviartThomasNodal(const unsigned int degree)
   : FE_PolyTensor<dim>(
-      PolynomialsVectorAnisotropic<dim>(degree + 1,
-                                        degree,
-                                        get_lexicographic_numbering(degree + 1,
-                                                                    degree)),
+      PolynomialsVectorAnisotropic<dim>(
+        degree + 1,
+        degree,
+        FE_RaviartThomas<dim>::get_lexicographic_numbering(degree)),
       FiniteElementData<dim>(get_rt_dpo_vector(dim, degree),
                              dim,
                              degree + 1,
@@ -85,7 +85,7 @@ FE_RaviartThomasNodal<dim>::FE_RaviartThomasNodal(const unsigned int degree)
   this->mapping_kind = {mapping_raviart_thomas};
 
   const std::vector<unsigned int> numbering =
-    get_lexicographic_numbering(degree + 1, degree);
+    FE_RaviartThomas<dim>::get_lexicographic_numbering(degree);
 
   // First, initialize the generalized support points and quadrature weights,
   // since they are required for interpolation.
@@ -242,71 +242,6 @@ FE_RaviartThomasNodal<dim>::has_support_on_face(
 
 
 
-template <int dim>
-std::vector<unsigned int>
-FE_RaviartThomasNodal<dim>::get_lexicographic_numbering(
-  const unsigned int normal_degree,
-  const unsigned int tangential_degree) const
-{
-  const unsigned int n_dofs_face =
-    Utilities::pow(tangential_degree + 1, dim - 1);
-  std::vector<unsigned int> lexicographic_numbering;
-  // component 1
-  for (unsigned int j = 0; j < n_dofs_face; ++j)
-    {
-      lexicographic_numbering.push_back(j);
-      if (normal_degree > 1)
-        for (unsigned int i = n_dofs_face * 2 * dim;
-             i < n_dofs_face * 2 * dim + normal_degree - 1;
-             ++i)
-          lexicographic_numbering.push_back(i + j * (normal_degree - 1));
-      lexicographic_numbering.push_back(n_dofs_face + j);
-    }
-
-  // component 2
-  unsigned int layers = (dim == 3) ? tangential_degree + 1 : 1;
-  for (unsigned int k = 0; k < layers; ++k)
-    {
-      unsigned int k_add = k * (tangential_degree + 1);
-      for (unsigned int j = n_dofs_face * 2;
-           j < n_dofs_face * 2 + tangential_degree + 1;
-           ++j)
-        lexicographic_numbering.push_back(j + k_add);
-
-      if (normal_degree > 1)
-        for (unsigned int i = n_dofs_face * (2 * dim + (normal_degree - 1));
-             i < n_dofs_face * (2 * dim + (normal_degree - 1)) +
-                   (normal_degree - 1) * (tangential_degree + 1);
-             ++i)
-          {
-            lexicographic_numbering.push_back(i + k_add * tangential_degree);
-          }
-      for (unsigned int j = n_dofs_face * 3;
-           j < n_dofs_face * 3 + tangential_degree + 1;
-           ++j)
-        lexicographic_numbering.push_back(j + k_add);
-    }
-
-  // component 3
-  if (dim == 3)
-    {
-      for (unsigned int i = 4 * n_dofs_face; i < 5 * n_dofs_face; ++i)
-        lexicographic_numbering.push_back(i);
-      if (normal_degree > 1)
-        for (unsigned int i =
-               6 * n_dofs_face + n_dofs_face * 2 * (normal_degree - 1);
-             i < 6 * n_dofs_face + n_dofs_face * 3 * (normal_degree - 1);
-             ++i)
-          lexicographic_numbering.push_back(i);
-      for (unsigned int i = 5 * n_dofs_face; i < 6 * n_dofs_face; ++i)
-        lexicographic_numbering.push_back(i);
-    }
-
-  return lexicographic_numbering;
-}
-
-
-
 template <int dim>
 void
 FE_RaviartThomasNodal<dim>::

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.