Actual source code: mgimpl.h

  1: /*
  2:       Data structure used for Multigrid preconditioner.
  3: */
 6:  #include private/pcimpl.h
 7:  #include petscmg.h
 8:  #include petscksp.h


 11: /*
 12:      Structure for abstract multigrid solver. 

 14:      Level (0) is always the coarsest level and Level (levels-1) is the finest.
 15: */
 16: typedef struct
 17: {
 18:   PCMGType   am;                           /* Multiplicative, additive or full */
 19:   int        cycles;                       /* Number cycles to run */
 20:   int        level;                        /* level = 0 coarsest level */
 21:   int        levels;                       /* number of active levels used */
 22:   int        maxlevels;                    /* total number of levels allocated */
 23:   PetscTruth galerkin;                     /* use Galerkin process to compute coarser matrices */
 24:   PetscTruth galerkinused;                 /* destroy the Mat created by the Galerkin process */
 25:   Vec        b;                            /* Right hand side */
 26:   Vec        x;                            /* Solution */
 27:   Vec        r;                            /* Residual */
 28:   PetscErrorCode (*residual)(Mat,Vec,Vec,Vec);
 29:   Mat        A;                            /* matrix used in forming residual*/
 30:   KSP        smoothd;                      /* pre smoother */
 31:   KSP        smoothu;                      /* post smoother */
 32:   Mat        interpolate;
 33:   Mat        restrct;                      /* restrict is a reserved word on the Cray!!!*/
 34:   int        default_smoothu;              /* number of smooths per level if not over-ridden */
 35:   int        default_smoothd;              /*  with calls to KSPSetTolerances() */
 36:   PetscReal  rtol,abstol,dtol,ttol;        /* tolerances for when running with PCApplyRichardson_MG */
 37:   PetscEvent eventsetup;                   /* if logging times for each level */
 38:   PetscEvent eventsolve;
 39: }  PC_MG;


 42: #endif