]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Avoid calling the deprecated version of project_to_subface(). 18140/head
authorDavid Wells <drwells@email.unc.edu>
Sat, 15 Feb 2025 16:37:59 +0000 (11:37 -0500)
committerDavid Wells <drwells@email.unc.edu>
Sat, 15 Feb 2025 20:45:24 +0000 (15:45 -0500)
include/deal.II/fe/fe_tools.templates.h
source/fe/fe_abf.cc
source/fe/fe_bernstein.cc
source/fe/fe_p1nc.cc
source/fe/fe_q_base.cc
source/fe/fe_raviart_thomas.cc
source/fe/fe_raviart_thomas_nodal.cc
source/fe/fe_simplex_p.cc

index 15431416f3bcdec33c8eb83ea0ebe3b29e2e89f4..77442f4897366a0e44000f367f26ae410cf5c731 100644 (file)
@@ -1968,7 +1968,12 @@ namespace FETools
          ++cell_number)
       {
         const Quadrature<dim> q_coarse = QProjector<dim>::project_to_subface(
-          fe.reference_cell(), q_gauss, face_coarse, cell_number);
+          fe.reference_cell(),
+          q_gauss,
+          face_coarse,
+          cell_number,
+          numbers::default_geometric_orientation,
+          RefinementCase<dim - 1>::isotropic_refinement);
         FEValues<dim> coarse(mapping, fe, q_coarse, update_values);
 
         typename Triangulation<dim, spacedim>::active_cell_iterator fine_cell =
index c93890a486067cddf54e95dde797556b77e3cca3..be10e3f6fc947b54c8c951bce4b6a41a9e1dc3f6 100644 (file)
@@ -382,7 +382,12 @@ FE_ABF<dim>::initialize_restriction()
           // evaluated on the subface
           // only.
           Quadrature<dim> q_sub = QProjector<dim>::project_to_subface(
-            this->reference_cell(), q_base, face, sub);
+            this->reference_cell(),
+            q_base,
+            face,
+            sub,
+            numbers::default_geometric_orientation,
+            RefinementCase<dim - 1>::isotropic_refinement);
           const unsigned int child = GeometryInfo<dim>::child_cell_on_face(
             RefinementCase<dim>::isotropic_refinement, face, sub);
 
index 145ffb59836542c9e1dd1e58044935a630a205bb..02764a88d9b64cf5b2a322ae09d4f8e7e6d9b393 100644 (file)
@@ -158,10 +158,13 @@ FE_Bernstein<dim, spacedim>::get_subface_interpolation_matrix(
             quad_face_support,
             0,
             numbers::default_geometric_orientation) :
-          QProjector<dim>::project_to_subface(this->reference_cell(),
-                                              quad_face_support,
-                                              0,
-                                              subface);
+          QProjector<dim>::project_to_subface(
+            this->reference_cell(),
+            quad_face_support,
+            0,
+            subface,
+            numbers::default_geometric_orientation,
+            RefinementCase<dim - 1>::isotropic_refinement);
 
       for (unsigned int i = 0; i < source_fe->n_dofs_per_face(face_no); ++i)
         {
index 212ca1ef4ab1696a5e94e67adbbc34c85619862e..6162aa925ae4de37dd3c03899d29b060923c98e2 100644 (file)
@@ -309,8 +309,13 @@ FE_P1NC::fill_fe_subface_values(
   const ndarray<double, 4, 3> coeffs = get_linear_shape_coefficients(cell);
 
   // compute on the subface
-  const Quadrature<2> quadrature_on_subface = QProjector<2>::project_to_subface(
-    this->reference_cell(), quadrature, face_no, sub_no);
+  const Quadrature<2> quadrature_on_subface =
+    QProjector<2>::project_to_subface(this->reference_cell(),
+                                      quadrature,
+                                      face_no,
+                                      sub_no,
+                                      numbers::default_geometric_orientation,
+                                      RefinementCase<1>::isotropic_refinement);
 
   if (flags & update_values)
     for (unsigned int i = 0; i < quadrature_on_subface.size(); ++i)
index 70cdd53d90e7241cababc160d008c597db063795..0084fa5e75cb6e12387e0332f15156750b5c9b4e 100644 (file)
@@ -247,38 +247,38 @@ struct FE_Q_Base<xdim, xspacedim>::Implementation
 
     if (q_deg > 1)
       {
-        const unsigned int          n    = q_deg - 1;
-        const double                step = 1. / q_deg;
-        std::vector<Point<dim - 2>> line_support_points(n);
+        const unsigned int    n    = q_deg - 1;
+        const double          step = 1. / q_deg;
+        std::vector<Point<1>> line_support_points(n);
         for (unsigned int i = 0; i < n; ++i)
           line_support_points[i][0] = (i + 1) * step;
-        const Quadrature<dim - 2> qline(line_support_points);
-
-        // auxiliary points in 2d
-        std::vector<Point<dim - 1>> p_line(n);
+        const Quadrature<1> qline(line_support_points);
 
         // Add nodes of lines interior in the "mother-face"
+        auto get_points = [&](const unsigned int face_no,
+                              const unsigned int subface_no) {
+          return QProjector<2>::project_to_subface(
+                   ReferenceCells::get_hypercube<2>(),
+                   qline,
+                   face_no,
+                   subface_no,
+                   numbers::default_geometric_orientation,
+                   RefinementCase<1>::isotropic_refinement)
+            .get_points();
+        };
 
         // line 5: use line 9
-        QProjector<dim - 1>::project_to_subface(
-          ReferenceCells::get_hypercube<dim - 1>(), qline, 0, 0, p_line);
-        for (unsigned int i = 0; i < n; ++i)
-          constraint_points.push_back(p_line[i] + Point<dim - 1>(0.5, 0));
+        for (const Point<dim - 1> &p : get_points(0, 0))
+          constraint_points.push_back(p + Point<dim - 1>(0.5, 0));
         // line 6: use line 10
-        QProjector<dim - 1>::project_to_subface(
-          ReferenceCells::get_hypercube<dim - 1>(), qline, 0, 1, p_line);
-        for (unsigned int i = 0; i < n; ++i)
-          constraint_points.push_back(p_line[i] + Point<dim - 1>(0.5, 0));
+        for (const Point<dim - 1> &p : get_points(0, 1))
+          constraint_points.push_back(p + Point<dim - 1>(0.5, 0));
         // line 7: use line 13
-        QProjector<dim - 1>::project_to_subface(
-          ReferenceCells::get_hypercube<dim - 1>(), qline, 2, 0, p_line);
-        for (unsigned int i = 0; i < n; ++i)
-          constraint_points.push_back(p_line[i] + Point<dim - 1>(0, 0.5));
+        for (const Point<dim - 1> &p : get_points(2, 0))
+          constraint_points.push_back(p + Point<dim - 1>(0, 0.5));
         // line 8: use line 14
-        QProjector<dim - 1>::project_to_subface(
-          ReferenceCells::get_hypercube<dim - 1>(), qline, 2, 1, p_line);
-        for (unsigned int i = 0; i < n; ++i)
-          constraint_points.push_back(p_line[i] + Point<dim - 1>(0, 0.5));
+        for (const Point<dim - 1> &p : get_points(2, 1))
+          constraint_points.push_back(p + Point<dim - 1>(0, 0.5));
 
         // DoFs on bordering lines lines 9-16
         for (unsigned int face = 0;
@@ -288,12 +288,7 @@ struct FE_Q_Base<xdim, xspacedim>::Implementation
                subface < GeometryInfo<dim - 1>::max_children_per_face;
                ++subface)
             {
-              QProjector<dim - 1>::project_to_subface(
-                ReferenceCells::get_hypercube<dim - 1>(),
-                qline,
-                face,
-                subface,
-                p_line);
+              const auto p_line = get_points(face, subface);
               constraint_points.insert(constraint_points.end(),
                                        p_line.begin(),
                                        p_line.end());
@@ -672,10 +667,13 @@ FE_Q_Base<dim, spacedim>::get_subface_interpolation_matrix(
             quad_face_support,
             0,
             numbers::default_geometric_orientation) :
-          QProjector<dim>::project_to_subface(this->reference_cell(),
-                                              quad_face_support,
-                                              0,
-                                              subface);
+          QProjector<dim>::project_to_subface(
+            this->reference_cell(),
+            quad_face_support,
+            0,
+            subface,
+            numbers::default_geometric_orientation,
+            RefinementCase<dim - 1>::isotropic_refinement);
       for (unsigned int i = 0; i < source_fe.n_dofs_per_face(face_no); ++i)
         {
           const Point<dim> &p = subface_quadrature.point(i);
index ebdec7a3c2b5387e1b62acc82e34a9038368e89c..8667afc104c7252eecf5abd1868b5137a8c05e3d 100644 (file)
@@ -490,7 +490,12 @@ FE_RaviartThomas<dim>::initialize_restriction()
           // evaluated on the subface
           // only.
           const Quadrature<dim> q_sub = QProjector<dim>::project_to_subface(
-            this->reference_cell(), q_base, face, sub);
+            this->reference_cell(),
+            q_base,
+            face,
+            sub,
+            numbers::default_geometric_orientation,
+            RefinementCase<dim - 1>::isotropic_refinement);
           const unsigned int child = GeometryInfo<dim>::child_cell_on_face(
             RefinementCase<dim>::isotropic_refinement, face, sub);
 
index d2abaf304e055b68f49aa37d2f52fe2cf7b13dc2..4e9dbfb28df27546b02ffa14f6c8632c90a4940f 100644 (file)
@@ -623,10 +623,13 @@ FE_RaviartThomasNodal<dim>::get_subface_interpolation_matrix(
   // compute the interpolation matrix by simply taking the value at the
   // support points.
   const Quadrature<dim> subface_projection =
-    QProjector<dim>::project_to_subface(this->reference_cell(),
-                                        quad_face_support,
-                                        0,
-                                        subface);
+    QProjector<dim>::project_to_subface(
+      this->reference_cell(),
+      quad_face_support,
+      0,
+      subface,
+      numbers::default_geometric_orientation,
+      RefinementCase<dim - 1>::isotropic_refinement);
 
   for (unsigned int i = 0; i < source_fe.n_dofs_per_face(face_no); ++i)
     {
index 7035f3a23a8fb371e795358a35b37477e68b5d95..8c5daebd42210b4ca3f7dee4f4ad805fac4b06e0 100644 (file)
@@ -652,20 +652,21 @@ FE_SimplexPoly<dim, spacedim>::get_subface_interpolation_matrix(
 
       const double eps = 2e-13 * this->degree * (dim - 1);
 
-      std::vector<Point<dim>> subface_quadrature_points(
-        quad_face_support.size());
-      QProjector<dim>::project_to_subface(this->reference_cell(),
-                                          quad_face_support,
-                                          face_no,
-                                          subface,
-                                          subface_quadrature_points);
+      const Quadrature<dim> subface_quadrature =
+        QProjector<dim>::project_to_subface(
+          this->reference_cell(),
+          quad_face_support,
+          face_no,
+          subface,
+          numbers::default_geometric_orientation,
+          RefinementCase<dim - 1>::isotropic_refinement);
 
       for (unsigned int i = 0; i < source_fe.n_dofs_per_face(face_no); ++i)
         for (unsigned int j = 0; j < this->n_dofs_per_face(face_no); ++j)
           {
             double matrix_entry =
               this->shape_value(this->face_to_cell_index(j, 0),
-                                subface_quadrature_points[i]);
+                                subface_quadrature.point(i));
 
             // Correct the interpolated value. I.e. if it is close to 1 or
             // 0, make it exactly 1 or 0. Unfortunately, this is required to

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.