From: Daniel Arndt Date: Mon, 4 Sep 2017 23:31:55 +0000 (+0200) Subject: Implement variadic template constructor and std::initializer_list constructor for... X-Git-Tag: v9.0.0-rc1~1079^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bf3a03c3b7611495ebcbc5c6c9cf5d25af4f1baa;p=dealii.git Implement variadic template constructor and std::initializer_list constructor for FESystem --- diff --git a/include/deal.II/fe/fe_system.h b/include/deal.II/fe/fe_system.h index e439fb71b1..30254613f8 100644 --- a/include/deal.II/fe/fe_system.h +++ b/include/deal.II/fe/fe_system.h @@ -23,6 +23,7 @@ #include #include #include +#include #include #include #include @@ -421,6 +422,48 @@ public: FESystem (const std::vector*> &fes, const std::vector &multiplicities); + /** + * Constructor taking an arbitrary number of parameters of type + * std::pair>, unsigned int>. + * In combination with FiniteElement::operator^, this allows to construct FESystem objects + * as follows: + * @code + * FiniteElementType1 fe_1; + * FiniteElementType1 fe_2; + * FESystem fe_system = ( fe_1^dim, fe_2^1 ); + * @endcode + * + * The FiniteElement objects are not actually used for anything other than creating a + * copy that will then be owned by the current object. In other words, it is + * completely fine to call this constructor with a temporary object for the + * finite element, as in this code snippet: + * @code + * FESystem fe (FE_Q(2)^2); + * @endcode + * Here, FE_Q@(2) constructs an unnamed, temporary object + * that is passed to the FESystem constructor to create a finite element + * that consists of two components, both of which are quadratic FE_Q + * elements. The temporary is destroyed again at the end of the code that + * corresponds to this line, but this does not matter because FESystem + * creates its own copy of the FE_Q object. + */ + template ::type..., + std::pair>, + unsigned int>>::value>::type> + FESystem (FESystems&& ... fe_systems); + + /** + * Same as above allowing the following syntax: + * @code + * FiniteElementType1 fe_1; + * FiniteElementType1 fe_2; + * FESystem fe_system = { fe_1^dim, fe_2^1 }; + * @endcode + */ + FESystem (const std::initializer_list>, + unsigned int>> &fe_systems); + /** * Copy constructor. This constructor is deleted, i.e., copying * FESystem objects is not allowed. @@ -430,7 +473,7 @@ public: /** * Destructor. */ - virtual ~FESystem (); + virtual ~FESystem () = default; /** * Return a string that uniquely identifies a finite element. This element @@ -1081,9 +1124,66 @@ private: friend class FE_Enriched; }; +//------------------------variadic template constructor------------------------ + +#ifndef DOXYGEN +namespace +{ + template + unsigned int count_nonzeros + (const std::initializer_list>, + unsigned int>> &fe_systems) + { + return std::count_if(fe_systems.begin(), + fe_systems.end(), + [](const std::pair>, + unsigned int> &fe_system) + { + return fe_system.second > 0; + }); + } +} + +// We are just forwarding/delegating to the constructor taking a std::initializer_list. +// If we decide to remove the deprecated constructors, we might just use the variadic +// constructor with a suitable static_assert instead of the std::enable_if. +template +template +FESystem::FESystem (FESystems &&... fe_systems) + : FESystem ({std::forward(fe_systems)...}) +{} + +template +FESystem::FESystem +(const std::initializer_list>, + unsigned int>> &fe_systems) + : + FiniteElement + (FETools::Compositing::multiply_dof_numbers(fe_systems), + FETools::Compositing::compute_restriction_is_additive_flags (fe_systems), + FETools::Compositing::compute_nonzero_components(fe_systems)), + base_elements(count_nonzeros(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); + }; + + for (const auto &p : fe_systems) + extract (p); + + initialize(fes, multiplicities); +} + +#endif //DOXYGEN DEAL_II_NAMESPACE_CLOSE -/*---------------------------- fe_system.h ---------------------------*/ #endif /*---------------------------- fe_system.h ---------------------------*/ diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 91e28b0812..c1f7312c91 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -981,12 +981,10 @@ namespace FETools * std::pair>, unsigned int> * and do_tensor_product = true. */ - template ::type..., - std::pair>, - unsigned int>>::value>::type> + template FiniteElementData - multiply_dof_numbers (FESystems &&... fe_systems); + multiply_dof_numbers + (const std::initializer_list>, unsigned int>> &fe_systems); /** * Same as above but for a specific number of sub-elements. @@ -1025,12 +1023,10 @@ namespace FETools * Same as above for an arbitrary number of parameters of type * std::pair>, unsigned int>. */ - template ::type..., - std::pair>, - unsigned int>>::value>::type> + template std::vector - compute_restriction_is_additive_flags (FESystems &&... fe_systems); + compute_restriction_is_additive_flags + (const std::initializer_list>, unsigned int>> &fe_systems); /** * Take a @p FiniteElement object and return a boolean vector @@ -1087,12 +1083,10 @@ namespace FETools * std::pair>, unsigned int> * and do_tensor_product = true. */ - template ::type..., - std::pair>, - unsigned int>>::value>::type> + template std::vector - compute_nonzero_components (FESystems &&... fe_systems); + compute_nonzero_components + (const std::initializer_list>, unsigned int>> &fe_systems); /** * Compute the non-zero vector components of a composed finite @@ -1382,9 +1376,10 @@ namespace FETools namespace Compositing { - template + template std::vector - compute_restriction_is_additive_flags (FESystems &&... fe_systems) + compute_restriction_is_additive_flags + (const std::initializer_list>, unsigned int>> &fe_systems) { std::vector*> fes; std::vector multiplicities; @@ -1397,18 +1392,18 @@ namespace FETools multiplicities.push_back(fe_system.second); }; - const auto fe_system_pointers = {&fe_systems...}; - for (const auto p : fe_system_pointers) - extract (*p); + for (const auto &p : fe_systems) + extract (p); return compute_restriction_is_additive_flags (fes, multiplicities); } - template + template FiniteElementData - multiply_dof_numbers (FESystems &&... fe_systems) + multiply_dof_numbers + (const std::initializer_list>, unsigned int>> &fe_systems) { std::vector*> fes; std::vector multiplicities; @@ -1421,18 +1416,18 @@ namespace FETools multiplicities.push_back(fe_system.second); }; - const auto fe_system_pointers = {&fe_systems...}; - for (const auto p : fe_system_pointers) - extract (*p); + for (const auto &p : fe_systems) + extract (p); return multiply_dof_numbers (fes, multiplicities, true); } - template + template std::vector - compute_nonzero_components (FESystems &&... fe_systems) + compute_nonzero_components + (const std::initializer_list>, unsigned int>> &fe_systems) { std::vector*> fes; std::vector multiplicities; @@ -1445,9 +1440,8 @@ namespace FETools multiplicities.push_back(fe_system.second); }; - const auto fe_system_pointers = {&fe_systems...}; - for (const auto p : fe_system_pointers) - extract (*p); + for (const auto &p : fe_systems) + extract (p); return compute_nonzero_components (fes, multiplicities, true); } diff --git a/source/fe/fe_system.cc b/source/fe/fe_system.cc index 5412942c81..97d58fa604 100644 --- a/source/fe/fe_system.cc +++ b/source/fe/fe_system.cc @@ -285,12 +285,6 @@ FESystem::FESystem ( -template -FESystem::~FESystem () -{} - - - template std::string FESystem::get_name () const