This class is a wrapper to be used when the library Petsc is installed and used with OFELI. More...
Public Member Functions | |
PETScWrapper (int argc, char **args, string help="") | |
Constructor with program arguments. More... | |
~PETScWrapper () | |
Destructor. More... | |
PetscErrorCode | getIntOption (string s, PetscInt &n, PetscBool &set) const |
Get an option as an integer number. More... | |
PetscErrorCode | getBoolOption (string s, PetscBool &b, PetscBool &set) const |
Get an option as a bool variable. More... | |
PetscMPIInt | size () const |
Return wrapper size, i.e. number of processors. | |
void | setMesh (Mesh &ms) |
Set mesh. More... | |
void | setPartition (Partition &p) |
Set mesh partition. More... | |
void | setMatrix (PETScMatrix< T_ > &A) |
Define problem matrix. More... | |
void | setLinearSystem (PETScMatrix< T_ > &A, PETScVect< T_ > &b, string s=KSPCG, string p=PCJACOBI, real_t tol=1.e-12, size_t max_it=1000) |
Set linear system features. More... | |
void | setPreconditioner (string p) |
Choose preconditioner for the iterative procedure. More... | |
void | setIterationParameters (real_t tol, size_t max_it) |
Choose iteration parameters. More... | |
void | setIterationMethod (string m) |
Choose the iterative method. More... | |
void | solve (PETScVect< T_ > &x) |
Solve the linear system. More... | |
void | solve (const PETScVect< T_ > &b, PETScVect< T_ > &x) |
Solve the linear system. More... | |
void | checkError (PETScVect< T_ > &u) const |
Check residual error. More... | |
int | getIterationNumber () const |
Return the number of iterations. | |
void | setLSTolerances (real_t rel_tol, real_t abs_tol, real_t div_tol=PETSC_DEFAULT, int max_it=PETSC_DEFAULT) const |
Set tolerance parameters for a linear system. More... | |
PetscMPIInt | getRank () const |
Return the rank of the current processor. | |
Friends | |
template<class S_ > | |
ostream & | operator<< (ostream &s, const PETScWrapper< S_ > &w) |
Output wrapper information. | |
Detailed Description
template<class T_>
class OFELI::PETScWrapper< T_ >
This class is a wrapper to be used when the library Petsc is installed and used with OFELI.
When Petsc is used, an instance of class PETScWrapper is to be declared. It initializes the use of Petsc and enables calling solver functions in Petsc.
- Template Parameters
-
T_ Data type (double, int, complex<double>, ...)
When a linear system is invoked, the choice of iterative solvers can be made among the following methods (see Petsc documentation for more details):
-
KSPRICHARDSON
: The Richardson iterative method (Default damping parameter is1.0
) -
KSPCHEBYSHEV
: The Chebyshev iterative method -
KSPCG
: The conjugate gradient method [Default] -
KSPCGNE
: The CG method for normal equations (without explicitly forming the productA^TA
-
KSPGMRES
: [Default] The GMRES iterative method (see A Generalized Minimal Residual Algorithm for Solving Nonsymmetric Linear Systems. Y. Saad and M. H. Schultz, SIAM J. Sci. Stat. Comput., Vo|. 7, No. 3, July 1986, pp. 856-869) -
KSPFGMRES
: The Flexible GMRES method (with restart) -
KSPLGMRES
: The 'augmented' standard GMRES method where the subspace uses approximations to the error from previous restart cycles -
KSPTCQMR
: A variant of QMR (quasi minimal residual) developed by Tony Chan -
KSPBCGS
: The BiCGStab (Stabilized version of BiConjugate Gradient Squared) method -
KSPIBCGS
: The IBiCGStab (Improved Stabilized version of BiConjugate Gradient Squared) method in an alternative form to have only a single global reduction operation instead of the usual 3 (or 4) -
KSPFBCGS
: The flexible BiCGStab method. -
KSPCGS
: The CGS (Conjugate Gradient Squared) method -
KSPTFQMR
: A transpose free QMR (quasi minimal residual) -
KSPCR
: The conjugate residuals method -
KSPLSQR
: The LSQR method -
KSPBICG
: The Biconjugate gradient method (similar to running the conjugate gradient on the normal equations) -
KSPMINRES
: The MINRES (Minimum Residual) method -
KSPSYMMLQ
: The SYMMLQ method -
KSPGCR
: The Generalized Conjugate Residual method
When a linear system is invoked, the choice of a preconditioner can be made among the following methods (see Petsc documentation for more details):
-
PCJACOBI
: [Default] Jacobi (i.e. diagonal scaling) preconditioning -
PCBJACOBI
: Block Jacobi preconditioning, each block is (approximately) solved with its own KSP object -
PCSOR
: (S)SOR (successive over relaxation, Gauss-Seidel) preconditioning -
PCEISENSTAT
: An implementation of SSOR (symmetric successive over relaxation, symmetric Gauss-Seidel) preconditioning that incorporates Eisenstat's trick to reduce the amount of computation needed -
PCICC
: Incomplete Cholesky factorization preconditioners -
PCILU
: Incomplete factorization preconditioners -
PCASM
: Use the (restricted) additive Schwarz method, each block is (approximately) solved with its own KSP object -
PCLU
: Uses a direct solver, based on LU factorization, as a preconditioner -
PCCHOLESKY
: Uses a direct solver, based on Cholesky factorization, as a preconditioner
Constructor & Destructor Documentation
PETScWrapper | ( | int | argc, |
char ** | args, | ||
string | help = "" |
||
) |
Constructor with program arguments.
- Parameters
-
[in] argc Count of number of command line arguments [in] args The command line arguments. Here is the list of arguments: -
-start_in_debugger [noxterm,dbx,xdb,gdb,...]
- Starts program in debugger
-
-on_error_attach_debugger [noxterm,dbx,xdb,gdb,...]
- Starts debugger when error detected
-
-on_error_emacs <machinename>
causes emacsclient to jump to error file- . -on_error_abort calls abort() when error detected (no traceback)
-
-on_error_mpiabort
calls MPI_abort() when error detected- . -error_output_stderr prints error messages to stderr instead of the default stdout
-
-error_output_none
does not print the error messages (but handles errors in the same way as if this was not called)- . -debugger_nodes [node1,node2,...] - Indicates nodes to start in debugger
-
-debugger_pause [sleeptime]
(in seconds)- Pauses debugger
-
-stop_for_debugger
- Print message on how to attach debugger manually to process and wait (
-debugger_pause
) seconds for attachment
- Print message on how to attach debugger manually to process and wait (
-
-malloc
- Indicates use of PETSc error-checking malloc (on by default for debug version of libraries)
-
-malloc no
- Indicates not to use error-checking malloc
-
-malloc_debug
- check for memory corruption at EVERY malloc or free
-
-malloc_dump
- prints a list of all unfreed memory at the end of the run
-
-malloc_test
- like
-malloc_dump
-malloc_debug
, but only active for debugging builds
- like
-
-fp_trap
- Stops on floating point exceptions (Note that on the IBM RS6000 this slows code by at least a factor of 10.)
-
-no_signal_handler
- Indicates not to trap error signals
-
-shared_tmp
- indicates
/tmp
directory is shared by all processors
- indicates
-
-not_shared_tmp
- each processor has own
/tmp
- each processor has own
-
-tmp
- alternative name of
/tmp
directory
- alternative name of
-
-get_total_flops
- returns total flops done by all processors
-
-memory_info
- Print memory usage at end of run
[in] help String that contains message to display when argument -v
is used -
- Warning
- This class is available only when OFELI has been installed with Petsc
~PETScWrapper | ( | ) |
Destructor.
Destroy the KSP context and release memory allocated by petsc
Member Function Documentation
PetscErrorCode getIntOption | ( | string | s, |
PetscInt & | n, | ||
PetscBool & | set | ||
) | const |
Get an option as an integer number.
- Parameters
-
[in] s String to preprend the name of the option [out] n Obtained integer value [out] set true
if found,false
if not.
PetscErrorCode getBoolOption | ( | string | s, |
PetscBool & | b, | ||
PetscBool & | set | ||
) | const |
Get an option as a bool variable.
- Parameters
-
[in] s String to preprend the name of the option [out] b Obtained boolean value [out] set true
if found,false
if not.
void setPartition | ( | Partition & | p | ) |
Set mesh partition.
This function is to be used for parallel computing
- Parameters
-
[in] p Partition instance
void setMatrix | ( | PETScMatrix< T_ > & | A | ) |
Define problem matrix.
- Parameters
-
[in] A PETScMatrix instance that contains matrix
void setLinearSystem | ( | PETScMatrix< T_ > & | A, |
PETScVect< T_ > & | b, | ||
string | s = KSPCG , |
||
string | p = PCJACOBI , |
||
real_t | tol = 1.e-12 , |
||
size_t | max_it = 1000 |
||
) |
Set linear system features.
- Parameters
-
[in] A PETScMatrix instance that contains matrix [in] b Vector containing the right-hand side [in] s Option to choose iterative solver. See the definition of the class for iterative methods [in] p Option to choose preconditioner. See the definition of the class for available preconditioners. [in] tol Tolerance for convergence of iteration process [Default: 1.e-12
][in] max_it Maximum number of linear solver iterations [Default: 1000
]
void setPreconditioner | ( | string | p | ) |
Choose preconditioner for the iterative procedure.
- Parameters
-
[in] p Option to choose preconditioner. See the definition of the class for available preconditioners.
void setIterationParameters | ( | real_t | tol, |
size_t | max_it | ||
) |
Choose iteration parameters.
- Parameters
-
[in] tol Tolerance for convergence of iteration process [in] max_it Maximum number of linear solver iterations
void setIterationMethod | ( | string | m | ) |
Choose the iterative method.
- Parameters
-
[in] m Option to choose iterative solver. See the definition of the class for available iterative solvers.
void solve | ( | PETScVect< T_ > & | x | ) |
Solve the linear system.
If the member functions setIterationMethod and setPreconditioner have not been used, default methods are used
- Parameters
-
[in,out] x Vector containing the initial guess on input and, if convergence is achieved, the solution on output
Solve the linear system.
If the member functions setIterationMethod and setPreconditioner have not been used, default methods are used
- Parameters
-
[in] b Vector containing the right-hand side [in,out] x Vector containing the initial guess on input and, if convergence is achieved, the solution on output
void checkError | ( | PETScVect< T_ > & | u | ) | const |
Check residual error.
This function computes the residual A*x - b and outputs the number of iterations
- Parameters
-
[out] u Residual vector
void setLSTolerances | ( | real_t | rel_tol, |
real_t | abs_tol, | ||
real_t | div_tol = PETSC_DEFAULT , |
||
int | max_it = PETSC_DEFAULT |
||
) | const |
Set tolerance parameters for a linear system.
- Parameters
-
[in] rel_tol Relative convergence tolerance, relative decrease in the preconditioned residual norm [in] abs_tol Absolute convergence tolerance of the preconditioned residual norm [in] div_tol Divergence tolerance: Amount preconditioned residual norm [in] max_it Maximum number of iterations