6.4.7. Multigrid Preconditioners

Up: Contents Next: SNES: Nonlinear Solvers and Unconstrained Minimization Previous: Combining Preconditioners

A large suite of routines is available for using multigrid as a preconditioner. In the PC framework the user is required to provide the coarse grid solver, smoothers, restriction, and interpolation, as well as the code to calculate residuals. The PC component allows all of that to be wrapped up into a PETSc compliant preconditioner. We fully support both matrix-free and matrix-based multigrid solvers.

A multigrid preconditioner is created with the four commands

   ierr = SLESCreate(MPI_Comm comm,SLES *sles); 
   ierr = SLESGetPC(SLES sles,PC *pc); 
   ierr = PCSetType(PC pc,PCMG); 
   ierr = MGSetLevels(pc,int levels); 
A large number of parameters affect the multigrid behavior. The command
   ierr = MGSetType(PC pc,MGType mode);  
indicates which form of multigrid to apply [19].

For standard V or W-cycle multigrids, one sets the mode to be MGMULTIPLICATIVE; for the additive form (which in certain cases reduces to the BPX method, or additive multilevel Schwarz, or multilevel diagonal scaling), one uses MGADDITIVE as the mode. For a variant of full multigrid, one can use MGFULL, and for the Kaskade algorithm MGKASKADE. For the multiplicative and full multigrid options, one can use a W-cycle by calling

   ierr = MGSetCycles(PC pc,int cycles); 
with a value of MG_W_CYCLE for cycles. The commands above can also be set from the options database. The option names are -pc_mg_type [multiplicative, additive, full, kaskade], and -pc_mg_cycles cycles.

The user can control the amount of pre- and postsmoothing by using either the options -pc_mg_smoothup m and -pc_mg_smoothdown n or the routines

   ierr = MGSetNumberSmoothUp(PC pc,int m); 
   ierr = MGSetNumberSmoothDown(PC pc,int n); 
Note that if the command MGSetSmoother() (discussed below) has been employed, the same amounts of pre- and postsmoothing will be used.

The remainder of the multigrid routines, which determine the solvers and interpolation/restriction operators that are used, are mandatory. To set the coarse grid solver, one must call

   ierr = MGGetCoarseSolve(PC pc,SLES *sles); 
and set the appropriate options in sles. Similarly, the smoothers are set by calling
   ierr = MGGetSmoother(PC pc,int level,SLES *sles); 
and setting the various options in sles. To use a different pre- and postsmoother, one should call the following routines instead operations to be matrix free (see Section Matrix-Free Matrices ), he or she should make sure that these operations are defined. Note that this system is arranged so that if the interpolation is the transpose of the restriction, the same mat argument can be passed to both MGSetRestriction() and MGSetInterpolation().

On each level except the coarsest, one must also set the routine to compute the residual. The following command suffices:

   MGSetResidual(PC pc,int level,int (*residual)(Mat,Vec,Vec,Vec),Mat mat); 
The residual() function can be set to be MGDefaultResidual() if one's operator is stored in a Mat format. In certain circumstances, where it is much cheaper to calculate the residual directly, rather than through the usual formula b - Ax, the user may wish to provide an alternative.

Finally, the user must provide three work vectors for each level (except on the finest, where only the residual work vector is required). The work vectors are set with the commands

   ierr = MGSetRhs(PC pc,int level,Vec b); 
   ierr = MGSetX(PC pc,int level,Vec x); 
   ierr = MGSetR(PC pc,int level,Vec r); 
The user is responsible for freeing these vectors once the iteration is complete.

One can control the KSP and PC options used on the various levels (as well as the coarse grid) using the prefix mg_levels_ ( mg_coarse_ for the coarse grid). For example,

  -mg_levels_ksp_type cg 
will cause the CG method to be used as the Krylov method for each level. Or
  -mg_levels_pc_type ilu -mg_levels_pc_ilu_levels 2 
will cause the the ILU preconditioner to be used on each level with two levels of fill in the incomplete factorization.


Up: Contents Next: SNES: Nonlinear Solvers and Unconstrained Minimization Previous: Combining Preconditioners