From 895d384dcccc0cf55ca2b9e1f22033ad27cf0393 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Mon, 5 Oct 2020 19:34:25 -0600 Subject: [PATCH] Clean up hp-coarsening code. --- doc/news/changes/minor/20201005Fehling | 4 + include/deal.II/dofs/dof_accessor.h | 42 +++++++++++ include/deal.II/dofs/dof_accessor.templates.h | 34 +++++++++ include/deal.II/hp/fe_collection.h | 10 --- source/distributed/cell_weights.cc | 30 ++------ source/distributed/error_predictor.cc | 21 ++---- source/distributed/solution_transfer.cc | 30 ++------ source/dofs/dof_handler.cc | 34 +++------ source/hp/refinement.cc | 28 ++----- source/numerics/solution_transfer.cc | 6 +- tests/grid/accessor_04.cc | 75 +++++++++++++++++++ tests/grid/accessor_04.output | 4 + 12 files changed, 203 insertions(+), 115 deletions(-) create mode 100644 doc/news/changes/minor/20201005Fehling create mode 100644 tests/grid/accessor_04.cc create mode 100644 tests/grid/accessor_04.output diff --git a/doc/news/changes/minor/20201005Fehling b/doc/news/changes/minor/20201005Fehling new file mode 100644 index 0000000000..a9d4679274 --- /dev/null +++ b/doc/news/changes/minor/20201005Fehling @@ -0,0 +1,4 @@ +New: Helper function DoFCellAccessor::dominated_future_fe_on_children() +to clean up code for hp-coarsening. +
+(Marc Fehling, 2020/10/05) diff --git a/include/deal.II/dofs/dof_accessor.h b/include/deal.II/dofs/dof_accessor.h index 14fe641575..75b130f0f9 100644 --- a/include/deal.II/dofs/dof_accessor.h +++ b/include/deal.II/dofs/dof_accessor.h @@ -2039,6 +2039,48 @@ public: */ void clear_future_fe_index() const; + + /** + * Return the index of the finite element from the entire hp::FECollection + * that is dominated by those assigned as future finite elements to the + * children of this cell. + * + * We find the corresponding finite element among the future finite elements + * on the children of this cell. If none of them qualify, we extend our search + * on the whole hp::FECollection, which is the element that describes the + * smallest finite element space that includes all future finite elements + * assigned to the children. If the function is not able to find a finite + * element at all, an assertion will be triggered. + * + * In this way, we determine the finite element of the parent cell in case of + * h-coarsening in the hp-context. + * + * @note This function can only be called on direct parent cells, i.e., + * non-active cells whose children are all active. + */ + unsigned int + dominated_future_fe_on_children() const; + /** + * @} + */ + + /** + * @name Exceptions + */ + /** + * @{ + */ + + /** + * Exception + * + * @ingroup Exceptions + */ + DeclExceptionMsg( + ExcNoDominatedFiniteElementOnChildren, + "No FiniteElement has been found in your FECollection that is " + "dominated by all children of a cell you are trying to coarsen!"); + /** * @} */ diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index 3284b305fb..40f1038ad0 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -2920,6 +2920,40 @@ DoFCellAccessor:: +template +inline unsigned int +DoFCellAccessor:: + dominated_future_fe_on_children() const +{ + Assert(!this->is_active(), + ExcMessage( + "You ask for information on children of this cell which is only " + "available for active cells. This cell has no children.")); + + std::set future_fe_indices_children; + for (const auto &child : this->child_iterators()) + { + Assert( + child->is_active(), + ExcMessage( + "You ask for information on children of this cell which is only " + "available for active cells. One of its children is not active.")); + future_fe_indices_children.insert(child->future_fe_index()); + } + Assert(!future_fe_indices_children.empty(), ExcInternalError()); + + const unsigned int future_fe_index = + this->dof_handler->get_fe_collection().find_dominated_fe_extended( + future_fe_indices_children, /*codim=*/0); + + Assert(future_fe_index != numbers::invalid_unsigned_int, + ExcNoDominatedFiniteElementOnChildren()); + + return future_fe_index; +} + + + template template inline void diff --git a/include/deal.II/hp/fe_collection.h b/include/deal.II/hp/fe_collection.h index c6e6f87956..5c319f378b 100644 --- a/include/deal.II/hp/fe_collection.h +++ b/include/deal.II/hp/fe_collection.h @@ -735,16 +735,6 @@ namespace hp */ DeclException0(ExcNoFiniteElements); - /** - * Exception - * - * @ingroup Exceptions - */ - DeclExceptionMsg( - ExcNoDominatedFiniteElementAmongstChildren, - "No FiniteElement has been found in your FECollection that is " - "dominated by all children of a cell you are trying to coarsen!"); - /** * @} */ diff --git a/source/distributed/cell_weights.cc b/source/distributed/cell_weights.cc index 72f0c0850f..8aaefbcd4e 100644 --- a/source/distributed/cell_weights.cc +++ b/source/distributed/cell_weights.cc @@ -192,28 +192,14 @@ namespace parallel break; case Triangulation::CELL_COARSEN: - { - std::set fe_indices_children; - for (unsigned int child_index = 0; child_index < cell->n_children(); - ++child_index) - { - const auto &child = cell->child(child_index); - Assert(child->is_active() && child->coarsen_flag_set(), - typename dealii::Triangulation< - dim>::ExcInconsistentCoarseningFlags()); - - fe_indices_children.insert(child->future_fe_index()); - } - Assert(!fe_indices_children.empty(), ExcInternalError()); - - fe_index = - dof_handler.get_fe_collection().find_dominated_fe_extended( - fe_indices_children, /*codim=*/0); - - Assert(fe_index != numbers::invalid_unsigned_int, - typename dealii::hp::FECollection< - dim>::ExcNoDominatedFiniteElementAmongstChildren()); - } +#ifdef DEBUG + for (const auto &child : cell->child_iterators()) + Assert(child->is_active() && child->coarsen_flag_set(), + typename dealii::Triangulation< + dim>::ExcInconsistentCoarseningFlags()); +#endif + + fe_index = cell->dominated_future_fe_on_children(); break; default: diff --git a/source/distributed/error_predictor.cc b/source/distributed/error_predictor.cc index ef12125b81..65d8bbba29 100644 --- a/source/distributed/error_predictor.cc +++ b/source/distributed/error_predictor.cc @@ -233,22 +233,15 @@ namespace parallel // parent cell after h-adaptation analogously to // dealii::internal::hp::DoFHandlerImplementation::Implementation:: // collect_fe_indices_on_cells_to_be_refined() - std::set fe_indices_children; - for (unsigned int child_index = 0; - child_index < cell->n_children(); - ++child_index) - { - const auto child = cell->child(child_index); - Assert(child->is_active() && child->coarsen_flag_set(), - typename dealii::Triangulation< - dim>::ExcInconsistentCoarseningFlags()); - - fe_indices_children.insert(child->future_fe_index()); - } +# ifdef DEBUG + for (const auto &child : cell->child_iterators()) + Assert(child->is_active() && child->coarsen_flag_set(), + typename dealii::Triangulation< + dim>::ExcInconsistentCoarseningFlags()); +# endif const unsigned int future_fe_index = - dof_handler->get_fe_collection().find_dominated_fe_extended( - fe_indices_children, /*codim=*/0); + cell->dominated_future_fe_on_children(); const unsigned int future_fe_degree = dof_handler->get_fe_collection()[future_fe_index].degree; diff --git a/source/distributed/solution_transfer.cc b/source/distributed/solution_transfer.cc index 481990c56d..c1862d0126 100644 --- a/source/distributed/solution_transfer.cc +++ b/source/distributed/solution_transfer.cc @@ -311,28 +311,14 @@ namespace parallel // In case of coarsening, we need to find a suitable fe index // for the parent cell. We choose the 'least dominant fe' // on all children from the associated FECollection. - std::set fe_indices_children; - for (unsigned int child_index = 0; - child_index < cell->n_children(); - ++child_index) - { - const auto &child = cell->child(child_index); - Assert(child->is_active() && child->coarsen_flag_set(), - typename dealii::Triangulation< - dim>::ExcInconsistentCoarseningFlags()); - - fe_indices_children.insert(child->future_fe_index()); - } - Assert(!fe_indices_children.empty(), ExcInternalError()); - - fe_index = - dof_handler->get_fe_collection().find_dominated_fe_extended( - fe_indices_children, /*codim=*/0); - - Assert(fe_index != numbers::invalid_unsigned_int, - typename dealii::hp::FECollection< - dim>::ExcNoDominatedFiniteElementAmongstChildren()); - +# ifdef DEBUG + for (const auto &child : cell->child_iterators()) + Assert(child->is_active() && child->coarsen_flag_set(), + typename dealii::Triangulation< + dim>::ExcInconsistentCoarseningFlags()); +# endif + + fe_index = cell->dominated_future_fe_on_children(); break; } diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index 22a0c26f5a..9f75ec85c4 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -1879,30 +1879,16 @@ namespace internal // based on the 'least dominant finite element' of its // children. Consider the childrens' hypothetical future // index when they have been flagged for p-refinement. - std::set fe_indices_children; - for (unsigned int child_index = 0; - child_index < parent->n_children(); - ++child_index) - { - const auto sibling = parent->child(child_index); - Assert(sibling->is_active() && - sibling->coarsen_flag_set(), - typename dealii::Triangulation< - dim>::ExcInconsistentCoarseningFlags()); - - fe_indices_children.insert( - sibling->future_fe_index()); - } - Assert(!fe_indices_children.empty(), - ExcInternalError()); +#ifdef DEBUG + for (const auto &child : parent->child_iterators()) + Assert(child->is_active() && + child->coarsen_flag_set(), + typename dealii::Triangulation< + dim>::ExcInconsistentCoarseningFlags()); +#endif const unsigned int fe_index = - dof_handler.fe_collection.find_dominated_fe_extended( - fe_indices_children, /*codim=*/0); - - Assert(fe_index != numbers::invalid_unsigned_int, - typename dealii::hp::FECollection:: - ExcNoDominatedFiniteElementAmongstChildren()); + parent->dominated_future_fe_on_children(); fe_transfer->coarsened_cells_fe_index.insert( {parent, fe_index}); @@ -1999,8 +1985,8 @@ namespace internal /*codim=*/0); Assert(dominated_fe_index != numbers::invalid_unsigned_int, - typename dealii::hp::FECollection< - dim>::ExcNoDominatedFiniteElementAmongstChildren()); + (typename dealii::DoFCellAccessor:: + ExcNoDominatedFiniteElementOnChildren())); return dominated_fe_index; } diff --git a/source/hp/refinement.cc b/source/hp/refinement.cc index 3e80a37954..7c9c02a494 100644 --- a/source/hp/refinement.cc +++ b/source/hp/refinement.cc @@ -542,29 +542,15 @@ namespace hp if (future_fe_indices_on_coarsened_cells.find(parent) == future_fe_indices_on_coarsened_cells.end()) { - std::set fe_indices_children; - for (unsigned int child_index = 0; - child_index < parent->n_children(); - ++child_index) - { - const auto &child = parent->child(child_index); - Assert(child->is_active() && child->coarsen_flag_set(), - typename dealii::Triangulation< - dim>::ExcInconsistentCoarseningFlags()); - - fe_indices_children.insert(child->future_fe_index()); - } - Assert(!fe_indices_children.empty(), ExcInternalError()); +#ifdef DEBUG + for (const auto &child : parent->child_iterators()) + Assert(child->is_active() && child->coarsen_flag_set(), + typename dealii::Triangulation< + dim>::ExcInconsistentCoarseningFlags()); +#endif parent_future_fe_index = - dof_handler.get_fe_collection() - .find_dominated_fe_extended(fe_indices_children, - /*codim=*/0); - - Assert( - parent_future_fe_index != numbers::invalid_unsigned_int, - typename dealii::hp::FECollection< - dim>::ExcNoDominatedFiniteElementAmongstChildren()); + parent->dominated_future_fe_on_children(); future_fe_indices_on_coarsened_cells.insert( {parent, parent_future_fe_index}); diff --git a/source/numerics/solution_transfer.cc b/source/numerics/solution_transfer.cc index 67c25059ae..4ea1fa8aa6 100644 --- a/source/numerics/solution_transfer.cc +++ b/source/numerics/solution_transfer.cc @@ -372,8 +372,10 @@ SolutionTransfer:: fe_indices_children, /*codim=*/0); Assert(target_fe_index != numbers::invalid_unsigned_int, - typename dealii::hp::FECollection< - dim>::ExcNoDominatedFiniteElementAmongstChildren()); + (typename dealii::DoFCellAccessor< + dim, + DoFHandlerType::space_dimension, + false>::ExcNoDominatedFiniteElementOnChildren())); const unsigned int dofs_per_cell = dof_handler->get_fe(target_fe_index).n_dofs_per_cell(); diff --git a/tests/grid/accessor_04.cc b/tests/grid/accessor_04.cc new file mode 100644 index 0000000000..578c2d953b --- /dev/null +++ b/tests/grid/accessor_04.cc @@ -0,0 +1,75 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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 functionality of +// DoFCellAccessor:dominated_future_fe_on_children() + +#include + +#include +#include + +#include + +#include +#include + +#include + +#include "../tests.h" + + + +template +void +test() +{ + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + + hp::FECollection fes; + fes.push_back(FE_Q(1)); + fes.push_back(FE_Q(2)); + + DoFHandler dofh(tria, /*hp_capability_enabled=*/true); + dofh.set_fe(fes); + dofh.begin_active()->set_active_fe_index(1); + + const auto & parent = dofh.begin(/*level=*/0); + const unsigned int parent_future_fe = + parent->dominated_future_fe_on_children(); + + deallog << parent_future_fe << std::endl; +} + + +int +main() +{ + initlog(); + + deallog.push("1d"); + test<1>(); + deallog.pop(); + deallog.push("2d"); + test<2>(); + deallog.pop(); + deallog.push("3d"); + test<3>(); + deallog.pop(); +} diff --git a/tests/grid/accessor_04.output b/tests/grid/accessor_04.output new file mode 100644 index 0000000000..8b5ba12227 --- /dev/null +++ b/tests/grid/accessor_04.output @@ -0,0 +1,4 @@ + +DEAL:1d::1 +DEAL:2d::1 +DEAL:3d::1 -- 2.39.5