]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Minor cleanups. 16064/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Thu, 28 Sep 2023 23:02:33 +0000 (17:02 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Fri, 29 Sep 2023 07:16:19 +0000 (01:16 -0600)
include/deal.II/multigrid/mg_constrained_dofs.h

index ee22cdb2aedb6ddbd6a318eecf3d79ca18e98bc5..dd3e9a2b394903a392f4a093c2319958450e5591 100644 (file)
@@ -271,6 +271,7 @@ private:
 };
 
 
+
 template <int dim, int spacedim>
 inline void
 MGConstrainedDoFs::initialize(
@@ -287,8 +288,8 @@ MGConstrainedDoFs::initialize(
   const unsigned int max_level = (level_relevant_dofs.max_level() == 0) ?
                                    nlevels - 1 :
                                    level_relevant_dofs.max_level();
-  const bool         user_level_dofs =
-    (level_relevant_dofs.max_level() == 0) ? false : true;
+  const bool         use_provided_level_relevant_dofs =
+    (level_relevant_dofs.max_level() > 0);
 
   // At this point level_constraint and refinement_edge_indices are empty.
   refinement_edge_indices.resize(nlevels);
@@ -296,7 +297,7 @@ MGConstrainedDoFs::initialize(
   user_constraints.resize(nlevels);
   for (unsigned int l = min_level; l <= max_level; ++l)
     {
-      if (user_level_dofs)
+      if (use_provided_level_relevant_dofs)
         {
           level_constraints[l].reinit(dof.locally_owned_mg_dofs(l),
                                       level_relevant_dofs[l]);
@@ -311,60 +312,54 @@ MGConstrainedDoFs::initialize(
 
       // Loop through relevant cells and faces finding those which are periodic
       // neighbors.
-      typename DoFHandler<dim, spacedim>::cell_iterator cell = dof.begin(l),
-                                                        endc = dof.end(l);
-      for (; cell != endc; ++cell)
-        if (cell->level_subdomain_id() != numbers::artificial_subdomain_id)
-          {
-            for (auto f : cell->face_indices())
-              if (cell->has_periodic_neighbor(f) &&
-                  cell->periodic_neighbor(f)->level() == cell->level())
-                {
-                  if (cell->is_locally_owned_on_level())
-                    {
-                      Assert(
-                        cell->periodic_neighbor(f)->level_subdomain_id() !=
-                          numbers::artificial_subdomain_id,
-                        ExcMessage(
-                          "Periodic neighbor of a locally owned cell must either be owned or ghost."));
-                    }
-                  // Cell is a level-ghost and its neighbor is a
-                  // level-artificial cell nothing to do here
-                  else if (cell->periodic_neighbor(f)->level_subdomain_id() ==
-                           numbers::artificial_subdomain_id)
+      for (const auto &cell : dof.cell_iterators_on_level(l))
+        if (cell->is_artificial_on_level() == false)
+          for (auto f : cell->face_indices())
+            if (cell->has_periodic_neighbor(f) &&
+                cell->periodic_neighbor(f)->level() == cell->level())
+              {
+                if (cell->is_locally_owned_on_level())
+                  {
+                    Assert(
+                      cell->periodic_neighbor(f)->level_subdomain_id() !=
+                        numbers::artificial_subdomain_id,
+                      ExcMessage(
+                        "Periodic neighbor of a locally owned cell must either be owned or ghost."));
+                  }
+                // Cell is a level-ghost and its neighbor is a
+                // level-artificial cell nothing to do here
+                else if (cell->periodic_neighbor(f)->level_subdomain_id() ==
+                         numbers::artificial_subdomain_id)
+                  {
+                    Assert(cell->is_locally_owned_on_level() == false,
+                           ExcInternalError());
+                    continue;
+                  }
+
+                const unsigned int dofs_per_face =
+                  dof.get_fe(0).n_dofs_per_face(f);
+                std::vector<types::global_dof_index> dofs_1(dofs_per_face);
+                std::vector<types::global_dof_index> dofs_2(dofs_per_face);
+
+                cell->periodic_neighbor(f)
+                  ->face(cell->periodic_neighbor_face_no(f))
+                  ->get_mg_dof_indices(l, dofs_1, 0);
+                cell->face(f)->get_mg_dof_indices(l, dofs_2, 0);
+                // Store periodicity information in the level
+                // AffineConstraints object. Skip DoFs for which we've
+                // previously entered periodicity constraints already; this
+                // can happen, for example, for a vertex dof at a periodic
+                // boundary that we visit from more than one cell
+                for (unsigned int i = 0; i < dofs_per_face; ++i)
+                  if (level_constraints[l].can_store_line(dofs_2[i]) &&
+                      level_constraints[l].can_store_line(dofs_1[i]) &&
+                      !level_constraints[l].is_constrained(dofs_2[i]) &&
+                      !level_constraints[l].is_constrained(dofs_1[i]))
                     {
-                      Assert(cell->is_locally_owned_on_level() == false,
-                             ExcInternalError());
-                      continue;
+                      level_constraints[l].add_line(dofs_2[i]);
+                      level_constraints[l].add_entry(dofs_2[i], dofs_1[i], 1.);
                     }
-
-                  const unsigned int dofs_per_face =
-                    dof.get_fe(0).n_dofs_per_face(f);
-                  std::vector<types::global_dof_index> dofs_1(dofs_per_face);
-                  std::vector<types::global_dof_index> dofs_2(dofs_per_face);
-
-                  cell->periodic_neighbor(f)
-                    ->face(cell->periodic_neighbor_face_no(f))
-                    ->get_mg_dof_indices(l, dofs_1, 0);
-                  cell->face(f)->get_mg_dof_indices(l, dofs_2, 0);
-                  // Store periodicity information in the level
-                  // AffineConstraints object. Skip DoFs for which we've
-                  // previously entered periodicity constraints already; this
-                  // can happen, for example, for a vertex dof at a periodic
-                  // boundary that we visit from more than one cell
-                  for (unsigned int i = 0; i < dofs_per_face; ++i)
-                    if (level_constraints[l].can_store_line(dofs_2[i]) &&
-                        level_constraints[l].can_store_line(dofs_1[i]) &&
-                        !level_constraints[l].is_constrained(dofs_2[i]) &&
-                        !level_constraints[l].is_constrained(dofs_1[i]))
-                      {
-                        level_constraints[l].add_line(dofs_2[i]);
-                        level_constraints[l].add_entry(dofs_2[i],
-                                                       dofs_1[i],
-                                                       1.);
-                      }
-                }
-          }
+              }
       level_constraints[l].close();
 
       // Initialize with empty IndexSet of correct size
@@ -375,6 +370,7 @@ MGConstrainedDoFs::initialize(
 }
 
 
+
 template <int dim, int spacedim>
 inline void
 MGConstrainedDoFs::make_zero_boundary_constraints(
@@ -470,6 +466,7 @@ MGConstrainedDoFs::make_no_normal_flux_constraints(
 }
 
 
+
 inline void
 MGConstrainedDoFs::add_user_constraints(
   const unsigned int               level,
@@ -521,6 +518,8 @@ MGConstrainedDoFs::is_boundary_index(const unsigned int            level,
   return boundary_indices[level].is_element(index);
 }
 
+
+
 inline bool
 MGConstrainedDoFs::at_refinement_edge(const unsigned int            level,
                                       const types::global_dof_index index) const
@@ -530,6 +529,8 @@ MGConstrainedDoFs::at_refinement_edge(const unsigned int            level,
   return refinement_edge_indices[level].is_element(index);
 }
 
+
+
 inline bool
 MGConstrainedDoFs::is_interface_matrix_entry(
   const unsigned int            level,
@@ -590,6 +591,7 @@ MGConstrainedDoFs::get_user_constraint_matrix(const unsigned int level) const
 }
 
 
+
 template <typename Number>
 inline void
 MGConstrainedDoFs::merge_constraints(AffineConstraints<Number> &constraints,

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.