From: Martin Kronbichler <martin.kronbichler@uni-a.de>
Date: Sun, 29 May 2022 20:29:11 +0000 (+0200)
Subject: Avoid some uses of Triangulation::user_flags
X-Git-Tag: v9.5.0-rc1~1166^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=38c1c5262afbf594f1fa9f9cf42d502461d6c060;p=dealii.git

Avoid some uses of Triangulation::user_flags
---

diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc
index c7aaa49023..5726e3f1d1 100644
--- a/source/dofs/dof_handler.cc
+++ b/source/dofs/dof_handler.cc
@@ -921,43 +921,14 @@ namespace internal
           // *have* to be different, and so we need to prepare for this
           // as well.
           //
-          // The way we do things is that we loop over all active
-          // cells (these are the only ones that have DoFs
-          // anyway) and all their faces. We note in the
-          // user flags whether we have previously visited a face and
-          // if so skip it (consequently, we have to save and later
-          // restore the face flags)
+          // The way we do things is that we loop over all active cells (these
+          // are the only ones that have DoFs anyway) and all their faces. We
+          // note in the vector face_touched whether we have previously
+          // visited a face and if so skip it
           {
-            std::vector<bool> saved_face_user_flags;
-            switch (dim)
-              {
-                case 2:
-                  {
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .save_user_flags_line(saved_face_user_flags);
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .clear_user_flags_line();
-
-                    break;
-                  }
-
-                case 3:
-                  {
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .save_user_flags_quad(saved_face_user_flags);
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .clear_user_flags_quad();
-
-                    break;
-                  }
-
-                default:
-                  Assert(false, ExcNotImplemented());
-              }
+            std::vector<bool> face_touched(dim == 2 ?
+                                             dof_handler.tria->n_raw_lines() :
+                                             dof_handler.tria->n_raw_quads());
 
             const unsigned int d = dim - 1;
             const unsigned int l = 0;
@@ -977,7 +948,7 @@ namespace internal
             for (const auto &cell : dof_handler.active_cell_iterators())
               if (!cell->is_artificial())
                 for (const auto face : cell->face_indices())
-                  if (cell->face(face)->user_flag_set() == false)
+                  if (!face_touched[cell->face(face)->index()])
                     {
                       unsigned int fe_slots_needed = 0;
 
@@ -1009,7 +980,7 @@ namespace internal
                         }
 
                       // mark this face as visited
-                      cell->face(face)->set_user_flag();
+                      face_touched[cell->face(face)->index()] = true;
 
                       dof_handler
                         .hp_object_fe_ptr[d][cell->face(face)->index() + 1] =
@@ -1033,34 +1004,12 @@ namespace internal
             // With the memory now allocated, loop over the
             // dof_handler cells again and prime the _offset values as
             // well as the fe_index fields
-            switch (dim)
-              {
-                case 2:
-                  {
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .clear_user_flags_line();
-
-                    break;
-                  }
-
-                case 3:
-                  {
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .clear_user_flags_quad();
-
-                    break;
-                  }
-
-                default:
-                  Assert(false, ExcNotImplemented());
-              }
+            face_touched = std::vector<bool>(face_touched.size());
 
             for (const auto &cell : dof_handler.active_cell_iterators())
               if (!cell->is_artificial())
                 for (const auto face : cell->face_indices())
-                  if (!cell->face(face)->user_flag_set())
+                  if (!face_touched[cell->face(face)->index()])
                     {
                       // Same decision tree as before
                       if (cell->at_boundary(face) ||
@@ -1135,7 +1084,7 @@ namespace internal
                         }
 
                       // mark this face as visited
-                      cell->face(face)->set_user_flag();
+                      face_touched[cell->face(face)->index()] = true;
                     }
 
             for (unsigned int i = 1;
@@ -1143,31 +1092,6 @@ namespace internal
                  i++)
               dof_handler.object_dof_ptr[l][d][i] +=
                 dof_handler.object_dof_ptr[l][d][i - 1];
-
-            // at the end, restore the user flags for the faces
-            switch (dim)
-              {
-                case 2:
-                  {
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .load_user_flags_line(saved_face_user_flags);
-
-                    break;
-                  }
-
-                case 3:
-                  {
-                    const_cast<dealii::Triangulation<dim, spacedim> &>(
-                      *dof_handler.tria)
-                      .load_user_flags_quad(saved_face_user_flags);
-
-                    break;
-                  }
-
-                default:
-                  Assert(false, ExcNotImplemented());
-              }
           }
         }
 
diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc
index 91334cf40d..da4b13eefc 100644
--- a/source/dofs/dof_handler_policy.cc
+++ b/source/dofs/dof_handler_policy.cc
@@ -446,14 +446,6 @@ namespace internal
           std::map<types::global_dof_index, types::global_dof_index>
             dof_identities;
 
-          // we will mark lines that we have already treated, so first save and
-          // clear the user flags on lines and later restore them
-          std::vector<bool> user_flags;
-          dof_handler.get_triangulation().save_user_flags_line(user_flags);
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .clear_user_flags_line();
-
           // An implementation of the algorithm described in the hp-paper,
           // including the modification mentioned later in the "complications in
           // 3-d" subsections
@@ -476,12 +468,14 @@ namespace internal
           dealii::Table<2, std::unique_ptr<DoFIdentities>> line_dof_identities(
             dof_handler.fe_collection.size(), dof_handler.fe_collection.size());
 
+          std::vector<bool> line_touched(
+            dof_handler.get_triangulation().n_raw_lines());
           for (const auto &cell : dof_handler.active_cell_iterators())
             for (const auto l : cell->line_indices())
-              if (cell->line(l)->user_flag_set() == false)
+              if (!line_touched[cell->line(l)->index()])
                 {
-                  const auto line = cell->line(l);
-                  line->set_user_flag();
+                  const auto line             = cell->line(l);
+                  line_touched[line->index()] = true;
 
                   unsigned int unique_sets_of_dofs =
                     line->n_active_fe_indices();
@@ -774,11 +768,6 @@ namespace internal
                     }
                 }
 
-          // finally restore the user flags
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .load_user_flags_line(user_flags);
-
           return dof_identities;
         }
 
@@ -819,16 +808,6 @@ namespace internal
           std::map<types::global_dof_index, types::global_dof_index>
             dof_identities;
 
-
-          // we will mark quads that we have already treated, so first
-          // save and clear the user flags on quads and later restore
-          // them
-          std::vector<bool> user_flags;
-          dof_handler.get_triangulation().save_user_flags_quad(user_flags);
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .clear_user_flags_quad();
-
           // An implementation of the algorithm described in the hp-
           // paper, including the modification mentioned later in the
           // "complications in 3-d" subsections
@@ -845,13 +824,15 @@ namespace internal
             dof_handler.fe_collection.size(),
             2 /*triangle (0) or quadrilateral (1)*/);
 
+          std::vector<bool> quad_touched(
+            dof_handler.get_triangulation().n_raw_quads());
           for (const auto &cell : dof_handler.active_cell_iterators())
             for (const auto q : cell->face_indices())
-              if ((cell->quad(q)->user_flag_set() == false) &&
+              if (!quad_touched[cell->quad(q)->index()] &&
                   (cell->quad(q)->n_active_fe_indices() == 2))
                 {
-                  const auto quad = cell->quad(q);
-                  quad->set_user_flag();
+                  const auto quad             = cell->quad(q);
+                  quad_touched[quad->index()] = true;
 
                   const std::set<unsigned int> fe_indices =
                     quad->get_active_fe_indices();
@@ -938,11 +919,6 @@ namespace internal
                     }
                 }
 
-          // finally restore the user flags
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .load_user_flags_quad(user_flags);
-
           return dof_identities;
         }
 
@@ -1129,7 +1105,8 @@ namespace internal
             vertex_dof_identities(dof_handler.get_fe_collection().size(),
                                   dof_handler.get_fe_collection().size());
 
-          // mark all vertices on ghost cells
+          // mark all vertices on ghost cells to identify those cells that we
+          // have already treated
           std::vector<bool> include_vertex(
             dof_handler.get_triangulation().n_vertices(), false);
           if (dynamic_cast<const dealii::parallel::
@@ -1290,19 +1267,13 @@ namespace internal
             dof_handler.hp_capability_enabled == true,
             (typename DoFHandler<dim, spacedim>::ExcOnlyAvailableWithHP()));
 
-          // we will mark lines that we have already treated, so first save and
-          // clear the user flags on lines and later restore them
-          std::vector<bool> user_flags;
-          dof_handler.get_triangulation().save_user_flags_line(user_flags);
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .clear_user_flags_line();
-
           // mark all lines on ghost cells
+          std::vector<bool> line_marked(
+            dof_handler.get_triangulation().n_raw_lines());
           for (const auto &cell : dof_handler.active_cell_iterators())
             if (cell->is_ghost())
               for (const auto l : cell->line_indices())
-                cell->line(l)->set_user_flag();
+                line_marked[cell->line(l)->index()] = true;
 
           // An implementation of the algorithm described in the hp-paper,
           // including the modification mentioned later in the "complications in
@@ -1329,10 +1300,10 @@ namespace internal
           for (const auto &cell : dof_handler.active_cell_iterators())
             for (const auto l : cell->line_indices())
               if ((cell->is_locally_owned()) &&
-                  (cell->line(l)->user_flag_set() == true))
+                  line_marked[cell->line(l)->index()])
                 {
-                  const auto line = cell->line(l);
-                  line->clear_user_flag();
+                  const auto line            = cell->line(l);
+                  line_marked[line->index()] = false;
 
                   unsigned int unique_sets_of_dofs =
                     line->n_active_fe_indices();
@@ -1531,11 +1502,6 @@ namespace internal
                         }
                     }
                 }
-
-          // finally restore the user flags
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .load_user_flags_line(user_flags);
         }
 
 
@@ -1571,20 +1537,13 @@ namespace internal
 
           const int dim = 3;
 
-          // we will mark quads that we have already treated, so first
-          // save and clear the user flags on quads and later restore
-          // them
-          std::vector<bool> user_flags;
-          dof_handler.get_triangulation().save_user_flags_quad(user_flags);
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .clear_user_flags_quad();
-
           // mark all quads on ghost cells
+          std::vector<bool> quad_marked(
+            dof_handler.get_triangulation().n_raw_quads());
           for (const auto &cell : dof_handler.active_cell_iterators())
             if (cell->is_ghost())
               for (const auto q : cell->face_indices())
-                cell->quad(q)->set_user_flag();
+                quad_marked[cell->quad(q)->index()] = true;
 
           // An implementation of the algorithm described in the hp-
           // paper, including the modification mentioned later in the
@@ -1605,11 +1564,11 @@ namespace internal
           for (const auto &cell : dof_handler.active_cell_iterators())
             for (const auto q : cell->face_indices())
               if ((cell->is_locally_owned()) &&
-                  (cell->quad(q)->user_flag_set() == true) &&
+                  quad_marked[cell->quad(q)->index()] &&
                   (cell->quad(q)->n_active_fe_indices() == 2))
                 {
-                  const auto quad = cell->quad(q);
-                  quad->clear_user_flag();
+                  const auto quad            = cell->quad(q);
+                  quad_marked[quad->index()] = false;
 
                   const std::set<unsigned int> fe_indices =
                     quad->get_active_fe_indices();
@@ -1685,11 +1644,6 @@ namespace internal
                           }
                     }
                 }
-
-          // finally restore the user flags
-          const_cast<dealii::Triangulation<dim, spacedim> &>(
-            dof_handler.get_triangulation())
-            .load_user_flags_quad(user_flags);
         }
 
 
@@ -2190,23 +2144,15 @@ namespace internal
 
           // deal with DoFs on lines
           {
-            // save user flags on lines so we can use them to mark lines
-            // we've already treated
-            std::vector<bool> saved_line_user_flags;
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .save_user_flags_line(saved_line_user_flags);
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .clear_user_flags_line();
-
+            std::vector<bool> line_touched(
+              dof_handler.get_triangulation().n_raw_lines());
             for (const auto &cell : dof_handler.active_cell_iterators())
               if (!cell->is_artificial())
                 for (const auto l : cell->line_indices())
-                  if (cell->line(l)->user_flag_set() == false)
+                  if (!line_touched[cell->line(l)->index()])
                     {
-                      const auto line = cell->line(l);
-                      line->set_user_flag();
+                      const auto line             = cell->line(l);
+                      line_touched[line->index()] = true;
 
                       const unsigned int n_active_fe_indices =
                         line->n_active_fe_indices();
@@ -2265,12 +2211,6 @@ namespace internal
                             }
                         }
                     }
-
-            // at the end, restore the user
-            // flags for the lines
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .load_user_flags_line(saved_line_user_flags);
           }
         }
 
@@ -2299,23 +2239,15 @@ namespace internal
 
           // deal with DoFs on lines
           {
-            // save user flags on lines so we can use them to mark lines
-            // we've already treated
-            std::vector<bool> saved_line_user_flags;
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .save_user_flags_line(saved_line_user_flags);
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .clear_user_flags_line();
-
+            std::vector<bool> line_touched(
+              dof_handler.get_triangulation().n_raw_lines());
             for (const auto &cell : dof_handler.active_cell_iterators())
               if (!cell->is_artificial())
                 for (const auto l : cell->line_indices())
-                  if (cell->line(l)->user_flag_set() == false)
+                  if (!line_touched[cell->line(l)->index()])
                     {
-                      const auto line = cell->line(l);
-                      line->set_user_flag();
+                      const auto line             = cell->line(l);
+                      line_touched[line->index()] = true;
 
                       const unsigned int n_active_fe_indices =
                         line->n_active_fe_indices();
@@ -2369,31 +2301,19 @@ namespace internal
                             }
                         }
                     }
-
-            // at the end, restore the user
-            // flags for the lines
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .load_user_flags_line(saved_line_user_flags);
           }
 
           // then deal with dofs on quads
           {
-            std::vector<bool> saved_quad_user_flags;
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .save_user_flags_quad(saved_quad_user_flags);
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .clear_user_flags_quad();
-
+            std::vector<bool> quad_touched(
+              dof_handler.get_triangulation().n_raw_quads());
             for (const auto &cell : dof_handler.active_cell_iterators())
               if (!cell->is_artificial())
                 for (const auto q : cell->face_indices())
-                  if (cell->quad(q)->user_flag_set() == false)
+                  if (!quad_touched[cell->quad(q)->index()])
                     {
-                      const auto quad = cell->quad(q);
-                      quad->set_user_flag();
+                      const auto quad             = cell->quad(q);
+                      quad_touched[quad->index()] = true;
 
                       const unsigned int n_active_fe_indices =
                         quad->n_active_fe_indices();
@@ -2452,11 +2372,6 @@ namespace internal
                             }
                         }
                     }
-
-            // at the end, restore the user flags for the quads
-            const_cast<dealii::Triangulation<dim, spacedim> &>(
-              dof_handler.get_triangulation())
-              .load_user_flags_quad(saved_quad_user_flags);
           }
         }
 
@@ -2637,28 +2552,23 @@ namespace internal
         {
           if (dof_handler.get_fe().n_dofs_per_line() > 0)
             {
-              // save user flags as they will be modified
-              std::vector<bool> user_flags;
-              dof_handler.get_triangulation().save_user_flags(user_flags);
-              const_cast<dealii::Triangulation<2, spacedim> &>(
-                dof_handler.get_triangulation())
-                .clear_user_flags();
-
-              // flag all lines adjacent to cells of the current
+              // mark all lines adjacent to cells of the current
               // level, as those lines logically belong to the same
               // level as the cell, at least for for isotropic
               // refinement
+              std::vector<bool> line_marked(
+                dof_handler.get_triangulation().n_raw_lines());
               for (const auto &cell :
                    dof_handler.cell_iterators_on_level(level))
                 if (cell->level_subdomain_id() !=
                     numbers::artificial_subdomain_id)
                   for (const unsigned int line : cell->face_indices())
-                    cell->face(line)->set_user_flag();
+                    line_marked[cell->face(line)->index()] = true;
 
               for (const auto &cell :
                    dof_handler.cell_iterators_on_level(level))
                 for (const auto l : cell->line_indices())
-                  if (cell->line(l)->user_flag_set())
+                  if (line_marked[cell->line(l)->index()])
                     {
                       for (unsigned int d = 0;
                            d < dof_handler.get_fe().n_dofs_per_line();
@@ -2679,12 +2589,8 @@ namespace internal
                                  new_numbers[indices_we_care_about
                                                .index_within_set(idx)]));
                         }
-                      cell->line(l)->clear_user_flag();
+                      line_marked[cell->line(l)->index()] = false;
                     }
-              // finally, restore user flags
-              const_cast<dealii::Triangulation<2, spacedim> &>(
-                dof_handler.get_triangulation())
-                .load_user_flags(user_flags);
             }
         }
 
@@ -2702,27 +2608,22 @@ namespace internal
           if (dof_handler.get_fe().n_dofs_per_line() > 0 ||
               dof_handler.get_fe().max_dofs_per_quad() > 0)
             {
-              // save user flags as they will be modified
-              std::vector<bool> user_flags;
-              dof_handler.get_triangulation().save_user_flags(user_flags);
-              const_cast<dealii::Triangulation<3, spacedim> &>(
-                dof_handler.get_triangulation())
-                .clear_user_flags();
-
               // flag all lines adjacent to cells of the current
               // level, as those lines logically belong to the same
               // level as the cell, at least for isotropic refinement
+              std::vector<bool> line_marked(
+                dof_handler.get_triangulation().n_raw_lines());
               for (const auto &cell :
                    dof_handler.cell_iterators_on_level(level))
                 if (cell->level_subdomain_id() !=
                     numbers::artificial_subdomain_id)
                   for (const auto line : cell->line_indices())
-                    cell->line(line)->set_user_flag();
+                    line_marked[cell->line(line)->index()] = true;
 
               for (const auto &cell :
                    dof_handler.cell_iterators_on_level(level))
                 for (const auto l : cell->line_indices())
-                  if (cell->line(l)->user_flag_set())
+                  if (line_marked[cell->line(l)->index()])
                     {
                       for (unsigned int d = 0;
                            d < dof_handler.get_fe().n_dofs_per_line();
@@ -2743,22 +2644,24 @@ namespace internal
                                  new_numbers[indices_we_care_about
                                                .index_within_set(idx)]));
                         }
-                      cell->line(l)->clear_user_flag();
+                      line_marked[cell->line(l)->index()] = false;
                     }
 
               // flag all quads adjacent to cells of the current level, as
               // those quads logically belong to the same level as the cell,
               // at least for isotropic refinement
+              std::vector<bool> quad_marked(
+                dof_handler.get_triangulation().n_raw_quads());
               for (const auto &cell :
                    dof_handler.cell_iterators_on_level(level))
                 if (cell->level_subdomain_id() !=
                     numbers::artificial_subdomain_id)
                   for (const auto quad : cell->face_indices())
-                    cell->quad(quad)->set_user_flag();
+                    quad_marked[cell->quad(quad)->index()] = true;
 
               for (const auto &cell : dof_handler.cell_iterators())
                 for (const auto l : cell->face_indices())
-                  if (cell->quad(l)->user_flag_set())
+                  if (quad_marked[cell->quad(l)->index()])
                     {
                       for (unsigned int d = 0;
                            d < dof_handler.get_fe().n_dofs_per_quad(l);
@@ -2779,13 +2682,8 @@ namespace internal
                                  new_numbers[indices_we_care_about
                                                .index_within_set(idx)]));
                         }
-                      cell->quad(l)->clear_user_flag();
+                      quad_marked[cell->quad(l)->index()] = false;
                     }
-
-              // finally, restore user flags
-              const_cast<dealii::Triangulation<3, spacedim> &>(
-                dof_handler.get_triangulation())
-                .load_user_flags(user_flags);
             }
         }
 
@@ -2873,13 +2771,6 @@ namespace internal
       std::vector<NumberCache>
       Sequential<dim, spacedim>::distribute_mg_dofs() const
       {
-        std::vector<bool> user_flags;
-        dof_handler->get_triangulation().save_user_flags(user_flags);
-
-        const_cast<dealii::Triangulation<dim, spacedim> &>(
-          dof_handler->get_triangulation())
-          .clear_user_flags();
-
         std::vector<NumberCache> number_caches;
         number_caches.reserve(dof_handler->get_triangulation().n_levels());
         for (unsigned int level = 0;
@@ -2895,10 +2786,6 @@ namespace internal
             number_caches.emplace_back(n_level_dofs);
           }
 
-        const_cast<dealii::Triangulation<dim, spacedim> &>(
-          dof_handler->get_triangulation())
-          .load_user_flags(user_flags);
-
         return number_caches;
       }
 
@@ -3608,13 +3495,12 @@ namespace internal
       {
         template <int dim, int spacedim>
         void
-        communicate_mg_ghost_cells(DoFHandler<dim, spacedim> &dof_handler)
+        communicate_mg_ghost_cells(DoFHandler<dim, spacedim> &     dof_handler,
+                                   std::vector<std::vector<bool>> &cell_marked)
         {
-          const auto pack = [&](const auto &cell) {
+          const auto pack = [](const auto &cell) {
             // why would somebody request a cell that is not ours?
-            Assert(cell->level_subdomain_id() ==
-                     dof_handler.get_triangulation().locally_owned_subdomain(),
-                   ExcInternalError());
+            Assert(cell->is_locally_owned_on_level(), ExcInternalError());
 
             std::vector<dealii::types::global_dof_index> data(
               cell->get_fe().n_dofs_per_cell());
@@ -3623,7 +3509,8 @@ namespace internal
             return data;
           };
 
-          const auto unpack = [](const auto &cell, const auto &dofs) {
+          const auto unpack = [&cell_marked](const auto &cell,
+                                             const auto &dofs) {
             Assert(cell->get_fe().n_dofs_per_cell() == dofs.size(),
                    ExcInternalError());
 
@@ -3654,17 +3541,17 @@ namespace internal
               true);
 
             if (!complete)
-              const_cast<
-                typename DoFHandler<dim, spacedim>::level_cell_iterator &>(cell)
-                ->set_user_flag();
+              {
+                // We should have the cell already marked
+                Assert(cell_marked[cell->level()][cell->index()],
+                       ExcInternalError());
+              }
             else
-              const_cast<
-                typename DoFHandler<dim, spacedim>::level_cell_iterator &>(cell)
-                ->clear_user_flag();
+              cell_marked[cell->level()][cell->index()] = false;
           };
 
-          const auto filter = [](const auto &cell) {
-            return cell->user_flag_set();
+          const auto filter = [&cell_marked](const auto &cell) {
+            return cell_marked[cell->level()][cell->index()];
           };
 
           GridTools::exchange_cell_data_to_level_ghosts<
@@ -3695,7 +3582,8 @@ namespace internal
         template <int dim, int spacedim>
         void
         communicate_dof_indices_on_marked_cells(
-          const DoFHandler<dim, spacedim> &dof_handler)
+          const DoFHandler<dim, spacedim> &dof_handler,
+          std::vector<bool> &              cell_marked)
         {
 #  ifndef DEAL_II_WITH_MPI
           (void)dof_handler;
@@ -3718,7 +3606,8 @@ namespace internal
             return data;
           };
 
-          const auto unpack = [](const auto &cell, const auto &dofs) {
+          const auto unpack = [&cell_marked](const auto &cell,
+                                             const auto &dofs) {
             Assert(cell->get_fe().n_dofs_per_cell() == dofs.size(),
                    ExcInternalError());
 
@@ -3750,19 +3639,17 @@ namespace internal
               false);
 
             if (!complete)
-              const_cast<
-                typename DoFHandler<dim, spacedim>::active_cell_iterator &>(
-                cell)
-                ->set_user_flag();
+              {
+                // We should have the cell already marked
+                Assert(cell_marked[cell->active_cell_index()],
+                       ExcInternalError());
+              }
             else
-              const_cast<
-                typename DoFHandler<dim, spacedim>::active_cell_iterator &>(
-                cell)
-                ->clear_user_flag();
+              cell_marked[cell->active_cell_index()] = false;
           };
 
-          const auto filter = [](const auto &cell) {
-            return cell->user_flag_set();
+          const auto filter = [&cell_marked](const auto &cell) {
+            return cell_marked[cell->active_cell_index()];
           };
 
           GridTools::exchange_cell_data_to_ghosts<
@@ -3918,16 +3805,13 @@ namespace internal
         //                    cells somewhere else, send our own DoF indices
         //                    to the appropriate set of other processors
         {
-          std::vector<bool> user_flags;
-          triangulation->save_user_flags(user_flags);
-          triangulation->clear_user_flags();
-
           // mark all cells that either have to send data (locally
           // owned cells that are adjacent to ghost neighbors in some
           // way) or receive data (all ghost cells) via the user flags
+          std::vector<bool> cell_marked(triangulation->n_active_cells());
           for (const auto &cell : dof_handler->active_cell_iterators())
             if (cell->is_ghost())
-              cell->set_user_flag();
+              cell_marked[cell->active_cell_index()] = true;
 
           // Send and receive cells. After this, only the local cells
           // are marked, that received new data. This has to be
@@ -3935,7 +3819,7 @@ namespace internal
           //
           // as explained in the 'distributed' paper, this has to be
           // done twice
-          communicate_dof_indices_on_marked_cells(*dof_handler);
+          communicate_dof_indices_on_marked_cells(*dof_handler, cell_marked);
 
           // If the DoFHandler has hp-capabilities enabled, then we may have
           // received valid indices of degrees of freedom that are dominated
@@ -3948,16 +3832,15 @@ namespace internal
           //                    DoF indices set. however, some ghost cells
           //                    may still have invalid ones. thus, exchange
           //                    one more time.
-          communicate_dof_indices_on_marked_cells(*dof_handler);
+          communicate_dof_indices_on_marked_cells(*dof_handler, cell_marked);
 
           // at this point, we must have taken care of the data transfer
           // on all cells we had previously marked. verify this
 #  ifdef DEBUG
           for (const auto &cell : dof_handler->active_cell_iterators())
-            Assert(cell->user_flag_set() == false, ExcInternalError());
+            Assert(cell_marked[cell->active_cell_index()] == false,
+                   ExcInternalError());
 #  endif
-
-          triangulation->load_user_flags(user_flags);
         }
 
         update_all_active_cell_dof_indices_caches(*this->dof_handler);
@@ -4145,38 +4028,28 @@ namespace internal
         // these cells from the corresponding owners. As this information
         // can be incomplete,
         {
-          std::vector<bool> user_flags;
-          triangulation->save_user_flags(user_flags);
-          triangulation->clear_user_flags();
-
-          // mark all ghost cells for transfer
-          {
-            for (const auto &cell : dof_handler->cell_iterators())
-              if (cell->is_ghost_on_level())
-                cell->set_user_flag();
-          }
+          std::vector<std::vector<bool>> cell_marked(triangulation->n_levels());
+          for (unsigned int l = 0; l < triangulation->n_levels(); ++l)
+            cell_marked[l].resize(triangulation->n_raw_cells(l));
+          for (const auto &cell : dof_handler->cell_iterators())
+            if (cell->is_ghost_on_level())
+              cell_marked[cell->level()][cell->index()] = true;
 
           // Phase 1. Request all marked cells from corresponding owners. If we
           // managed to get every DoF, remove the user_flag, otherwise we
           // will request them again in the step below.
-          communicate_mg_ghost_cells(*dof_handler);
+          communicate_mg_ghost_cells(*dof_handler, cell_marked);
 
           // Phase 2, only request the cells that were not completed
           // in Phase 1.
-          communicate_mg_ghost_cells(*dof_handler);
+          communicate_mg_ghost_cells(*dof_handler, cell_marked);
 
 #  ifdef DEBUG
-          // make sure we have removed all flags:
-          {
-            for (const auto &cell : dof_handler->cell_iterators())
-              if (cell->level_subdomain_id() !=
-                    dealii::numbers::artificial_subdomain_id &&
-                  !cell->is_locally_owned_on_level())
-                Assert(cell->user_flag_set() == false, ExcInternalError());
-          }
+          // make sure we have finished all cells:
+          for (const auto &cell : dof_handler->cell_iterators())
+            Assert(cell_marked[cell->level()][cell->index()] == false,
+                   ExcInternalError());
 #  endif
-
-          triangulation->load_user_flags(user_flags);
         }
 
 
@@ -4307,15 +4180,11 @@ namespace internal
         // taking into account that we have to unify a few DoFs in between
         // then communication phases if we do hp-numbering
         {
-          std::vector<bool> user_flags;
-          triangulation->save_user_flags(user_flags);
-          triangulation->clear_user_flags();
-
-          // mark all own cells for transfer
+          // mark all ghost cells for transfer
+          std::vector<bool> cell_marked(triangulation->n_active_cells());
           for (const auto &cell : dof_handler->active_cell_iterators())
             if (cell->is_ghost())
-              cell->set_user_flag();
-
+              cell_marked[cell->active_cell_index()] = true;
 
           // Send and receive cells. After this, only the local cells
           // are marked, that received new data. This has to be
@@ -4323,7 +4192,7 @@ namespace internal
           //
           // as explained in the 'distributed' paper, this has to be
           // done twice
-          communicate_dof_indices_on_marked_cells(*dof_handler);
+          communicate_dof_indices_on_marked_cells(*dof_handler, cell_marked);
 
           // if the DoFHandler has hp-capabilities then we may have
           // received valid indices of degrees of freedom that are
@@ -4333,9 +4202,8 @@ namespace internal
           Implementation::merge_invalid_dof_indices_on_ghost_interfaces(
             *dof_handler);
 
-          communicate_dof_indices_on_marked_cells(*dof_handler);
+          communicate_dof_indices_on_marked_cells(*dof_handler, cell_marked);
 
-          triangulation->load_user_flags(user_flags);
           update_all_active_cell_dof_indices_caches(*this->dof_handler);
         }
 
@@ -4430,20 +4298,16 @@ namespace internal
 
         // communicate newly assigned DoF indices with other processors
         {
-          std::vector<bool> user_flags;
-          triangulation->save_user_flags(user_flags);
-          triangulation->clear_user_flags();
-
-          // mark only cells on the chosen level
+          std::vector<std::vector<bool>> cell_marked(triangulation->n_levels());
+          for (unsigned int l = 0; l < triangulation->n_levels(); ++l)
+            cell_marked[l].resize(triangulation->n_raw_cells(l));
           for (const auto &cell : dof_handler->cell_iterators_on_level(level))
             if (cell->is_ghost_on_level())
-              cell->set_user_flag();
-
-          communicate_mg_ghost_cells(*dof_handler);
+              cell_marked[cell->level()][cell->index()] = true;
 
-          communicate_mg_ghost_cells(*dof_handler);
+          communicate_mg_ghost_cells(*dof_handler, cell_marked);
 
-          triangulation->load_user_flags(user_flags);
+          communicate_mg_ghost_cells(*dof_handler, cell_marked);
         }
 
         NumberCache number_cache;
diff --git a/source/grid/tria_description.cc b/source/grid/tria_description.cc
index 2d3d6334f6..4de96bb503 100644
--- a/source/grid/tria_description.cc
+++ b/source/grid/tria_description.cc
@@ -385,12 +385,13 @@ namespace TriangulationDescription
        */
       template <int dim, int spacedim>
       void
-      set_user_flag_and_of_its_parents(
-        const TriaIterator<CellAccessor<dim, spacedim>> &cell)
+      mark_cell_and_its_parents(
+        const TriaIterator<CellAccessor<dim, spacedim>> &cell,
+        std::vector<std::vector<bool>> &                 cell_marked)
       {
-        cell->set_user_flag();
+        cell_marked[cell->level()][cell->index()] = true;
         if (cell->level() != 0)
-          set_user_flag_and_of_its_parents(cell->parent());
+          mark_cell_and_its_parents(cell->parent(), cell_marked);
       }
 
       /**
@@ -466,16 +467,12 @@ namespace TriangulationDescription
                   }
               };
 
-          // 1) collect locally relevant cells (set user_flag)
-          std::vector<bool> old_user_flags;
-          tria.save_user_flags(old_user_flags);
+          // 1) loop over levels (from fine to coarse) and mark on each level
+          //    the locally relevant cells
+          std::vector<std::vector<bool>> cell_marked(tria.n_levels());
+          for (unsigned int l = 0; l < tria.n_levels(); ++l)
+            cell_marked[l].resize(tria.n_raw_cells(l));
 
-          // 1a) clear user_flags
-          const_cast<dealii::Triangulation<dim, spacedim> &>(tria)
-            .clear_user_flags();
-
-          // 1b) loop over levels (from fine to coarse) and mark on each level
-          //     the locally relevant cells
           for (int level = tria.get_triangulation().n_global_levels() - 1;
                level >= 0;
                --level)
@@ -509,7 +506,7 @@ namespace TriangulationDescription
               // mark all locally relevant cells
               for (const auto &cell : tria.cell_iterators_on_level(level))
                 if (is_locally_relevant_on_level(cell))
-                  set_user_flag_and_of_its_parents(cell);
+                  mark_cell_and_its_parents(cell, cell_marked);
             }
 
           // 2) set_up coarse-grid triangulation
@@ -520,7 +517,7 @@ namespace TriangulationDescription
             // a) loop over all cells
             for (const auto &cell : tria.cell_iterators_on_level(0))
               {
-                if (!cell->user_flag_set())
+                if (!cell_marked[cell->level()][cell->index()])
                   continue;
 
                 // extract cell definition (with old numbering of vertices)
@@ -597,7 +594,7 @@ namespace TriangulationDescription
               for (const auto &cell : tria.cell_iterators_on_level(level))
                 {
                   // check if cell is locally relevant
-                  if (!(cell->user_flag_set()))
+                  if (!cell_marked[cell->level()][cell->index()])
                     continue;
 
                   CellData<dim> cell_info;
@@ -658,9 +655,6 @@ namespace TriangulationDescription
                 }
             }
 
-          const_cast<dealii::Triangulation<dim, spacedim> &>(tria)
-            .load_user_flags(old_user_flags);
-
           return construction_data;
         }
 
diff --git a/source/multigrid/mg_tools.cc b/source/multigrid/mg_tools.cc
index 9737241408..34cfe0ff07 100644
--- a/source/multigrid/mg_tools.cc
+++ b/source/multigrid/mg_tools.cc
@@ -115,17 +115,10 @@ namespace MGTools
     // Function starts here by
     // resetting the counters.
     std::fill(row_lengths.begin(), row_lengths.end(), 0);
-    // We need the user flags, so we
-    // save them for later restoration
-    std::vector<bool> old_flags;
-    // We need a non-constant
-    // triangulation for the user
-    // flags. Since we restore them in
-    // the end, this cast is safe.
-    Triangulation<dim, spacedim> &user_flags_triangulation =
-      const_cast<Triangulation<dim, spacedim> &>(dofs.get_triangulation());
-    user_flags_triangulation.save_user_flags(old_flags);
-    user_flags_triangulation.clear_user_flags();
+
+    std::vector<bool> face_touched(dim == 2 ?
+                                     dofs.get_triangulation().n_raw_lines() :
+                                     dofs.get_triangulation().n_raw_quads());
 
     std::vector<types::global_dof_index> cell_indices;
     std::vector<types::global_dof_index> neighbor_indices;
@@ -262,9 +255,10 @@ namespace MGTools
 
             // Do this only once per
             // face.
-            if (face->user_flag_set())
+            if (face_touched[face->index()])
               continue;
-            face->set_user_flag();
+            face_touched[face->index()] = true;
+
             // At this point, we assume
             // that each cell added its
             // dofs minus the face to
@@ -289,7 +283,6 @@ namespace MGTools
                 fe.n_dofs_per_face(face_no);
           }
       }
-    user_flags_triangulation.load_user_flags(old_flags);
   }
 
 
@@ -308,17 +301,10 @@ namespace MGTools
     // Function starts here by
     // resetting the counters.
     std::fill(row_lengths.begin(), row_lengths.end(), 0);
-    // We need the user flags, so we
-    // save them for later restoration
-    std::vector<bool> old_flags;
-    // We need a non-constant
-    // triangulation for the user
-    // flags. Since we restore them in
-    // the end, this cast is safe.
-    Triangulation<dim, spacedim> &user_flags_triangulation =
-      const_cast<Triangulation<dim, spacedim> &>(dofs.get_triangulation());
-    user_flags_triangulation.save_user_flags(old_flags);
-    user_flags_triangulation.clear_user_flags();
+
+    std::vector<bool> face_touched(dim == 2 ?
+                                     dofs.get_triangulation().n_raw_lines() :
+                                     dofs.get_triangulation().n_raw_quads());
 
     std::vector<types::global_dof_index> cell_indices;
     std::vector<types::global_dof_index> neighbor_indices;
@@ -530,12 +516,11 @@ namespace MGTools
                       row_lengths[cell_indices[local_dof]] += dof_increment;
                     }
 
-            // Do this only once per
-            // face and not on the
-            // hanging faces.
-            if (face->user_flag_set())
+            // Do this only once per face and not on the hanging faces.
+            if (face_touched[face->index()])
               continue;
-            face->set_user_flag();
+            face_touched[face->index()] = true;
+
             // At this point, we assume
             // that each cell added its
             // dofs minus the face to
@@ -582,7 +567,6 @@ namespace MGTools
                       fe.base_element(base).n_dofs_per_face(face_no);
           }
       }
-    user_flags_triangulation.load_user_flags(old_flags);
   }
 
 
@@ -642,9 +626,7 @@ namespace MGTools
     const unsigned int dofs_per_cell = dof.get_fe().n_dofs_per_cell();
     std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
     std::vector<types::global_dof_index> dofs_on_other_cell(dofs_per_cell);
-    typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(level),
-                                                      endc = dof.end(level);
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof.cell_iterators_on_level(level))
       {
         if (!cell->is_locally_owned_on_level())
           continue;
@@ -723,9 +705,7 @@ namespace MGTools
     const unsigned int dofs_per_cell = dof.get_fe().n_dofs_per_cell();
     std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
     std::vector<types::global_dof_index> dofs_on_other_cell(dofs_per_cell);
-    typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(level),
-                                                      endc = dof.end(level);
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof.cell_iterators_on_level(level))
       {
         if (!cell->is_locally_owned_on_level())
           continue;
@@ -801,9 +781,6 @@ namespace MGTools
     Table<2, bool>                       support_on_face(total_dofs,
                                    GeometryInfo<dim>::faces_per_cell);
 
-    typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(level),
-                                                      endc = dof.end(level);
-
     const Table<2, DoFTools::Coupling>
       int_dof_mask =
         DoFTools::dof_couplings_from_component_couplings(fe, int_mask),
@@ -814,20 +791,11 @@ namespace MGTools
       for (auto f : GeometryInfo<dim>::face_indices())
         support_on_face(i, f) = fe.has_support_on_face(i, f);
 
-    // 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 state as it was at the
-    // beginning of this function.
-    std::vector<bool> user_flags;
-    dof.get_triangulation().save_user_flags(user_flags);
-    const_cast<Triangulation<dim, spacedim> &>(dof.get_triangulation())
-      .clear_user_flags();
-
-    for (; cell != endc; ++cell)
+    std::vector<bool> face_touched(dim == 2 ?
+                                     dof.get_triangulation().n_raw_lines() :
+                                     dof.get_triangulation().n_raw_quads());
+
+    for (const auto &cell : dof.cell_iterators_on_level(level))
       {
         if (!cell->is_locally_owned_on_level())
           continue;
@@ -844,7 +812,7 @@ namespace MGTools
           {
             typename DoFHandler<dim, spacedim>::face_iterator cell_face =
               cell->face(face);
-            if (cell_face->user_flag_set())
+            if (face_touched[cell_face->index()])
               continue;
 
             if (cell->at_boundary(face) && !cell->has_periodic_neighbor(face))
@@ -944,14 +912,10 @@ namespace MGTools
                           }
                       }
                   }
-                neighbor->face(neighbor_face)->set_user_flag();
+                face_touched[neighbor->face(neighbor_face)->index()] = true;
               }
           }
       }
-
-    // finally restore the user flags
-    const_cast<Triangulation<dim, spacedim> &>(dof.get_triangulation())
-      .load_user_flags(user_flags);
   }
 
 
@@ -992,9 +956,6 @@ namespace MGTools
     Table<2, bool>                       support_on_face(dofs_per_cell,
                                    GeometryInfo<dim>::faces_per_cell);
 
-    typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(level),
-                                                      endc = dof.end(level);
-
     const Table<2, DoFTools::Coupling> flux_dof_mask =
       DoFTools::dof_couplings_from_component_couplings(fe, flux_mask);
 
@@ -1002,7 +963,7 @@ namespace MGTools
       for (auto f : GeometryInfo<dim>::face_indices())
         support_on_face(i, f) = fe.has_support_on_face(i, f);
 
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof.cell_iterators_on_level(level))
       {
         if (!cell->is_locally_owned_on_level())
           continue;
@@ -1065,9 +1026,7 @@ namespace MGTools
 
     const unsigned int dofs_per_cell = dof.get_fe().n_dofs_per_cell();
     std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
-    typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(level),
-                                                      endc = dof.end(level);
-    for (; cell != endc; ++cell)
+    for (const auto &cell : dof.cell_iterators_on_level(level))
       if (cell->is_locally_owned_on_level())
         {
           cell->get_mg_dof_indices(dofs_on_this_cell);
@@ -1494,10 +1453,7 @@ namespace MGTools
 
     std::vector<bool> cell_dofs(dofs_per_cell, false);
 
-    typename DoFHandler<dim>::cell_iterator cell = mg_dof_handler.begin(),
-                                            endc = mg_dof_handler.end();
-
-    for (; cell != endc; ++cell)
+    for (const auto &cell : mg_dof_handler.cell_iterators())
       {
         // Do not look at artificial level cells (in a serial computation we
         // need to ignore the level_subdomain_id() because it is never set).
@@ -1573,10 +1529,7 @@ namespace MGTools
     // the first locally owned cell we find will have
     // the lowest level in the particular subdomain.
     unsigned int min_level = tria.n_global_levels();
-    typename Triangulation<dim, spacedim>::active_cell_iterator
-      cell = tria.begin_active(),
-      endc = tria.end();
-    for (; cell != endc; ++cell)
+    for (const auto &cell : tria.active_cell_iterators())
       if (cell->is_locally_owned())
         {
           min_level = cell->level();