Actual source code: petscmg.h

  1: /*
  2:       Structure used for Multigrid preconditioners 
  3: */
 6:  #include petscksp.h

  9: /*E
 10:     MGType - Determines the type of multigrid method that is run.

 12:    Level: beginner

 14:    Values:
 15: +  MGMULTIPLICATIVE (default) - traditional V or W cycle as determined by MGSetCycles()
 16: .  MGADDITIVE - the additive multigrid preconditioner where all levels are
 17:                 smoothed before updating the residual
 18: .  MGFULL - same as multiplicative except one also performs grid sequencing, 
 19:             that is starts on the coarsest grid, performs a cycle, interpolates
 20:             to the next, performs a cycle etc
 21: -  MGKASKADE - like full multigrid except one never goes back to a coarser level
 22:                from a finer

 24: .seealso: MGSetType()

 26: E*/
 27: typedef enum { MGMULTIPLICATIVE,MGADDITIVE,MGFULL,MGKASKADE } MGType;
 28: #define MGCASCADE MGKASKADE;

 30: #define MG_V_CYCLE     1
 31: #define MG_W_CYCLE     2

 33: EXTERN PetscErrorCode MGSetType(PC,MGType);
 34: EXTERN PetscErrorCode MGCheck(PC);
 35: EXTERN PetscErrorCode MGSetLevels(PC,PetscInt,MPI_Comm*);
 36: EXTERN PetscErrorCode MGGetLevels(PC,PetscInt*);

 38: EXTERN PetscErrorCode MGSetNumberSmoothUp(PC,PetscInt);
 39: EXTERN PetscErrorCode MGSetNumberSmoothDown(PC,PetscInt);
 40: EXTERN PetscErrorCode MGSetCycles(PC,PetscInt);
 41: EXTERN PetscErrorCode MGSetCyclesOnLevel(PC,PetscInt,PetscInt);

 43: EXTERN PetscErrorCode MGGetSmoother(PC,PetscInt,KSP*);
 44: EXTERN PetscErrorCode MGGetSmootherDown(PC,PetscInt,KSP*);
 45: EXTERN PetscErrorCode MGGetSmootherUp(PC,PetscInt,KSP*);
 46: EXTERN PetscErrorCode MGGetCoarseSolve(PC,KSP*);

 48: EXTERN PetscErrorCode MGSetRhs(PC,PetscInt,Vec);
 49: EXTERN PetscErrorCode MGSetX(PC,PetscInt,Vec);
 50: EXTERN PetscErrorCode MGSetR(PC,PetscInt,Vec);

 52: EXTERN PetscErrorCode MGSetRestriction(PC,PetscInt,Mat);
 53: EXTERN PetscErrorCode MGSetInterpolate(PC,PetscInt,Mat);
 54: EXTERN PetscErrorCode MGSetResidual(PC,PetscInt,PetscErrorCode (*)(Mat,Vec,Vec,Vec),Mat);
 55: EXTERN PetscErrorCode MGDefaultResidual(Mat,Vec,Vec,Vec);


 59: #endif