Ifpack2 Templated Preconditioning Package
Version 1.0
|
Diagonal preconditioner. More...
#include <Ifpack2_Diagonal_decl.hpp>
Public Types | |
typedef Tpetra::RowMatrix < scalar_type, local_ordinal_type, global_ordinal_type, node_type > | row_matrix_type |
Tpetra::RowMatrix specialization used by this class. More... | |
typedef Tpetra::CrsMatrix < scalar_type, local_ordinal_type, global_ordinal_type, node_type > | crs_matrix_type |
Tpetra::CrsMatrix specialization used by this class. More... | |
typedef Tpetra::Vector < scalar_type, local_ordinal_type, global_ordinal_type, node_type > | vector_type |
Tpetra::Vector specialization used by this class. More... | |
typedef Tpetra::Map < local_ordinal_type, global_ordinal_type, node_type > | map_type |
Tpetra::Map specialization used by this class. More... | |
![]() | |
typedef Teuchos::ScalarTraits < MatrixType::scalar_type > ::magnitudeType | magnitude_type |
The type of the magnitude (absolute value) of a matrix entry. More... | |
Public Member Functions | |
Diagonal (const Teuchos::RCP< const row_matrix_type > &A) | |
Constructor that takes a Tpetra::RowMatrix. More... | |
Diagonal (const Teuchos::RCP< const crs_matrix_type > &A_in) | |
Constructor that takes a Tpetra::CrsMatrix. More... | |
Diagonal (const Teuchos::RCP< const vector_type > &diag) | |
Constructor that accepts a Tpetra::Vector of inverse diagonal entries. More... | |
virtual | ~Diagonal () |
Destructor. More... | |
void | setParameters (const Teuchos::ParameterList ¶ms) |
Sets parameters on this object. More... | |
void | initialize () |
Initialize. More... | |
bool | isInitialized () const |
Returns true if the preconditioner has been successfully initialized. More... | |
void | compute () |
Compute the preconditioner. 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... | |
Implementation of Tpetra::Operator | |
void | apply (const Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &X, Tpetra::MultiVector< scalar_type, local_ordinal_type, global_ordinal_type, node_type > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const |
Apply the preconditioner to X, putting the result in Y. More... | |
Teuchos::RCP< const map_type > | getDomainMap () const |
The Tpetra::Map representing this operator's domain. More... | |
Teuchos::RCP< const map_type > | getRangeMap () const |
The Tpetra::Map representing this operator's range. More... | |
Attribute accessor methods | |
Teuchos::RCP< const row_matrix_type > | getMatrix () const |
Return the communicator associated with this matrix operator. More... | |
double | getComputeFlops () const |
Return the number of flops in the computation phase. More... | |
double | getApplyFlops () const |
Return the number of flops for the application of the preconditioner. More... | |
int | getNumInitialize () const |
Return the number of calls to initialize(). More... | |
int | getNumCompute () const |
Return the number of calls to compute(). More... | |
int | getNumApply () const |
Return the number of calls to apply(). More... | |
double | getInitializeTime () const |
Return the time spent in initialize(). More... | |
double | getComputeTime () const |
Return the time spent in compute(). More... | |
double | getApplyTime () const |
Return the time spent in apply(). More... | |
Implementation of Teuchos::Describable | |
std::string | description () const |
Return 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... | |
![]() | |
virtual | ~Preconditioner () |
Destructor. More... | |
virtual void | setZeroStartingSolution (bool zeroStartingSolution) |
Set this preconditioner's parameters. More... | |
![]() | |
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... | |
Diagonal preconditioner.
MatrixType | A specialization of Tpetra::RowMatrix. |
This class wraps a Tpetra::Vector as a diagonal preconditioner. The preconditioner is defined as
\[ z_i = D_i r_i, \]
where \(r\) is the Vector to be preconditioned, \(z\) is the Vector result of applying the preconditioner to \(r\), and \(D_i\) is the i-th element of the scaling vector.
When Diagonal is constructed with a matrix, \(D\) contains the inverted diagonal elements from the matrix. When Diagonal is constructed with a Tpetra::Vector, \(D\) is the caller-supplied vector.
typedef Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Diagonal< MatrixType >::row_matrix_type |
Tpetra::RowMatrix specialization used by this class.
typedef Tpetra::CrsMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Diagonal< MatrixType >::crs_matrix_type |
Tpetra::CrsMatrix specialization used by this class.
typedef Tpetra::Vector<scalar_type, local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Diagonal< MatrixType >::vector_type |
Tpetra::Vector specialization used by this class.
typedef Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Diagonal< MatrixType >::map_type |
Tpetra::Map specialization used by this class.
Ifpack2::Diagonal< MatrixType >::Diagonal | ( | const Teuchos::RCP< const row_matrix_type > & | A | ) |
Constructor that takes a Tpetra::RowMatrix.
A_in | [in] The input matrix. |
Ifpack2::Diagonal< MatrixType >::Diagonal | ( | const Teuchos::RCP< const crs_matrix_type > & | A_in | ) |
Constructor that takes a Tpetra::CrsMatrix.
This constructor exists to avoid "ambiguous constructor" warnings. It does the same thing as the constructor that takes a Tpetra::RowMatrix.
A_in | [in] The input matrix. |
Ifpack2::Diagonal< MatrixType >::Diagonal | ( | const Teuchos::RCP< const vector_type > & | diag | ) |
Constructor that accepts a Tpetra::Vector of inverse diagonal entries.
diag | [in] Vector of inverse diagonal entries. |
If your compiler complains about this constructor being ambigous with the other constructor overload, instead call the free-standing function Ifpack2::createDiagonalPreconditioner which is located at the bottom of this header file. (This issue may arise if this constructor is called with a Teuchos::RCP<Tpetra::Vector>
that isn't const-qualified exactly as declared here.)
|
virtual |
Destructor.
|
virtual |
Sets parameters on this object.
Currently doesn't need any parameters.
|
virtual |
|
inlinevirtual |
Returns true
if the preconditioner has been successfully initialized.
|
virtual |
Compute the preconditioner.
|
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.
|
virtual |
Apply the preconditioner to X, putting the result in Y.
If the result of applying this preconditioner to a vector X is \(F \cdot X\), then this method computes \(\beta Y + \alpha F \cdot X\). The typical case is \(\beta = 0\) and \(\alpha = 1\).
|
virtual |
The Tpetra::Map representing this operator's domain.
|
virtual |
The Tpetra::Map representing this operator's range.
|
inlinevirtual |
Return the communicator associated with this matrix operator.
The original input matrix to be preconditioned.
This could be null, for example if the user created this object using the constructor that takes a Tpetra::Vector, or if the user called setMatrix() with a null input.
double Ifpack2::Diagonal< MatrixType >::getComputeFlops | ( | ) | const |
Return the number of flops in the computation phase.
double Ifpack2::Diagonal< MatrixType >::getApplyFlops | ( | ) | const |
Return the number of flops for the application of the preconditioner.
|
virtual |
Return the number of calls to initialize().
|
virtual |
Return the number of calls to compute().
|
virtual |
Return the number of calls to apply().
|
virtual |
Return the time spent in initialize().
|
virtual |
Return the time spent in compute().
|
virtual |
Return the time spent in apply().
std::string Ifpack2::Diagonal< MatrixType >::description | ( | ) | const |
Return a one-line description of this object.
void Ifpack2::Diagonal< MatrixType >::describe | ( | Teuchos::FancyOStream & | out, |
const Teuchos::EVerbosityLevel | verbLevel = Teuchos::Describable::verbLevel_default |
||
) | const |
Print the object with some verbosity level to an FancyOStream object.