From: Matthias Maier Date: Thu, 7 Sep 2017 22:14:48 +0000 (-0500) Subject: Add FETools::convert_generalized_support_point_values_to_dof_values X-Git-Tag: v9.0.0-rc1~1055^2~9 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=00314125da1e30f8eebe27512eda86af3685479a;p=dealii.git Add FETools::convert_generalized_support_point_values_to_dof_values --- diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 202fec1088..1d28e058bc 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -557,13 +557,31 @@ namespace FETools * information. */ template - void - compute_projection_from_face_quadrature_points_matrix (const FiniteElement &fe, - const Quadrature &lhs_quadrature, - const Quadrature &rhs_quadrature, - const typename DoFHandler::active_cell_iterator &cell, - const unsigned int face, - FullMatrix &X); + void compute_projection_from_face_quadrature_points_matrix( + const FiniteElement &fe, + const Quadrature &lhs_quadrature, + const Quadrature &rhs_quadrature, + const typename DoFHandler::active_cell_iterator &cell, + const unsigned int face, + FullMatrix &X); + + + + /** + * Wrapper around + * FiniteElement::convert_generalized_support_point_values_to_dof_values() + * that works with arbitrary number types. + * + * If the number type is double this function simply calls the + * FiniteElement::convert_generalized_support_point_values_to_dof_values() + * directly. Otherwise, temporary (per-thread storage) is used to call to + * above function. + */ + template + void convert_generalized_support_point_values_to_dof_values( + const FiniteElement &finite_element, + const std::vector > &support_point_values, + std::vector &nodal_values); diff --git a/include/deal.II/fe/fe_tools.templates.h b/include/deal.II/fe/fe_tools.templates.h index b10f0efd72..cd064e4b85 100644 --- a/include/deal.II/fe/fe_tools.templates.h +++ b/include/deal.II/fe/fe_tools.templates.h @@ -2910,6 +2910,79 @@ namespace FETools + namespace + { + // internal helper class to partially specialize the number type + + template + struct ConvertHelper + { + static void convert_generalized_support_point_values_to_dof_values( + const FiniteElement &finite_element, + const std::vector > &support_point_values, + std::vector &dof_values) + { + static Threads::ThreadLocalStorage > > + double_support_point_values; + static Threads::ThreadLocalStorage > + double_dof_values; + + double_support_point_values.get().resize(support_point_values.size()); + double_dof_values.get().resize(dof_values.size()); + + for (unsigned int i = 0; i < support_point_values.size(); ++i) + { + double_support_point_values.get()[i].reinit( + finite_element.n_components(), false); + std::copy(std::begin(support_point_values[i]), + std::end(support_point_values[i]), + std::begin(double_support_point_values.get()[i])); + } + + finite_element.convert_generalized_support_point_values_to_dof_values( + double_support_point_values.get(), double_dof_values.get()); + + std::copy(std::begin(double_dof_values.get()), + std::end(double_dof_values.get()), + std::begin(dof_values)); + } + }; + + // In case of double simply call the finite element directly + template + struct ConvertHelper + { + static void convert_generalized_support_point_values_to_dof_values( + const FiniteElement &finite_element, + const std::vector > &support_point_values, + std::vector &dof_values) + { + finite_element.convert_generalized_support_point_values_to_dof_values( + support_point_values, dof_values); + } + }; + + } /* anonymous namespace */ + + + + template + void convert_generalized_support_point_values_to_dof_values( + const FiniteElement &finite_element, + const std::vector > &support_point_values, + std::vector &dof_values) + { + AssertDimension(support_point_values.size(), + finite_element.get_generalized_support_points().size()); + AssertDimension(dof_values.size(), finite_element.dofs_per_cell); + + ConvertHelper:: + convert_generalized_support_point_values_to_dof_values( + finite_element, support_point_values, dof_values); + } + + + template void hierarchic_to_lexicographic_numbering (unsigned int degree, std::vector &h2l) diff --git a/source/fe/fe_tools.inst.in b/source/fe/fe_tools.inst.in index 53815d8f45..53b7632fec 100644 --- a/source/fe/fe_tools.inst.in +++ b/source/fe/fe_tools.inst.in @@ -258,3 +258,21 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS #endif \} } + + +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; number : REAL_SCALARS) +{ +#if deal_II_dimension <= deal_II_space_dimension + namespace FETools + \{ + + template + void + convert_generalized_support_point_values_to_dof_values ( + const FiniteElement &, + const std::vector > &, + std::vector &); + + \} +#endif +}