]> https://gitweb.dealii.org/ - dealii.git/commitdiff
temporary fix of vector_tools.templates.h
authorDenis Davydov <davydden@gmail.com>
Sun, 1 Mar 2015 15:37:59 +0000 (16:37 +0100)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 8 May 2015 16:23:50 +0000 (11:23 -0500)
include/deal.II/numerics/vector_tools.templates.h

index b71d1ac8094ed5bfc4255003cc11fcc24123d832..0a2bed01e0d4333456d51fd44e61db05e78f3383 100644 (file)
@@ -6344,6 +6344,7 @@ namespace VectorTools
                              const Function<spacedim>   *weight,
                              const double           exponent_1)
     {
+      typedef typename InVector::value_type Number;
       // 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:
@@ -6406,6 +6407,11 @@ namespace VectorTools
       dealii::hp::FECollection<dim,spacedim> fe_collection (dof.get_fe());
       IDScratchData<dim,spacedim> data(mapping, fe_collection, q, update_flags);
 
+      //FIXME
+      // temporary vectors of consistent with InVector type
+      std::vector<dealii::Vector<Number>> function_values;
+      std::vector<std::vector<Tensor<1,spacedim,Number> >> function_grads;
+
       // loop over all cells
       for (typename DH::active_cell_iterator cell = dof.begin_active();
            cell != dof.end(); ++cell)
@@ -6418,10 +6424,24 @@ namespace VectorTools
             const unsigned int   n_q_points = fe_values.n_quadrature_points;
             data.resize_vectors (n_q_points, n_components);
 
+            //FIXME:
+            function_values.resize (n_q_points,
+                                    dealii::Vector<Number>(n_components));
+            function_grads.resize (n_q_points,
+                                   std::vector<Tensor<1,spacedim,Number> >(n_components));
+
             if (update_flags & update_values)
-              fe_values.get_function_values (fe_function, data.function_values);
+              fe_values.get_function_values (fe_function, function_values);
             if (update_flags & update_gradients)
-              fe_values.get_function_gradients (fe_function, data.function_grads);
+              fe_values.get_function_gradients (fe_function, function_grads);
+
+            // FIXME
+            for (unsigned int q = 0; q < n_q_points; q++)
+              for (unsigned int c = 0; c < n_components; c++)
+                {
+                  data.function_values[q][c] = function_values[q][c];
+                  data.function_grads[q][c] = function_grads[q][c];
+                }
 
             difference(cell->active_cell_index()) =
               integrate_difference_inner (exact_solution, norm, weight,
@@ -6545,6 +6565,7 @@ namespace VectorTools
                     Vector<double>        &difference,
                     const Point<spacedim>      &point)
   {
+    typedef typename InVector::value_type Number;
     const FiniteElement<dim> &fe = dof.get_fe();
 
     Assert(difference.size() == fe.n_components(),
@@ -6568,7 +6589,7 @@ namespace VectorTools
 
     // then use this to get at the values of
     // the given fe_function at this point
-    std::vector<Vector<double> > u_value(1, Vector<double> (fe.n_components()));
+    std::vector<Vector<Number> > u_value(1, Vector<Number> (fe.n_components()));
     fe_values.get_function_values(fe_function, u_value);
 
     if (fe.n_components() == 1)
@@ -6646,6 +6667,7 @@ namespace VectorTools
                const Point<spacedim>      &point,
                Vector<double>        &value)
   {
+    typedef typename InVector::value_type Number;
     const FiniteElement<dim> &fe = dof.get_fe();
 
     Assert(value.size() == fe.n_components(),
@@ -6671,7 +6693,7 @@ namespace VectorTools
 
     // then use this to get at the values of
     // the given fe_function at this point
-    std::vector<Vector<double> > u_value(1, Vector<double> (fe.n_components()));
+    std::vector<Vector<Number> > u_value(1, Vector<Number> (fe.n_components()));
     fe_values.get_function_values(fe_function, u_value);
 
     value = u_value[0];
@@ -6686,6 +6708,7 @@ namespace VectorTools
                const Point<spacedim>      &point,
                Vector<double>        &value)
   {
+    typedef typename InVector::value_type Number;
     const hp::FECollection<dim, spacedim> &fe = dof.get_fe();
 
     Assert(value.size() == fe.n_components(),
@@ -6710,7 +6733,7 @@ namespace VectorTools
 
     // then use this to get at the values of
     // the given fe_function at this point
-    std::vector<Vector<double> > u_value(1, Vector<double> (fe.n_components()));
+    std::vector<Vector<Number> > u_value(1, Vector<Number> (fe.n_components()));
     fe_values.get_function_values(fe_function, u_value);
 
     value = u_value[0];
@@ -6842,7 +6865,9 @@ namespace VectorTools
 
     // then use this to get the gradients of
     // the given fe_function at this point
-    std::vector<std::vector<Tensor<1, dim> > > u_gradient(1, std::vector<Tensor<1, dim> > (fe.n_components()));
+    typedef typename InVector::value_type Number;
+    std::vector<std::vector<Tensor<1, dim, Number> > >
+      u_gradient(1, std::vector<Tensor<1, dim, Number> > (fe.n_components()));
     fe_values.get_function_gradients(fe_function, u_gradient);
 
     gradient = u_gradient[0];
@@ -6857,6 +6882,7 @@ namespace VectorTools
                   const Point<spacedim>      &point,
                   std::vector<Tensor<1, spacedim> > &gradient)
   {
+    typedef typename InVector::value_type Number;
     const hp::FECollection<dim, spacedim> &fe = dof.get_fe();
 
     Assert(gradient.size() == fe.n_components(),
@@ -6881,7 +6907,9 @@ namespace VectorTools
 
     // then use this to get the gradients of
     // the given fe_function at this point
-    std::vector<std::vector<Tensor<1, dim> > > u_gradient(1, std::vector<Tensor<1, dim> > (fe.n_components()));
+    typedef typename InVector::value_type Number;
+    std::vector<std::vector<Tensor<1, dim, Number> > >
+      u_gradient(1, std::vector<Tensor<1, dim, Number> > (fe.n_components()));
     fe_values.get_function_gradients(fe_function, u_gradient);
 
     gradient = u_gradient[0];
@@ -6980,6 +7008,7 @@ namespace VectorTools
                       const InVector        &v,
                       const unsigned int     component)
   {
+    typedef typename InVector::value_type Number;
     Assert (v.size() == dof.n_dofs(),
             ExcDimensionMismatch (v.size(), dof.n_dofs()));
     Assert (component < dof.get_fe().n_components(),
@@ -6990,10 +7019,10 @@ namespace VectorTools
                                           | update_values));
 
     typename DoFHandler<dim,spacedim>::active_cell_iterator cell;
-    std::vector<Vector<double> > values(quadrature.size(),
-                                        Vector<double> (dof.get_fe().n_components()));
+    std::vector<Vector<Number> > values(quadrature.size(),
+                                        Vector<Number> (dof.get_fe().n_components()));
 
-    double mean = 0.;
+    Number mean = Number();
     double area = 0.;
     // Compute mean value
     for (cell = dof.begin_active(); cell != dof.end(); ++cell)
@@ -7015,7 +7044,8 @@ namespace VectorTools
         p_d_triangulation
         = dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim> *>(&dof.get_tria()))
       {
-        double my_values[2] = { mean, area };
+        double mean_double = mean;
+        double my_values[2] = { mean_double, area };
         double global_values[2];
 
         MPI_Allreduce (&my_values, &global_values, 2, MPI_DOUBLE,

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.