From: Daniel Arndt Date: Mon, 4 Sep 2017 23:27:23 +0000 (+0200) Subject: Variadic template versions of FETools::Compositing functions X-Git-Tag: v9.0.0-rc1~1079^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7e4267be82e2a28414215a4b8ab939ede96359bc;p=dealii.git Variadic template versions of FETools::Compositing functions --- diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 8e8228eda1..91e28b0812 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -976,6 +976,18 @@ namespace FETools const std::vector &multiplicities, const bool do_tensor_product = true); + /** + * Same as above for an arbitrary number of parameters of type + * std::pair>, unsigned int> + * and do_tensor_product = true. + */ + template ::type..., + std::pair>, + unsigned int>>::value>::type> + FiniteElementData + multiply_dof_numbers (FESystems &&... fe_systems); + /** * Same as above but for a specific number of sub-elements. */ @@ -1009,6 +1021,17 @@ namespace FETools compute_restriction_is_additive_flags (const std::vector*> &fes, const std::vector &multiplicities); + /** + * Same as above for an arbitrary number of parameters of type + * std::pair>, unsigned int>. + */ + template ::type..., + std::pair>, + unsigned int>>::value>::type> + std::vector + compute_restriction_is_additive_flags (FESystems &&... fe_systems); + /** * Take a @p FiniteElement object and return a boolean vector * describing the @p restriction_is_additive_flags (see the @@ -1036,6 +1059,7 @@ namespace FETools const FiniteElement *fe5=nullptr, const unsigned int N5=0); + /** * Compute the nonzero components for each shape function of a * composed finite element described by a list of finite elements @@ -1058,6 +1082,18 @@ namespace FETools const std::vector &multiplicities, const bool do_tensor_product = true); + /** + * Same as above for an arbitrary number of parameters of type + * std::pair>, unsigned int> + * and do_tensor_product = true. + */ + template ::type..., + std::pair>, + unsigned int>>::value>::type> + std::vector + compute_nonzero_components (FESystems &&... fe_systems); + /** * Compute the non-zero vector components of a composed finite * element. This function is similar to the previous one, except @@ -1343,7 +1379,80 @@ namespace FETools { return new FE(degree); } -} + + namespace Compositing + { + template + std::vector + compute_restriction_is_additive_flags (FESystems &&... fe_systems) + { + std::vector*> fes; + std::vector multiplicities; + + const auto extract = + [&fes, &multiplicities] + (const std::pair>, unsigned int> &fe_system) + { + fes.push_back(fe_system.first.get()); + multiplicities.push_back(fe_system.second); + }; + + const auto fe_system_pointers = {&fe_systems...}; + for (const auto p : fe_system_pointers) + extract (*p); + + return compute_restriction_is_additive_flags (fes, multiplicities); + } + + + + template + FiniteElementData + multiply_dof_numbers (FESystems &&... fe_systems) + { + std::vector*> fes; + std::vector multiplicities; + + const auto extract = + [&fes, &multiplicities] + (const std::pair>, unsigned int> &fe_system) + { + fes.push_back(fe_system.first.get()); + multiplicities.push_back(fe_system.second); + }; + + const auto fe_system_pointers = {&fe_systems...}; + for (const auto p : fe_system_pointers) + extract (*p); + + return multiply_dof_numbers (fes, multiplicities, true); + } + + + + template + std::vector + compute_nonzero_components (FESystems &&... fe_systems) + { + std::vector*> fes; + std::vector multiplicities; + + const auto extract = + [&fes, &multiplicities] + (const std::pair>, unsigned int> &fe_system) + { + fes.push_back(fe_system.first.get()); + multiplicities.push_back(fe_system.second); + }; + + const auto fe_system_pointers = {&fe_systems...}; + for (const auto p : fe_system_pointers) + extract (*p); + + return compute_nonzero_components (fes, multiplicities, true); + } + } //Compositing +} //FETools #endif