]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added hierarchy to hp::FECollection.
authorMarc Fehling <marc.fehling@gmx.net>
Tue, 26 Feb 2019 15:34:09 +0000 (16:34 +0100)
committerMarc Fehling <marc.fehling@gmx.net>
Mon, 4 Mar 2019 20:45:15 +0000 (21:45 +0100)
doc/doxygen/headers/hp.h
doc/news/changes/minor/20190228Fehling [new file with mode: 0644]
include/deal.II/hp/fe_collection.h
source/hp/fe_collection.cc
tests/hp/fe_hierarchy.cc [new file with mode: 0644]
tests/hp/fe_hierarchy.output [new file with mode: 0644]

index 5fceed9e248cc658019970d9f475814cc515543d..b9c69531681f735c7d8b71dad8bc1dea746d8e42 100644 (file)
  * collection of mappings is used, the same holds for hp::MappingCollection
  * objects as well.
  *
+ * Whenever p adaptivity is considered in an hp finite element program,
+ * a hierarchy of finite elements needs to be established to determine
+ * succeeding finite elements for refinement and preceding ones for coarsening.
+ * Typically, this hierarchy considers how finite element spaces are nested:
+ * for example, a $Q_1$ element describes a sub-space of a $Q_2$ element,
+ * and so doing $p$ refinement usually means using a larger (more accurate)
+ * finite element space. In other words, the hierarchy of finite elements is built
+ * by considering whether some elements of the collection are sub- or
+ * super-spaces of others.
+ *
+ * By default, we assume that finite elements are stored in an ascending order
+ * based on their polynomial degree. If the order of elements differs,
+ * a corresponding hierarchy needs to be supplied to the collection via the
+ * hp::FECollection::set_hierarchy() member function.
+ *
  * @ingroup hp
  */
 
diff --git a/doc/news/changes/minor/20190228Fehling b/doc/news/changes/minor/20190228Fehling
new file mode 100644 (file)
index 0000000..85c0089
--- /dev/null
@@ -0,0 +1,7 @@
+New: Hierarchy of finite elements in hp::FECollection objects. Get succeeding
+and preceding indices via hp::FECollection::next_in_hierarchy() and
+hp::FECollection::prev_in_hierarchy(). By default, a hierarchy corresponding
+to indices is established. Hierarchy can be overridden via
+hp::FECollection::set_hierarchy().
+<br>
+(Marc Fehling, 2019/02/28)
index a58bef7ea2dcccd4b59e08db768a1d1ae1024190..c4668b89dcc5e99f06b678ac4cb5af0ac9b9fd28 100644 (file)
@@ -54,11 +54,53 @@ namespace hp
   class FECollection : public Subscriptor
   {
   public:
+    /**
+     * Whenever p adaptivity is considered in an hp finite element program,
+     * a hierarchy of finite elements needs to be established to determine
+     * succeeding finite elements for refinement and preceding ones for
+     * coarsening.
+     *
+     * In this struct, we supply a hierachy that is imposed on all FECollection
+     * objects by default.
+     */
+    struct DefaultHierarchy
+    {
+      /**
+       * Return the index succeeding @p fe_index in the @p fe_collection.
+       *
+       * Once the last element of the @p fe_collection is reached, there is no element on a higher level in
+       * the hierarchy and thus we return the last value.
+       */
+      static unsigned int
+      next_index(const typename hp::FECollection<dim, spacedim> &fe_collection,
+                 const unsigned int                              fe_index)
+      {
+        return ((fe_index + 1) < fe_collection.size()) ? fe_index + 1 :
+                                                         fe_index;
+      }
+
+      /**
+       * Return the index preceding @p fe_index in the @p fe_collection.
+       *
+       * Once the first element of the @p fe_collection is reached, there is no element on a lower level in
+       * the hierarchy and thus we return the first value.
+       */
+      static unsigned int
+      previous_index(
+        const typename hp::FECollection<dim, spacedim> &fe_collection,
+        const unsigned int                              fe_index)
+      {
+        (void)fe_collection;
+        return (fe_index > 0) ? fe_index - 1 : fe_index;
+      }
+    };
+
     /**
      * Default constructor. Leads to an empty collection that can later be
-     * filled using push_back().
+     * filled using push_back(). Establishes a hierarchy of finite elements
+     * corresponding to their index in the collection.
      */
-    FECollection() = default;
+    FECollection();
 
     /**
      * Conversion constructor. This constructor creates a FECollection from a
@@ -94,7 +136,7 @@ namespace hp
     /**
      * Move constructor.
      */
-    FECollection(FECollection<dim, spacedim> &&) noexcept = default;
+    FECollection(FECollection<dim, spacedim> &&) = default;
 
     /**
      * Move assignment operator.
@@ -371,6 +413,60 @@ namespace hp
     find_dominating_fe_in_subset(const std::set<unsigned int> &fes,
                                  const unsigned int            codim = 0) const;
 
+    /**
+     * Set functions determining the hierarchy of finite elements, i.e. a
+     * function @p next that returns the index of the finite element following
+     * the given one, and a function @p prev returning the preceding one.
+     *
+     * Both functions expect an hp::FECollection to be passed along with a
+     * finite element index, on whose basis the new index will be found and
+     * returned.
+     *
+     * @note Both passed and returned indices have to be valid within the index
+     * range of this collection, i.e. within [0, size()).
+     */
+    void
+    set_hierarchy(const std::function<unsigned int(
+                    const typename hp::FECollection<dim, spacedim> &,
+                    const unsigned int)> &next,
+                  const std::function<unsigned int(
+                    const typename hp::FECollection<dim, spacedim> &,
+                    const unsigned int)> &prev);
+
+    /**
+     * Set the default hierarchy corresponding to the index of each finite
+     * element in the collection.
+     *
+     * This default hierarchy is established with functions
+     * DefaultHierarchy::next_index() and DefaultHierarchy::previous_index().
+     */
+    void
+    set_default_hierarchy();
+
+    /**
+     * Function returning the index of the finite element following the given
+     * @p fe_index in hierarchy.
+     *
+     * By default, the index succeeding @p fe_index will be returned. If @p fe_index
+     * already corresponds to the last index, the last index will be returned.
+     * A custom hierarchy can be supplied via the member function
+     * set_hierachy().
+     */
+    unsigned int
+    next_in_hierarchy(const unsigned int fe_index) const;
+
+    /**
+     * Function returning the index of the finite element preceding the given
+     * @p fe_index in hierarchy.
+     *
+     * By default, the index preceding @p fe_index will be returned. If @p fe_index
+     * already corresponds to the first index, the first index will be returned.
+     * A custom hierarchy can be supplied via the member function
+     * set_hierachy().
+     */
+    unsigned int
+    previous_in_hierarchy(const unsigned int fe_index) const;
+
     /**
      * Return a component mask with as many elements as this object has vector
      * components and of which exactly the one component is true that
@@ -576,6 +672,22 @@ namespace hp
      */
     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.
+     */
+    std::function<unsigned int(const typename hp::FECollection<dim, spacedim> &,
+                               const unsigned int)>
+      hierarchy_next;
+
+    /**
+     * Function returning the index of the finite element preceding the given
+     * one in hierarchy.
+     */
+    std::function<unsigned int(const typename hp::FECollection<dim, spacedim> &,
+                               const unsigned int)>
+      hierarchy_prev;
   };
 
 
index 616359767d63148407ccad28f7656864853bdba3..cc3c26634a8408b4ab09d2a2f0a7e4609f90b631 100644 (file)
@@ -160,9 +160,18 @@ namespace hp
 
 
 
+  template <int dim, int spacedim>
+  FECollection<dim, spacedim>::FECollection()
+  {
+    set_default_hierarchy();
+  }
+
+
+
   template <int dim, int spacedim>
   FECollection<dim, spacedim>::FECollection(
     const FiniteElement<dim, spacedim> &fe)
+    : FECollection()
   {
     push_back(fe);
   }
@@ -172,6 +181,7 @@ namespace hp
   template <int dim, int spacedim>
   FECollection<dim, spacedim>::FECollection(
     const std::vector<const FiniteElement<dim, spacedim> *> &fes)
+    : FECollection()
   {
     Assert(fes.size() > 0,
            ExcMessage("Need to pass at least one finite element."));
@@ -202,6 +212,64 @@ namespace hp
 
 
 
+  template <int dim, int spacedim>
+  void
+  FECollection<dim, spacedim>::set_hierarchy(
+    const std::function<
+      unsigned int(const typename hp::FECollection<dim, spacedim> &,
+                   const unsigned int)> &next,
+    const std::function<
+      unsigned int(const typename hp::FECollection<dim, spacedim> &,
+                   const unsigned int)> &prev)
+  {
+    // copy hierarchy functions
+    hierarchy_next = next;
+    hierarchy_prev = prev;
+  }
+
+
+
+  template <int dim, int spacedim>
+  void
+  FECollection<dim, spacedim>::set_default_hierarchy()
+  {
+    // establish hierarchy corresponding to order of indices
+    set_hierarchy(&DefaultHierarchy::next_index,
+                  &DefaultHierarchy::previous_index);
+  }
+
+
+
+  template <int dim, int spacedim>
+  unsigned int
+  FECollection<dim, spacedim>::next_in_hierarchy(
+    const unsigned int fe_index) const
+  {
+    Assert(fe_index < size(), ExcIndexRange(fe_index, 0, size()));
+
+    const unsigned int new_fe_index = hierarchy_next(*this, fe_index);
+    Assert(new_fe_index < size(), ExcIndexRange(new_fe_index, 0, size()));
+
+    return new_fe_index;
+  }
+
+
+
+  template <int dim, int spacedim>
+  unsigned int
+  FECollection<dim, spacedim>::previous_in_hierarchy(
+    const unsigned int fe_index) const
+  {
+    Assert(fe_index < size(), ExcIndexRange(fe_index, 0, size()));
+
+    const unsigned int new_fe_index = hierarchy_prev(*this, fe_index);
+    Assert(new_fe_index < size(), ExcIndexRange(new_fe_index, 0, size()));
+
+    return new_fe_index;
+  }
+
+
+
   template <int dim, int spacedim>
   ComponentMask
   FECollection<dim, spacedim>::component_mask(
diff --git a/tests/hp/fe_hierarchy.cc b/tests/hp/fe_hierarchy.cc
new file mode 100644 (file)
index 0000000..185421a
--- /dev/null
@@ -0,0 +1,61 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2005 - 2019 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Test default hierarchy of any hp::FECollection object.
+// In the default implementation, the new indices are either the
+// succeeding or preceding indices until limits are reached.
+
+
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/hp/fe_collection.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  hp::FECollection<dim> fe_collection;
+
+  while (fe_collection.size() < 3)
+    {
+      // add dummy fe to collection
+      fe_collection.push_back(FE_Q<dim>(1));
+      deallog << "size:" << fe_collection.size() << std::endl;
+
+      for (unsigned int fe_index = 0; fe_index < fe_collection.size();
+           ++fe_index)
+        deallog << " idx:" << fe_index
+                << " next:" << fe_collection.next_in_hierarchy(fe_index)
+                << " prev:" << fe_collection.previous_in_hierarchy(fe_index)
+                << std::endl;
+    }
+}
+
+
+
+int
+main()
+{
+  initlog();
+
+  test<1>();
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/hp/fe_hierarchy.output b/tests/hp/fe_hierarchy.output
new file mode 100644 (file)
index 0000000..64e52b6
--- /dev/null
@@ -0,0 +1,11 @@
+
+DEAL::size:1
+DEAL:: idx:0 next:0 prev:0
+DEAL::size:2
+DEAL:: idx:0 next:1 prev:0
+DEAL:: idx:1 next:1 prev:0
+DEAL::size:3
+DEAL:: idx:0 next:1 prev:0
+DEAL:: idx:1 next:2 prev:0
+DEAL:: idx:2 next:2 prev:1
+DEAL::OK

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.