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::Details::TriDiSolver< MatrixType, false > Class Template Reference

Partial specialization for stub=false (the default). More...

#include <Ifpack2_Details_TriDiSolver_decl.hpp>

Inheritance diagram for Ifpack2::Details::TriDiSolver< MatrixType, false >:
Inheritance graph
[legend]

Public Types

Public typedefs
typedef MatrixType matrix_type
 The first template parameter of this class. More...
 
typedef MatrixType::scalar_type scalar_type
 The type of entries in the input (global) matrix. More...
 
typedef
MatrixType::local_ordinal_type 
local_ordinal_type
 The type of local indices in the input (global) matrix. More...
 
typedef
MatrixType::global_ordinal_type 
global_ordinal_type
 The type of global indices in the input (global) matrix. More...
 
typedef MatrixType::node_type node_type
 The Node type of the input (global) matrix. More...
 
typedef Teuchos::ScalarTraits
< scalar_type >::magnitudeType 
magnitude_type
 The type of the absolute value (magnitude) of a scalar_type. More...
 
typedef Tpetra::RowMatrix
< scalar_type,
local_ordinal_type,
global_ordinal_type, node_type
row_matrix_type
 Specialization of Tpetra::RowMatrix used by this class. More...
 
typedef Tpetra::Map
< local_ordinal_type,
global_ordinal_type, node_type
map_type
 Specialization of Tpetra::Map 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

void setParameters (const Teuchos::ParameterList &params)
 Set the solvers's parameters. More...
 
void initialize ()
 Set up the graph structure of this preconditioner. More...
 
bool isInitialized () const
 True if the preconditioner has been successfully initialized, else false. More...
 
void compute ()
 Set up the numerical values in this preconditioner. More...
 
bool isComputed () const
 True if the preconditioner has been successfully computed, else false. More...
 
Teuchos::RCP< const
row_matrix_type
getMatrix () const
 The input matrix given to the constructor. More...
 
void setMatrix (const Teuchos::RCP< const row_matrix_type > &A)
 Change the matrix to precondition. More...
 
int getNumInitialize () const
 The number of calls to initialize(). More...
 
int getNumCompute () const
 The number of calls to compute(). More...
 
int getNumApply () const
 The number of calls to apply(). More...
 
double getInitializeTime () const
 The total time (in seconds) spent in initialize(). More...
 
double getComputeTime () const
 The total time (in seconds) spent in compute(). More...
 
double getApplyTime () const
 The total time (in seconds) spent in apply(). More...
 
Constructor and destructor
 TriDiSolver (const Teuchos::RCP< const row_matrix_type > &matrix)
 Constructor. More...
 
virtual ~TriDiSolver ()
 Destructor (declared virtual for memory safety of derived classes). More...
 
Teuchos::RCP< const map_typegetDomainMap () const
 Implementation of Tpetra::Operator. More...
 
Teuchos::RCP< const map_typegetRangeMap () const
 The range Map of this operator. More...
 
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...
 
Implementation of Teuchos::Describable
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 the given FancyOStream. More...
 
void describeLocal (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
 
- 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...
 
- 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::Details::TriDiSolver< MatrixType, false >

Partial specialization for stub=false (the default).

Member Typedef Documentation

template<class MatrixType >
typedef MatrixType Ifpack2::Details::TriDiSolver< MatrixType, false >::matrix_type

The first template parameter of this class.

This must be a Tpetra::RowMatrix specialization.

template<class MatrixType >
typedef MatrixType::scalar_type Ifpack2::Details::TriDiSolver< MatrixType, false >::scalar_type

The type of entries in the input (global) matrix.

template<class MatrixType >
typedef MatrixType::local_ordinal_type Ifpack2::Details::TriDiSolver< MatrixType, false >::local_ordinal_type

The type of local indices in the input (global) matrix.

template<class MatrixType >
typedef MatrixType::global_ordinal_type Ifpack2::Details::TriDiSolver< MatrixType, false >::global_ordinal_type

The type of global indices in the input (global) matrix.

template<class MatrixType >
typedef MatrixType::node_type Ifpack2::Details::TriDiSolver< MatrixType, false >::node_type

The Node type of the input (global) matrix.

template<class MatrixType >
typedef Teuchos::ScalarTraits<scalar_type>::magnitudeType Ifpack2::Details::TriDiSolver< MatrixType, false >::magnitude_type

The type of the absolute value (magnitude) of a scalar_type.

template<class MatrixType >
typedef Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Details::TriDiSolver< MatrixType, false >::row_matrix_type

Specialization of Tpetra::RowMatrix used by this class.

template<class MatrixType >
typedef Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Details::TriDiSolver< MatrixType, false >::map_type

Specialization of Tpetra::Map used by this class.

Constructor & Destructor Documentation

template<class MatrixType >
Ifpack2::Details::TriDiSolver< MatrixType, false >::TriDiSolver ( const Teuchos::RCP< const row_matrix_type > &  matrix)

Constructor.

Parameters
matrix[in] The original input matrix.
template<class MatrixType >
Ifpack2::Details::TriDiSolver< MatrixType, false >::~TriDiSolver ( )
virtual

Destructor (declared virtual for memory safety of derived classes).

Member Function Documentation

template<class MatrixType >
Teuchos::RCP< const typename TriDiSolver< MatrixType, false >::map_type > Ifpack2::Details::TriDiSolver< MatrixType, false >::getDomainMap ( ) const
virtual

Implementation of Tpetra::Operator.

The domain Map of this operator.

The domain Map describes the distribution of valid input vectors X to the apply() method.

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

template<class MatrixType >
Teuchos::RCP< const typename TriDiSolver< MatrixType, false >::map_type > Ifpack2::Details::TriDiSolver< MatrixType, false >::getRangeMap ( ) const
virtual

The range Map of this operator.

The range Map describes the distribution of valid output vectors Y to the apply() method.

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

template<class MatrixType >
void Ifpack2::Details::TriDiSolver< MatrixType, false >::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

\(M^{-1} \cdot X$, then this method computes \) Y + M^{-1} X \(. The typical case is \) = 0 \( and \) = 1 \(. 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; //@} //! Set the solvers's parameters. void setParameters (const Teuchos::ParameterList& params); \brief Set up the graph structure of this preconditioner. If the graph structure of the constructor's input matrix has changed, or if you have not yet called initialize(), you must call initialize() before you may call compute() or apply(). Thus, initialize() corresponds to the "symbolic factorization" step of a sparse factorization, whether or not the specific preconditioner actually does a sparse factorization. void initialize (); //! True if the preconditioner has been successfully initialized, else false. bool isInitialized () const; \brief Set up the numerical values in this preconditioner. If the values of the constructor's input matrix have changed, or if you have not yet called compute(), you must call compute() before you may call apply(). Thus, compute() corresponds to the "numeric factorization" step of a sparse factorization, whether or not the specific preconditioner actually does a sparse factorization. void compute (); //! True if the preconditioner has been successfully computed, else false. bool isComputed () const; //! The input matrix given to the constructor. Teuchos::RCP<const row_matrix_type> getMatrix () const; //! Change the matrix to precondition. void setMatrix (const Teuchos::RCP<const row_matrix_type>& A); //! The number of calls to initialize(). int getNumInitialize () const; //! The number of calls to compute(). int getNumCompute () const; //! The number of calls to apply(). int getNumApply() const; //! The total time (in seconds) spent in initialize(). double getInitializeTime() const; //! The total time (in seconds) spent in compute(). double getComputeTime() const; //! The total time (in seconds) spent in apply(). double getApplyTime() const; //@} //! @name Implementation of Teuchos::Describable //@{ //! A one-line description of this object. std::string description () const; //! Print the object with some verbosity level to the given FancyOStream. void describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel = Teuchos::Describable::verbLevel_default) const; //@} private: //! Implementation of describe() for "local" (per process) data. void describeLocal (Teuchos::FancyOStream& out, const Teuchos::EVerbosityLevel verbLevel) const; //! Reset stored data. void reset (); \brief Extract the local tridi matrix from the local sparse matrix. \param A_local_tridi [out] The local tridi matrix. \param A_local [out] The local sparse matrix; the result of applying a LocalFilter to the original input matrix A. (If A is already "local," then this need not necessarily be the result of a LocalFilter.) static void extract (Teuchos::SerialTriDiMatrix<int, scalar_type>& A_local_tridi, const row_matrix_type& A_local); \brief Factor the local tridi matrix A in place. \param A [in/out] On input: The local tridi matrix to factor. On output: The resulting LU factorization. \param ipiv [out] The pivot array from the LU factorization (with partial pivoting). Call this after calling extract(). static void factor (Teuchos::SerialTriDiMatrix<int, scalar_type>& A, const Teuchos::ArrayView<int>& ipiv); //! Specialization of Tpetra::MultiVector used in implementation. typedef Tpetra::MultiVector<scalar_type, local_ordinal_type, global_ordinal_type, node_type> MV; //! Specialization of Tpetra::Import used in implementation. typedef Tpetra::Import<local_ordinal_type, global_ordinal_type, node_type> import_type; //! Specialization of Tpetra::Export used in implementation. typedef Tpetra::Export<local_ordinal_type, global_ordinal_type, node_type> export_type; //! Specialization of Teuchos::ScalarTraits used in implementation. typedef Teuchos::ScalarTraits<scalar_type> STS; \brief Post-permutation, post-view version of apply(). apply() first does any necessary subset permutation and view creation (or copying data), then calls this method to solve the linear system with the diagonal block. \param X [in] Subset permutation of the input X of apply(). \param Y [in] Subset permutation of the input/output Y of apply(). void applyImpl (const MV& X, MV& Y, const Teuchos::ETransp mode, const scalar_type alpha, const scalar_type beta) const; //! The (original) input matrix. Teuchos::RCP<const row_matrix_type> A_; //! The result of a LocalFilter on A_. Teuchos::RCP<const row_matrix_type> A_local_; //! The local (tridi) matrix, which compute() extracts and factor() factors in place. Teuchos::SerialTriDiMatrix<int, scalar_type> A_local_tridi_; //! Permutation array from LAPACK (GETRF). Teuchos::Array<int> ipiv_; //! The total time (in seconds) spent in initialize(). double initializeTime_; //! The total time (in seconds) spent in compute(). double computeTime_; //! The total time (in seconds) spent in apply(). mutable double applyTime_; //! Total number of successful initialize() calls. int numInitialize_; //! Total number of successful compute() calls. int numCompute_; //! Total number of successful apply() calls. mutable int numApply_; //! If \c true, the container has been successfully initialized. bool isInitialized_; //! If \c true, the container has been successfully computed. bool isComputed_; }; //! Partial specialization for stub=true. template<class MatrixType> class TriDiSolver<MatrixType, true> : public Ifpack2::Preconditioner<typename MatrixType::scalar_type, typename MatrixType::local_ordinal_type, typename MatrixType::global_ordinal_type, typename MatrixType::node_type>, virtual public Ifpack2::Details::CanChangeMatrix<Tpetra::RowMatrix<typename MatrixType::scalar_type, typename MatrixType::local_ordinal_type, typename MatrixType::global_ordinal_type, typename MatrixType::node_type> > { public: //! \name Public typedefs //@{ \brief The first template parameter of this class. This must be a Tpetra::RowMatrix specialization. typedef MatrixType matrix_type; //! The type of entries in the input (global) matrix. typedef typename MatrixType::scalar_type scalar_type; //! The type of local indices in the input (global) matrix. typedef typename MatrixType::local_ordinal_type local_ordinal_type; //! The type of global indices in the input (global) matrix. typedef typename MatrixType::global_ordinal_type global_ordinal_type; //! The Node type of the input (global) matrix. typedef typename MatrixType::node_type node_type; //! The type of the absolute value (magnitude) of a \c scalar_type. typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitude_type; //! Specialization of Tpetra::RowMatrix used by this class. typedef Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> row_matrix_type; static_assert(std::is_same<MatrixType, row_matrix_type>::value, "Ifpack2::Details::TriDiSolver: The template parameter MatrixType must be a Tpetra::RowMatrix specialization. Please don't use Tpetra::CrsMatrix (a subclass of Tpetra::RowMatrix) here anymore. The constructor can take either a RowMatrix or a CrsMatrix just fine."); //! Specialization of Tpetra::Map used by this class. typedef Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> map_type; //@} //! \name Constructor and destructor //@{ \brief Constructor. \param matrix [in] The original input matrix. TriDiSolver (const Teuchos::RCP<const row_matrix_type>& matrix); //! Destructor (declared virtual for memory safety of derived classes). virtual ~TriDiSolver (); //@} //! Implementation of Tpetra::Operator //@{ \brief The domain Map of this operator. The domain Map describes the distribution of valid input vectors X to the apply() method. Teuchos::RCP<const map_type> getDomainMap () const; \brief The range Map of this operator. The range Map describes the distribution of valid output vectors Y to the apply() method. Teuchos::RCP<const map_type> getRangeMap () const; \brief Apply the preconditioner to X, putting the result in Y. If the result of applying this preconditioner to a vector X is \)^{-1} X$, then this method computes \(\beta Y + \alpha M^{-1} \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 >
void Ifpack2::Details::TriDiSolver< MatrixType, false >::setParameters ( const Teuchos::ParameterList params)
virtual
template<class MatrixType >
void Ifpack2::Details::TriDiSolver< MatrixType, false >::initialize ( )
virtual

Set up the graph structure of this preconditioner.

If the graph structure of the constructor's input matrix has changed, or if you have not yet called initialize(), you must call initialize() before you may call compute() or apply().

Thus, initialize() corresponds to the "symbolic factorization" step of a sparse factorization, whether or not the specific preconditioner actually does a sparse factorization.

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

template<class MatrixType >
bool Ifpack2::Details::TriDiSolver< MatrixType, false >::isInitialized ( ) const
virtual
template<class MatrixType >
void Ifpack2::Details::TriDiSolver< MatrixType, false >::compute ( )
virtual

Set up the numerical values in this preconditioner.

If the values of the constructor's input matrix have changed, or if you have not yet called compute(), you must call compute() before you may call apply().

Thus, compute() corresponds to the "numeric factorization" step of a sparse factorization, whether or not the specific preconditioner actually does a sparse factorization.

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

template<class MatrixType >
bool Ifpack2::Details::TriDiSolver< MatrixType, false >::isComputed ( ) const
virtual
template<class MatrixType >
Teuchos::RCP< const typename TriDiSolver< MatrixType, false >::row_matrix_type > Ifpack2::Details::TriDiSolver< MatrixType, false >::getMatrix ( ) const
virtual
template<class MatrixType >
void Ifpack2::Details::TriDiSolver< MatrixType, false >::setMatrix ( const Teuchos::RCP< const row_matrix_type > &  A)

Change the matrix to precondition.

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

A one-line description of this object.

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

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


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