6.3. Krylov Methods

Up: Contents Next: Preconditioning within KSP Previous: Solving Successive Linear Systems

The Krylov subspace methods accept a number of options, many of which are discussed below. First, to set the Krylov subspace method that is to be used, one calls the command

   ierr = KSPSetType(KSP ksp,KSPType method); 
The type can be one of KSPRICHARDSON, KSPCHEBYCHEV, KSPCG, KSPGMRES, KSPTCQMR, KSPBCGS, KSPCGS, KSPTFQMR, KSPCR, KSPLSQR, KSPBICG, or KSPPREONLY. The KSP method can also be set with the options database command -ksp_type, followed by one of the options richardson, chebychev, cg, gmres, tcqmr, bcgs, cgs, tfqmr, cr, lsqr, bicg, or preonly. There are method-specific options for the Richardson, Chebychev, and GMRES methods.
   ierr = KSPRichardsonSetScale(KSP ksp,double damping_factor); 
   ierr = KSPChebychevSetEigenvalues(KSP ksp,double emax,double emin); 
   ierr = KSPGMRESSetRestart(KSP ksp,int max_steps); 
The default parameter values are damping_factor=1.0, emax=0.01, emin=100.0, and max_steps=30. The GMRES restart and Richardson damping factor can also be set with the options -ksp_gmres_restart <n> and -ksp_richardson_scale <factor>.

The default technique for orthogonalization of the Hessenberg matrix in GMRES is the iterative refinement Gram-Schmidt method. This can be set by using the command line option -ksp_gmres_irorthog. Or via

   ierr = KSPGMRESSetOrthogonalization(KSP ksp, 
                               KSPGMRESModifiedGramSchmidtOrthogonalization); 
A slightly faster approachis to use the unmodified (classical) Gram-Schmidt method, which can be set with
   ierr = KSPGMRESSetOrthogonalization(KSP ksp, 
                               KSPGMRESUnmodifiedGramSchmidtOrthogonalization); 
or the options database command -ksp_gmres_unmodifiedgramschmidt. Note that this algorithm is numerically unstable, but may deliver slightly better speed performance. One can also use modifed Gram-Schmidt, by setting the orthogonalization routine, KSPGMRESModifiedGramSchmidtOrthogonalization(), by using the command line option -ksp_gmres_modifiedgramschmidt.

For the conjugate gradient method with complex numbers, there are two slightly different algorithms depending on whether the matrix is Hermitian symmetric or truly symmetric (the default is to assume that it is Hermitian symmetric). To indicate that it is symmetric, one uses the command

   ierr = KSPCGSetType(KSP ksp,KSPCGType KSP_CG_SYMMETRIC); 
Note that this option is not valid for all matrices.

The LSQR algorithm does not involve a preconditioner, any preconditioner set to work with the KSP object is ignored if LSQR was selected.

By default, KSP assumes an initial guess of zero by zeroing the initial value for the solution vector that is given; this zeroing is done at the call to SLESSolve() (or KSPSolve()). To use a nonzero initial guess, the user must call

   ierr = KSPSetInitialGuessNonzero(KSP ksp); 


Up: Contents Next: Preconditioning within KSP Previous: Solving Successive Linear Systems