]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added hp::FECollection::get_hierarchy_sequence(). 11793/head
authorMarc Fehling <mafehling.git@gmail.com>
Mon, 22 Feb 2021 19:54:21 +0000 (12:54 -0700)
committerMarc Fehling <mafehling.git@gmail.com>
Wed, 24 Feb 2021 05:48:17 +0000 (22:48 -0700)
doc/news/changes/minor/20210222Fehling [new file with mode: 0644]
include/deal.II/hp/fe_collection.h
source/hp/fe_collection.cc
tests/hp/fe_hierarchy_sequence.cc [new file with mode: 0644]
tests/hp/fe_hierarchy_sequence.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20210222Fehling b/doc/news/changes/minor/20210222Fehling
new file mode 100644 (file)
index 0000000..fab75bc
--- /dev/null
@@ -0,0 +1,5 @@
+New: Member function hp::FECollection::get_hierarchy_sequence() returning
+the sequence of finite element indices that correspond to the registered
+hierarchy. 
+<br>
+(Marc Fehling, 2021/02/22)
index 07a32368ef2ddc315b2918e110ae4af313e785de..3ff44f2b2fafadff9e2faab926f6e3c7ed8f5a6b 100644 (file)
@@ -505,6 +505,28 @@ namespace hp
     void
     set_default_hierarchy();
 
+    /**
+     * Returns a sequence of FE indices that corresponds to the registered
+     * hierarchy in ascending order, i.e., FE indices are sorted from lowest to
+     * highest level.
+     *
+     * Multiple sequences of FE indices are possible with a single custom
+     * hierarchy that can be registered with set_hierarchy(). This function
+     * will return the sequence that contains the user-provided index
+     * @p fe_index which could be located anywhere inside the sequence. The
+     * default hierarchy set via set_default_hierarchy(), which corresponds to
+     * FE indices in ascending order, consists of only one sequence.
+     *
+     * This function can be used, for example, to verify that your provided
+     * hierarchy covers all elements in the desired order.
+     *
+     * Only one sequence of FE indices exists if the size of the returned
+     * container equals the number of elements of this object, i.e.,
+     * FECollection::size().
+     */
+    std::vector<unsigned int>
+    get_hierarchy_sequence(const unsigned int fe_index = 0) const;
+
     /**
      * Function returning the index of the finite element following the given
      * @p fe_index in hierarchy.
index c3f84a30fb8ed74629ff60b028816d5e6deda965..fbee5a2718c20a03c4bcdfe312468d8d3da2c581 100644 (file)
@@ -314,6 +314,52 @@ namespace hp
 
 
 
+  template <int dim, int spacedim>
+  std::vector<unsigned int>
+  FECollection<dim, spacedim>::get_hierarchy_sequence(
+    const unsigned int fe_index) const
+  {
+    AssertIndexRange(fe_index, size());
+
+    std::deque<unsigned int> sequence = {fe_index};
+
+    // get predecessors
+    {
+      unsigned int front = sequence.front();
+      unsigned int previous;
+      while ((previous = previous_in_hierarchy(front)) != front)
+        {
+          sequence.push_front(previous);
+          front = previous;
+
+          Assert(sequence.size() <= finite_elements.size(),
+                 ExcMessage(
+                   "The registered hierarchy is not terminated: "
+                   "previous_in_hierarchy() does not stop at a final index."));
+        }
+    }
+
+    // get successors
+    {
+      unsigned int back = sequence.back();
+      unsigned int next;
+      while ((next = next_in_hierarchy(back)) != back)
+        {
+          sequence.push_back(next);
+          back = next;
+
+          Assert(sequence.size() <= finite_elements.size(),
+                 ExcMessage(
+                   "The registered hierarchy is not terminated: "
+                   "next_in_hierarchy() does not stop at a final index."));
+        }
+    }
+
+    return {sequence.begin(), sequence.end()};
+  }
+
+
+
   template <int dim, int spacedim>
   unsigned int
   FECollection<dim, spacedim>::next_in_hierarchy(
diff --git a/tests/hp/fe_hierarchy_sequence.cc b/tests/hp/fe_hierarchy_sequence.cc
new file mode 100644 (file)
index 0000000..d0922db
--- /dev/null
@@ -0,0 +1,71 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// Verify sequences for a custom hierarchy registered on a hp::FECollection.
+
+
+#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;
+
+  // register custom hierarchy with two sequences: odd and even indices
+  fe_collection.set_hierarchy(
+    [](const hp::FECollection<dim> &fe_collection,
+       const unsigned int           fe_index) {
+      return ((fe_index + 2) < fe_collection.size()) ? fe_index + 2 : fe_index;
+    },
+    [](const hp::FECollection<dim> &fe_collection,
+       const unsigned int           fe_index) {
+      return (fe_index > 1) ? fe_index - 2 : fe_index;
+    });
+
+  // add dummy FEs to collection
+  while (fe_collection.size() < 6)
+    fe_collection.push_back(FE_Q<dim>(1));
+  deallog << "size: " << fe_collection.size() << std::endl;
+
+  // verify sequence for each FE index
+  for (unsigned int fe_index = 0; fe_index < fe_collection.size(); ++fe_index)
+    {
+      const auto sequence = fe_collection.get_hierarchy_sequence(fe_index);
+
+      deallog << " idx: " << fe_index << ", sequence:";
+      for (const auto index : sequence)
+        deallog << " " << index;
+      deallog << std::endl;
+    }
+}
+
+
+int
+main()
+{
+  initlog();
+
+  test<1>();
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/hp/fe_hierarchy_sequence.output b/tests/hp/fe_hierarchy_sequence.output
new file mode 100644 (file)
index 0000000..70a66da
--- /dev/null
@@ -0,0 +1,9 @@
+
+DEAL::size: 6
+DEAL:: idx: 0, sequence: 0 2 4
+DEAL:: idx: 1, sequence: 1 3 5
+DEAL:: idx: 2, sequence: 0 2 4
+DEAL:: idx: 3, sequence: 1 3 5
+DEAL:: idx: 4, sequence: 0 2 4
+DEAL:: idx: 5, sequence: 1 3 5
+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.