From 7ade5f019f77c27378ac6ef4fe04acf9d539bdc9 Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Fri, 12 Feb 2021 14:58:58 -0700 Subject: [PATCH] New class TemporarilyRestoreSubdomainIds. Temporarily changes subdomain IDs of cells on p:s:T to their true owner. --- include/deal.II/base/aligned_vector.h | 9 +- include/deal.II/distributed/shared_tria.h | 67 ++++++++++ source/distributed/shared_tria.cc | 53 ++++++++ source/distributed/shared_tria.inst.in | 24 +++- source/dofs/dof_handler.cc | 66 ++++----- source/dofs/dof_handler_policy.cc | 59 ++------ source/numerics/solution_transfer.cc | 156 ++++++---------------- 7 files changed, 222 insertions(+), 212 deletions(-) diff --git a/include/deal.II/base/aligned_vector.h b/include/deal.II/base/aligned_vector.h index 09e23fdcd4..7ad95848ee 100644 --- a/include/deal.II/base/aligned_vector.h +++ b/include/deal.II/base/aligned_vector.h @@ -387,7 +387,7 @@ namespace internal * @relatesalso AlignedVector */ template - class AlignedVectorCopy : private parallel::ParallelForInteger + class AlignedVectorCopy : private dealii::parallel::ParallelForInteger { static const std::size_t minimum_parallel_grain_size = 160000 / sizeof(T) + 1; @@ -454,7 +454,7 @@ namespace internal * @relatesalso AlignedVector */ template - class AlignedVectorMove : private parallel::ParallelForInteger + class AlignedVectorMove : private dealii::parallel::ParallelForInteger { static const std::size_t minimum_parallel_grain_size = 160000 / sizeof(T) + 1; @@ -532,7 +532,7 @@ namespace internal * @relatesalso AlignedVector */ template - class AlignedVectorSet : private parallel::ParallelForInteger + class AlignedVectorSet : private dealii::parallel::ParallelForInteger { static const std::size_t minimum_parallel_grain_size = 160000 / sizeof(T) + 1; @@ -634,7 +634,8 @@ namespace internal * @relatesalso AlignedVector */ template - class AlignedVectorDefaultInitialize : private parallel::ParallelForInteger + class AlignedVectorDefaultInitialize + : private dealii::parallel::ParallelForInteger { static const std::size_t minimum_parallel_grain_size = 160000 / sizeof(T) + 1; diff --git a/include/deal.II/distributed/shared_tria.h b/include/deal.II/distributed/shared_tria.h index da4bdd967d..0d4efafe5e 100644 --- a/include/deal.II/distributed/shared_tria.h +++ b/include/deal.II/distributed/shared_tria.h @@ -460,6 +460,73 @@ namespace parallel #endif } // namespace parallel + +namespace internal +{ + namespace parallel + { + namespace shared + { + /** + * This class temporarily modifies the subdomain ID of all active cells to + * their respective "true" owner. + * + * The modification only happens on parallel::shared::Triangulation + * objects with artificial cells, and persists for the lifetime of an + * instantiation of this class. + * + * The TemporarilyRestoreSubdomainIds class should only be used for + * temporary read-only purposes. For example, whenever your implementation + * requires to treat artificial cells temporarily as locally relevant to + * access their dof indices. + * + * This class has effect only if artificial cells are allowed. Without + * artificial cells, the current subdomain IDs already correspond to the + * true subdomain IDs. See the @ref GlossArtificialCell "glossary" + * for more information about artificial cells. + */ + template + class TemporarilyRestoreSubdomainIds : public Subscriptor + { + public: + /** + * Constructor. + * + * Stores the subdomain ID of all active cells if the provided + * Triangulation is of type parallel::shared::Triangulation. + * + * Replaces them by their true subdomain ID equivalent. + */ + TemporarilyRestoreSubdomainIds( + const Triangulation &tria); + + /** + * Destructor. + * + * Returns the subdomain ID of all active cells on the + * parallel::shared::Triangulation into their previous state. + */ + ~TemporarilyRestoreSubdomainIds(); + + private: + /** + * The modified parallel::shared::Triangulation. + */ + SmartPointer< + const dealii::parallel::shared::Triangulation> + shared_tria; + + /** + * A vector that temporarily stores the subdomain IDs on all active + * cells before they have been modified on the + * parallel::shared::Triangulation. + */ + std::vector saved_subdomain_ids; + }; + } // namespace shared + } // namespace parallel +} // namespace internal + DEAL_II_NAMESPACE_CLOSE #endif diff --git a/source/distributed/shared_tria.cc b/source/distributed/shared_tria.cc index 69a90bfe3e..ac6a11d26b 100644 --- a/source/distributed/shared_tria.cc +++ b/source/distributed/shared_tria.cc @@ -468,6 +468,59 @@ namespace parallel #endif + +namespace internal +{ + namespace parallel + { + namespace shared + { + template + TemporarilyRestoreSubdomainIds:: + TemporarilyRestoreSubdomainIds(const Triangulation &tria) + : shared_tria( + dynamic_cast< + const dealii::parallel::shared::Triangulation *>( + &tria)) + { + if (shared_tria && shared_tria->with_artificial_cells()) + { + // Save the current set of subdomain IDs, and set subdomain IDs + // to the "true" owner of each cell. + const std::vector &true_subdomain_ids = + shared_tria->get_true_subdomain_ids_of_cells(); + + saved_subdomain_ids.resize(shared_tria->n_active_cells()); + for (const auto &cell : shared_tria->active_cell_iterators()) + { + const unsigned int index = cell->active_cell_index(); + saved_subdomain_ids[index] = cell->subdomain_id(); + cell->set_subdomain_id(true_subdomain_ids[index]); + } + } + } + + + + template + TemporarilyRestoreSubdomainIds:: + ~TemporarilyRestoreSubdomainIds() + { + if (shared_tria && shared_tria->with_artificial_cells()) + { + // Undo the subdomain modification. + for (const auto &cell : shared_tria->active_cell_iterators()) + { + const unsigned int index = cell->active_cell_index(); + cell->set_subdomain_id(saved_subdomain_ids[index]); + } + } + } + } // namespace shared + } // namespace parallel +} // namespace internal + + /*-------------- Explicit Instantiations -------------------------------*/ #include "shared_tria.inst" diff --git a/source/distributed/shared_tria.inst.in b/source/distributed/shared_tria.inst.in index 32100be239..4c8bd1661f 100644 --- a/source/distributed/shared_tria.inst.in +++ b/source/distributed/shared_tria.inst.in @@ -15,19 +15,31 @@ -for (deal_II_dimension : DIMENSIONS) +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) { namespace parallel \{ namespace shared \{ - 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 shared + \{ +#if deal_II_dimension <= deal_II_space_dimension + template class TemporarilyRestoreSubdomainIds< + deal_II_dimension, + deal_II_space_dimension>; #endif + \} \} \} } diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index 6abec40208..a5cd353aca 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -2559,47 +2559,31 @@ DoFHandler::distribute_dofs( #endif } - // We would like to enumerate all dofs for shared::Triangulations. If an - // underlying shared::Tria allows artificial cells, we need to restore the - // true cell owners temporarily. We save the current set of subdomain ids, set - // subdomain ids to the "true" owner of each cell, and later restore these - // flags. - std::vector saved_subdomain_ids; - const parallel::shared::Triangulation *shared_tria = - (dynamic_cast *>( - &this->get_triangulation())); - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - { - saved_subdomain_ids.resize(shared_tria->n_active_cells()); - - const std::vector &true_subdomain_ids = - shared_tria->get_true_subdomain_ids_of_cells(); - - for (const auto &cell : shared_tria->active_cell_iterators()) - { - const unsigned int index = cell->active_cell_index(); - saved_subdomain_ids[index] = cell->subdomain_id(); - cell->set_subdomain_id(true_subdomain_ids[index]); - } - } - - // Adjust size of levels to the triangulation. Note that we still have to - // allocate space for all degrees of freedom on this mesh (including ghost - // and cells that are entirely stored on different processors), though we - // may not assign numbers to some of them (i.e. they will remain at - // invalid_dof_index). We need to allocate the space because we will want - // to be able to query the dof_indices on each cell, and simply be told - // that we don't know them on some cell (i.e. get back invalid_dof_index) - if (hp_capability_enabled) - internal::hp::DoFHandlerImplementation::Implementation::reserve_space( - *this); - else - internal::DoFHandlerImplementation::Implementation::reserve_space(*this); - - // undo the subdomain modification - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - for (const auto &cell : shared_tria->active_cell_iterators()) - cell->set_subdomain_id(saved_subdomain_ids[cell->active_cell_index()]); + { + // We would like to enumerate all dofs for shared::Triangulations. If an + // underlying shared::Tria allows artificial cells, we need to restore the + // true cell owners temporarily. + // We use the TemporarilyRestoreSubdomainIds class for this purpose: we save + // the current set of subdomain ids, set subdomain ids to the "true" owner + // of each cell upon construction of the TemporarilyRestoreSubdomainIds + // object, and later restore these flags when it is destroyed. + const internal::parallel::shared::TemporarilyRestoreSubdomainIds + subdomain_modifier(this->get_triangulation()); + + // Adjust size of levels to the triangulation. Note that we still have to + // allocate space for all degrees of freedom on this mesh (including ghost + // and cells that are entirely stored on different processors), though we + // may not assign numbers to some of them (i.e. they will remain at + // invalid_dof_index). We need to allocate the space because we will want + // to be able to query the dof_indices on each cell, and simply be told + // that we don't know them on some cell (i.e. get back invalid_dof_index) + if (hp_capability_enabled) + internal::hp::DoFHandlerImplementation::Implementation::reserve_space( + *this); + else + internal::DoFHandlerImplementation::Implementation::reserve_space(*this); + } // hand the actual work over to the policy this->number_cache = this->policy->distribute_dofs(); diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index 76de280b8a..e9f07561cc 100644 --- a/source/dofs/dof_handler_policy.cc +++ b/source/dofs/dof_handler_policy.cc @@ -3022,25 +3022,16 @@ namespace internal const unsigned int n_procs = Utilities::MPI::n_mpi_processes(tr->get_communicator()); - // If the underlying shared::Tria allows artificial cells, - // then save the current set of subdomain ids, and set - // subdomain ids to the "true" owner of each cell. we later - // restore these flags - std::vector saved_subdomain_ids; - if (tr->with_artificial_cells()) - { - saved_subdomain_ids.resize(tr->n_active_cells()); - - const std::vector &true_subdomain_ids = - tr->get_true_subdomain_ids_of_cells(); - - for (const auto &cell : tr->active_cell_iterators()) - { - const unsigned int index = cell->active_cell_index(); - saved_subdomain_ids[index] = cell->subdomain_id(); - cell->set_subdomain_id(true_subdomain_ids[index]); - } - } + // If an underlying shared::Tria allows artificial cells, we need to + // restore the true cell owners temporarily. + // We use the TemporarilyRestoreSubdomainIds class for this purpose: we + // save the current set of subdomain ids, set subdomain ids to the + // "true" owner of each cell upon construction of the + // TemporarilyRestoreSubdomainIds object, and later restore these flags + // when it is destroyed. + const internal::parallel::shared:: + TemporarilyRestoreSubdomainIds + subdomain_modifier(*tr); // first let the sequential algorithm do its magic. it is going to // enumerate DoFs on all cells, regardless of owner @@ -3154,12 +3145,6 @@ namespace internal } } - // finally, restore current subdomain ids - if (tr->with_artificial_cells()) - for (const auto &cell : tr->active_cell_iterators()) - cell->set_subdomain_id( - saved_subdomain_ids[cell->active_cell_index()]); - // return a NumberCache object made up from the sets of locally // owned DoFs return NumberCache( @@ -3375,20 +3360,10 @@ namespace internal &this->dof_handler->get_triangulation())); Assert(tr != nullptr, ExcInternalError()); - typename dealii::parallel::shared::Triangulation:: - active_cell_iterator - cell = this->dof_handler->get_triangulation().begin_active(), - endc = this->dof_handler->get_triangulation().end(); - std::vector current_subdomain_ids( - tr->n_active_cells()); - const std::vector &true_subdomain_ids = - tr->get_true_subdomain_ids_of_cells(); - if (tr->with_artificial_cells()) - for (unsigned int index = 0; cell != endc; cell++, index++) - { - current_subdomain_ids[index] = cell->subdomain_id(); - cell->set_subdomain_id(true_subdomain_ids[index]); - } + // Set subdomain IDs to the "true" owner of each cell. + const internal::parallel::shared:: + TemporarilyRestoreSubdomainIds + subdomain_modifier(*tr); std::vector global_gathered_numbers( this->dof_handler->n_dofs(), 0); @@ -3511,12 +3486,6 @@ namespace internal DoFTools::locally_owned_dofs_per_subdomain(*this->dof_handler), this->dof_handler->get_triangulation().locally_owned_subdomain()); - // restore artificial cells - cell = tr->begin_active(); - if (tr->with_artificial_cells()) - for (unsigned int index = 0; cell != endc; cell++, index++) - cell->set_subdomain_id(current_subdomain_ids[index]); - return number_cache; #endif } diff --git a/source/numerics/solution_transfer.cc b/source/numerics/solution_transfer.cc index dacfb341ba..c67ae4dc41 100644 --- a/source/numerics/solution_transfer.cc +++ b/source/numerics/solution_transfer.cc @@ -94,30 +94,16 @@ SolutionTransfer::prepare_for_pure_refinement() clear(); // We need to access dof indices on the entire domain. For - // shared::Triangulations, ownership of cells might change. If they allow the - // use of artificial cells, we need to restore the true cell owners - // temporarily. We save the current set of subdomain ids, set subdomain ids to - // the "true" owner of each cell, and later restore these flags. - std::vector saved_subdomain_ids; - const parallel::shared::Triangulation - *shared_tria = - (dynamic_cast *>( - &dof_handler->get_triangulation())); - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - { - saved_subdomain_ids.resize(shared_tria->n_active_cells()); - - const std::vector &true_subdomain_ids = - shared_tria->get_true_subdomain_ids_of_cells(); - - for (const auto &cell : shared_tria->active_cell_iterators()) - { - const unsigned int index = cell->active_cell_index(); - saved_subdomain_ids[index] = cell->subdomain_id(); - cell->set_subdomain_id(true_subdomain_ids[index]); - } - } + // parallel::shared::Triangulations, ownership of cells might change. If they + // allow artificial cells, we need to restore the "true" cell owners + // temporarily. + // We use the TemporarilyRestoreSubdomainIds class for this purpose: we save + // the current set of subdomain ids, set subdomain ids to the "true" owner of + // each cell upon construction of the TemporarilyRestoreSubdomainIds object, + // and later restore these flags when it is destroyed. + const internal::parallel::shared:: + TemporarilyRestoreSubdomainIds + subdomain_modifier(dof_handler->get_triangulation()); const unsigned int n_active_cells = dof_handler->get_triangulation().n_active_cells(); @@ -142,11 +128,6 @@ SolutionTransfer::prepare_for_pure_refinement() Pointerstruct(&indices_on_cell[i], cell->active_fe_index()); } prepared_for = pure_refinement; - - // Undo the subdomain modification. - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - for (const auto &cell : shared_tria->active_cell_iterators()) - cell->set_subdomain_id(saved_subdomain_ids[cell->active_cell_index()]); } @@ -166,30 +147,16 @@ SolutionTransfer::refine_interpolate( " at the same time!")); // We need to access dof indices on the entire domain. For - // shared::Triangulations, ownership of cells might change. If they allow the - // use of artificial cells, we need to restore the true cell owners - // temporarily. We save the current set of subdomain ids, set subdomain ids to - // the "true" owner of each cell, and later restore these flags. - std::vector saved_subdomain_ids; - const parallel::shared::Triangulation - *shared_tria = - (dynamic_cast *>( - &dof_handler->get_triangulation())); - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - { - saved_subdomain_ids.resize(shared_tria->n_active_cells()); - - const std::vector &true_subdomain_ids = - shared_tria->get_true_subdomain_ids_of_cells(); - - for (const auto &cell : shared_tria->active_cell_iterators()) - { - const unsigned int index = cell->active_cell_index(); - saved_subdomain_ids[index] = cell->subdomain_id(); - cell->set_subdomain_id(true_subdomain_ids[index]); - } - } + // parallel::shared::Triangulations, ownership of cells might change. If they + // allow artificial cells, we need to restore the "true" cell owners + // temporarily. + // We use the TemporarilyRestoreSubdomainIds class for this purpose: we save + // the current set of subdomain ids, set subdomain ids to the "true" owner of + // each cell upon construction of the TemporarilyRestoreSubdomainIds object, + // and later restore these flags when it is destroyed. + const internal::parallel::shared:: + TemporarilyRestoreSubdomainIds + subdomain_modifier(dof_handler->get_triangulation()); Vector local_values(0); @@ -229,11 +196,6 @@ SolutionTransfer::refine_interpolate( this_fe_index); } } - - // Undo the subdomain modification. - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - for (const auto &cell : shared_tria->active_cell_iterators()) - cell->set_subdomain_id(saved_subdomain_ids[cell->active_cell_index()]); } @@ -338,30 +300,16 @@ SolutionTransfer:: #endif // We need to access dof indices on the entire domain. For - // shared::Triangulations, ownership of cells might change. If they allow the - // use of artificial cells, we need to restore the true cell owners - // temporarily. We save the current set of subdomain ids, set subdomain ids to - // the "true" owner of each cell, and later restore these flags. - std::vector saved_subdomain_ids; - const parallel::shared::Triangulation - *shared_tria = - (dynamic_cast *>( - &dof_handler->get_triangulation())); - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - { - saved_subdomain_ids.resize(shared_tria->n_active_cells()); - - const std::vector &true_subdomain_ids = - shared_tria->get_true_subdomain_ids_of_cells(); - - for (const auto &cell : shared_tria->active_cell_iterators()) - { - const unsigned int index = cell->active_cell_index(); - saved_subdomain_ids[index] = cell->subdomain_id(); - cell->set_subdomain_id(true_subdomain_ids[index]); - } - } + // parallel::shared::Triangulations, ownership of cells might change. If they + // allow artificial cells, we need to restore the "true" cell owners + // temporarily. + // We use the TemporarilyRestoreSubdomainIds class for this purpose: we save + // the current set of subdomain ids, set subdomain ids to the "true" owner of + // each cell upon construction of the TemporarilyRestoreSubdomainIds object, + // and later restore these flags when it is destroyed. + const internal::parallel::shared:: + TemporarilyRestoreSubdomainIds + subdomain_modifier(dof_handler->get_triangulation()); // first count the number // of cells that will be coarsened @@ -480,11 +428,6 @@ SolutionTransfer:: Assert(n_cf == n_coarsen_fathers, ExcInternalError()); prepared_for = coarsening_and_refinement; - - // Undo the subdomain modification. - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - for (const auto &cell : shared_tria->active_cell_iterators()) - cell->set_subdomain_id(saved_subdomain_ids[cell->active_cell_index()]); } @@ -524,30 +467,16 @@ SolutionTransfer::interpolate( #endif // We need to access dof indices on the entire domain. For - // shared::Triangulations, ownership of cells might change. If they allow the - // use of artificial cells, we need to restore the true cell owners - // temporarily. We save the current set of subdomain ids, set subdomain ids to - // the "true" owner of each cell, and later restore these flags. - std::vector saved_subdomain_ids; - const parallel::shared::Triangulation - *shared_tria = - (dynamic_cast *>( - &dof_handler->get_triangulation())); - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - { - saved_subdomain_ids.resize(shared_tria->n_active_cells()); - - const std::vector &true_subdomain_ids = - shared_tria->get_true_subdomain_ids_of_cells(); - - for (const auto &cell : shared_tria->active_cell_iterators()) - { - const unsigned int index = cell->active_cell_index(); - saved_subdomain_ids[index] = cell->subdomain_id(); - cell->set_subdomain_id(true_subdomain_ids[index]); - } - } + // parallel::shared::Triangulations, ownership of cells might change. If they + // allow artificial cells, we need to restore the "true" cell owners + // temporarily. + // We use the TemporarilyRestoreSubdomainIds class for this purpose: we save + // the current set of subdomain ids, set subdomain ids to the "true" owner of + // each cell upon construction of the TemporarilyRestoreSubdomainIds object, + // and later restore these flags when it is destroyed. + const internal::parallel::shared:: + TemporarilyRestoreSubdomainIds + subdomain_modifier(dof_handler->get_triangulation()); Vector local_values; std::vector dofs; @@ -654,11 +583,6 @@ SolutionTransfer::interpolate( Assert(false, ExcInternalError()); } } - - // Undo the subdomain modification. - if (shared_tria != nullptr && shared_tria->with_artificial_cells()) - for (const auto &cell : shared_tria->active_cell_iterators()) - cell->set_subdomain_id(saved_subdomain_ids[cell->active_cell_index()]); } -- 2.39.5