From 123844140974125ff9cd6db156a09a1e4e055b9d Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 16 Apr 2013 01:23:18 +0000 Subject: [PATCH] New version of SparseMatrix::copy_from for TrilinosWrappers::SparseMatrix. git-svn-id: https://svn.dealii.org/trunk@29291 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 6 ++ deal.II/include/deal.II/lac/sparse_matrix.h | 24 ++++++- .../deal.II/lac/sparse_matrix.templates.h | 49 ++++++++++++- tests/trilinos/sparse_matrix_copy_from_01.cc | 69 +++++++++++++++++++ .../sparse_matrix_copy_from_01/cmp/generic | 17 +++++ 5 files changed, 163 insertions(+), 2 deletions(-) create mode 100644 tests/trilinos/sparse_matrix_copy_from_01.cc create mode 100644 tests/trilinos/sparse_matrix_copy_from_01/cmp/generic diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 9bdbf2bdc2..f691adcc08 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -99,6 +99,12 @@ this function.

Specific improvements

    +
  1. New: There is now a version of SparseMatrix::copy_from that can copy +from TrilinosWrappers::SparseMatrix. +
    +(Wolfgang Bangerth, Jörg Frohne, 2013/04/15) +
  2. +
  3. Improved: The SolverCG implementation now uses only three auxiliary vectors, down from previously four. Also, there are some shortcuts in case PreconditionIdentity is used that improve the solver's performance. diff --git a/deal.II/include/deal.II/lac/sparse_matrix.h b/deal.II/include/deal.II/lac/sparse_matrix.h index d09bd61aac..1eacd52948 100644 --- a/deal.II/include/deal.II/lac/sparse_matrix.h +++ b/deal.II/include/deal.II/lac/sparse_matrix.h @@ -27,6 +27,14 @@ template class FullMatrix; template class BlockMatrixBase; template class SparseILU; + +#ifdef DEAL_II_WITH_TRILINOS +namespace TrilinosWrappers +{ + class SparseMatrix; +} +#endif + /** * @addtogroup Matrix1 * @{ @@ -862,7 +870,7 @@ public: void symmetrize (); /** - * Copy the given matrix to this one. The operation throws an error if the + * Copy the given matrix to this one. The operation triggers an assertion if the * sparsity patterns of the two involved matrices do not point to the same * object, since in this case the copy operation is cheaper. Since this * operation is notheless not for free, we do not make it available through @@ -906,6 +914,20 @@ public: template void copy_from (const FullMatrix &matrix); +#ifdef DEAL_II_WITH_TRILINOS + /** + * Copy the given Trilinos matrix to this one. The operation triggers an + * assertion if the sparsity patterns of the current object does not contain + * the location of a non-zero entry of the given argument. + * + * This function assumes that the two matrices have the same sizes. + * + * The function returns a reference to *this. + */ + SparseMatrix & + copy_from (const TrilinosWrappers::SparseMatrix &matrix); +#endif + /** * Add matrix scaled by factor to this matrix, i.e. the * matrix factor*matrix is added to this. This function diff --git a/deal.II/include/deal.II/lac/sparse_matrix.templates.h b/deal.II/include/deal.II/lac/sparse_matrix.templates.h index 8e75a8ad9f..2dbf6e0ee9 100644 --- a/deal.II/include/deal.II/lac/sparse_matrix.templates.h +++ b/deal.II/include/deal.II/lac/sparse_matrix.templates.h @@ -1,7 +1,7 @@ //--------------------------------------------------------------------------- // $Id$ // -// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011, 2012 by the deal.II authors +// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 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 @@ -21,6 +21,7 @@ #include #include #include +#include #include #include #include @@ -350,6 +351,52 @@ SparseMatrix::copy_from (const FullMatrix &matrix) +#ifdef DEAL_II_WITH_TRILINOS + +template +SparseMatrix & +SparseMatrix::copy_from (const TrilinosWrappers::SparseMatrix &matrix) +{ + Assert (m() == matrix.m(), ExcDimensionMismatch(m(), matrix.m())); + Assert (n() == matrix.n(), ExcDimensionMismatch(n(), matrix.n())); + + // first delete previous content + *this = 0; + + std::vector < TrilinosScalar > value_cache; + std::vector colnum_cache; + + for (unsigned int row = 0; row < matrix.m(); ++row) + { + value_cache.resize(matrix.n()); + colnum_cache.resize(matrix.n()); + + // copy column indices and values and at the same time enquire about the + // length of the row + int ncols; + int ierr + = matrix.trilinos_matrix().ExtractGlobalRowCopy + (row, matrix.row_length(row), ncols, + &(value_cache[0]), + reinterpret_cast(&(colnum_cache[0]))); + Assert (ierr==0, ExcTrilinosError(ierr)); + + // resize arrays to the size actually used + value_cache.resize(ncols); + colnum_cache.resize(ncols); + + // then copy everything + for (int i = 0; i < ncols; ++i) + this->set(row, colnum_cache[i], + value_cache[i]); + } + + return *this; +} + +#endif + + template template void diff --git a/tests/trilinos/sparse_matrix_copy_from_01.cc b/tests/trilinos/sparse_matrix_copy_from_01.cc new file mode 100644 index 0000000000..7e305b72cf --- /dev/null +++ b/tests/trilinos/sparse_matrix_copy_from_01.cc @@ -0,0 +1,69 @@ +//----------------- sparse_matrix_copy_from_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. +// +//----------------- sparse_matrix_copy_from_01.cc ------------------------- + + +// test SparseMatrix::copy_from from a TrilinosWrappers::SparseMatrix + +#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_copy_from_01/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); + + // now copy things back into a SparseMatrix + SparseMatrix copy (sparsity); + copy.copy_from (tmatrix); + + // print some output + deallog << "Copy with all values:" << std::endl; + matrix.print (deallog.get_file_stream()); + + // also compare for equality with the original + for (SparsityPattern::const_iterator + p = sparsity.begin(); p != sparsity.end(); ++p) + Assert (copy(p->row(), p->column()) == matrix(p->row(), p->column()), + ExcInternalError()); +} diff --git a/tests/trilinos/sparse_matrix_copy_from_01/cmp/generic b/tests/trilinos/sparse_matrix_copy_from_01/cmp/generic new file mode 100644 index 0000000000..6a95ff5cf1 --- /dev/null +++ b/tests/trilinos/sparse_matrix_copy_from_01/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,4) 8.00000 +(4,3) 9.00000 -- 2.39.5