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

Diagonal preconditioner. More...

#include <Ifpack2_Diagonal_decl.hpp>

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

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...
 
- 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...
 

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 &params)
 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_typegetDomainMap () const
 The Tpetra::Map representing this operator's domain. More...
 
Teuchos::RCP< const map_typegetRangeMap () 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...
 
- 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 setZeroStartingSolution (bool zeroStartingSolution)
 Set this preconditioner's parameters. 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...
 

Detailed Description

template<class MatrixType>
class Ifpack2::Diagonal< MatrixType >

Diagonal preconditioner.

Template Parameters
MatrixTypeA 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.

Member Typedef Documentation

template<class MatrixType >
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.

template<class MatrixType >
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.

template<class MatrixType >
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.

template<class MatrixType >
typedef Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Diagonal< MatrixType >::map_type

Tpetra::Map specialization used by this class.

Constructor & Destructor Documentation

template<class MatrixType >
Ifpack2::Diagonal< MatrixType >::Diagonal ( const Teuchos::RCP< const row_matrix_type > &  A)

Constructor that takes a Tpetra::RowMatrix.

Parameters
A_in[in] The input matrix.
template<class MatrixType >
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.

Parameters
A_in[in] The input matrix.
template<class MatrixType >
Ifpack2::Diagonal< MatrixType >::Diagonal ( const Teuchos::RCP< const vector_type > &  diag)

Constructor that accepts a Tpetra::Vector of inverse diagonal entries.

Parameters
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.)

template<class MatrixType >
Ifpack2::Diagonal< MatrixType >::~Diagonal ( )
virtual

Destructor.

Member Function Documentation

template<class MatrixType >
void Ifpack2::Diagonal< MatrixType >::setParameters ( const Teuchos::ParameterList params)
virtual
template<class MatrixType >
void Ifpack2::Diagonal< MatrixType >::initialize ( )
virtual
template<class MatrixType >
bool Ifpack2::Diagonal< MatrixType >::isInitialized ( ) const
inlinevirtual
template<class MatrixType >
void Ifpack2::Diagonal< MatrixType >::compute ( )
virtual
template<class MatrixType >
bool Ifpack2::Diagonal< MatrixType >::isComputed ( ) const
inlinevirtual
template<class MatrixType >
void Ifpack2::Diagonal< MatrixType >::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 >
void Ifpack2::Diagonal< MatrixType >::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
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\).

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

template<class MatrixType >
Teuchos::RCP< const typename Diagonal< MatrixType >::map_type > Ifpack2::Diagonal< MatrixType >::getDomainMap ( ) const
virtual
template<class MatrixType >
Teuchos::RCP< const typename Diagonal< MatrixType >::map_type > Ifpack2::Diagonal< MatrixType >::getRangeMap ( ) const
virtual
template<class MatrixType >
Teuchos::RCP<const row_matrix_type> Ifpack2::Diagonal< MatrixType >::getMatrix ( ) const
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.

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

template<class MatrixType >
double Ifpack2::Diagonal< MatrixType >::getComputeFlops ( ) const

Return the number of flops in the computation phase.

template<class MatrixType >
double Ifpack2::Diagonal< MatrixType >::getApplyFlops ( ) const

Return the number of flops for the application of the preconditioner.

template<class MatrixType >
int Ifpack2::Diagonal< MatrixType >::getNumInitialize ( ) const
virtual
template<class MatrixType >
int Ifpack2::Diagonal< MatrixType >::getNumCompute ( ) const
virtual
template<class MatrixType >
int Ifpack2::Diagonal< MatrixType >::getNumApply ( ) const
virtual
template<class MatrixType >
double Ifpack2::Diagonal< MatrixType >::getInitializeTime ( ) const
virtual
template<class MatrixType >
double Ifpack2::Diagonal< MatrixType >::getComputeTime ( ) const
virtual
template<class MatrixType >
double Ifpack2::Diagonal< MatrixType >::getApplyTime ( ) const
virtual
template<class MatrixType >
std::string Ifpack2::Diagonal< MatrixType >::description ( ) const

Return a one-line description of this object.

template<class MatrixType >
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.


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