]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Test ConstraintMatrix::merge() for constraint matrices localized by an IndexSet ...
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 27 Jun 2011 16:32:10 +0000 (16:32 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 27 Jun 2011 16:32:10 +0000 (16:32 +0000)
git-svn-id: https://svn.dealii.org/trunk@23868 0785d39b-7218-0410-832d-ea1e28bc413d

tests/deal.II/constraints_merge_08.cc [new file with mode: 0644]
tests/deal.II/constraints_merge_08/cmp/generic [new file with mode: 0644]
tests/deal.II/constraints_merge_09.cc [new file with mode: 0644]
tests/deal.II/constraints_merge_09/cmp/generic [new file with mode: 0644]

diff --git a/tests/deal.II/constraints_merge_08.cc b/tests/deal.II/constraints_merge_08.cc
new file mode 100644 (file)
index 0000000..3b5d317
--- /dev/null
@@ -0,0 +1,119 @@
+//----------------------------  constraints_merge_08.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2003, 2004, 2005, 2008, 2010, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  constraints_merge_08.cc  ---------------------------
+
+
+// merge and print a bunch of ConstrainMatrices. test the case that we have
+// inhomogeneities in the constraints and the constraint matrix is constructed
+// based on an IndexSet to transform large global indices into local ones
+
+#include "../tests.h"
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/base/logstream.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+std::ofstream logfile("constraints_merge_08/output");
+
+
+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_lines (100000000);
+  local_lines.add_range (99999890, 99999900);
+  local_lines.add_range (99999990,100000000);
+  local_lines.compress();
+
+                               // the test is the same as
+                               // constraints_merge_02, but we add very large
+                               // indices here
+  const unsigned int index_0  = local_lines.nth_index_in_set(0);
+  const unsigned int index_1  = local_lines.nth_index_in_set(1);
+  const unsigned int index_3  = local_lines.nth_index_in_set(3);
+  const unsigned int index_4  = local_lines.nth_index_in_set(4);
+  const unsigned int index_10 = local_lines.nth_index_in_set(10);
+  const unsigned int index_11 = local_lines.nth_index_in_set(11);
+  const unsigned int index_12 = local_lines.nth_index_in_set(12);
+  const unsigned int index_13 = local_lines.nth_index_in_set(13);
+
+  deallog << "Number of local lines: "
+         << local_lines.n_elements() << std::endl;
+
+                                  // 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_lines), c2 (local_lines);
+
+                                      // 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 (logfile);
+    };
+}
+
+
+int main ()
+{
+  deallog << std::setprecision (2);
+  logfile << std::setprecision (2);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  merge_check ();
+}
+
diff --git a/tests/deal.II/constraints_merge_08/cmp/generic b/tests/deal.II/constraints_merge_08/cmp/generic
new file mode 100644 (file)
index 0000000..3d7fa02
--- /dev/null
@@ -0,0 +1,27 @@
+
+DEAL::Checking ConstraintMatrix::merge with localized lines
+DEAL::Number of local lines: 20
+DEAL::Checking with open objects
+    99999890 99999991:  1.0
+    99999890: 42.
+    99999891 99999992:  0.12
+    99999891 99999993:  0.38
+    99999891 99999894:  0.50
+    99999891: 2.2e+02
+    99999990 99999991:  1.0
+    99999990: 1.4e+02
+    99999893 99999992:  0.25
+    99999893 99999993:  0.75
+    99999893: 2.4e+02
+DEAL::Checking with closed objects
+    99999890 99999991:  1.0
+    99999890: 42.
+    99999891 99999894:  0.50
+    99999891 99999992:  0.12
+    99999891 99999993:  0.38
+    99999891: 2.2e+02
+    99999893 99999992:  0.25
+    99999893 99999993:  0.75
+    99999893: 2.4e+02
+    99999990 99999991:  1.0
+    99999990: 1.4e+02
diff --git a/tests/deal.II/constraints_merge_09.cc b/tests/deal.II/constraints_merge_09.cc
new file mode 100644 (file)
index 0000000..d0c184a
--- /dev/null
@@ -0,0 +1,114 @@
+//----------------------------  constraints_merge_09.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2003, 2004, 2005, 2008, 2010, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  constraints_merge_09.cc  ---------------------------
+
+
+// merge and print a bunch of ConstrainMatrices. test the case where we have
+// conflicting constraints and the right object wins (winning constraint
+// introduces a chain)
+//
+// like _06, but using localized lines
+
+#include "../tests.h"
+#include <deal.II/lac/constraint_matrix.h>
+#include <deal.II/base/logstream.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+std::ofstream logfile("constraints_merge_09/output");
+
+
+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_lines (100000000);
+  local_lines.add_range (99999890, 99999900);
+  local_lines.add_range (99999990,100000000);
+  local_lines.compress();
+
+                               // the test is the same as
+                               // constraints_merge_02, but we add very large
+                               // indices here
+  const unsigned int index_0  = local_lines.nth_index_in_set(0);
+  const unsigned int index_2  = local_lines.nth_index_in_set(2);
+  const unsigned int index_11 = local_lines.nth_index_in_set(11);
+  const unsigned int index_13 = local_lines.nth_index_in_set(13);
+
+  deallog << "Number of local lines: "
+         << local_lines.n_elements() << std::endl;
+
+                                  // 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_lines), c2 (local_lines);
+
+                                      // enter simple line
+      c1.add_line          (index_0);
+      c1.add_entry         (index_0, index_11, 1.);
+      c1.set_inhomogeneity (index_0, 42);
+
+      c1.add_line          (index_13);
+      c1.add_entry         (index_13, index_2, 0.5);
+      c1.set_inhomogeneity (index_13, 2);
+
+                                      // fill second constraints
+                                      // object that has a conflict
+      c2.add_line          (index_0);
+      c2.add_entry         (index_0, index_13, 2.);
+      c2.set_inhomogeneity (index_0, 142);
+                                      // in one of the two runs,
+                                      // close the objects
+      if (run == 1)
+       {
+         c1.close ();
+         c2.close ();
+       };
+
+                                      // now merge the two and print the
+                                      // results
+      try
+       {
+         c1.merge (c2, ConstraintMatrix::right_object_wins);
+       }
+      catch (...)
+       {
+         Assert (false, ExcInternalError());
+       }
+
+      c1.print (logfile);
+    }
+}
+
+
+int main ()
+{
+  deallog << std::setprecision (2);
+  logfile << std::setprecision (2);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  merge_check ();
+}
+
diff --git a/tests/deal.II/constraints_merge_09/cmp/generic b/tests/deal.II/constraints_merge_09/cmp/generic
new file mode 100644 (file)
index 0000000..477b6c0
--- /dev/null
@@ -0,0 +1,13 @@
+
+DEAL::Checking ConstraintMatrix::merge with localized lines
+DEAL::Number of local lines: 20
+DEAL::Checking with open objects
+    99999890 99999993:  2.0
+    99999890: 1.4e+02
+    99999993 99999892:  0.50
+    99999993: 2.0
+DEAL::Checking with closed objects
+    99999890 99999892:  1.0
+    99999890: 1.5e+02
+    99999993 99999892:  0.50
+    99999993: 2.0

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.