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;
18: /* restrict the RHS through all levels to coarsest. */
19: for (i=l-1; i>0; i--){
20: MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);
21: }
22:
23: /* work our way up through the levels */
24: VecSet(mg[0]->x,0.0);
25: for (i=0; i<l-1; i++) {
26: PCMGMCycle_Private(&mg[i],PETSC_NULL);
27: MatInterpolate(mg[i+1]->interpolate,mg[i]->x,mg[i+1]->x);
28: }
29: PCMGMCycle_Private(&mg[l-1],PETSC_NULL);
30: return(0);
31: }
35: PetscErrorCode PCMGKCycle_Private(PC_MG **mg)
36: {
38: PetscInt i,l = mg[0]->levels;
41: /* restrict the RHS through all levels to coarsest. */
42: for (i=l-1; i>0; i--){
43: MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);
44: }
45:
46: /* work our way up through the levels */
47: VecSet(mg[0]->x,0.0);
48: for (i=0; i<l-1; i++) {
49: if (mg[i]->eventsolve) {PetscLogEventBegin(mg[i]->eventsolve,0,0,0,0);}
50: KSPSolve(mg[i]->smoothd,mg[i]->b,mg[i]->x);
51: if (mg[i]->eventsolve) {PetscLogEventEnd(mg[i]->eventsolve,0,0,0,0);}
52: MatInterpolate(mg[i+1]->interpolate,mg[i]->x,mg[i+1]->x);
53: }
54: if (mg[l-1]->eventsolve) {PetscLogEventBegin(mg[l-1]->eventsolve,0,0,0,0);}
55: KSPSolve(mg[l-1]->smoothd,mg[l-1]->b,mg[l-1]->x);
56: if (mg[l-1]->eventsolve) {PetscLogEventEnd(mg[l-1]->eventsolve,0,0,0,0);}
58: return(0);
59: }