]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Introduce MappingCollection in MatrixFree 11209/head
authorPeter Munch <peterrmuench@gmail.com>
Thu, 19 Nov 2020 20:39:36 +0000 (21:39 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Fri, 20 Nov 2020 14:17:20 +0000 (15:17 +0100)
include/deal.II/matrix_free/mapping_info.h
include/deal.II/matrix_free/mapping_info.templates.h
include/deal.II/matrix_free/matrix_free.h
include/deal.II/matrix_free/matrix_free.templates.h
source/matrix_free/matrix_free.inst.in

index d717e3fc1790dfbdf3f2582d92e5c892c39f913e..941d530556351b7809cfcf01eba7f483081223d0 100644 (file)
@@ -27,6 +27,7 @@
 #include <deal.II/fe/fe.h>
 #include <deal.II/fe/mapping.h>
 
+#include <deal.II/hp/mapping_collection.h>
 #include <deal.II/hp/q_collection.h>
 
 #include <deal.II/matrix_free/face_info.h>
@@ -335,10 +336,10 @@ namespace internal
         const dealii::Triangulation<dim> &                        tria,
         const std::vector<std::pair<unsigned int, unsigned int>> &cells,
         const FaceInfo<VectorizedArrayType::size()> &             faces,
-        const std::vector<unsigned int> &                active_fe_index,
-        const Mapping<dim> &                             mapping,
-        const std::vector<dealii::hp::QCollection<dim>> &quad,
-        const UpdateFlags                                update_flags_cells,
+        const std::vector<unsigned int> &active_fe_index,
+        const std::shared_ptr<dealii::hp::MappingCollection<dim>> &mapping,
+        const std::vector<dealii::hp::QCollection<dim>> &          quad,
+        const UpdateFlags update_flags_cells,
         const UpdateFlags update_flags_boundary_faces,
         const UpdateFlags update_flags_inner_faces,
         const UpdateFlags update_flags_faces_by_cells);
@@ -355,7 +356,7 @@ namespace internal
         const std::vector<std::pair<unsigned int, unsigned int>> &cells,
         const FaceInfo<VectorizedArrayType::size()> &             faces,
         const std::vector<unsigned int> &active_fe_index,
-        const Mapping<dim> &             mapping);
+        const std::shared_ptr<dealii::hp::MappingCollection<dim>> &mapping);
 
       /**
        * Return the type of a given cell as detected during initialization.
@@ -444,7 +445,12 @@ namespace internal
         face_data_by_cells;
 
       /**
-       * The pointer to the underlying Mapping object.
+       * The pointer to the underlying hp::MappingCollection object.
+       */
+      std::shared_ptr<dealii::hp::MappingCollection<dim>> mapping_collection;
+
+      /**
+       * The pointer to the first entry of mapping_collection.
        */
       SmartPointer<const Mapping<dim>> mapping;
 
@@ -484,8 +490,8 @@ namespace internal
       initialize_cells(
         const dealii::Triangulation<dim> &                        tria,
         const std::vector<std::pair<unsigned int, unsigned int>> &cells,
-        const std::vector<unsigned int> &active_fe_index,
-        const Mapping<dim> &             mapping);
+        const std::vector<unsigned int> &         active_fe_index,
+        const dealii::hp::MappingCollection<dim> &mapping);
 
       /**
        * Computes the information in the given faces, called within
@@ -496,8 +502,8 @@ namespace internal
         const dealii::Triangulation<dim> &                        tria,
         const std::vector<std::pair<unsigned int, unsigned int>> &cells,
         const std::vector<FaceToCellTopology<VectorizedArrayType::size()>>
-          &                 faces,
-        const Mapping<dim> &mapping);
+          &                                       faces,
+        const dealii::hp::MappingCollection<dim> &mapping);
 
       /**
        * Computes the information in the given faces, called within
@@ -507,7 +513,7 @@ namespace internal
       initialize_faces_by_cells(
         const dealii::Triangulation<dim> &                        tria,
         const std::vector<std::pair<unsigned int, unsigned int>> &cells,
-        const Mapping<dim> &                                      mapping);
+        const dealii::hp::MappingCollection<dim> &                mapping);
 
       /**
        * Helper function to determine which update flags must be set in the
index e81fcb0fe8595bcc10449747f3ccecaccba49106..4c127012408fe74139395ec59c7ac4fa1efa5d8c 100644 (file)
@@ -272,7 +272,8 @@ namespace internal
       face_data_by_cells.clear();
       cell_type.clear();
       face_type.clear();
-      mapping = nullptr;
+      mapping_collection = nullptr;
+      mapping            = nullptr;
     }
 
 
@@ -331,15 +332,16 @@ namespace internal
       const std::vector<std::pair<unsigned int, unsigned int>> &cells,
       const FaceInfo<VectorizedArrayType::size()> &             face_info,
       const std::vector<unsigned int> &                         active_fe_index,
-      const Mapping<dim> &                                      mapping,
-      const std::vector<dealii::hp::QCollection<dim>> &         quad,
+      const std::shared_ptr<dealii::hp::MappingCollection<dim>> &mapping,
+      const std::vector<dealii::hp::QCollection<dim>> &          quad,
       const UpdateFlags update_flags_cells,
       const UpdateFlags update_flags_boundary_faces,
       const UpdateFlags update_flags_inner_faces,
       const UpdateFlags update_flags_faces_by_cells)
     {
       clear();
-      this->mapping = &mapping;
+      this->mapping_collection = mapping;
+      this->mapping            = &mapping->operator[](0);
 
       cell_data.resize(quad.size());
       face_data.resize(quad.size());
@@ -435,16 +437,16 @@ namespace internal
       // In case we have no hp adaptivity (active_fe_index is empty), we have
       // cells, and the mapping is MappingQGeneric or a derived class, we can
       // use the fast method.
-      if (active_fe_index.empty() && !cells.empty() &&
-          dynamic_cast<const MappingQGeneric<dim> *>(&mapping))
+      if (active_fe_index.empty() && !cells.empty() && mapping->size() == 1 &&
+          dynamic_cast<const MappingQGeneric<dim> *>(&mapping->operator[](0)))
         compute_mapping_q(tria, cells, face_info.faces);
       else
         {
           // Could call these functions in parallel, but not useful because
           // the work inside is nicely split up already
-          initialize_cells(tria, cells, active_fe_index, mapping);
-          initialize_faces(tria, cells, face_info.faces, mapping);
-          initialize_faces_by_cells(tria, cells, mapping);
+          initialize_cells(tria, cells, active_fe_index, *mapping);
+          initialize_faces(tria, cells, face_info.faces, *mapping);
+          initialize_faces_by_cells(tria, cells, *mapping);
         }
     }
 
@@ -457,7 +459,7 @@ namespace internal
       const std::vector<std::pair<unsigned int, unsigned int>> &cells,
       const FaceInfo<VectorizedArrayType::size()> &             face_info,
       const std::vector<unsigned int> &                         active_fe_index,
-      const Mapping<dim> &                                      mapping)
+      const std::shared_ptr<dealii::hp::MappingCollection<dim>> &mapping)
     {
       AssertDimension(cells.size() / VectorizedArrayType::size(),
                       cell_type.size());
@@ -469,18 +471,19 @@ namespace internal
       for (auto &data : face_data_by_cells)
         data.clear_data_fields();
 
-      this->mapping = &mapping;
+      this->mapping_collection = mapping;
+      this->mapping            = &mapping->operator[](0);
 
-      if (active_fe_index.empty() && !cells.empty() &&
-          dynamic_cast<const MappingQGeneric<dim> *>(&mapping))
+      if (active_fe_index.empty() && !cells.empty() && mapping->size() == 1 &&
+          dynamic_cast<const MappingQGeneric<dim> *>(&mapping->operator[](0)))
         compute_mapping_q(tria, cells, face_info.faces);
       else
         {
           // Could call these functions in parallel, but not useful because
           // the work inside is nicely split up already
-          initialize_cells(tria, cells, active_fe_index, mapping);
-          initialize_faces(tria, cells, face_info.faces, mapping);
-          initialize_faces_by_cells(tria, cells, mapping);
+          initialize_cells(tria, cells, active_fe_index, *mapping);
+          initialize_faces(tria, cells, face_info.faces, *mapping);
+          initialize_faces_by_cells(tria, cells, *mapping);
         }
     }
 
@@ -843,12 +846,15 @@ namespace internal
         const dealii::Triangulation<dim> &                        tria,
         const std::vector<std::pair<unsigned int, unsigned int>> &cells,
         const std::vector<unsigned int> &              active_fe_index,
-        const Mapping<dim> &                           mapping,
+        const dealii::hp::MappingCollection<dim> &     mapping_in,
         MappingInfo<dim, Number, VectorizedArrayType> &mapping_info,
         std::pair<std::vector<
                     MappingInfoStorage<dim, dim, Number, VectorizedArrayType>>,
                   CompressedCellData<dim, Number, VectorizedArrayType>> &data)
       {
+        AssertDimension(mapping_in.size(), 1);
+        const auto &mapping = mapping_in[0];
+
         FE_Nothing<dim> dummy_fe;
 
         // when we make comparisons about the size of Jacobians we need to
@@ -1515,7 +1521,7 @@ namespace internal
       const dealii::Triangulation<dim> &                        tria,
       const std::vector<std::pair<unsigned int, unsigned int>> &cells,
       const std::vector<unsigned int> &                         active_fe_index,
-      const Mapping<dim> &                                      mapping)
+      const dealii::hp::MappingCollection<dim> &                mapping)
     {
       const unsigned int n_cells = cells.size();
       const unsigned int n_lanes = VectorizedArrayType::size();
@@ -1748,13 +1754,16 @@ namespace internal
         const std::vector<std::pair<unsigned int, unsigned int>> &cells,
         const std::vector<FaceToCellTopology<VectorizedArrayType::size()>>
           &                                            faces,
-        const Mapping<dim> &                           mapping,
+        const dealii::hp::MappingCollection<dim> &     mapping_in,
         MappingInfo<dim, Number, VectorizedArrayType> &mapping_info,
         std::pair<
           std::vector<
             MappingInfoStorage<dim - 1, dim, Number, VectorizedArrayType>>,
           CompressedFaceData<dim, Number, VectorizedArrayType>> &data)
       {
+        AssertDimension(mapping_in.size(), 1);
+        const auto &mapping = mapping_in[0];
+
         FE_Nothing<dim> dummy_fe;
 
         std::vector<std::vector<std::shared_ptr<FEFaceValues<dim>>>>
@@ -2392,7 +2401,7 @@ namespace internal
       const dealii::Triangulation<dim> &                                  tria,
       const std::vector<std::pair<unsigned int, unsigned int>> &          cells,
       const std::vector<FaceToCellTopology<VectorizedArrayType::size()>> &faces,
-      const Mapping<dim> &mapping)
+      const dealii::hp::MappingCollection<dim> &mapping)
     {
       face_type.resize(faces.size(), general);
 
@@ -2591,8 +2600,11 @@ namespace internal
     {
       // step 1: extract quadrature point data with the data appropriate for
       // MappingQGeneric
+      AssertDimension(this->mapping_collection->size(), 1);
+
       const MappingQGeneric<dim> *mapping_q =
-        dynamic_cast<const MappingQGeneric<dim> *>(&*this->mapping);
+        dynamic_cast<const MappingQGeneric<dim> *>(
+          &this->mapping_collection->operator[](0));
       Assert(mapping_q != nullptr, ExcInternalError());
 
       const unsigned int mapping_degree = mapping_q->get_degree();
@@ -2937,7 +2949,7 @@ namespace internal
       // transitioned to extracting the information from cell quadrature
       // points but we need to figure out the correct indices of neighbors
       // within the list of arrays still
-      initialize_faces_by_cells(tria, cell_array, *this->mapping);
+      initialize_faces_by_cells(tria, cell_array, *this->mapping_collection);
     }
 
 
@@ -2947,11 +2959,15 @@ namespace internal
     MappingInfo<dim, Number, VectorizedArrayType>::initialize_faces_by_cells(
       const dealii::Triangulation<dim> &                        tria,
       const std::vector<std::pair<unsigned int, unsigned int>> &cells,
-      const Mapping<dim> &                                      mapping)
+      const dealii::hp::MappingCollection<dim> &                mapping_in)
     {
       if (update_flags_faces_by_cells == update_default)
         return;
 
+
+      AssertDimension(mapping_in.size(), 1);
+      const auto &mapping = mapping_in[0];
+
       const unsigned int n_quads = face_data_by_cells.size();
       const unsigned int n_lanes = VectorizedArrayType::size();
       UpdateFlags        update_flags =
index 165e5a23eecd7ff47d4c347259fcc61f5b4aa420..43d3f882619d8cb328149eb321d37d2791e0d0dc 100644 (file)
@@ -35,6 +35,7 @@
 #include <deal.II/grid/grid_tools.h>
 
 #include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/mapping_collection.h>
 #include <deal.II/hp/q_collection.h>
 
 #include <deal.II/lac/affine_constraints.h>
@@ -580,9 +581,9 @@ public:
    * not allowed. In that case, use the initialization function with
    * several DoFHandler arguments.
    */
-  template <typename QuadratureType, typename number2>
+  template <typename QuadratureType, typename number2, typename MappingType>
   void
-  reinit(const Mapping<dim> &              mapping,
+  reinit(const MappingType &               mapping,
          const DoFHandler<dim> &           dof_handler,
          const AffineConstraints<number2> &constraint,
          const QuadratureType &            quad,
@@ -622,9 +623,9 @@ public:
    * can be sets independently from the number of DoFHandlers, when several
    * elements are always integrated with the same quadrature formula.
    */
-  template <typename QuadratureType, typename number2>
+  template <typename QuadratureType, typename number2, typename MappingType>
   void
-  reinit(const Mapping<dim> &                                   mapping,
+  reinit(const MappingType &                                    mapping,
          const std::vector<const DoFHandler<dim> *> &           dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const std::vector<QuadratureType> &                    quad,
@@ -635,9 +636,12 @@ public:
    *
    * @deprecated Use the overload taking a DoFHandler object instead.
    */
-  template <typename QuadratureType, typename number2, typename DoFHandlerType>
+  template <typename QuadratureType,
+            typename number2,
+            typename DoFHandlerType,
+            typename MappingType>
   DEAL_II_DEPRECATED void
-  reinit(const Mapping<dim> &                                   mapping,
+  reinit(const MappingType &                                    mapping,
          const std::vector<const DoFHandlerType *> &            dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const std::vector<QuadratureType> &                    quad,
@@ -675,9 +679,9 @@ public:
    * as might be necessary when several components in a vector-valued problem
    * are integrated together based on the same quadrature formula.
    */
-  template <typename QuadratureType, typename number2>
+  template <typename QuadratureType, typename number2, typename MappingType>
   void
-  reinit(const Mapping<dim> &                                   mapping,
+  reinit(const MappingType &                                    mapping,
          const std::vector<const DoFHandler<dim> *> &           dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const QuadratureType &                                 quad,
@@ -688,9 +692,12 @@ public:
    *
    * @deprecated Use the overload taking a DoFHandler object instead.
    */
-  template <typename QuadratureType, typename number2, typename DoFHandlerType>
+  template <typename QuadratureType,
+            typename number2,
+            typename DoFHandlerType,
+            typename MappingType>
   DEAL_II_DEPRECATED void
-  reinit(const Mapping<dim> &                                   mapping,
+  reinit(const MappingType &                                    mapping,
          const std::vector<const DoFHandlerType *> &            dof_handler,
          const std::vector<const AffineConstraints<number2> *> &constraint,
          const QuadratureType &                                 quad,
@@ -742,6 +749,12 @@ public:
   void
   update_mapping(const Mapping<dim> &mapping);
 
+  /**
+   * Same as above but with hp::MappingCollection.
+   */
+  void
+  update_mapping(const std::shared_ptr<hp::MappingCollection<dim>> &mapping);
+
   /**
    * Clear all data fields and brings the class into a condition similar to
    * after having called the default constructor.
@@ -2030,7 +2043,7 @@ private:
   template <typename number2, int q_dim>
   void
   internal_reinit(
-    const Mapping<dim> &                                   mapping,
+    const std::shared_ptr<hp::MappingCollection<dim>> &    mapping,
     const std::vector<const DoFHandler<dim, dim> *> &      dof_handlers,
     const std::vector<const AffineConstraints<number2> *> &constraint,
     const std::vector<IndexSet> &                          locally_owned_set,
@@ -2853,7 +2866,8 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
   std::vector<hp::QCollection<dim>> quad_hp;
   quad_hp.emplace_back(quad);
 
-  internal_reinit(StaticMappingQ1<dim>::mapping,
+  internal_reinit(std::make_shared<hp::MappingCollection<dim>>(
+                    StaticMappingQ1<dim>::mapping),
                   dof_handlers,
                   constraints,
                   locally_owned_sets,
@@ -2864,10 +2878,10 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2>
+template <typename QuadratureType, typename number2, typename MappingType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
-  const Mapping<dim> &              mapping,
+  const MappingType &               mapping,
   const DoFHandler<dim> &           dof_handler,
   const AffineConstraints<number2> &constraints_in,
   const QuadratureType &            quad,
@@ -2887,7 +2901,7 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
   std::vector<hp::QCollection<dim>> quad_hp;
   quad_hp.emplace_back(quad);
 
-  internal_reinit(mapping,
+  internal_reinit(std::make_shared<hp::MappingCollection<dim>>(mapping),
                   dof_handlers,
                   constraints,
                   locally_owned_sets,
@@ -2913,7 +2927,9 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
   std::vector<hp::QCollection<dim>> quad_hp;
   for (unsigned int q = 0; q < quad.size(); ++q)
     quad_hp.emplace_back(quad[q]);
-  internal_reinit(StaticMappingQ1<dim>::mapping,
+
+  internal_reinit(std::make_shared<hp::MappingCollection<dim>>(
+                    StaticMappingQ1<dim>::mapping),
                   dof_handler,
                   constraint,
                   locally_owned_set,
@@ -2924,10 +2940,13 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2, typename DoFHandlerType>
+template <typename QuadratureType,
+          typename number2,
+          typename DoFHandlerType,
+          typename MappingType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
-  const Mapping<dim> &                                   mapping,
+  const MappingType &                                    mapping,
   const std::vector<const DoFHandlerType *> &            dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const std::vector<QuadratureType> &                    quad,
@@ -2961,7 +2980,9 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
       dof_handler, additional_data.mg_level);
   std::vector<hp::QCollection<dim>> quad_hp;
   quad_hp.emplace_back(quad);
-  internal_reinit(StaticMappingQ1<dim>::mapping,
+
+  internal_reinit(std::make_shared<hp::MappingCollection<dim>>(
+                    StaticMappingQ1<dim>::mapping),
                   dof_handler,
                   constraint,
                   locally_owned_set,
@@ -2994,10 +3015,10 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2>
+template <typename QuadratureType, typename number2, typename MappingType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
-  const Mapping<dim> &                                   mapping,
+  const MappingType &                                    mapping,
   const std::vector<const DoFHandler<dim> *> &           dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const QuadratureType &                                 quad,
@@ -3009,7 +3030,8 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
       dof_handler, additional_data.mg_level);
   std::vector<hp::QCollection<dim>> quad_hp;
   quad_hp.emplace_back(quad);
-  internal_reinit(mapping,
+
+  internal_reinit(std::make_shared<hp::MappingCollection<dim>>(mapping),
                   dof_handler,
                   constraint,
                   locally_owned_set,
@@ -3020,10 +3042,10 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2>
+template <typename QuadratureType, typename number2, typename MappingType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
-  const Mapping<dim> &                                   mapping,
+  const MappingType &                                    mapping,
   const std::vector<const DoFHandler<dim> *> &           dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const std::vector<QuadratureType> &                    quad,
@@ -3036,7 +3058,8 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
   std::vector<hp::QCollection<dim>> quad_hp;
   for (unsigned int q = 0; q < quad.size(); ++q)
     quad_hp.emplace_back(quad[q]);
-  internal_reinit(mapping,
+
+  internal_reinit(std::make_shared<hp::MappingCollection<dim>>(mapping),
                   dof_handler,
                   constraint,
                   locally_owned_set,
@@ -3047,10 +3070,13 @@ MatrixFree<dim, Number, VectorizedArrayType>::reinit(
 
 
 template <int dim, typename Number, typename VectorizedArrayType>
-template <typename QuadratureType, typename number2, typename DoFHandlerType>
+template <typename QuadratureType,
+          typename number2,
+          typename DoFHandlerType,
+          typename MappingType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::reinit(
-  const Mapping<dim> &                                   mapping,
+  const MappingType &                                    mapping,
   const std::vector<const DoFHandlerType *> &            dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const QuadratureType &                                 quad,
index 3b5643dd84660a4baae75a8875978b6501ad1d09..01809a8ea98ec3ec510a1fb003902ebbb9bcdf99 100644 (file)
@@ -290,7 +290,7 @@ template <int dim, typename Number, typename VectorizedArrayType>
 template <typename number2, int q_dim>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::internal_reinit(
-  const Mapping<dim> &                                   mapping,
+  const std::shared_ptr<hp::MappingCollection<dim>> &    mapping,
   const std::vector<const DoFHandler<dim, dim> *> &      dof_handler,
   const std::vector<const AffineConstraints<number2> *> &constraint,
   const std::vector<IndexSet> &                          locally_owned_dofs,
@@ -466,6 +466,16 @@ template <int dim, typename Number, typename VectorizedArrayType>
 void
 MatrixFree<dim, Number, VectorizedArrayType>::update_mapping(
   const Mapping<dim> &mapping)
+{
+  update_mapping(std::make_shared<hp::MappingCollection<dim>>(mapping));
+}
+
+
+
+template <int dim, typename Number, typename VectorizedArrayType>
+void
+MatrixFree<dim, Number, VectorizedArrayType>::update_mapping(
+  const std::shared_ptr<hp::MappingCollection<dim>> &mapping)
 {
   AssertDimension(shape_info.size(1), mapping_info.cell_data.size());
   mapping_info.update_mapping(dof_handlers[0]->get_triangulation(),
index aaacd0001c7697d771d5e2c4e3bc30c1d7c57f99..0ae3ae5d0dc2514186d5ca1b2aa505d7ef693bbc 100644 (file)
@@ -35,7 +35,7 @@ for (deal_II_dimension : DIMENSIONS;
                              deal_II_scalar_vectorized::value_type,
                              deal_II_scalar_vectorized>::
       internal_reinit<deal_II_scalar_vectorized::value_type>(
-        const Mapping<deal_II_dimension> &,
+        const std::shared_ptr<hp::MappingCollection<deal_II_dimension>> &,
         const std::vector<const DoFHandler<deal_II_dimension> *> &,
         const std::vector<
           const AffineConstraints<deal_II_scalar_vectorized::value_type> *> &,
@@ -52,7 +52,7 @@ for (deal_II_dimension : DIMENSIONS;
                              deal_II_float_vectorized::value_type,
                              deal_II_float_vectorized>::
       internal_reinit<double>(
-        const Mapping<deal_II_dimension> &,
+        const std::shared_ptr<hp::MappingCollection<deal_II_dimension>> &,
         const std::vector<const DoFHandler<deal_II_dimension> *> &,
         const std::vector<const AffineConstraints<double> *> &,
         const std::vector<IndexSet> &,

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.