From: David Wells <drwells@email.unc.edu>
Date: Tue, 6 Dec 2022 20:55:04 +0000 (-0500)
Subject: MappingRelatedData: move to separate header and source files.
X-Git-Tag: v9.5.0-rc1~725^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a5385921a5dbe712ace3438344b1ffc54e0a055a;p=dealii.git

MappingRelatedData: move to separate header and source files.
---

diff --git a/include/deal.II/fe/fe.h b/include/deal.II/fe/fe.h
index ac2b887bab..dcbc1e6a07 100644
--- a/include/deal.II/fe/fe.h
+++ b/include/deal.II/fe/fe.h
@@ -24,6 +24,7 @@
 #include <deal.II/fe/fe_update_flags.h>
 #include <deal.II/fe/fe_values_extractors.h>
 #include <deal.II/fe/mapping.h>
+#include <deal.II/fe/mapping_related_data.h>
 
 #include <deal.II/lac/full_matrix.h>
 #include <deal.II/lac/vector.h>
diff --git a/include/deal.II/fe/fe_update_flags.h b/include/deal.II/fe/fe_update_flags.h
index b17d55053d..14ce6435f8 100644
--- a/include/deal.II/fe/fe_update_flags.h
+++ b/include/deal.II/fe/fe_update_flags.h
@@ -19,8 +19,6 @@
 
 #include <deal.II/base/config.h>
 
-#include <deal.II/base/derivative_form.h>
-#include <deal.II/base/point.h>
 #include <deal.II/base/table.h>
 #include <deal.II/base/tensor.h>
 
@@ -387,116 +385,6 @@ namespace internal
 {
   namespace FEValuesImplementation
   {
-    /**
-     * A class that stores all of the mapping related data used in
-     * dealii::FEValues, dealii::FEFaceValues, and dealii::FESubfaceValues
-     * objects. Objects of this kind will be given as <i>output</i> argument
-     * when dealii::FEValues::reinit() calls Mapping::fill_fe_values() for a
-     * given cell, face, or subface.
-     *
-     * The data herein will then be provided as <i>input</i> argument in the
-     * following call to FiniteElement::fill_fe_values().
-     *
-     * @ingroup feaccess
-     */
-    template <int dim, int spacedim = dim>
-    class MappingRelatedData
-    {
-    public:
-      /**
-       * Initialize all vectors to correct size.
-       */
-      void
-      initialize(const unsigned int n_quadrature_points,
-                 const UpdateFlags  flags);
-
-      /**
-       * Compute and return an estimate for the memory consumption (in bytes)
-       * of this object.
-       */
-      std::size_t
-      memory_consumption() const;
-
-      /**
-       * Store an array of weights times the Jacobi determinant at the
-       * quadrature points. This function is reset each time reinit() is
-       * called. The Jacobi determinant is actually the reciprocal value of
-       * the Jacobi matrices stored in this class, see the general
-       * documentation of this class for more information.
-       *
-       * However, if this object refers to an FEFaceValues or FESubfaceValues
-       * object, then the JxW_values correspond to the Jacobian of the
-       * transformation of the face, not the cell, i.e. the dimensionality is
-       * that of a surface measure, not of a volume measure. In this case, it
-       * is computed from the boundary forms, rather than the Jacobian matrix.
-       */
-      std::vector<double> JxW_values;
-
-      /**
-       * Array of the Jacobian matrices at the quadrature points.
-       */
-      std::vector<DerivativeForm<1, dim, spacedim>> jacobians;
-
-      /**
-       * Array of the derivatives of the Jacobian matrices at the quadrature
-       * points.
-       */
-      std::vector<DerivativeForm<2, dim, spacedim>> jacobian_grads;
-
-      /**
-       * Array of the inverse Jacobian matrices at the quadrature points.
-       */
-      std::vector<DerivativeForm<1, spacedim, dim>> inverse_jacobians;
-
-      /**
-       * Array of the derivatives of the Jacobian matrices at the quadrature
-       * points, pushed forward to the real cell coordinates.
-       */
-      std::vector<Tensor<3, spacedim>> jacobian_pushed_forward_grads;
-
-      /**
-       * Array of the second derivatives of the Jacobian matrices at the
-       * quadrature points.
-       */
-      std::vector<DerivativeForm<3, dim, spacedim>> jacobian_2nd_derivatives;
-
-      /**
-       * Array of the  second derivatives of the Jacobian matrices at the
-       * quadrature points, pushed forward to the real cell coordinates.
-       */
-      std::vector<Tensor<4, spacedim>> jacobian_pushed_forward_2nd_derivatives;
-
-      /**
-       * Array of the  third derivatives of the Jacobian matrices at the
-       * quadrature points.
-       */
-      std::vector<DerivativeForm<4, dim, spacedim>> jacobian_3rd_derivatives;
-
-      /**
-       * Array of the  third derivatives of the Jacobian matrices at the
-       * quadrature points, pushed forward to the real cell coordinates.
-       */
-      std::vector<Tensor<5, spacedim>> jacobian_pushed_forward_3rd_derivatives;
-
-      /**
-       * Array of quadrature points. This array is set up upon calling
-       * reinit() and contains the quadrature points on the real element,
-       * rather than on the reference element.
-       */
-      std::vector<Point<spacedim>> quadrature_points;
-
-      /**
-       * List of outward normal vectors at the quadrature points.
-       */
-      std::vector<Tensor<1, spacedim>> normal_vectors;
-
-      /**
-       * List of boundary forms at the quadrature points.
-       */
-      std::vector<Tensor<1, spacedim>> boundary_forms;
-    };
-
-
     /**
      * A class that stores all of the shape function related data used in
      * dealii::FEValues, dealii::FEFaceValues, and dealii::FESubfaceValues
diff --git a/include/deal.II/fe/fe_values.h b/include/deal.II/fe/fe_values.h
index 89279498fc..683a802a80 100644
--- a/include/deal.II/fe/fe_values.h
+++ b/include/deal.II/fe/fe_values.h
@@ -34,6 +34,7 @@
 #include <deal.II/fe/fe_update_flags.h>
 #include <deal.II/fe/fe_values_extractors.h>
 #include <deal.II/fe/mapping.h>
+#include <deal.II/fe/mapping_related_data.h>
 
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_iterator.h>
diff --git a/include/deal.II/fe/mapping.h b/include/deal.II/fe/mapping.h
index ef42a7d6d8..9c2568abc6 100644
--- a/include/deal.II/fe/mapping.h
+++ b/include/deal.II/fe/mapping.h
@@ -23,6 +23,7 @@
 #include <deal.II/base/derivative_form.h>
 
 #include <deal.II/fe/fe_update_flags.h>
+#include <deal.II/fe/mapping_related_data.h>
 
 #include <deal.II/grid/tria.h>
 
diff --git a/include/deal.II/fe/mapping_related_data.h b/include/deal.II/fe/mapping_related_data.h
new file mode 100644
index 0000000000..ab525b9070
--- /dev/null
+++ b/include/deal.II/fe/mapping_related_data.h
@@ -0,0 +1,157 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 - 2022 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_fe_mapping_related_data_h
+#define dealii_fe_mapping_related_data_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/derivative_form.h>
+#include <deal.II/base/point.h>
+#include <deal.II/base/tensor.h>
+
+#include <deal.II/fe/fe_update_flags.h>
+
+#include <vector>
+
+
+DEAL_II_NAMESPACE_OPEN
+
+/** @addtogroup feaccess */
+/** @{ */
+
+namespace internal
+{
+  namespace FEValuesImplementation
+  {
+    /**
+     * A class that stores all of the mapping related data used in
+     * dealii::FEValues, dealii::FEFaceValues, and dealii::FESubfaceValues
+     * objects. Objects of this kind will be given as <i>output</i> argument
+     * when dealii::FEValues::reinit() calls Mapping::fill_fe_values() for a
+     * given cell, face, or subface.
+     *
+     * The data herein will then be provided as <i>input</i> argument in the
+     * following call to FiniteElement::fill_fe_values().
+     *
+     * @ingroup feaccess
+     */
+    template <int dim, int spacedim = dim>
+    class MappingRelatedData
+    {
+    public:
+      /**
+       * Initialize all vectors to correct size.
+       */
+      void
+      initialize(const unsigned int n_quadrature_points,
+                 const UpdateFlags  flags);
+
+      /**
+       * Compute and return an estimate for the memory consumption (in bytes)
+       * of this object.
+       */
+      std::size_t
+      memory_consumption() const;
+
+      /**
+       * Store an array of weights times the Jacobi determinant at the
+       * quadrature points. This function is reset each time reinit() is
+       * called. The Jacobi determinant is actually the reciprocal value of
+       * the Jacobi matrices stored in this class, see the general
+       * documentation of this class for more information.
+       *
+       * However, if this object refers to an FEFaceValues or FESubfaceValues
+       * object, then the JxW_values correspond to the Jacobian of the
+       * transformation of the face, not the cell, i.e. the dimensionality is
+       * that of a surface measure, not of a volume measure. In this case, it
+       * is computed from the boundary forms, rather than the Jacobian matrix.
+       */
+      std::vector<double> JxW_values;
+
+      /**
+       * Array of the Jacobian matrices at the quadrature points.
+       */
+      std::vector<DerivativeForm<1, dim, spacedim>> jacobians;
+
+      /**
+       * Array of the derivatives of the Jacobian matrices at the quadrature
+       * points.
+       */
+      std::vector<DerivativeForm<2, dim, spacedim>> jacobian_grads;
+
+      /**
+       * Array of the inverse Jacobian matrices at the quadrature points.
+       */
+      std::vector<DerivativeForm<1, spacedim, dim>> inverse_jacobians;
+
+      /**
+       * Array of the derivatives of the Jacobian matrices at the quadrature
+       * points, pushed forward to the real cell coordinates.
+       */
+      std::vector<Tensor<3, spacedim>> jacobian_pushed_forward_grads;
+
+      /**
+       * Array of the second derivatives of the Jacobian matrices at the
+       * quadrature points.
+       */
+      std::vector<DerivativeForm<3, dim, spacedim>> jacobian_2nd_derivatives;
+
+      /**
+       * Array of the  second derivatives of the Jacobian matrices at the
+       * quadrature points, pushed forward to the real cell coordinates.
+       */
+      std::vector<Tensor<4, spacedim>> jacobian_pushed_forward_2nd_derivatives;
+
+      /**
+       * Array of the  third derivatives of the Jacobian matrices at the
+       * quadrature points.
+       */
+      std::vector<DerivativeForm<4, dim, spacedim>> jacobian_3rd_derivatives;
+
+      /**
+       * Array of the  third derivatives of the Jacobian matrices at the
+       * quadrature points, pushed forward to the real cell coordinates.
+       */
+      std::vector<Tensor<5, spacedim>> jacobian_pushed_forward_3rd_derivatives;
+
+      /**
+       * Array of quadrature points. This array is set up upon calling
+       * reinit() and contains the quadrature points on the real element,
+       * rather than on the reference element.
+       */
+      std::vector<Point<spacedim>> quadrature_points;
+
+      /**
+       * List of outward normal vectors at the quadrature points.
+       */
+      std::vector<Tensor<1, spacedim>> normal_vectors;
+
+      /**
+       * List of boundary forms at the quadrature points.
+       */
+      std::vector<Tensor<1, spacedim>> boundary_forms;
+    };
+  } // namespace FEValuesImplementation
+} // namespace internal
+
+
+/** @} */
+
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
diff --git a/include/deal.II/non_matching/mapping_info.h b/include/deal.II/non_matching/mapping_info.h
index 1a70680da9..241622d326 100644
--- a/include/deal.II/non_matching/mapping_info.h
+++ b/include/deal.II/non_matching/mapping_info.h
@@ -28,6 +28,7 @@
 #include <deal.II/fe/fe_update_flags.h>
 #include <deal.II/fe/mapping_cartesian.h>
 #include <deal.II/fe/mapping_q.h>
+#include <deal.II/fe/mapping_related_data.h>
 
 #include <memory>
 
diff --git a/source/fe/CMakeLists.txt b/source/fe/CMakeLists.txt
index 177a8c2df8..776c82fe45 100644
--- a/source/fe/CMakeLists.txt
+++ b/source/fe/CMakeLists.txt
@@ -65,6 +65,7 @@ set(_unity_include_src
   mapping_q.cc
   mapping_q_cache.cc
   mapping_manifold.cc
+  mapping_related_data.cc
   )
 
 set(_separate_src
@@ -142,6 +143,7 @@ set(_inst
   mapping_q_eulerian.inst.in
   mapping_q.inst.in
   mapping_manifold.inst.in
+  mapping_related_data.inst.in
   )
 
 file(GLOB _header
diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc
index d742414a27..cff0ae89ef 100644
--- a/source/fe/fe_values.cc
+++ b/source/fe/fe_values.cc
@@ -2699,90 +2699,6 @@ namespace internal
 {
   namespace FEValuesImplementation
   {
-    template <int dim, int spacedim>
-    void
-    MappingRelatedData<dim, spacedim>::initialize(
-      const unsigned int n_quadrature_points,
-      const UpdateFlags  flags)
-    {
-      if (flags & update_quadrature_points)
-        this->quadrature_points.resize(
-          n_quadrature_points,
-          Point<spacedim>(numbers::signaling_nan<Tensor<1, spacedim>>()));
-
-      if (flags & update_JxW_values)
-        this->JxW_values.resize(n_quadrature_points,
-                                numbers::signaling_nan<double>());
-
-      if (flags & update_jacobians)
-        this->jacobians.resize(
-          n_quadrature_points,
-          numbers::signaling_nan<DerivativeForm<1, dim, spacedim>>());
-
-      if (flags & update_jacobian_grads)
-        this->jacobian_grads.resize(
-          n_quadrature_points,
-          numbers::signaling_nan<DerivativeForm<2, dim, spacedim>>());
-
-      if (flags & update_jacobian_pushed_forward_grads)
-        this->jacobian_pushed_forward_grads.resize(
-          n_quadrature_points, numbers::signaling_nan<Tensor<3, spacedim>>());
-
-      if (flags & update_jacobian_2nd_derivatives)
-        this->jacobian_2nd_derivatives.resize(
-          n_quadrature_points,
-          numbers::signaling_nan<DerivativeForm<3, dim, spacedim>>());
-
-      if (flags & update_jacobian_pushed_forward_2nd_derivatives)
-        this->jacobian_pushed_forward_2nd_derivatives.resize(
-          n_quadrature_points, numbers::signaling_nan<Tensor<4, spacedim>>());
-
-      if (flags & update_jacobian_3rd_derivatives)
-        this->jacobian_3rd_derivatives.resize(n_quadrature_points);
-
-      if (flags & update_jacobian_pushed_forward_3rd_derivatives)
-        this->jacobian_pushed_forward_3rd_derivatives.resize(
-          n_quadrature_points, numbers::signaling_nan<Tensor<5, spacedim>>());
-
-      if (flags & update_inverse_jacobians)
-        this->inverse_jacobians.resize(
-          n_quadrature_points,
-          numbers::signaling_nan<DerivativeForm<1, spacedim, dim>>());
-
-      if (flags & update_boundary_forms)
-        this->boundary_forms.resize(
-          n_quadrature_points, numbers::signaling_nan<Tensor<1, spacedim>>());
-
-      if (flags & update_normal_vectors)
-        this->normal_vectors.resize(
-          n_quadrature_points, numbers::signaling_nan<Tensor<1, spacedim>>());
-    }
-
-
-
-    template <int dim, int spacedim>
-    std::size_t
-    MappingRelatedData<dim, spacedim>::memory_consumption() const
-    {
-      return (
-        MemoryConsumption::memory_consumption(JxW_values) +
-        MemoryConsumption::memory_consumption(jacobians) +
-        MemoryConsumption::memory_consumption(jacobian_grads) +
-        MemoryConsumption::memory_consumption(jacobian_pushed_forward_grads) +
-        MemoryConsumption::memory_consumption(jacobian_2nd_derivatives) +
-        MemoryConsumption::memory_consumption(
-          jacobian_pushed_forward_2nd_derivatives) +
-        MemoryConsumption::memory_consumption(jacobian_3rd_derivatives) +
-        MemoryConsumption::memory_consumption(
-          jacobian_pushed_forward_3rd_derivatives) +
-        MemoryConsumption::memory_consumption(inverse_jacobians) +
-        MemoryConsumption::memory_consumption(quadrature_points) +
-        MemoryConsumption::memory_consumption(normal_vectors) +
-        MemoryConsumption::memory_consumption(boundary_forms));
-    }
-
-
-
     template <int dim, int spacedim>
     void
     FiniteElementRelatedData<dim, spacedim>::initialize(
diff --git a/source/fe/fe_values.inst.in b/source/fe/fe_values.inst.in
index 8d5cd661f1..63d6b4362b 100644
--- a/source/fe/fe_values.inst.in
+++ b/source/fe/fe_values.inst.in
@@ -40,8 +40,6 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
     \{
       namespace FEValuesImplementation
       \{
-        template class MappingRelatedData<deal_II_dimension,
-                                          deal_II_space_dimension>;
         template class FiniteElementRelatedData<deal_II_dimension,
                                                 deal_II_space_dimension>;
       \}
diff --git a/source/fe/mapping_related_data.cc b/source/fe/mapping_related_data.cc
new file mode 100644
index 0000000000..9c9dfd1875
--- /dev/null
+++ b/source/fe/mapping_related_data.cc
@@ -0,0 +1,114 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 1998 - 2022 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/memory_consumption.h>
+#include <deal.II/base/numbers.h>
+#include <deal.II/base/signaling_nan.h>
+
+#include <deal.II/fe/mapping_related_data.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace internal
+{
+  namespace FEValuesImplementation
+  {
+    template <int dim, int spacedim>
+    void
+    MappingRelatedData<dim, spacedim>::initialize(
+      const unsigned int n_quadrature_points,
+      const UpdateFlags  flags)
+    {
+      if (flags & update_quadrature_points)
+        this->quadrature_points.resize(
+          n_quadrature_points,
+          Point<spacedim>(numbers::signaling_nan<Tensor<1, spacedim>>()));
+
+      if (flags & update_JxW_values)
+        this->JxW_values.resize(n_quadrature_points,
+                                numbers::signaling_nan<double>());
+
+      if (flags & update_jacobians)
+        this->jacobians.resize(
+          n_quadrature_points,
+          numbers::signaling_nan<DerivativeForm<1, dim, spacedim>>());
+
+      if (flags & update_jacobian_grads)
+        this->jacobian_grads.resize(
+          n_quadrature_points,
+          numbers::signaling_nan<DerivativeForm<2, dim, spacedim>>());
+
+      if (flags & update_jacobian_pushed_forward_grads)
+        this->jacobian_pushed_forward_grads.resize(
+          n_quadrature_points, numbers::signaling_nan<Tensor<3, spacedim>>());
+
+      if (flags & update_jacobian_2nd_derivatives)
+        this->jacobian_2nd_derivatives.resize(
+          n_quadrature_points,
+          numbers::signaling_nan<DerivativeForm<3, dim, spacedim>>());
+
+      if (flags & update_jacobian_pushed_forward_2nd_derivatives)
+        this->jacobian_pushed_forward_2nd_derivatives.resize(
+          n_quadrature_points, numbers::signaling_nan<Tensor<4, spacedim>>());
+
+      if (flags & update_jacobian_3rd_derivatives)
+        this->jacobian_3rd_derivatives.resize(n_quadrature_points);
+
+      if (flags & update_jacobian_pushed_forward_3rd_derivatives)
+        this->jacobian_pushed_forward_3rd_derivatives.resize(
+          n_quadrature_points, numbers::signaling_nan<Tensor<5, spacedim>>());
+
+      if (flags & update_inverse_jacobians)
+        this->inverse_jacobians.resize(
+          n_quadrature_points,
+          numbers::signaling_nan<DerivativeForm<1, spacedim, dim>>());
+
+      if (flags & update_boundary_forms)
+        this->boundary_forms.resize(
+          n_quadrature_points, numbers::signaling_nan<Tensor<1, spacedim>>());
+
+      if (flags & update_normal_vectors)
+        this->normal_vectors.resize(
+          n_quadrature_points, numbers::signaling_nan<Tensor<1, spacedim>>());
+    }
+
+
+
+    template <int dim, int spacedim>
+    std::size_t
+    MappingRelatedData<dim, spacedim>::memory_consumption() const
+    {
+      return (
+        MemoryConsumption::memory_consumption(JxW_values) +
+        MemoryConsumption::memory_consumption(jacobians) +
+        MemoryConsumption::memory_consumption(jacobian_grads) +
+        MemoryConsumption::memory_consumption(jacobian_pushed_forward_grads) +
+        MemoryConsumption::memory_consumption(jacobian_2nd_derivatives) +
+        MemoryConsumption::memory_consumption(
+          jacobian_pushed_forward_2nd_derivatives) +
+        MemoryConsumption::memory_consumption(jacobian_3rd_derivatives) +
+        MemoryConsumption::memory_consumption(
+          jacobian_pushed_forward_3rd_derivatives) +
+        MemoryConsumption::memory_consumption(inverse_jacobians) +
+        MemoryConsumption::memory_consumption(quadrature_points) +
+        MemoryConsumption::memory_consumption(normal_vectors) +
+        MemoryConsumption::memory_consumption(boundary_forms));
+    }
+  } // namespace FEValuesImplementation
+} // namespace internal
+
+#include "mapping_related_data.inst"
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/fe/mapping_related_data.inst.in b/source/fe/mapping_related_data.inst.in
new file mode 100644
index 0000000000..6d586bb650
--- /dev/null
+++ b/source/fe/mapping_related_data.inst.in
@@ -0,0 +1,31 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 - 2022 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef DOXYGEN
+
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
+  {
+#  if deal_II_dimension <= deal_II_space_dimension
+    namespace internal
+    \{
+      namespace FEValuesImplementation
+      \{
+        template class MappingRelatedData<deal_II_dimension,
+                                          deal_II_space_dimension>;
+      \}
+    \}
+#  endif
+  }
+#endif