Actual source code: fmg.c
1: /*$Id: fmg.c,v 1.26 2001/08/21 21:03:20 bsmith Exp $*/
2: /*
3: Full multigrid using either additive or multiplicative V or W cycle
4: */
5: #include src/ksp/pc/impls/mg/mgimpl.h
7: EXTERN int MGMCycle_Private(MG *,PetscTruth*);
9: /*
10: MGFCycle_Private - Given an MG structure created with MGCreate() runs
11: full multigrid.
13: Iput Parameters:
14: . mg - structure created with MGCreate().
16: Note: This may not be what others call full multigrid. What we
17: do is restrict the rhs to all levels, then starting
18: on the coarsest level work our way up generating
19: initial guess for the next level. This provides an
20: improved preconditioner but not a great improvement.
21: */
24: int MGFCycle_Private(MG *mg)
25: {
26: int i,l = mg[0]->levels,ierr;
27: PetscScalar zero = 0.0;
30: /* restrict the RHS through all levels to coarsest. */
31: for (i=l-1; i>0; i--){
32: MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);
33: }
34:
35: /* work our way up through the levels */
36: VecSet(&zero,mg[0]->x);
37: for (i=0; i<l-1; i++) {
38: MGMCycle_Private(&mg[i],PETSC_NULL);
39: MatInterpolate(mg[i+1]->interpolate,mg[i]->x,mg[i+1]->x);
40: }
41: MGMCycle_Private(&mg[l-1],PETSC_NULL);
42: return(0);
43: }
45: /*
46: MGKCycle_Private - Given an MG structure created with MGCreate() runs
47: full Kascade MG solve.
49: Iput Parameters:
50: . mg - structure created with MGCreate().
52: Note: This may not be what others call Kascadic MG.
53: */
56: int MGKCycle_Private(MG *mg)
57: {
58: int i,l = mg[0]->levels,ierr;
59: PetscScalar zero = 0.0;
62: /* restrict the RHS through all levels to coarsest. */
63: for (i=l-1; i>0; i--){
64: MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);
65: }
66:
67: /* work our way up through the levels */
68: VecSet(&zero,mg[0]->x);
69: for (i=0; i<l-1; i++) {
70: if (mg[i]->eventsolve) {PetscLogEventBegin(mg[i]->eventsolve,0,0,0,0);}
71: KSPSetRhs(mg[i]->smoothd,mg[i]->b);
72: KSPSetSolution(mg[i]->smoothd,mg[i]->x);
73: KSPSolve(mg[i]->smoothd);
74: if (mg[i]->eventsolve) {PetscLogEventEnd(mg[i]->eventsolve,0,0,0,0);}
75: MatInterpolate(mg[i+1]->interpolate,mg[i]->x,mg[i+1]->x);
76: }
77: if (mg[l-1]->eventsolve) {PetscLogEventBegin(mg[l-1]->eventsolve,0,0,0,0);}
78: KSPSetRhs(mg[l-1]->smoothd,mg[l-1]->b);
79: KSPSetSolution(mg[l-1]->smoothd,mg[l-1]->x);
80: KSPSolve(mg[l-1]->smoothd);
81: if (mg[l-1]->eventsolve) {PetscLogEventEnd(mg[l-1]->eventsolve,0,0,0,0);}
83: return(0);
84: }