*
* 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,
*
* 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,
*
* 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,
*
* 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,
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
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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();
+ }
+}