From 3322611456ed2598f80a9f661ef242c48bd087b2 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Mon, 22 Feb 2021 12:54:21 -0700 Subject: [PATCH] Added hp::FECollection::get_hierarchy_sequence(). --- doc/news/changes/minor/20210222Fehling | 5 ++ include/deal.II/hp/fe_collection.h | 22 ++++++++ source/hp/fe_collection.cc | 46 +++++++++++++++++ tests/hp/fe_hierarchy_sequence.cc | 71 ++++++++++++++++++++++++++ tests/hp/fe_hierarchy_sequence.output | 9 ++++ 5 files changed, 153 insertions(+) create mode 100644 doc/news/changes/minor/20210222Fehling create mode 100644 tests/hp/fe_hierarchy_sequence.cc create mode 100644 tests/hp/fe_hierarchy_sequence.output diff --git a/doc/news/changes/minor/20210222Fehling b/doc/news/changes/minor/20210222Fehling new file mode 100644 index 0000000000..fab75bc877 --- /dev/null +++ b/doc/news/changes/minor/20210222Fehling @@ -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. +
+(Marc Fehling, 2021/02/22) diff --git a/include/deal.II/hp/fe_collection.h b/include/deal.II/hp/fe_collection.h index 07a32368ef..3ff44f2b2f 100644 --- a/include/deal.II/hp/fe_collection.h +++ b/include/deal.II/hp/fe_collection.h @@ -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 + 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. diff --git a/source/hp/fe_collection.cc b/source/hp/fe_collection.cc index c3f84a30fb..fbee5a2718 100644 --- a/source/hp/fe_collection.cc +++ b/source/hp/fe_collection.cc @@ -314,6 +314,52 @@ namespace hp + template + std::vector + FECollection::get_hierarchy_sequence( + const unsigned int fe_index) const + { + AssertIndexRange(fe_index, size()); + + std::deque 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 unsigned int FECollection::next_in_hierarchy( diff --git a/tests/hp/fe_hierarchy_sequence.cc b/tests/hp/fe_hierarchy_sequence.cc new file mode 100644 index 0000000000..d0922dba40 --- /dev/null +++ b/tests/hp/fe_hierarchy_sequence.cc @@ -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 + +#include + +#include "../tests.h" + + +template +void +test() +{ + hp::FECollection fe_collection; + + // register custom hierarchy with two sequences: odd and even indices + fe_collection.set_hierarchy( + [](const hp::FECollection &fe_collection, + const unsigned int fe_index) { + return ((fe_index + 2) < fe_collection.size()) ? fe_index + 2 : fe_index; + }, + [](const hp::FECollection &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(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 index 0000000000..70a66da985 --- /dev/null +++ b/tests/hp/fe_hierarchy_sequence.output @@ -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 -- 2.39.5