]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Improve documentation of TrilinosWrappers::SparseMatrix::set
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 5 Aug 2015 06:52:51 +0000 (08:52 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Wed, 5 Aug 2015 08:45:24 +0000 (10:45 +0200)
include/deal.II/lac/trilinos_sparse_matrix.h
tests/trilinos/sparse_matrix_set_03.cc [new file with mode: 0644]
tests/trilinos/sparse_matrix_set_03.mpirun=1.output [new file with mode: 0644]
tests/trilinos/sparse_matrix_set_03.mpirun=2.output [new file with mode: 0644]
tests/trilinos/trilinos_sparsity_pattern_04.cc

index 6abe7d3cfb9cd2ca9588a0c544d671039f0187e4..bfc51149a1eb2942e4181133dffb8c73281540a7 100644 (file)
@@ -1112,15 +1112,21 @@ namespace TrilinosWrappers
      *
      * This function is able to insert new elements into the matrix as long as
      * compress() has not been called, so the sparsity pattern will be
-     * extended. When compress() is called for the first time, then this is no
-     * longer possible and an insertion of elements at positions which have
-     * not been initialized will throw an exception. Note that in case
-     * elements need to be inserted, it is mandatory that elements are
-     * inserted only once. Otherwise, the elements will actually be added in
-     * the end (since it is not possible to efficiently find values to the
-     * same entry before compress() has been called). In the case that an
-     * element is set more than once, initialize the matrix with a sparsity
-     * pattern first.
+     * extended. When compress() is called for the first time (or in case the
+     * matrix is initialized from a sparsity pattern), no new elements can be
+     * added and an insertion of elements at positions which have not been
+     * initialized will throw an exception.
+     *
+     * For the case that the matrix is constructed without a sparsity pattern
+     * and new matrix entries are added on demand, please note the following
+     * behavior imposed by the underlying Epetra_FECrsMatrix data
+     * structure: If the same matrix entry is inserted more than once, the
+     * matrix entries will be added upon calling compress() (since Epetra does
+     * not track values to the same entry before the final compress() is
+     * called), even if VectorOperation::insert is specified as argument to
+     * compress(). In the case you cannot make sure that matrix entries are
+     * only set once, initialize the matrix with a sparsity pattern to fix the
+     * matrix structure before inserting elements.
      */
     void set (const size_type i,
               const size_type j,
@@ -1137,14 +1143,26 @@ namespace TrilinosWrappers
      *
      * This function is able to insert new elements into the matrix as long as
      * compress() has not been called, so the sparsity pattern will be
-     * extended. When compress() is called for the first time, then this is no
-     * longer possible and an insertion of elements at positions which have
-     * not been initialized will throw an exception.
+     * extended. After compress() has been called for the first time or the
+     * matrix has been initialized from a sparsity pattern, extending the
+     * sparsity pattern is no longer possible and an insertion of elements at
+     * positions which have not been initialized will throw an exception.
      *
      * The optional parameter <tt>elide_zero_values</tt> can be used to
      * specify whether zero values should be inserted anyway or they should be
      * filtered away. The default value is <tt>false</tt>, i.e., even zero
      * values are inserted/replaced.
+     *
+     * For the case that the matrix is constructed without a sparsity pattern
+     * and new matrix entries are added on demand, please note the following
+     * behavior imposed by the underlying Epetra_FECrsMatrix data
+     * structure: If the same matrix entry is inserted more than once, the
+     * matrix entries will be added upon calling compress() (since Epetra does
+     * not track values to the same entry before the final compress() is
+     * called), even if VectorOperation::insert is specified as argument to
+     * compress(). In the case you cannot make sure that matrix entries are
+     * only set once, initialize the matrix with a sparsity pattern to fix the
+     * matrix structure before inserting elements.
      */
     void set (const std::vector<size_type>     &indices,
               const FullMatrix<TrilinosScalar> &full_matrix,
@@ -1166,14 +1184,26 @@ namespace TrilinosWrappers
      *
      * This function is able to insert new elements into the matrix as long as
      * compress() has not been called, so the sparsity pattern will be
-     * extended. When compress() is called for the first time, then this is no
-     * longer possible and an insertion of elements at positions which have
-     * not been initialized will throw an exception.
+     * extended. After compress() has been called for the first time or the
+     * matrix has been initialized from a sparsity pattern, extending the
+     * sparsity pattern is no longer possible and an insertion of elements at
+     * positions which have not been initialized will throw an exception.
      *
      * The optional parameter <tt>elide_zero_values</tt> can be used to
      * specify whether zero values should be inserted anyway or they should be
      * filtered away. The default value is <tt>false</tt>, i.e., even zero
      * values are inserted/replaced.
+     *
+     * For the case that the matrix is constructed without a sparsity pattern
+     * and new matrix entries are added on demand, please note the following
+     * behavior imposed by the underlying Epetra_FECrsMatrix data
+     * structure: If the same matrix entry is inserted more than once, the
+     * matrix entries will be added upon calling compress() (since Epetra does
+     * not track values to the same entry before the final compress() is
+     * called), even if VectorOperation::insert is specified as argument to
+     * compress(). In the case you cannot make sure that matrix entries are
+     * only set once, initialize the matrix with a sparsity pattern to fix the
+     * matrix structure before inserting elements.
      */
     void set (const size_type                    row,
               const std::vector<size_type>      &col_indices,
@@ -1186,14 +1216,26 @@ namespace TrilinosWrappers
      *
      * This function is able to insert new elements into the matrix as long as
      * compress() has not been called, so the sparsity pattern will be
-     * extended. When compress() is called for the first time, then this is no
-     * longer possible and an insertion of elements at positions which have
-     * not been initialized will throw an exception.
+     * extended. After compress() has been called for the first time or the
+     * matrix has been initialized from a sparsity pattern, extending the
+     * sparsity pattern is no longer possible and an insertion of elements at
+     * positions which have not been initialized will throw an exception.
      *
      * The optional parameter <tt>elide_zero_values</tt> can be used to
      * specify whether zero values should be inserted anyway or they should be
      * filtered away. The default value is <tt>false</tt>, i.e., even zero
      * values are inserted/replaced.
+     *
+     * For the case that the matrix is constructed without a sparsity pattern
+     * and new matrix entries are added on demand, please note the following
+     * behavior imposed by the underlying Epetra_FECrsMatrix data
+     * structure: If the same matrix entry is inserted more than once, the
+     * matrix entries will be added upon calling compress() (since Epetra does
+     * not track values to the same entry before the final compress() is
+     * called), even if VectorOperation::insert is specified as argument to
+     * compress(). In the case you cannot make sure that matrix entries are
+     * only set once, initialize the matrix with a sparsity pattern to fix the
+     * matrix structure before inserting elements.
      */
     void set (const size_type       row,
               const size_type       n_cols,
@@ -2509,10 +2551,11 @@ namespace TrilinosWrappers
   SparseMatrix::size_type
   SparseMatrix::n () const
   {
+    Assert(column_space_map.get() != 0, ExcInternalError());
 #ifndef DEAL_II_WITH_64BIT_INDICES
-    return matrix->NumGlobalCols();
+    return column_space_map->NumGlobalElements();
 #else
-    return matrix->NumGlobalCols64();
+    return column_space_map->NumGlobalElements64();
 #endif
   }
 
diff --git a/tests/trilinos/sparse_matrix_set_03.cc b/tests/trilinos/sparse_matrix_set_03.cc
new file mode 100644 (file)
index 0000000..0df5d18
--- /dev/null
@@ -0,0 +1,154 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2004 - 2015 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+// check setting off-processor entries of Epetra_CrsMatrix. It turns out that
+// the underlying Epetra data structures actually add the entries even though
+// we want to insert them only.
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+#include <deal.II/lac/trilinos_sparse_matrix.h>
+#include <fstream>
+#include <iostream>
+
+
+void test ()
+{
+  const unsigned int n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
+  const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+
+  IndexSet rows(5*n_procs);
+  rows.add_range(5*myid, 5*(myid+1));
+  rows.compress();
+  IndexSet columns(5);
+  columns.add_range(5*myid/n_procs, 5*(myid+1)/n_procs);
+  columns.compress();
+
+  {
+    TrilinosWrappers::SparseMatrix m(rows, columns, MPI_COMM_WORLD);
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j*.5+.5);
+
+    m.compress (VectorOperation::insert);
+    deallog << "Matrix norm set local once: " << m.frobenius_norm() << std::endl;
+  }
+
+  {
+    TrilinosWrappers::SparseMatrix m(rows, columns, MPI_COMM_WORLD);
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j*.5+.5);
+
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j+1.);
+
+    m.compress (VectorOperation::insert);
+    deallog << "Matrix norm set local twice: " << m.frobenius_norm() << std::endl;
+  }
+
+  {
+    TrilinosWrappers::SparseMatrix m(rows, columns, MPI_COMM_WORLD);
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j*.5+.5);
+
+    if (myid == 1)
+      m.set(1, 3, 10.);
+
+    m.compress (VectorOperation::insert);
+    deallog << "Matrix norm set non-local once: " << m.frobenius_norm() << std::endl;
+  }
+
+  {
+    TrilinosWrappers::SparseMatrix m(rows, columns, MPI_COMM_WORLD);
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j*.5+.5);
+
+    m.set(1, 3, 10.);
+
+    m.compress (VectorOperation::insert);
+    deallog << "Matrix norm set non-local twice: " << m.frobenius_norm() << std::endl;
+  }
+
+  {
+    TrilinosWrappers::SparseMatrix m(rows, columns, MPI_COMM_WORLD);
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j*.5+.5);
+
+    m.set(1, 3, 10.);
+    m.set(2, 3, 2*3*0.5+0.5);
+
+    m.compress (VectorOperation::insert);
+    deallog << "Matrix norm set non-local twice: " << m.frobenius_norm() << std::endl;
+
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j*.5+.5);
+
+    m.set(1, 3, 10.);
+    m.set(2, 3, 2*3*0.5+0.5);
+
+    m.compress (VectorOperation::insert);
+    deallog << "Matrix norm set twice, 2nd attempt: " << m.frobenius_norm() << std::endl;
+
+    for (unsigned int i=5*myid; i<5*(myid+1); ++i)
+      for (unsigned int j=0; j<m.n(); ++j)
+        if ((i+2*j+1) % 3 == 0)
+          m.set (i,j, i*j*.5+.5);
+
+    m.set(1, 3, 10.);
+    m.set(2, 3, 2*3*0.5+0.5);
+    m.set(2, 3, 2*3*0.5+0.5);
+
+    m.compress (VectorOperation::insert);
+    deallog << "Matrix norm set twice-twice: " << m.frobenius_norm() << std::endl;
+  }
+
+  deallog << "OK" << std::endl;
+}
+
+
+
+int main (int argc,char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
+
+  if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+    {
+      std::ofstream logfile ("output");
+      deallog.attach (logfile);
+      deallog.depth_console (0);
+      deallog.threshold_double (1.e-10);
+
+      test();
+    }
+  else
+    {
+      test();
+    }
+}
diff --git a/tests/trilinos/sparse_matrix_set_03.mpirun=1.output b/tests/trilinos/sparse_matrix_set_03.mpirun=1.output
new file mode 100644 (file)
index 0000000..b28c15d
--- /dev/null
@@ -0,0 +1,9 @@
+
+DEAL::Matrix norm set local once: 9.04157
+DEAL::Matrix norm set local twice: 27.1247
+DEAL::Matrix norm set non-local once: 9.04157
+DEAL::Matrix norm set non-local twice: 13.4815
+DEAL::Matrix norm set non-local twice: 14.7817
+DEAL::Matrix norm set twice, 2nd attempt: 13.4815
+DEAL::Matrix norm set twice-twice: 13.4815
+DEAL::OK
diff --git a/tests/trilinos/sparse_matrix_set_03.mpirun=2.output b/tests/trilinos/sparse_matrix_set_03.mpirun=2.output
new file mode 100644 (file)
index 0000000..f4a36ad
--- /dev/null
@@ -0,0 +1,9 @@
+
+DEAL::Matrix norm set local once: 29.9082
+DEAL::Matrix norm set local twice: 89.7246
+DEAL::Matrix norm set non-local once: 31.5357
+DEAL::Matrix norm set non-local twice: 35.9792
+DEAL::Matrix norm set non-local twice: 37.3162
+DEAL::Matrix norm set twice, 2nd attempt: 31.5357
+DEAL::Matrix norm set twice-twice: 31.5357
+DEAL::OK
index e184f49af737f453dcd252ad4e1372492a9f7b1f..89d2fc6ff2d5a334c0c9dff1f41c60ef3ccb7130 100644 (file)
@@ -196,11 +196,6 @@ void test ()
 
 int main (int argc, char **argv)
 {
-  std::ofstream logfile("output");
-  deallog.attach(logfile);
-  deallog.depth_console(0);
-  deallog.threshold_double(1.e-10);
-
   Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, testing_max_num_threads());
 
   if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)

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.