const TrilinosScalar *values,
const bool elide_zero_values = true,
const bool col_indices_are_sorted = false);
+
+ /**
+ * Set the element (<i>i,j</i>) to @p value.
+ *
+ * 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 (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 Tpetra::CrsMatrix data structure:
+ * If the same matrix entry is inserted more than once, the matrix entries
+ * will be added upon calling compress() (since Tpetra 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, const Number value);
+
+ /**
+ * Set all elements given in a FullMatrix<double> into the sparse matrix
+ * locations given by <tt>indices</tt>. In other words, this function
+ * writes the elements in <tt>full_matrix</tt> into the calling matrix,
+ * using the local-to-global indexing specified by <tt>indices</tt> for
+ * both the rows and the columns of the matrix. This function assumes a
+ * quadratic sparse matrix and a quadratic full_matrix, the usual
+ * situation in FE calculations.
+ *
+ * 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. 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 Tpetra::CrsMatrix 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<Number> &full_matrix,
+ const bool elide_zero_values = false);
+
+ /**
+ * Same function as before, but now including the possibility to use
+ * rectangular full_matrices and different local-to-global indexing on
+ * rows and columns, respectively.
+ */
+ void
+ set(const std::vector<size_type> &row_indices,
+ const std::vector<size_type> &col_indices,
+ const FullMatrix<Number> &full_matrix,
+ const bool elide_zero_values = false);
+
+ /**
+ * Set several elements in the specified row of the matrix with column
+ * indices as given by <tt>col_indices</tt> to the respective value.
+ *
+ * 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. 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 Tpetra::CrsMatrix 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,
+ const std::vector<Number> &values,
+ const bool elide_zero_values = false);
+
+ /**
+ * Set several elements to values given by <tt>values</tt> in a given row
+ * in columns given by col_indices into the sparse 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. 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 Tpetra::CrsMatrix 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,
+ const size_type *col_indices,
+ const Number *values,
+ const bool elide_zero_values = false);
+
/** @} */
/**
#ifdef DEAL_II_TRILINOS_WITH_TPETRA
# include <deal.II/lac/dynamic_sparsity_pattern.h>
+# include <deal.II/lac/full_matrix.h>
# include <deal.II/lac/trilinos_tpetra_sparse_matrix.h>
# include <deal.II/lac/trilinos_tpetra_sparsity_pattern.h>
+ template <typename Number, typename MemorySpace>
+ void
+ SparseMatrix<Number, MemorySpace>::set(const size_type row,
+ const size_type n_cols,
+ const size_type *col_indices,
+ const Number *values,
+ const bool elide_zero_values)
+ {
+ AssertIndexRange(row, this->m());
+ const types::signed_global_dof_index *col_index_ptr;
+ const Number *col_value_ptr;
+ const types::signed_global_dof_index trilinos_row = row;
+ types::signed_global_dof_index n_columns;
+
+ boost::container::small_vector<Number, 200> local_value_array(
+ elide_zero_values ? n_cols : 0);
+ boost::container::small_vector<types::signed_global_dof_index, 200>
+ local_index_array(elide_zero_values ? n_cols : 0);
+
+ // If we don't elide zeros, the pointers are already available... need to
+ // cast to non-const pointers as that is the format taken by Trilinos (but
+ // we will not modify const data)
+ if (elide_zero_values == false)
+ {
+ col_index_ptr =
+ reinterpret_cast<const dealii::types::signed_global_dof_index *>(
+ col_indices);
+ col_value_ptr = values;
+ n_columns = n_cols;
+ }
+ else
+ {
+ // Otherwise, extract nonzero values in each row and get the
+ // respective indices.
+ col_index_ptr = local_index_array.data();
+ col_value_ptr = local_value_array.data();
+
+ n_columns = 0;
+ for (size_type j = 0; j < n_cols; ++j)
+ {
+ const double value = values[j];
+ AssertIsFinite(value);
+ if (value != 0)
+ {
+ local_index_array[n_columns] = col_indices[j];
+ local_value_array[n_columns] = value;
+ ++n_columns;
+ }
+ }
+
+ AssertIndexRange(n_columns, n_cols + 1);
+ }
+
+ // We distinguish between two cases: the first one is when the matrix is
+ // not filled (i.e., it is possible to add new elements to the sparsity
+ // pattern), and the second one is when the pattern is already fixed. In
+ // the former case, we add the possibility to insert new values, and in
+ // the second we just replace data.
+ if (compressed || matrix->isFillComplete())
+ {
+ matrix->resumeFill();
+ compressed = false;
+ }
+
+ if (!matrix->isStaticGraph())
+ matrix->insertGlobalValues(trilinos_row,
+ n_columns,
+ col_value_ptr,
+ col_index_ptr);
+ else
+ matrix->replaceGlobalValues(trilinos_row,
+ n_columns,
+ col_value_ptr,
+ col_index_ptr);
+ }
+
+
+
+ template <typename Number, typename MemorySpace>
+ inline void
+ SparseMatrix<Number, MemorySpace>::set(const size_type i,
+ const size_type j,
+ const Number value)
+ {
+ AssertIsFinite(value);
+
+ set(i, 1, &j, &value, false);
+ }
+
+
+
+ template <typename Number, typename MemorySpace>
+ inline void
+ SparseMatrix<Number, MemorySpace>::set(
+ const std::vector<size_type> &indices,
+ const FullMatrix<Number> &values,
+ const bool elide_zero_values)
+ {
+ Assert(indices.size() == values.m(),
+ ExcDimensionMismatch(indices.size(), values.m()));
+ Assert(values.m() == values.n(), ExcNotQuadratic());
+
+ for (size_type i = 0; i < indices.size(); ++i)
+ set(indices[i],
+ indices.size(),
+ indices.data(),
+ &values(i, 0),
+ elide_zero_values);
+ }
+
+
+
template <typename Number, typename MemorySpace>
void
SparseMatrix<Number, MemorySpace>::add(
#include <deal.II/lac/sparse_matrix.h>
#include <deal.II/lac/trilinos_sparse_matrix.h>
#include <deal.II/lac/trilinos_vector.h>
+#ifdef DEAL_II_TRILINOS_WITH_TPETRA
+# include <deal.II/lac/trilinos_tpetra_sparse_matrix.h>
+# include <deal.II/lac/trilinos_tpetra_vector.h>
+#endif
#include <deal.II/lac/vector.h>
#include "../tests.h"
u = lo_A_plus_id_2 * f;
}
+#ifdef DEAL_II_TRILINOS_WITH_TPETRA
+ {
+ DynamicSparsityPattern csp(dim, dim);
+ testproblem.five_point_structure(csp);
+
+ LinearAlgebra::TpetraWrappers::SparseMatrix<double> A;
+ A.reinit(csp);
+ testproblem.five_point(A);
+
+ LinearAlgebra::TpetraWrappers::Vector<double> f;
+ f.reinit(complete_index_set(dim));
+ LinearAlgebra::TpetraWrappers::Vector<double> u;
+ u.reinit(complete_index_set(dim));
+
+ A.compress(VectorOperation::insert);
+ f.compress(VectorOperation::insert);
+ u.compress(VectorOperation::insert);
+
+ const auto lo_A =
+ linear_operator<LinearAlgebra::TpetraWrappers::Vector<double>>(A);
+ const auto lo_id_1 =
+ identity_operator<LinearAlgebra::TpetraWrappers::Vector<double>>(
+ lo_A.reinit_range_vector);
+ const auto lo_id_2 =
+ identity_operator<LinearAlgebra::TpetraWrappers::Vector<double>>(lo_A);
+ const auto lo_A_plus_id_1 = lo_A + lo_id_1; // Not a good idea. See below.
+ const auto lo_A_plus_id_2 = lo_A + lo_id_2;
+
+ // u = lo_A_plus_id_1 * f; // This is not allowed -- different payloads
+ u = lo_A_plus_id_2 * f;
+ }
+#endif
+
deallog << "OK" << std::endl;
return 0;
projects / dealii.git / commitdiff