From b2792d9f1d01603df69933ac0bc6b7609d0024e1 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Wed, 24 Mar 2021 16:22:12 -0600 Subject: [PATCH] Added `p:d:TemporarilyMatchRefineFlags`. --- include/deal.II/distributed/tria.h | 91 +++++++++++++++++++++++- source/distributed/tria.cc | 86 ++++++++++++++++++++++ source/distributed/tria.inst.in | 23 ++++-- tests/mpi/limit_p_level_difference_01.cc | 17 +++-- tests/mpi/limit_p_level_difference_02.cc | 18 +++-- tests/mpi/limit_p_level_difference_03.cc | 17 +++-- tests/mpi/limit_p_level_difference_04.cc | 19 +++-- 7 files changed, 243 insertions(+), 28 deletions(-) diff --git a/include/deal.II/distributed/tria.h b/include/deal.II/distributed/tria.h index 6deb092a97..3727241811 100644 --- a/include/deal.II/distributed/tria.h +++ b/include/deal.II/distributed/tria.h @@ -73,6 +73,18 @@ namespace GridTools template struct PeriodicFacePair; } + +namespace internal +{ + namespace parallel + { + namespace distributed + { + template + class TemporarilyMatchRefineFlags; + } + } // namespace parallel +} // namespace internal # endif namespace parallel @@ -803,6 +815,10 @@ namespace parallel template friend class dealii::FETools::internal::ExtrapolateImplementation; + + template + friend class dealii::internal::parallel::distributed:: + TemporarilyMatchRefineFlags; }; @@ -910,6 +926,10 @@ namespace parallel virtual types::coarse_cell_id coarse_cell_index_to_coarse_cell_id( const unsigned int coarse_cell_index) const override; + + template + friend class dealii::internal::parallel::distributed:: + TemporarilyMatchRefineFlags; }; } // namespace distributed } // namespace parallel @@ -934,7 +954,7 @@ namespace parallel */ template class Triangulation - : public dealii::parallel::TriangulationBase + : public dealii::parallel::DistributedTriangulationBase { public: /** @@ -950,6 +970,75 @@ namespace parallel #endif +namespace internal +{ + namespace parallel + { + namespace distributed + { + /** + * This class temporarily modifies the refine and coarsen flags of all + * active cells to match the p4est oracle. + * + * The modification only happens on parallel::distributed::Triangulation + * objects, and persists for the lifetime of an instantiation of this + * class. + * + * The TemporarilyMatchRefineFlags class should only be used in + * combination with the Triangulation::Signals::post_p4est_refinement + * signal. At this stage, the p4est orcale already has been refined, but + * the triangulation is still unchanged. After the modification, all + * refine and coarsen flags describe how the traingulation will acutally + * be refined. + */ + template + class TemporarilyMatchRefineFlags : public Subscriptor + { + public: + /** + * Constructor. + * + * Stores the refine and coarsen flags of all active cells if the + * provided Triangulation is of type + * parallel::distributed::Triangulation. + * + * Adjusts them to be consistent with the p4est oracle. + */ + TemporarilyMatchRefineFlags(Triangulation &tria); + + /** + * Destructor. + * + * Returns the refine and coarsen flags of all active cells on the + * parallel::distributed::Triangulation into their previous state. + */ + ~TemporarilyMatchRefineFlags(); + + private: + /** + * The modified parallel::distributed::Triangulation. + */ + const SmartPointer< + dealii::parallel::distributed::Triangulation> + distributed_tria; + + /** + * A vector that temporarily stores the refine flags before they have + * been modified on the parallel::distributed::Triangulation. + */ + std::vector saved_refine_flags; + + /** + * A vector that temporarily stores the coarsen flags before they have + * been modified on the parallel::distributed::Triangulation. + */ + std::vector saved_coarsen_flags; + }; + } // namespace distributed + } // namespace parallel +} // namespace internal + + DEAL_II_NAMESPACE_CLOSE #endif diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index 65b34777b0..09c21bf910 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -3559,6 +3559,92 @@ namespace parallel +namespace internal +{ + namespace parallel + { + namespace distributed + { + template + TemporarilyMatchRefineFlags::TemporarilyMatchRefineFlags( + Triangulation &tria) + : distributed_tria( + dynamic_cast< + dealii::parallel::distributed::Triangulation *>( + &tria)) + { +#ifdef DEAL_II_WITH_P4EST + if (distributed_tria != nullptr) + { + // Save the current set of refinement flags, and adjust the + // refinement flags to be consistent with the p4est oracle. + distributed_tria->save_coarsen_flags(saved_coarsen_flags); + distributed_tria->save_refine_flags(saved_refine_flags); + + for (const auto &pair : distributed_tria->local_cell_relations) + { + const auto &cell = pair.first; + const auto &status = pair.second; + + switch (status) + { + case dealii::Triangulation::CELL_PERSIST: + // cell remains unchanged + cell->clear_refine_flag(); + cell->clear_coarsen_flag(); + break; + + case dealii::Triangulation::CELL_REFINE: + // cell will be refined + cell->clear_coarsen_flag(); + cell->set_refine_flag(); + break; + + case dealii::Triangulation::CELL_COARSEN: + // children of this cell will be coarsened + for (const auto &child : cell->child_iterators()) + { + child->clear_refine_flag(); + child->set_coarsen_flag(); + } + break; + + case dealii::Triangulation::CELL_INVALID: + // do nothing as cell does not exist yet + break; + + default: + Assert(false, ExcInternalError()); + break; + } + } + } +#endif + } + + + + template + TemporarilyMatchRefineFlags::~TemporarilyMatchRefineFlags() + { +#ifdef DEAL_II_WITH_P4EST + if (distributed_tria) + { + // Undo the refinement flags modification. + distributed_tria->load_coarsen_flags(saved_coarsen_flags); + distributed_tria->load_refine_flags(saved_refine_flags); + } +#else + // pretend that this destructor does something to silence clang-tidy + (void)distributed_tria; +#endif + } + } // namespace distributed + } // namespace parallel +} // namespace internal + + + /*-------------- Explicit Instantiations -------------------------------*/ #include "tria.inst" diff --git a/source/distributed/tria.inst.in b/source/distributed/tria.inst.in index 141e63fc17..c0e4321358 100644 --- a/source/distributed/tria.inst.in +++ b/source/distributed/tria.inst.in @@ -15,19 +15,30 @@ -for (deal_II_dimension : DIMENSIONS) +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { namespace parallel \{ namespace distributed \{ - template class Triangulation; -#if deal_II_dimension < 3 - template class Triangulation; +#if deal_II_dimension <= deal_II_space_dimension + template class Triangulation; #endif -#if deal_II_dimension < 2 - template class Triangulation; + \} + \} + + namespace internal + \{ + namespace parallel + \{ + namespace distributed + \{ +#if deal_II_dimension <= deal_II_space_dimension + template class TemporarilyMatchRefineFlags; #endif + \} \} \} } diff --git a/tests/mpi/limit_p_level_difference_01.cc b/tests/mpi/limit_p_level_difference_01.cc index 9ca79b1e5f..a21913a9bd 100644 --- a/tests/mpi/limit_p_level_difference_01.cc +++ b/tests/mpi/limit_p_level_difference_01.cc @@ -75,13 +75,20 @@ test(const unsigned int fes_size, const unsigned int max_difference) if (cell->is_locally_owned() && cell->center() == Point()) cell->set_active_fe_index(sequence.back()); - const bool fe_indices_changed = - hp::Refinement::limit_p_level_difference(dofh, - max_difference, - contains_fe_index); + bool fe_indices_changed = false; + tria.signals.post_p4est_refinement.connect( + [&]() { + const internal::parallel::distributed::TemporarilyMatchRefineFlags + refine_modifier(tria); + fe_indices_changed = + hp::Refinement::limit_p_level_difference(dofh, + max_difference, + contains_fe_index); + }, + boost::signals2::at_front); + tria.execute_coarsening_and_refinement(); - (void)fe_indices_changed; Assert(fe_indices_changed, ExcInternalError()); // display number of cells for each FE index diff --git a/tests/mpi/limit_p_level_difference_02.cc b/tests/mpi/limit_p_level_difference_02.cc index f87bad8190..5c7cf6c1af 100644 --- a/tests/mpi/limit_p_level_difference_02.cc +++ b/tests/mpi/limit_p_level_difference_02.cc @@ -70,6 +70,18 @@ test(const unsigned int fes_size, const unsigned int max_difference) DoFHandler dofh(tria); dofh.distribute_dofs(fes); + bool fe_indices_changed; + tria.signals.post_p4est_refinement.connect( + [&]() { + const internal::parallel::distributed::TemporarilyMatchRefineFlags + refine_modifier(tria); + fe_indices_changed = + hp::Refinement::limit_p_level_difference(dofh, + max_difference, + contains_fe_index); + }, + boost::signals2::at_front); + // increase p-level in center subsequently for (unsigned int i = 0; i < sequence.size() - 1; ++i) { @@ -79,13 +91,9 @@ test(const unsigned int fes_size, const unsigned int max_difference) cell->set_future_fe_index( fes.next_in_hierarchy(cell->active_fe_index())); - const bool fe_indices_changed = - hp::Refinement::limit_p_level_difference(dofh, - max_difference, - contains_fe_index); + fe_indices_changed = false; tria.execute_coarsening_and_refinement(); - (void)fe_indices_changed; if (i >= max_difference) Assert(fe_indices_changed, ExcInternalError()); diff --git a/tests/mpi/limit_p_level_difference_03.cc b/tests/mpi/limit_p_level_difference_03.cc index 127e339b2a..8ac8ed9308 100644 --- a/tests/mpi/limit_p_level_difference_03.cc +++ b/tests/mpi/limit_p_level_difference_03.cc @@ -90,13 +90,20 @@ test(const unsigned int fes_size, const unsigned int max_difference) cell->set_active_fe_index(sequence.back()); dofh.distribute_dofs(fes); - const bool fe_indices_changed = - hp::Refinement::limit_p_level_difference(dofh, - max_difference, - contains_fe_index); + bool fe_indices_changed = false; + tria.signals.post_p4est_refinement.connect( + [&]() { + const internal::parallel::distributed::TemporarilyMatchRefineFlags + refine_modifier(tria); + fe_indices_changed = + hp::Refinement::limit_p_level_difference(dofh, + max_difference, + contains_fe_index); + }, + boost::signals2::at_front); + tria.execute_coarsening_and_refinement(); - (void)fe_indices_changed; Assert(fe_indices_changed, ExcInternalError()); #ifdef DEBUG diff --git a/tests/mpi/limit_p_level_difference_04.cc b/tests/mpi/limit_p_level_difference_04.cc index f3e2d1c5f4..07546f3b62 100644 --- a/tests/mpi/limit_p_level_difference_04.cc +++ b/tests/mpi/limit_p_level_difference_04.cc @@ -82,15 +82,22 @@ test(const unsigned int max_difference) cell->set_active_fe_index(sequence.back()); dofh.distribute_dofs(fes); - const bool fe_indices_changed = - hp::Refinement::limit_p_level_difference(dofh, - max_difference, - contains_fe_index); - (void)fe_indices_changed; - Assert(fe_indices_changed, ExcInternalError()); + bool fe_indices_changed = false; + tria.signals.post_p4est_refinement.connect( + [&]() { + const internal::parallel::distributed::TemporarilyMatchRefineFlags + refine_modifier(tria); + fe_indices_changed = + hp::Refinement::limit_p_level_difference(dofh, + max_difference, + contains_fe_index); + }, + boost::signals2::at_front); tria.execute_coarsening_and_refinement(); + Assert(fe_indices_changed, ExcInternalError()); + deallog << "active FE indices after adaptation:" << std::endl; for (const auto &cell : dofh.active_cell_iterators()) if (cell->is_locally_owned()) -- 2.39.5