]> https://gitweb.dealii.org/ - dealii.git/commitdiff
New class TemporarilyRestoreSubdomainIds. 11741/head
authorMarc Fehling <mafehling.git@gmail.com>
Fri, 12 Feb 2021 21:58:58 +0000 (14:58 -0700)
committerMarc Fehling <mafehling.git@gmail.com>
Wed, 24 Feb 2021 23:10:42 +0000 (16:10 -0700)
Temporarily changes subdomain IDs of cells on p:s:T to their true owner.

include/deal.II/base/aligned_vector.h
include/deal.II/distributed/shared_tria.h
source/distributed/shared_tria.cc
source/distributed/shared_tria.inst.in
source/dofs/dof_handler.cc
source/dofs/dof_handler_policy.cc
source/numerics/solution_transfer.cc

index 09e23fdcd404c1a58192809f12c4f5b9918fae1a..7ad95848ee8213571487be3eada840e6a1fbe968 100644 (file)
@@ -387,7 +387,7 @@ namespace internal
    * @relatesalso AlignedVector
    */
   template <typename T>
-  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 <typename T>
-  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 <typename T, bool initialize_memory>
-  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 <typename T, bool initialize_memory>
-  class AlignedVectorDefaultInitialize : private parallel::ParallelForInteger
+  class AlignedVectorDefaultInitialize
+    : private dealii::parallel::ParallelForInteger
   {
     static const std::size_t minimum_parallel_grain_size =
       160000 / sizeof(T) + 1;
index da4bdd967d7cd98c178441978920c5ff2aa0eaea..0d4efafe5e2442e11bd0121c8c9c3e26aad01462 100644 (file)
@@ -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 <int dim, int spacedim>
+      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<dim, spacedim> &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<dim, spacedim>>
+          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<unsigned int> saved_subdomain_ids;
+      };
+    } // namespace shared
+  }   // namespace parallel
+} // namespace internal
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 69a90bfe3e59ca219ab610738e57003675afddad..ac6a11d26b6f489a49e3da86b4748eec50033953 100644 (file)
@@ -468,6 +468,59 @@ namespace parallel
 #endif
 
 
+
+namespace internal
+{
+  namespace parallel
+  {
+    namespace shared
+    {
+      template <int dim, int spacedim>
+      TemporarilyRestoreSubdomainIds<dim, spacedim>::
+        TemporarilyRestoreSubdomainIds(const Triangulation<dim, spacedim> &tria)
+        : shared_tria(
+            dynamic_cast<
+              const dealii::parallel::shared::Triangulation<dim, spacedim> *>(
+              &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<types::subdomain_id> &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 <int dim, int spacedim>
+      TemporarilyRestoreSubdomainIds<dim, spacedim>::
+        ~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"
 
index 32100be23954133f142bf2a17c069e06924f39fe..4c8bd1661f60bcc4f3754984631801623c5ed9b3 100644 (file)
 
 
 
-for (deal_II_dimension : DIMENSIONS)
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
   {
     namespace parallel
     \{
       namespace shared
       \{
-        template class Triangulation<deal_II_dimension>;
-#if deal_II_dimension < 3
-        template class Triangulation<deal_II_dimension, deal_II_dimension + 1>;
+#if deal_II_dimension <= deal_II_space_dimension
+        template class Triangulation<deal_II_dimension,
+                                     deal_II_space_dimension>;
 #endif
-#if deal_II_dimension < 2
-        template class Triangulation<deal_II_dimension, deal_II_dimension + 2>;
+      \}
+    \}
+
+    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
+        \}
       \}
     \}
   }
index 6abec402085590be61589de92c94c7f91bf56a0d..a5cd353acaebb5a8026a3e2f0f8a915f3d778ef0 100644 (file)
@@ -2559,47 +2559,31 @@ DoFHandler<dim, spacedim>::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<types::subdomain_id>                      saved_subdomain_ids;
-  const parallel::shared::Triangulation<dim, spacedim> *shared_tria =
-    (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim> *>(
-      &this->get_triangulation()));
-  if (shared_tria != nullptr && shared_tria->with_artificial_cells())
-    {
-      saved_subdomain_ids.resize(shared_tria->n_active_cells());
-
-      const std::vector<types::subdomain_id> &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<dim,
+                                                                     spacedim>
+      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();
index 76de280b8a4934dc31744b8a0f9a11b7090f4320..e9f07561ccca071cfd467c5967bbd4e02395a356 100644 (file)
@@ -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<types::subdomain_id> saved_subdomain_ids;
-        if (tr->with_artificial_cells())
-          {
-            saved_subdomain_ids.resize(tr->n_active_cells());
-
-            const std::vector<types::subdomain_id> &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<dim, spacedim>
+            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<dim, spacedim>::
-          active_cell_iterator
-            cell = this->dof_handler->get_triangulation().begin_active(),
-            endc = this->dof_handler->get_triangulation().end();
-        std::vector<types::subdomain_id> current_subdomain_ids(
-          tr->n_active_cells());
-        const std::vector<types::subdomain_id> &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<dim, spacedim>
+            subdomain_modifier(*tr);
 
         std::vector<types::global_dof_index> 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
       }
index dacfb341baad14ade745b09155949c46e108d57f..c67ae4dc419e8ce770e9a1d417008ffbcd0c98e3 100644 (file)
@@ -94,30 +94,16 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::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<types::subdomain_id> saved_subdomain_ids;
-  const parallel::shared::Triangulation<dim, DoFHandlerType::space_dimension>
-    *shared_tria =
-      (dynamic_cast<const parallel::shared::
-                      Triangulation<dim, DoFHandlerType::space_dimension> *>(
-        &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<types::subdomain_id> &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<dim, DoFHandlerType::space_dimension>
+      subdomain_modifier(dof_handler->get_triangulation());
 
   const unsigned int n_active_cells =
     dof_handler->get_triangulation().n_active_cells();
@@ -142,11 +128,6 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::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<dim, VectorType, DoFHandlerType>::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<types::subdomain_id> saved_subdomain_ids;
-  const parallel::shared::Triangulation<dim, DoFHandlerType::space_dimension>
-    *shared_tria =
-      (dynamic_cast<const parallel::shared::
-                      Triangulation<dim, DoFHandlerType::space_dimension> *>(
-        &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<types::subdomain_id> &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<dim, DoFHandlerType::space_dimension>
+      subdomain_modifier(dof_handler->get_triangulation());
 
   Vector<typename VectorType::value_type> local_values(0);
 
@@ -229,11 +196,6 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::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<dim, VectorType, DoFHandlerType>::
 #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<types::subdomain_id> saved_subdomain_ids;
-  const parallel::shared::Triangulation<dim, DoFHandlerType::space_dimension>
-    *shared_tria =
-      (dynamic_cast<const parallel::shared::
-                      Triangulation<dim, DoFHandlerType::space_dimension> *>(
-        &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<types::subdomain_id> &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<dim, DoFHandlerType::space_dimension>
+      subdomain_modifier(dof_handler->get_triangulation());
 
   // first count the number
   // of cells that will be coarsened
@@ -480,11 +428,6 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::
   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<dim, VectorType, DoFHandlerType>::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<types::subdomain_id> saved_subdomain_ids;
-  const parallel::shared::Triangulation<dim, DoFHandlerType::space_dimension>
-    *shared_tria =
-      (dynamic_cast<const parallel::shared::
-                      Triangulation<dim, DoFHandlerType::space_dimension> *>(
-        &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<types::subdomain_id> &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<dim, DoFHandlerType::space_dimension>
+      subdomain_modifier(dof_handler->get_triangulation());
 
   Vector<typename VectorType::value_type> local_values;
   std::vector<types::global_dof_index>    dofs;
@@ -654,11 +583,6 @@ SolutionTransfer<dim, VectorType, DoFHandlerType>::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()]);
 }
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.