--- /dev/null
+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)
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.
+ 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(
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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;
+}
--- /dev/null
+
+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