]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add std::complex support to integrate_difference()
authorDaniel Garcia-Sanchez <daniel.garcia-sanchez@insp.upmc.fr>
Wed, 9 Oct 2019 07:43:41 +0000 (09:43 +0200)
committerDaniel Garcia-Sanchez <daniel.garcia-sanchez@insp.upmc.fr>
Wed, 9 Oct 2019 07:43:41 +0000 (09:43 +0200)
The type of exact_solution was double and fe_function could have any numerical
type (double, float, std::complex<double> or std::complex<float>). This would
lead in certain cases to unnecessary casts, for example from double to float.
In addition it was not possible to compare a std::complex exact_solution to a
std::complex fe_function.

Now the types of exact_solution and fe_function must be the same, as a result
integrate_difference() can be used to compare a std::complex exact_solution to
a std::complex fe_function.

This change is not backward-compatible. The old version of the function has
been deprecated.

include/deal.II/numerics/vector_tools.h
include/deal.II/numerics/vector_tools.templates.h
source/numerics/vector_tools_integrate_difference.inst.in

index 986e3125eff45acbd23cca32a8b6339d16f704a6..52ae5023a1a119a80f8488bed8ecb40db8d645b5 100644 (file)
@@ -2516,6 +2516,82 @@ namespace VectorTools
    */
   template <int dim, class InVector, class OutVector, int spacedim>
   void
+  integrate_difference(
+    const Mapping<dim, spacedim> &                           mapping,
+    const DoFHandler<dim, spacedim> &                        dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const Quadrature<dim> &                                  q,
+    const NormType &                                         norm,
+    const Function<spacedim, double> *                       weight   = nullptr,
+    const double                                             exponent = 2.);
+
+  /**
+   * Call the integrate_difference() function, see above, with
+   * <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
+   */
+  template <int dim, class InVector, class OutVector, int spacedim>
+  void
+  integrate_difference(
+    const DoFHandler<dim, spacedim> &                        dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const Quadrature<dim> &                                  q,
+    const NormType &                                         norm,
+    const Function<spacedim, double> *                       weight   = nullptr,
+    const double                                             exponent = 2.);
+
+  /**
+   * Same as above for hp.
+   */
+  template <int dim, class InVector, class OutVector, int spacedim>
+  void
+  integrate_difference(
+    const hp::MappingCollection<dim, spacedim> &             mapping,
+    const hp::DoFHandler<dim, spacedim> &                    dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const hp::QCollection<dim> &                             q,
+    const NormType &                                         norm,
+    const Function<spacedim, double> *                       weight   = nullptr,
+    const double                                             exponent = 2.);
+
+  /**
+   * Call the integrate_difference() function, see above, with
+   * <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
+   */
+  template <int dim, class InVector, class OutVector, int spacedim>
+  void
+  integrate_difference(
+    const hp::DoFHandler<dim, spacedim> &                    dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const hp::QCollection<dim> &                             q,
+    const NormType &                                         norm,
+    const Function<spacedim, double> *                       weight   = nullptr,
+    const double                                             exponent = 2.);
+
+  /**
+   * Compute the cellwise error of the finite element solution.  Integrate the
+   * difference between a reference function which is given as a continuous
+   * function object, and a finite element function. The result of this
+   * function is the vector @p difference that contains one value per active
+   * cell $K$ of the triangulation. Each of the values of this vector $d$
+   * equals
+   * @f{align*}{
+   * d_K = \| u-u_h \|_X
+   * @f}
+   * where $X$ denotes the norm chosen and $u$ represents the exact solution.
+   *
+   * @deprecated Use integrate_difference(const Mapping<dim, spacedim> &, const DoFHandler<dim, spacedim> &, const InVector &, const Function<spacedim, typename InVector::value_type> &, OutVector &, const Quadrature<dim> &, const NormType &, const Function<spacedim, double> *, const double) instead.
+   */
+  template <int dim, class InVector, class OutVector, int spacedim>
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(const Mapping<dim, spacedim> &    mapping,
                        const DoFHandler<dim, spacedim> & dof,
                        const InVector &                  fe_function,
@@ -2529,9 +2605,12 @@ namespace VectorTools
   /**
    * Call the integrate_difference() function, see above, with
    * <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
+   *
+   * @deprecated Use integrate_difference(const DoFHandler<dim, spacedim> &, const InVector &, const Function<spacedim, typename InVector::value_type> &exact_solution, OutVector &, const Quadrature<dim> &, const NormType &, const Function<spacedim, double> *, const double) instead.
    */
   template <int dim, class InVector, class OutVector, int spacedim>
-  void
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(const DoFHandler<dim, spacedim> & dof,
                        const InVector &                  fe_function,
                        const Function<spacedim, double> &exact_solution,
@@ -2543,9 +2622,12 @@ namespace VectorTools
 
   /**
    * Same as above for hp.
+   *
+   * @deprecated Use integrate_difference(const hp::MappingCollection<dim, spacedim> &, const hp::DoFHandler<dim, spacedim> &, const InVector &, const Function<spacedim, typename InVector::value_type> &, OutVector &, const hp::QCollection<dim> &, const NormType &, const Function<spacedim, double> *, const double) instead.
    */
   template <int dim, class InVector, class OutVector, int spacedim>
-  void
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(const hp::MappingCollection<dim, spacedim> &mapping,
                        const hp::DoFHandler<dim, spacedim> &       dof,
                        const InVector &                            fe_function,
@@ -2559,9 +2641,12 @@ namespace VectorTools
   /**
    * Call the integrate_difference() function, see above, with
    * <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
+   *
+   * @deprecated Use integrate_difference(const hp::DoFHandler<dim, spacedim> &, const InVector &, const Function<spacedim, typename InVector::value_type> &, OutVector &, const hp::QCollection<dim> &, const NormType &, const Function<spacedim, double> *, const double) instead.
    */
   template <int dim, class InVector, class OutVector, int spacedim>
-  void
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(const hp::DoFHandler<dim, spacedim> &dof,
                        const InVector &                     fe_function,
                        const Function<spacedim, double> &   exact_solution,
index 61a9b3588e8960216e3eecfbb639be5933cc06f0..35b085efd897ab61ed5ffdf5582fc72a3536461f 100644 (file)
@@ -7819,10 +7819,11 @@ namespace VectorTools
       std::vector<std::vector<Tensor<1, spacedim, Number>>> psi_grads;
       std::vector<Number>                                   psi_scalar;
 
-      std::vector<double>                           tmp_values;
-      std::vector<Vector<double>>                   tmp_vector_values;
-      std::vector<Tensor<1, spacedim>>              tmp_gradients;
-      std::vector<std::vector<Tensor<1, spacedim>>> tmp_vector_gradients;
+      std::vector<Number>                      tmp_values;
+      std::vector<Vector<Number>>              tmp_vector_values;
+      std::vector<Tensor<1, spacedim, Number>> tmp_gradients;
+      std::vector<std::vector<Tensor<1, spacedim, Number>>>
+        tmp_vector_gradients;
 
       dealii::hp::FEValues<dim, spacedim> x_fe_values;
     };
@@ -7859,55 +7860,420 @@ namespace VectorTools
       weight_values.resize(n_q_points);
       weight_vectors.resize(n_q_points, Vector<double>(n_components));
 
-      psi_values.resize(n_q_points, Vector<Number>(n_components));
-      psi_grads.resize(n_q_points,
-                       std::vector<Tensor<1, spacedim, Number>>(n_components));
-      psi_scalar.resize(n_q_points);
+      psi_values.resize(n_q_points, Vector<Number>(n_components));
+      psi_grads.resize(n_q_points,
+                       std::vector<Tensor<1, spacedim, Number>>(n_components));
+      psi_scalar.resize(n_q_points);
+
+      tmp_values.resize(n_q_points);
+      tmp_vector_values.resize(n_q_points, Vector<Number>(n_components));
+      tmp_gradients.resize(n_q_points);
+      tmp_vector_gradients.resize(
+        n_q_points, std::vector<Tensor<1, spacedim, Number>>(n_components));
+    }
+
+    template <int dim, int spacedim, typename Number>
+    struct DEAL_II_DEPRECATED DeprecatedIDScratchData
+    {
+      DeprecatedIDScratchData(
+        const dealii::hp::MappingCollection<dim, spacedim> &mapping,
+        const dealii::hp::FECollection<dim, spacedim> &     fe,
+        const dealii::hp::QCollection<dim> &                q,
+        const UpdateFlags                                   update_flags);
+
+      DeprecatedIDScratchData(const DeprecatedIDScratchData &data);
+
+      void
+      resize_vectors(const unsigned int n_q_points,
+                     const unsigned int n_components);
+
+      std::vector<Vector<Number>>                           function_values;
+      std::vector<std::vector<Tensor<1, spacedim, Number>>> function_grads;
+      std::vector<double>                                   weight_values;
+      std::vector<Vector<double>>                           weight_vectors;
+
+      std::vector<Vector<Number>>                           psi_values;
+      std::vector<std::vector<Tensor<1, spacedim, Number>>> psi_grads;
+      std::vector<Number>                                   psi_scalar;
+
+      std::vector<double>                           tmp_values;
+      std::vector<Vector<double>>                   tmp_vector_values;
+      std::vector<Tensor<1, spacedim>>              tmp_gradients;
+      std::vector<std::vector<Tensor<1, spacedim>>> tmp_vector_gradients;
+
+      dealii::hp::FEValues<dim, spacedim> x_fe_values;
+    };
+
+
+    template <int dim, int spacedim, typename Number>
+    DeprecatedIDScratchData<dim, spacedim, Number>::DeprecatedIDScratchData(
+      const dealii::hp::MappingCollection<dim, spacedim> &mapping,
+      const dealii::hp::FECollection<dim, spacedim> &     fe,
+      const dealii::hp::QCollection<dim> &                q,
+      const UpdateFlags                                   update_flags)
+      : x_fe_values(mapping, fe, q, update_flags)
+    {}
+
+    template <int dim, int spacedim, typename Number>
+    DeprecatedIDScratchData<dim, spacedim, Number>::DeprecatedIDScratchData(
+      const DeprecatedIDScratchData &data)
+      : x_fe_values(data.x_fe_values.get_mapping_collection(),
+                    data.x_fe_values.get_fe_collection(),
+                    data.x_fe_values.get_quadrature_collection(),
+                    data.x_fe_values.get_update_flags())
+    {}
+
+    template <int dim, int spacedim, typename Number>
+    void
+    DeprecatedIDScratchData<dim, spacedim, Number>::resize_vectors(
+      const unsigned int n_q_points,
+      const unsigned int n_components)
+    {
+      function_values.resize(n_q_points, Vector<Number>(n_components));
+      function_grads.resize(
+        n_q_points, std::vector<Tensor<1, spacedim, Number>>(n_components));
+
+      weight_values.resize(n_q_points);
+      weight_vectors.resize(n_q_points, Vector<double>(n_components));
+
+      psi_values.resize(n_q_points, Vector<Number>(n_components));
+      psi_grads.resize(n_q_points,
+                       std::vector<Tensor<1, spacedim, Number>>(n_components));
+      psi_scalar.resize(n_q_points);
+
+      tmp_values.resize(n_q_points);
+      tmp_vector_values.resize(n_q_points, Vector<double>(n_components));
+      tmp_gradients.resize(n_q_points);
+      tmp_vector_gradients.resize(
+        n_q_points, std::vector<Tensor<1, spacedim>>(n_components));
+    }
+
+    namespace internal
+    {
+      template <typename number>
+      double
+      mean_to_double(const number &mean_value)
+      {
+        return mean_value;
+      }
+
+      template <typename number>
+      double
+      mean_to_double(const std::complex<number> &mean_value)
+      {
+        // we need to return double as a norm, but mean value is a complex
+        // number. Panic and return real-part while warning the user that
+        // he shall never do that.
+        Assert(
+          false,
+          ExcMessage(
+            "Mean value norm is not implemented for complex-valued vectors"));
+        return mean_value.real();
+      }
+    } // namespace internal
+
+
+    // avoid compiling inner function for many vector types when we always
+    // really do the same thing by putting the main work into this helper
+    // function
+    template <int dim, int spacedim, typename Number>
+    double
+    integrate_difference_inner(const Function<spacedim, Number> &exact_solution,
+                               const NormType &                  norm,
+                               const Function<spacedim> *        weight,
+                               const UpdateFlags                 update_flags,
+                               const double                      exponent,
+                               const unsigned int                n_components,
+                               IDScratchData<dim, spacedim, Number> &data)
+    {
+      const bool                             fe_is_system = (n_components != 1);
+      const dealii::FEValues<dim, spacedim> &fe_values =
+        data.x_fe_values.get_present_fe_values();
+      const unsigned int n_q_points = fe_values.n_quadrature_points;
+
+      if (weight != nullptr)
+        {
+          if (weight->n_components > 1)
+            weight->vector_value_list(fe_values.get_quadrature_points(),
+                                      data.weight_vectors);
+          else
+            {
+              weight->value_list(fe_values.get_quadrature_points(),
+                                 data.weight_values);
+              for (unsigned int k = 0; k < n_q_points; ++k)
+                data.weight_vectors[k] = data.weight_values[k];
+            }
+        }
+      else
+        {
+          for (unsigned int k = 0; k < n_q_points; ++k)
+            data.weight_vectors[k] = 1.;
+        }
+
+
+      if (update_flags & update_values)
+        {
+          // first compute the exact solution (vectors) at the quadrature
+          // points. try to do this as efficient as possible by avoiding a
+          // second virtual function call in case the function really has only
+          // one component
+          //
+          // TODO: we have to work a bit here because the Function<dim,double>
+          //   interface of the argument denoting the exact function only
+          //   provides us with double/Tensor<1,dim> values, rather than
+          //   with the correct data type. so evaluate into a temp
+          //   object, then copy around
+          if (fe_is_system)
+            {
+              exact_solution.vector_value_list(
+                fe_values.get_quadrature_points(), data.tmp_vector_values);
+              for (unsigned int i = 0; i < n_q_points; ++i)
+                data.psi_values[i] = data.tmp_vector_values[i];
+            }
+          else
+            {
+              exact_solution.value_list(fe_values.get_quadrature_points(),
+                                        data.tmp_values);
+              for (unsigned int i = 0; i < n_q_points; ++i)
+                data.psi_values[i](0) = data.tmp_values[i];
+            }
+
+          // then subtract finite element fe_function
+          for (unsigned int q = 0; q < n_q_points; ++q)
+            for (unsigned int i = 0; i < data.psi_values[q].size(); ++i)
+              data.psi_values[q][i] -= data.function_values[q][i];
+        }
+
+      // Do the same for gradients, if required
+      if (update_flags & update_gradients)
+        {
+          // try to be a little clever to avoid recursive virtual function
+          // calls when calling gradient_list for functions that are really
+          // scalar functions
+          if (fe_is_system)
+            {
+              exact_solution.vector_gradient_list(
+                fe_values.get_quadrature_points(), data.tmp_vector_gradients);
+              for (unsigned int i = 0; i < n_q_points; ++i)
+                for (unsigned int comp = 0; comp < data.psi_grads[i].size();
+                     ++comp)
+                  data.psi_grads[i][comp] = data.tmp_vector_gradients[i][comp];
+            }
+          else
+            {
+              exact_solution.gradient_list(fe_values.get_quadrature_points(),
+                                           data.tmp_gradients);
+              for (unsigned int i = 0; i < n_q_points; ++i)
+                data.psi_grads[i][0] = data.tmp_gradients[i];
+            }
+
+          // then subtract finite element function_grads. We need to be
+          // careful in the codimension one case, since there we only have
+          // tangential gradients in the finite element function, not the full
+          // gradient. This is taken care of, by subtracting the normal
+          // component of the gradient from the exact function.
+          if (update_flags & update_normal_vectors)
+            for (unsigned int k = 0; k < n_components; ++k)
+              for (unsigned int q = 0; q < n_q_points; ++q)
+                {
+                  // compute (f.n) n
+                  const typename ProductType<Number, double>::type f_dot_n =
+                    data.psi_grads[q][k] * fe_values.normal_vector(q);
+                  const Tensor<1, spacedim, Number> f_dot_n_times_n =
+                    f_dot_n * fe_values.normal_vector(q);
+
+                  data.psi_grads[q][k] -=
+                    (data.function_grads[q][k] + f_dot_n_times_n);
+                }
+          else
+            for (unsigned int k = 0; k < n_components; ++k)
+              for (unsigned int q = 0; q < n_q_points; ++q)
+                for (unsigned int d = 0; d < spacedim; ++d)
+                  data.psi_grads[q][k][d] -= data.function_grads[q][k][d];
+        }
+
+      double diff      = 0;
+      Number diff_mean = 0;
+
+      // First work on function values:
+      switch (norm)
+        {
+          case mean:
+            // Compute values in quadrature points and integrate
+            for (unsigned int q = 0; q < n_q_points; ++q)
+              {
+                Number sum = 0;
+                for (unsigned int k = 0; k < n_components; ++k)
+                  if (data.weight_vectors[q](k) != 0)
+                    sum += data.psi_values[q](k) * data.weight_vectors[q](k);
+                diff_mean += sum * fe_values.JxW(q);
+              }
+            break;
+
+          case Lp_norm:
+          case L1_norm:
+          case W1p_norm:
+            // Compute values in quadrature points and integrate
+            for (unsigned int q = 0; q < n_q_points; ++q)
+              {
+                double sum = 0;
+                for (unsigned int k = 0; k < n_components; ++k)
+                  if (data.weight_vectors[q](k) != 0)
+                    sum += std::pow(static_cast<double>(
+                                      numbers::NumberTraits<Number>::abs_square(
+                                        data.psi_values[q](k))),
+                                    exponent / 2.) *
+                           data.weight_vectors[q](k);
+                diff += sum * fe_values.JxW(q);
+              }
+
+            // Compute the root only if no derivative values are added later
+            if (!(update_flags & update_gradients))
+              diff = std::pow(diff, 1. / exponent);
+            break;
+
+          case L2_norm:
+          case H1_norm:
+            // Compute values in quadrature points and integrate
+            for (unsigned int q = 0; q < n_q_points; ++q)
+              {
+                double sum = 0;
+                for (unsigned int k = 0; k < n_components; ++k)
+                  if (data.weight_vectors[q](k) != 0)
+                    sum += numbers::NumberTraits<Number>::abs_square(
+                             data.psi_values[q](k)) *
+                           data.weight_vectors[q](k);
+                diff += sum * fe_values.JxW(q);
+              }
+            // Compute the root only, if no derivative values are added later
+            if (norm == L2_norm)
+              diff = std::sqrt(diff);
+            break;
+
+          case Linfty_norm:
+          case W1infty_norm:
+            for (unsigned int q = 0; q < n_q_points; ++q)
+              for (unsigned int k = 0; k < n_components; ++k)
+                if (data.weight_vectors[q](k) != 0)
+                  diff = std::max(diff,
+                                  double(std::abs(data.psi_values[q](k) *
+                                                  data.weight_vectors[q](k))));
+            break;
+
+          case H1_seminorm:
+          case Hdiv_seminorm:
+          case W1p_seminorm:
+          case W1infty_seminorm:
+            // function values are not used for these norms
+            break;
+
+          default:
+            Assert(false, ExcNotImplemented());
+            break;
+        }
+
+      // Now compute terms depending on derivatives:
+      switch (norm)
+        {
+          case W1p_seminorm:
+          case W1p_norm:
+            for (unsigned int q = 0; q < n_q_points; ++q)
+              {
+                double sum = 0;
+                for (unsigned int k = 0; k < n_components; ++k)
+                  if (data.weight_vectors[q](k) != 0)
+                    sum += std::pow(data.psi_grads[q][k].norm_square(),
+                                    exponent / 2.) *
+                           data.weight_vectors[q](k);
+                diff += sum * fe_values.JxW(q);
+              }
+            diff = std::pow(diff, 1. / exponent);
+            break;
+
+          case H1_seminorm:
+          case H1_norm:
+            for (unsigned int q = 0; q < n_q_points; ++q)
+              {
+                double sum = 0;
+                for (unsigned int k = 0; k < n_components; ++k)
+                  if (data.weight_vectors[q](k) != 0)
+                    sum += data.psi_grads[q][k].norm_square() *
+                           data.weight_vectors[q](k);
+                diff += sum * fe_values.JxW(q);
+              }
+            diff = std::sqrt(diff);
+            break;
+
+          case Hdiv_seminorm:
+            for (unsigned int q = 0; q < n_q_points; ++q)
+              {
+                unsigned int idx = 0;
+                if (weight != nullptr)
+                  for (; idx < n_components; ++idx)
+                    if (data.weight_vectors[0](idx) > 0)
+                      break;
+
+                Assert(
+                  n_components >= idx + dim,
+                  ExcMessage(
+                    "You can only ask for the Hdiv norm for a finite element "
+                    "with at least 'dim' components. In that case, this function "
+                    "will find the index of the first non-zero weight and take "
+                    "the divergence of the 'dim' components that follow it."));
+
+                Number sum = 0;
+                // take the trace of the derivatives scaled by the weight and
+                // square it
+                for (unsigned int k = idx; k < idx + dim; ++k)
+                  if (data.weight_vectors[q](k) != 0)
+                    sum += data.psi_grads[q][k][k - idx] *
+                           std::sqrt(data.weight_vectors[q](k));
+                diff += numbers::NumberTraits<Number>::abs_square(sum) *
+                        fe_values.JxW(q);
+              }
+            diff = std::sqrt(diff);
+            break;
 
-      tmp_values.resize(n_q_points);
-      tmp_vector_values.resize(n_q_points, Vector<double>(n_components));
-      tmp_gradients.resize(n_q_points);
-      tmp_vector_gradients.resize(
-        n_q_points, std::vector<Tensor<1, spacedim>>(n_components));
-    }
+          case W1infty_seminorm:
+          case W1infty_norm:
+            {
+              double t = 0;
+              for (unsigned int q = 0; q < n_q_points; ++q)
+                for (unsigned int k = 0; k < n_components; ++k)
+                  if (data.weight_vectors[q](k) != 0)
+                    for (unsigned int d = 0; d < dim; ++d)
+                      t = std::max(t,
+                                   double(std::abs(data.psi_grads[q][k][d]) *
+                                          data.weight_vectors[q](k)));
 
-    namespace internal
-    {
-      template <typename number>
-      double
-      mean_to_double(const number &mean_value)
-      {
-        return mean_value;
-      }
+              // then add seminorm to norm if that had previously been computed
+              diff += t;
+            }
+            break;
+          default:
+            break;
+        }
 
-      template <typename number>
-      double
-      mean_to_double(const std::complex<number> &mean_value)
-      {
-        // we need to return double as a norm, but mean value is a complex
-        // number. Panic and return real-part while warning the user that
-        // he shall never do that.
-        Assert(
-          false,
-          ExcMessage(
-            "Mean value norm is not implemented for complex-valued vectors"));
-        return mean_value.real();
-      }
-    } // namespace internal
+      if (norm == mean)
+        diff = internal::mean_to_double(diff_mean);
 
+      // append result of this cell to the end of the vector
+      AssertIsFinite(diff);
+      return diff;
+    }
 
-    // avoid compiling inner function for many vector types when we always
-    // really do the same thing by putting the main work into this helper
-    // function
     template <int dim, int spacedim, typename Number>
-    double
-    integrate_difference_inner(const Function<spacedim> &exact_solution,
-                               const NormType &          norm,
-                               const Function<spacedim> *weight,
-                               const UpdateFlags         update_flags,
-                               const double              exponent,
-                               const unsigned int        n_components,
-                               IDScratchData<dim, spacedim, Number> &data)
+    DEAL_II_DEPRECATED
+      typename std::enable_if<!std::is_same<Number, double>::value,
+                              double>::type
+      integrate_difference_inner(
+        const Function<spacedim> &                      exact_solution,
+        const NormType &                                norm,
+        const Function<spacedim> *                      weight,
+        const UpdateFlags                               update_flags,
+        const double                                    exponent,
+        const unsigned int                              n_components,
+        DeprecatedIDScratchData<dim, spacedim, Number> &data)
     {
       const bool                             fe_is_system = (n_components != 1);
       const dealii::FEValues<dim, spacedim> &fe_values =
@@ -8193,6 +8559,131 @@ namespace VectorTools
               typename DoFHandlerType,
               int spacedim>
     static void
+    do_integrate_difference(
+      const dealii::hp::MappingCollection<dim, spacedim> &     mapping,
+      const DoFHandlerType &                                   dof,
+      const InVector &                                         fe_function,
+      const Function<spacedim, typename InVector::value_type> &exact_solution,
+      OutVector &                                              difference,
+      const dealii::hp::QCollection<dim> &                     q,
+      const NormType &                                         norm,
+      const Function<spacedim> *                               weight,
+      const double                                             exponent_1)
+    {
+      using Number = typename InVector::value_type;
+      // we mark the "exponent" parameter to this function "const" since it is
+      // strictly incoming, but we need to set it to something different later
+      // on, if necessary, so have a read-write version of it:
+      double exponent = exponent_1;
+
+      const unsigned int n_components = dof.get_fe(0).n_components();
+
+      Assert(exact_solution.n_components == n_components,
+             ExcDimensionMismatch(exact_solution.n_components, n_components));
+
+      if (weight != nullptr)
+        {
+          Assert((weight->n_components == 1) ||
+                   (weight->n_components == n_components),
+                 ExcDimensionMismatch(weight->n_components, n_components));
+        }
+
+      difference.reinit(dof.get_triangulation().n_active_cells());
+
+      switch (norm)
+        {
+          case L2_norm:
+          case H1_seminorm:
+          case H1_norm:
+          case Hdiv_seminorm:
+            exponent = 2.;
+            break;
+
+          case L1_norm:
+            exponent = 1.;
+            break;
+
+          default:
+            break;
+        }
+
+      UpdateFlags update_flags =
+        UpdateFlags(update_quadrature_points | update_JxW_values);
+      switch (norm)
+        {
+          case H1_seminorm:
+          case Hdiv_seminorm:
+          case W1p_seminorm:
+          case W1infty_seminorm:
+            update_flags |= UpdateFlags(update_gradients);
+            if (spacedim == dim + 1)
+              update_flags |= UpdateFlags(update_normal_vectors);
+
+            break;
+
+          case H1_norm:
+          case W1p_norm:
+          case W1infty_norm:
+            update_flags |= UpdateFlags(update_gradients);
+            if (spacedim == dim + 1)
+              update_flags |= UpdateFlags(update_normal_vectors);
+            DEAL_II_FALLTHROUGH;
+
+          default:
+            update_flags |= UpdateFlags(update_values);
+            break;
+        }
+
+      const dealii::hp::FECollection<dim, spacedim> &fe_collection =
+        dof.get_fe_collection();
+      IDScratchData<dim, spacedim, Number> data(mapping,
+                                                fe_collection,
+                                                q,
+                                                update_flags);
+
+      // loop over all cells
+      for (typename DoFHandlerType::active_cell_iterator cell =
+             dof.begin_active();
+           cell != dof.end();
+           ++cell)
+        if (cell->is_locally_owned())
+          {
+            // initialize for this cell
+            data.x_fe_values.reinit(cell);
+
+            const dealii::FEValues<dim, spacedim> &fe_values =
+              data.x_fe_values.get_present_fe_values();
+            const unsigned int n_q_points = fe_values.n_quadrature_points;
+            data.resize_vectors(n_q_points, n_components);
+
+            if (update_flags & update_values)
+              fe_values.get_function_values(fe_function, data.function_values);
+            if (update_flags & update_gradients)
+              fe_values.get_function_gradients(fe_function,
+                                               data.function_grads);
+
+            difference(cell->active_cell_index()) =
+              integrate_difference_inner<dim, spacedim, Number>(exact_solution,
+                                                                norm,
+                                                                weight,
+                                                                update_flags,
+                                                                exponent,
+                                                                n_components,
+                                                                data);
+          }
+        else
+          // the cell is a ghost cell or is artificial. write a zero into the
+          // corresponding value of the returned vector
+          difference(cell->active_cell_index()) = 0;
+    }
+
+    template <int dim,
+              class InVector,
+              class OutVector,
+              typename DoFHandlerType,
+              int spacedim>
+    DEAL_II_DEPRECATED static typename std::enable_if<
+      !std::is_same<typename InVector::value_type, double>::value>::type
     do_integrate_difference(
       const dealii::hp::MappingCollection<dim, spacedim> &mapping,
       const DoFHandlerType &                              dof,
@@ -8270,10 +8761,10 @@ namespace VectorTools
 
       const dealii::hp::FECollection<dim, spacedim> &fe_collection =
         dof.get_fe_collection();
-      IDScratchData<dim, spacedim, Number> data(mapping,
-                                                fe_collection,
-                                                q,
-                                                update_flags);
+      DeprecatedIDScratchData<dim, spacedim, Number> data(mapping,
+                                                          fe_collection,
+                                                          q,
+                                                          update_flags);
 
       // loop over all cells
       for (typename DoFHandlerType::active_cell_iterator cell =
@@ -8317,6 +8808,32 @@ namespace VectorTools
 
   template <int dim, class InVector, class OutVector, int spacedim>
   void
+  integrate_difference(
+    const Mapping<dim, spacedim> &                           mapping,
+    const DoFHandler<dim, spacedim> &                        dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const Quadrature<dim> &                                  q,
+    const NormType &                                         norm,
+    const Function<spacedim> *                               weight,
+    const double                                             exponent)
+  {
+    internal ::do_integrate_difference(hp::MappingCollection<dim, spacedim>(
+                                         mapping),
+                                       dof,
+                                       fe_function,
+                                       exact_solution,
+                                       difference,
+                                       hp::QCollection<dim>(q),
+                                       norm,
+                                       weight,
+                                       exponent);
+  }
+
+  template <int dim, class InVector, class OutVector, int spacedim>
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(const Mapping<dim, spacedim> &   mapping,
                        const DoFHandler<dim, spacedim> &dof,
                        const InVector &                 fe_function,
@@ -8342,6 +8859,32 @@ namespace VectorTools
 
   template <int dim, class InVector, class OutVector, int spacedim>
   void
+  integrate_difference(
+    const DoFHandler<dim, spacedim> &                        dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const Quadrature<dim> &                                  q,
+    const NormType &                                         norm,
+    const Function<spacedim> *                               weight,
+    const double                                             exponent)
+  {
+    internal ::do_integrate_difference(
+      hp::StaticMappingQ1<dim, spacedim>::mapping_collection,
+      dof,
+      fe_function,
+      exact_solution,
+      difference,
+      hp::QCollection<dim>(q),
+      norm,
+      weight,
+      exponent);
+  }
+
+
+  template <int dim, class InVector, class OutVector, int spacedim>
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(const DoFHandler<dim, spacedim> &dof,
                        const InVector &                 fe_function,
                        const Function<spacedim> &       exact_solution,
@@ -8367,6 +8910,32 @@ namespace VectorTools
 
   template <int dim, class InVector, class OutVector, int spacedim>
   void
+  integrate_difference(
+    const dealii::hp::MappingCollection<dim, spacedim> &     mapping,
+    const dealii::hp::DoFHandler<dim, spacedim> &            dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const dealii::hp::QCollection<dim> &                     q,
+    const NormType &                                         norm,
+    const Function<spacedim> *                               weight,
+    const double                                             exponent)
+  {
+    internal ::do_integrate_difference(hp::MappingCollection<dim, spacedim>(
+                                         mapping),
+                                       dof,
+                                       fe_function,
+                                       exact_solution,
+                                       difference,
+                                       q,
+                                       norm,
+                                       weight,
+                                       exponent);
+  }
+
+  template <int dim, class InVector, class OutVector, int spacedim>
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(
     const dealii::hp::MappingCollection<dim, spacedim> &mapping,
     const dealii::hp::DoFHandler<dim, spacedim> &       dof,
@@ -8393,6 +8962,31 @@ namespace VectorTools
 
   template <int dim, class InVector, class OutVector, int spacedim>
   void
+  integrate_difference(
+    const dealii::hp::DoFHandler<dim, spacedim> &            dof,
+    const InVector &                                         fe_function,
+    const Function<spacedim, typename InVector::value_type> &exact_solution,
+    OutVector &                                              difference,
+    const dealii::hp::QCollection<dim> &                     q,
+    const NormType &                                         norm,
+    const Function<spacedim> *                               weight,
+    const double                                             exponent)
+  {
+    internal ::do_integrate_difference(
+      hp::StaticMappingQ1<dim, spacedim>::mapping_collection,
+      dof,
+      fe_function,
+      exact_solution,
+      difference,
+      q,
+      norm,
+      weight,
+      exponent);
+  }
+
+  template <int dim, class InVector, class OutVector, int spacedim>
+  DEAL_II_DEPRECATED typename std::enable_if<
+    !std::is_same<typename InVector::value_type, double>::value>::type
   integrate_difference(const dealii::hp::DoFHandler<dim, spacedim> &dof,
                        const InVector &                             fe_function,
                        const Function<spacedim> &          exact_solution,
index 6b05ae0d7ba0785a91052e575ebfd702fedc3791..62b17d219775b8c890969c1af39e12b9d1d9503e 100644 (file)
@@ -29,7 +29,7 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
         const Mapping<deal_II_dimension, deal_II_space_dimension> &,
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
-        const Function<deal_II_space_dimension> &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
         Vector<float> &,
         const Quadrature<deal_II_dimension> &,
         const NormType &,
@@ -43,7 +43,7 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
                            deal_II_space_dimension>(
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
-        const Function<deal_II_space_dimension> &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
         Vector<float> &,
         const Quadrature<deal_II_dimension> &,
         const NormType &,
@@ -58,7 +58,7 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
         const Mapping<deal_II_dimension, deal_II_space_dimension> &,
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
-        const Function<deal_II_space_dimension> &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
         Vector<double> &,
         const Quadrature<deal_II_dimension> &,
         const NormType &,
@@ -72,7 +72,7 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
                            deal_II_space_dimension>(
         const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
-        const Function<deal_II_space_dimension> &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
         Vector<double> &,
         const Quadrature<deal_II_dimension> &,
         const NormType &,
@@ -88,7 +88,7 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
           &,
         const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
-        const Function<deal_II_space_dimension> &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
         Vector<double> &,
         const hp::QCollection<deal_II_dimension> &,
         const NormType &,
@@ -102,7 +102,7 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
                            deal_II_space_dimension>(
         const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
-        const Function<deal_II_space_dimension> &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
         Vector<double> &,
         const hp::QCollection<deal_II_dimension> &,
         const NormType &,
@@ -118,7 +118,7 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
           &,
         const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
-        const Function<deal_II_space_dimension> &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
         Vector<float> &,
         const hp::QCollection<deal_II_dimension> &,
         const NormType &,
@@ -132,12 +132,660 @@ for (VEC : VECTOR_TYPES; deal_II_dimension : DIMENSIONS;
                            deal_II_space_dimension>(
         const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
         const VEC &,
+        const Function<deal_II_space_dimension, typename VEC::value_type> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+    \}
+#endif
+  }
+
+// All the functions in this for loop are deprecated
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
+  {
+#if deal_II_dimension <= deal_II_space_dimension
+    namespace VectorTools
+    \{
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+#  endif
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float>,
+          &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+#  endif
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const Mapping<deal_II_dimension, deal_II_space_dimension> &,
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+#  endif
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const Quadrature<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+#  endif
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+#  endif
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<double>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<double> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+#  endif
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension>
+          &,
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+#  endif
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::Vector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::distributed::BlockVector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::distributed::BlockVector<float> &,
+        const Function<deal_II_space_dimension> &,
+        Vector<float> &,
+        const hp::QCollection<deal_II_dimension> &,
+        const NormType &,
+        const Function<deal_II_space_dimension> *,
+        const double);
+
+#  ifdef DEAL_II_TRILINOS_WITH_TPETRA
+      template void
+      integrate_difference<deal_II_dimension,
+                           LinearAlgebra::TpetraWrappers::Vector<float>,
+                           Vector<float>,
+                           deal_II_space_dimension>(
+        const hp::DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+        const LinearAlgebra::TpetraWrappers::Vector<float> &,
         const Function<deal_II_space_dimension> &,
         Vector<float> &,
         const hp::QCollection<deal_II_dimension> &,
         const NormType &,
         const Function<deal_II_space_dimension> *,
         const double);
+#  endif
 
     \}
 #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.