]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Eliminate FEEvaluationData::get_all_face_numbers () and get_all_face_orientations () 13078/head
authorPeter Munch <peterrmuench@gmail.com>
Wed, 15 Dec 2021 08:21:15 +0000 (09:21 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 15 Dec 2021 12:09:15 +0000 (13:09 +0100)
include/deal.II/matrix_free/evaluation_kernels.h
include/deal.II/matrix_free/fe_evaluation.h
include/deal.II/matrix_free/fe_evaluation_data.h

index c85113a8678aee28b57ad823fb52dc2df7f54ea8..ea64b69fc95677114bce3b6f4fc7e63cef460aff 100644 (file)
@@ -3192,10 +3192,8 @@ namespace internal
     bool all_faces_are_same = n_filled_lanes == n_lanes;
     if (n_face_orientations == n_lanes)
       for (unsigned int v = 1; v < n_lanes; ++v)
-        if (fe_eval.get_all_face_numbers()[v] !=
-              fe_eval.get_all_face_numbers()[0] ||
-            fe_eval.get_all_face_orientations()[v] !=
-              fe_eval.get_all_face_orientations()[0])
+        if (fe_eval.get_face_no(v) != fe_eval.get_face_no(0) ||
+            fe_eval.get_face_orientation(v) != fe_eval.get_face_orientation(0))
           {
             all_faces_are_same = false;
             break;
@@ -3224,9 +3222,9 @@ namespace internal
             continue;
 
           if (shape_data.nodal_at_cell_boundaries &&
-              fe_eval.get_all_face_orientations()[v] != 0)
+              fe_eval.get_face_orientation(v) != 0)
             orientation[v] = &fe_eval.get_shape_info().face_orientations_dofs(
-              fe_eval.get_all_face_orientations()[v], 0);
+              fe_eval.get_face_orientation(v), 0);
         }
     else if (fe_eval.get_face_orientation() != 0)
       orientation[0] = &fe_eval.get_shape_info().face_orientations_dofs(
@@ -3242,20 +3240,22 @@ namespace internal
             &fe_eval.get_shape_info().face_to_cell_index_hermite(face_no, 0);
         else
           {
-            const auto &face_nos = fe_eval.get_all_face_numbers();
             for (unsigned int v = 0; v < n_lanes; ++v)
               {
                 if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int)
                   continue;
 
+                const auto face_no = fe_eval.get_face_no(v);
+
                 grad_weight[v] =
-                  shape_data.shape_data_on_face
-                    [0][fe_degree + (integrate ? (2 - (face_nos[v] % 2)) :
-                                                 (1 + (face_nos[v] % 2)))][0];
+                  shape_data.shape_data_on_face[0][fe_degree +
+                                                   (integrate ?
+                                                      (2 - (face_no % 2)) :
+                                                      (1 + (face_no % 2)))][0];
 
                 index_array_hermite[v] =
-                  &fe_eval.get_shape_info().face_to_cell_index_hermite(
-                    face_nos[v], 0);
+                  &fe_eval.get_shape_info().face_to_cell_index_hermite(face_no,
+                                                                       0);
               }
           }
       }
@@ -3270,15 +3270,16 @@ namespace internal
             &fe_eval.get_shape_info().face_to_cell_index_nodal(face_no, 0);
         else
           {
-            const auto &face_nos = fe_eval.get_all_face_numbers();
             for (unsigned int v = 0; v < n_lanes; ++v)
               {
                 if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int)
                   continue;
 
+                const auto face_no = fe_eval.get_face_no(v);
+
                 index_array_nodal[v] =
-                  &fe_eval.get_shape_info().face_to_cell_index_nodal(
-                    face_nos[v], 0);
+                  &fe_eval.get_shape_info().face_to_cell_index_nodal(face_no,
+                                                                     0);
               }
           }
       }
@@ -3766,14 +3767,14 @@ namespace internal
                       numbers::invalid_unsigned_int)
                     continue;
 
-                  if (fe_eval.get_all_face_orientations()[v] != 0)
+                  if (fe_eval.get_face_orientation(v) != 0)
                     adjust_for_face_orientation_per_lane(
                       dim,
                       n_components,
                       v,
                       evaluation_flag,
                       &fe_eval.get_shape_info().face_orientations_quad(
-                        fe_eval.get_all_face_orientations()[v], 0),
+                        fe_eval.get_face_orientation(v), 0),
                       false,
                       Utilities::pow(n_q_points_1d, dim - 1),
                       &temp[0][0],
@@ -3944,14 +3945,14 @@ namespace internal
                 if (fe_eval.get_cell_ids()[v] == numbers::invalid_unsigned_int)
                   continue;
 
-                if (fe_eval.get_all_face_orientations()[v] != 0)
+                if (fe_eval.get_face_orientation(v) != 0)
                   adjust_for_face_orientation_per_lane(
                     dim,
                     n_components,
                     v,
                     integration_flag,
                     &fe_eval.get_shape_info().face_orientations_quad(
-                      fe_eval.get_all_face_orientations()[v], 0),
+                      fe_eval.get_face_orientation(v), 0),
                     true,
                     Utilities::pow(n_q_points_1d, dim - 1),
                     &temp[0][0],
index 484049329233ffe97ec681c31034aca94deb8255..dfeb9506a9757f78f87f7b5ebe4a1d8faa411406 100644 (file)
@@ -7658,9 +7658,9 @@ FEFaceEvaluation<dim,
 
           if (face_index == numbers::invalid_unsigned_int)
             {
-              this->cell_ids[i]              = numbers::invalid_unsigned_int;
-              this->all_face_numbers[i]      = static_cast<std::uint8_t>(-1);
-              this->all_face_orientations[i] = static_cast<std::uint8_t>(-1);
+              this->cell_ids[i]          = numbers::invalid_unsigned_int;
+              this->face_numbers[i]      = static_cast<std::uint8_t>(-1);
+              this->face_orientations[i] = static_cast<std::uint8_t>(-1);
               continue; // invalid face ID: no neighbor on boundary
             }
 
@@ -7678,13 +7678,13 @@ FEFaceEvaluation<dim,
           // compare the IDs with the given cell ID
           if (face_identifies_as_interior)
             {
-              this->cell_ids[i]         = cell_m; // neighbor has the other ID
-              this->all_face_numbers[i] = faces.interior_face_no;
+              this->cell_ids[i]     = cell_m; // neighbor has the other ID
+              this->face_numbers[i] = faces.interior_face_no;
             }
           else
             {
-              this->cell_ids[i]         = cell_p;
-              this->all_face_numbers[i] = faces.exterior_face_no;
+              this->cell_ids[i]     = cell_p;
+              this->face_numbers[i] = faces.exterior_face_no;
             }
 
           const bool   orientation_interior_face = faces.face_orientation >= 8;
@@ -7695,16 +7695,13 @@ FEFaceEvaluation<dim,
                 {0, 1, 2, 3, 6, 5, 4, 7}};
               face_orientation = table[face_orientation];
             }
-          this->all_face_orientations[i] = face_orientation;
+          this->face_orientations[i] = face_orientation;
         }
-
-      this->face_no          = this->all_face_numbers[0];
-      this->face_orientation = this->all_face_orientations[0];
     }
   else
     {
-      this->face_orientation = 0;
-      this->face_no          = face_number;
+      this->face_orientations[0] = 0;
+      this->face_numbers[0]      = face_number;
       for (unsigned int i = 0; i < n_lanes; ++i)
         this->cell_ids[i] = cell_index * n_lanes + i;
     }
@@ -8250,7 +8247,7 @@ FEFaceEvaluation<dim,
              internal::ExcMatrixFreeAccessToUninitializedMappingField(
                "update_quadrature_points"));
       const unsigned int index =
-        this->cell * GeometryInfo<dim>::faces_per_cell + this->face_no;
+        this->cell * GeometryInfo<dim>::faces_per_cell + this->face_numbers[0];
       AssertIndexRange(index,
                        this->matrix_free->get_mapping_info()
                          .face_data_by_cells[this->quad_no]
index cd0b0c6226188a50779bb7083c96ff213f4b2e57..de3cb3753535a152fa04ddc3e7f504d6bd589c5e 100644 (file)
@@ -406,16 +406,19 @@ public:
   get_first_selected_component() const;
 
   /**
-   * If is_face is true, this function returns the face number within a cell for
-   * the face this object was initialized to. On cells where the face makes no
-   * sense, an exception is thrown.
+   * If is_face is true, this function returns the face number of the given lane
+   * within a cell for the face this object was initialized to. On cells where
+   * the face makes no sense, an exception is thrown.
+   *
+   * @note Only available for `dof_access_index == dof_access_cell` and
+   * `is_interior_face == false`.
    *
    * @note This function depends on the internal representation of data, which
    * is not stable between releases of deal.II, and is hence mostly for
    * internal use.
    */
-  unsigned int
-  get_face_no() const;
+  std::uint8_t
+  get_face_no(const unsigned int v = 0) const;
 
   /**
    * If is_face is true, this function returns the index of a subface along a
@@ -430,16 +433,15 @@ public:
   get_subface_index() const;
 
   /**
-   * If is_face is true, this function returns the orientation index within an
-   * array of orientations as stored in ShapeInfo for unknowns and quadrature
-   * points.
+   * If is_face is true, this function returns for a given lane the
+   * orientation index within an array of orientations as stored in
+   * ShapeInfo for unknowns and quadrature points.
    *
-   * @note This function depends on the internal representation of data, which
-   * is not stable between releases of deal.II, and is hence mostly for
-   * internal use.
+   * @note Only available for `dof_access_index == dof_access_cell` and
+   * `is_interior_face == false`.
    */
-  unsigned int
-  get_face_orientation() const;
+  std::uint8_t
+  get_face_orientation(const unsigned int v = 0) const;
 
   /**
    * Return the current index in the access to compressed indices.
@@ -502,58 +504,6 @@ public:
       return cell_or_face_ids;
   }
 
-  /**
-   * Return the (non-vectorized) number of faces within cells in case of ECL
-   * for the exterior cells where the single number accessible via
-   * get_face_no() is not enough to cover all cases.
-   *
-   * @note Only available for `dof_access_index == dof_access_cell` and
-   * `is_interior_face == false`.
-   *
-   * @note This function depends on the internal representation of data, which
-   * is not stable between releases of deal.II, and is hence mostly for
-   * internal use.
-   */
-  const std::array<std::uint8_t, n_lanes> &
-  get_all_face_numbers() const
-  {
-    // implemented inline to avoid compilation problems on Windows
-    Assert(cell != numbers::invalid_unsigned_int, ExcNotInitialized());
-    Assert(is_face &&
-             dof_access_index ==
-               internal::MatrixFreeFunctions::DoFInfo::dof_access_cell &&
-             is_interior_face == false,
-           ExcMessage(
-             "All face numbers can only be queried for ECL at exterior "
-             "faces. Use get_face_no() in other cases."));
-
-    return all_face_numbers;
-  }
-
-  /**
-   * Store the orientation of the neighbor's faces with respect to the current
-   * cell for the case of exterior faces on ECL with possibly different
-   * orientations behind different cells.
-   *
-   * @note Only available for `dof_access_index == dof_access_cell` and
-   * `is_interior_face == false`.
-   */
-  const std::array<std::uint8_t, n_lanes> &
-  get_all_face_orientations() const
-  {
-    // implemented inline to avoid compilation problems on Windows
-    Assert(cell != numbers::invalid_unsigned_int, ExcNotInitialized());
-    Assert(is_face &&
-             dof_access_index ==
-               internal::MatrixFreeFunctions::DoFInfo::dof_access_cell &&
-             is_interior_face == false,
-           ExcMessage(
-             "All face numbers can only be queried for ECL at exterior "
-             "faces. Use get_face_no() in other cases."));
-
-    return all_face_orientations;
-  }
-
   //@}
 
   /**
@@ -825,16 +775,24 @@ protected:
   internal::MatrixFreeFunctions::DoFInfo::DoFAccessIndex dof_access_index;
 
   /**
-   * Stores the current number of a face within the given cell in case
-   * `is_face==true`, using values between `0` and `2*dim`.
+   * Stores the (non-vectorized) number of faces within cells in case of ECL
+   * for the exterior cells where the single number `face_no` is not enough to
+   * cover all cases.
+   *
+   * @note Only available for `dof_access_index == dof_access_cell` and
+   * `is_interior_face == false`.
    */
-  unsigned int face_no;
+  std::array<std::uint8_t, n_lanes> face_numbers;
 
   /**
-   * Stores the orientation of the given face with respect to the standard
-   * orientation, 0 if in standard orientation.
+   * Store the orientation of the neighbor's faces with respect to the current
+   * cell for the case of exterior faces on ECL with possibly different
+   * orientations behind different cells.
+   *
+   * @note Only available for `dof_access_index == dof_access_cell` and
+   * `is_interior_face == false`.
    */
-  unsigned int face_orientation;
+  std::array<std::uint8_t, n_lanes> face_orientations;
 
   /**
    * Stores the subface index of the given face. Usually, this variable takes
@@ -865,26 +823,6 @@ protected:
    */
   std::array<unsigned int, n_lanes> cell_or_face_ids;
 
-  /**
-   * Stores the (non-vectorized) number of faces within cells in case of ECL
-   * for the exterior cells where the single number `face_no` is not enough to
-   * cover all cases.
-   *
-   * @note Only available for `dof_access_index == dof_access_cell` and
-   * `is_interior_face == false`.
-   */
-  std::array<std::uint8_t, n_lanes> all_face_numbers;
-
-  /**
-   * Store the orientation of the neighbor's faces with respect to the current
-   * cell for the case of exterior faces on ECL with possibly different
-   * orientations behind different cells.
-   *
-   * @note Only available for `dof_access_index == dof_access_cell` and
-   * `is_interior_face == false`.
-   */
-  std::array<std::uint8_t, n_lanes> all_face_orientations;
-
   /**
    * Geometry data that can be generated FEValues on the fly with the
    * respective constructor, as an alternative to the entry point with
@@ -972,8 +910,6 @@ inline FEEvaluationData<dim, Number, is_face>::FEEvaluationData(
            internal::MatrixFreeFunctions::DoFInfo::dof_access_face_interior :
            internal::MatrixFreeFunctions::DoFInfo::dof_access_face_exterior) :
         internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
-  , face_no(0)
-  , face_orientation(0)
   , subface_index(0)
   , cell_type(internal::MatrixFreeFunctions::general)
 {}
@@ -1054,10 +990,10 @@ FEEvaluationData<dim, Number, is_face>::operator=(const FEEvaluationData &other)
          internal::MatrixFreeFunctions::DoFInfo::dof_access_face_interior :
          internal::MatrixFreeFunctions::DoFInfo::dof_access_face_exterior) :
       internal::MatrixFreeFunctions::DoFInfo::dof_access_cell;
-  face_no          = 0;
-  face_orientation = 0;
-  subface_index    = 0;
-  cell_type        = internal::MatrixFreeFunctions::general;
+  face_numbers[0]      = 0;
+  face_orientations[0] = 0;
+  subface_index        = 0;
+  cell_type            = internal::MatrixFreeFunctions::general;
 
   return *this;
 }
@@ -1113,7 +1049,8 @@ FEEvaluationData<dim, Number, is_face>::reinit_face(
   Assert(is_face == true,
          ExcMessage("Faces can only be set if the is_face template parameter "
                     "is true"));
-  face_no = (is_interior_face ? face.interior_face_no : face.exterior_face_no);
+  face_numbers[0] =
+    (is_interior_face ? face.interior_face_no : face.exterior_face_no);
   subface_index = is_interior_face == true ?
                     GeometryInfo<dim>::max_children_per_cell :
                     face.subface_index;
@@ -1123,9 +1060,9 @@ FEEvaluationData<dim, Number, is_face>::reinit_face(
   // standard-orientation else copy the first three bits
   // (which is equivalent to modulo 8). See also the documentation of
   // internal::MatrixFreeFunctions::FaceToCellTopology::face_orientation.
-  face_orientation = (is_interior_face == (face.face_orientation >= 8)) ?
-                       (face.face_orientation % 8) :
-                       0;
+  face_orientations[0] = (is_interior_face == (face.face_orientation >= 8)) ?
+                           (face.face_orientation % 8) :
+                           0;
 
   if (is_interior_face)
     cell_ids = face.cells_interior;
@@ -1357,7 +1294,7 @@ FEEvaluationData<dim, Number, is_face>::get_current_cell_index() const
 {
   if (is_face && dof_access_index ==
                    internal::MatrixFreeFunctions::DoFInfo::dof_access_cell)
-    return cell * GeometryInfo<dim>::faces_per_cell + face_no;
+    return cell * GeometryInfo<dim>::faces_per_cell + face_numbers[0];
   else
     return cell;
 }
@@ -1401,10 +1338,18 @@ FEEvaluationData<dim, Number, is_face>::get_scratch_data() const
 
 
 template <int dim, typename Number, bool is_face>
-inline unsigned int
-FEEvaluationData<dim, Number, is_face>::get_face_no() const
+std::uint8_t
+FEEvaluationData<dim, Number, is_face>::get_face_no(const unsigned int v) const
 {
-  return face_no;
+  Assert(is_face, ExcNotInitialized());
+  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
+                    dof_access_index ==
+                      internal::MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+                    is_interior_face == false),
+         ExcMessage("All face numbers can only be queried for ECL at exterior "
+                    "faces. Use get_face_no() in other cases."));
+
+  return face_numbers[v];
 }
 
 
@@ -1419,10 +1364,19 @@ FEEvaluationData<dim, Number, is_face>::get_subface_index() const
 
 
 template <int dim, typename Number, bool is_face>
-inline unsigned int
-FEEvaluationData<dim, Number, is_face>::get_face_orientation() const
+std::uint8_t
+FEEvaluationData<dim, Number, is_face>::get_face_orientation(
+  const unsigned int v) const
 {
-  return face_orientation;
+  Assert(is_face, ExcNotInitialized());
+  Assert(v == 0 || (cell != numbers::invalid_unsigned_int &&
+                    dof_access_index ==
+                      internal::MatrixFreeFunctions::DoFInfo::dof_access_cell &&
+                    is_interior_face == false),
+         ExcMessage("All face numbers can only be queried for ECL at exterior "
+                    "faces. Use get_face_no() in other cases."));
+
+  return face_orientations[v];
 }
 
 

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.