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