]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added a test case for the Tvmult method of the SparseILU class. It is
authoroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 15 Sep 2005 10:17:19 +0000 (10:17 +0000)
committeroliver <oliver@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 15 Sep 2005 10:17:19 +0000 (10:17 +0000)
essentially the same code as sparse_ilu.cc with the difference that
the Tvmult method is tested instead of the vmult method.

git-svn-id: https://svn.dealii.org/trunk@11445 0785d39b-7218-0410-832d-ea1e28bc413d

tests/lac/Makefile
tests/lac/sparse_ilu_t.cc [new file with mode: 0644]

index d1f09eba66f0ad6b95f09f0feb438c5e1e88478b..3cae089c98bf0dc9e9c8ca5f5057563b70a7be8f 100644 (file)
@@ -26,7 +26,7 @@ tests_x = vector-vector \
        block_vector block_vector_iterator block_matrices \
        matrix_lib matrix_out \
        solver eigen \
-       sparse_ilu lapack
+       sparse_ilu sparse_ilu_t lapack
 
 
 # from above list of regular expressions, generate the real set of
diff --git a/tests/lac/sparse_ilu_t.cc b/tests/lac/sparse_ilu_t.cc
new file mode 100644 (file)
index 0000000..61ece33
--- /dev/null
@@ -0,0 +1,125 @@
+//----------------------------  sparse_ilu_t.cc  -------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2001, 2002, 2003, 2004, 2005 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.
+//
+//----------------------------  sparse_ilu_t.cc  -------------------------
+
+
+// make sure that the SparseILU applied with infinite fill-in
+// generates the exact inverse matrix
+
+#include "../tests.h"
+#include <cmath>
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include <cstdlib>
+#include "testmatrix.h"
+#include <base/logstream.h>
+#include <lac/sparse_matrix.h>
+#include <lac/sparse_ilu.h>
+#include <lac/vector.h>
+
+//TODO:[WB] find test that is less sensitive to floating point accuracy
+
+int main()
+{
+  std::ofstream logfile("sparse_ilu_t.output");
+  logfile.setf(std::ios::fixed);
+  logfile.precision(3);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+  
+
+  for (unsigned int size=4; size <= 16; size *= 2)
+    {
+      unsigned int dim = (size-1)*(size-1);
+
+      deallog << "Size " << size << " Unknowns " << dim << std::endl;
+      
+                                      // Make matrix
+      FDMatrix testproblem(size, size);
+      SparsityPattern structure(dim, dim, 5);
+      testproblem.five_point_structure(structure);
+      structure.compress();
+      SparseMatrix<double>  A(structure);
+      testproblem.five_point(A);
+
+      
+      for (unsigned int test=0; test<2; ++test)
+       {
+         deallog << "Test " << test << std::endl;
+         
+                                          // generate sparse ILU.
+                                          //
+                                          // for test 1, test with
+                                          // full pattern.  for test
+                                          // 2, test with same
+                                          // pattern as A
+         SparsityPattern ilu_pattern (dim, dim,
+                                      (test==0 ? dim : 5));
+         switch (test)
+           {
+             case 0:
+                   for (unsigned int i=0; i<dim; ++i)
+                     for (unsigned int j=0; j<dim; ++j)
+                       ilu_pattern.add(i,j);
+                   break;
+
+             case 1:
+                   for (unsigned int i=0; i<dim; ++i)
+                     for (unsigned int j=0; j<dim; ++j)
+                       if (structure(i,j) != SparsityPattern::invalid_entry)
+                         ilu_pattern.add(i,j);
+                   break;
+
+             default:
+                   Assert (false, ExcNotImplemented());
+           };
+         ilu_pattern.compress();
+         SparseILU<double> ilu (ilu_pattern);
+         ilu.decompose (A);
+         
+                                          // now for three test vectors v
+                                          // determine norm of
+                                          // (I-B'A')v, where B' is the ILU 
+                                          // of A transposed.
+                                          // Since matrix is symmetric,
+                                          // likewise test for right
+                                          // preconditioner
+         Vector<double> v(dim);
+         Vector<double> tmp1(dim), tmp2(dim);
+         for (unsigned int i=0; i<3; ++i)
+           {
+             for (unsigned int j=0; j<dim; ++j)
+               v(j) = 1. * std::rand()/RAND_MAX;
+             
+             A.Tvmult (tmp1, v);
+             ilu.Tvmult (tmp2, tmp1);
+             tmp2 -= v;
+             const double left_residual = tmp2.l2_norm();
+             
+             ilu.Tvmult (tmp1, v);
+             A.Tvmult (tmp2, tmp1);
+             tmp2 -= v;
+             const double right_residual = tmp2.l2_norm();
+             
+             
+             deallog << "Residual with test vector " << i << ":  "
+                     << " left=" << left_residual
+                     << ", right=" << right_residual
+                     << std::endl;
+           };
+       };
+      
+    };
+}
+

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.