From: Denis Davydov Date: Sun, 1 Mar 2015 15:37:59 +0000 (+0100) Subject: temporary fix of vector_tools.templates.h X-Git-Tag: v8.3.0-rc1~184^2~18 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b165e17b98cf2c151f3ae747899974db7f341c7a;p=dealii.git temporary fix of vector_tools.templates.h --- diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index b71d1ac809..0a2bed01e0 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -6344,6 +6344,7 @@ namespace VectorTools const Function *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 fe_collection (dof.get_fe()); IDScratchData data(mapping, fe_collection, q, update_flags); + //FIXME + // temporary vectors of consistent with InVector type + std::vector> function_values; + std::vector >> 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(n_components)); + function_grads.resize (n_q_points, + std::vector >(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 &difference, const Point &point) { + typedef typename InVector::value_type Number; const FiniteElement &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 > u_value(1, Vector (fe.n_components())); + std::vector > u_value(1, Vector (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 &point, Vector &value) { + typedef typename InVector::value_type Number; const FiniteElement &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 > u_value(1, Vector (fe.n_components())); + std::vector > u_value(1, Vector (fe.n_components())); fe_values.get_function_values(fe_function, u_value); value = u_value[0]; @@ -6686,6 +6708,7 @@ namespace VectorTools const Point &point, Vector &value) { + typedef typename InVector::value_type Number; const hp::FECollection &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 > u_value(1, Vector (fe.n_components())); + std::vector > u_value(1, Vector (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 > > u_gradient(1, std::vector > (fe.n_components())); + typedef typename InVector::value_type Number; + std::vector > > + u_gradient(1, std::vector > (fe.n_components())); fe_values.get_function_gradients(fe_function, u_gradient); gradient = u_gradient[0]; @@ -6857,6 +6882,7 @@ namespace VectorTools const Point &point, std::vector > &gradient) { + typedef typename InVector::value_type Number; const hp::FECollection &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 > > u_gradient(1, std::vector > (fe.n_components())); + typedef typename InVector::value_type Number; + std::vector > > + u_gradient(1, std::vector > (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::active_cell_iterator cell; - std::vector > values(quadrature.size(), - Vector (dof.get_fe().n_components())); + std::vector > values(quadrature.size(), + Vector (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 *>(&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,