Ifpack_AdditiveSchwarz: a class to define Additive Schwarz preconditioners of Epetra_RowMatrix's. More...
#include <Ifpack_AdditiveSchwarz.h>
Public Member Functions | |
virtual bool | IsInitialized () const |
Returns true if the preconditioner has been successfully initialized. | |
virtual bool | IsComputed () const |
Returns true if the preconditioner has been successfully computed. | |
virtual int | SetParameters (Teuchos::ParameterList &List) |
Sets the parameters. More... | |
Ifpack_AdditiveSchwarz (Epetra_RowMatrix *Matrix_in, int OverlapLevel_in=0) | |
Ifpack_AdditiveSchwarz constructor with given Epetra_RowMatrix. More... | |
virtual | ~Ifpack_AdditiveSchwarz () |
Destructor. | |
virtual int | SetUseTranspose (bool UseTranspose_in) |
If set true, transpose of this operator will be applied (not implemented). More... | |
virtual int | Apply (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Applies the matrix to X, returns the result in Y. More... | |
virtual int | ApplyInverse (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const |
Applies the preconditioner to X, returns the result in Y. More... | |
virtual double | NormInf () const |
Returns the infinity norm of the global matrix (not implemented) | |
virtual const char * | Label () const |
Returns a character string describing the operator. | |
virtual bool | UseTranspose () const |
Returns the current UseTranspose setting. | |
virtual bool | HasNormInf () const |
Returns true if the this object can provide an approximate Inf-norm, false otherwise. | |
virtual const Epetra_Comm & | Comm () const |
Returns a pointer to the Epetra_Comm communicator associated with this operator. | |
virtual const Epetra_Map & | OperatorDomainMap () const |
Returns the Epetra_Map object associated with the domain of this operator. | |
virtual const Epetra_Map & | OperatorRangeMap () const |
Returns the Epetra_Map object associated with the range of this operator. | |
virtual int | Initialize () |
Initialized the preconditioner. | |
virtual int | Compute () |
Computes the preconditioner. | |
virtual double | Condest (const Ifpack_CondestType CT=Ifpack_Cheap, const int MaxIters=1550, const double Tol=1e-9, Epetra_RowMatrix *Matrix_in=0) |
Computes the estimated condition number and returns its value. | |
virtual double | Condest () const |
Returns the estimated condition number, or -1.0 if not computed. | |
virtual const Epetra_RowMatrix & | Matrix () const |
Returns a refernence to the internally stored matrix. | |
virtual bool | IsOverlapping () const |
Returns true is an overlapping matrix is present. | |
virtual std::ostream & | Print (std::ostream &) const |
Prints major information about this preconditioner. | |
virtual const T * | Inverse () const |
virtual int | NumInitialize () const |
Returns the number of calls to Initialize(). | |
virtual int | NumCompute () const |
Returns the number of calls to Compute(). | |
virtual int | NumApplyInverse () const |
Returns the number of calls to ApplyInverse(). | |
virtual double | InitializeTime () const |
Returns the time spent in Initialize(). | |
virtual double | ComputeTime () const |
Returns the time spent in Compute(). | |
virtual double | ApplyInverseTime () const |
Returns the time spent in ApplyInverse(). | |
virtual double | InitializeFlops () const |
Returns the number of flops in the initialization phase. | |
virtual double | ComputeFlops () const |
Returns the number of flops in the computation phase. | |
virtual double | ApplyInverseFlops () const |
Returns the number of flops in the application of the preconditioner. | |
virtual int | OverlapLevel () const |
Returns the level of overlap. | |
virtual const Teuchos::ParameterList & | List () const |
Returns a reference to the internally stored list. | |
Protected Member Functions | |
Ifpack_AdditiveSchwarz (const Ifpack_AdditiveSchwarz &RHS) | |
Copy constructor (should never be used) | |
int | Setup () |
Sets up the localized matrix and the singleton filter. | |
Protected Attributes | |
Teuchos::RefCountPtr< const Epetra_RowMatrix > | Matrix_ |
Pointers to the matrix to be preconditioned. | |
Teuchos::RefCountPtr < Ifpack_OverlappingRowMatrix > | OverlappingMatrix_ |
Pointers to the overlapping matrix. | |
Teuchos::RefCountPtr < Ifpack_LocalFilter > | LocalizedMatrix_ |
Localized version of Matrix_ or OverlappingMatrix_. | |
std::string | Label_ |
Contains the label of this object. | |
bool | IsInitialized_ |
If true, the preconditioner has been successfully initialized. | |
bool | IsComputed_ |
If true, the preconditioner has been successfully computed. | |
bool | UseTranspose_ |
If true , solve with the transpose (not supported by all solvers). | |
bool | IsOverlapping_ |
If true, overlapping is used. | |
int | OverlapLevel_ |
Level of overlap among the processors. | |
Teuchos::ParameterList | List_ |
Stores a copy of the list given in SetParameters() | |
Epetra_CombineMode | CombineMode_ |
Combine mode for off-process elements (only if overlap is used) | |
double | Condest_ |
Contains the estimated condition number. | |
bool | ComputeCondest_ |
If true , compute the condition number estimate each time Compute() is called. | |
bool | UseReordering_ |
If true , reorder the local matrix. | |
std::string | ReorderingType_ |
Type of reordering of the local matrix. | |
Teuchos::RefCountPtr < Ifpack_Reordering > | Reordering_ |
Pointer to a reordering object. | |
Teuchos::RefCountPtr < Ifpack_ReorderFilter > | ReorderedLocalizedMatrix_ |
Pointer to the reorderd matrix. | |
bool | FilterSingletons_ |
Filter for singletons. | |
Teuchos::RefCountPtr < Ifpack_SingletonFilter > | SingletonFilter_ |
filtering object. | |
int | NumInitialize_ |
Contains the number of successful calls to Initialize(). | |
int | NumCompute_ |
Contains the number of successful call to Compute(). | |
int | NumApplyInverse_ |
Contains the number of successful call to ApplyInverse(). | |
double | InitializeTime_ |
Contains the time for all successful calls to Initialize(). | |
double | ComputeTime_ |
Contains the time for all successful calls to Compute(). | |
double | ApplyInverseTime_ |
Contains the time for all successful calls to ApplyInverse(). | |
double | InitializeFlops_ |
Contains the number of flops for Initialize(). | |
double | ComputeFlops_ |
Contains the number of flops for Compute(). | |
double | ApplyInverseFlops_ |
Contain sthe number of flops for ApplyInverse(). | |
Teuchos::RefCountPtr< Epetra_Time > | Time_ |
Object used for timing purposes. | |
Teuchos::RefCountPtr< T > | Inverse_ |
Pointer to the local solver. | |
Ifpack_AdditiveSchwarz: a class to define Additive Schwarz preconditioners of Epetra_RowMatrix's.
Class Ifpack_AdditiveSchwarz enables the construction of Additive Schwarz (one-level overlapping domain decomposition) preconditioners, for a given Epetra_RowMatrix. Ifpack_AdditiveSchwarz is derived from Ifpack_Preconditioner, itself derived from Epetra_Operator. An application of the Additive Schwarz preconditioner can be obtained by calling method ApplyInverse().
One-level overlapping domain decomposition preconditioners use local solvers, of Dirichlet type. This means that the inverse of the local matrix (with minimal or wider overlap) is applied to the residual to be preconditioned.
The preconditioner can be written as:
\[ P_{AS}^{-1} = \sum_{i=1}^M P_i A_i^{-1} R_i , \]
where \(M\) is the number of subdomains (that is, the number of processors in the computation), \(R_i\) is an operator that restricts the global vector to the vector lying on subdomain \(i\), \(P_i\) is the prolongator operator, and
\[ A_i = R_i A P_i. \]
The construction of Schwarz preconditioners is mainly composed by two steps:
The local matrix \(A_i\) can be filtered, to eliminate singletons, and reordered. At the present time, RCM and METIS can be used to reorder the local matrix.
The complete list of supported parameters is reported in page List of Supported Parameters.
Definition at line 142 of file Ifpack_AdditiveSchwarz.h.
Ifpack_AdditiveSchwarz< T >::Ifpack_AdditiveSchwarz | ( | Epetra_RowMatrix * | Matrix_in, |
int | OverlapLevel_in = 0 |
||
) |
Ifpack_AdditiveSchwarz constructor with given Epetra_RowMatrix.
Creates an Ifpack_AdditiveSchwarz preconditioner with overlap. To use minimal-overlap, OverlappingMatrix is omitted (as defaulted to 0).
Matrix | - (In) Pointer to matrix to be preconditioned |
OverlappingMatrix | - (In) Pointer to the matrix extended with the desired level of overlap. |
Definition at line 509 of file Ifpack_AdditiveSchwarz.h.
References Ifpack_AdditiveSchwarz< T >::IsOverlapping_, Ifpack_AdditiveSchwarz< T >::Matrix_, Ifpack_AdditiveSchwarz< T >::OverlapLevel_, Ifpack_AdditiveSchwarz< T >::SetParameters(), and Zero.
|
virtual |
Applies the matrix to X, returns the result in Y.
X | - (In) A Epetra_MultiVector of dimension NumVectors to multiply with matrix. |
Y | -(Out) A Epetra_MultiVector of dimension NumVectors containing the result. |
Implements Epetra_Operator.
Definition at line 946 of file Ifpack_AdditiveSchwarz.h.
|
virtual |
Applies the preconditioner to X, returns the result in Y.
X | - (In) A Epetra_MultiVector of dimension NumVectors to be preconditioned. |
Y | -(Out) A Epetra_MultiVector of dimension NumVectors containing result. |
Implements Ifpack_Preconditioner.
Definition at line 1004 of file Ifpack_AdditiveSchwarz.h.
References Insert.
|
virtual |
Sets the parameters.
Sets the parameter for the additive Schwarz preconditioner, as well as for all the preconditioners that may need to be defined on each subblock. Parameters accepted by List are:
"schwarz: combine mode"
: It must be an Epetra_CombineMode. Default: Zero
. It Can be assume of the following values:"schwarz: compute condest"
: if true
, Compute()
will estimate the condition number of the preconditioner. Default: true
. Implements Ifpack_Preconditioner.
Definition at line 719 of file Ifpack_AdditiveSchwarz.h.
References AbsMax, Add, Average, Insert, InsertAdd, and Zero.
Referenced by Ifpack_AdditiveSchwarz< T >::Ifpack_AdditiveSchwarz().
|
virtual |
If set true, transpose of this operator will be applied (not implemented).
This flag allows the transpose of the given operator to be used
implicitly.
UseTranspose_in | - (In) If true, multiply by the transpose of operator, otherwise just use operator. |
Implements Epetra_Operator.
Definition at line 931 of file Ifpack_AdditiveSchwarz.h.