Classes | |
class | ExcIncompatibleColNumbers |
class | ExcIncompatibleRowNumbers |
Public Types | |
typedef BlockMatrixBase < SparseMatrix > | BaseClass |
typedef BaseClass::BlockType | BlockType |
typedef BaseClass::value_type | value_type |
typedef BaseClass::pointer | pointer |
typedef BaseClass::const_pointer | const_pointer |
typedef BaseClass::reference | reference |
typedef BaseClass::const_reference | const_reference |
typedef BaseClass::size_type | size_type |
typedef BaseClass::iterator | iterator |
typedef BaseClass::const_iterator | const_iterator |
Public Member Functions | |
BlockSparseMatrix () | |
~BlockSparseMatrix () | |
BlockSparseMatrix & | operator= (const BlockSparseMatrix &) |
BlockSparseMatrix & | operator= (const double d) |
void | reinit (const unsigned int n_block_rows, const unsigned int n_block_columns) |
template<typename BlockSparsityType > | |
void | reinit (const std::vector< Epetra_Map > &input_maps, const BlockSparsityType &block_sparsity_pattern) |
template<typename BlockSparsityType > | |
void | reinit (const BlockSparsityType &block_sparsity_pattern) |
void | reinit (const std::vector< Epetra_Map > &input_maps, const ::::BlockSparseMatrix< double > &deal_ii_sparse_matrix, const double drop_tolerance=1e-13) |
void | reinit (const ::::BlockSparseMatrix< double > &deal_ii_sparse_matrix, const double drop_tolerance=1e-13) |
void | compress () |
bool | is_compressed () const |
void | collect_sizes () |
unsigned int | n_nonzero_elements () const |
void | vmult (MPI::BlockVector &dst, const MPI::BlockVector &src) const |
void | vmult (BlockVector &dst, const BlockVector &src) const |
void | vmult (MPI::BlockVector &dst, const MPI::Vector &src) const |
void | vmult (BlockVector &dst, const Vector &src) const |
void | vmult (MPI::Vector &dst, const MPI::BlockVector &src) const |
void | vmult (Vector &dst, const BlockVector &src) const |
void | vmult (VectorBase &dst, const VectorBase &src) const |
void | Tvmult (MPI::BlockVector &dst, const MPI::BlockVector &src) const |
void | Tvmult (BlockVector &dst, const BlockVector &src) const |
void | Tvmult (MPI::BlockVector &dst, const MPI::Vector &src) const |
void | Tvmult (BlockVector &dst, const Vector &src) const |
void | Tvmult (MPI::Vector &dst, const MPI::BlockVector &src) const |
void | Tvmult (Vector &dst, const BlockVector &src) const |
void | Tvmult (VectorBase &dst, const VectorBase &src) const |
TrilinosScalar | residual (MPI::BlockVector &dst, const MPI::BlockVector &x, const MPI::BlockVector &b) const |
TrilinosScalar | residual (BlockVector &dst, const BlockVector &x, const BlockVector &b) const |
TrilinosScalar | residual (MPI::BlockVector &dst, const MPI::Vector &x, const MPI::BlockVector &b) const |
TrilinosScalar | residual (BlockVector &dst, const Vector &x, const BlockVector &b) const |
TrilinosScalar | residual (MPI::Vector &dst, const MPI::BlockVector &x, const MPI::Vector &b) const |
TrilinosScalar | residual (Vector &dst, const BlockVector &x, const Vector &b) const |
TrilinosScalar | residual (VectorBase &dst, const VectorBase &x, const VectorBase &b) const |
In contrast to the deal.II-type SparseMatrix class, the Trilinos matrices do not have external objects for the sparsity patterns. Thus, one does not determine the size of the individual blocks of a block matrix of this type by attaching a block sparsity pattern, but by calling reinit() to set the number of blocks and then by setting the size of each block separately. In order to fix the data structures of the block matrix, it is then necessary to let it know that we have changed the sizes of the underlying matrices. For this, one has to call the collect_sizes() function, for much the same reason as is documented with the BlockSparsityPattern class.
Typedef the base class for simpler access to its own typedefs.
Typedef the type of the underlying matrix.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Import the typedefs from the base class.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Reimplemented from BlockMatrixBase< SparseMatrix >.
TrilinosWrappers::BlockSparseMatrix::BlockSparseMatrix | ( | ) |
Constructor; initializes the matrix to be empty, without any structure, i.e. the matrix is not usable at all. This constructor is therefore only useful for matrices which are members of a class. All other matrices should be created at a point in the data flow where all necessary information is available.
You have to initialize the matrix before usage with reinit(BlockSparsityPattern). The number of blocks per row and column are then determined by that function.
TrilinosWrappers::BlockSparseMatrix::~BlockSparseMatrix | ( | ) |
Destructor.
BlockSparseMatrix& TrilinosWrappers::BlockSparseMatrix::operator= | ( | const BlockSparseMatrix & | ) |
Pseudo copy operator only copying empty objects. The sizes of the block matrices need to be the same.
BlockSparseMatrix& TrilinosWrappers::BlockSparseMatrix::operator= | ( | const double | d | ) |
This operator assigns a scalar to a matrix. Since this does usually not make much sense (should we set all matrix entries to this value? Only the nonzero entries of the sparsity pattern?), this operation is only allowed if the actual value to be assigned is zero. This operator only exists to allow for the obvious notation matrix=0
, which sets all elements of the matrix to zero, but keep the sparsity pattern previously used.
void TrilinosWrappers::BlockSparseMatrix::reinit | ( | const unsigned int | n_block_rows, | |
const unsigned int | n_block_columns | |||
) |
Resize the matrix, by setting the number of block rows and columns. This deletes all blocks and replaces them by unitialized ones, i.e. ones for which also the sizes are not yet set. You have to do that by calling the reinit
functions of the blocks themselves. Do not forget to call collect_sizes() after that on this object.
The reason that you have to set sizes of the blocks yourself is that the sizes may be varying, the maximum number of elements per row may be varying, etc. It is simpler not to reproduce the interface of the SparsityPattern
class here but rather let the user call whatever function she desires.
void TrilinosWrappers::BlockSparseMatrix::reinit | ( | const std::vector< Epetra_Map > & | input_maps, | |
const BlockSparsityType & | block_sparsity_pattern | |||
) | [inline] |
Resize the matrix, by using an array of Epetra maps to determine the distribution of the individual matrices. This function assumes that a quadratic block matrix is generated.
void TrilinosWrappers::BlockSparseMatrix::reinit | ( | const BlockSparsityType & | block_sparsity_pattern | ) | [inline] |
Resize the matrix and initialize it by the given sparsity pattern. Since no distribution map is given, the result is a block matrix for which all elements are stored locally.
void TrilinosWrappers::BlockSparseMatrix::reinit | ( | const std::vector< Epetra_Map > & | input_maps, | |
const ::::BlockSparseMatrix< double > & | deal_ii_sparse_matrix, | |||
const double | drop_tolerance = 1e-13 | |||
) |
This function initializes the Trilinos matrix using the deal.II sparse matrix and the entries stored therein. It uses a threshold to copy only elements whose modulus is larger than the threshold (so zeros in the deal.II matrix can be filtered away).
void TrilinosWrappers::BlockSparseMatrix::reinit | ( | const ::::BlockSparseMatrix< double > & | deal_ii_sparse_matrix, | |
const double | drop_tolerance = 1e-13 | |||
) |
This function initializes the Trilinos matrix using the deal.II sparse matrix and the entries stored therein. It uses a threshold to copy only elements whose modulus is larger than the threshold (so zeros in the deal.II matrix can be filtered away). Since no Epetra_Map is given, all the elements will be locally stored.
void TrilinosWrappers::BlockSparseMatrix::compress | ( | ) |
This function calls the compress() command of all matrices after the assembly is completed. Note that all MPI processes need to call this command (whereas the individual assembly routines will most probably only be called on each processor individually) before any can complete it.
Reimplemented from BlockMatrixBase< SparseMatrix >.
bool BlockSparseMatrix< number >::is_compressed | ( | ) | const [inline] |
Returns the state of the matrix, i.e., whether compress() needs to be called after an operation requiring data exchange. Does only return non-true values when used in debug
mode, since it is quite expensive to keep track of all operations that lead to the need for compress().
References BlockMatrixBase< SparseMatrix >::block(), BlockMatrixBase< SparseMatrix >::n_block_cols(), and BlockMatrixBase< SparseMatrix >::n_block_rows().
void TrilinosWrappers::BlockSparseMatrix::collect_sizes | ( | ) |
This function collects the sizes of the sub-objects and stores them in internal arrays, in order to be able to relay global indices into the matrix to indices into the subobjects. You *must* call this function each time after you have changed the size of the sub-objects. Note that this is a collective operation, i.e., it needs to be called on all MPI processes. This command internally calls the method compress()
, so you don't need to call that function in case you use collect_sizes()
.
Reimplemented from BlockMatrixBase< SparseMatrix >.
Return the number of nonzero elements of this matrix.
void BlockSparseMatrix< number >::vmult | ( | MPI::BlockVector & | dst, | |
const MPI::BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication: let with
being this matrix.
References BlockMatrixBase< SparseMatrix >::vmult_block_block().
void BlockSparseMatrix< number >::vmult | ( | BlockVector & | dst, | |
const BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication: let with
being this matrix, now applied to localized block vectors (works only when run on one processor).
References BlockMatrixBase< SparseMatrix >::vmult_block_block().
void BlockSparseMatrix< number >::vmult | ( | MPI::BlockVector & | dst, | |
const MPI::Vector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block column.
References BlockMatrixBase< SparseMatrix >::vmult_block_nonblock().
void BlockSparseMatrix< number >::vmult | ( | BlockVector & | dst, | |
const Vector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block column, now applied to localized vectors (works only when run on one processor).
References BlockMatrixBase< SparseMatrix >::vmult_block_nonblock().
void BlockSparseMatrix< number >::vmult | ( | MPI::Vector & | dst, | |
const MPI::BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block row.
References BlockMatrixBase< SparseMatrix >::vmult_nonblock_block().
void BlockSparseMatrix< number >::vmult | ( | Vector & | dst, | |
const BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block row, now applied to localized vectors (works only when run on one processor).
References BlockMatrixBase< SparseMatrix >::vmult_nonblock_block().
void BlockSparseMatrix< number >::vmult | ( | VectorBase & | dst, | |
const VectorBase & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block.
References BlockMatrixBase< SparseMatrix >::vmult_nonblock_nonblock().
void BlockSparseMatrix< number >::Tvmult | ( | MPI::BlockVector & | dst, | |
const MPI::BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication: let with
being this matrix. This function does the same as vmult() but takes the transposed matrix.
References BlockMatrixBase< SparseMatrix >::Tvmult_block_block().
void BlockSparseMatrix< number >::Tvmult | ( | BlockVector & | dst, | |
const BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication: let with
being this matrix. This function does the same as vmult() but takes the transposed matrix, now applied to localized Trilinos vectors (works only when run on one processor).
References BlockMatrixBase< SparseMatrix >::Tvmult_block_block().
void BlockSparseMatrix< number >::Tvmult | ( | MPI::BlockVector & | dst, | |
const MPI::Vector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block row.
References BlockMatrixBase< SparseMatrix >::Tvmult_block_nonblock().
void BlockSparseMatrix< number >::Tvmult | ( | BlockVector & | dst, | |
const Vector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block row, now applied to localized Trilinos vectors (works only when run on one processor).
References BlockMatrixBase< SparseMatrix >::Tvmult_block_nonblock().
void BlockSparseMatrix< number >::Tvmult | ( | MPI::Vector & | dst, | |
const MPI::BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block column.
References BlockMatrixBase< SparseMatrix >::Tvmult_nonblock_block().
void BlockSparseMatrix< number >::Tvmult | ( | Vector & | dst, | |
const BlockVector & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block column, now applied to localized Trilinos vectors (works only when run on one processor).
References BlockMatrixBase< SparseMatrix >::Tvmult_nonblock_block().
void BlockSparseMatrix< number >::Tvmult | ( | VectorBase & | dst, | |
const VectorBase & | src | |||
) | const [inline] |
Matrix-vector multiplication. Just like the previous function, but only applicable if the matrix has only one block.
References BlockMatrixBase< SparseMatrix >::Tvmult_nonblock_nonblock().
TrilinosScalar TrilinosWrappers::BlockSparseMatrix::residual | ( | MPI::BlockVector & | dst, | |
const MPI::BlockVector & | x, | |||
const MPI::BlockVector & | b | |||
) | const |
Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst
. The l2 norm of the residual vector is returned.
Source x and destination dst must not be the same vector.
Note that both vectors have to be distributed vectors generated using the same Map as was used for the matrix in case you work on a distributed memory architecture, using the interface in the TrilinosWrappers::MPI::BlockVector class.
TrilinosScalar TrilinosWrappers::BlockSparseMatrix::residual | ( | BlockVector & | dst, | |
const BlockVector & | x, | |||
const BlockVector & | b | |||
) | const |
Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst
. The l2 norm of the residual vector is returned.
Source x and destination dst must not be the same vector.
Note that both vectors have to be distributed vectors generated using the same Map as was used for the matrix in case you work on a distributed memory architecture, using the interface in the TrilinosWrappers::BlockVector class. Since the block matrix is in general distributed among processes, this function only works when running the program on one processor.
TrilinosScalar TrilinosWrappers::BlockSparseMatrix::residual | ( | MPI::BlockVector & | dst, | |
const MPI::Vector & | x, | |||
const MPI::BlockVector & | b | |||
) | const |
Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst
. The l2 norm of the residual vector is returned. Just like the previous function, but only applicable if the matrix only has one block row.
TrilinosScalar TrilinosWrappers::BlockSparseMatrix::residual | ( | BlockVector & | dst, | |
const Vector & | x, | |||
const BlockVector & | b | |||
) | const |
Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst
. The l2 norm of the residual vector is returned. Just like the previous function, but only applicable if the matrix only has one block row.
TrilinosScalar TrilinosWrappers::BlockSparseMatrix::residual | ( | MPI::Vector & | dst, | |
const MPI::BlockVector & | x, | |||
const MPI::Vector & | b | |||
) | const |
Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst
. The l2 norm of the residual vector is returned. Just like the previous function, but only applicable if the matrix only has one block column.
TrilinosScalar TrilinosWrappers::BlockSparseMatrix::residual | ( | Vector & | dst, | |
const BlockVector & | x, | |||
const Vector & | b | |||
) | const |
Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst
. The l2 norm of the residual vector is returned. Just like the previous function, but only applicable if the matrix only has one block column.
TrilinosScalar TrilinosWrappers::BlockSparseMatrix::residual | ( | VectorBase & | dst, | |
const VectorBase & | x, | |||
const VectorBase & | b | |||
) | const |
Compute the residual of an equation Mx=b, where the residual is defined to be r=b-Mx. Write the residual into dst
. The l2 norm of the residual vector is returned. Just like the previous function, but only applicable if the matrix only has one block.