Actual source code: fmg.c
1: #define PETSCKSP_DLL
3: /*
4: Full multigrid using either additive or multiplicative V or W cycle
5: */
6: #include src/ksp/pc/impls/mg/mgimpl.h
8: EXTERN PetscErrorCode PCMGMCycle_Private(PC_MG **,PetscTruth*);
12: PetscErrorCode PCMGFCycle_Private(PC_MG **mg)
13: {
15: PetscInt i,l = mg[0]->levels;
16: PetscScalar zero = 0.0;
19: /* restrict the RHS through all levels to coarsest. */
20: for (i=l-1; i>0; i--){
21: MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);
22: }
23:
24: /* work our way up through the levels */
25: VecSet(mg[0]->x,zero);
26: for (i=0; i<l-1; i++) {
27: PCMGMCycle_Private(&mg[i],PETSC_NULL);
28: MatInterpolate(mg[i+1]->interpolate,mg[i]->x,mg[i+1]->x);
29: }
30: PCMGMCycle_Private(&mg[l-1],PETSC_NULL);
31: return(0);
32: }
36: PetscErrorCode PCMGKCycle_Private(PC_MG **mg)
37: {
39: PetscInt i,l = mg[0]->levels;
40: PetscScalar zero = 0.0;
43: /* restrict the RHS through all levels to coarsest. */
44: for (i=l-1; i>0; i--){
45: MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);
46: }
47:
48: /* work our way up through the levels */
49: VecSet(mg[0]->x,zero);
50: for (i=0; i<l-1; i++) {
51: if (mg[i]->eventsolve) {PetscLogEventBegin(mg[i]->eventsolve,0,0,0,0);}
52: KSPSolve(mg[i]->smoothd,mg[i]->b,mg[i]->x);
53: if (mg[i]->eventsolve) {PetscLogEventEnd(mg[i]->eventsolve,0,0,0,0);}
54: MatInterpolate(mg[i+1]->interpolate,mg[i]->x,mg[i+1]->x);
55: }
56: if (mg[l-1]->eventsolve) {PetscLogEventBegin(mg[l-1]->eventsolve,0,0,0,0);}
57: KSPSolve(mg[l-1]->smoothd,mg[l-1]->b,mg[l-1]->x);
58: if (mg[l-1]->eventsolve) {PetscLogEventEnd(mg[l-1]->eventsolve,0,0,0,0);}
60: return(0);
61: }