Restore true cell owners before setting up DoF data structures.

This is necessary when using a parallel::shared::Triangulation because
there, we want to enumerate DoFs on *all* cells, not just some, and so
we need to know who owns them etc and what the active_fe_index is.

diff --git a/source/hp/dof_handler.cc b/source/hp/dof_handler.cc
index a1d12159a4746056a68d2348910431e7f5f44520..1a18ad6e0c6a3347a9118454f7f363633fc00a2c 100644 (file)
--- a/source/hp/dof_handler.cc
+++ b/source/hp/dof_handler.cc
@@ -1468,7 +1468,32 @@ namespace hp
     // active_fe_indices on both its own cells and all ghost cells
     communicate_active_fe_indices (*this);
 
-    // This call ensures that the active_fe_indices vectors are
+    // If an 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 (const parallel::shared::Triangulation<dim, spacedim> *shared_tria =
+          (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim>*> (&get_triangulation())))
+      if (shared_tria->with_artificial_cells())
+        {
+          saved_subdomain_ids.resize (shared_tria->n_active_cells());
+
+          typename parallel::shared::Triangulation<dim,spacedim>::active_cell_iterator
+          cell = get_triangulation().begin_active(),
+          endc = get_triangulation().end();
+
+          const std::vector<types::subdomain_id> &true_subdomain_ids
+            = shared_tria->get_true_subdomain_ids_of_cells();
+
+          for (unsigned int index=0; cell != endc; ++cell, ++index)
+            {
+              saved_subdomain_ids[index] = cell->subdomain_id();
+              cell->set_subdomain_id(true_subdomain_ids[index]);
+            }
+        }
+
+    // ensure that the active_fe_indices vectors are
     // initialized correctly.
     create_active_fe_table ();
 
@@ -1483,6 +1508,21 @@ namespace hp
     // then allocate space for all the other tables
     dealii::internal::hp::DoFHandler::Implementation::reserve_space (*this);
 
+
+    // now undo the subdomain modification
+    if (const parallel::shared::Triangulation<dim, spacedim> *shared_tria =
+          (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim>*> (&get_triangulation())))
+      if (shared_tria->with_artificial_cells())
+        {
+          typename parallel::shared::Triangulation<dim,spacedim>::active_cell_iterator
+          cell = get_triangulation().begin_active(),
+          endc = get_triangulation().end();
+
+          for (unsigned int index=0; cell != endc; ++cell, ++index)
+            cell->set_subdomain_id(saved_subdomain_ids[index]);
+        }
+
+
     // Clear user flags because we will need them. But first we save
     // them and make sure that we restore them later such that at the
     // end of this function the Triangulation will be in the same
@@ -1499,12 +1539,7 @@ namespace hp
 
     /////////////////////////////////
 
-
-    Assert ((dynamic_cast<const parallel::distributed::Triangulation< dim, spacedim >*>
-             (&this->get_triangulation())
-             == nullptr),
-            ExcNotImplemented());
-
+    // do some housekeeping: compress indices
     {
       Threads::TaskGroup<> tg;
       for (int level=levels.size()-1; level>=0; --level)