KSPGCR#
Implements the preconditioned flexible Generalized Conjugate Residual method. Flexible Krylov Methods,
Options Database Key#
-ksp_gcr_restart
- the number of stored vectors to orthogonalize against
Notes#
The GCR Krylov method supports non-symmetric matrices and permits the use of a preconditioner which may vary from one iteration to the next.
Users can can define a method to vary the
preconditioner between iterates via KSPGCRSetModifyPC()
.
Restarts are solves with x0 not equal to zero. When a restart occurs, the initial starting solution is given by the current estimate for x which was obtained by the last restart iterations of the GCR algorithm.
Unlike KSPGMRES
and KSPFGMRES
, when using GCR, the solution and residual vector can be directly accessed at any iterate,
with zero computational cost, via a call to KSPBuildSolution()
and KSPBuildResidual()
respectively.
This implementation of GCR will only apply the stopping condition test whenever ksp->its > ksp->chknorm,
where ksp->chknorm is specified via the command line argument -ksp_check_norm_iteration or via
the function KSPSetCheckNormIteration()
. Hence the residual norm reported by the monitor and stored
in the residual history will be listed as 0.0 before this iteration. It is actually not 0.0; just not calculated.
The method implemented requires the storage of 2 x restart + 1 vectors, twice as much as KSPGMRES
.
Support only for right preconditioning.
Contributed by#
Dave May
References#
**** -*** S. C. Eisenstat, H. C. Elman, and H. C. Schultz. Variational iterative methods for nonsymmetric systems of linear equations. SIAM J. Numer. Anal., 20, 1983
See Also#
KSP: Linear System Solvers, Flexible Krylov Methods, KSPCreate()
, KSPSetType()
, KSPType
, KSP
, KSPGCRSetRestart()
, KSPGCRGetRestart()
,
KSPGCRSetRestart()
, KSPGCRSetModifyPC()
, KSPGMRES
, KSPFGMRES
Level#
beginner
Location#
Examples#
src/dm/impls/stag/tutorials/ex4.c
Index of all KSP routines
Table of Contents for all manual pages
Index of all manual pages