1 /* 2 Full multigrid using either additive or multiplicative V or W cycle 3 */ 4 #include "../src/ksp/pc/impls/mg/mgimpl.h" 5 6 EXTERN PetscErrorCode PCMGMCycle_Private(PC,PC_MG_Levels **,PCRichardsonConvergedReason*); 7 8 #undef __FUNCT__ 9 #define __FUNCT__ "PCMGFCycle_Private" 10 PetscErrorCode PCMGFCycle_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 /* restrict the RHS through all levels to coarsest. */ 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 24 /* work our way up through the levels */ 25 ierr = VecSet(mglevels[0]->x,0.0);CHKERRQ(ierr); 26 for (i=0; i<l-1; i++) { 27 ierr = PCMGMCycle_Private(pc,&mglevels[i],PETSC_NULL);CHKERRQ(ierr); 28 if (mg->eventinterprestrict) {ierr = PetscLogEventBegin(mg->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 29 ierr = MatInterpolate(mglevels[i+1]->interpolate,mglevels[i]->x,mglevels[i+1]->x);CHKERRQ(ierr); 30 if (mg->eventinterprestrict) {ierr = PetscLogEventEnd(mg->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 31 } 32 ierr = PCMGMCycle_Private(pc,&mglevels[l-1],PETSC_NULL);CHKERRQ(ierr); 33 PetscFunctionReturn(0); 34 } 35 36 #undef __FUNCT__ 37 #define __FUNCT__ "PCMGKCycle_Private" 38 PetscErrorCode PCMGKCycle_Private(PC pc,PC_MG_Levels **mglevels) 39 { 40 PC_MG *mg = (PC_MG*)pc->data; 41 PetscErrorCode ierr; 42 PetscInt i,l = mglevels[0]->levels; 43 44 PetscFunctionBegin; 45 /* restrict the RHS through all levels to coarsest. */ 46 for (i=l-1; i>0; i--){ 47 if (mg->eventinterprestrict) {ierr = PetscLogEventBegin(mg->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 48 ierr = MatRestrict(mglevels[i]->restrct,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr); 49 if (mg->eventinterprestrict) {ierr = PetscLogEventEnd(mg->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 50 } 51 52 /* work our way up through the levels */ 53 ierr = VecSet(mglevels[0]->x,0.0);CHKERRQ(ierr); 54 for (i=0; i<l-1; i++) { 55 if (mg->eventsmoothsolve) {ierr = PetscLogEventBegin(mg->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 56 ierr = KSPSolve(mglevels[i]->smoothd,mglevels[i]->b,mglevels[i]->x);CHKERRQ(ierr); 57 if (mg->eventsmoothsolve) {ierr = PetscLogEventEnd(mg->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 58 if (mg->eventinterprestrict) {ierr = PetscLogEventBegin(mg->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 59 ierr = MatInterpolate(mglevels[i+1]->interpolate,mglevels[i]->x,mglevels[i+1]->x);CHKERRQ(ierr); 60 if (mg->eventinterprestrict) {ierr = PetscLogEventEnd(mg->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 61 } 62 if (mg->eventsmoothsolve) {ierr = PetscLogEventBegin(mg->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 63 ierr = KSPSolve(mglevels[l-1]->smoothd,mglevels[l-1]->b,mglevels[l-1]->x);CHKERRQ(ierr); 64 if (mg->eventsmoothsolve) {ierr = PetscLogEventEnd(mg->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 65 66 PetscFunctionReturn(0); 67 } 68 69 70