Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Public Types | Public Member Functions | Static Public Member Functions | Private Member Functions | List of all members

#include <deal.II/lac/trilinos_block_sparse_matrix.h>

Inheritance diagram for TrilinosWrappers::BlockSparseMatrix:
[legend]

Public Types

using BaseClass = BlockMatrixBase< SparseMatrix >
 
using BlockType = BaseClass::BlockType
 
using value_type = BaseClass::value_type
 
using pointer = BaseClass::pointer
 
using const_pointer = BaseClass::const_pointer
 
using reference = BaseClass::reference
 
using const_reference = BaseClass::const_reference
 
using size_type = BaseClass::size_type
 
using iterator = BaseClass::iterator
 
using const_iterator = BaseClass::const_iterator
 
- Public Types inherited from BlockMatrixBase< SparseMatrix >
using BlockType = SparseMatrix
 
using value_type = typename BlockType::value_type
 
using real_type = typename numbers::NumberTraits< value_type >::real_type
 
using pointer = value_type *
 
using const_pointer = const value_type *
 
using reference = value_type &
 
using const_reference = const value_type &
 
using size_type = types::global_dof_index
 
using iterator = MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, false > >
 
using const_iterator = MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, true > >
 

Public Member Functions

 BlockSparseMatrix ()=default
 
 ~BlockSparseMatrix () override
 
BlockSparseMatrixoperator= (const BlockSparseMatrix &)=default
 
BlockSparseMatrixoperator= (const double d)
 
void reinit (const size_type n_block_rows, const size_type n_block_columns)
 
template<typename BlockSparsityPatternType >
void reinit (const std::vector< IndexSet > &input_maps, const BlockSparsityPatternType &block_sparsity_pattern, const MPI_Comm &communicator=MPI_COMM_WORLD, const bool exchange_data=false)
 
template<typename BlockSparsityPatternType >
void reinit (const BlockSparsityPatternType &block_sparsity_pattern)
 
void reinit (const std::vector< IndexSet > &parallel_partitioning, const ::BlockSparseMatrix< double > &dealii_block_sparse_matrix, const MPI_Comm &communicator=MPI_COMM_WORLD, const double drop_tolerance=1e-13)
 
void reinit (const ::BlockSparseMatrix< double > &deal_ii_sparse_matrix, const double drop_tolerance=1e-13)
 
bool is_compressed () const
 
void collect_sizes ()
 
size_type n_nonzero_elements () const
 
MPI_Comm get_mpi_communicator () const
 
std::vector< IndexSetlocally_owned_domain_indices () const
 
std::vector< IndexSetlocally_owned_range_indices () const
 
template<typename VectorType1 , typename VectorType2 >
void vmult (VectorType1 &dst, const VectorType2 &src) const
 
template<typename VectorType1 , typename VectorType2 >
void Tvmult (VectorType1 &dst, const VectorType2 &src) const
 
TrilinosScalar residual (MPI::BlockVector &dst, const MPI::BlockVector &x, const MPI::BlockVector &b) const
 
TrilinosScalar residual (MPI::BlockVector &dst, const MPI::Vector &x, const MPI::BlockVector &b) const
 
TrilinosScalar residual (MPI::Vector &dst, const MPI::BlockVector &x, const MPI::Vector &b) const
 
TrilinosScalar residual (MPI::Vector &dst, const MPI::Vector &x, const MPI::Vector &b) const
 
template<>
void reinit (const BlockSparsityPattern &block_sparsity_pattern)
 
- Public Member Functions inherited from BlockMatrixBase< SparseMatrix >
 BlockMatrixBase ()=default
 
 ~BlockMatrixBase () override
 
BlockMatrixBasecopy_from (const BlockMatrixType &source)
 
BlockTypeblock (const unsigned int row, const unsigned int column)
 
const BlockTypeblock (const unsigned int row, const unsigned int column) const
 
size_type m () const
 
size_type n () const
 
unsigned int n_block_rows () const
 
unsigned int n_block_cols () const
 
void set (const size_type i, const size_type j, const value_type value)
 
void set (const std::vector< size_type > &indices, const FullMatrix< number > &full_matrix, const bool elide_zero_values=false)
 
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)
 
void set (const size_type row, const std::vector< size_type > &col_indices, const std::vector< number > &values, const bool elide_zero_values=false)
 
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)
 
void add (const size_type i, const size_type j, const value_type value)
 
void add (const std::vector< size_type > &indices, const FullMatrix< number > &full_matrix, const bool elide_zero_values=true)
 
void add (const std::vector< size_type > &row_indices, const std::vector< size_type > &col_indices, const FullMatrix< number > &full_matrix, const bool elide_zero_values=true)
 
void add (const size_type row, const std::vector< size_type > &col_indices, const std::vector< number > &values, const bool elide_zero_values=true)
 
void add (const size_type row, const size_type n_cols, const size_type *col_indices, const number *values, const bool elide_zero_values=true, const bool col_indices_are_sorted=false)
 
void add (const value_type factor, const BlockMatrixBase< SparseMatrix > &matrix)
 
value_type operator() (const size_type i, const size_type j) const
 
value_type el (const size_type i, const size_type j) const
 
value_type diag_element (const size_type i) const
 
void compress (::VectorOperation::values operation)
 
BlockMatrixBaseoperator*= (const value_type factor)
 
BlockMatrixBaseoperator/= (const value_type factor)
 
void vmult_add (BlockVectorType &dst, const BlockVectorType &src) const
 
void Tvmult_add (BlockVectorType &dst, const BlockVectorType &src) const
 
value_type matrix_norm_square (const BlockVectorType &v) const
 
real_type frobenius_norm () const
 
value_type matrix_scalar_product (const BlockVectorType &u, const BlockVectorType &v) const
 
value_type residual (BlockVectorType &dst, const BlockVectorType &x, const BlockVectorType &b) const
 
void print (std::ostream &out, const bool alternative_output=false) const
 
iterator begin ()
 
iterator begin (const size_type r)
 
const_iterator begin () const
 
const_iterator begin (const size_type r) const
 
iterator end ()
 
iterator end (const size_type r)
 
const_iterator end () const
 
const_iterator end (const size_type r) const
 
const BlockIndicesget_row_indices () const
 
const BlockIndicesget_column_indices () const
 
std::size_t memory_consumption () const
 

Static Public Member Functions

static ::ExceptionBaseExcIncompatibleRowNumbers (int arg1, int arg2, int arg3, int arg4)
 
static ::ExceptionBaseExcIncompatibleColNumbers (int arg1, int arg2, int arg3, int arg4)
 
- Static Public Member Functions inherited from BlockMatrixBase< SparseMatrix >
static ::ExceptionBaseExcIncompatibleRowNumbers (int arg1, int arg2, int arg3, int arg4)
 
static ::ExceptionBaseExcIncompatibleColNumbers (int arg1, int arg2, int arg3, int arg4)
 

Private Member Functions

template<typename VectorType1 , typename VectorType2 >
void vmult (VectorType1 &dst, const VectorType2 &src, const bool transpose, const std::integral_constant< bool, true >, const std::integral_constant< bool, true >) const
 
template<typename VectorType1 , typename VectorType2 >
void vmult (VectorType1 &dst, const VectorType2 &src, const bool transpose, const std::integral_constant< bool, false >, const std::integral_constant< bool, true >) const
 
template<typename VectorType1 , typename VectorType2 >
void vmult (VectorType1 &dst, const VectorType2 &src, const bool transpose, const std::integral_constant< bool, true >, const std::integral_constant< bool, false >) const
 
template<typename VectorType1 , typename VectorType2 >
void vmult (VectorType1 &dst, const VectorType2 &src, const bool transpose, const std::integral_constant< bool, false >, const std::integral_constant< bool, false >) const
 

Additional Inherited Members

- Protected Member Functions inherited from BlockMatrixBase< SparseMatrix >
void clear ()
 
void collect_sizes ()
 
void vmult_block_block (BlockVectorType &dst, const BlockVectorType &src) const
 
void vmult_block_nonblock (BlockVectorType &dst, const VectorType &src) const
 
void vmult_nonblock_block (VectorType &dst, const BlockVectorType &src) const
 
void vmult_nonblock_nonblock (VectorType &dst, const VectorType &src) const
 
void Tvmult_block_block (BlockVectorType &dst, const BlockVectorType &src) const
 
void Tvmult_block_nonblock (BlockVectorType &dst, const VectorType &src) const
 
void Tvmult_nonblock_block (VectorType &dst, const BlockVectorType &src) const
 
void Tvmult_nonblock_nonblock (VectorType &dst, const VectorType &src) const
 
void prepare_add_operation ()
 
void prepare_set_operation ()
 
- Protected Attributes inherited from BlockMatrixBase< SparseMatrix >
BlockIndices row_block_indices
 
BlockIndices column_block_indices
 
Table< 2, SmartPointer< BlockType, BlockMatrixBase< SparseMatrix > > > sub_objects
 

Detailed Description

Blocked sparse matrix based on the TrilinosWrappers::SparseMatrix class. This class implements the functions that are specific to the Trilinos SparseMatrix base objects for a blocked sparse matrix, and leaves the actual work relaying most of the calls to the individual blocks to the functions implemented in the base class. See there also for a description of when this class is useful.

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.

@ Block (linear algebra)

Author
Martin Kronbichler, Wolfgang Bangerth, 2008

Definition at line 72 of file trilinos_block_sparse_matrix.h.

Member Function Documentation

◆ reinit()

template<>
void TrilinosWrappers::BlockSparseMatrix::reinit ( const BlockSparsityPattern block_sparsity_pattern)

Definition at line 142 of file trilinos_block_sparse_matrix.cc.


The documentation for this class was generated from the following files: