Ifpack2 Templated Preconditioning Package
Version 1.0
|
Block relaxation preconditioners (or smoothers) for Tpetra::RowMatrix and Tpetra::CrsMatrix sparse matrices. More...
#include <Ifpack2_BlockRelaxation_decl.hpp>
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 ¶ms) |
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... | |
Block relaxation preconditioners (or smoothers) for Tpetra::RowMatrix and Tpetra::CrsMatrix sparse matrices.
MatrixType | A specialization of Tpetra::RowMatrix. |
ContainerType | DO 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().
typedef MatrixType::scalar_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::scalar_type |
The type of the entries of the input MatrixType.
typedef MatrixType::local_ordinal_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::local_ordinal_type |
The type of local indices in the input MatrixType.
typedef MatrixType::global_ordinal_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::global_ordinal_type |
The type of global indices in the input MatrixType.
typedef MatrixType::node_type Ifpack2::BlockRelaxation< MatrixType, ContainerType >::node_type |
Node type of the input MatrixType.
typedef Teuchos::ScalarTraits<scalar_type>::magnitudeType Ifpack2::BlockRelaxation< MatrixType, ContainerType >::magnitude_type |
The type of the magnitude (absolute value) of a matrix entry.
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
.
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.
|
explicit |
Constructor.
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):
but you may not do this:
|
virtual |
Destructor.
|
virtual |
Sets all the parameters for the preconditioner.
Valid parameters are:
|
inlinevirtual |
Set this preconditioner's parameters.
Reimplemented from Ifpack2::Preconditioner< MatrixType::scalar_type, MatrixType::local_ordinal_type, MatrixType::global_ordinal_type, MatrixType::node_type >.
Teuchos::RCP< const Teuchos::ParameterList > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getValidParameters | ( | ) | const |
Return a list of all the parameters that this class accepts.
|
virtual |
|
inlinevirtual |
Returns true
if the preconditioner has been successfully initialized.
|
virtual |
compute the preconditioner for the specified matrix, diagonal perturbation thresholds and relaxation parameters.
|
inlinevirtual |
Return true if compute() has been called.
|
virtual |
Change the matrix to be preconditioned.
A | [in] The new matrix. |
! 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.
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.
X | - (In) A Tpetra::MultiVector of dimension NumVectors to be preconditioned. |
Y | - (InOut) A Tpetra::MultiVector of dimension NumVectors containing result. |
|
virtual |
Returns the Tpetra::Map object associated with the domain of this operator.
|
virtual |
Returns the Tpetra::Map object associated with the range of this operator.
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.
X | - (In) A Tpetra::MultiVector of dimension NumVectors to multiply with matrix. |
Y | - (Out) A Tpetra::MultiVector of dimension NumVectors containing the result. |
Teuchos::RCP< const Teuchos::Comm< int > > Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getComm | ( | ) | const |
The communicator over which the input matrix is distributed.
|
virtual |
The input matrix of this preconditioner's constructor.
|
virtual |
Returns the number of calls to initialize().
|
virtual |
Returns the number of calls to compute().
|
virtual |
Returns the number of calls to apply().
|
virtual |
Returns the time spent in initialize().
|
virtual |
Returns the time spent in compute().
|
virtual |
Returns the time spent in apply().
size_t Ifpack2::BlockRelaxation< MatrixType, ContainerType >::getNodeSmootherComplexity | ( | ) | const |
Get a rough estimate of cost per iteration.
std::string Ifpack2::BlockRelaxation< MatrixType, ContainerType >::description | ( | ) | const |
A one-line description of this object.
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.
|
inline |
For diagnostic purposes.