1 /*
2 Data structure used for Multigrid preconditioner.
3 */
4 #if !defined(__MG_IMPL)
5 #define __MG_IMPL
6 #include <petsc/private/pcimpl.h>
7 #include <petscksp.h>
8 #define PETSC_MG_MAXLEVELS 10
9 /*
10 Each level has its own copy of this data.
11 Level (0) is always the coarsest level and Level (levels-1) is the finest.
12 */
13 typedef struct {
14 PetscInt cycles; /* Type of cycle to run: 1 V 2 W */
15 PetscInt level; /* level = 0 coarsest level */
16 PetscInt levels; /* number of active levels used */
17 Vec b; /* Right hand side */
18 Vec x; /* Solution */
19 Vec r; /* Residual */
20 Vec *coarseSpace; /* A vector space which should be accurately captured by the next coarser mesh,
21 and thus accurately interpolated. This array should have the same size on each
22 level, and the vectors should correspond to the same function discretized in
23 the sequence of spaces. */
24
25 PetscErrorCode (*residual)(Mat,Vec,Vec,Vec);
26
27 Mat A; /* matrix used in forming residual*/
28 KSP smoothd; /* pre smoother */
29 KSP smoothu; /* post smoother */
30 Mat interpolate;
31 Mat restrct; /* restrict is a reserved word in C99 and on Cray */
32 Mat inject; /* Used for moving state if provided. */
33 Vec rscale; /* scaling of restriction matrix */
34 PetscLogEvent eventsmoothsetup; /* if logging times for each level */
35 PetscLogEvent eventsmoothsolve;
36 PetscLogEvent eventresidual;
37 PetscLogEvent eventinterprestrict;
38 } PC_MG_Levels;
39
40 /*
41 This data structure is shared by all the levels.
42 */
43 typedef struct {
44 PCMGType am; /* Multiplicative, additive or full */
45 PetscInt cyclesperpcapply; /* Number of cycles to use in each PCApply(), multiplicative only*/
46 PetscInt maxlevels; /* total number of levels allocated */
47 PCMGGalerkinType galerkin; /* use Galerkin process to compute coarser matrices */
48 PetscBool usedmfornumberoflevels; /* sets the number of levels by getting this information out of the DM */
49
50 PetscBool adaptInterpolation; /* flag to adapt the interpolator based upon the coarseSpace */
51 PCMGCoarseSpaceType coarseSpaceType; /* Type of coarse space: polynomials, harmonics, eigenvectors, ... */
52 PetscInt Nc; /* The number of vectors in coarseSpace */
53 PetscInt eigenvalue; /* Key for storing the eigenvalue as a scalar in the eigenvector Vec */
54 PetscBool mespMonitor; /* flag to monitor the multilevel eigensolver */
55
56 PetscInt nlevels;
57 PC_MG_Levels **levels;
58 PetscInt default_smoothu; /* number of smooths per level if not over-ridden */
59 PetscInt default_smoothd; /* with calls to KSPSetTolerances() */
60 PetscReal rtol,abstol,dtol,ttol; /* tolerances for when running with PCApplyRichardson_MG */
61
62 void *innerctx; /* optional data for preconditioner, like PCEXOTIC that inherits off of PCMG */
63 PetscLogStage stageApply;
64 PetscErrorCode (*view)(PC,PetscViewer); /* GAMG and other objects that use PCMG can set their own viewer here */
65 PetscReal min_eigen_DinvA[PETSC_MG_MAXLEVELS];
66 PetscReal max_eigen_DinvA[PETSC_MG_MAXLEVELS];
67 } PC_MG;
68
69 PETSC_INTERN PetscErrorCode PCSetUp_MG(PC);
70 PETSC_INTERN PetscErrorCode PCDestroy_MG(PC);
71 PETSC_INTERN PetscErrorCode PCSetFromOptions_MG(PetscOptionItems *PetscOptionsObject,PC);
72 PETSC_INTERN PetscErrorCode PCView_MG(PC,PetscViewer);
73 PETSC_INTERN PetscErrorCode PCMGGetLevels_MG(PC,PetscInt *);
74 PETSC_INTERN PetscErrorCode PCMGSetLevels_MG(PC,PetscInt,MPI_Comm *);
PCMGResidual_Default(Mat A,Vec b,Vec x,Vec r)75 PETSC_DEPRECATED_FUNCTION("Use PCMGResidualDefault() (since version 3.5)") PETSC_STATIC_INLINE PetscErrorCode PCMGResidual_Default(Mat A,Vec b,Vec x,Vec r) {
76 return PCMGResidualDefault(A,b,x,r);
77 }
78
79 PETSC_INTERN PetscErrorCode DMSetBasisFunction_Internal(PetscInt, PetscBool, PetscInt, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal[], PetscInt, PetscScalar *, void *));
80 PETSC_INTERN PetscErrorCode PCMGComputeCoarseSpace_Internal(PC, PetscInt, PCMGCoarseSpaceType, PetscInt, const Vec[], Vec *[]);
81 PETSC_INTERN PetscErrorCode PCMGAdaptInterpolator_Internal(PC, PetscInt, KSP, KSP, PetscInt, Vec[], Vec[]);
82 PETSC_INTERN PetscErrorCode PCMGRecomputeLevelOperators_Internal(PC, PetscInt);
83
84
85 #endif
86