]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove DoFHandlerType from numerics/derivative_approximation 10598/head
authorPeter Munch <peterrmuench@gmail.com>
Thu, 25 Jun 2020 05:59:02 +0000 (07:59 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Thu, 25 Jun 2020 08:02:14 +0000 (10:02 +0200)
include/deal.II/numerics/derivative_approximation.h
source/numerics/derivative_approximation.cc
source/numerics/derivative_approximation.inst.in

index 732b6e52ff33a535e9c1094f6885d27c1dd5c45d..9474ba5d8b907c6bea63a52ceebd365d57dec6e8 100644 (file)
@@ -21,6 +21,8 @@
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/synchronous_iterator.h>
 
+#include <deal.II/dofs/dof_handler.h>
+
 #include <deal.II/fe/fe_update_flags.h>
 #include <deal.II/fe/mapping.h>
 
@@ -74,7 +76,7 @@ DEAL_II_NAMESPACE_OPEN
  * of file <source/numerics/derivative_approximation.cc> in function
  *     void DerivativeApproximation::approximate(...)
  * [with DerivativeDescription = DerivativeApproximation::Gradient<3>, int
- * dim = 3, DoFHandlerType = DoFHandler, InputVector = Vector<double>]
+ * dim = 3, InputVector = Vector<double>, spacedim = 3]
  * The violated condition was:
  *     determinant(Y) != 0
  * The name and call sequence of the exception was:
@@ -173,29 +175,23 @@ namespace DerivativeApproximation
    * In a parallel computation the @p solution vector needs to contain the
    * locally relevant unknowns.
    */
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_gradient(const Mapping<dim, spacedim> &       mapping,
-                       const DoFHandlerType<dim, spacedim> &dof,
-                       const InputVector &                  solution,
-                       Vector<float> &                      derivative_norm,
-                       const unsigned int                   component = 0);
+  approximate_gradient(const Mapping<dim, spacedim> &   mapping,
+                       const DoFHandler<dim, spacedim> &dof,
+                       const InputVector &              solution,
+                       Vector<float> &                  derivative_norm,
+                       const unsigned int               component = 0);
 
   /**
    * Call the function above with <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
    */
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_gradient(const DoFHandlerType<dim, spacedim> &dof,
-                       const InputVector &                  solution,
-                       Vector<float> &                      derivative_norm,
-                       const unsigned int                   component = 0);
+  approximate_gradient(const DoFHandler<dim, spacedim> &dof,
+                       const InputVector &              solution,
+                       Vector<float> &                  derivative_norm,
+                       const unsigned int               component = 0);
 
   /**
    * This function is the analogue to the one above, computing finite
@@ -214,27 +210,21 @@ namespace DerivativeApproximation
    * In a parallel computation the @p solution vector needs to contain the
    * locally relevant unknowns.
    */
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_second_derivative(const Mapping<dim, spacedim> &       mapping,
-                                const DoFHandlerType<dim, spacedim> &dof,
-                                const InputVector &                  solution,
+  approximate_second_derivative(const Mapping<dim, spacedim> &   mapping,
+                                const DoFHandler<dim, spacedim> &dof,
+                                const InputVector &              solution,
                                 Vector<float> &    derivative_norm,
                                 const unsigned int component = 0);
 
   /**
    * Call the function above with <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
    */
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_second_derivative(const DoFHandlerType<dim, spacedim> &dof,
-                                const InputVector &                  solution,
+  approximate_second_derivative(const DoFHandler<dim, spacedim> &dof,
+                                const InputVector &              solution,
                                 Vector<float> &    derivative_norm,
                                 const unsigned int component = 0);
 
@@ -251,42 +241,37 @@ namespace DerivativeApproximation
    * In a parallel computation the @p solution vector needs to contain the
    * locally relevant unknowns.
    */
-  template <typename DoFHandlerType, class InputVector, int order>
+  template <int dim, int spacedim, class InputVector, int order>
   void
   approximate_derivative_tensor(
-    const Mapping<DoFHandlerType::dimension, DoFHandlerType::space_dimension>
-      &                   mapping,
-    const DoFHandlerType &dof,
-    const InputVector &   solution,
+    const Mapping<dim, spacedim> &   mapping,
+    const DoFHandler<dim, spacedim> &dof,
+    const InputVector &              solution,
 #ifndef _MSC_VER
-    const typename DoFHandlerType::active_cell_iterator &cell,
+    const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
 #else
-    const TriaActiveIterator<
-      dealii::DoFCellAccessor<DoFHandlerType::dimension,
-                              DoFHandlerType::space_dimension,
-                              false>> &cell,
+    const TriaActiveIterator<dealii::DoFCellAccessor<dim, spacedim, false>>
+      &cell,
 #endif
-    Tensor<order, DoFHandlerType::dimension> &derivative,
-    const unsigned int                        component = 0);
+    Tensor<order, dim> &derivative,
+    const unsigned int  component = 0);
 
   /**
    * Same as above, with <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
    */
-  template <typename DoFHandlerType, class InputVector, int order>
+  template <int dim, int spacedim, class InputVector, int order>
   void
   approximate_derivative_tensor(
-    const DoFHandlerType &dof,
-    const InputVector &   solution,
+    const DoFHandler<dim, spacedim> &dof,
+    const InputVector &              solution,
 #ifndef _MSC_VER
-    const typename DoFHandlerType::active_cell_iterator &cell,
+    const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
 #else
-    const TriaActiveIterator<
-      dealii::DoFCellAccessor<DoFHandlerType::dimension,
-                              DoFHandlerType::space_dimension,
-                              false>> &cell,
+    const TriaActiveIterator<dealii::DoFCellAccessor<dim, spacedim, false>>
+      &cell,
 #endif
-    Tensor<order, DoFHandlerType::dimension> &derivative,
-    const unsigned int                        component = 0);
+    Tensor<order, dim> &derivative,
+    const unsigned int  component = 0);
 
   /**
    * Return the norm of the derivative.
index 3c6aa413c971a32d1593ca74a7a502580d16b9ce..38ad27d01b98ddee235c4f334398acf50fcc16a3 100644 (file)
@@ -736,16 +736,15 @@ namespace DerivativeApproximation
      */
     template <class DerivativeDescription,
               int dim,
-              template <int, int> class DoFHandlerType,
               class InputVector,
               int spacedim>
     void
     approximate_cell(
-      const Mapping<dim, spacedim> &       mapping,
-      const DoFHandlerType<dim, spacedim> &dof_handler,
-      const InputVector &                  solution,
-      const unsigned int                   component,
-      const typename DoFHandlerType<dim, spacedim>::active_cell_iterator &cell,
+      const Mapping<dim, spacedim> &   mapping,
+      const DoFHandler<dim, spacedim> &dof_handler,
+      const InputVector &              solution,
+      const unsigned int               component,
+      const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
       typename DerivativeDescription::Derivative &derivative)
     {
       QMidpoint<dim> midpoint_rule;
@@ -776,7 +775,7 @@ namespace DerivativeApproximation
       // active neighbors of a cell
       // reserve the maximal number of
       // active neighbors
-      std::vector<typename DoFHandlerType<dim, spacedim>::active_cell_iterator>
+      std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>
         active_neighbors;
 
       active_neighbors.reserve(GeometryInfo<dim>::faces_per_cell *
@@ -815,7 +814,7 @@ namespace DerivativeApproximation
       // first collect all neighbor
       // cells in a vector, and then
       // collect the data from them
-      GridTools::get_active_neighbors<DoFHandlerType<dim, spacedim>>(
+      GridTools::get_active_neighbors<DoFHandler<dim, spacedim>>(
         cell, active_neighbors);
 
       // now loop over all active
@@ -907,16 +906,15 @@ namespace DerivativeApproximation
      */
     template <class DerivativeDescription,
               int dim,
-              template <int, int> class DoFHandlerType,
               class InputVector,
               int spacedim>
     void
     approximate(
       SynchronousIterators<
-        std::tuple<typename DoFHandlerType<dim, spacedim>::active_cell_iterator,
+        std::tuple<typename DoFHandler<dim, spacedim>::active_cell_iterator,
                    Vector<float>::iterator>> const &cell,
       const Mapping<dim, spacedim> &                mapping,
-      const DoFHandlerType<dim, spacedim> &         dof_handler,
+      const DoFHandler<dim, spacedim> &             dof_handler,
       const InputVector &                           solution,
       const unsigned int                            component)
     {
@@ -928,16 +926,13 @@ namespace DerivativeApproximation
           typename DerivativeDescription::Derivative derivative;
           // call the function doing the actual
           // work on this cell
-          approximate_cell<DerivativeDescription,
-                           dim,
-                           DoFHandlerType,
-                           InputVector,
-                           spacedim>(mapping,
-                                     dof_handler,
-                                     solution,
-                                     component,
-                                     std::get<0>(*cell),
-                                     derivative);
+          approximate_cell<DerivativeDescription, dim, InputVector, spacedim>(
+            mapping,
+            dof_handler,
+            solution,
+            component,
+            std::get<0>(*cell),
+            derivative);
 
           // evaluate the norm and fill the vector
           //*derivative_norm_on_this_cell
@@ -959,15 +954,14 @@ namespace DerivativeApproximation
      */
     template <class DerivativeDescription,
               int dim,
-              template <int, int> class DoFHandlerType,
               class InputVector,
               int spacedim>
     void
-    approximate_derivative(const Mapping<dim, spacedim> &       mapping,
-                           const DoFHandlerType<dim, spacedim> &dof_handler,
-                           const InputVector &                  solution,
-                           const unsigned int                   component,
-                           Vector<float> &                      derivative_norm)
+    approximate_derivative(const Mapping<dim, spacedim> &   mapping,
+                           const DoFHandler<dim, spacedim> &dof_handler,
+                           const InputVector &              solution,
+                           const unsigned int               component,
+                           Vector<float> &                  derivative_norm)
     {
       Assert(derivative_norm.size() ==
                dof_handler.get_triangulation().n_active_cells(),
@@ -977,7 +971,7 @@ namespace DerivativeApproximation
       AssertIndexRange(component, dof_handler.get_fe(0).n_components());
 
       using Iterators =
-        std::tuple<typename DoFHandlerType<dim, spacedim>::active_cell_iterator,
+        std::tuple<typename DoFHandler<dim, spacedim>::active_cell_iterator,
                    Vector<float>::iterator>;
       SynchronousIterators<Iterators> begin(
         Iterators(dof_handler.begin_active(), derivative_norm.begin())),
@@ -993,11 +987,7 @@ namespace DerivativeApproximation
           SynchronousIterators<Iterators> const &cell,
           Assembler::Scratch const &,
           Assembler::CopyData &) {
-          approximate<DerivativeDescription,
-                      dim,
-                      DoFHandlerType,
-                      InputVector,
-                      spacedim>(
+          approximate<DerivativeDescription, dim, InputVector, spacedim>(
             cell, mapping, dof_handler, solution, component);
         },
         std::function<void(internal::Assembler::CopyData const &)>(),
@@ -1014,31 +1004,25 @@ namespace DerivativeApproximation
 
 namespace DerivativeApproximation
 {
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_gradient(const Mapping<dim, spacedim> &       mapping,
-                       const DoFHandlerType<dim, spacedim> &dof_handler,
-                       const InputVector &                  solution,
-                       Vector<float> &                      derivative_norm,
-                       const unsigned int                   component)
+  approximate_gradient(const Mapping<dim, spacedim> &   mapping,
+                       const DoFHandler<dim, spacedim> &dof_handler,
+                       const InputVector &              solution,
+                       Vector<float> &                  derivative_norm,
+                       const unsigned int               component)
   {
     internal::approximate_derivative<internal::Gradient<dim>, dim>(
       mapping, dof_handler, solution, component, derivative_norm);
   }
 
 
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_gradient(const DoFHandlerType<dim, spacedim> &dof_handler,
-                       const InputVector &                  solution,
-                       Vector<float> &                      derivative_norm,
-                       const unsigned int                   component)
+  approximate_gradient(const DoFHandler<dim, spacedim> &dof_handler,
+                       const InputVector &              solution,
+                       Vector<float> &                  derivative_norm,
+                       const unsigned int               component)
   {
     internal::approximate_derivative<internal::Gradient<dim>, dim>(
       StaticMappingQ1<dim>::mapping,
@@ -1049,33 +1033,25 @@ namespace DerivativeApproximation
   }
 
 
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_second_derivative(
-    const Mapping<dim, spacedim> &       mapping,
-    const DoFHandlerType<dim, spacedim> &dof_handler,
-    const InputVector &                  solution,
-    Vector<float> &                      derivative_norm,
-    const unsigned int                   component)
+  approximate_second_derivative(const Mapping<dim, spacedim> &   mapping,
+                                const DoFHandler<dim, spacedim> &dof_handler,
+                                const InputVector &              solution,
+                                Vector<float> &    derivative_norm,
+                                const unsigned int component)
   {
     internal::approximate_derivative<internal::SecondDerivative<dim>, dim>(
       mapping, dof_handler, solution, component, derivative_norm);
   }
 
 
-  template <int dim,
-            template <int, int> class DoFHandlerType,
-            class InputVector,
-            int spacedim>
+  template <int dim, class InputVector, int spacedim>
   void
-  approximate_second_derivative(
-    const DoFHandlerType<dim, spacedim> &dof_handler,
-    const InputVector &                  solution,
-    Vector<float> &                      derivative_norm,
-    const unsigned int                   component)
+  approximate_second_derivative(const DoFHandler<dim, spacedim> &dof_handler,
+                                const InputVector &              solution,
+                                Vector<float> &    derivative_norm,
+                                const unsigned int component)
   {
     internal::approximate_derivative<internal::SecondDerivative<dim>, dim>(
       StaticMappingQ1<dim>::mapping,
@@ -1086,56 +1062,49 @@ namespace DerivativeApproximation
   }
 
 
-  template <typename DoFHandlerType, class InputVector, int order>
+  template <int dim, int spacedim, class InputVector, int order>
   void
   approximate_derivative_tensor(
-    const Mapping<DoFHandlerType::dimension, DoFHandlerType::space_dimension>
-      &                   mapping,
-    const DoFHandlerType &dof,
-    const InputVector &   solution,
+    const Mapping<dim, spacedim> &   mapping,
+    const DoFHandler<dim, spacedim> &dof,
+    const InputVector &              solution,
 #ifndef _MSC_VER
-    const typename DoFHandlerType::active_cell_iterator &cell,
+    const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
 #else
-    const TriaActiveIterator<
-      dealii::DoFCellAccessor<DoFHandlerType::dimension,
-                              DoFHandlerType::space_dimension,
-                              false>> &cell,
+    const TriaActiveIterator<dealii::DoFCellAccessor<dim, spacedim, false>>
+      &cell,
 #endif
-    Tensor<order, DoFHandlerType::dimension> &derivative,
-    const unsigned int                        component)
+    Tensor<order, dim> &derivative,
+    const unsigned int  component)
   {
     internal::approximate_cell<
-      typename internal::DerivativeSelector<order, DoFHandlerType::dimension>::
-        DerivDescr>(mapping, dof, solution, component, cell, derivative);
+      typename internal::DerivativeSelector<order, dim>::DerivDescr>(
+      mapping, dof, solution, component, cell, derivative);
   }
 
 
 
-  template <typename DoFHandlerType, class InputVector, int order>
+  template <int dim, int spacedim, class InputVector, int order>
   void
   approximate_derivative_tensor(
-    const DoFHandlerType &dof,
-    const InputVector &   solution,
+    const DoFHandler<dim, spacedim> &dof,
+    const InputVector &              solution,
 #ifndef _MSC_VER
-    const typename DoFHandlerType::active_cell_iterator &cell,
+    const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
 #else
-    const TriaActiveIterator<
-      dealii::DoFCellAccessor<DoFHandlerType::dimension,
-                              DoFHandlerType::space_dimension,
-                              false>> &cell,
+    const TriaActiveIterator<dealii::DoFCellAccessor<dim, spacedim, false>>
+      &cell,
 #endif
-    Tensor<order, DoFHandlerType::dimension> &derivative,
-    const unsigned int                        component)
+    Tensor<order, dim> &derivative,
+    const unsigned int  component)
   {
     // just call the respective function with Q1 mapping
-    approximate_derivative_tensor(
-      StaticMappingQ1<DoFHandlerType::dimension,
-                      DoFHandlerType::space_dimension>::mapping,
-      dof,
-      solution,
-      cell,
-      derivative,
-      component);
+    approximate_derivative_tensor(StaticMappingQ1<dim, spacedim>::mapping,
+                                  dof,
+                                  solution,
+                                  cell,
+                                  derivative,
+                                  component);
   }
 
 
index e34c0b85ff982be9634e224249ee228dd0895942..0eff551f4bd48fa2c1d95f10465d0ab2d23891f5 100644 (file)
 // ---------------------------------------------------------------------
 
 
-for (deal_II_dimension : DIMENSIONS; VEC : REAL_VECTOR_TYPES;
-     DH : DOFHANDLER_TEMPLATES)
+for (deal_II_dimension : DIMENSIONS; VEC : REAL_VECTOR_TYPES)
   {
     namespace DerivativeApproximation
     \{
       template void
       approximate_gradient<deal_II_dimension>(
-        const Mapping<deal_II_dimension> &mapping,
-        const DH<deal_II_dimension> &     dof_handler,
-        const VEC &                       solution,
-        Vector<float> &                   derivative_norm,
-        const unsigned int                component);
+        const Mapping<deal_II_dimension> &   mapping,
+        const DoFHandler<deal_II_dimension> &dof_handler,
+        const VEC &                          solution,
+        Vector<float> &                      derivative_norm,
+        const unsigned int                   component);
 
       template void
       approximate_gradient<deal_II_dimension>(
-        const DH<deal_II_dimension> &dof_handler,
-        const VEC &                  solution,
-        Vector<float> &              derivative_norm,
-        const unsigned int           component);
+        const DoFHandler<deal_II_dimension> &dof_handler,
+        const VEC &                          solution,
+        Vector<float> &                      derivative_norm,
+        const unsigned int                   component);
 
       template void
       approximate_second_derivative<deal_II_dimension>(
-        const Mapping<deal_II_dimension> &mapping,
-        const DH<deal_II_dimension> &     dof_handler,
-        const VEC &                       solution,
-        Vector<float> &                   derivative_norm,
-        const unsigned int                component);
+        const Mapping<deal_II_dimension> &   mapping,
+        const DoFHandler<deal_II_dimension> &dof_handler,
+        const VEC &                          solution,
+        Vector<float> &                      derivative_norm,
+        const unsigned int                   component);
 
       template void
       approximate_second_derivative<deal_II_dimension>(
-        const DH<deal_II_dimension> &dof_handler,
-        const VEC &                  solution,
-        Vector<float> &              derivative_norm,
-        const unsigned int           component);
+        const DoFHandler<deal_II_dimension> &dof_handler,
+        const VEC &                          solution,
+        Vector<float> &                      derivative_norm,
+        const unsigned int                   component);
 
       template void
-      approximate_derivative_tensor<DH<deal_II_dimension>, VEC, 1>(
-        const Mapping<deal_II_dimension> &                 mapping,
-        const DH<deal_II_dimension> &                      dof_handler,
-        const VEC &                                        solution,
-        const DH<deal_II_dimension>::active_cell_iterator &cell,
-        Tensor<1, deal_II_dimension> &                     derivative,
-        const unsigned int                                 component);
+      approximate_derivative_tensor<deal_II_dimension,
+                                    deal_II_dimension,
+                                    VEC,
+                                    1>(
+        const Mapping<deal_II_dimension> &                         mapping,
+        const DoFHandler<deal_II_dimension> &                      dof_handler,
+        const VEC &                                                solution,
+        const DoFHandler<deal_II_dimension>::active_cell_iterator &cell,
+        Tensor<1, deal_II_dimension> &                             derivative,
+        const unsigned int                                         component);
 
       template void
-      approximate_derivative_tensor<DH<deal_II_dimension>, VEC, 2>(
-        const Mapping<deal_II_dimension> &                 mapping,
-        const DH<deal_II_dimension> &                      dof_handler,
-        const VEC &                                        solution,
-        const DH<deal_II_dimension>::active_cell_iterator &cell,
-        Tensor<2, deal_II_dimension> &                     derivative,
-        const unsigned int                                 component);
+      approximate_derivative_tensor<deal_II_dimension,
+                                    deal_II_dimension,
+                                    VEC,
+                                    2>(
+        const Mapping<deal_II_dimension> &                         mapping,
+        const DoFHandler<deal_II_dimension> &                      dof_handler,
+        const VEC &                                                solution,
+        const DoFHandler<deal_II_dimension>::active_cell_iterator &cell,
+        Tensor<2, deal_II_dimension> &                             derivative,
+        const unsigned int                                         component);
 
       template void
-      approximate_derivative_tensor<DH<deal_II_dimension>, VEC, 3>(
-        const Mapping<deal_II_dimension> &                 mapping,
-        const DH<deal_II_dimension> &                      dof_handler,
-        const VEC &                                        solution,
-        const DH<deal_II_dimension>::active_cell_iterator &cell,
-        Tensor<3, deal_II_dimension> &                     derivative,
-        const unsigned int                                 component);
+      approximate_derivative_tensor<deal_II_dimension,
+                                    deal_II_dimension,
+                                    VEC,
+                                    3>(
+        const Mapping<deal_II_dimension> &                         mapping,
+        const DoFHandler<deal_II_dimension> &                      dof_handler,
+        const VEC &                                                solution,
+        const DoFHandler<deal_II_dimension>::active_cell_iterator &cell,
+        Tensor<3, deal_II_dimension> &                             derivative,
+        const unsigned int                                         component);
 
 
       template void
-      approximate_derivative_tensor<DH<deal_II_dimension>, VEC, 1>(
-        const DH<deal_II_dimension> &                      dof_handler,
-        const VEC &                                        solution,
-        const DH<deal_II_dimension>::active_cell_iterator &cell,
-        Tensor<1, deal_II_dimension> &                     derivative,
-        const unsigned int                                 component);
+      approximate_derivative_tensor<deal_II_dimension,
+                                    deal_II_dimension,
+                                    VEC,
+                                    1>(
+        const DoFHandler<deal_II_dimension> &                      dof_handler,
+        const VEC &                                                solution,
+        const DoFHandler<deal_II_dimension>::active_cell_iterator &cell,
+        Tensor<1, deal_II_dimension> &                             derivative,
+        const unsigned int                                         component);
 
       template void
-      approximate_derivative_tensor<DH<deal_II_dimension>, VEC, 2>(
-        const DH<deal_II_dimension> &                      dof_handler,
-        const VEC &                                        solution,
-        const DH<deal_II_dimension>::active_cell_iterator &cell,
-        Tensor<2, deal_II_dimension> &                     derivative,
-        const unsigned int                                 component);
+      approximate_derivative_tensor<deal_II_dimension,
+                                    deal_II_dimension,
+                                    VEC,
+                                    2>(
+        const DoFHandler<deal_II_dimension> &                      dof_handler,
+        const VEC &                                                solution,
+        const DoFHandler<deal_II_dimension>::active_cell_iterator &cell,
+        Tensor<2, deal_II_dimension> &                             derivative,
+        const unsigned int                                         component);
 
       template void
-      approximate_derivative_tensor<DH<deal_II_dimension>, VEC, 3>(
-        const DH<deal_II_dimension> &                      dof_handler,
-        const VEC &                                        solution,
-        const DH<deal_II_dimension>::active_cell_iterator &cell,
-        Tensor<3, deal_II_dimension> &                     derivative,
-        const unsigned int                                 component);
+      approximate_derivative_tensor<deal_II_dimension,
+                                    deal_II_dimension,
+                                    VEC,
+                                    3>(
+        const DoFHandler<deal_II_dimension> &                      dof_handler,
+        const VEC &                                                solution,
+        const DoFHandler<deal_II_dimension>::active_cell_iterator &cell,
+        Tensor<3, deal_II_dimension> &                             derivative,
+        const unsigned int                                         component);
 
     \}
   }

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.