]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement variadic template constructor and std::initializer_list constructor for...
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 4 Sep 2017 23:31:55 +0000 (01:31 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 11 Sep 2017 16:21:06 +0000 (18:21 +0200)
include/deal.II/fe/fe_system.h
include/deal.II/fe/fe_tools.h
source/fe/fe_system.cc

index e439fb71b1ac00df02c8cba965a4272e92134a31..30254613f87d506c78946ea68a0d0b02e505f08d 100644 (file)
@@ -23,6 +23,7 @@
 #include <deal.II/base/config.h>
 #include <deal.II/base/thread_management.h>
 #include <deal.II/fe/fe.h>
+#include <deal.II/fe/fe_tools.h>
 #include <vector>
 #include <memory>
 #include <utility>
@@ -421,6 +422,48 @@ public:
   FESystem (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
             const std::vector<unsigned int>                   &multiplicities);
 
+  /**
+   * Constructor taking an arbitrary number of parameters of type
+   * <code>std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int></code>.
+   * In combination with FiniteElement::operator^, this allows to construct FESystem objects
+   * as follows:
+   * @code
+   *   FiniteElementType1<dim,spacedim> fe_1;
+   *   FiniteElementType1<dim,spacedim> fe_2;
+   *   FESystem<dim,spacedim> 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<dim> fe (FE_Q<dim>(2)^2);
+   * @endcode
+   * Here, <code>FE_Q@<dim@>(2)</code> 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 <class... FESystems,
+            typename = typename std::enable_if<all_same_as<typename std::decay<FESystems>::type...,
+                                                           std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
+                                                               unsigned int>>::value>::type>
+  FESystem (FESystems&& ... fe_systems);
+
+  /**
+   * Same as above allowing the following syntax:
+   * @code
+   *   FiniteElementType1<dim,spacedim> fe_1;
+   *   FiniteElementType1<dim,spacedim> fe_2;
+   *   FESystem<dim,spacedim> fe_system = { fe_1^dim, fe_2^1 };
+   * @endcode
+   */
+  FESystem (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
+            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<dim,spacedim>;
 };
 
+//------------------------variadic template constructor------------------------
+
+#ifndef DOXYGEN
+namespace
+{
+  template <int dim, int spacedim>
+  unsigned int count_nonzeros
+  (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
+   unsigned int>> &fe_systems)
+  {
+    return std::count_if(fe_systems.begin(),
+                         fe_systems.end(),
+                         [](const std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
+                            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<int dim, int spacedim>
+template <class... FESystems, typename>
+FESystem<dim,spacedim>::FESystem (FESystems &&... fe_systems)
+  : FESystem<dim, spacedim> ({std::forward<FESystems>(fe_systems)...})
+{}
+
+template <int dim, int spacedim>
+FESystem<dim,spacedim>::FESystem
+(const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
+ unsigned int>> &fe_systems)
+  :
+  FiniteElement<dim,spacedim>
+  (FETools::Compositing::multiply_dof_numbers<dim,spacedim>(fe_systems),
+   FETools::Compositing::compute_restriction_is_additive_flags<dim,spacedim> (fe_systems),
+   FETools::Compositing::compute_nonzero_components<dim,spacedim>(fe_systems)),
+  base_elements(count_nonzeros(fe_systems))
+{
+  std::vector<const FiniteElement<dim,spacedim>*> fes;
+  std::vector<unsigned int> multiplicities;
+
+  const auto extract
+    = [&fes, &multiplicities]
+      (const std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, 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  ---------------------------*/
index 91e28b0812550bd53235c7f86f81edd5ab725304..c1f7312c914cc0453c8498948d43a3455f4b4888 100644 (file)
@@ -981,12 +981,10 @@ namespace FETools
      * <code>std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int></code>
      * and <code>do_tensor_product = true</code>.
      */
-    template <int dim, int spacedim, class... FESystems,
-              typename = typename std::enable_if<all_same_as<typename std::decay<FESystems>::type...,
-                                                             std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
-                                                                 unsigned int>>::value>::type>
+    template <int dim, int spacedim>
     FiniteElementData<dim>
-    multiply_dof_numbers (FESystems &&... fe_systems);
+    multiply_dof_numbers
+    (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, 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
      * <code>std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int></code>.
      */
-    template <int dim, int spacedim, class... FESystems,
-              typename = typename std::enable_if<all_same_as<typename std::decay<FESystems>::type...,
-                                                             std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
-                                                                 unsigned int>>::value>::type>
+    template <int dim, int spacedim>
     std::vector<bool>
-    compute_restriction_is_additive_flags (FESystems &&... fe_systems);
+    compute_restriction_is_additive_flags
+    (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>> &fe_systems);
 
     /**
      * Take a @p FiniteElement object and return a boolean vector
@@ -1087,12 +1083,10 @@ namespace FETools
      * <code>std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int></code>
      * and <code>do_tensor_product = true</code>.
      */
-    template <int dim, int spacedim, class... FESystems,
-              typename = typename std::enable_if<all_same_as<typename std::decay<FESystems>::type...,
-                                                             std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>,
-                                                                 unsigned int>>::value>::type>
+    template <int dim, int spacedim>
     std::vector<ComponentMask>
-    compute_nonzero_components (FESystems &&... fe_systems);
+    compute_nonzero_components
+    (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>> &fe_systems);
 
     /**
      * Compute the non-zero vector components of a composed finite
@@ -1382,9 +1376,10 @@ namespace FETools
 
   namespace Compositing
   {
-    template <int dim, int spacedim, class... FESystems, typename>
+    template <int dim, int spacedim>
     std::vector<bool>
-    compute_restriction_is_additive_flags (FESystems &&... fe_systems)
+    compute_restriction_is_additive_flags
+    (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>> &fe_systems)
     {
       std::vector<const FiniteElement<dim,spacedim>*> fes;
       std::vector<unsigned int> 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 <int dim, int spacedim, class... FESystems, typename>
+    template <int dim, int spacedim>
     FiniteElementData<dim>
-    multiply_dof_numbers (FESystems &&... fe_systems)
+    multiply_dof_numbers
+    (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>> &fe_systems)
     {
       std::vector<const FiniteElement<dim,spacedim>*> fes;
       std::vector<unsigned int> 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 <int dim, int spacedim, class... FESystems, typename>
+    template <int dim, int spacedim>
     std::vector<ComponentMask>
-    compute_nonzero_components (FESystems &&... fe_systems)
+    compute_nonzero_components
+    (const std::initializer_list<std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>> &fe_systems)
     {
       std::vector<const FiniteElement<dim,spacedim>*> fes;
       std::vector<unsigned int> 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);
     }
index 5412942c81139cf1a4cda8c29d6e8aef4343f95e..97d58fa60440cb88376cc5d1389b2d6cb161faf2 100644 (file)
@@ -285,12 +285,6 @@ FESystem<dim,spacedim>::FESystem (
 
 
 
-template <int dim, int spacedim>
-FESystem<dim,spacedim>::~FESystem ()
-{}
-
-
-
 template <int dim, int spacedim>
 std::string
 FESystem<dim,spacedim>::get_name () const

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.