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

Diagonally scaled Chebyshev iteration for Tpetra sparse matrices. More...

#include <Ifpack2_Chebyshev_decl.hpp>

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

Public Types

Typedefs
typedef MatrixType matrix_type
 The template parameter of this class. More...
 
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::device_type 
device_type
 The Kokkos::Device specialization used by the input MatrixType. More...
 
typedef MatrixType::node_type node_type
 The Node type used by 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
 The Tpetra::RowMatrix specialization matching MatrixType. More...
 
typedef Tpetra::Map
< local_ordinal_type,
global_ordinal_type, node_type
map_type
 The Tpetra::Map specialization matching MatrixType. More...
 
typedef Tpetra::Vector
< scalar_type,
local_ordinal_type,
global_ordinal_type, node_type
vector_type
 The Tpetra::Vector specialization matching MatrixType. 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

 Chebyshev (const Teuchos::RCP< const row_matrix_type > &A)
 Constructor. More...
 
virtual ~Chebyshev ()
 Destructor. More...
 
Preconditioner computation methods
void setParameters (const Teuchos::ParameterList &params)
 Set (or reset) parameters. More...
 
void initialize ()
 Initialize the preconditioner. More...
 
bool isInitialized () const
 
void compute ()
 (Re)compute the left scaling, and (if applicable) estimate max and min eigenvalues of D_inv * A. More...
 
bool isComputed () const
 
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, returning the result in Y. More...
 
Teuchos::RCP< const map_typegetDomainMap () const
 The Tpetra::Map representing the domain of this operator. More...
 
Teuchos::RCP< const map_typegetRangeMap () const
 The Tpetra::Map representing the range of this operator. More...
 
bool hasTransposeApply () const
 Whether it's possible to apply the transpose of this operator. More...
 
void applyMat (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) const
 Compute Y = Op(A)*X, where Op(A) is either A, \(A^T\), or \(A^H\). More...
 
Attribute accessor methods
Teuchos::RCP< const
Teuchos::Comm< int > > 
getComm () const
 The communicator over which the matrix is distributed. More...
 
Teuchos::RCP< const
row_matrix_type
getMatrix () const
 The matrix for which this is a preconditioner. More...
 
Teuchos::RCP< const
Tpetra::CrsMatrix< scalar_type,
local_ordinal_type,
global_ordinal_type, node_type > > 
getCrsMatrix () const
 Attempt to return the matrix A as a Tpetra::CrsMatrix. More...
 
double getComputeFlops () const
 The total number of floating-point operations taken by all calls to compute(). More...
 
double getApplyFlops () const
 The total number of floating-point operations taken by all calls to apply(). More...
 
int getNumInitialize () const
 The total number of successful calls to initialize(). More...
 
int getNumCompute () const
 The total number of successful calls to compute(). More...
 
int getNumApply () const
 The total number of successful calls to apply(). More...
 
double getInitializeTime () const
 The total time spent in all calls to initialize(). More...
 
double getComputeTime () const
 The total time spent in all calls to compute(). More...
 
double getApplyTime () const
 The total time spent in all calls to apply(). More...
 
size_t getNodeSmootherComplexity () const
 Get a rough estimate of cost per iteration. More...
 
MatrixType::scalar_type getLambdaMaxForApply () const
 The estimate of the maximum eigenvalue used in the apply(). More...
 
Implementation of Teuchos::Describable
std::string description () const
 A simple 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 a Teuchos::FancyOStream. 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...
 
- 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::Chebyshev< MatrixType >

Diagonally scaled Chebyshev iteration for Tpetra sparse matrices.

Template Parameters
MatrixTypeA specialization of Tpetra::RowMatrix.

Summary

This class implements a Chebyshev polynomial preconditioner or smoother for a Tpetra sparse matrix. Given a matrix A, it applies Chebyshev iteration to the left-scaled matrix \(D^{-1} A\), where D = diag(A) is the matrix of the diagonal entries of A. This class' constructor accepts a Tpetra::RowMatrix or any subclass thereof (including Tpetra::CrsMatrix). Its template parameter may be a specialization of either class.

Chebyshev is derived from Preconditioner, which itself is derived from Tpetra::Operator. Therefore, a Chebyshev instance may be used as an operator in any code that invokes the operator as apply().

Warning
Our implementation currently only works with a real symmetric positive definite (SPD) matrix. Results for matrices that are not SPD, or for complex-valued Scalar types, are not defined.

Algorithm

Given a matrix A, a right-hand side X, and an initial guess Y, this class computes an approximate solution to \(AY=X\) via Chebyshev iteration using the left-scaled matrix \(D^{-1} A\), where D is the matrix of the diagonal elements of A. (You may control left scaling yourself if you wish, by providing an optional vector of the entries of \(D^{-1}\).) While Chebyshev iteration works for any matrix, we have chosen only to allow real-valued, symmetric positive definite matrices.

Chebyshev iteration was originally intended as an iterative solver for linear systems. See the following publication (the spelling of "Chebyshev" in Latin characters differs in some publications):

T. Manteuffel, "The Tchebychev iteration for nonsymmetric linear systems," Numer. Math., 28 (1977), pp. 307-327.

It also works as a smoother for algebraic multigrid, which is the target use case of this implementation.

Eigenvalue bounds

We require that the input matrix A be real-valued and symmetric positive definite. Thus, all of its eigenvalues must lie in a positive interval on the real line. Furthermore, if D is the matrix of the diagonal elements of A, then the same is true of \(D^{-1} A\).

Suppose \([\lambda_{min}, \lambda_{max}]\) is the interval of the eigenvalues of \(D^{-1} A\). Users may either give us an estimate of the maximum eigenvalue \(\lambda_{max}\), or let us compute it (which we do with a few power iterations). They may optionally also give us the (estimated) ratio \(\eta = \lambda_{max} / \lambda_{min}\), or (an estimate of) the minimum eigenvalue \(\lambda_{min}\). The \(\eta\) parameter corresponds to the "smoother: Chebyshev alpha" parameter of ML. (We use "eta" instead of "alpha" to avoid confusion with the "alpha" argument of the apply() method of Tpetra::Operator.)

When using Chebyshev iteration by itself to solve linear systems, it is important to have good estimates of both the minimum and maximum eigenvalues. However, when using a small number of Chebyshev iterations as a smoother in multigrid, the maximum eigenvalue estimate is more important. (The point of a smoother is to smooth out the high-frequency components of the error, that is, those that correspond to the largest eigenvalues. The coarser grids below the current grid will take care of the lower-frequency components of the error.) This is why we use a ratio \(\eta = \lambda_{max} / \lambda_{min}\), rather than requiring a guess for \(\lambda_{min}\). In fact, we only use \(\lambda_{min}\) for error checking, not when determining the Chebyshev coefficients. Often, if users give us \(\lambda_{max}\), our default value of \(\eta\) suffices.

Underestimating \(\lambda_{min}\) may make Chebyshev fail to converge, or fail to reduce the highest-frequency components of the error, if used as a smoother. Thus, we always multiply the given \(\lambda_{min}\) by a small factor (1.1). This heuristic accounts for the fact that typical methods for estimating extremal eigenvalues (like Lanczos or CG) underestimate them.

If you do not give us an estimate for the maximum eigenvalue, we estimate it using a few iterations of the power method in the compute() method. We do not attempt to refine the eigenvalue bounds over Chebyshev iterations, as the typical smoother case does not use very many iterations. For an example of a Chebyshev implementation that updates eigenvalue bound estimates, see Steve Ashby's CHEBYCODE:

S. ASHBY, "CHEBYCODE: A Fortran implementation of Manteuffel's adaptive Chebyshev algorithm," Tech. Rep. UIUCDCS-R-85-1203, University of Illinois, 1985.

Setting parameters

Call the setParameters() method to give this instance your eigenvalue bound estimates (if you have them), as well as to set other options controlling the behavior of Chebyshev iteration. The documentation of setParameters() lists all the parameters that this class accepts. Where possible, we list comparable parameters in the Ifpack package and the ML multigrid package.

Performance

Chebyshev should spend most of its time in Tpetra's native sparse matrix-vector multiply kernel. This should give good performance, since we have spent a lot of effort tuning that kernel. Depending on the Node type of your Tpetra matrix, the kernel may also exploit threads for additional parallelism within each MPI process ("hybrid parallelism" a.k.a. "MPI + X"). If your application depends on hybrid parallelism for performance, you should favor Chebyshev smoothers whenever possible over "serial within a process" smoothers like Gauss-Seidel or SOR (Symmetric Over-Relaxation).

History

The original implementation of this class was an adaptation of ML's ML_Cheby routine. The original author was Ulrich Hetmaniuk, a Sandia employee in what was then (2006) Org 1416. Ifpack2 has seen significant development since then.

Member Typedef Documentation

template<class MatrixType>
typedef MatrixType Ifpack2::Chebyshev< MatrixType >::matrix_type

The template parameter of this class.

template<class MatrixType>
typedef MatrixType::scalar_type Ifpack2::Chebyshev< MatrixType >::scalar_type

The type of the entries of the input MatrixType.

template<class MatrixType>
typedef MatrixType::local_ordinal_type Ifpack2::Chebyshev< MatrixType >::local_ordinal_type

The type of local indices in the input MatrixType.

template<class MatrixType>
typedef MatrixType::global_ordinal_type Ifpack2::Chebyshev< MatrixType >::global_ordinal_type

The type of global indices in the input MatrixType.

template<class MatrixType>
typedef MatrixType::node_type::device_type Ifpack2::Chebyshev< MatrixType >::device_type

The Kokkos::Device specialization used by the input MatrixType.

template<class MatrixType>
typedef MatrixType::node_type Ifpack2::Chebyshev< MatrixType >::node_type

The Node type used by the input MatrixType.

template<class MatrixType>
typedef Teuchos::ScalarTraits<scalar_type>::magnitudeType Ifpack2::Chebyshev< MatrixType >::magnitude_type

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

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

The Tpetra::RowMatrix specialization matching MatrixType.

MatrixType must be a Tpetra::RowMatrix specialization. This typedef will always be a Tpetra::RowMatrix specialization.

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

The Tpetra::Map specialization matching MatrixType.

template<class MatrixType>
typedef Tpetra::Vector<scalar_type, local_ordinal_type, global_ordinal_type, node_type> Ifpack2::Chebyshev< MatrixType >::vector_type

The Tpetra::Vector specialization matching MatrixType.

If you wish to supply setParameters() a precomputed vector of diagonal entries of the matrix, use a pointer to an object of this type.

Constructor & Destructor Documentation

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

Constructor.

Parameters
[in]AThe sparse matrix to which to apply Chebyshev iteration. The matrix A must be square, and its domain Map and range Map must be the same. The latter means that the vectors x and y in the sparse matrix-vector product y = A*x must both have the same distribution over process(es).

We do not require that the row Map and the range Map of A be the same. However, set-up will take less time if they are identical (in terms of pointer equality). This is because we have to extract the diagonal entries of A as a row Map vector: if the row and range Maps are not identical, we have to redistribute the vector from the row Map to the range Map.

The constructor will only check the requirements on the various Maps of A if the CMake configuration option Teuchos_ENABLE_DEBUG was set to ON before building Trilinos. The checks require \(O(1)\) global reductions over all processes in A's communicator, so we prefer to avoid them if we can.

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

Destructor.

Member Function Documentation

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

Set (or reset) parameters.

This method fills in the input ParameterList with missing parameters set to their default values. You may call this method as many times as you want. On each call, the input ParameterList is treated as a complete list of the desired parameters, not as a "delta" or change list from the current set of parameters. (That is, if you remove parameters from the list that were there in the last call to setParameters() and call setParameters() again with the revised list, this method will use default values for the removed parameters, rather than letting the current settings remain.) However, since the method fills in missing parameters, you may keep calling it with the ParameterList used in the previous call in order to get the same behavior as before.

List of parameters

Parameters that govern spectral bounds of the matrix:

  • "chebyshev: max eigenvalue" (ScalarType): lambdaMax, an upper bound of the bounding ellipse of the eigenvalues of the matrix A. If you do not set this parameter, we will compute an approximation. See "Parameters that govern eigenvalue analysis" to control this approximation process.
  • "chebyshev: ratio eigenvalue" (ScalarType): eigRatio, the ratio of lambdaMax to the lower bound of the bounding ellipse of the eigenvalues of A. We use lambdaMax and eigRatio to determine the Chebyshev iteration coefficients. This parameter is optional and defaults to 30.
  • "chebyshev: min eigenvalue" (ScalarType): lambdaMin, a lower bound of real part of bounding ellipse of eigenvalues of the matrix A. This parameter is optional and only used for a quick check if the matrix is the identity matrix (if lambdaMax == lambdaMin == 1).

Parameters that govern the number of Chebyshev iterations:

  • "chebyshev: degree" (int): numIters, the number of iterations. This overrides "relaxation: sweeps" and "smoother: sweeps" (see below).
  • "relaxation: sweeps" (int): numIters, the number of iterations. We include this for compatibility with Ifpack. This overrides "smoother: sweeps" (see below).
  • "smoother: sweeps" (int): numIters, as above. We include this for compatibility with ML.

Parameters that govern eigenvalue analysis:

  • "chebyshev: eigenvalue max iterations" (int): eigMaxIters, the number of power method iterations used to compute the maximum eigenvalue. This overrides "eigen-analysis: iterations" (see below).
  • "eigen-analysis: iterations" (int): eigMaxIters, as above. We include this parameter for compatibility with ML.
  • "eigen-analysis: type" (std::string): The algorithm to use for estimating the max eigenvalue. This parameter is optional. Currently, we only support "power-method" (or "power method"), which is what Ifpack::Chebyshev uses for eigenanalysis. We include this parameter for compatibility with ML.

Parameters that govern other algorithmic details:

  • "chebyshev: assume matrix does not change": Whether compute() should always assume that the matrix has not changed since the last call to compute(). The default is false. If true, compute() will not recompute the inverse diagonal or the estimates of the max and min eigenvalues. compute() will always compute any quantity which the user did not provide and which we have not yet computed before.
  • "chebyshev: operator inv diagonal" (RCP<const V> or const V*): If nonnull, we will use a deep copy of this vector for left scaling as the inverse diagonal of the matrix A, instead of computing the inverse diagonal ourselves. We will make a copy every time you call setParameters(). If you ever call setParameters() without this parameter, we will clear our copy and compute the inverse diagonal ourselves again. If you choose to provide this parameter, you are responsible for updating this if the matrix has changed.
  • "chebyshev: min diagonal value" (ST): minDiagVal. If any entry of the diagonal of the matrix is less than this in magnitude, it will be replaced with this value in the inverse diagonal used for left scaling.
  • "chebyshev: zero starting solution" (bool): If true, then always use the zero vector(s) as the initial guess(es). If false, then apply() will use X on input as the initial guess(es).
Precondition
lambdaMin, lambdaMax, and eigRatio are real
0 < lambdaMin <= lambdaMax
numIters >= 0
eigMaxIters >= 0

Note on compatibility with Ifpack and ML

Both the Ifpack and ML packages implement a Chebyshev smoother. We accept Ifpack and ML names for parameters whenever Ifpack2 has an equivalent parameter. Default settings for parameters relating to spectral bounds come from Ifpack.

The following list maps from an ML parameter to its corresponding Ifpack2 parameter.

  • "smoother: Chebyshev alpha": "chebyshev: ratio eigenvalue"
  • "smoother: sweeps": "chebyshev: degree"

ML does not have a parameter corresponding to "chebyshev: max eigenvalue", because ML estimates the spectral radius automatically. Ifpack and Ifpack2 both can estimate this automatically, but also let the user provide an estimate. Similarly, ML does not have a parameter corresponding to "chebyshev: min eigenvalue".

The following list maps from an Ifpack parameter to its corresponding Ifpack2 parameter. Many of the parameters have the same names, in which case we simply write same.

  • "chebyshev: max eigenvalue": same
  • "chebyshev: ratio eigenvalue": same
  • "chebyshev: min eigenvalue": same
  • "chebyshev: degree": same
  • "relaxation: sweeps": "chebyshev: degree"
  • "chebyshev: min diagonal value": same
  • "relaxation: min diagonal value": "chebyshev: min diagonal value"
  • "chebyshev: zero starting solution": same
  • "relaxation: zero starting solution": "chebyshev: zero starting solution"
  • "chebyshev: operator inv diagonal": same

Details on parameters

The optional user-provided vector of diagonal entries of the matrix may have any distribution for which an Export to the range Map of the matrix is legal. However, if the vector is already distributed according to the range Map, that saves us the communication cost of an Export. We also avoid the Export in case the row Map and the range Map of the matrix are the same. If they are not the same, and if the vector is distributed according to the row Map, we will reuse the Export from the matrix. Otherwise, we have to make a fresh Export object, which is more expensive. To avoid this cost, you should always provide a row Map or range Map vector for this parameter.

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

template<class MatrixType >
void Ifpack2::Chebyshev< MatrixType >::initialize ( )
virtual

Initialize the preconditioner.

The compute() method will call initialize() automatically if it has not yet been called, so you do not normally need to call this. However, it is correct to call initialize() yourself, and compute() will not call it again if it already has been called.

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

template<class MatrixType>
bool Ifpack2::Chebyshev< MatrixType >::isInitialized ( ) const
inlinevirtual
template<class MatrixType >
void Ifpack2::Chebyshev< MatrixType >::compute ( )
virtual

(Re)compute the left scaling, and (if applicable) estimate max and min eigenvalues of D_inv * A.

You must call this method before calling apply(),

  • if you have not yet called this method,
  • if the matrix (either its values or its structure) has changed, or
  • any time after you call setParameters().

Users have the option to supply the left scaling vector D_inv and estimates of the min and max eigenvalues of D_inv * A as parameters to setParameters(). If users did not supply a left scaling, then this method will compute it by default (if assumeMatrixUnchanged is false). Likewise, if users did not supply at least an estimate of the max eigenvalue, this method will estimate it by default. If estimation of the eigenvalues is required, this method may take as long as several Chebyshev iterations.

Advanced users may avoid recomputing the left scaling vector and eigenvalue estimates by setting the "chebyshev: assume matrix does not change" parameter of setParameters() to true. The left scaling vector and eigenvalue estimates will always be computed if the user did not provide them and we have not yet computed them. Any changes to parameters that affect computation of the inverse diagonal or estimation of the eigenvalue bounds will not affect subsequent apply() operations, until the "chebyshev: assume matrix does not change" parameter is set back to false (its default value).

This method will call initialize() if it has not already been called. However, you may call initialize() before calling this method if you wish.

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

template<class MatrixType>
bool Ifpack2::Chebyshev< MatrixType >::isComputed ( ) const
inlinevirtual

Whether compute() has been called at least once.

Note that you must always call compute() if the matrix has changed, if you have called setParameters(), or if you have not yet called compute(). This method only tells you if compute() has been called at least once, not if you need to call compute(). Ifpack2 doesn't have an efficient way to tell if the matrix has changed, so we ask users to tell Ifpack2 if the matrix has changed.

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

template<class MatrixType >
void Ifpack2::Chebyshev< MatrixType >::setMatrix ( const Teuchos::RCP< const row_matrix_type > &  A)
virtual

Change the matrix to be preconditioned.

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

Calling this method 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::Chebyshev< 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, returning the result in Y.

This method actually computes Y = beta*Y + alpha*(M*X), where M*X represents the result of Chebyshev iteration on X, using the matrix Op(A). Op(A) is either A itself, its transpose \(A^T\), or its Hermitian transpose \(A^H\), depending on the mode argument. Since this class currently requires A to be real and symmetric positive definite, it should always be the case that \(A = A^T = A^H\), but we will still respect the mode argument.

Warning
If you did not set the "chebyshev: zero starting solution" parameter to true, then this method will use X as the starting guess for Chebyshev iteration. If you did not initialize X before calling this method, then the resulting solution will be undefined, since it will be computed using uninitialized data.
Parameters
[in]XA (multi)vector to which to apply the preconditioner.
[in,out]YA (multi)vector containing the result of applying the preconditioner to X.
[in]modeIf Teuchos::NO_TRANS, apply the matrix A. If mode is Teuchos::NO_TRANS, apply its transpose \(A^T\). If Teuchos::CONJ_TRANS, apply its Hermitian transpose \(A^H\).
[in]alphaScaling factor for the result of Chebyshev iteration. The default is 1.
[in]betaScaling factor for Y. The default is 0.

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

template<class MatrixType >
Teuchos::RCP< const typename Chebyshev< MatrixType >::map_type > Ifpack2::Chebyshev< MatrixType >::getDomainMap ( ) const
virtual
template<class MatrixType >
Teuchos::RCP< const typename Chebyshev< MatrixType >::map_type > Ifpack2::Chebyshev< MatrixType >::getRangeMap ( ) const
virtual
template<class MatrixType >
bool Ifpack2::Chebyshev< MatrixType >::hasTransposeApply ( ) const

Whether it's possible to apply the transpose of this operator.

template<class MatrixType >
void Ifpack2::Chebyshev< MatrixType >::applyMat ( 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 
) const

Compute Y = Op(A)*X, where Op(A) is either A, \(A^T\), or \(A^H\).

Parameters
[in]XInput (multi)vector of sparse matrix-vector multiply. If mode == Teuchos::NO_TRANS, X must be in the domain Map of the matrix A. Otherwise, X must be in the range Map of A.
[out]YOutput (multi)vector of sparse matrix-vector multiply. If mode == Teuchos::NO_TRANS, Y must be in the range Map of the matrix A. Otherwise, Y must be in the domain Map of A.
[in]modeWhether to apply the matrix A, its transpose \(A^T\), or its conjugate transpose \(A^H\). This method applies A if mode is Teuchos::NO_TRANS, \(A^T\) if mode is Teuchos::TRANS, and \(A^H\) (the Hermitian transpose) if mode is Teuchos::CONJ_TRANS.

Since this class currently requires A to be real and symmetric positive definite, setting mode should not affect the result.

template<class MatrixType >
Teuchos::RCP< const Teuchos::Comm< int > > Ifpack2::Chebyshev< MatrixType >::getComm ( ) const

The communicator over which the matrix is distributed.

template<class MatrixType >
Teuchos::RCP< const typename Chebyshev< MatrixType >::row_matrix_type > Ifpack2::Chebyshev< MatrixType >::getMatrix ( ) const
virtual
template<class MatrixType >
Teuchos::RCP< const Tpetra::CrsMatrix< typename MatrixType::scalar_type, typename MatrixType::local_ordinal_type, typename MatrixType::global_ordinal_type, typename MatrixType::node_type > > Ifpack2::Chebyshev< MatrixType >::getCrsMatrix ( ) const

Attempt to return the matrix A as a Tpetra::CrsMatrix.

This class does not require that A be a Tpetra::CrsMatrix. If it is NOT, this method will return Teuchos::null.

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

The total number of floating-point operations taken by all calls to compute().

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

The total number of floating-point operations taken by all calls to apply().

template<class MatrixType >
int Ifpack2::Chebyshev< MatrixType >::getNumInitialize ( ) const
virtual
template<class MatrixType >
int Ifpack2::Chebyshev< MatrixType >::getNumCompute ( ) const
virtual
template<class MatrixType >
int Ifpack2::Chebyshev< MatrixType >::getNumApply ( ) const
virtual
template<class MatrixType >
double Ifpack2::Chebyshev< MatrixType >::getInitializeTime ( ) const
virtual
template<class MatrixType >
double Ifpack2::Chebyshev< MatrixType >::getComputeTime ( ) const
virtual
template<class MatrixType >
double Ifpack2::Chebyshev< MatrixType >::getApplyTime ( ) const
virtual
template<class MatrixType >
size_t Ifpack2::Chebyshev< MatrixType >::getNodeSmootherComplexity ( ) const

Get a rough estimate of cost per iteration.

template<class MatrixType >
MatrixType::scalar_type Ifpack2::Chebyshev< MatrixType >::getLambdaMaxForApply ( ) const

The estimate of the maximum eigenvalue used in the apply().

template<class MatrixType >
std::string Ifpack2::Chebyshev< MatrixType >::description ( ) const

A simple one-line description of this object.

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

Print the object with some verbosity level to a Teuchos::FancyOStream.


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