1 /* 2 Additive Multigrid V Cycle routine 3 */ 4 #include "src/ksp/pc/impls/mg/mgimpl.h" 5 6 /* 7 MGACycle_Private - Given an MG structure created with MGCreate() runs 8 one cycle down through the levels and back up. Applys 9 the smoothers in an additive manner. 10 11 Iput Parameters: 12 . mg - structure created with MGCreate(). 13 14 */ 15 #undef __FUNCT__ 16 #define __FUNCT__ "MGACycle_Private" 17 PetscErrorCode MGACycle_Private(MG *mg) 18 { 19 int i,l = mg[0]->levels,ierr; 20 PetscScalar zero = 0.0; 21 22 PetscFunctionBegin; 23 /* compute RHS on each level */ 24 for (i=l-1; i>0; i--) { 25 ierr = MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);CHKERRQ(ierr); 26 } 27 /* solve seperately on each level */ 28 for (i=0; i<l; i++) { 29 ierr = VecSet(&zero,mg[i]->x);CHKERRQ(ierr); 30 if (mg[i]->eventsolve) {ierr = PetscLogEventBegin(mg[i]->eventsolve,0,0,0,0);CHKERRQ(ierr);} 31 ierr = KSPSolve(mg[i]->smoothd,mg[i]->b,mg[i]->x);CHKERRQ(ierr); 32 if (mg[i]->eventsolve) {ierr = PetscLogEventEnd(mg[i]->eventsolve,0,0,0,0);CHKERRQ(ierr);} 33 } 34 for (i=1; i<l; i++) { 35 ierr = MatInterpolateAdd(mg[i]->interpolate,mg[i-1]->x,mg[i]->x,mg[i]->x);CHKERRQ(ierr); 36 } 37 PetscFunctionReturn(0); 38 } 39