IFPACK  Development
 All Classes Namespaces Files Functions Variables Enumerations Friends Pages
Public Member Functions | List of all members
Ifpack_AdditiveSchwarz< T > Class Template Reference

Ifpack_AdditiveSchwarz: a class to define Additive Schwarz preconditioners of Epetra_RowMatrix's. More...

#include <Ifpack_AdditiveSchwarz.h>

Inheritance diagram for Ifpack_AdditiveSchwarz< T >:
Inheritance graph
[legend]
Collaboration diagram for Ifpack_AdditiveSchwarz< T >:
Collaboration graph
[legend]

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_CommComm () const
 Returns a pointer to the Epetra_Comm communicator associated with this operator.
 
virtual const Epetra_MapOperatorDomainMap () const
 Returns the Epetra_Map object associated with the domain of this operator.
 
virtual const Epetra_MapOperatorRangeMap () 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_RowMatrixMatrix () 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_TimeTime_
 Object used for timing purposes.
 
Teuchos::RefCountPtr< T > Inverse_
 Pointer to the local solver.
 

Detailed Description

template<typename T>
class Ifpack_AdditiveSchwarz< T >

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.

Author
Marzio Sala, SNL 9214.
Date
Last modified on 22-Jan-05.

Definition at line 142 of file Ifpack_AdditiveSchwarz.h.

Constructor & Destructor Documentation

template<typename T >
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).

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

Member Function Documentation

template<typename T >
int Ifpack_AdditiveSchwarz< T >::Apply ( const Epetra_MultiVector X,
Epetra_MultiVector Y 
) const
virtual

Applies the matrix to X, returns the result in Y.

Parameters
X- (In) A Epetra_MultiVector of dimension NumVectors to multiply with matrix.
Y-(Out) A Epetra_MultiVector of dimension NumVectors containing the result.
Returns
Integer error code, set to 0 if successful.

Implements Epetra_Operator.

Definition at line 946 of file Ifpack_AdditiveSchwarz.h.

template<typename T >
int Ifpack_AdditiveSchwarz< T >::ApplyInverse ( const Epetra_MultiVector X,
Epetra_MultiVector Y 
) const
virtual

Applies the preconditioner to X, returns the result in Y.

Parameters
X- (In) A Epetra_MultiVector of dimension NumVectors to be preconditioned.
Y-(Out) A Epetra_MultiVector of dimension NumVectors containing result.
Returns
Integer error code, set to 0 if successful.
Warning
In order to work with AztecOO, any implementation of this method must support the case where X and Y are the same object.

Implements Ifpack_Preconditioner.

Definition at line 1004 of file Ifpack_AdditiveSchwarz.h.

References Insert.

template<typename T >
int Ifpack_AdditiveSchwarz< T >::SetParameters ( Teuchos::ParameterList &  List)
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:
    • Add: Components on the receiving processor will be added together;
    • Zero: Off-processor components will be ignored;
    • Insert: Off-processor components will be inserted into locations on receiving processor replacing existing values.
    • Average: Off-processor components will be averaged with existing;
    • AbsMax: Magnitudes of Off-processor components will be maxed with magnitudes of existing components on the receiving processor.
  • "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().

template<typename T >
int Ifpack_AdditiveSchwarz< T >::SetUseTranspose ( bool  UseTranspose_in)
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.

Parameters
UseTranspose_in- (In) If true, multiply by the transpose of operator, otherwise just use operator.
Returns
Integer error code, set to 0 if successful. Set to -1 if this implementation does not support transpose.

Implements Epetra_Operator.

Definition at line 931 of file Ifpack_AdditiveSchwarz.h.


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