1 2 /* 3 Additive Multigrid V Cycle routine 4 */ 5 #include <petsc/private/pcmgimpl.h> 6 7 #undef __FUNCT__ 8 #define __FUNCT__ "PCMGACycle_Private" 9 PetscErrorCode PCMGACycle_Private(PC pc,PC_MG_Levels **mglevels) 10 { 11 PetscErrorCode ierr; 12 PetscInt i,l = mglevels[0]->levels; 13 14 PetscFunctionBegin; 15 /* compute RHS on each level */ 16 for (i=l-1; i>0; i--) { 17 if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 18 ierr = MatRestrict(mglevels[i]->restrct,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr); 19 if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 20 } 21 /* solve separately on each level */ 22 for (i=0; i<l; i++) { 23 ierr = VecSet(mglevels[i]->x,0.0);CHKERRQ(ierr); 24 if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 25 ierr = KSPSolve(mglevels[i]->smoothd,mglevels[i]->b,mglevels[i]->x);CHKERRQ(ierr); 26 if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 27 } 28 for (i=1; i<l; i++) { 29 if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 30 ierr = MatInterpolateAdd(mglevels[i]->interpolate,mglevels[i-1]->x,mglevels[i]->x,mglevels[i]->x);CHKERRQ(ierr); 31 if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 32 } 33 PetscFunctionReturn(0); 34 } 35