]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Allow merging of ConstraintMatrix objects even if local_line is not the same
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Fri, 13 Jan 2017 17:04:41 +0000 (18:04 +0100)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sat, 14 Jan 2017 11:35:30 +0000 (12:35 +0100)
doc/news/changes/minor/20170113DanielArndt [new file with mode: 0644]
source/lac/constraint_matrix.cc
tests/lac/constraints_merge_10.cc [new file with mode: 0644]
tests/lac/constraints_merge_10.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20170113DanielArndt b/doc/news/changes/minor/20170113DanielArndt
new file mode 100644 (file)
index 0000000..bf908e2
--- /dev/null
@@ -0,0 +1,4 @@
+New: Two ConstraintMatrix objects can be merged even if their
+local_lines are not the same.
+<br>
+(Daniel Arndt, 2017/01/13)
index 90860a2c69f8e302d868203e4aaaf8ab36be1d8c..29f61125cb8a73c5eb35be0d538127b56dd45249 100644 (file)
@@ -477,17 +477,10 @@ void
 ConstraintMatrix::merge (const ConstraintMatrix &other_constraints,
                          const MergeConflictBehavior merge_conflict_behavior)
 {
-  AssertThrow(local_lines == other_constraints.local_lines,
-              ExcNotImplemented());
-
   // store the previous state with respect to sorting
   const bool object_was_sorted = sorted;
   sorted = false;
 
-  if (other_constraints.lines_cache.size() > lines_cache.size())
-    lines_cache.resize(other_constraints.lines_cache.size(),
-                       numbers::invalid_size_type);
-
   // first action is to fold into the present object possible constraints
   // in the second object. we don't strictly need to do this any more since
   // the ConstraintMatrix has learned to deal with chains of constraints in
@@ -503,15 +496,16 @@ ConstraintMatrix::merge (const ConstraintMatrix &other_constraints,
       tmp.clear ();
       for (size_type i=0; i<line->entries.size(); ++i)
         {
-          // if the present dof is not constrained, or if we won't take the
+          // if the present dof is not stored, or not constrained, or if we won't take the
           // constraint from the other object, then simply copy it over
-          if (other_constraints.is_constrained(line->entries[i].first) == false
+          if ((other_constraints.local_lines.size() != 0
+               && other_constraints.local_lines.is_element(line->entries[i].first) == false)
+              ||
+              other_constraints.is_constrained(line->entries[i].first) == false
               ||
               ((merge_conflict_behavior != right_object_wins)
-               &&
-               other_constraints.is_constrained(line->entries[i].first)
-               &&
-               this->is_constrained(line->entries[i].first)))
+               && other_constraints.is_constrained(line->entries[i].first)
+               && this->is_constrained(line->entries[i].first)))
             tmp.push_back(line->entries[i]);
           else
             // otherwise resolve further constraints by replacing the old
@@ -530,56 +524,85 @@ ConstraintMatrix::merge (const ConstraintMatrix &other_constraints,
                 tmp.push_back (std::pair<size_type,double>(j->first,
                                                            j->second*weight));
 
-              line->inhomogeneity += other_constraints.get_inhomogeneity(line->entries[i].first) *
-                                     weight;
+              line->inhomogeneity
+              += other_constraints.get_inhomogeneity(line->entries[i].first) *
+                 weight;
             }
         }
       // finally exchange old and newly resolved line
       line->entries.swap (tmp);
     }
 
-
-
-  // next action: append those lines at the end that we want to add
-  for (std::vector<ConstraintLine>::const_iterator
-       line=other_constraints.lines.begin();
-       line!=other_constraints.lines.end(); ++line)
-    if (is_constrained(line->line) == false)
-      lines.push_back (*line);
+  {
+    size_type new_size;
+    if (local_lines.size() != 0)
+      {
+        Assert (other_constraints.local_lines.size() !=0, ExcNotImplemented());
+        // The last stored index determines the new size of lines_cache
+        const size_type last_stored_index
+          = std::max(local_lines.nth_index_in_set(lines_cache.size()-1),
+                     other_constraints.local_lines.nth_index_in_set(other_constraints.lines_cache.size()-1));
+        // we don't need the old local_lines anymore, so update it
+        local_lines.add_indices(other_constraints.local_lines);
+        new_size = local_lines.index_within_set (last_stored_index)+1;
+      }
     else
       {
-        // the constrained dof we want to copy from the other object is
-        // also constrained here. let's see what we should do with that
-        switch (merge_conflict_behavior)
+        Assert (other_constraints.local_lines.size() == 0, ExcNotImplemented());
+        new_size = std::max(lines_cache.size(), other_constraints.lines_cache.size());
+      }
+    // we don't need the old lines_cache anymore
+    lines_cache.assign(new_size, numbers::invalid_size_type);
+  }
+
+  {
+    size_type index = 0;
+    // reset lines_cache for our own constraints
+    for (std::vector<ConstraintLine>::const_iterator line = lines.begin();
+         line != lines.end(); ++line)
+      lines_cache[calculate_line_index(line->line)] = index++;
+
+    // Then consider other_constraints
+    for (std::vector<ConstraintLine>::const_iterator line = other_constraints.lines.begin();
+         line != other_constraints.lines.end(); ++line)
+      {
+        const size_type local_line_no = calculate_line_index(line->line);
+        if (lines_cache[local_line_no] == numbers::invalid_size_type)
           {
-          case no_conflicts_allowed:
-            AssertThrow (false,
-                         ExcDoFIsConstrainedFromBothObjects (line->line));
-            break;
+            // there are no constraints for that line yet
+            lines.push_back(*line);
+            lines_cache[local_line_no] = index++;
+          }
+        else
+          {
+            // we already store that line
+            switch (merge_conflict_behavior)
+              {
+              case no_conflicts_allowed:
+                AssertThrow (false,
+                             ExcDoFIsConstrainedFromBothObjects (line->line));
+                break;
 
-          case left_object_wins:
-            // ignore this constraint
-            break;
+              case left_object_wins:
+                // ignore this constraint
+                break;
 
-          case right_object_wins:
-            // we need to replace the existing constraint by the one from
-            // the other object
-            lines[lines_cache[calculate_line_index(line->line)]].entries
-              = line->entries;
-            lines[lines_cache[calculate_line_index(line->line)]].inhomogeneity
-              = line->inhomogeneity;
-            break;
+              case right_object_wins:
+                lines[lines_cache[local_line_no]] = *line;
+                break;
 
-          default:
-            Assert (false, ExcNotImplemented());
+              default:
+                Assert (false, ExcNotImplemented());
+              }
           }
       }
 
-  // update the lines cache
-  size_type counter = 0;
-  for (std::vector<ConstraintLine>::const_iterator line=lines.begin();
-       line!=lines.end(); ++line, ++counter)
-    lines_cache[calculate_line_index(line->line)] = counter;
+    // check that we set the pointers correctly
+    for (size_type i=0; i<lines_cache.size(); ++i)
+      if (lines_cache[i] != numbers::invalid_size_type)
+        Assert (i == calculate_line_index(lines[lines_cache[i]].line),
+                ExcInternalError());
+  }
 
   // if the object was sorted before, then make sure it is so afterward as
   // well. otherwise leave everything in the unsorted state
diff --git a/tests/lac/constraints_merge_10.cc b/tests/lac/constraints_merge_10.cc
new file mode 100644 (file)
index 0000000..7d0fc2d
--- /dev/null
@@ -0,0 +1,118 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// This is the same as constraints_merge_08 but local_lines differs for
+// the objects to be merged.
+
+#include "../tests.h"
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/base/logstream.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+void merge_check ()
+{
+  deallog << "Checking ConstraintMatrix::merge with localized lines" << std::endl;
+
+  // set local lines to a very large range that
+  // surely triggers an error if the
+  // implementation is wrong
+  IndexSet local_lines1 (100000000);
+  local_lines1.add_range (99999890, 99999900);
+  local_lines1.add_range (99999990, 99999992);
+  local_lines1.compress();
+  local_lines1.print(std::cout);
+
+  IndexSet local_lines2 (100000000);
+  local_lines2.add_range (99999893, 99999900);
+  local_lines2.add_range (99999990, 100000000);
+  local_lines2.compress();
+  local_lines2.print(std::cout);
+
+  // the test is the same as
+  // constraints_merge_02, but we add very large
+  // indices here
+  const unsigned int index_0  = 99999890;
+  const unsigned int index_1  = 99999891;
+  const unsigned int index_3  = 99999893;
+  const unsigned int index_4  = 99999894;
+  const unsigned int index_10 = 99999990;
+  const unsigned int index_11 = 99999991;
+  const unsigned int index_12 = 99999992;
+  const unsigned int index_13 = 99999993;
+
+  // check twice, once with closed
+  // objects, once with open ones
+  for (unsigned int run=0; run<2; ++run)
+    {
+      deallog << "Checking with " << (run == 0 ? "open" : "closed")
+              << " objects" << std::endl;
+
+      // check that the `merge' function
+      // works correctly
+      ConstraintMatrix c1 (local_lines1), c2 (local_lines2);
+
+      // enter simple line
+      c1.add_line          (index_0);
+      c1.add_entry         (index_0, index_11, 1.);
+      c1.set_inhomogeneity (index_0, 42);
+
+      // add more complex line
+      c1.add_line          (index_1);
+      c1.add_entry         (index_1, index_3,  0.5);
+      c1.add_entry         (index_1, index_4,  0.5);
+      c1.set_inhomogeneity (index_1, 100);
+
+      // fill second constraints
+      // object with one trivial line
+      // and one which further
+      // constrains one of the
+      // entries in the first object
+      c2.add_line          (index_10);
+      c2.add_entry         (index_10, index_11, 1.);
+      c2.set_inhomogeneity (index_10, 142);
+
+      c2.add_line          (index_3);
+      c2.add_entry         (index_3, index_12, 0.25);
+      c2.add_entry         (index_3, index_13, 0.75);
+      c2.set_inhomogeneity (index_3, 242);
+
+      // in one of the two runs,
+      // close the objects
+      if (run == 1)
+        {
+          c1.close ();
+          c2.close ();
+        };
+
+      // now merge the two and print the
+      // results
+      c1.merge (c2);
+      c1.print(deallog.get_file_stream());
+    };
+}
+
+
+int main ()
+{
+  initlog();
+
+  merge_check ();
+}
+
diff --git a/tests/lac/constraints_merge_10.output b/tests/lac/constraints_merge_10.output
new file mode 100644 (file)
index 0000000..981ff5c
--- /dev/null
@@ -0,0 +1,26 @@
+
+DEAL::Checking ConstraintMatrix::merge with localized lines
+DEAL::Checking with open objects
+    99999890 99999991:  1.00000
+    99999890: 42.0000
+    99999891 99999992:  0.125000
+    99999891 99999993:  0.375000
+    99999891 99999894:  0.500000
+    99999891: 221.000
+    99999990 99999991:  1.00000
+    99999990: 142.000
+    99999893 99999992:  0.250000
+    99999893 99999993:  0.750000
+    99999893: 242.000
+DEAL::Checking with closed objects
+    99999890 99999991:  1.00000
+    99999890: 42.0000
+    99999891 99999894:  0.500000
+    99999891 99999992:  0.125000
+    99999891 99999993:  0.375000
+    99999891: 221.000
+    99999893 99999992:  0.250000
+    99999893 99999993:  0.750000
+    99999893: 242.000
+    99999990 99999991:  1.00000
+    99999990: 142.000

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.