]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Keep track of constraints we have already determined are finalized.
authorWolfgang Bangerth <bangerth@colostate.edu>
Tue, 11 Jul 2023 05:08:16 +0000 (23:08 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Wed, 12 Jul 2023 15:45:31 +0000 (09:45 -0600)
include/deal.II/lac/affine_constraints.templates.h

index 834384fa668a7d23b04bcac71ac93db08df4d4c5..f9a6e01eb65d64f1ba5782562668d3bf76cf3413 100644 (file)
@@ -712,139 +712,147 @@ AffineConstraints<number>::close()
 #ifdef DEBUG
   size_type iteration = 0;
 #endif
-  bool chained_constraint_replaced = false;
+  bool              chained_constraint_replaced = false;
+  std::vector<bool> line_finalized(lines.size(), false);
   do
     {
       chained_constraint_replaced = false;
-      for (ConstraintLine &line : lines)
-        {
+      for (unsigned int line_index = 0; line_index < lines.size(); ++line_index)
+        if (line_finalized[line_index] == false)
+          {
+            ConstraintLine &line = lines[line_index];
+
 #ifdef DEBUG
-          // we need to keep track of how many replacements we do in this line,
-          // because we can end up in a cycle A->B->C->A without the number of
-          // entries growing.
-          size_type n_replacements = 0;
+            // we need to keep track of how many replacements we do in this
+            // line, because we can end up in a cycle A->B->C->A without the
+            // number of entries growing.
+            size_type n_replacements = 0;
 #endif
 
-          // loop over all entries of this line (including ones that we
-          // have appended in this go around) and see whether they are
-          // further constrained. ignore elements that we don't store on
-          // the current processor.
-          //
-          // If we find that one of the entries in the current line is indeed
-          // constrained, then store the address of the line it corresponds to
-          // for further processing; this avoids having to look it up again,
-          // which requires an IndexSet lookup. If not, store a nullptr.
-          std::vector<const ConstraintLine *> sub_constraints(
-            line.entries.size(), nullptr);
-          for (unsigned int entry = 0; entry < line.entries.size(); ++entry)
-            if (((local_lines.size() == 0) ||
-                 (local_lines.is_element(line.entries[entry].first))) &&
-                is_constrained(line.entries[entry].first))
+            // loop over all entries of this line (including ones that we
+            // have appended in this go around) and see whether they are
+            // further constrained. ignore elements that we don't store on
+            // the current processor.
+            //
+            // If we find that one of the entries in the current line is indeed
+            // constrained, then store the address of the line it corresponds to
+            // for further processing; this avoids having to look it up again,
+            // which requires an IndexSet lookup. If not, store a nullptr.
+            std::vector<const ConstraintLine *> sub_constraints(
+              line.entries.size(), nullptr);
+            for (unsigned int entry = 0; entry < line.entries.size(); ++entry)
+              if (((local_lines.size() == 0) ||
+                   (local_lines.is_element(line.entries[entry].first))) &&
+                  is_constrained(line.entries[entry].first))
+                {
+                  const size_type dof_index = line.entries[entry].first;
+                  sub_constraints[entry] =
+                    &lines[lines_cache[calculate_line_index(dof_index)]];
+                }
+
+            // If none of the entries in the current line refer to DoFs that are
+            // themselves constrained, then we can move on:
+            if (std::none_of(sub_constraints.begin(),
+                             sub_constraints.end(),
+                             [](const auto p) { return (p != nullptr); }))
               {
-                const size_type dof_index = line.entries[entry].first;
-                sub_constraints[entry] =
-                  &lines[lines_cache[calculate_line_index(dof_index)]];
+                line_finalized[line_index] = true;
+                continue;
               }
 
-          // If none of the entries in the current line refer to DoFs that are
-          // themselves constrained, then we can move on:
-          if (std::none_of(sub_constraints.begin(),
-                           sub_constraints.end(),
-                           [](const auto p) { return (p != nullptr); }))
-            continue;
+            // Now walk through the original entries. We may replace some of
+            // them, and we may add some, but we only walk through the original
+            // entries (and also don't touch the replacements) even though in
+            // principle we could also do the ones we have added or replaced. We
+            // don't because we don't have information for those in the
+            // sub_constraints array above. But there is no harm in doing so:
+            // we will simply treat those in the next round around of the outer
+            // iteration.
+            //
+            // Since we want to delete some entries (see below) but don't want
+            // to disturb the order of the correspondence between lines.entries
+            // and sub_constraints, we store up which entries we will later
+            // have to delete.
+            const unsigned int     n_original_entries = line.entries.size();
+            std::set<unsigned int> entries_to_delete;
+            for (unsigned int entry = 0; entry < n_original_entries; ++entry)
+              if (sub_constraints[entry] != nullptr)
+                {
+                  // ok, this entry is further constrained:
+                  chained_constraint_replaced = true;
 
-          // Now walk through the original entries. We may replace some of them,
-          // and we may add some, but we only walk through the original entries
-          // (and also don't touch the replacements) even though in principle
-          // we could also do the ones we have added or replaced.
-          // We don't because we don't have information for those in the
-          // sub_constraints array above. But there is no harm in doing so:
-          // we will simply treat those in the next round around of the outer
-          // iteration.
-          //
-          // Since we want to delete some entries (see below) but don't want
-          // to disturb the order of the correspondence between lines.entries
-          // and sub_constraints, we store up which entries we will later
-          // have to delete.
-          const unsigned int     n_original_entries = line.entries.size();
-          std::set<unsigned int> entries_to_delete;
-          for (unsigned int entry = 0; entry < n_original_entries; ++entry)
-            if (sub_constraints[entry] != nullptr)
-              {
-                // ok, this entry is further constrained:
-                chained_constraint_replaced = true;
-
-                // look up the chain of constraints for this entry
-                const size_type dof_index = line.entries[entry].first;
-                const number    weight    = line.entries[entry].second;
-
-                Assert(dof_index != line.index,
-                       ExcMessage("Cycle in constraints detected!"));
-
-                const ConstraintLine &constrained_line =
-                  *sub_constraints[entry];
-                Assert(constrained_line.index == dof_index, ExcInternalError());
-
-                // now we have to replace an entry by its expansion. we do
-                // that by overwriting the entry by the first entry of the
-                // expansion and adding the remaining ones to the end,
-                // where we will later process them once more
-                //
-                // we can of course only do that if the DoF that we are
-                // currently handling is constrained by a linear combination
-                // of other dofs:
-                if (constrained_line.entries.size() > 0)
-                  {
-                    for (size_type i = 0; i < constrained_line.entries.size();
-                         ++i)
-                      Assert(dof_index != constrained_line.entries[i].first,
-                             ExcMessage("Cycle in constraints detected!"));
+                  // look up the chain of constraints for this entry
+                  const size_type dof_index = line.entries[entry].first;
+                  const number    weight    = line.entries[entry].second;
 
-                    // replace first entry, then tack the rest to the end
-                    // of the list
-                    line.entries[entry] = std::pair<size_type, number>(
-                      constrained_line.entries[0].first,
-                      constrained_line.entries[0].second * weight);
+                  Assert(dof_index != line.index,
+                         ExcMessage("Cycle in constraints detected!"));
 
-                    for (size_type i = 1; i < constrained_line.entries.size();
-                         ++i)
-                      line.entries.emplace_back(
-                        constrained_line.entries[i].first,
-                        constrained_line.entries[i].second * weight);
+                  const ConstraintLine &constrained_line =
+                    *sub_constraints[entry];
+                  Assert(constrained_line.index == dof_index,
+                         ExcInternalError());
+
+                  // now we have to replace an entry by its expansion. we do
+                  // that by overwriting the entry by the first entry of the
+                  // expansion and adding the remaining ones to the end,
+                  // where we will later process them once more
+                  //
+                  // we can of course only do that if the DoF that we are
+                  // currently handling is constrained by a linear combination
+                  // of other dofs:
+                  if (constrained_line.entries.size() > 0)
+                    {
+                      for (size_type i = 0; i < constrained_line.entries.size();
+                           ++i)
+                        Assert(dof_index != constrained_line.entries[i].first,
+                               ExcMessage("Cycle in constraints detected!"));
+
+                      // replace first entry, then tack the rest to the end
+                      // of the list
+                      line.entries[entry] = std::pair<size_type, number>(
+                        constrained_line.entries[0].first,
+                        constrained_line.entries[0].second * weight);
+
+                      for (size_type i = 1; i < constrained_line.entries.size();
+                           ++i)
+                        line.entries.emplace_back(
+                          constrained_line.entries[i].first,
+                          constrained_line.entries[i].second * weight);
 
 #ifdef DEBUG
-                    // keep track of how many entries we replace in this
-                    // line. If we do more than there are constraints or
-                    // dofs in our system, we must have a cycle.
-                    ++n_replacements;
-                    Assert(n_replacements / 2 < largest_idx,
-                           ExcMessage("Cycle in constraints detected!"));
+                      // keep track of how many entries we replace in this
+                      // line. If we do more than there are constraints or
+                      // dofs in our system, we must have a cycle.
+                      ++n_replacements;
+                      Assert(n_replacements / 2 < largest_idx,
+                             ExcMessage("Cycle in constraints detected!"));
 #endif
-                  }
-                else
-                  // the DoF that we encountered is not constrained by a
-                  // linear combination of other dofs but is equal to just
-                  // the inhomogeneity (i.e. its chain of entries is
-                  // empty). in that case, we can't just overwrite the
-                  // current entry, but we have to actually eliminate it
-                  {
-                    entries_to_delete.insert(entry);
-                  }
+                    }
+                  else
+                    // the DoF that we encountered is not constrained by a
+                    // linear combination of other dofs but is equal to just
+                    // the inhomogeneity (i.e. its chain of entries is
+                    // empty). in that case, we can't just overwrite the
+                    // current entry, but we have to actually eliminate it
+                    {
+                      entries_to_delete.insert(entry);
+                    }
 
-                line.inhomogeneity += constrained_line.inhomogeneity * weight;
-              }
+                  line.inhomogeneity += constrained_line.inhomogeneity * weight;
+                }
 
-          // Now delete the elements we have marked for deletion. Do
-          // so in reverse order so that we compress the array walking
-          // backward from the end without having to keep track that
-          // we have already erased earlier elements (as we would have to
-          // do if we walked the list of entries to delete in forward
-          // order).
-          for (auto it = entries_to_delete.rbegin();
-               it != entries_to_delete.rend();
-               ++it)
-            line.entries.erase(line.entries.begin() + *it);
-        }
+            // Now delete the elements we have marked for deletion. Do
+            // so in reverse order so that we compress the array walking
+            // backward from the end without having to keep track that
+            // we have already erased earlier elements (as we would have to
+            // do if we walked the list of entries to delete in forward
+            // order).
+            for (auto it = entries_to_delete.rbegin();
+                 it != entries_to_delete.rend();
+                 ++it)
+              line.entries.erase(line.entries.begin() + *it);
+          }
 
 #ifdef DEBUG
       // increase iteration count. note that we should not iterate more

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.