]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove DoFHandlerType from internal::DataOutImplementation::{DataEntryBase|DataEntry... 11182/head
authorMarc Fehling <mafehling.git@gmail.com>
Sat, 14 Nov 2020 02:31:36 +0000 (19:31 -0700)
committerMarc Fehling <mafehling.git@gmail.com>
Sat, 14 Nov 2020 02:35:28 +0000 (19:35 -0700)
include/deal.II/numerics/data_out_dof_data.h
include/deal.II/numerics/data_out_dof_data.templates.h
source/numerics/data_out_dof_data.inst.in
source/numerics/data_out_dof_data_codim.inst.in

index c326b475fa61a549dc496605fad4a9d149c3287c..86fe1fcea07583fdd7f04a572b76385a6153ec30 100644 (file)
@@ -218,7 +218,7 @@ namespace internal
      * <a href="https://www.artima.com/cppsource/type_erasure.html">type
      * erasure</a> design pattern.
      */
-    template <typename DoFHandlerType>
+    template <int dim, int spacedim>
     class DataEntryBase
     {
     public:
@@ -227,8 +227,8 @@ namespace internal
        * the vector and their interpretation as scalar or vector data. This
        * constructor assumes that no postprocessor is going to be used.
        */
-      DataEntryBase(const DoFHandlerType *          dofs,
-                    const std::vector<std::string> &names,
+      DataEntryBase(const DoFHandler<dim, spacedim> *dofs,
+                    const std::vector<std::string> & names,
                     const std::vector<
                       DataComponentInterpretation::DataComponentInterpretation>
                       &data_component_interpretation);
@@ -238,9 +238,8 @@ namespace internal
        * case, the names and vector declarations are going to be acquired from
        * the postprocessor.
        */
-      DataEntryBase(const DoFHandlerType *dofs,
-                    const DataPostprocessor<DoFHandlerType::space_dimension>
-                      *data_postprocessor);
+      DataEntryBase(const DoFHandler<dim, spacedim> *  dofs,
+                    const DataPostprocessor<spacedim> *data_postprocessor);
 
       /**
        * Destructor made virtual.
@@ -260,11 +259,9 @@ namespace internal
        * from the vector we actually store.
        */
       virtual void
-      get_function_values(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<double> &patch_values) const = 0;
+      get_function_values(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                          const ComponentExtractor           extract_component,
+                          std::vector<double> &patch_values) const = 0;
 
       /**
        * Given a FEValuesBase object, extract the values on the present cell
@@ -273,9 +270,8 @@ namespace internal
        */
       virtual void
       get_function_values(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
+        const FEValuesBase<dim, spacedim> &  fe_patch_values,
+        const ComponentExtractor             extract_component,
         std::vector<dealii::Vector<double>> &patch_values_system) const = 0;
 
       /**
@@ -284,11 +280,9 @@ namespace internal
        */
       virtual void
       get_function_gradients(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<Tensor<1, DoFHandlerType::space_dimension>>
-          &patch_gradients) const = 0;
+        const FEValuesBase<dim, spacedim> &fe_patch_values,
+        const ComponentExtractor           extract_component,
+        std::vector<Tensor<1, spacedim>> & patch_gradients) const = 0;
 
       /**
        * Given a FEValuesBase object, extract the gradients on the present
@@ -296,12 +290,10 @@ namespace internal
        * as the one above but for vector-valued finite elements.
        */
       virtual void
-      get_function_gradients(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<std::vector<Tensor<1, DoFHandlerType::space_dimension>>>
-          &patch_gradients_system) const = 0;
+      get_function_gradients(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                             const ComponentExtractor extract_component,
+                             std::vector<std::vector<Tensor<1, spacedim>>>
+                               &patch_gradients_system) const = 0;
 
       /**
        * Given a FEValuesBase object, extract the second derivatives on the
@@ -309,11 +301,9 @@ namespace internal
        */
       virtual void
       get_function_hessians(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<Tensor<2, DoFHandlerType::space_dimension>> &patch_hessians)
-        const = 0;
+        const FEValuesBase<dim, spacedim> &fe_patch_values,
+        const ComponentExtractor           extract_component,
+        std::vector<Tensor<2, spacedim>> & patch_hessians) const = 0;
 
       /**
        * Given a FEValuesBase object, extract the second derivatives on the
@@ -321,12 +311,10 @@ namespace internal
        * the same as the one above but for vector-valued finite elements.
        */
       virtual void
-      get_function_hessians(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<std::vector<Tensor<2, DoFHandlerType::space_dimension>>>
-          &patch_hessians_system) const = 0;
+      get_function_hessians(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                            const ComponentExtractor extract_component,
+                            std::vector<std::vector<Tensor<2, spacedim>>>
+                              &patch_hessians_system) const = 0;
 
       /**
        * Return whether the data represented by (a derived class of) this object
@@ -351,7 +339,7 @@ namespace internal
       /**
        * Pointer to the DoFHandler object that the vector is based on.
        */
-      SmartPointer<const DoFHandlerType> dof_handler;
+      SmartPointer<const DoFHandler<dim, spacedim>> dof_handler;
 
       /**
        * Names of the components of this data vector.
@@ -371,9 +359,7 @@ namespace internal
        * Pointer to a DataPostprocessing object which shall be applied to this
        * data vector.
        */
-      SmartPointer<
-        const dealii::DataPostprocessor<DoFHandlerType::space_dimension>>
-        postprocessor;
+      SmartPointer<const dealii::DataPostprocessor<spacedim>> postprocessor;
 
       /**
        * Number of output variables this dataset provides (either number of
@@ -428,10 +414,9 @@ namespace internal
 
       ParallelDataBase(const ParallelDataBase &data);
 
-      template <typename DoFHandlerType>
       void
       reinit_all_fe_values(
-        std::vector<std::shared_ptr<DataEntryBase<DoFHandlerType>>> &dof_data,
+        std::vector<std::shared_ptr<DataEntryBase<dim, spacedim>>> &dof_data,
         const typename dealii::Triangulation<dim, spacedim>::cell_iterator
           &                cell,
         const unsigned int face = numbers::invalid_unsigned_int);
@@ -994,15 +979,17 @@ protected:
   /**
    * List of data elements with vectors of values for each degree of freedom.
    */
-  std::vector<std::shared_ptr<
-    internal::DataOutImplementation::DataEntryBase<DoFHandlerType>>>
+  std::vector<std::shared_ptr<internal::DataOutImplementation::DataEntryBase<
+    DoFHandlerType::dimension,
+    DoFHandlerType::space_dimension>>>
     dof_data;
 
   /**
    * List of data elements with vectors of values for each cell.
    */
-  std::vector<std::shared_ptr<
-    internal::DataOutImplementation::DataEntryBase<DoFHandlerType>>>
+  std::vector<std::shared_ptr<internal::DataOutImplementation::DataEntryBase<
+    DoFHandlerType::dimension,
+    DoFHandlerType::space_dimension>>>
     cell_data;
 
   /**
index 19332c8938f9e340d98d18ec234131e25ead0327..f7b9b6fbe2b9e2aaefbdac97ea9dc23a358b6f87 100644 (file)
@@ -347,10 +347,9 @@ namespace internal
 
 
     template <int dim, int spacedim>
-    template <typename DoFHandlerType>
     void
     ParallelDataBase<dim, spacedim>::reinit_all_fe_values(
-      std::vector<std::shared_ptr<DataEntryBase<DoFHandlerType>>> &dof_data,
+      std::vector<std::shared_ptr<DataEntryBase<dim, spacedim>>> &dof_data,
       const typename dealii::Triangulation<dim, spacedim>::cell_iterator &cell,
       const unsigned int                                                  face)
     {
@@ -365,11 +364,11 @@ namespace internal
             {
               if (cell->active())
                 {
-                  typename DoFHandlerType::active_cell_iterator dh_cell(
-                    &cell->get_triangulation(),
-                    cell->level(),
-                    cell->index(),
-                    dof_data[dataset]->dof_handler);
+                  typename DoFHandler<dim, spacedim>::active_cell_iterator
+                    dh_cell(&cell->get_triangulation(),
+                            cell->level(),
+                            cell->index(),
+                            dof_data[dataset]->dof_handler);
                   if (x_fe_values.empty())
                     {
                       AssertIndexRange(face, GeometryInfo<dim>::faces_per_cell);
@@ -596,19 +595,18 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType>
-    DataEntryBase<DoFHandlerType>::DataEntryBase(
-      const DoFHandlerType *          dofs,
-      const std::vector<std::string> &names_in,
+    template <int dim, int spacedim>
+    DataEntryBase<dim, spacedim>::DataEntryBase(
+      const DoFHandler<dim, spacedim> *dofs,
+      const std::vector<std::string> & names_in,
       const std::vector<
         DataComponentInterpretation::DataComponentInterpretation>
         &data_component_interpretation)
-      : dof_handler(dofs,
-                    typeid(
-                      dealii::DataOut_DoFData<DoFHandlerType,
-                                              DoFHandlerType::dimension,
-                                              DoFHandlerType::space_dimension>)
-                      .name())
+      : dof_handler(
+          dofs,
+          typeid(
+            dealii::DataOut_DoFData<DoFHandler<dim, spacedim>, dim, spacedim>)
+            .name())
       , names(names_in)
       , data_component_interpretation(data_component_interpretation)
       , postprocessor(nullptr, typeid(*this).name())
@@ -635,17 +633,15 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType>
-    DataEntryBase<DoFHandlerType>::DataEntryBase(
-      const DoFHandlerType *dofs,
-      const DataPostprocessor<DoFHandlerType::space_dimension>
-        *data_postprocessor)
-      : dof_handler(dofs,
-                    typeid(
-                      dealii::DataOut_DoFData<DoFHandlerType,
-                                              DoFHandlerType::dimension,
-                                              DoFHandlerType::space_dimension>)
-                      .name())
+    template <int dim, int spacedim>
+    DataEntryBase<dim, spacedim>::DataEntryBase(
+      const DoFHandler<dim, spacedim> *  dofs,
+      const DataPostprocessor<spacedim> *data_postprocessor)
+      : dof_handler(
+          dofs,
+          typeid(
+            dealii::DataOut_DoFData<DoFHandler<dim, spacedim>, dim, spacedim>)
+            .name())
       , names(data_postprocessor->get_names())
       , data_component_interpretation(
           data_postprocessor->get_data_component_interpretation())
@@ -679,8 +675,8 @@ namespace internal
      * Class that stores a pointer to a vector of type equal to the template
      * argument, and provides the functions to extract data from it.
      */
-    template <typename DoFHandlerType, typename VectorType>
-    class DataEntry : public DataEntryBase<DoFHandlerType>
+    template <int dim, int spacedim, typename VectorType>
+    class DataEntry : public DataEntryBase<dim, spacedim>
     {
     public:
       /**
@@ -688,9 +684,9 @@ namespace internal
        * the vector and their interpretation as scalar or vector data. This
        * constructor assumes that no postprocessor is going to be used.
        */
-      DataEntry(const DoFHandlerType *          dofs,
-                const VectorType *              data,
-                const std::vector<std::string> &names,
+      DataEntry(const DoFHandler<dim, spacedim> *dofs,
+                const VectorType *               data,
+                const std::vector<std::string> & names,
                 const std::vector<
                   DataComponentInterpretation::DataComponentInterpretation>
                   &data_component_interpretation);
@@ -700,10 +696,9 @@ namespace internal
        * case, the names and vector declarations are going to be acquired from
        * the postprocessor.
        */
-      DataEntry(const DoFHandlerType *dofs,
-                const VectorType *    data,
-                const DataPostprocessor<DoFHandlerType::space_dimension>
-                  *data_postprocessor);
+      DataEntry(const DoFHandler<dim, spacedim> *  dofs,
+                const VectorType *                 data,
+                const DataPostprocessor<spacedim> *data_postprocessor);
 
       /**
        * Assuming that the stored vector is a cell vector, extract the given
@@ -719,11 +714,9 @@ namespace internal
        * from the vector we actually store.
        */
       virtual void
-      get_function_values(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<double> &patch_values) const override;
+      get_function_values(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                          const ComponentExtractor           extract_component,
+                          std::vector<double> &patch_values) const override;
 
       /**
        * Given a FEValuesBase object, extract the values on the present cell
@@ -731,12 +724,10 @@ namespace internal
        * one above but for vector-valued finite elements.
        */
       virtual void
-      get_function_values(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<dealii::Vector<double>> &patch_values_system)
-        const override;
+      get_function_values(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                          const ComponentExtractor           extract_component,
+                          std::vector<dealii::Vector<double>>
+                            &patch_values_system) const override;
 
       /**
        * Given a FEValuesBase object, extract the gradients on the present
@@ -744,11 +735,9 @@ namespace internal
        */
       virtual void
       get_function_gradients(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<Tensor<1, DoFHandlerType::space_dimension>>
-          &patch_gradients) const override;
+        const FEValuesBase<dim, spacedim> &fe_patch_values,
+        const ComponentExtractor           extract_component,
+        std::vector<Tensor<1, spacedim>> & patch_gradients) const override;
 
       /**
        * Given a FEValuesBase object, extract the gradients on the present
@@ -756,12 +745,10 @@ namespace internal
        * as the one above but for vector-valued finite elements.
        */
       virtual void
-      get_function_gradients(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<std::vector<Tensor<1, DoFHandlerType::space_dimension>>>
-          &patch_gradients_system) const override;
+      get_function_gradients(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                             const ComponentExtractor extract_component,
+                             std::vector<std::vector<Tensor<1, spacedim>>>
+                               &patch_gradients_system) const override;
 
       /**
        * Given a FEValuesBase object, extract the second derivatives on the
@@ -769,11 +756,9 @@ namespace internal
        */
       virtual void
       get_function_hessians(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<Tensor<2, DoFHandlerType::space_dimension>> &patch_hessians)
-        const override;
+        const FEValuesBase<dim, spacedim> &fe_patch_values,
+        const ComponentExtractor           extract_component,
+        std::vector<Tensor<2, spacedim>> & patch_hessians) const override;
 
       /**
        * Given a FEValuesBase object, extract the second derivatives on the
@@ -781,12 +766,10 @@ namespace internal
        * the same as the one above but for vector-valued finite elements.
        */
       virtual void
-      get_function_hessians(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<std::vector<Tensor<2, DoFHandlerType::space_dimension>>>
-          &patch_hessians_system) const override;
+      get_function_hessians(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                            const ComponentExtractor extract_component,
+                            std::vector<std::vector<Tensor<2, spacedim>>>
+                              &patch_hessians_system) const override;
 
       /**
        * Return whether the data represented by (a derived class of) this object
@@ -818,37 +801,34 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
-    DataEntry<DoFHandlerType, VectorType>::DataEntry(
-      const DoFHandlerType *          dofs,
-      const VectorType *              data,
-      const std::vector<std::string> &names,
+    template <int dim, int spacedim, typename VectorType>
+    DataEntry<dim, spacedim, VectorType>::DataEntry(
+      const DoFHandler<dim, spacedim> *dofs,
+      const VectorType *               data,
+      const std::vector<std::string> & names,
       const std::vector<
         DataComponentInterpretation::DataComponentInterpretation>
         &data_component_interpretation)
-      : DataEntryBase<DoFHandlerType>(dofs,
-                                      names,
-                                      data_component_interpretation)
+      : DataEntryBase<dim, spacedim>(dofs, names, data_component_interpretation)
       , vector(data)
     {}
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
-    DataEntry<DoFHandlerType, VectorType>::DataEntry(
-      const DoFHandlerType *dofs,
-      const VectorType *    data,
-      const DataPostprocessor<DoFHandlerType::space_dimension>
-        *data_postprocessor)
-      : DataEntryBase<DoFHandlerType>(dofs, data_postprocessor)
+    template <int dim, int spacedim, typename VectorType>
+    DataEntry<dim, spacedim, VectorType>::DataEntry(
+      const DoFHandler<dim, spacedim> *  dofs,
+      const VectorType *                 data,
+      const DataPostprocessor<spacedim> *data_postprocessor)
+      : DataEntryBase<dim, spacedim>(dofs, data_postprocessor)
       , vector(data)
     {}
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     double
-    DataEntry<DoFHandlerType, VectorType>::get_cell_data_value(
+    DataEntry<dim, spacedim, VectorType>::get_cell_data_value(
       const unsigned int       cell_number,
       const ComponentExtractor extract_component) const
     {
@@ -859,12 +839,11 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    DataEntry<DoFHandlerType, VectorType>::get_function_values(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
+    DataEntry<dim, spacedim, VectorType>::get_function_values(
+      const FEValuesBase<dim, spacedim> &  fe_patch_values,
+      const ComponentExtractor             extract_component,
       std::vector<dealii::Vector<double>> &patch_values_system) const
     {
       if (typeid(typename VectorType::value_type) == typeid(double))
@@ -916,13 +895,12 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    DataEntry<DoFHandlerType, VectorType>::get_function_values(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
-      std::vector<double> &                                patch_values) const
+    DataEntry<dim, spacedim, VectorType>::get_function_values(
+      const FEValuesBase<dim, spacedim> &fe_patch_values,
+      const ComponentExtractor           extract_component,
+      std::vector<double> &              patch_values) const
     {
       if (typeid(typename VectorType::value_type) == typeid(double))
         {
@@ -952,14 +930,13 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    DataEntry<DoFHandlerType, VectorType>::get_function_gradients(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
-      std::vector<std::vector<Tensor<1, DoFHandlerType::space_dimension>>>
-        &patch_gradients_system) const
+    DataEntry<dim, spacedim, VectorType>::get_function_gradients(
+      const FEValuesBase<dim, spacedim> &            fe_patch_values,
+      const ComponentExtractor                       extract_component,
+      std::vector<std::vector<Tensor<1, spacedim>>> &patch_gradients_system)
+      const
     {
       if (typeid(typename VectorType::value_type) == typeid(double))
         {
@@ -973,10 +950,8 @@ namespace internal
             // above, this is the identity cast whenever the code is executed,
             // but the cast is necessary to allow compilation even if we don't
             // get here
-            reinterpret_cast<std::vector<
-              std::vector<Tensor<1,
-                                 DoFHandlerType::space_dimension,
-                                 typename VectorType::value_type>>> &>(
+            reinterpret_cast<std::vector<std::vector<
+              Tensor<1, spacedim, typename VectorType::value_type>>> &>(
               patch_gradients_system));
         }
       else
@@ -991,9 +966,8 @@ namespace internal
           const unsigned int n_eval_points =
             fe_patch_values.n_quadrature_points;
 
-          std::vector<std::vector<Tensor<1,
-                                         DoFHandlerType::space_dimension,
-                                         typename VectorType::value_type>>>
+          std::vector<
+            std::vector<Tensor<1, spacedim, typename VectorType::value_type>>>
             tmp(n_eval_points);
           for (unsigned int i = 0; i < n_eval_points; i++)
             tmp[i].resize(n_components);
@@ -1014,14 +988,12 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    DataEntry<DoFHandlerType, VectorType>::get_function_gradients(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
-      std::vector<Tensor<1, DoFHandlerType::space_dimension>> &patch_gradients)
-      const
+    DataEntry<dim, spacedim, VectorType>::get_function_gradients(
+      const FEValuesBase<dim, spacedim> &fe_patch_values,
+      const ComponentExtractor           extract_component,
+      std::vector<Tensor<1, spacedim>> & patch_gradients) const
     {
       if (typeid(typename VectorType::value_type) == typeid(double))
         {
@@ -1035,18 +1007,13 @@ namespace internal
             // above, this is the identity cast whenever the code is executed,
             // but the cast is necessary to allow compilation even if we don't
             // get here
-            reinterpret_cast<
-              std::vector<Tensor<1,
-                                 DoFHandlerType::space_dimension,
-                                 typename VectorType::value_type>> &>(
+            reinterpret_cast<std::vector<
+              Tensor<1, spacedim, typename VectorType::value_type>> &>(
               patch_gradients));
         }
       else
         {
-          std::vector<Tensor<1,
-                             DoFHandlerType::space_dimension,
-                             typename VectorType::value_type>>
-            tmp;
+          std::vector<Tensor<1, spacedim, typename VectorType::value_type>> tmp;
           tmp.resize(patch_gradients.size());
 
           fe_patch_values.get_function_gradients(*vector, tmp);
@@ -1058,14 +1025,13 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    DataEntry<DoFHandlerType, VectorType>::get_function_hessians(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
-      std::vector<std::vector<Tensor<2, DoFHandlerType::space_dimension>>>
-        &patch_hessians_system) const
+    DataEntry<dim, spacedim, VectorType>::get_function_hessians(
+      const FEValuesBase<dim, spacedim> &            fe_patch_values,
+      const ComponentExtractor                       extract_component,
+      std::vector<std::vector<Tensor<2, spacedim>>> &patch_hessians_system)
+      const
     {
       if (typeid(typename VectorType::value_type) == typeid(double))
         {
@@ -1079,10 +1045,8 @@ namespace internal
             // above, this is the identity cast whenever the code is executed,
             // but the cast is necessary to allow compilation even if we don't
             // get here
-            reinterpret_cast<std::vector<
-              std::vector<Tensor<2,
-                                 DoFHandlerType::space_dimension,
-                                 typename VectorType::value_type>>> &>(
+            reinterpret_cast<std::vector<std::vector<
+              Tensor<2, spacedim, typename VectorType::value_type>>> &>(
               patch_hessians_system));
         }
       else
@@ -1097,9 +1061,8 @@ namespace internal
           const unsigned int n_eval_points =
             fe_patch_values.n_quadrature_points;
 
-          std::vector<std::vector<Tensor<2,
-                                         DoFHandlerType::space_dimension,
-                                         typename VectorType::value_type>>>
+          std::vector<
+            std::vector<Tensor<2, spacedim, typename VectorType::value_type>>>
             tmp(n_eval_points);
           for (unsigned int i = 0; i < n_eval_points; i++)
             tmp[i].resize(n_components);
@@ -1120,14 +1083,12 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    DataEntry<DoFHandlerType, VectorType>::get_function_hessians(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
-      std::vector<Tensor<2, DoFHandlerType::space_dimension>> &patch_hessians)
-      const
+    DataEntry<dim, spacedim, VectorType>::get_function_hessians(
+      const FEValuesBase<dim, spacedim> &fe_patch_values,
+      const ComponentExtractor           extract_component,
+      std::vector<Tensor<2, spacedim>> & patch_hessians) const
     {
       if (typeid(typename VectorType::value_type) == typeid(double))
         {
@@ -1141,18 +1102,14 @@ namespace internal
             // above, this is the identity cast whenever the code is executed,
             // but the cast is necessary to allow compilation even if we don't
             // get here
-            reinterpret_cast<
-              std::vector<Tensor<2,
-                                 DoFHandlerType::space_dimension,
-                                 typename VectorType::value_type>> &>(
+            reinterpret_cast<std::vector<
+              Tensor<2, spacedim, typename VectorType::value_type>> &>(
               patch_hessians));
         }
       else
         {
-          std::vector<Tensor<2,
-                             DoFHandlerType::space_dimension,
-                             typename VectorType::value_type>>
-            tmp(patch_hessians.size());
+          std::vector<Tensor<2, spacedim, typename VectorType::value_type>> tmp(
+            patch_hessians.size());
 
           fe_patch_values.get_function_hessians(*vector, tmp);
 
@@ -1163,18 +1120,18 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     bool
-    DataEntry<DoFHandlerType, VectorType>::is_complex_valued() const
+    DataEntry<dim, spacedim, VectorType>::is_complex_valued() const
     {
       return numbers::NumberTraits<typename VectorType::value_type>::is_complex;
     }
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     std::size_t
-    DataEntry<DoFHandlerType, VectorType>::memory_consumption() const
+    DataEntry<dim, spacedim, VectorType>::memory_consumption() const
     {
       return (sizeof(vector) +
               MemoryConsumption::memory_consumption(this->names));
@@ -1182,9 +1139,9 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    DataEntry<DoFHandlerType, VectorType>::clear()
+    DataEntry<dim, spacedim, VectorType>::clear()
     {
       vector            = nullptr;
       this->dof_handler = nullptr;
@@ -1198,19 +1155,19 @@ namespace internal
      * MGLevelObject<VectorType> given on the specific level instead of
      * interpolating data to coarser cells.
      */
-    template <typename DoFHandlerType, typename VectorType>
-    class MGDataEntry : public DataEntryBase<DoFHandlerType>
+    template <int dim, int spacedim, typename VectorType>
+    class MGDataEntry : public DataEntryBase<dim, spacedim>
     {
     public:
-      MGDataEntry(const DoFHandlerType *           dofs,
+      MGDataEntry(const DoFHandler<dim, spacedim> *dofs,
                   const MGLevelObject<VectorType> *vectors,
                   const std::vector<std::string> & names,
                   const std::vector<
                     DataComponentInterpretation::DataComponentInterpretation>
                     &data_component_interpretation)
-        : DataEntryBase<DoFHandlerType>(dofs,
-                                        names,
-                                        data_component_interpretation)
+        : DataEntryBase<dim, spacedim>(dofs,
+                                       names,
+                                       data_component_interpretation)
         , vectors(vectors)
       {}
 
@@ -1220,11 +1177,9 @@ namespace internal
         const ComponentExtractor extract_component) const override;
 
       virtual void
-      get_function_values(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<double> &patch_values) const override;
+      get_function_values(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                          const ComponentExtractor           extract_component,
+                          std::vector<double> &patch_values) const override;
 
       /**
        * Given a FEValuesBase object, extract the values on the present cell
@@ -1232,12 +1187,10 @@ namespace internal
        * one above but for vector-valued finite elements.
        */
       virtual void
-      get_function_values(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension> &fe_patch_values,
-        const ComponentExtractor                             extract_component,
-        std::vector<dealii::Vector<double>> &patch_values_system)
-        const override;
+      get_function_values(const FEValuesBase<dim, spacedim> &fe_patch_values,
+                          const ComponentExtractor           extract_component,
+                          std::vector<dealii::Vector<double>>
+                            &patch_values_system) const override;
 
       /**
        * Given a FEValuesBase object, extract the gradients on the present
@@ -1245,12 +1198,9 @@ namespace internal
        */
       virtual void
       get_function_gradients(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension>
-          & /*fe_patch_values*/,
+        const FEValuesBase<dim, spacedim> & /*fe_patch_values*/,
         const ComponentExtractor /*extract_component*/,
-        std::vector<Tensor<1, DoFHandlerType::space_dimension>>
-          & /*patch_gradients*/) const override
+        std::vector<Tensor<1, spacedim>> & /*patch_gradients*/) const override
       {
         Assert(false, ExcNotImplemented());
       }
@@ -1262,11 +1212,9 @@ namespace internal
        */
       virtual void
       get_function_gradients(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension>
-          & /*fe_patch_values*/,
+        const FEValuesBase<dim, spacedim> & /*fe_patch_values*/,
         const ComponentExtractor /*extract_component*/,
-        std::vector<std::vector<Tensor<1, DoFHandlerType::space_dimension>>>
+        std::vector<std::vector<Tensor<1, spacedim>>>
           & /*patch_gradients_system*/) const override
       {
         Assert(false, ExcNotImplemented());
@@ -1279,12 +1227,9 @@ namespace internal
        */
       virtual void
       get_function_hessians(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension>
-          & /*fe_patch_values*/,
+        const FEValuesBase<dim, spacedim> & /*fe_patch_values*/,
         const ComponentExtractor /*extract_component*/,
-        std::vector<Tensor<2, DoFHandlerType::space_dimension>>
-          & /*patch_hessians*/) const override
+        std::vector<Tensor<2, spacedim>> & /*patch_hessians*/) const override
       {
         Assert(false, ExcNotImplemented());
       }
@@ -1296,11 +1241,9 @@ namespace internal
        */
       virtual void
       get_function_hessians(
-        const FEValuesBase<DoFHandlerType::dimension,
-                           DoFHandlerType::space_dimension>
-          & /*fe_patch_values*/,
+        const FEValuesBase<dim, spacedim> & /*fe_patch_values*/,
         const ComponentExtractor /*extract_component*/,
-        std::vector<std::vector<Tensor<2, DoFHandlerType::space_dimension>>>
+        std::vector<std::vector<Tensor<2, spacedim>>>
           & /*patch_hessians_system*/) const override
       {
         Assert(false, ExcNotImplemented());
@@ -1346,9 +1289,9 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     double
-    MGDataEntry<DoFHandlerType, VectorType>::get_cell_data_value(
+    MGDataEntry<dim, spacedim, VectorType>::get_cell_data_value(
       const unsigned int       cell_number,
       const ComponentExtractor extract_component) const
     {
@@ -1361,18 +1304,17 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    MGDataEntry<DoFHandlerType, VectorType>::get_function_values(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
-      std::vector<double> &                                patch_values) const
+    MGDataEntry<dim, spacedim, VectorType>::get_function_values(
+      const FEValuesBase<dim, spacedim> &fe_patch_values,
+      const ComponentExtractor           extract_component,
+      std::vector<double> &              patch_values) const
     {
       Assert(extract_component == ComponentExtractor::real_part,
              ExcNotImplemented());
 
-      const typename DoFHandlerType::level_cell_iterator dof_cell(
+      const typename DoFHandler<dim, spacedim>::level_cell_iterator dof_cell(
         &fe_patch_values.get_cell()->get_triangulation(),
         fe_patch_values.get_cell()->level(),
         fe_patch_values.get_cell()->index(),
@@ -1400,18 +1342,17 @@ namespace internal
 
 
 
-    template <typename DoFHandlerType, typename VectorType>
+    template <int dim, int spacedim, typename VectorType>
     void
-    MGDataEntry<DoFHandlerType, VectorType>::get_function_values(
-      const FEValuesBase<DoFHandlerType::dimension,
-                         DoFHandlerType::space_dimension> &fe_patch_values,
-      const ComponentExtractor                             extract_component,
+    MGDataEntry<dim, spacedim, VectorType>::get_function_values(
+      const FEValuesBase<dim, spacedim> &  fe_patch_values,
+      const ComponentExtractor             extract_component,
       std::vector<dealii::Vector<double>> &patch_values_system) const
     {
       Assert(extract_component == ComponentExtractor::real_part,
              ExcNotImplemented());
 
-      typename DoFHandlerType::level_cell_iterator dof_cell(
+      typename DoFHandler<dim, spacedim>::level_cell_iterator dof_cell(
         &fe_patch_values.get_cell()->get_triangulation(),
         fe_patch_values.get_cell()->level(),
         fe_patch_values.get_cell()->index(),
@@ -1543,7 +1484,9 @@ DataOut_DoFData<DoFHandlerType, patch_dim, patch_space_dim>::add_data_vector(
 
 
   auto new_entry = std::make_unique<
-    internal::DataOutImplementation::DataEntry<DoFHandlerType, VectorType>>(
+    internal::DataOutImplementation::DataEntry<DoFHandlerType::dimension,
+                                               DoFHandlerType::space_dimension,
+                                               VectorType>>(
     &dof_handler, &vec, &data_postprocessor);
   dof_data.emplace_back(std::move(new_entry));
 }
@@ -1665,7 +1608,9 @@ DataOut_DoFData<DoFHandlerType, patch_dim, patch_space_dim>::
 
   // finally, add the data vector:
   auto new_entry = std::make_unique<
-    internal::DataOutImplementation::DataEntry<DoFHandlerType, VectorType>>(
+    internal::DataOutImplementation::DataEntry<DoFHandlerType::dimension,
+                                               DoFHandlerType::space_dimension,
+                                               VectorType>>(
     dof_handler, &data_vector, deduced_names, data_component_interpretation);
 
   if (actual_type == type_dof_data)
@@ -1737,9 +1682,14 @@ DataOut_DoFData<DoFHandlerType, patch_dim, patch_space_dim>::add_mg_data_vector(
          ExcMessage(
            "Invalid number of entries in data_component_interpretation."));
 
-  auto new_entry = std::make_unique<
-    internal::DataOutImplementation::MGDataEntry<DoFHandlerType, VectorType>>(
-    &dof_handler, &data, deduced_names, data_component_interpretation);
+  auto new_entry =
+    std::make_unique<internal::DataOutImplementation::MGDataEntry<
+      DoFHandlerType::dimension,
+      DoFHandlerType::space_dimension,
+      VectorType>>(&dof_handler,
+                   &data,
+                   deduced_names,
+                   data_component_interpretation);
   dof_data.emplace_back(std::move(new_entry));
 }
 
index 6e4b8d8307bd3189cdf3b597d305de720c0b8d4c..559cf27b97d9b05871ccc0537a43acd33c9c2037 100644 (file)
@@ -145,24 +145,3 @@ for (deal_II_dimension : DIMENSIONS)
       \}
     \}
   }
-
-
-for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS)
-  {
-    namespace internal
-    \{
-      namespace DataOutImplementation
-      \{
-        template void
-        ParallelDataBase<deal_II_dimension, deal_II_dimension>::
-          reinit_all_fe_values<
-            dealii::DH<deal_II_dimension, deal_II_dimension>>(
-            std::vector<std::shared_ptr<
-              DataEntryBase<dealii::DH<deal_II_dimension, deal_II_dimension>>>>
-              &dof_data,
-            const dealii::Triangulation<deal_II_dimension,
-                                        deal_II_dimension>::cell_iterator &cell,
-            const unsigned int face);
-      \}
-    \}
-  }
index 5915471d972c3f7bf50e79bd001f753c0652a55f..cf5e763c8ac6e4cbc32848dfa6a946797b9da5b8 100644 (file)
@@ -82,28 +82,3 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : DIMENSIONS)
       \}
     \}
   }
-
-
-for (DH : DOFHANDLER_TEMPLATES; deal_II_dimension : DIMENSIONS;
-     deal_II_space_dimension : DIMENSIONS)
-  {
-    namespace internal
-    \{
-      namespace DataOutImplementation
-      \{
-#if deal_II_dimension < deal_II_space_dimension
-        template void
-        ParallelDataBase<deal_II_dimension, deal_II_space_dimension>::
-          reinit_all_fe_values<
-            dealii::DH<deal_II_dimension, deal_II_space_dimension>>(
-            std::vector<std::shared_ptr<DataEntryBase<
-              dealii::DH<deal_II_dimension, deal_II_space_dimension>>>>
-              &dof_data,
-            const dealii::Triangulation<deal_II_dimension,
-                                        deal_II_space_dimension>::cell_iterator
-              &                cell,
-            const unsigned int face);
-#endif
-      \}
-    \}
-  }

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.