From: Wolfgang Bangerth <bangerth@colostate.edu>
Date: Mon, 16 Oct 2023 22:31:26 +0000 (-0600)
Subject: Avoid constraining DoFs more than once.
X-Git-Tag: relicensing~382^2~1
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d6603822b81a3772b8c1e147a496c11498b406be;p=dealii.git

Avoid constraining DoFs more than once.
---

diff --git a/source/dofs/dof_tools_constraints.cc b/source/dofs/dof_tools_constraints.cc
index 6ec936f7b2..ca27c13dac 100644
--- a/source/dofs/dof_tools_constraints.cc
+++ b/source/dofs/dof_tools_constraints.cc
@@ -1022,20 +1022,27 @@ namespace DoFTools
                 Assert(dofs_on_children.size() <= n_dofs_on_children,
                        ExcInternalError());
 
-                // for each row in the AffineConstraints object for this line:
+                // For each row in the AffineConstraints object for
+                // this line, add the constraint. Ignore rows that
+                // have already been added (e.g., in 3d degrees of
+                // freedom on edges with hanging nodes will be visited
+                // more than once).
                 for (unsigned int row = 0; row != dofs_on_children.size();
                      ++row)
-                  {
-                    constraint_entries.clear();
-                    constraint_entries.reserve(dofs_on_mother.size());
-                    for (unsigned int i = 0; i != dofs_on_mother.size(); ++i)
-                      constraint_entries.emplace_back(dofs_on_mother[i],
-                                                      fe.constraints()(row, i));
-
-                    constraints.add_constraint(dofs_on_children[row],
-                                               constraint_entries,
-                                               0.);
-                  }
+                  if (constraints.is_constrained(dofs_on_children[row]) ==
+                      false)
+                    {
+                      constraint_entries.clear();
+                      constraint_entries.reserve(dofs_on_mother.size());
+                      for (unsigned int i = 0; i != dofs_on_mother.size(); ++i)
+                        constraint_entries.emplace_back(dofs_on_mother[i],
+                                                        fe.constraints()(row,
+                                                                         i));
+
+                      constraints.add_constraint(dofs_on_children[row],
+                                                 constraint_entries,
+                                                 0.);
+                    }
               }
             else
               {
@@ -3576,9 +3583,22 @@ namespace DoFTools
                             }
 
                         if (nonzero)
-                          zero_boundary_constraints.add_constraint(face_dof,
-                                                                   {},
-                                                                   0.);
+                          {
+                            // Check that either (i) the DoF is not
+                            // yet constrained, or (ii) if it is, its
+                            // inhomogeneity is zero:
+                            if (zero_boundary_constraints.is_constrained(
+                                  face_dof) == false)
+                              zero_boundary_constraints.add_constraint(face_dof,
+                                                                       {},
+                                                                       0.);
+                            else
+                              Assert(zero_boundary_constraints
+                                         .is_inhomogeneously_constrained(
+                                           face_dof) == false,
+                                     ExcInternalError());
+                          }
+
                         // We already dealt with this DoF. Make sure we
                         // don't touch it again.
                         dofs_already_treated.insert(face_dof);
diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc
index b62a661b90..8f0a628d0e 100644
--- a/source/grid/grid_tools.cc
+++ b/source/grid/grid_tools.cc
@@ -2281,11 +2281,13 @@ namespace GridTools
     typename std::map<unsigned int, Point<dim>>::const_iterator map_end =
       new_points.end();
 
-    // fill these maps using the data given by new_points
+    // Fill these maps using the data given by new_points
     for (const auto &cell : dof_handler.active_cell_iterators())
       {
-        // loop over all vertices of the cell and see if it is listed in the map
-        // given as first argument of the function
+        // Loop over all vertices of the cell and see if it is listed in the map
+        // given as first argument of the function. We visit vertices multiple
+        // times, so also check that if we have already added a constraint, we
+        // don't do it a second time again.
         for (const unsigned int vertex_no : cell->vertex_indices())
           {
             const unsigned int vertex_index = cell->vertex_index(vertex_no);
@@ -2296,14 +2298,16 @@ namespace GridTools
 
             if (map_iter != map_end)
               for (unsigned int i = 0; i < dim; ++i)
-                {
-                  constraints[i].add_constraint(
-                    cell->vertex_dof_index(vertex_no, 0),
-                    {},
-                    (solve_for_absolute_positions ?
-                       map_iter->second(i) :
-                       map_iter->second(i) - vertex_point[i]));
-                }
+                if (constraints[i].is_constrained(
+                      cell->vertex_dof_index(vertex_no, 0)) == false)
+                  {
+                    constraints[i].add_constraint(
+                      cell->vertex_dof_index(vertex_no, 0),
+                      {},
+                      (solve_for_absolute_positions ?
+                         map_iter->second(i) :
+                         map_iter->second(i) - vertex_point[i]));
+                  }
           }
       }