]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Introduce hp::Collection 11934/head
authorPeter Munch <peterrmuench@gmail.com>
Wed, 17 Mar 2021 14:24:20 +0000 (15:24 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Wed, 17 Mar 2021 15:06:20 +0000 (16:06 +0100)
include/deal.II/hp/collection.h [new file with mode: 0644]
include/deal.II/hp/fe_collection.h
include/deal.II/hp/mapping_collection.h
include/deal.II/hp/q_collection.h
source/hp/fe_collection.cc
source/hp/mapping_collection.cc

diff --git a/include/deal.II/hp/collection.h b/include/deal.II/hp/collection.h
new file mode 100644 (file)
index 0000000..47e0bc5
--- /dev/null
@@ -0,0 +1,128 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_hp_collection_h
+#define dealii_hp_collection_h
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/base/memory_consumption.h>
+#include <deal.II/base/subscriptor.h>
+
+#include <memory>
+#include <vector>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace hp
+{
+  /**
+   * This class implements a collection of objects.
+   *
+   * It implements the concepts stated in the @ref hpcollection
+   * module described in the doxygen documentation.
+   *
+   * @ingroup hp hpcollection
+   */
+  template <typename T>
+  class Collection : public Subscriptor
+  {
+  public:
+    /**
+     * Default constructor. Leads to an empty collection that can later be
+     * filled using push_back().
+     */
+    Collection() = default;
+
+    /**
+     * Add a new object.
+     */
+    void
+    push_back(const std::shared_ptr<const T> &new_entry);
+
+    /**
+     * Return the object which was specified by the user for the
+     * active FE index which is provided as a parameter to this method.
+     *
+     * @pre @p index must be between zero and the number of elements of the
+     * collection.
+     */
+    const T &operator[](const unsigned int index) const;
+
+    /**
+     * Return the number of objects stored in this container.
+     */
+    unsigned int
+    size() const;
+
+    /**
+     * Determine an estimate for the memory consumption (in bytes) of this
+     * object.
+     */
+    std::size_t
+    memory_consumption() const;
+
+  private:
+    /**
+     * The real container, which stores pointers to the different objects.
+     */
+    std::vector<std::shared_ptr<const T>> entries;
+  };
+
+
+  /* --------------- inline functions ------------------- */
+
+
+
+  template <typename T>
+  std::size_t
+  Collection<T>::memory_consumption() const
+  {
+    return (sizeof(*this) + MemoryConsumption::memory_consumption(entries));
+  }
+
+
+
+  template <typename T>
+  void
+  Collection<T>::push_back(const std::shared_ptr<const T> &new_entry)
+  {
+    entries.push_back(new_entry);
+  }
+
+
+
+  template <typename T>
+  inline unsigned int
+  Collection<T>::size() const
+  {
+    return entries.size();
+  }
+
+
+
+  template <typename T>
+  inline const T &Collection<T>::operator[](const unsigned int index) const
+  {
+    AssertIndexRange(index, entries.size());
+    return *entries[index];
+  }
+
+} // namespace hp
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 142611c502d974869f9d1e042b033e0c4c19c43e..c022a82a44a375d13077cb17fe1bc1e703fc06b9 100644 (file)
@@ -22,6 +22,8 @@
 #include <deal.II/fe/fe.h>
 #include <deal.II/fe/fe_values_extractors.h>
 
+#include <deal.II/hp/collection.h>
+
 #include <memory>
 
 DEAL_II_NAMESPACE_OPEN
@@ -48,7 +50,7 @@ namespace hp
    * @ingroup hp hpcollection
    */
   template <int dim, int spacedim = dim>
-  class FECollection : public Subscriptor
+  class FECollection : public Collection<FiniteElement<dim, spacedim>>
   {
   public:
     /**
@@ -178,21 +180,6 @@ namespace hp
     void
     push_back(const FiniteElement<dim, spacedim> &new_fe);
 
-    /**
-     * Get a reference to the given element in this collection.
-     *
-     * @pre @p index must be between zero and the number of elements of the
-     * collection.
-     */
-    const FiniteElement<dim, spacedim> &
-    operator[](const unsigned int index) const;
-
-    /**
-     * Return the number of finite element objects stored in this collection.
-     */
-    unsigned int
-    size() const;
-
     /**
      * Return the number of vector components of the finite elements in this
      * collection.  This number must be the same for all elements in the
@@ -273,12 +260,6 @@ namespace hp
     unsigned int
     max_dofs_per_cell() const;
 
-    /**
-     * Return an estimate for the memory allocated for this object.
-     */
-    std::size_t
-    memory_consumption() const;
-
 
     /**
      * Return whether all elements in this collection implement the hanging
@@ -761,12 +742,6 @@ namespace hp
      */
 
   private:
-    /**
-     * Array of pointers to the finite elements stored by this collection.
-     */
-    std::vector<std::shared_ptr<const FiniteElement<dim, spacedim>>>
-      finite_elements;
-
     /**
      * %Function returning the index of the finite element following the given
      * one in hierarchy.
@@ -807,19 +782,12 @@ namespace hp
   }
 
 
-  template <int dim, int spacedim>
-  inline unsigned int
-  FECollection<dim, spacedim>::size() const
-  {
-    return finite_elements.size();
-  }
-
 
   template <int dim, int spacedim>
   inline unsigned int
   FECollection<dim, spacedim>::n_components() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     // note that there is no need
     // here to enforce that indeed
@@ -829,7 +797,7 @@ namespace hp
     // adding a new element to the
     // collection.
 
-    return finite_elements[0]->n_components();
+    return this->operator[](0).n_components();
   }
 
 
@@ -839,12 +807,12 @@ namespace hp
   FECollection<dim, spacedim>::
   operator==(const FECollection<dim, spacedim> &fe_collection) const
   {
-    const unsigned int n_elements = size();
+    const unsigned int n_elements = this->size();
     if (n_elements != fe_collection.size())
       return false;
 
     for (unsigned int i = 0; i < n_elements; ++i)
-      if (!(*finite_elements[i] == fe_collection[i]))
+      if (!(this->operator[](i) == fe_collection[i]))
         return false;
 
     return true;
@@ -862,25 +830,15 @@ namespace hp
 
 
 
-  template <int dim, int spacedim>
-  inline const FiniteElement<dim, spacedim> &FECollection<dim, spacedim>::
-                                             operator[](const unsigned int index) const
-  {
-    AssertIndexRange(index, finite_elements.size());
-    return *finite_elements[index];
-  }
-
-
-
   template <int dim, int spacedim>
   unsigned int
   FECollection<dim, spacedim>::max_degree() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     unsigned int max = 0;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      max = std::max(max, finite_elements[i]->degree);
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).degree);
 
     return max;
   }
@@ -891,11 +849,11 @@ namespace hp
   unsigned int
   FECollection<dim, spacedim>::max_dofs_per_vertex() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     unsigned int max = 0;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      max = std::max(max, finite_elements[i]->n_dofs_per_vertex());
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).n_dofs_per_vertex());
 
     return max;
   }
@@ -906,11 +864,11 @@ namespace hp
   unsigned int
   FECollection<dim, spacedim>::max_dofs_per_line() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     unsigned int max = 0;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      max = std::max(max, finite_elements[i]->n_dofs_per_line());
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).n_dofs_per_line());
 
     return max;
   }
@@ -921,11 +879,11 @@ namespace hp
   unsigned int
   FECollection<dim, spacedim>::max_dofs_per_quad() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     unsigned int max = 0;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      max = std::max(max, finite_elements[i]->max_dofs_per_quad());
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).max_dofs_per_quad());
 
     return max;
   }
@@ -936,11 +894,11 @@ namespace hp
   unsigned int
   FECollection<dim, spacedim>::max_dofs_per_hex() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     unsigned int max = 0;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      max = std::max(max, finite_elements[i]->n_dofs_per_hex());
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).n_dofs_per_hex());
 
     return max;
   }
@@ -951,11 +909,11 @@ namespace hp
   unsigned int
   FECollection<dim, spacedim>::max_dofs_per_face() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     unsigned int max = 0;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      max = std::max(max, finite_elements[i]->max_dofs_per_face());
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).max_dofs_per_face());
 
     return max;
   }
@@ -966,11 +924,11 @@ namespace hp
   unsigned int
   FECollection<dim, spacedim>::max_dofs_per_cell() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
     unsigned int max = 0;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      max = std::max(max, finite_elements[i]->n_dofs_per_cell());
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).n_dofs_per_cell());
 
     return max;
   }
@@ -980,13 +938,13 @@ namespace hp
   bool
   FECollection<dim, spacedim>::hp_constraints_are_implemented() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
-    return std::all_of(
-      finite_elements.cbegin(),
-      finite_elements.cend(),
-      [](const std::shared_ptr<const FiniteElement<dim, spacedim>> &fe) {
-        return fe->hp_constraints_are_implemented();
-      });
+    Assert(this->size() > 0, ExcNoFiniteElements());
+
+    for (unsigned int i = 0; i < this->size(); ++i)
+      if (this->operator[](i).hp_constraints_are_implemented() == false)
+        return false;
+
+    return true;
   }
 
 
index e0f321582cbb4f395aa6f9fecfe7b48fb4e3ffd4..118b82c4dec8c2ee66ea06382aceb83c868ce88c 100644 (file)
@@ -23,6 +23,8 @@
 #include <deal.II/fe/fe.h>
 #include <deal.II/fe/mapping_q1.h>
 
+#include <deal.II/hp/collection.h>
+
 #include <memory>
 #include <vector>
 
@@ -51,7 +53,7 @@ namespace hp
    * @ingroup hp hpcollection
    */
   template <int dim, int spacedim = dim>
-  class MappingCollection : public Subscriptor
+  class MappingCollection : public Collection<Mapping<dim, spacedim>>
   {
   public:
     /**
@@ -97,35 +99,6 @@ namespace hp
      */
     void
     push_back(const Mapping<dim, spacedim> &new_mapping);
-
-    /**
-     * Return the mapping object which was specified by the user for the
-     * active FE index which is provided as a parameter to this method.
-     *
-     * @pre @p index must be between zero and the number of elements of the
-     * collection.
-     */
-    const Mapping<dim, spacedim> &operator[](const unsigned int index) const;
-
-    /**
-     * Return the number of mapping objects stored in this container.
-     */
-    unsigned int
-    size() const;
-
-    /**
-     * Determine an estimate for the memory consumption (in bytes) of this
-     * object.
-     */
-    std::size_t
-    memory_consumption() const;
-
-  private:
-    /**
-     * The real container, which stores pointers to the different Mapping
-     * objects.
-     */
-    std::vector<std::shared_ptr<const Mapping<dim, spacedim>>> mappings;
   };
 
 
@@ -177,23 +150,6 @@ namespace hp
       push_back(*p);
   }
 
-  template <int dim, int spacedim>
-  inline unsigned int
-  MappingCollection<dim, spacedim>::size() const
-  {
-    return mappings.size();
-  }
-
-
-
-  template <int dim, int spacedim>
-  inline const Mapping<dim, spacedim> &MappingCollection<dim, spacedim>::
-                                       operator[](const unsigned int index) const
-  {
-    AssertIndexRange(index, mappings.size());
-    return *mappings[index];
-  }
-
 } // namespace hp
 
 
index 46339178cd55570bf51f744a7f540e007813f513..692582d62d64a87f19959787bf6c8c88d7a672e8 100644 (file)
@@ -22,6 +22,8 @@
 #include <deal.II/base/quadrature.h>
 #include <deal.II/base/subscriptor.h>
 
+#include <deal.II/hp/collection.h>
+
 #include <memory>
 #include <vector>
 
@@ -41,7 +43,7 @@ namespace hp
    * @ingroup hp hpcollection
    */
   template <int dim>
-  class QCollection : public Subscriptor
+  class QCollection : public Collection<Quadrature<dim>>
   {
   public:
     /**
@@ -102,14 +104,6 @@ namespace hp
     void
     push_back(const Quadrature<dim_in> &new_quadrature);
 
-    /**
-     * Return a reference to the quadrature rule specified by the argument.
-     *
-     * @pre @p index must be between zero and the number of elements of the
-     * collection.
-     */
-    const Quadrature<dim> &operator[](const unsigned int index) const;
-
     /**
      * Equality comparison operator. All stored Quadrature objects are compared
      * in order.
@@ -117,12 +111,6 @@ namespace hp
     bool
     operator==(const QCollection<dim> &q_collection) const;
 
-    /**
-     * Return the number of quadrature pointers stored in this object.
-     */
-    unsigned int
-    size() const;
-
     /**
      * Return the maximum number of quadrature points over all the elements of
      * the collection. This is mostly useful to initialize arrays to allocate
@@ -132,13 +120,6 @@ namespace hp
     unsigned int
     max_n_quadrature_points() const;
 
-    /**
-     * Determine an estimate for the memory consumption (in bytes) of this
-     * object.
-     */
-    std::size_t
-    memory_consumption() const;
-
     /**
      * Exception
      */
@@ -195,38 +176,18 @@ namespace hp
 
 
 
-  template <int dim>
-  inline unsigned int
-  QCollection<dim>::size() const
-  {
-    return quadratures.size();
-  }
-
-
-
   template <int dim>
   inline unsigned int
   QCollection<dim>::max_n_quadrature_points() const
   {
-    Assert(quadratures.size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("You can't call this function for an empty collection"));
 
-    unsigned int m = 0;
-    for (unsigned int i = 0; i < quadratures.size(); ++i)
-      if (quadratures[i]->size() > m)
-        m = quadratures[i]->size();
-
-    return m;
-  }
-
-
+    unsigned int max = 0;
+    for (unsigned int i = 0; i < this->size(); ++i)
+      max = std::max(max, this->operator[](i).size());
 
-  template <int dim>
-  inline const Quadrature<dim> &QCollection<dim>::
-                                operator[](const unsigned int index) const
-  {
-    AssertIndexRange(index, quadratures.size());
-    return *quadratures[index];
+    return max;
   }
 
 
@@ -235,12 +196,12 @@ namespace hp
   inline bool
   QCollection<dim>::operator==(const QCollection<dim> &q_collection) const
   {
-    const unsigned int n_quadratures = size();
+    const unsigned int n_quadratures = this->size();
     if (n_quadratures != q_collection.size())
       return false;
 
     for (unsigned int i = 0; i < n_quadratures; ++i)
-      if (!(*quadratures[i] == q_collection[i]))
+      if ((this->operator[](i) == q_collection[i]) == false)
         return false;
 
     return true;
@@ -252,16 +213,7 @@ namespace hp
   template <int dim_in>
   inline QCollection<dim>::QCollection(const Quadrature<dim_in> &quadrature)
   {
-    quadratures.push_back(std::make_shared<const Quadrature<dim>>(quadrature));
-  }
-
-
-
-  template <int dim>
-  inline std::size_t
-  QCollection<dim>::memory_consumption() const
-  {
-    return (sizeof(*this) + MemoryConsumption::memory_consumption(quadratures));
+    this->push_back(quadrature);
   }
 
 
@@ -270,7 +222,7 @@ namespace hp
   inline void
   QCollection<dim>::push_back(const Quadrature<dim_in> &new_quadrature)
   {
-    quadratures.push_back(
+    Collection<Quadrature<dim>>::push_back(
       std::make_shared<const Quadrature<dim>>(new_quadrature));
   }
 
index fbee5a2718c20a03c4bcdfe312468d8d3da2c581..6b2c39215a3a0418b1677c161ea9271aedf6e072 100644 (file)
@@ -31,25 +31,25 @@ namespace hp
 #ifdef DEBUG
     // Validate user inputs.
     Assert(codim <= dim, ExcImpossibleInDim(dim));
-    Assert(size() > 0, ExcEmptyObject());
+    Assert(this->size() > 0, ExcEmptyObject());
     for (const auto &fe : fes)
-      AssertIndexRange(fe, finite_elements.size());
+      AssertIndexRange(fe, this->size());
 #endif
 
     // Check if any element of this FECollection is able to dominate all
     // elements of @p fes. If one was found, we add it to the set of
     // dominating elements.
     std::set<unsigned int> dominating_fes;
-    for (unsigned int current_fe = 0; current_fe < finite_elements.size();
-         ++current_fe)
+    for (unsigned int current_fe = 0; current_fe < this->size(); ++current_fe)
       {
         // Check if current_fe can dominate all elements in @p fes.
         FiniteElementDomination::Domination domination =
           FiniteElementDomination::no_requirements;
         for (const auto &other_fe : fes)
           domination =
-            domination & finite_elements[current_fe]->compare_for_domination(
-                           *finite_elements[other_fe], codim);
+            domination &
+            this->operator[](current_fe)
+              .compare_for_domination(this->operator[](other_fe), codim);
 
         // If current_fe dominates, add it to the set.
         if ((domination == FiniteElementDomination::this_element_dominates) ||
@@ -71,25 +71,25 @@ namespace hp
 #ifdef DEBUG
     // Validate user inputs.
     Assert(codim <= dim, ExcImpossibleInDim(dim));
-    Assert(size() > 0, ExcEmptyObject());
+    Assert(this->size() > 0, ExcEmptyObject());
     for (const auto &fe : fes)
-      AssertIndexRange(fe, finite_elements.size());
+      AssertIndexRange(fe, this->size());
 #endif
 
     // Check if any element of this FECollection is dominated by all
     // elements of @p fes. If one was found, we add it to the set of
     // dominated elements.
     std::set<unsigned int> dominated_fes;
-    for (unsigned int current_fe = 0; current_fe < finite_elements.size();
-         ++current_fe)
+    for (unsigned int current_fe = 0; current_fe < this->size(); ++current_fe)
       {
         // Check if current_fe is dominated by all other elements in @p fes.
         FiniteElementDomination::Domination domination =
           FiniteElementDomination::no_requirements;
         for (const auto &other_fe : fes)
           domination =
-            domination & finite_elements[current_fe]->compare_for_domination(
-                           *finite_elements[other_fe], codim);
+            domination &
+            this->operator[](current_fe)
+              .compare_for_domination(this->operator[](other_fe), codim);
 
         // If current_fe is dominated, add it to the set.
         if ((domination == FiniteElementDomination::other_element_dominates) ||
@@ -116,9 +116,9 @@ namespace hp
 #ifdef DEBUG
     // Validate user inputs.
     Assert(codim <= dim, ExcImpossibleInDim(dim));
-    Assert(size() > 0, ExcEmptyObject());
+    Assert(this->size() > 0, ExcEmptyObject());
     for (const auto &fe : fes)
-      AssertIndexRange(fe, finite_elements.size());
+      AssertIndexRange(fe, this->size());
 #endif
 
     // There may also be others, in which case we'll check if any of these
@@ -132,8 +132,9 @@ namespace hp
         for (const auto &other_fe : fes)
           if (current_fe != other_fe)
             domination =
-              domination & finite_elements[current_fe]->compare_for_domination(
-                             *finite_elements[other_fe], codim);
+              domination &
+              this->operator[](current_fe)
+                .compare_for_domination(this->operator[](other_fe), codim);
 
         // If current_fe dominates, return its index.
         if ((domination == FiniteElementDomination::this_element_dominates) ||
@@ -162,9 +163,9 @@ namespace hp
 #ifdef DEBUG
     // Validate user inputs.
     Assert(codim <= dim, ExcImpossibleInDim(dim));
-    Assert(size() > 0, ExcEmptyObject());
+    Assert(this->size() > 0, ExcEmptyObject());
     for (const auto &fe : fes)
-      AssertIndexRange(fe, finite_elements.size());
+      AssertIndexRange(fe, this->size());
 #endif
 
     // There may also be others, in which case we'll check if any of these
@@ -178,8 +179,9 @@ namespace hp
         for (const auto &other_fe : fes)
           if (current_fe != other_fe)
             domination =
-              domination & finite_elements[current_fe]->compare_for_domination(
-                             *finite_elements[other_fe], codim);
+              domination &
+              this->operator[](current_fe)
+                .compare_for_domination(this->operator[](other_fe), codim);
 
         // If current_fe is dominated, return its index.
         if ((domination == FiniteElementDomination::other_element_dominates) ||
@@ -271,17 +273,15 @@ namespace hp
   FECollection<dim, spacedim>::push_back(
     const FiniteElement<dim, spacedim> &new_fe)
   {
-    // check that the new element has the right
-    // number of components. only check with
-    // the first element, since all the other
-    // elements have already passed the test
-    // against the first element
-    if (finite_elements.size() != 0)
-      Assert(new_fe.n_components() == finite_elements[0]->n_components(),
-             ExcMessage("All elements inside a collection need to have the "
-                        "same number of vector components!"));
-
-    finite_elements.push_back(new_fe.clone());
+    // check that the new element has the right number of components. only check
+    // with the first element, since all the other elements have already passed
+    // the test against the first element
+    Assert(this->size() == 0 ||
+             new_fe.n_components() == this->operator[](0).n_components(),
+           ExcMessage("All elements inside a collection need to have the "
+                      "same number of vector components!"));
+
+    Collection<FiniteElement<dim, spacedim>>::push_back(new_fe.clone());
   }
 
 
@@ -319,7 +319,7 @@ namespace hp
   FECollection<dim, spacedim>::get_hierarchy_sequence(
     const unsigned int fe_index) const
   {
-    AssertIndexRange(fe_index, size());
+    AssertIndexRange(fe_index, this->size());
 
     std::deque<unsigned int> sequence = {fe_index};
 
@@ -332,7 +332,7 @@ namespace hp
           sequence.push_front(previous);
           front = previous;
 
-          Assert(sequence.size() <= finite_elements.size(),
+          Assert(sequence.size() <= this->size(),
                  ExcMessage(
                    "The registered hierarchy is not terminated: "
                    "previous_in_hierarchy() does not stop at a final index."));
@@ -348,7 +348,7 @@ namespace hp
           sequence.push_back(next);
           back = next;
 
-          Assert(sequence.size() <= finite_elements.size(),
+          Assert(sequence.size() <= this->size(),
                  ExcMessage(
                    "The registered hierarchy is not terminated: "
                    "next_in_hierarchy() does not stop at a final index."));
@@ -365,10 +365,10 @@ namespace hp
   FECollection<dim, spacedim>::next_in_hierarchy(
     const unsigned int fe_index) const
   {
-    AssertIndexRange(fe_index, size());
+    AssertIndexRange(fe_index, this->size());
 
     const unsigned int new_fe_index = hierarchy_next(*this, fe_index);
-    AssertIndexRange(new_fe_index, size());
+    AssertIndexRange(new_fe_index, this->size());
 
     return new_fe_index;
   }
@@ -380,10 +380,10 @@ namespace hp
   FECollection<dim, spacedim>::previous_in_hierarchy(
     const unsigned int fe_index) const
   {
-    AssertIndexRange(fe_index, size());
+    AssertIndexRange(fe_index, this->size());
 
     const unsigned int new_fe_index = hierarchy_prev(*this, fe_index);
-    AssertIndexRange(new_fe_index, size());
+    AssertIndexRange(new_fe_index, this->size());
 
     return new_fe_index;
   }
@@ -395,7 +395,7 @@ namespace hp
   FECollection<dim, spacedim>::component_mask(
     const FEValuesExtractors::Scalar &scalar) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -403,7 +403,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].component_mask(scalar), ExcInternalError());
 
     return mask;
@@ -415,7 +415,7 @@ namespace hp
   FECollection<dim, spacedim>::component_mask(
     const FEValuesExtractors::Vector &vector) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -423,7 +423,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].component_mask(vector), ExcInternalError());
 
     return mask;
@@ -435,7 +435,7 @@ namespace hp
   FECollection<dim, spacedim>::component_mask(
     const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -443,7 +443,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].component_mask(sym_tensor), ExcInternalError());
 
     return mask;
@@ -454,7 +454,7 @@ namespace hp
   ComponentMask
   FECollection<dim, spacedim>::component_mask(const BlockMask &block_mask) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -462,7 +462,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].component_mask(block_mask),
              ExcMessage("Not all elements of this collection agree on what "
                         "the appropriate mask should be."));
@@ -476,7 +476,7 @@ namespace hp
   FECollection<dim, spacedim>::block_mask(
     const FEValuesExtractors::Scalar &scalar) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -484,7 +484,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].block_mask(scalar),
              ExcMessage("Not all elements of this collection agree on what "
                         "the appropriate mask should be."));
@@ -498,7 +498,7 @@ namespace hp
   FECollection<dim, spacedim>::block_mask(
     const FEValuesExtractors::Vector &vector) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -506,7 +506,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].block_mask(vector),
              ExcMessage("Not all elements of this collection agree on what "
                         "the appropriate mask should be."));
@@ -520,7 +520,7 @@ namespace hp
   FECollection<dim, spacedim>::block_mask(
     const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -528,7 +528,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].block_mask(sym_tensor),
              ExcMessage("Not all elements of this collection agree on what "
                         "the appropriate mask should be."));
@@ -543,7 +543,7 @@ namespace hp
   FECollection<dim, spacedim>::block_mask(
     const ComponentMask &component_mask) const
   {
-    Assert(size() > 0,
+    Assert(this->size() > 0,
            ExcMessage("This collection contains no finite element."));
 
     // get the mask from the first element of the collection
@@ -551,7 +551,7 @@ namespace hp
 
     // but then also verify that the other elements of the collection
     // would return the same mask
-    for (unsigned int c = 1; c < size(); ++c)
+    for (unsigned int c = 1; c < this->size(); ++c)
       Assert(mask == (*this)[c].block_mask(component_mask),
              ExcMessage("Not all elements of this collection agree on what "
                         "the appropriate mask should be."));
@@ -565,30 +565,16 @@ namespace hp
   unsigned int
   FECollection<dim, spacedim>::n_blocks() const
   {
-    Assert(finite_elements.size() > 0, ExcNoFiniteElements());
+    Assert(this->size() > 0, ExcNoFiniteElements());
 
-    const unsigned int nb = finite_elements[0]->n_blocks();
-    for (unsigned int i = 1; i < finite_elements.size(); ++i)
-      Assert(finite_elements[i]->n_blocks() == nb,
+    const unsigned int nb = this->operator[](0).n_blocks();
+    for (unsigned int i = 1; i < this->size(); ++i)
+      Assert(this->operator[](i).n_blocks() == nb,
              ExcMessage("Not all finite elements in this collection have "
                         "the same number of components."));
 
     return nb;
   }
-
-
-
-  template <int dim, int spacedim>
-  std::size_t
-  FECollection<dim, spacedim>::memory_consumption() const
-  {
-    std::size_t mem =
-      (sizeof(*this) + MemoryConsumption::memory_consumption(finite_elements));
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      mem += finite_elements[i]->memory_consumption();
-
-    return mem;
-  }
 } // namespace hp
 
 
index 9b7b275ba20f3303ecabce9e1dbcdc8bc12ef8ff..e257249e0408cf579b08586ab12ccc91f9d34346 100644 (file)
@@ -27,38 +27,17 @@ namespace hp
   MappingCollection<dim, spacedim>::MappingCollection(
     const Mapping<dim, spacedim> &mapping)
   {
-    mappings.push_back(
-      std::shared_ptr<const Mapping<dim, spacedim>>(mapping.clone()));
+    this->push_back(mapping);
   }
 
 
 
   template <int dim, int spacedim>
   MappingCollection<dim, spacedim>::MappingCollection(
-    const MappingCollection<dim, spacedim> &mapping_collection)
-    : Subscriptor()
-    ,
-    // copy the array
-    // of shared
-    // pointers. nothing
-    // bad should
-    // happen -- they
-    // simply all point
-    // to the same
-    // objects, and the
-    // last one to die
-    // will delete the
-    // mappings
-    mappings(mapping_collection.mappings)
-  {}
-
-
-
-  template <int dim, int spacedim>
-  std::size_t
-  MappingCollection<dim, spacedim>::memory_consumption() const
+    const MappingCollection<dim, spacedim> &other)
   {
-    return (sizeof(*this) + MemoryConsumption::memory_consumption(mappings));
+    for (unsigned int i = 0; i < other.size(); ++i)
+      push_back(other[i]);
   }
 
 
@@ -68,7 +47,7 @@ namespace hp
   MappingCollection<dim, spacedim>::push_back(
     const Mapping<dim, spacedim> &new_mapping)
   {
-    mappings.push_back(
+    Collection<Mapping<dim, spacedim>>::push_back(
       std::shared_ptr<const Mapping<dim, spacedim>>(new_mapping.clone()));
   }
 

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.