]> https://gitweb.dealii.org/ - dealii.git/commitdiff
QProjector: deprecate the other project_to_subface() functions.
authorDavid Wells <drwells@email.unc.edu>
Sat, 15 Feb 2025 16:43:03 +0000 (11:43 -0500)
committerDavid Wells <drwells@email.unc.edu>
Sat, 15 Feb 2025 20:45:24 +0000 (15:45 -0500)
include/deal.II/base/qprojector.h
source/base/qprojector.cc

index 7aa720229ebeb8eb6b893b8ca68c762f2283d4bc..760c30e5c083107f1afe2cc3645b8a6dbe62b1b7 100644 (file)
@@ -159,7 +159,13 @@ public:
    *
    * @note Only the points are transformed. The quadrature weights are the
    * same as those of the original rule.
+   *
+   * @deprecated Use the version of project_to_face() which takes an
+   * orientation argument instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Use the version of project_to_subface() which takes an orientation "
+    "argument instead.")
   static void
   project_to_subface(const ReferenceCell           &reference_cell,
                      const SubQuadrature           &quadrature,
@@ -181,7 +187,13 @@ public:
    * @note This function is deprecated since it makes an implicit assumption
    * that the cell is a line (1D), a quad (2d), or a hex (3d). Use the other
    * version of this function that takes the reference cell type instead.
+   *
+   * @deprecated Use the version of project_to_face() which takes an
+   * orientation argument instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Use the version of project_to_subface() which takes an orientation "
+    "argument instead.")
   static Quadrature<dim>
   project_to_subface(const ReferenceCell           &reference_cell,
                      const SubQuadrature           &quadrature,
@@ -198,7 +210,13 @@ public:
    *
    * @note Only the points are transformed. The quadrature weights are the
    * same as those of the original rule.
+   *
+   * @deprecated Use the version of project_to_face() which takes an
+   * orientation argument instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Use the version of project_to_subface() which takes an orientation "
+    "argument instead.")
   static Quadrature<dim>
   project_to_oriented_subface(const ReferenceCell             &reference_cell,
                               const SubQuadrature             &quadrature,
@@ -217,7 +235,13 @@ public:
    *
    * @note Only the points are transformed. The quadrature weights are the
    * same as those of the original rule.
+   *
+   * @deprecated Use the version of project_to_face() which takes an
+   * orientation argument instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Use the version of project_to_subface() which takes an orientation "
+    "argument instead.")
   static Quadrature<dim>
   project_to_subface(const ReferenceCell               &reference_cell,
                      const SubQuadrature               &quadrature,
index 1af849c4b420f094536e32a70c798759d2085c43..ee9ed6aac17ac632fca9461553930fc2aa7b564f 100644 (file)
@@ -520,94 +520,22 @@ QProjector<1>::project_to_subface(const ReferenceCell &reference_cell,
 
 template <>
 void
-QProjector<2>::project_to_subface(const ReferenceCell   &reference_cell,
-                                  const Quadrature<1>   &quadrature,
-                                  const unsigned int     face_no,
-                                  const unsigned int     subface_no,
-                                  std::vector<Point<2>> &q_points,
-                                  const RefinementCase<1> &)
+QProjector<2>::project_to_subface(const ReferenceCell     &reference_cell,
+                                  const Quadrature<1>     &quadrature,
+                                  const unsigned int       face_no,
+                                  const unsigned int       subface_no,
+                                  std::vector<Point<2>>   &q_points,
+                                  const RefinementCase<1> &ref_case)
 {
-  const unsigned int dim = 2;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  AssertIndexRange(subface_no, GeometryInfo<dim>::max_children_per_face);
-
-  Assert(q_points.size() == quadrature.size(),
-         ExcDimensionMismatch(q_points.size(), quadrature.size()));
-
-  q_points.resize(0);
-  if (reference_cell == ReferenceCells::Triangle)
-    {
-      // use linear polynomial to map the reference quadrature points correctly
-      // on faces, i.e., BarycentricPolynomials<1>(1)
-      for (unsigned int p = 0; p < quadrature.size(); ++p)
-        {
-          if (face_no == 0)
-            {
-              if (subface_no == 0)
-                q_points.emplace_back(quadrature.point(p)[0] / 2, 0);
-              else
-                q_points.emplace_back(0.5 + quadrature.point(p)[0] / 2, 0);
-            }
-          else if (face_no == 1)
-            {
-              if (subface_no == 0)
-                q_points.emplace_back(1 - quadrature.point(p)[0] / 2,
-                                      quadrature.point(p)[0] / 2);
-              else
-                q_points.emplace_back(0.5 - quadrature.point(p)[0] / 2,
-                                      0.5 + quadrature.point(p)[0] / 2);
-            }
-          else if (face_no == 2)
-            {
-              if (subface_no == 0)
-                q_points.emplace_back(0, 1 - quadrature.point(p)[0] / 2);
-              else
-                q_points.emplace_back(0, 0.5 - quadrature.point(p)[0] / 2);
-            }
-          else
-            DEAL_II_ASSERT_UNREACHABLE();
-        }
-    }
-  else if (reference_cell == ReferenceCells::Quadrilateral)
-    {
-      for (unsigned int p = 0; p < quadrature.size(); ++p)
-        {
-          if (face_no == 0)
-            {
-              if (subface_no == 0)
-                q_points.emplace_back(0, quadrature.point(p)[0] / 2);
-              else
-                q_points.emplace_back(0, quadrature.point(p)[0] / 2 + 0.5);
-            }
-          else if (face_no == 1)
-            {
-              if (subface_no == 0)
-                q_points.emplace_back(1, quadrature.point(p)[0] / 2);
-              else
-                q_points.emplace_back(1, quadrature.point(p)[0] / 2 + 0.5);
-            }
-          else if (face_no == 2)
-            {
-              if (subface_no == 0)
-                q_points.emplace_back(quadrature.point(p)[0] / 2, 0);
-              else
-                q_points.emplace_back(quadrature.point(p)[0] / 2 + 0.5, 0);
-            }
-          else if (face_no == 3)
-            {
-              if (subface_no == 0)
-                q_points.emplace_back(quadrature.point(p)[0] / 2, 1);
-              else
-                q_points.emplace_back(quadrature.point(p)[0] / 2 + 0.5, 1);
-            }
-          else
-            DEAL_II_ASSERT_UNREACHABLE();
-        }
-    }
-  else
-    {
-      DEAL_II_ASSERT_UNREACHABLE();
-    }
+  AssertDimension(quadrature.size(), q_points.size());
+  const auto face_quadrature =
+    project_to_subface(reference_cell,
+                       quadrature,
+                       face_no,
+                       subface_no,
+                       numbers::default_geometric_orientation,
+                       ref_case);
+  q_points = face_quadrature.get_points();
 }
 
 
@@ -623,15 +551,15 @@ QProjector<3>::project_to_subface(const ReferenceCell     &reference_cell,
 {
   Assert(reference_cell == ReferenceCells::Hexahedron, ExcNotImplemented());
   (void)reference_cell;
-
-  AssertIndexRange(face_no, GeometryInfo<3>::faces_per_cell);
-  AssertIndexRange(subface_no, GeometryInfo<3>::max_children_per_face);
-  Assert(q_points.size() == quadrature.size(),
-         ExcDimensionMismatch(q_points.size(), quadrature.size()));
-
-  q_points.clear();
-  internal::QProjector::project_to_hex_face_and_append(
-    quadrature.get_points(), face_no, q_points, ref_case, subface_no);
+  AssertDimension(quadrature.size(), q_points.size());
+  const auto face_quadrature =
+    project_to_subface(reference_cell,
+                       quadrature,
+                       face_no,
+                       subface_no,
+                       numbers::default_geometric_orientation,
+                       ref_case);
+  q_points = face_quadrature.get_points();
 }
 
 

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.