]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New version of SparseMatrix::copy_from for TrilinosWrappers::SparseMatrix.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Apr 2013 01:23:18 +0000 (01:23 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 16 Apr 2013 01:23:18 +0000 (01:23 +0000)
git-svn-id: https://svn.dealii.org/trunk@29291 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/lac/sparse_matrix.h
deal.II/include/deal.II/lac/sparse_matrix.templates.h
tests/trilinos/sparse_matrix_copy_from_01.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_copy_from_01/cmp/generic [new file with mode: 0644]

index 9bdbf2bdc291339ee04cc4df1ad264bc765a21a3..f691adcc0822a8c8d0d0e7973f55d7b871b96751 100644 (file)
@@ -99,6 +99,12 @@ this function.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: There is now a version of SparseMatrix::copy_from that can copy
+from TrilinosWrappers::SparseMatrix.
+<br>
+(Wolfgang Bangerth, J&ouml;rg Frohne, 2013/04/15)
+</li>
+
 <li> 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.
index d09bd61aacd325dc09a0a9a0cff73b2d312414b7..1eacd52948841ae410d2e8faf288d0e1c89ad014 100644 (file)
@@ -27,6 +27,14 @@ template <typename number> class FullMatrix;
 template <typename Matrix> class BlockMatrixBase;
 template <typename number> 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 <typename somenumber>
   void copy_from (const FullMatrix<somenumber> &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 <tt>*this</tt>.
+   */
+  SparseMatrix<number> &
+  copy_from (const TrilinosWrappers::SparseMatrix &matrix);
+#endif
+  
   /**
    * Add <tt>matrix</tt> scaled by <tt>factor</tt> to this matrix, i.e. the
    * matrix <tt>factor*matrix</tt> is added to <tt>this</tt>. This function
index 8e75a8ad9f97277f106abc0c91cca034343fafd0..2dbf6e0ee9bc68a762ffc5f5340cdd9420147fb3 100644 (file)
@@ -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 <deal.II/base/multithread_info.h>
 #include <deal.II/base/utilities.h>
 #include <deal.II/lac/sparse_matrix.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
 #include <deal.II/lac/vector.h>
 #include <deal.II/lac/full_matrix.h>
 #include <deal.II/lac/compressed_simple_sparsity_pattern.h>
@@ -350,6 +351,52 @@ SparseMatrix<number>::copy_from (const FullMatrix<somenumber> &matrix)
 
 
 
+#ifdef DEAL_II_WITH_TRILINOS
+
+template <typename number>
+SparseMatrix<number> &
+SparseMatrix<number>::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<unsigned int> 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<int*>(&(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 <typename number>
 template <typename somenumber>
 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 (file)
index 0000000..7e305b7
--- /dev/null
@@ -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 <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_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<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);
+
+  // now copy things back into a SparseMatrix
+  SparseMatrix<double> 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 (file)
index 0000000..6a95ff5
--- /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,4) 8.00000
+(4,3) 9.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.