Ifpack2 Templated Preconditioning Package  Version 1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Member Functions | Static Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | List of all members
Ifpack2::Container< MatrixType > Class Template Referenceabstract

Interface for creating and solving a set of local linear problems. More...

#include <Ifpack2_Container_decl.hpp>

Inheritance diagram for Ifpack2::Container< MatrixType >:
Inheritance graph
[legend]

Public Member Functions

 Container (const Teuchos::RCP< const row_matrix_type > &matrix, const Teuchos::Array< Teuchos::Array< LO > > &partitions, bool pointIndexed)
 Constructor. More...
 
virtual ~Container ()
 Destructor. More...
 
Teuchos::ArrayView< const LO > getBlockRows (int blockIndex) const
 Local indices of the rows of the input matrix that belong to this block. More...
 
virtual void initialize ()=0
 Do all set-up operations that only require matrix structure. More...
 
void setBlockSizes (const Teuchos::Array< Teuchos::Array< LO > > &partitions)
 Initialize arrays with information about block sizes. More...
 
bool isInitialized () const
 Whether the container has been successfully initialized. More...
 
bool isComputed () const
 Whether the container has been successfully computed. More...
 
virtual void compute ()=0
 Extract the local diagonal blocks and prepare the solver. More...
 
virtual void setParameters (const Teuchos::ParameterList &List)=0
 Set parameters, if any. More...
 
virtual void apply (HostView X, HostView Y, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const =0
 Compute Y := alpha * M^{-1} X + beta*Y. More...
 
virtual void weightedApply (HostView X, HostView Y, HostView D, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const =0
 Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y. More...
 
virtual void applyInverseJacobi (const mv_type &, mv_type &, SC dampingFactor, bool, int) const =0
 Compute Y := (1 - a) Y + a D^{-1} (X - R*Y). More...
 
virtual void applyMV (mv_type &X, mv_type &Y) const
 Wrapper for apply with MultiVector. More...
 
virtual void weightedApplyMV (mv_type &X, mv_type &Y, vector_type &W) const
 Wrapper for weightedApply with MultiVector. More...
 
virtual std::ostream & print (std::ostream &os) const =0
 Print basic information about the container to os. More...
 

Static Public Member Functions

static std::string getName ()
 

Protected Types

using ISC = typename Kokkos::Details::ArithTraits< SC >::val_type
 Internal representation of Scalar in Kokkos::View. More...
 
using HostView = typename mv_type::dual_view_type::t_host
 

Protected Member Functions

virtual void DoGSBlock (HostView X, HostView Y, HostView Y2, HostView Resid, SC dampingFactor, LO i) const
 Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS) More...
 

Protected Attributes

Teuchos::RCP< const
row_matrix_type > 
inputMatrix_
 The input matrix to the constructor. More...
 
Teuchos::RCP< const
crs_matrix_type > 
inputCrsMatrix_
 The input matrix, dynamic cast to CrsMatrix. May be null. More...
 
Teuchos::RCP< const
block_crs_matrix_type > 
inputBlockMatrix_
 The input matrix, dynamic cast to BlockCrsMatrix. May be null. More...
 
int numBlocks_
 The number of blocks (partitions) in the container. More...
 
Teuchos::Array< LO > blockRows_
 Local indices of the rows of the input matrix that belong to this block. More...
 
Teuchos::Array< LO > blockSizes_
 Number of rows in each block. More...
 
Teuchos::Array< LO > blockOffsets_
 Starting index in blockRows_ of local row indices for each block. More...
 
Teuchos::RCP< vector_type > Diag_
 Diagonal elements. More...
 
bool IsParallel_
 Whether the problem is distributed across multiple MPI processes. More...
 
LO NumLocalRows_
 Number of local rows in input matrix. More...
 
GO NumGlobalRows_
 Number of global rows in input matrix. More...
 
GO NumGlobalNonzeros_
 Number of nonzeros in input matrix. More...
 
bool hasBlockCrs_
 Whether the input matrix is a BlockCRS matrix. More...
 
int bcrsBlockSize_
 If hasBlockCrs_, the number of DOFs per vertex. Otherwise 1. More...
 
bool pointIndexed_
 (If hasBlockCrs_) Whether the blocks are described using sub-block row indices instead of full block rows. More...
 
LO scalarsPerRow_
 
bool IsInitialized_
 If true, the container has been successfully initialized. More...
 
bool IsComputed_
 If true, the container has been successfully computed. More...
 

Detailed Description

template<class MatrixType>
class Ifpack2::Container< MatrixType >

Interface for creating and solving a set of local linear problems.

Template Parameters
MatrixTypeA specialization of Tpetra::RowMatrix.

This class is mainly useful for the implementation of BlockRelaxation, and other preconditioners that need to solve linear systems with diagonal blocks of a sparse matrix.

Users of BlockRelaxation (and any analogous preconditioners) do not normally need to interact with the Container interface. However, they do need to specify a specific Container subclass to use, for example as the second template parameter (ContainerType) of BlockRelaxation. Implementations of Container specify

For example, the SparseContainer subclass uses sparse matrices (Tpetra::CrsMatrix) to store each diagonal block, and can use any given Ifpack2 Preconditioner subclass to solve linear systems.

A Container can create, populate, and solve local linear systems. A local linear system matrix, B, is a submatrix of the local components of the distributed matrix, A. The idea of Container is to specify the rows of A that are contained in each B, then extract B from A, and compute all it is necessary to solve a linear system in B. Then, set the initial guess (if necessary) and right-hand sides for B, and solve the linear systems in B.

If you are writing a class (comparable to BlockRelaxation) that uses Container, you should use it in the following way:

  1. Create a Container object, specifying the global matrix A and the indices of the local rows of A that are contained in each B. The latter indices come from a Partitioner object.
  2. Optionally, set linear solve parameters using setParameters().
  3. Initialize the container by calling initialize().
  4. Prepare the linear system solver using compute().
  5. Solve a linear system using apply() for each block index (from 0 to numBlocks_) as needed.

For an example of Steps 1-5 above, see ExtractSubmatrices() and ApplyInverseGS() in Ifpack2_BlockRelaxation_def.hpp.

Member Typedef Documentation

template<class MatrixType>
using Ifpack2::Container< MatrixType >::ISC = typename Kokkos::Details::ArithTraits<SC>::val_type
protected

Internal representation of Scalar in Kokkos::View.

template<class MatrixType>
using Ifpack2::Container< MatrixType >::HostView = typename mv_type::dual_view_type::t_host
protected

HostView (the host-space internal representation for Tpetra::Multivector) is the type of the vector arguments of DoJacobi, DoGaussSeidel, and DoSGS.

Constructor & Destructor Documentation

template<class MatrixType >
Ifpack2::Container< MatrixType >::Container ( const Teuchos::RCP< const row_matrix_type > &  matrix,
const Teuchos::Array< Teuchos::Array< LO > > &  partitions,
bool  pointIndexed 
)

Constructor.

Parameters
matrix[in] The original input matrix. This Container will construct local diagonal blocks from its rows according to partitions.
partitioner[in] The Partitioner object that assigns local rows of the input matrix to blocks.
pointIndexed[in] If the input matrix is a Tpetra::BlockCrsMatrix, whether elements of partitions[k] identify rows within blocks (true) or whole blocks (false).
template<class MatrixType >
Ifpack2::Container< MatrixType >::~Container ( )
virtual

Destructor.

Member Function Documentation

template<class MatrixType >
Teuchos::ArrayView< const typename MatrixType::local_ordinal_type > Ifpack2::Container< MatrixType >::getBlockRows ( int  blockIndex) const

Local indices of the rows of the input matrix that belong to this block.

The set of (local) rows assigned to this Container is defined by passing in a set of indices blockRows[i] = j to the constructor, where i (from 0 to getNumRows() - 1) indicates the Container's row, and j indicates the local row in the calling process. Subclasses must always pass along these indices to the base class.

The indices are usually used to reorder the local row index (on the calling process) of the i-th row in the Container.

For an example of how to use these indices, see the implementation of BlockRelaxation::ExtractSubmatrices() in Ifpack2_BlockRelaxation_def.hpp.

template<class MatrixType>
virtual void Ifpack2::Container< MatrixType >::initialize ( )
pure virtual
template<class MatrixType >
void Ifpack2::Container< MatrixType >::setBlockSizes ( const Teuchos::Array< Teuchos::Array< LO > > &  partitions)

Initialize arrays with information about block sizes.

template<class MatrixType >
bool Ifpack2::Container< MatrixType >::isInitialized ( ) const

Whether the container has been successfully initialized.

template<class MatrixType >
bool Ifpack2::Container< MatrixType >::isComputed ( ) const

Whether the container has been successfully computed.

template<class MatrixType>
virtual void Ifpack2::Container< MatrixType >::compute ( )
pure virtual
template<class MatrixType>
virtual void Ifpack2::Container< MatrixType >::setParameters ( const Teuchos::ParameterList List)
pure virtual
template<class MatrixType>
virtual void Ifpack2::Container< MatrixType >::apply ( HostView  X,
HostView  Y,
int  blockIndex,
Teuchos::ETransp  mode = Teuchos::NO_TRANS,
SC  alpha = Teuchos::ScalarTraits< SC >::one(),
SC  beta = Teuchos::ScalarTraits< SC >::zero() 
) const
pure virtual

Compute Y := alpha * M^{-1} X + beta*Y.

X is in the domain Map of the original matrix (the argument to compute()), and Y is in the range Map of the original matrix. This method only reads resp. modifies the permuted subset of entries of X resp. Y related to the diagonal block M. That permuted subset is defined by the indices passed into the constructor.

This method is marked const for compatibility with Tpetra::Operator's method of the same name. This might require subclasses to mark some of their instance data as mutable.

Implemented in Ifpack2::ContainerImpl< MatrixType, LocalScalarType >, Ifpack2::ContainerImpl< MatrixType, InverseType::scalar_type >, and Ifpack2::SparseContainer< MatrixType, InverseType >.

template<class MatrixType>
virtual void Ifpack2::Container< MatrixType >::weightedApply ( HostView  X,
HostView  Y,
HostView  D,
int  blockIndex,
Teuchos::ETransp  mode = Teuchos::NO_TRANS,
SC  alpha = Teuchos::ScalarTraits< SC >::one(),
SC  beta = Teuchos::ScalarTraits< SC >::zero() 
) const
pure virtual
template<class MatrixType>
virtual void Ifpack2::Container< MatrixType >::applyInverseJacobi ( const mv_type &  ,
mv_type &  ,
SC  dampingFactor,
bool  ,
int   
) const
pure virtual
template<class MatrixType >
void Ifpack2::Container< MatrixType >::applyMV ( mv_type &  X,
mv_type &  Y 
) const
virtual
template<class MatrixType >
void Ifpack2::Container< MatrixType >::weightedApplyMV ( mv_type &  X,
mv_type &  Y,
vector_type &  W 
) const
virtual
template<class MatrixType>
virtual std::ostream& Ifpack2::Container< MatrixType >::print ( std::ostream &  os) const
pure virtual
template<class MatrixType >
std::string Ifpack2::Container< MatrixType >::getName ( )
static

Returns string describing the container. See Details::ContainerFactory.

template<class MatrixType >
void Ifpack2::Container< MatrixType >::DoGSBlock ( HostView  X,
HostView  Y,
HostView  Y2,
HostView  Resid,
SC  dampingFactor,
LO  i 
) const
protectedvirtual

Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS)

Reimplemented in Ifpack2::ContainerImpl< MatrixType, LocalScalarType >, and Ifpack2::ContainerImpl< MatrixType, InverseType::scalar_type >.

Member Data Documentation

template<class MatrixType>
Teuchos::RCP<const row_matrix_type> Ifpack2::Container< MatrixType >::inputMatrix_
protected

The input matrix to the constructor.

template<class MatrixType>
Teuchos::RCP<const crs_matrix_type> Ifpack2::Container< MatrixType >::inputCrsMatrix_
protected

The input matrix, dynamic cast to CrsMatrix. May be null.

template<class MatrixType>
Teuchos::RCP<const block_crs_matrix_type> Ifpack2::Container< MatrixType >::inputBlockMatrix_
protected

The input matrix, dynamic cast to BlockCrsMatrix. May be null.

template<class MatrixType>
int Ifpack2::Container< MatrixType >::numBlocks_
protected

The number of blocks (partitions) in the container.

template<class MatrixType>
Teuchos::Array<LO> Ifpack2::Container< MatrixType >::blockRows_
protected

Local indices of the rows of the input matrix that belong to this block.

template<class MatrixType>
Teuchos::Array<LO> Ifpack2::Container< MatrixType >::blockSizes_
protected

Number of rows in each block.

template<class MatrixType>
Teuchos::Array<LO> Ifpack2::Container< MatrixType >::blockOffsets_
protected

Starting index in blockRows_ of local row indices for each block.

template<class MatrixType>
Teuchos::RCP<vector_type> Ifpack2::Container< MatrixType >::Diag_
mutableprotected

Diagonal elements.

template<class MatrixType>
bool Ifpack2::Container< MatrixType >::IsParallel_
protected

Whether the problem is distributed across multiple MPI processes.

template<class MatrixType>
LO Ifpack2::Container< MatrixType >::NumLocalRows_
protected

Number of local rows in input matrix.

template<class MatrixType>
GO Ifpack2::Container< MatrixType >::NumGlobalRows_
protected

Number of global rows in input matrix.

template<class MatrixType>
GO Ifpack2::Container< MatrixType >::NumGlobalNonzeros_
protected

Number of nonzeros in input matrix.

template<class MatrixType>
bool Ifpack2::Container< MatrixType >::hasBlockCrs_
protected

Whether the input matrix is a BlockCRS matrix.

template<class MatrixType>
int Ifpack2::Container< MatrixType >::bcrsBlockSize_
protected

If hasBlockCrs_, the number of DOFs per vertex. Otherwise 1.

template<class MatrixType>
bool Ifpack2::Container< MatrixType >::pointIndexed_
protected

(If hasBlockCrs_) Whether the blocks are described using sub-block row indices instead of full block rows.

template<class MatrixType>
LO Ifpack2::Container< MatrixType >::scalarsPerRow_
protected

Number of scalars corresponding to one element of blockRows_ If pointIndexed_, always 1. Otherwise if hasBlockCrs_, then bcrsBlockSize_. Otherwise 1.

template<class MatrixType>
bool Ifpack2::Container< MatrixType >::IsInitialized_
protected

If true, the container has been successfully initialized.

template<class MatrixType>
bool Ifpack2::Container< MatrixType >::IsComputed_
protected

If true, the container has been successfully computed.


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