]> https://gitweb.dealii.org/ - dealii.git/commitdiff
extend VectorTools to complex-valued PETSc vector 2217/head
authorDenis Davydov <davydden@gmail.com>
Mon, 22 Feb 2016 12:37:20 +0000 (13:37 +0100)
committerDenis Davydov <davydden@gmail.com>
Mon, 22 Feb 2016 17:00:33 +0000 (18:00 +0100)
doc/news/changes.h
include/deal.II/numerics/vector_tools.h
include/deal.II/numerics/vector_tools.templates.h
source/numerics/vector_tools_interpolate.inst.in
source/numerics/vector_tools_mean_value.inst.in
source/numerics/vector_tools_point_value.inst.in

index 6bff97718f5497811f7c829575c45f81c29e3231..1f84ee4d2df5be52f87b9c6bbd0a4e2fb066ec25 100644 (file)
@@ -38,8 +38,9 @@ inconvenience this causes.
 </p>
 
 <ol>
-  <li> Changed: ConstraintMatrix::distribute_local_to_global() now requires
-  matching data types. This is done to correctly handle complex-valued case.
+  <li> Changed: ConstraintMatrix::distribute_local_to_global() and numerous
+  functions in VectorTools namespace now require matching data types.
+  This is done to correctly handle complex-valued case.
   <br>
   (Denis Davydov, 2016/02/22)
   </li>
index 635ee356d2b4b5ded407a605507646008a2ff02f..64f046f0c4486f7a6eab79369893b6bea590196b 100644 (file)
@@ -403,7 +403,7 @@ namespace VectorTools
   template <typename VectorType, int dim, int spacedim, template <int, int> class DoFHandlerType>
   void interpolate (const Mapping<dim,spacedim>        &mapping,
                     const DoFHandlerType<dim,spacedim> &dof,
-                    const Function<spacedim,double>    &function,
+                    const Function<spacedim,typename VectorType::value_type>    &function,
                     VectorType                         &vec);
 
   /**
@@ -412,7 +412,7 @@ namespace VectorTools
    */
   template <typename VectorType, typename DoFHandlerType>
   void interpolate (const DoFHandlerType                                   &dof,
-                    const Function<DoFHandlerType::space_dimension,double> &function,
+                    const Function<DoFHandlerType::space_dimension,typename VectorType::value_type> &function,
                     VectorType                                             &vec);
 
   /**
@@ -489,7 +489,7 @@ namespace VectorTools
   interpolate_based_on_material_id
   (const Mapping<DoFHandlerType::dimension, DoFHandlerType::space_dimension> &mapping,
    const DoFHandlerType                                                  &dof_handler,
-   const std::map<types::material_id, const Function<DoFHandlerType::space_dimension, double> *> &function_map,
+   const std::map<types::material_id, const Function<DoFHandlerType::space_dimension, typename VectorType::value_type> *> &function_map,
    VectorType                                                            &dst,
    const ComponentMask                                                   &component_mask = ComponentMask());
 
@@ -592,7 +592,7 @@ namespace VectorTools
                 const DoFHandler<dim,spacedim>  &dof,
                 const ConstraintMatrix          &constraints,
                 const Quadrature<dim>           &quadrature,
-                const Function<spacedim,double> &function,
+                const Function<spacedim,typename VectorType::value_type> &function,
                 VectorType                      &vec,
                 const bool                      enforce_zero_boundary = false,
                 const Quadrature<dim-1>         &q_boundary = (dim > 1 ?
@@ -608,7 +608,7 @@ namespace VectorTools
   void project (const DoFHandler<dim,spacedim>  &dof,
                 const ConstraintMatrix          &constraints,
                 const Quadrature<dim>           &quadrature,
-                const Function<spacedim,double> &function,
+                const Function<spacedim,typename VectorType::value_type> &function,
                 VectorType                      &vec,
                 const bool                      enforce_zero_boundary = false,
                 const Quadrature<dim-1>         &q_boundary = (dim > 1 ?
@@ -625,7 +625,7 @@ namespace VectorTools
                 const hp::DoFHandler<dim,spacedim>         &dof,
                 const ConstraintMatrix                     &constraints,
                 const hp::QCollection<dim>                 &quadrature,
-                const Function<spacedim,double>            &function,
+                const Function<spacedim,typename VectorType::value_type>            &function,
                 VectorType                                 &vec,
                 const bool                                 enforce_zero_boundary = false,
                 const hp::QCollection<dim-1> &q_boundary = hp::QCollection<dim-1>(dim > 1 ?
@@ -641,7 +641,7 @@ namespace VectorTools
   void project (const hp::DoFHandler<dim,spacedim> &dof,
                 const ConstraintMatrix             &constraints,
                 const hp::QCollection<dim>         &quadrature,
-                const Function<spacedim,double>    &function,
+                const Function<spacedim,typename VectorType::value_type>    &function,
                 VectorType                         &vec,
                 const bool                         enforce_zero_boundary = false,
                 const hp::QCollection<dim-1>       &q_boundary = hp::QCollection<dim-1>(dim > 1 ?
@@ -1989,8 +1989,8 @@ namespace VectorTools
   template <int dim, typename VectorType, int spacedim>
   void point_difference (const DoFHandler<dim,spacedim>  &dof,
                          const VectorType                &fe_function,
-                         const Function<spacedim,double> &exact_solution,
-                         Vector<double>                  &difference,
+                         const Function<spacedim,typename VectorType::value_type> &exact_solution,
+                         Vector<typename VectorType::value_type>                  &difference,
                          const Point<spacedim>           &point);
 
   /**
@@ -2009,8 +2009,8 @@ namespace VectorTools
   void point_difference (const Mapping<dim, spacedim>    &mapping,
                          const DoFHandler<dim,spacedim>  &dof,
                          const VectorType                &fe_function,
-                         const Function<spacedim,double> &exact_solution,
-                         Vector<double>                  &difference,
+                         const Function<spacedim,typename VectorType::value_type> &exact_solution,
+                         Vector<typename VectorType::value_type>                  &difference,
                          const Point<spacedim>           &point);
 
   /**
@@ -2029,7 +2029,7 @@ namespace VectorTools
   point_value (const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
                const Point<spacedim>          &point,
-               Vector<double>                 &value);
+               Vector<typename VectorType::value_type>                 &value);
 
   /**
    * Same as above for hp.
@@ -2042,7 +2042,7 @@ namespace VectorTools
   point_value (const hp::DoFHandler<dim,spacedim> &dof,
                const VectorType                   &fe_function,
                const Point<spacedim>              &point,
-               Vector<double>                     &value);
+               Vector<typename VectorType::value_type>                     &value);
 
   /**
    * Evaluate a scalar finite element function defined by the given DoFHandler
@@ -2060,7 +2060,7 @@ namespace VectorTools
    * exception of type VectorTools::ExcPointNotAvailableHere is thrown.
    */
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
                const Point<spacedim>          &point);
@@ -2072,7 +2072,7 @@ namespace VectorTools
    * exception of type VectorTools::ExcPointNotAvailableHere is thrown.
    */
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const hp::DoFHandler<dim,spacedim> &dof,
                const VectorType                   &fe_function,
                const Point<spacedim>              &point);
@@ -2094,7 +2094,7 @@ namespace VectorTools
                const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
                const Point<spacedim>          &point,
-               Vector<double>                 &value);
+               Vector<typename VectorType::value_type>                 &value);
 
   /**
    * Same as above for hp.
@@ -2108,7 +2108,7 @@ namespace VectorTools
                const hp::DoFHandler<dim,spacedim>         &dof,
                const VectorType                           &fe_function,
                const Point<spacedim>                      &point,
-               Vector<double>                             &value);
+               Vector<typename VectorType::value_type>                             &value);
 
   /**
    * Evaluate a scalar finite element function defined by the given DoFHandler
@@ -2122,7 +2122,7 @@ namespace VectorTools
    * exception of type VectorTools::ExcPointNotAvailableHere is thrown.
    */
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const Mapping<dim,spacedim>    &mapping,
                const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
@@ -2135,7 +2135,7 @@ namespace VectorTools
    * exception of type VectorTools::ExcPointNotAvailableHere is thrown.
    */
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const hp::MappingCollection<dim,spacedim> &mapping,
                const hp::DoFHandler<dim,spacedim>        &dof,
                const VectorType                          &fe_function,
@@ -2346,21 +2346,23 @@ namespace VectorTools
    * mean value and subtract it right inside the evaluation routine.
    */
   template <int dim, typename VectorType, int spacedim>
-  double compute_mean_value (const Mapping<dim, spacedim>   &mapping,
-                             const DoFHandler<dim,spacedim> &dof,
-                             const Quadrature<dim>          &quadrature,
-                             const VectorType               &v,
-                             const unsigned int             component);
+  typename VectorType::value_type
+  compute_mean_value (const Mapping<dim, spacedim>   &mapping,
+                      const DoFHandler<dim,spacedim> &dof,
+                      const Quadrature<dim>          &quadrature,
+                      const VectorType               &v,
+                      const unsigned int             component);
 
   /**
    * Calls the other compute_mean_value() function, see above, with
    * <tt>mapping=MappingQGeneric@<dim@>(1)</tt>.
    */
   template <int dim, typename VectorType, int spacedim>
-  double compute_mean_value (const DoFHandler<dim,spacedim> &dof,
-                             const Quadrature<dim>          &quadrature,
-                             const VectorType               &v,
-                             const unsigned int             component);
+  typename VectorType::value_type
+  compute_mean_value (const DoFHandler<dim,spacedim> &dof,
+                      const Quadrature<dim>          &quadrature,
+                      const VectorType               &v,
+                      const unsigned int             component);
   //@}
   /**
    * Geometrical interpolation
index e6889a36b415401d333bc4168cf0070d087ee451..6ed3e4d05bece179d4f92a87b3a30b7d75a81814 100644 (file)
@@ -77,9 +77,10 @@ namespace VectorTools
             template <int, int> class DoFHandlerType>
   void interpolate (const Mapping<dim,spacedim>        &mapping,
                     const DoFHandlerType<dim,spacedim> &dof,
-                    const Function<spacedim>           &function,
+                    const Function<spacedim, typename VectorType::value_type>           &function,
                     VectorType                         &vec)
   {
+    typedef typename VectorType::value_type number;
     Assert (vec.size() == dof.n_dofs(),
             ExcDimensionMismatch (vec.size(), dof.n_dofs()));
     Assert (dof.get_fe().n_components() == function.n_components,
@@ -187,8 +188,8 @@ namespace VectorTools
     // have two versions, one for system fe
     // and one for scalar ones, to take the
     // more efficient one respectively
-    std::vector<std::vector<double> >         function_values_scalar(fe.size());
-    std::vector<std::vector<Vector<double> > > function_values_system(fe.size());
+    std::vector<std::vector<number> >         function_values_scalar(fe.size());
+    std::vector<std::vector<Vector<number> > > function_values_system(fe.size());
 
     // Make a quadrature rule from support points
     // to feed it into FEValues
@@ -234,7 +235,7 @@ namespace VectorTools
                   // all components
                   function_values_system[fe_index]
                   .resize (n_rep_points[fe_index],
-                           Vector<double> (fe[fe_index].n_components()));
+                           Vector<number> (fe[fe_index].n_components()));
                   function.vector_value_list (rep_points,
                                               function_values_system[fe_index]);
                   // distribute the function
@@ -273,7 +274,7 @@ namespace VectorTools
 
   template <typename VectorType, typename DoFHandlerType>
   void interpolate (const DoFHandlerType                            &dof,
-                    const Function<DoFHandlerType::space_dimension> &function,
+                    const Function<DoFHandlerType::space_dimension,typename VectorType::value_type> &function,
                     VectorType                                      &vec)
   {
     interpolate(StaticMappingQ1<DoFHandlerType::dimension,
@@ -294,8 +295,9 @@ namespace VectorTools
                const InVector                  &data_1,
                OutVector                       &data_2)
   {
-    Vector<double> cell_data_1(dof_1.get_fe().dofs_per_cell);
-    Vector<double> cell_data_2(dof_2.get_fe().dofs_per_cell);
+    typedef typename OutVector::value_type number;
+    Vector<number> cell_data_1(dof_1.get_fe().dofs_per_cell);
+    Vector<number> cell_data_2(dof_2.get_fe().dofs_per_cell);
 
     std::vector<short unsigned int> touch_count (dof_2.n_dofs(), 0); //TODO: check on datatype... kinda strange (UK)
     std::vector<types::global_dof_index>       local_dof_indices (dof_2.get_fe().dofs_per_cell);
@@ -6102,6 +6104,25 @@ namespace VectorTools
                                    std::vector<Tensor<1,spacedim> >(n_components));
     }
 
+    namespace
+    {
+      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. Panick 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();
+      }
+    }
+
 
     // avoid compiling inner function for many vector types when we always
     // really do the same thing by putting the main work into this helper
@@ -6219,6 +6240,7 @@ namespace VectorTools
         }
 
       double diff = 0;
+      Number diff_mean = 0;
 
       // First work on function values:
       switch (norm)
@@ -6227,10 +6249,10 @@ namespace VectorTools
           // Compute values in quadrature points and integrate
           for (unsigned int q=0; q<n_q_points; ++q)
             {
-              double sum = 0;
+              Number sum = 0;
               for (unsigned int k=0; k<n_components; ++k)
                 sum += data.psi_values[q](k) * data.weight_vectors[q](k);
-              diff += sum * fe_values.JxW(q);
+              diff_mean += sum * fe_values.JxW(q);
             }
           break;
 
@@ -6243,7 +6265,7 @@ namespace VectorTools
               double sum = 0;
               for (unsigned int k=0; k<n_components; ++k)
                 sum += std::pow(
-                         static_cast<double>(data.psi_values[q](k)*data.psi_values[q](k)),
+                         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);
             }
@@ -6260,7 +6282,7 @@ namespace VectorTools
             {
               double sum = 0;
               for (unsigned int k=0; k<n_components; ++k)
-                sum += data.psi_values[q](k) * data.psi_values[q](k) *
+                sum += numbers::NumberTraits<Number>::abs_square(data.psi_values[q](k)) *
                        data.weight_vectors[q](k);
               diff += sum * fe_values.JxW(q);
             }
@@ -6299,7 +6321,7 @@ namespace VectorTools
               double sum = 0;
               for (unsigned int k=0; k<n_components; ++k)
                 sum += std::pow(
-                         static_cast<double>(data.psi_grads[q][k]*data.psi_grads[q][k]),
+                         data.psi_grads[q][k].norm_square(),
                          exponent/2.) * data.weight_vectors[q](k);
               diff += sum * fe_values.JxW(q);
             }
@@ -6312,7 +6334,7 @@ namespace VectorTools
             {
               double sum = 0;
               for (unsigned int k=0; k<n_components; ++k)
-                sum += (data.psi_grads[q][k] * data.psi_grads[q][k]) *
+                sum += data.psi_grads[q][k].norm_square() *
                        data.weight_vectors[q](k);
               diff += sum * fe_values.JxW(q);
             }
@@ -6326,11 +6348,11 @@ namespace VectorTools
                       ExcMessage ("You can only ask for the Hdiv norm for a finite element "
                                   "with at least 'dim' components. In that case, this function "
                                   "will take the divergence of the first 'dim' components."));
-              double sum = 0;
+              Number sum = 0;
               // take the trace of the derivatives scaled by the weight and square it
               for (unsigned int k=0; k<dim; ++k)
                 sum += data.psi_grads[q][k][k] * std::sqrt(data.weight_vectors[q](k));
-              diff += sum * sum * fe_values.JxW(q);
+              diff += numbers::NumberTraits<Number>::abs_square(sum) * fe_values.JxW(q);
             }
           diff = std::sqrt(diff);
           break;
@@ -6354,6 +6376,9 @@ namespace VectorTools
           break;
         }
 
+      if (norm == mean)
+        diff = mean_to_double(diff_mean);
+
       // append result of this cell to the end of the vector
       AssertIsFinite(diff);
       return diff;
@@ -6554,8 +6579,8 @@ namespace VectorTools
   void
   point_difference (const DoFHandler<dim,spacedim> &dof,
                     const VectorType               &fe_function,
-                    const Function<spacedim>       &exact_function,
-                    Vector<double>                 &difference,
+                    const Function<spacedim, typename VectorType::value_type>       &exact_function,
+                    Vector<typename VectorType::value_type>                 &difference,
                     const Point<spacedim>          &point)
   {
     point_difference(StaticMappingQ1<dim>::mapping,
@@ -6572,8 +6597,8 @@ namespace VectorTools
   point_difference (const Mapping<dim, spacedim>   &mapping,
                     const DoFHandler<dim,spacedim> &dof,
                     const VectorType               &fe_function,
-                    const Function<spacedim>       &exact_function,
-                    Vector<double>                 &difference,
+                    const Function<spacedim, typename VectorType::value_type>       &exact_function,
+                    Vector<typename VectorType::value_type>                 &difference,
                     const Point<spacedim>          &point)
   {
     typedef typename VectorType::value_type Number;
@@ -6618,7 +6643,7 @@ namespace VectorTools
   point_value (const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
                const Point<spacedim>          &point,
-               Vector<double>                 &value)
+               Vector<typename VectorType::value_type>                 &value)
   {
 
     point_value (StaticMappingQ1<dim,spacedim>::mapping,
@@ -6634,7 +6659,7 @@ namespace VectorTools
   point_value (const hp::DoFHandler<dim,spacedim> &dof,
                const VectorType                   &fe_function,
                const Point<spacedim>              &point,
-               Vector<double>                     &value)
+               Vector<typename VectorType::value_type>                     &value)
   {
     point_value(hp::StaticMappingQ1<dim,spacedim>::mapping_collection,
                 dof,
@@ -6645,7 +6670,7 @@ namespace VectorTools
 
 
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
                const Point<spacedim>          &point)
@@ -6658,7 +6683,7 @@ namespace VectorTools
 
 
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const hp::DoFHandler<dim,spacedim> &dof,
                const VectorType                   &fe_function,
                const Point<spacedim>              &point)
@@ -6676,7 +6701,7 @@ namespace VectorTools
                const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
                const Point<spacedim>          &point,
-               Vector<double>                 &value)
+               Vector<typename VectorType::value_type>                 &value)
   {
     typedef typename VectorType::value_type Number;
     const FiniteElement<dim> &fe = dof.get_fe();
@@ -6717,7 +6742,7 @@ namespace VectorTools
                const hp::DoFHandler<dim,spacedim>         &dof,
                const VectorType                           &fe_function,
                const Point<spacedim>                      &point,
-               Vector<double>                             &value)
+               Vector<typename VectorType::value_type>                             &value)
   {
     typedef typename VectorType::value_type Number;
     const hp::FECollection<dim, spacedim> &fe = dof.get_fe();
@@ -6752,7 +6777,7 @@ namespace VectorTools
 
 
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const Mapping<dim, spacedim>   &mapping,
                const DoFHandler<dim,spacedim> &dof,
                const VectorType               &fe_function,
@@ -6761,7 +6786,7 @@ namespace VectorTools
     Assert(dof.get_fe().n_components() == 1,
            ExcMessage ("Finite element is not scalar as is necessary for this function"));
 
-    Vector<double> value(1);
+    Vector<typename VectorType::value_type> value(1);
     point_value(mapping, dof, fe_function, point, value);
 
     return value(0);
@@ -6769,7 +6794,7 @@ namespace VectorTools
 
 
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   point_value (const hp::MappingCollection<dim, spacedim> &mapping,
                const hp::DoFHandler<dim,spacedim>         &dof,
                const VectorType                           &fe_function,
@@ -6778,7 +6803,7 @@ namespace VectorTools
     Assert(dof.get_fe().n_components() == 1,
            ExcMessage ("Finite element is not scalar as is necessary for this function"));
 
-    Vector<double> value(1);
+    Vector<typename VectorType::value_type> value(1);
     point_value(mapping, dof, fe_function, point, value);
 
     return value(0);
@@ -6992,7 +7017,8 @@ namespace VectorTools
         for (unsigned int i=0; i<n; ++i)
           if (p_select[i])
             {
-              s += v(i);
+              typename VectorType::value_type vi = v(i);
+              s += vi;
               ++counter;
             }
         // Error out if we have not constrained anything. Note that in this
@@ -7030,7 +7056,7 @@ namespace VectorTools
 
 
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   compute_mean_value (const Mapping<dim, spacedim>   &mapping,
                       const DoFHandler<dim,spacedim> &dof,
                       const Quadrature<dim>          &quadrature,
@@ -7097,7 +7123,7 @@ namespace VectorTools
 
 
   template <int dim, typename VectorType, int spacedim>
-  double
+  typename VectorType::value_type
   compute_mean_value (const DoFHandler<dim,spacedim> &dof,
                       const Quadrature<dim>          &quadrature,
                       const VectorType               &v,
index bd96acf66fe37cf2d40606ea3f7daa4ca348d464..454b944ccf2c5517e86306f34c790832c8aa1278 100644 (file)
@@ -23,13 +23,13 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
         void interpolate
         (const Mapping<deal_II_dimension,deal_II_space_dimension> &,
          const DoFHandler<deal_II_dimension,deal_II_space_dimension>&,
-         const Function<deal_II_space_dimension>&,
+         const Function<deal_II_space_dimension, VEC::value_type>&,
          VEC&);
 
       template
         void interpolate
         (const DoFHandler<deal_II_dimension,deal_II_space_dimension>&,
-         const Function<deal_II_space_dimension>&,
+         const Function<deal_II_space_dimension, VEC::value_type>&,
          VEC&);
 
       template
@@ -64,12 +64,12 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
         void interpolate
         (const Mapping<deal_II_dimension>&,
          const hp::DoFHandler<deal_II_dimension>&,
-         const Function<deal_II_dimension>&,
+         const Function<deal_II_dimension, VEC::value_type>&,
          VEC&);
       template
         void interpolate
         (const hp::DoFHandler<deal_II_dimension>&,
-         const Function<deal_II_dimension>&,
+         const Function<deal_II_dimension, VEC::value_type>&,
          VEC&);
          
      template
index 6892ce087c8d6835a99d0205c2b5e1c1eb152c18..06d304eb707bd1bbff0685c447a6fab5b4817dfa 100644 (file)
@@ -20,7 +20,7 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
     namespace VectorTools \{
 
       template
-        double compute_mean_value<deal_II_dimension>
+        VEC::value_type compute_mean_value<deal_II_dimension>
         (const Mapping<deal_II_dimension,deal_II_space_dimension>&,
          const DoFHandler<deal_II_dimension,deal_II_space_dimension>&,
          const Quadrature<deal_II_dimension>&,
@@ -28,7 +28,7 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
          const unsigned int);
 
       template
-        double compute_mean_value<deal_II_dimension>
+        VEC::value_type compute_mean_value<deal_II_dimension>
         (const DoFHandler<deal_II_dimension,deal_II_space_dimension>&,
          const Quadrature<deal_II_dimension>&,
          const VEC&,
index c5ee68d2aa59c77fcb389778bb11980c9264b1b9..4c77fab7470cd98713a51ffb924503a90e8cd465 100644 (file)
@@ -27,10 +27,10 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
       const hp::DoFHandler<deal_II_dimension>&,
       const VEC&,
       const Point<deal_II_dimension>&,
-      Vector<double>&);
+      Vector<VEC::value_type>&);
 
   template
-    double point_value<deal_II_dimension> (
+    VEC::value_type point_value<deal_II_dimension> (
       const hp::DoFHandler<deal_II_dimension>&,
       const VEC&,
       const Point<deal_II_dimension>&);
@@ -41,10 +41,10 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
       const hp::DoFHandler<deal_II_dimension>&,
       const VEC&,
       const Point<deal_II_dimension>&,
-      Vector<double>&);
+      Vector<VEC::value_type>&);
 
   template
-    double point_value<deal_II_dimension> (
+    VEC::value_type point_value<deal_II_dimension> (
       const hp::MappingCollection<deal_II_dimension>&,
       const hp::DoFHandler<deal_II_dimension>&,
       const VEC&,
@@ -54,8 +54,8 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
         void point_difference<deal_II_dimension> (
           const DoFHandler<deal_II_dimension>&,
           const VEC&,
-          const Function<deal_II_dimension>&,
-          Vector<double>&,
+          const Function<deal_II_dimension, VEC::value_type>&,
+          Vector<VEC::value_type>&,
           const Point<deal_II_dimension>&);
 
       template
@@ -63,8 +63,8 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
           const Mapping<deal_II_dimension>&,
           const DoFHandler<deal_II_dimension>&,
           const VEC&,
-          const Function<deal_II_dimension>&,
-          Vector<double>&,
+          const Function<deal_II_dimension,VEC::value_type>&,
+          Vector<VEC::value_type>&,
           const Point<deal_II_dimension>&);
 
       template
@@ -72,10 +72,10 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
           const DoFHandler<deal_II_dimension>&,
           const VEC&,
           const Point<deal_II_dimension>&,
-          Vector<double>&);
+          Vector<VEC::value_type>&);
 
       template
-        double point_value<deal_II_dimension> (
+        VEC::value_type point_value<deal_II_dimension> (
           const DoFHandler<deal_II_dimension>&,
           const VEC&,
           const Point<deal_II_dimension>&);
@@ -86,10 +86,10 @@ for (VEC : SERIAL_VECTORS ; deal_II_dimension : DIMENSIONS; deal_II_space_dimens
           const DoFHandler<deal_II_dimension>&,
           const VEC&,
           const Point<deal_II_dimension>&,
-          Vector<double>&);
+          Vector<VEC::value_type>&);
 
       template
-        double point_value<deal_II_dimension> (
+        VEC::value_type point_value<deal_II_dimension> (
           const Mapping<deal_II_dimension>&,
           const DoFHandler<deal_II_dimension>&,
           const VEC&,

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.