From: bangerth Date: Thu, 17 Jan 2013 15:50:11 +0000 (+0000) Subject: Exhaustively test one function that previously wasn't tested. As always if one does... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=352fd9a0ec8fcb44d1447e95e50afdf7c870132b;p=dealii-svn.git Exhaustively test one function that previously wasn't tested. As always if one does so there are bugs to be found. git-svn-id: https://svn.dealii.org/trunk@28099 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index e3c0fcd829..8406fdc8d3 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -177,6 +177,14 @@ DoFHandler, in particular removal of specializations.

Specific improvements

    +
  1. Fixed: Various variants of the TrilinosWrappers::SparseMatrix::reinit +functions take a parameter drop_tolerance 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. +
    +(Wolfgang Bangerth, 2013/01/17) +
  2. New: All vector classes should now have a in_local_range() function indicating whether a given index is locally stored or not.
    diff --git a/deal.II/source/lac/trilinos_sparse_matrix.cc b/deal.II/source/lac/trilinos_sparse_matrix.cc index f7bc653e70..f9ee921ec8 100644 --- a/deal.II/source/lac/trilinos_sparse_matrix.cc +++ b/deal.II/source/lac/trilinos_sparse_matrix.cc @@ -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 index 0000000000..269ce40575 --- /dev/null +++ b/tests/trilinos/sparse_matrix_02.cc @@ -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 +#include +#include +#include +#include +#include + + +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 matrix(sparsity); + { + double value = 1; + for (SparseMatrix::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 index 0000000000..abdb4e34ea --- /dev/null +++ b/tests/trilinos/sparse_matrix_02/cmp/generic @@ -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 index 0000000000..f37d844bd3 --- /dev/null +++ b/tests/trilinos/sparse_matrix_03.cc @@ -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 +#include +#include +#include +#include +#include + + +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 matrix(sparsity); + { + double value = 1; + for (SparseMatrix::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 index 0000000000..2ec3dc60db --- /dev/null +++ b/tests/trilinos/sparse_matrix_03/cmp/generic @@ -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 index 0000000000..9b03e2da7a --- /dev/null +++ b/tests/trilinos/sparse_matrix_04.cc @@ -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 +#include +#include +#include +#include +#include + + +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 matrix(sparsity); + { + double value = 1; + for (SparseMatrix::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 index 0000000000..add5135bf7 --- /dev/null +++ b/tests/trilinos/sparse_matrix_04/cmp/generic @@ -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 index 0000000000..eb3849910a --- /dev/null +++ b/tests/trilinos/sparse_matrix_05.cc @@ -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 +#include +#include +#include +#include +#include + + +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 matrix(sparsity); + { + double value = 1; + for (SparseMatrix::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 index 0000000000..832242a4cc --- /dev/null +++ b/tests/trilinos/sparse_matrix_05/cmp/generic @@ -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 index 0000000000..73f32bc704 --- /dev/null +++ b/tests/trilinos/sparse_matrix_06.cc @@ -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 +#include +#include +#include +#include +#include + + +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 matrix(sparsity); + { + double value = 1; + for (SparseMatrix::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 index 0000000000..98e27e0d30 --- /dev/null +++ b/tests/trilinos/sparse_matrix_06/cmp/generic @@ -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 index 0000000000..7c3754e75c --- /dev/null +++ b/tests/trilinos/sparse_matrix_07.cc @@ -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 +#include +#include +#include +#include +#include + + +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 matrix(sparsity); + { + double value = 1; + for (SparseMatrix::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 index 0000000000..372283ed0a --- /dev/null +++ b/tests/trilinos/sparse_matrix_07/cmp/generic @@ -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