From 4a3ff764789b05b1ae009b87dae051663cc0c92f Mon Sep 17 00:00:00 2001 From: David Wells Date: Sun, 18 Oct 2015 22:48:08 -0400 Subject: [PATCH] Use VectorType, not InVector, as a template type. Some functions use a template for the input and output types, but other functions use InVector. This commit switches the second kind to use VectorType, which is more consistent with the usage of both VectorType and InVector. --- include/deal.II/numerics/vector_tools.h | 118 +++++----- .../deal.II/numerics/vector_tools.templates.h | 206 +++++++++--------- 2 files changed, 162 insertions(+), 162 deletions(-) diff --git a/include/deal.II/numerics/vector_tools.h b/include/deal.II/numerics/vector_tools.h index 0c893bfdc1..1e02f20c98 100644 --- a/include/deal.II/numerics/vector_tools.h +++ b/include/deal.II/numerics/vector_tools.h @@ -1963,12 +1963,12 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template - void point_difference (const DoFHandler &dof, - const InVector &fe_function, - const Function &exact_solution, - Vector &difference, - const Point &point); + template + void point_difference (const DoFHandler &dof, + const VectorType &fe_function, + const Function &exact_solution, + Vector &difference, + const Point &point); /** * Point error evaluation. Find the first cell containing the given point @@ -1982,13 +1982,13 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template - void point_difference (const Mapping &mapping, - const DoFHandler &dof, - const InVector &fe_function, - const Function &exact_solution, - Vector &difference, - const Point &point); + template + void point_difference (const Mapping &mapping, + const DoFHandler &dof, + const VectorType &fe_function, + const Function &exact_solution, + Vector &difference, + const Point &point); /** * Evaluate a possibly vector-valued finite element function defined by the @@ -2001,10 +2001,10 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_value (const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, Vector &value); @@ -2014,10 +2014,10 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_value (const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, Vector &value); @@ -2036,10 +2036,10 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template double point_value (const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); /** @@ -2048,10 +2048,10 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template double point_value (const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); /** @@ -2065,11 +2065,11 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_value (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, Vector &value); @@ -2079,11 +2079,11 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_value (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, Vector &value); @@ -2098,11 +2098,11 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template double point_value (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); /** @@ -2111,11 +2111,11 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template double point_value (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); /** @@ -2129,12 +2129,12 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_gradient (const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, - std::vector > &value); + std::vector > &value); /** * Same as above for hp. @@ -2142,12 +2142,12 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_gradient (const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, - std::vector > &value); + std::vector > &value); /** * Evaluate a scalar finite element function defined by the given DoFHandler @@ -2160,10 +2160,10 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template - Tensor<1, spacedim, typename InVector::value_type> + template + Tensor<1, spacedim, typename VectorType::value_type> point_gradient (const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); /** @@ -2172,10 +2172,10 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template - Tensor<1, spacedim, typename InVector::value_type> + template + Tensor<1, spacedim, typename VectorType::value_type> point_gradient (const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); /** @@ -2189,13 +2189,13 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_gradient (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, - std::vector > &value); + std::vector > &value); /** * Same as above for hp. @@ -2203,13 +2203,13 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template + template void point_gradient (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point, - std::vector > &value); + std::vector > &value); /** * Evaluate a scalar finite element function defined by the given DoFHandler @@ -2222,11 +2222,11 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template - Tensor<1, spacedim, typename InVector::value_type> + template + Tensor<1, spacedim, typename VectorType::value_type> point_gradient (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); /** @@ -2235,11 +2235,11 @@ namespace VectorTools * @note If the cell in which the point is found is not locally owned, an * exception of type VectorTools::ExcPointNotAvailableHere is thrown. */ - template - Tensor<1, spacedim, typename InVector::value_type> + template + Tensor<1, spacedim, typename VectorType::value_type> point_gradient (const hp::MappingCollection &mapping, const hp::DoFHandler &dof, - const InVector &fe_function, + const VectorType &fe_function, const Point &point); //@} @@ -2321,22 +2321,22 @@ namespace VectorTools * Lagrangian elements. For all other elements, you will need to compute the * mean value and subtract it right inside the evaluation routine. */ - template + template double compute_mean_value (const Mapping &mapping, const DoFHandler &dof, const Quadrature &quadrature, - const InVector &v, - const unsigned int component); + const VectorType &v, + const unsigned int component); /** * Calls the other compute_mean_value() function, see above, with * mapping=MappingQGeneric@(1). */ - template + template double compute_mean_value (const DoFHandler &dof, const Quadrature &quadrature, - const InVector &v, - const unsigned int component); + const VectorType &v, + const unsigned int component); //@} /** * Geometrical interpolation diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index eabe985aa0..2a0a2d6643 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -6541,13 +6541,13 @@ namespace VectorTools - template + template void point_difference (const DoFHandler &dof, - const InVector &fe_function, - const Function &exact_function, - Vector &difference, - const Point &point) + const VectorType &fe_function, + const Function &exact_function, + Vector &difference, + const Point &point) { point_difference(StaticMappingQ1::mapping, dof, @@ -6558,16 +6558,16 @@ namespace VectorTools } - template + template void - point_difference (const Mapping &mapping, + point_difference (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, - const Function &exact_function, - Vector &difference, - const Point &point) + const VectorType &fe_function, + const Function &exact_function, + Vector &difference, + const Point &point) { - typedef typename InVector::value_type Number; + typedef typename VectorType::value_type Number; const FiniteElement &fe = dof.get_fe(); Assert(difference.size() == fe.n_components(), @@ -6604,12 +6604,12 @@ namespace VectorTools } - template + template void point_value (const DoFHandler &dof, - const InVector &fe_function, - const Point &point, - Vector &value) + const VectorType &fe_function, + const Point &point, + Vector &value) { point_value (StaticMappingQ1::mapping, @@ -6620,12 +6620,12 @@ namespace VectorTools } - template + template void point_value (const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point, - Vector &value) + const VectorType &fe_function, + const Point &point, + Vector &value) { point_value(hp::StaticMappingQ1::mapping_collection, dof, @@ -6635,11 +6635,11 @@ namespace VectorTools } - template + template double point_value (const DoFHandler &dof, - const InVector &fe_function, - const Point &point) + const VectorType &fe_function, + const Point &point) { return point_value (StaticMappingQ1::mapping, dof, @@ -6648,11 +6648,11 @@ namespace VectorTools } - template + template double point_value (const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point) + const VectorType &fe_function, + const Point &point) { return point_value(hp::StaticMappingQ1::mapping_collection, dof, @@ -6661,15 +6661,15 @@ namespace VectorTools } - template + template void - point_value (const Mapping &mapping, + point_value (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, - const Point &point, - Vector &value) + const VectorType &fe_function, + const Point &point, + Vector &value) { - typedef typename InVector::value_type Number; + typedef typename VectorType::value_type Number; const FiniteElement &fe = dof.get_fe(); Assert(value.size() == fe.n_components(), @@ -6702,15 +6702,15 @@ namespace VectorTools } - template + template void - point_value (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point, - Vector &value) + point_value (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const VectorType &fe_function, + const Point &point, + Vector &value) { - typedef typename InVector::value_type Number; + typedef typename VectorType::value_type Number; const hp::FECollection &fe = dof.get_fe(); Assert(value.size() == fe.n_components(), @@ -6742,12 +6742,12 @@ namespace VectorTools } - template + template double - point_value (const Mapping &mapping, + point_value (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, - const Point &point) + const VectorType &fe_function, + const Point &point) { Assert(dof.get_fe().n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); @@ -6759,12 +6759,12 @@ namespace VectorTools } - template + template double - point_value (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point) + point_value (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const VectorType &fe_function, + const Point &point) { Assert(dof.get_fe().n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); @@ -6777,12 +6777,12 @@ namespace VectorTools - template + template void - point_gradient (const DoFHandler &dof, - const InVector &fe_function, - const Point &point, - std::vector > &gradients) + point_gradient (const DoFHandler &dof, + const VectorType &fe_function, + const Point &point, + std::vector > &gradients) { point_gradient (StaticMappingQ1::mapping, @@ -6793,12 +6793,12 @@ namespace VectorTools } - template + template void point_gradient (const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point, - std::vector > &gradients) + const VectorType &fe_function, + const Point &point, + std::vector > &gradients) { point_gradient(hp::StaticMappingQ1::mapping_collection, dof, @@ -6808,11 +6808,11 @@ namespace VectorTools } - template - Tensor<1, spacedim, typename InVector::value_type> + template + Tensor<1, spacedim, typename VectorType::value_type> point_gradient (const DoFHandler &dof, - const InVector &fe_function, - const Point &point) + const VectorType &fe_function, + const Point &point) { return point_gradient (StaticMappingQ1::mapping, dof, @@ -6821,11 +6821,11 @@ namespace VectorTools } - template - Tensor<1, spacedim, typename InVector::value_type> + template + Tensor<1, spacedim, typename VectorType::value_type> point_gradient (const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point) + const VectorType &fe_function, + const Point &point) { return point_gradient(hp::StaticMappingQ1::mapping_collection, dof, @@ -6834,13 +6834,13 @@ namespace VectorTools } - template + template void - point_gradient (const Mapping &mapping, + point_gradient (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, - const Point &point, - std::vector > &gradient) + const VectorType &fe_function, + const Point &point, + std::vector > &gradient) { const FiniteElement &fe = dof.get_fe(); @@ -6867,7 +6867,7 @@ namespace VectorTools // then use this to get the gradients of // the given fe_function at this point - typedef typename InVector::value_type Number; + typedef typename VectorType::value_type Number; std::vector > > u_gradient(1, std::vector > (fe.n_components())); fe_values.get_function_gradients(fe_function, u_gradient); @@ -6876,15 +6876,15 @@ namespace VectorTools } - template + template void - point_gradient (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point, - std::vector > &gradient) + point_gradient (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const VectorType &fe_function, + const Point &point, + std::vector > &gradient) { - typedef typename InVector::value_type Number; + typedef typename VectorType::value_type Number; const hp::FECollection &fe = dof.get_fe(); Assert(gradient.size() == fe.n_components(), @@ -6915,17 +6915,17 @@ namespace VectorTools } - template - Tensor<1, spacedim, typename InVector::value_type> - point_gradient (const Mapping &mapping, + template + Tensor<1, spacedim, typename VectorType::value_type> + point_gradient (const Mapping &mapping, const DoFHandler &dof, - const InVector &fe_function, - const Point &point) + const VectorType &fe_function, + const Point &point) { Assert(dof.get_fe().n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); - std::vector > gradient(1); + std::vector > gradient(1); point_gradient (mapping, dof, fe_function, point, gradient); return gradient[0]; @@ -6933,17 +6933,17 @@ namespace VectorTools - template - Tensor<1, spacedim, typename InVector::value_type> - point_gradient (const hp::MappingCollection &mapping, - const hp::DoFHandler &dof, - const InVector &fe_function, - const Point &point) + template + Tensor<1, spacedim, typename VectorType::value_type> + point_gradient (const hp::MappingCollection &mapping, + const hp::DoFHandler &dof, + const VectorType &fe_function, + const Point &point) { Assert(dof.get_fe().n_components() == 1, ExcMessage ("Finite element is not scalar as is necessary for this function")); - std::vector > gradient(1); + std::vector > gradient(1); point_gradient (mapping, dof, fe_function, point, gradient); return gradient[0]; @@ -7020,15 +7020,15 @@ namespace VectorTools } - template + template double - compute_mean_value (const Mapping &mapping, + compute_mean_value (const Mapping &mapping, const DoFHandler &dof, - const Quadrature &quadrature, - const InVector &v, - const unsigned int component) + const Quadrature &quadrature, + const VectorType &v, + const unsigned int component) { - typedef typename InVector::value_type Number; + typedef typename VectorType::value_type Number; Assert (v.size() == dof.n_dofs(), ExcDimensionMismatch (v.size(), dof.n_dofs())); Assert (component < dof.get_fe().n_components(), @@ -7087,20 +7087,20 @@ namespace VectorTools } - template + template double compute_mean_value (const DoFHandler &dof, - const Quadrature &quadrature, - const InVector &v, - const unsigned int component) + const Quadrature &quadrature, + const VectorType &v, + const unsigned int component) { return compute_mean_value(StaticMappingQ1::mapping, dof, quadrature, v, component); } template - void get_position_vector(const DH &dh, - VECTOR &vector, + void get_position_vector(const DH &dh, + VECTOR &vector, const ComponentMask &mask) { AssertDimension(vector.size(), dh.n_dofs()); -- 2.39.5