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_MG **mg) 11 { 12 PetscErrorCode ierr; 13 PetscInt i,l = mg[0]->levels; 14 15 PetscFunctionBegin; 16 /* compute RHS on each level */ 17 for (i=l-1; i>0; i--) { 18 if (mg[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mg[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 19 ierr = MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);CHKERRQ(ierr); 20 if (mg[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mg[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 21 } 22 /* solve separately on each level */ 23 for (i=0; i<l; i++) { 24 ierr = VecSet(mg[i]->x,0.0);CHKERRQ(ierr); 25 if (mg[i]->eventsmoothsolve) {ierr = PetscLogEventBegin(mg[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 26 ierr = KSPSolve(mg[i]->smoothd,mg[i]->b,mg[i]->x);CHKERRQ(ierr); 27 if (mg[i]->eventsmoothsolve) {ierr = PetscLogEventEnd(mg[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 28 } 29 for (i=1; i<l; i++) { 30 if (mg[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mg[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 31 ierr = MatInterpolateAdd(mg[i]->interpolate,mg[i-1]->x,mg[i]->x,mg[i]->x);CHKERRQ(ierr); 32 if (mg[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mg[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 33 } 34 PetscFunctionReturn(0); 35 } 36