]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Exhaustively test one function that previously wasn't tested. As always if one does...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 17 Jan 2013 15:50:11 +0000 (15:50 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 17 Jan 2013 15:50:11 +0000 (15:50 +0000)
git-svn-id: https://svn.dealii.org/trunk@28099 0785d39b-7218-0410-832d-ea1e28bc413d

14 files changed:
deal.II/doc/news/changes.h
deal.II/source/lac/trilinos_sparse_matrix.cc
tests/trilinos/sparse_matrix_02.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_02/cmp/generic [new file with mode: 0644]
tests/trilinos/sparse_matrix_03.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_03/cmp/generic [new file with mode: 0644]
tests/trilinos/sparse_matrix_04.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_04/cmp/generic [new file with mode: 0644]
tests/trilinos/sparse_matrix_05.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_05/cmp/generic [new file with mode: 0644]
tests/trilinos/sparse_matrix_06.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_06/cmp/generic [new file with mode: 0644]
tests/trilinos/sparse_matrix_07.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_07/cmp/generic [new file with mode: 0644]

index e3c0fcd829fe3d0fb1d89a05f1aa62f6089aea23..8406fdc8d38d3f746c022a4f6d8537d45088d0d4 100644 (file)
@@ -177,6 +177,14 @@ DoFHandler, in particular removal of specializations.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> Fixed: Various variants of the TrilinosWrappers::SparseMatrix::reinit
+functions take a parameter <code>drop_tolerance</code> that allows to remove
+small values from the matrix and replace them by zero instead. This was not
+enforced for values on the diagonal of the matrix but only for off-diagonal
+ones. This is now fixed.
+<br>
+(Wolfgang Bangerth, 2013/01/17)
+
 <li> New: All vector classes should now have a in_local_range()
 function indicating whether a given index is locally stored or not.
 <br>
index f7bc653e70578c845edef9e48239ea37a92dd37b..f9ee921ec88acc41d0f8feb8856127db48879e4e 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2008, 2009, 2010, 2011, 2012 by the deal.II authors
+//    Copyright (C) 2008, 2009, 2010, 2011, 2012, 2013 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -623,8 +623,11 @@ set_matrix_values:
                                 "diagonal, the deal.II matrix must optimize it,"
                                 " too. And vice versa."));
               AssertDimension(it->column(), row);
-              values[col] = it->value();
-              row_indices[col++] = it->column();
+              if (std::fabs(it->value()) > drop_tolerance)
+               {
+                 values[col] = it->value();
+                 row_indices[col++] = it->column();
+               }
               ++select_index;
               ++it;
             }
diff --git a/tests/trilinos/sparse_matrix_02.cc b/tests/trilinos/sparse_matrix_02.cc
new file mode 100644 (file)
index 0000000..269ce40
--- /dev/null
@@ -0,0 +1,59 @@
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2011, 2013 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.
+//
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+
+
+// test TrilinosWrappers::SparseMatrix::reinit with a dealii::SparseMatrix
+// and where we copy all values
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <fstream>
+#include <iostream>
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv);
+
+  std::ofstream logfile("sparse_matrix_02/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  SparsityPattern sparsity (5,5,5);
+  sparsity.add (1,2);
+  sparsity.add (2,3);
+  sparsity.add (3,4);
+  sparsity.add (4,3);
+  sparsity.compress();
+  SparseMatrix<double> matrix(sparsity);
+  {
+    double value = 1;
+    for (SparseMatrix<double>::iterator p=matrix.begin();
+        p != matrix.end(); ++p, ++value)
+      p->value() = value;
+  }
+  deallog << "Original:" << std::endl;
+  matrix.print_formatted (deallog.get_file_stream());
+
+  // now copy everything into a Trilinos matrix
+  Epetra_Map map(5,5,0,Utilities::Trilinos::comm_world());
+  TrilinosWrappers::SparseMatrix tmatrix;
+  tmatrix.reinit (map, map, matrix);
+
+  deallog << "Copy with all values:" << std::endl;
+  tmatrix.print (deallog.get_file_stream());
+}
diff --git a/tests/trilinos/sparse_matrix_02/cmp/generic b/tests/trilinos/sparse_matrix_02/cmp/generic
new file mode 100644 (file)
index 0000000..abdb4e3
--- /dev/null
@@ -0,0 +1,17 @@
+
+DEAL::Original:
+1.000e+00                                              
+           2.000e+00  3.000e+00                        
+                      4.000e+00  5.000e+00             
+                                 6.000e+00  7.000e+00  
+                                 9.000e+00  8.000e+00  
+DEAL::Copy with all values:
+(0,0) 1.00000
+(1,1) 2.00000
+(1,2) 3.00000
+(2,2) 4.00000
+(2,3) 5.00000
+(3,3) 6.00000
+(3,4) 7.00000
+(4,3) 9.00000
+(4,4) 8.00000
diff --git a/tests/trilinos/sparse_matrix_03.cc b/tests/trilinos/sparse_matrix_03.cc
new file mode 100644 (file)
index 0000000..f37d844
--- /dev/null
@@ -0,0 +1,59 @@
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2011, 2013 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.
+//
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+
+
+// test TrilinosWrappers::SparseMatrix::reinit with a dealii::SparseMatrix
+// with a drop tolerance
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <fstream>
+#include <iostream>
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv);
+
+  std::ofstream logfile("sparse_matrix_03/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  SparsityPattern sparsity (5,5,5);
+  sparsity.add (1,2);
+  sparsity.add (2,3);
+  sparsity.add (3,4);
+  sparsity.add (4,3);
+  sparsity.compress();
+  SparseMatrix<double> matrix(sparsity);
+  {
+    double value = 1;
+    for (SparseMatrix<double>::iterator p=matrix.begin();
+        p != matrix.end(); ++p, ++value)
+      p->value() = value;
+  }
+  deallog << "Original:" << std::endl;
+  matrix.print_formatted (deallog.get_file_stream());
+
+  // now copy everything into a Trilinos matrix
+  Epetra_Map map(5,5,0,Utilities::Trilinos::comm_world());
+  TrilinosWrappers::SparseMatrix tmatrix;
+  tmatrix.reinit (map, map, matrix, 4.);
+
+  deallog << "Copy with a drop tolerance of 4:" << std::endl;
+  tmatrix.print (deallog.get_file_stream());
+}
diff --git a/tests/trilinos/sparse_matrix_03/cmp/generic b/tests/trilinos/sparse_matrix_03/cmp/generic
new file mode 100644 (file)
index 0000000..2ec3dc6
--- /dev/null
@@ -0,0 +1,17 @@
+
+DEAL::Original:
+1.000e+00                                              
+           2.000e+00  3.000e+00                        
+                      4.000e+00  5.000e+00             
+                                 6.000e+00  7.000e+00  
+                                 9.000e+00  8.000e+00  
+DEAL::Copy with a drop tolerance of 4:
+(0,0) 0.00000
+(1,1) 0.00000
+(1,2) 0.00000
+(2,2) 0.00000
+(2,3) 5.00000
+(3,3) 6.00000
+(3,4) 7.00000
+(4,3) 9.00000
+(4,4) 8.00000
diff --git a/tests/trilinos/sparse_matrix_04.cc b/tests/trilinos/sparse_matrix_04.cc
new file mode 100644 (file)
index 0000000..9b03e2d
--- /dev/null
@@ -0,0 +1,59 @@
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2011, 2013 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.
+//
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+
+
+// test TrilinosWrappers::SparseMatrix::reinit with a dealii::SparseMatrix
+// without copying values
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <fstream>
+#include <iostream>
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv);
+
+  std::ofstream logfile("sparse_matrix_04/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  SparsityPattern sparsity (5,5,5);
+  sparsity.add (1,2);
+  sparsity.add (2,3);
+  sparsity.add (3,4);
+  sparsity.add (4,3);
+  sparsity.compress();
+  SparseMatrix<double> matrix(sparsity);
+  {
+    double value = 1;
+    for (SparseMatrix<double>::iterator p=matrix.begin();
+        p != matrix.end(); ++p, ++value)
+      p->value() = value;
+  }
+  deallog << "Original:" << std::endl;
+  matrix.print_formatted (deallog.get_file_stream());
+
+  // now copy everything into a Trilinos matrix
+  Epetra_Map map(5,5,0,Utilities::Trilinos::comm_world());
+  TrilinosWrappers::SparseMatrix tmatrix;
+  tmatrix.reinit (map, map, matrix, 0, false);
+
+  deallog << "Copy structure only:" << std::endl;
+  tmatrix.print (deallog.get_file_stream());
+}
diff --git a/tests/trilinos/sparse_matrix_04/cmp/generic b/tests/trilinos/sparse_matrix_04/cmp/generic
new file mode 100644 (file)
index 0000000..add5135
--- /dev/null
@@ -0,0 +1,17 @@
+
+DEAL::Original:
+1.000e+00                                              
+           2.000e+00  3.000e+00                        
+                      4.000e+00  5.000e+00             
+                                 6.000e+00  7.000e+00  
+                                 9.000e+00  8.000e+00  
+DEAL::Copy structure only:
+(0,0) 0.00000
+(1,1) 0.00000
+(1,2) 0.00000
+(2,2) 0.00000
+(2,3) 0.00000
+(3,3) 0.00000
+(3,4) 0.00000
+(4,3) 0.00000
+(4,4) 0.00000
diff --git a/tests/trilinos/sparse_matrix_05.cc b/tests/trilinos/sparse_matrix_05.cc
new file mode 100644 (file)
index 0000000..eb38499
--- /dev/null
@@ -0,0 +1,66 @@
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2011, 2013 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.
+//
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+
+
+// test TrilinosWrappers::SparseMatrix::reinit with a dealii::SparseMatrix
+// with a separate sparsity pattern that is a subset
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <fstream>
+#include <iostream>
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv);
+
+  std::ofstream logfile("sparse_matrix_05/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  SparsityPattern sparsity (5,5,5);
+  sparsity.add (1,2);
+  sparsity.add (2,3);
+  sparsity.add (3,4);
+  sparsity.add (4,3);
+  sparsity.compress();
+  SparseMatrix<double> matrix(sparsity);
+  {
+    double value = 1;
+    for (SparseMatrix<double>::iterator p=matrix.begin();
+        p != matrix.end(); ++p, ++value)
+      p->value() = value;
+  }
+  deallog << "Original:" << std::endl;
+  matrix.print_formatted (deallog.get_file_stream());
+
+  // create a separate sparsity pattern to use
+  SparsityPattern xsparsity (5,5,5);
+  xsparsity.add (1,2);
+  xsparsity.add (2,3);
+  xsparsity.compress();
+
+
+  // now copy everything into a Trilinos matrix
+  Epetra_Map map(5,5,0,Utilities::Trilinos::comm_world());
+  TrilinosWrappers::SparseMatrix tmatrix;
+  tmatrix.reinit (map, map, matrix, 0, true, &xsparsity);
+
+  deallog << "Copy structure only:" << std::endl;
+  tmatrix.print (deallog.get_file_stream());
+}
diff --git a/tests/trilinos/sparse_matrix_05/cmp/generic b/tests/trilinos/sparse_matrix_05/cmp/generic
new file mode 100644 (file)
index 0000000..832242a
--- /dev/null
@@ -0,0 +1,15 @@
+
+DEAL::Original:
+1.000e+00                                              
+           2.000e+00  3.000e+00                        
+                      4.000e+00  5.000e+00             
+                                 6.000e+00  7.000e+00  
+                                 9.000e+00  8.000e+00  
+DEAL::Copy structure only:
+(0,0) 1.00000
+(1,1) 2.00000
+(1,2) 3.00000
+(2,2) 4.00000
+(2,3) 5.00000
+(3,3) 6.00000
+(4,4) 8.00000
diff --git a/tests/trilinos/sparse_matrix_06.cc b/tests/trilinos/sparse_matrix_06.cc
new file mode 100644 (file)
index 0000000..73f32bc
--- /dev/null
@@ -0,0 +1,67 @@
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2011, 2013 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.
+//
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+
+
+// test TrilinosWrappers::SparseMatrix::reinit with a dealii::SparseMatrix
+// with a separate sparsity pattern that is a subset and for which we don't
+// even store the diagonal
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <fstream>
+#include <iostream>
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv);
+
+  std::ofstream logfile("sparse_matrix_06/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  SparsityPattern sparsity (5,5,5);
+  sparsity.add (1,2);
+  sparsity.add (2,3);
+  sparsity.add (3,4);
+  sparsity.add (4,3);
+  sparsity.compress();
+  SparseMatrix<double> matrix(sparsity);
+  {
+    double value = 1;
+    for (SparseMatrix<double>::iterator p=matrix.begin();
+        p != matrix.end(); ++p, ++value)
+      p->value() = value;
+  }
+  deallog << "Original:" << std::endl;
+  matrix.print_formatted (deallog.get_file_stream());
+
+  // create a separate sparsity pattern to use
+  SparsityPattern xsparsity (5,5,5,/*optimize_diagonal=*/false);
+  xsparsity.add (1,2);
+  xsparsity.add (2,3);
+  xsparsity.compress();
+
+
+  // now copy everything into a Trilinos matrix
+  Epetra_Map map(5,5,0,Utilities::Trilinos::comm_world());
+  TrilinosWrappers::SparseMatrix tmatrix;
+  tmatrix.reinit (map, map, matrix, 0, true, &xsparsity);
+
+  deallog << "Copy structure only:" << std::endl;
+  tmatrix.print (deallog.get_file_stream());
+}
diff --git a/tests/trilinos/sparse_matrix_06/cmp/generic b/tests/trilinos/sparse_matrix_06/cmp/generic
new file mode 100644 (file)
index 0000000..98e27e0
--- /dev/null
@@ -0,0 +1,10 @@
+
+DEAL::Original:
+1.000e+00                                              
+           2.000e+00  3.000e+00                        
+                      4.000e+00  5.000e+00             
+                                 6.000e+00  7.000e+00  
+                                 9.000e+00  8.000e+00  
+DEAL::Copy structure only:
+(1,2) 3.00000
+(2,3) 5.00000
diff --git a/tests/trilinos/sparse_matrix_07.cc b/tests/trilinos/sparse_matrix_07.cc
new file mode 100644 (file)
index 0000000..7c3754e
--- /dev/null
@@ -0,0 +1,68 @@
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2011, 2013 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.
+//
+//-----------------  trilinos_sparse_matrix_01.cc  -------------------------
+
+
+// test TrilinosWrappers::SparseMatrix::reinit with a dealii::SparseMatrix
+// with a separate sparsity pattern that is partly subset and partly superset
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparsity_pattern.h>
+#include <deal.II/lac/sparse_matrix.h>
+#include <fstream>
+#include <iostream>
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv);
+
+  std::ofstream logfile("sparse_matrix_07/output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  SparsityPattern sparsity (5,5,5);
+  sparsity.add (1,2);
+  sparsity.add (2,3);
+  sparsity.add (3,4);
+  sparsity.add (4,3);
+  sparsity.compress();
+  SparseMatrix<double> matrix(sparsity);
+  {
+    double value = 1;
+    for (SparseMatrix<double>::iterator p=matrix.begin();
+        p != matrix.end(); ++p, ++value)
+      p->value() = value;
+  }
+  deallog << "Original:" << std::endl;
+  matrix.print_formatted (deallog.get_file_stream());
+
+  // create a separate sparsity pattern to use
+  SparsityPattern xsparsity (5,5,5,/*optimize_diagonal=*/false);
+  xsparsity.add (1,2);
+  xsparsity.add (2,3);
+  xsparsity.add (2,4);
+  xsparsity.add (2,1);
+  xsparsity.compress();
+
+
+  // now copy everything into a Trilinos matrix
+  Epetra_Map map(5,5,0,Utilities::Trilinos::comm_world());
+  TrilinosWrappers::SparseMatrix tmatrix;
+  tmatrix.reinit (map, map, matrix, 0, true, &xsparsity);
+
+  deallog << "Copy structure only:" << std::endl;
+  tmatrix.print (deallog.get_file_stream());
+}
diff --git a/tests/trilinos/sparse_matrix_07/cmp/generic b/tests/trilinos/sparse_matrix_07/cmp/generic
new file mode 100644 (file)
index 0000000..372283e
--- /dev/null
@@ -0,0 +1,12 @@
+
+DEAL::Original:
+1.000e+00                                              
+           2.000e+00  3.000e+00                        
+                      4.000e+00  5.000e+00             
+                                 6.000e+00  7.000e+00  
+                                 9.000e+00  8.000e+00  
+DEAL::Copy structure only:
+(1,2) 3.00000
+(2,1) 0.00000
+(2,3) 5.00000
+(2,4) 0.00000

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.