]> https://gitweb.dealii.org/ - dealii.git/commitdiff
QProjector: deprecate the other project_to_face() functions.
authorDavid Wells <drwells@email.unc.edu>
Sun, 9 Feb 2025 02:54:08 +0000 (21:54 -0500)
committerDavid Wells <drwells@email.unc.edu>
Sun, 9 Feb 2025 20:05:43 +0000 (15:05 -0500)
include/deal.II/base/qprojector.h
source/base/qprojector.cc

index 37e16f9d4d202ac492bee6714b8e655acd5a86b0..07643e71ab266e7e1add4f867c3ce6d1150f22b7 100644 (file)
@@ -91,7 +91,13 @@ public:
    * Compute the quadrature points on the cell if the given quadrature formula
    * is used on face <tt>face_no</tt>. For further details, see the general
    * doc for this class.
+   *
+   * @deprecated Use the version of this function which takes a
+   * combined_orientation argument instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Use the version of this function which takes a combined_orientation "
+    "argument instead.")
   static void
   project_to_face(const ReferenceCell     &reference_cell,
                   const SubQuadrature     &quadrature,
@@ -102,7 +108,13 @@ public:
    * Compute the cell quadrature formula corresponding to using
    * <tt>quadrature</tt> on face <tt>face_no</tt>. For further details, see
    * the general doc for this class.
+   *
+   * @deprecated Use the version of this function which takes a
+   * combined_orientation argument instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Use the version of this function which takes a combined_orientation "
+    "argument instead.")
   static Quadrature<dim>
   project_to_face(const ReferenceCell &reference_cell,
                   const SubQuadrature &quadrature,
@@ -113,7 +125,13 @@ public:
    * <tt>quadrature</tt> on face <tt>face_no</tt> taking into account the
    * orientation of the face. For further details, see the general doc for this
    * class.
+   *
+   * @deprecated Use the version of project_to_face() which takes a
+   * combined_orientation argument instead.
    */
+  DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+    "Use the version of project_to_face() which takes a combined_orientation "
+    "argument instead.")
   static Quadrature<dim>
   project_to_oriented_face(const ReferenceCell &reference_cell,
                            const SubQuadrature &quadrature,
@@ -128,9 +146,9 @@ public:
    * the general doc for this class.
    */
   static Quadrature<dim>
-  project_to_face(const ReferenceCell &reference_cell,
-                  const SubQuadrature &quadrature,
-                  const unsigned int   face_no,
+  project_to_face(const ReferenceCell               &reference_cell,
+                  const SubQuadrature               &quadrature,
+                  const unsigned int                 face_no,
                   const types::geometric_orientation orientation);
 
   /**
index a6a897d566884c1eeb32398444db048e2f4084b1..127352e97a5e4016cbe06106cf0a015446f6874b 100644 (file)
@@ -341,19 +341,18 @@ namespace internal
 
 template <>
 void
-QProjector<1>::project_to_face(const ReferenceCell &reference_cell,
-                               const Quadrature<0> &,
+QProjector<1>::project_to_face(const ReferenceCell   &reference_cell,
+                               const Quadrature<0>   &quadrature,
                                const unsigned int     face_no,
                                std::vector<Point<1>> &q_points)
 {
-  Assert(reference_cell == ReferenceCells::Line, ExcNotImplemented());
-  (void)reference_cell;
-
-  const unsigned int dim = 1;
-  AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
-  AssertDimension(q_points.size(), 1);
-
-  q_points[0] = Point<dim>(static_cast<double>(face_no));
+  AssertDimension(quadrature.size(), q_points.size());
+  const auto face_quadrature =
+    QProjector<1>::project_to_face(reference_cell,
+                                   quadrature,
+                                   face_no,
+                                   numbers::default_geometric_orientation);
+  q_points = face_quadrature.get_points();
 }
 
 
@@ -365,43 +364,13 @@ QProjector<2>::project_to_face(const ReferenceCell   &reference_cell,
                                const unsigned int     face_no,
                                std::vector<Point<2>> &q_points)
 {
-  AssertIndexRange(face_no, reference_cell.n_faces());
-  AssertDimension(reference_cell.get_dimension(), 2);
-  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)
-          q_points.emplace_back(quadrature.point(p)[0], 0);
-        else if (face_no == 1)
-          q_points.emplace_back(1.0 - quadrature.point(p)[0],
-                                quadrature.point(p)[0]);
-        else if (face_no == 2)
-          q_points.emplace_back(0, 1.0 - quadrature.point(p)[0]);
-        else
-          DEAL_II_ASSERT_UNREACHABLE();
-      }
-  else if (reference_cell == ReferenceCells::Quadrilateral)
-    for (unsigned int p = 0; p < quadrature.size(); ++p)
-      {
-        if (face_no == 0)
-          q_points.emplace_back(0, quadrature.point(p)[0]);
-        else if (face_no == 1)
-          q_points.emplace_back(1, quadrature.point(p)[0]);
-        else if (face_no == 2)
-          q_points.emplace_back(quadrature.point(p)[0], 0);
-        else if (face_no == 3)
-          q_points.emplace_back(quadrature.point(p)[0], 1);
-        else
-          DEAL_II_ASSERT_UNREACHABLE();
-      }
-  else
-    DEAL_II_ASSERT_UNREACHABLE();
+  AssertDimension(quadrature.size(), q_points.size());
+  const auto face_quadrature =
+    QProjector<2>::project_to_face(reference_cell,
+                                   quadrature,
+                                   face_no,
+                                   numbers::default_geometric_orientation);
+  q_points = face_quadrature.get_points();
 }
 
 
@@ -426,38 +395,23 @@ QProjector<3>::project_to_face(const ReferenceCell   &reference_cell,
 }
 
 
+
 template <int dim>
 Quadrature<dim>
 QProjector<dim>::project_to_oriented_face(const ReferenceCell &reference_cell,
                                           const Quadrature<dim - 1> &quadrature,
                                           const unsigned int         face_no,
-                                          const bool,
-                                          const bool,
-                                          const bool)
-{
-  return QProjector<dim>::project_to_face(reference_cell, quadrature, face_no);
-}
-
-
-
-template <>
-Quadrature<3>
-QProjector<3>::project_to_oriented_face(const ReferenceCell &reference_cell,
-                                        const Quadrature<2> &quadrature,
-                                        const unsigned int   face_no,
-                                        const bool           face_orientation,
-                                        const bool           face_flip,
-                                        const bool           face_rotation)
+                                          const bool face_orientation,
+                                          const bool face_flip,
+                                          const bool face_rotation)
 {
-  Assert(reference_cell == ReferenceCells::Hexahedron, ExcNotImplemented());
-
-  const Quadrature<2> mutation = internal::QProjector::mutate_quadrature(
+  return QProjector<dim>::project_to_face(
+    reference_cell,
     quadrature,
+    face_no,
     internal::combined_face_orientation(face_orientation,
                                         face_rotation,
                                         face_flip));
-
-  return QProjector<3>::project_to_face(reference_cell, mutation, face_no);
 }
 
 
@@ -470,7 +424,9 @@ QProjector<dim>::project_to_face(const ReferenceCell       &reference_cell,
                                  const types::geometric_orientation orientation)
 {
   AssertIndexRange(face_no, reference_cell.n_faces());
-  AssertIndexRange(orientation, reference_cell.n_face_orientations(face_no));
+  // TODO once the default orientation is zero we can remove this special case
+  AssertIndexRange(dim == 1 ? 1 - orientation : orientation,
+                   reference_cell.n_face_orientations(face_no));
   AssertDimension(reference_cell.get_dimension(), dim);
 
   std::vector<Point<dim>> q_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.