Actual source code: pcdmmg.c
1: #define PETSCKSP_DLL
3: /*
4: Defines a multigrid preconditioner that is built from a DM
6: */
7: #include petscmg.h
8: #include petscda.h
13: /*@
14: PCDMMGSetDM - Sets the coarsest DM that is to be used to define the interpolation/restriction
15: for the multigrid preconditioner.
17: Not Collective
19: Input Parameter:
20: + pc - the preconditioner context
21: - dm - the coarsest dm
23: Level: intermediate
25: .keywords: MG, get, levels, multigrid
27: .seealso: PCMG, PCMGSetLevels()
28: @*/
29: PetscErrorCode PCDMMGSetDM(PC pc,DM dm)
30: {
32: PetscInt i,nlevels;
33: Mat R;
34: DM dmf;
35: MPI_Comm comm;
36: PetscTruth flg;
39: PetscTypeCompare((PetscObject)pc,PCMG,&flg);
40: if (!flg) return(0);
42: PetscObjectGetComm((PetscObject)pc,&comm);
43: PCMGGetLevels(pc,&nlevels);
44: PCMGSetGalerkin(pc);
45: PetscObjectReference((PetscObject)dm);
47: /* refine the DM nlevels - 1 and use it to fill up the PCMG restrictions/interpolations */
48: for (i=1; i<nlevels; i++) {
49: DMRefine(dm,comm,&dmf);
50: DMGetInterpolation(dm,dmf,&R,PETSC_NULL);
51: PCMGSetInterpolation(pc,i,R);
52: DMDestroy(dm);
53: dm = dmf;
54: }
55: DMDestroy(dm);
56: return(0);
57: }