Ifpack2 Templated Preconditioning Package  Version 1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Member Functions | List of all members
Ifpack2::BlockRelaxation< MatrixType, ContainerType > Class Template Reference

Block relaxation preconditioners (or smoothers) for Tpetra::RowMatrix and Tpetra::CrsMatrix sparse matrices. More...

#include <Ifpack2_BlockRelaxation_decl.hpp>

Inheritance diagram for Ifpack2::BlockRelaxation< MatrixType, ContainerType >:
Inheritance graph
[legend]

Public Member Functions

Teuchos::RCP
< Ifpack2::Partitioner
< Tpetra::RowGraph
< local_ordinal_type,
global_ordinal_type, node_type > > > 
getPartitioner ()
 For diagnostic purposes. More...
 
Preconditioner computation methods
void setParameters (const Teuchos::ParameterList &params)
 Sets all the parameters for the preconditioner. More...
 
bool supportsZeroStartingSolution ()
 
void setZeroStartingSolution (bool zeroStartingSolution)
 Set this preconditioner's parameters. More...
 
Teuchos::RCP< const
Teuchos::ParameterList
getValidParameters () const
 Return a list of all the parameters that this class accepts. More...
 
void initialize ()
 Initialize. More...
 
bool isInitialized () const
 Returns true if the preconditioner has been successfully initialized. More...
 
void compute ()
 compute the preconditioner for the specified matrix, diagonal perturbation thresholds and relaxation parameters. More...
 
bool isComputed () const
 Return true if compute() has been called. More...
 
Implementation of Ifpack2::Details::CanChangeMatrix
virtual void setMatrix (const Teuchos::RCP< const row_matrix_type > &A)
 Change the matrix to be preconditioned. More...
 
Methods implementing the Tpetra::Operator interface.
void apply (const MV &X, MV &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const
 Applies the preconditioner to X, returns the result in Y. More...
 
Teuchos::RCP< const map_type > getDomainMap () const
 Returns the Tpetra::Map object associated with the domain of this operator. More...
 
Teuchos::RCP< const map_type > getRangeMap () const
 Returns the Tpetra::Map object associated with the range of this operator. More...
 
bool hasTransposeApply () const
 
void applyMat (const MV &X, MV &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS) const
 Applies the matrix to a Tpetra::MultiVector. More...
 
Attribute accessor methods
Teuchos::RCP< const
Teuchos::Comm< int > > 
getComm () const
 The communicator over which the input matrix is distributed. More...
 
Teuchos::RCP< const
row_matrix_type
getMatrix () const
 The input matrix of this preconditioner's constructor. More...
 
int getNumInitialize () const
 Returns the number of calls to initialize(). More...
 
int getNumCompute () const
 Returns the number of calls to compute(). More...
 
int getNumApply () const
 Returns the number of calls to apply(). More...
 
double getInitializeTime () const
 Returns the time spent in initialize(). More...
 
double getComputeTime () const
 Returns the time spent in compute(). More...
 
double getApplyTime () const
 Returns the time spent in apply(). More...
 
size_t getNodeSmootherComplexity () const
 Get a rough estimate of cost per iteration. More...
 
Implementation of the Teuchos::Describable interface
std::string description () const
 A one-line description of this object. More...
 
void describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
 Print the object with some verbosity level to an FancyOStream object. More...
 
- Public Member Functions inherited from Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type >
virtual ~Preconditioner ()
 Destructor. More...
 
virtual void apply (const Tpetra::MultiVector< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > &X, Tpetra::MultiVector< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, MatrixType::scalar_typealpha=Teuchos::ScalarTraits< MatrixType::scalar_type >::one(), MatrixType::scalar_typebeta=Teuchos::ScalarTraits< MatrixType::scalar_type >::zero()) const =0
 Apply the preconditioner to X, putting the result in Y. More...
 
- Public Member Functions inherited from Ifpack2::Details::CanChangeMatrix< Tpetra::RowMatrix< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > >
virtual void setMatrix (const Teuchos::RCP< const Tpetra::RowMatrix< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type > > &A)=0
 Set the new matrix. More...
 
virtual ~CanChangeMatrix ()
 Destructor. More...
 

Typedefs

typedef MatrixType::scalar_type scalar_type
 The type of the entries of the input MatrixType. More...
 
typedef
MatrixType::local_ordinal_type 
local_ordinal_type
 The type of local indices in the input MatrixType. More...
 
typedef
MatrixType::global_ordinal_type 
global_ordinal_type
 The type of global indices in the input MatrixType. More...
 
typedef MatrixType::node_type node_type
 Node type of the input MatrixType. More...
 
typedef Teuchos::ScalarTraits
< scalar_type >::magnitudeType 
magnitude_type
 The type of the magnitude (absolute value) of a matrix entry. More...
 
typedef Tpetra::RowMatrix
< scalar_type,
local_ordinal_type,
global_ordinal_type, node_type
row_matrix_type
 Tpetra::RowMatrix specialization corresponding to MatrixType. More...
 
typedef Tpetra::Import
< local_ordinal_type,
global_ordinal_type, node_type
import_type
 Tpetra::Importer specialization for use with MatrixType and compatible MultiVectors. More...
 

Internal typedefs (handy for brevity and code clarity)

 BlockRelaxation (const Teuchos::RCP< const row_matrix_type > &Matrix)
 Constructor. More...
 
virtual ~BlockRelaxation ()
 Destructor. More...
 

Additional Inherited Members

- Public Types inherited from Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type >
typedef Teuchos::ScalarTraits
< MatrixType::scalar_type >
::magnitudeType 
magnitude_type
 The type of the magnitude (absolute value) of a matrix entry. More...
 

Detailed Description

template<class MatrixType, class ContainerType = Container<MatrixType>>
class Ifpack2::BlockRelaxation< MatrixType, ContainerType >

Block relaxation preconditioners (or smoothers) for Tpetra::RowMatrix and Tpetra::CrsMatrix sparse matrices.

Template Parameters
MatrixTypeA specialization of Tpetra::RowMatrix.
ContainerTypeDO NOT SPECIFY THIS EXPLICITLY. This exists only for backwards compatibility.

This class implements the construction and application of block relaxation preconditioners and smoothers, for sparse matrices represented as Tpetra::RowMatrix or Tpetra::CrsMatrix. This class implements Tpetra::Operator, and its apply() method applies the block relaxation.

BlockRelaxation implements block variants of the following relaxations:

For a list of supported parameters, please refer to the documentation of setParameters().

Member Typedef Documentation

template<class MatrixType, class ContainerType = Container<MatrixType>>
typedef MatrixType::scalar_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::scalar_type

The type of the entries of the input MatrixType.

template<class MatrixType, class ContainerType = Container<MatrixType>>
typedef MatrixType::local_ordinal_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::local_ordinal_type

The type of local indices in the input MatrixType.

template<class MatrixType, class ContainerType = Container<MatrixType>>
typedef MatrixType::global_ordinal_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::global_ordinal_type

The type of global indices in the input MatrixType.

template<class MatrixType, class ContainerType = Container<MatrixType>>
typedef MatrixType::node_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::node_type

Node type of the input MatrixType.

template<class MatrixType, class ContainerType = Container<MatrixType>>
typedef Teuchos::ScalarTraits<scalar_type>::magnitudeType Ifpack2::BlockRelaxation< MatrixType, ContainerType >::magnitude_type

The type of the magnitude (absolute value) of a matrix entry.

template<class MatrixType, class ContainerType = Container<MatrixType>>
typedef Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> Ifpack2::BlockRelaxation< MatrixType, ContainerType >::row_matrix_type

Tpetra::RowMatrix specialization corresponding to MatrixType.

template<class MatrixType, class ContainerType = Container<MatrixType>>
typedef Tpetra::Import<local_ordinal_type, global_ordinal_type, node_type> Ifpack2::BlockRelaxation< MatrixType, ContainerType >::import_type

Tpetra::Importer specialization for use with MatrixType and compatible MultiVectors.

Constructor & Destructor Documentation

template<class MatrixType , class ContainerType >
Ifpack2::BlockRelaxation< MatrixType, ContainerType >::BlockRelaxation ( const Teuchos::RCP< const row_matrix_type > &  Matrix)
explicit

Constructor.

Parameters
Matrix[in] The matrix for which to make the constructor. Tpetra::RowMatrix is the base class of Tpetra::CrsMatrix, so you may give either a Tpetra::RowMatrix or a Tpetra::CrsMatrix here.

The results of apply() are undefined if you change the sparse matrix after invoking this constructor, without first calling initialize() and compute() (in that order) to reinitialize the preconditioner.

The "explicit" keyword just means that you must invoke the Relaxation constructor explicitly; you aren't allowed to use it as an implicit conversion ("cast"). For example, you may do this (namespaces and Tpetra template parameters omitted for brevity):

RCP<const CrsMatrix<...> > A = ...;
BlockRelaxation<RowMatrix<...> > R (A);

but you may not do this:

// Declaration of some user-defined function.
void foo (const BlockRelaxation<RowMatrix<...> >& R);
RCP<const CrsMatrix<...> > A = ...;
foo (A);
template<class MatrixType , class ContainerType >
Ifpack2::BlockRelaxation< MatrixType, ContainerType >::~BlockRelaxation ( )
virtual

Destructor.

Member Function Documentation

template<class MatrixType , class ContainerType >
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::setParameters ( const Teuchos::ParameterList params)
virtual

Sets all the parameters for the preconditioner.

Valid parameters are:

  • "relaxation: type"
    Valid values (string):
    • "Jacobi"
    • "Gauss-Seidel"
    • "Symmetric Gauss-Seidel"
  • "relaxation: sweeps" (int)
  • "relaxation: damping factor" (scalar)
  • "relaxation: zero starting solution" (bool)
  • "relaxation: backward mode" (bool)
  • "partitioner: type"
    Valid values (string):
    • "linear"
    • "line"
    • "user"
  • "partitioner: local parts" (local ordinal)
  • "partitioner: overlap" (int)
See Also
Ifpack2::Details::UserPartitioner.

Implements Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type >.

template<class MatrixType, class ContainerType = Container<MatrixType>>
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::setZeroStartingSolution ( bool  zeroStartingSolution)
inlinevirtual
template<class MatrixType , class ContainerType >
Teuchos::RCP< const Teuchos::ParameterList > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getValidParameters ( ) const

Return a list of all the parameters that this class accepts.

template<class MatrixType , class ContainerType >
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::initialize ( )
virtual
template<class MatrixType, class ContainerType = Container<MatrixType>>
bool Ifpack2::BlockRelaxation< MatrixType, ContainerType >::isInitialized ( ) const
inlinevirtual
template<class MatrixType , class ContainerType >
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::compute ( )
virtual

compute the preconditioner for the specified matrix, diagonal perturbation thresholds and relaxation parameters.

Implements Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type >.

template<class MatrixType, class ContainerType = Container<MatrixType>>
bool Ifpack2::BlockRelaxation< MatrixType, ContainerType >::isComputed ( ) const
inlinevirtual
template<class MatrixType , class ContainerType >
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::setMatrix ( const Teuchos::RCP< const row_matrix_type > &  A)
virtual

Change the matrix to be preconditioned.

Parameters
A[in] The new matrix.
Postcondition
! isInitialized ()
! isComputed ()

Calling this method with a matrix different than the current matrix resets the preconditioner's state. After calling this method with a nonnull input, you must first call initialize() and compute() (in that order) before you may call apply().

You may call this method with a null input. If A is null, then you may not call initialize() or compute() until you first call this method again with a nonnull input. This method invalidates any previous factorization whether or not A is null, so calling setMatrix() with a null input is one way to clear the preconditioner's state (and free any memory that it may be using).

The new matrix A need not necessarily have the same Maps or even the same communicator as the original matrix.

template<class MatrixType, class ContainerType = Container<MatrixType>>
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::apply ( const MV &  X,
MV &  Y,
Teuchos::ETransp  mode = Teuchos::NO_TRANS,
scalar_type  alpha = Teuchos::ScalarTraits<scalar_type>::one(),
scalar_type  beta = Teuchos::ScalarTraits<scalar_type>::zero() 
) const

Applies the preconditioner to X, returns the result in Y.

Parameters
X- (In) A Tpetra::MultiVector of dimension NumVectors to be preconditioned.
Y- (InOut) A Tpetra::MultiVector of dimension NumVectors containing result.
Returns
Integer error code, set to 0 if successful.
Warning
This routine is NOT AztecOO compliant.
template<class MatrixType , class ContainerType >
Teuchos::RCP< const Tpetra::Map< typename MatrixType::local_ordinal_type, typename MatrixType::global_ordinal_type, typename MatrixType::node_type > > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getDomainMap ( ) const
virtual
template<class MatrixType , class ContainerType >
Teuchos::RCP< const Tpetra::Map< typename MatrixType::local_ordinal_type, typename MatrixType::global_ordinal_type, typename MatrixType::node_type > > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getRangeMap ( ) const
virtual
template<class MatrixType, class ContainerType = Container<MatrixType>>
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::applyMat ( const MV &  X,
MV &  Y,
Teuchos::ETransp  mode = Teuchos::NO_TRANS 
) const

Applies the matrix to a Tpetra::MultiVector.

Parameters
X- (In) A Tpetra::MultiVector of dimension NumVectors to multiply with matrix.
Y- (Out) A Tpetra::MultiVector of dimension NumVectors containing the result.
template<class MatrixType , class ContainerType >
Teuchos::RCP< const Teuchos::Comm< int > > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getComm ( ) const

The communicator over which the input matrix is distributed.

template<class MatrixType , class ContainerType >
Teuchos::RCP< const Tpetra::RowMatrix< typename MatrixType::scalar_type, typename MatrixType::local_ordinal_type, typename MatrixType::global_ordinal_type, typename MatrixType::node_type > > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getMatrix ( ) const
virtual
template<class MatrixType , class ContainerType >
int Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getNumInitialize ( ) const
virtual
template<class MatrixType , class ContainerType >
int Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getNumCompute ( ) const
virtual
template<class MatrixType , class ContainerType >
int Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getNumApply ( ) const
virtual
template<class MatrixType , class ContainerType >
double Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getInitializeTime ( ) const
virtual
template<class MatrixType , class ContainerType >
double Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getComputeTime ( ) const
virtual
template<class MatrixType , class ContainerType >
double Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getApplyTime ( ) const
virtual
template<class MatrixType , class ContainerType >
size_t Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getNodeSmootherComplexity ( ) const

Get a rough estimate of cost per iteration.

template<class MatrixType , class ContainerType >
std::string Ifpack2::BlockRelaxation< MatrixType, ContainerType >::description ( ) const

A one-line description of this object.

template<class MatrixType , class ContainerType >
void Ifpack2::BlockRelaxation< MatrixType, ContainerType >::describe ( Teuchos::FancyOStream out,
const Teuchos::EVerbosityLevel  verbLevel = Teuchos::Describable::verbLevel_default 
) const

Print the object with some verbosity level to an FancyOStream object.

template<class MatrixType, class ContainerType = Container<MatrixType>>
Teuchos::RCP<Ifpack2::Partitioner<Tpetra::RowGraph<local_ordinal_type,global_ordinal_type,node_type> > > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getPartitioner ( )
inline

For diagnostic purposes.


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