From 0a2e293538f9f96ff25fde18c51be8a9625578fc Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Tue, 8 Aug 2017 18:36:17 -0600 Subject: [PATCH] Added Scalar::get_function_values_from_local_dof_values --- include/deal.II/fe/fe_values.h | 62 ++++++++++++++++++++++++++++++++++ source/fe/fe_values.cc | 22 +++++++++++- source/fe/fe_values.inst.in | 13 +++++++ 3 files changed, 96 insertions(+), 1 deletion(-) diff --git a/include/deal.II/fe/fe_values.h b/include/deal.II/fe/fe_values.h index d49cf35727..3e747687d2 100644 --- a/include/deal.II/fe/fe_values.h +++ b/include/deal.II/fe/fe_values.h @@ -165,6 +165,44 @@ namespace FEValuesViews */ typedef dealii::Tensor<3,spacedim> third_derivative_type; + /** + * A struct that provides the output type for the product of the value + * and derivatives of basis functions of the Scalar view and any @p Number type. + */ + template + struct OutputType + { + /** + * A typedef for the data type of the product of a @p Number and the + * values of the view the Scalar class. + */ + typedef typename ProductType::type, typename Scalar::value_type>::type value_type; + + /** + * A typedef for the data type of the product of a @p Number and the + * gradients of the view the Scalar class. + */ + typedef typename ProductType::type, typename Scalar::gradient_type>::type gradient_type; + + /** + * A typedef for the data type of the product of a @p Number and the + * laplacians of the view the Scalar class. + */ + typedef typename ProductType::type, typename Scalar::value_type>::type laplacian_type; + + /** + * A typedef for the data type of the product of a @p Number and the + * hessians of the view the Scalar class. + */ + typedef typename ProductType::type, typename Scalar::hessian_type>::type hessian_type; + + /** + * A typedef for the data type of the product of a @p Number and the + * third derivatives of the view the Scalar class. + */ + typedef typename ProductType::type, typename Scalar::third_derivative_type>::type third_derivative_type; + }; + /** * A structure where for each shape function we pre-compute a bunch of * data that will make later accesses much cheaper. @@ -291,6 +329,30 @@ namespace FEValuesViews void get_function_values (const InputVector &fe_function, std::vector::type> &values) const; + /** + * Same as above, but using a vector of local degree-of-freedom values. + * + * The @p dof_values vector must have a length equal to number of DoFs on + * a cell, and each entry @p dof_values[i] is the value of the local DoF + * @p i. The fundamental prerequisite for the @p InputVector is that it must + * be possible to create an ArrayView from it; this is satisfied by the + * @p std::vector class. + * + * The DoF values typically would be obtained in the following way: + * @code + * Vector local_dof_values(cell->get_fe().dofs_per_cell); + * cell->get_dof_values(solution, local_dof_values); + * @endcode + * or, for a generic @p Number type, + * @code + * std::vector local_dof_values(cell->get_fe().dofs_per_cell); + * cell->get_dof_values(solution, local_dof_values.begin(), local_dof_values.end()); + * @endcode + */ + template + void get_function_values_from_local_dof_values (const InputVector &dof_values, + std::vector::value_type> &values) const; + /** * Return the gradients of the selected scalar component of the finite * element function characterized by fe_function at the diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc index 2d1ab1755f..039a7cec74 100644 --- a/source/fe/fe_values.cc +++ b/source/fe/fe_values.cc @@ -39,6 +39,8 @@ #include #include +#include + DEAL_II_NAMESPACE_OPEN @@ -459,7 +461,7 @@ namespace FEValuesViews do_function_values (const ArrayView &dof_values, const Table<2,double> &shape_values, const std::vector::ShapeFunctionData> &shape_function_data, - std::vector::type> &values) + std::vector::type,double>::type> &values) { const unsigned int dofs_per_cell = dof_values.size(); const unsigned int n_quadrature_points = dofs_per_cell > 0 ? @@ -1286,6 +1288,24 @@ namespace FEValuesViews fe_values->finite_element_output.shape_values, shape_function_data, values); } + template + template + void + Scalar:: + get_function_values_from_local_dof_values (const InputVector &dof_values, + std::vector::value_type> &values) const + { + Assert (fe_values->update_flags & update_values, + (typename FEValuesBase::ExcAccessToUninitializedField("update_values"))); + Assert (fe_values->present_cell.get() != nullptr, + ExcMessage ("FEValues object is not reinit'ed to any cell")); + AssertDimension (dof_values.size(), fe_values->dofs_per_cell); + + internal::do_function_values + (make_array_view(dof_values.begin(), dof_values.end()), + fe_values->finite_element_output.shape_values, shape_function_data, values); + } + template diff --git a/source/fe/fe_values.inst.in b/source/fe/fe_values.inst.in index 4526591997..90f24d926a 100644 --- a/source/fe/fe_values.inst.in +++ b/source/fe/fe_values.inst.in @@ -139,6 +139,19 @@ for (VEC : SERIAL_VECTORS; deal_II_dimension : DIMENSIONS; deal_II_space_dimensi } +for (VEC : GENERAL_CONTAINER_TYPES; + Number : ALL_SCALAR_TYPES; + deal_II_dimension : DIMENSIONS; + deal_II_space_dimension : SPACE_DIMENSIONS) +{ +#if deal_II_dimension <= deal_II_space_dimension + template + void FEValuesViews::Scalar + ::get_function_values_from_local_dof_values> + (const VEC&, std::vector::value_type>&) const; +#endif +} + for (VEC : SERIAL_VECTORS; deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { -- 2.39.5