14b9ad928SBarry Smith /* 24b9ad928SBarry Smith Full multigrid using either additive or multiplicative V or W cycle 34b9ad928SBarry Smith */ 4af0996ceSBarry Smith #include <petsc/private/pcmgimpl.h> 54b9ad928SBarry Smith 609573ac7SBarry Smith extern PetscErrorCode PCMGMCycle_Private(PC,PC_MG_Levels**,PCRichardsonConvergedReason*); 74b9ad928SBarry Smith 8f3fbd535SBarry Smith PetscErrorCode PCMGFCycle_Private(PC pc,PC_MG_Levels **mglevels) 94b9ad928SBarry Smith { 106849ba73SBarry Smith PetscErrorCode ierr; 11f3fbd535SBarry Smith PetscInt i,l = mglevels[0]->levels; 124b9ad928SBarry Smith 134b9ad928SBarry Smith PetscFunctionBegin; 144b9ad928SBarry Smith /* restrict the RHS through all levels to coarsest. */ 154b9ad928SBarry Smith for (i=l-1; i>0; i--) { 1663e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 17f3fbd535SBarry Smith ierr = MatRestrict(mglevels[i]->restrct,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr); 1863e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 194b9ad928SBarry Smith } 204b9ad928SBarry Smith 214b9ad928SBarry Smith /* work our way up through the levels */ 22f3fbd535SBarry Smith ierr = VecSet(mglevels[0]->x,0.0);CHKERRQ(ierr); 234b9ad928SBarry Smith for (i=0; i<l-1; i++) { 240298fd71SBarry Smith ierr = PCMGMCycle_Private(pc,&mglevels[i],NULL);CHKERRQ(ierr); 2563e6d426SJed Brown if (mglevels[i+1]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i+1]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 26f3fbd535SBarry Smith ierr = MatInterpolate(mglevels[i+1]->interpolate,mglevels[i]->x,mglevels[i+1]->x);CHKERRQ(ierr); 2763e6d426SJed Brown if (mglevels[i+1]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i+1]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 284b9ad928SBarry Smith } 290298fd71SBarry Smith ierr = PCMGMCycle_Private(pc,&mglevels[l-1],NULL);CHKERRQ(ierr); 304b9ad928SBarry Smith PetscFunctionReturn(0); 314b9ad928SBarry Smith } 324b9ad928SBarry Smith 3331567311SBarry Smith PetscErrorCode PCMGKCycle_Private(PC pc,PC_MG_Levels **mglevels) 344b9ad928SBarry Smith { 356849ba73SBarry Smith PetscErrorCode ierr; 36f3fbd535SBarry Smith PetscInt i,l = mglevels[0]->levels; 374b9ad928SBarry Smith 384b9ad928SBarry Smith PetscFunctionBegin; 394b9ad928SBarry Smith /* restrict the RHS through all levels to coarsest. */ 404b9ad928SBarry Smith for (i=l-1; i>0; i--) { 4163e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 42f3fbd535SBarry Smith ierr = MatRestrict(mglevels[i]->restrct,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr); 4363e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 444b9ad928SBarry Smith } 454b9ad928SBarry Smith 464b9ad928SBarry Smith /* work our way up through the levels */ 47f3fbd535SBarry Smith ierr = VecSet(mglevels[0]->x,0.0);CHKERRQ(ierr); 484b9ad928SBarry Smith for (i=0; i<l-1; i++) { 4963e6d426SJed Brown if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 50f3fbd535SBarry Smith ierr = KSPSolve(mglevels[i]->smoothd,mglevels[i]->b,mglevels[i]->x);CHKERRQ(ierr); 51*c0decd05SBarry Smith ierr = KSPCheckSolve(mglevels[i]->smoothd,pc,mglevels[i]->x);CHKERRQ(ierr); 5263e6d426SJed Brown if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 5363e6d426SJed Brown if (mglevels[i+1]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i+1]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 54f3fbd535SBarry Smith ierr = MatInterpolate(mglevels[i+1]->interpolate,mglevels[i]->x,mglevels[i+1]->x);CHKERRQ(ierr); 5563e6d426SJed Brown if (mglevels[i+1]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i+1]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 564b9ad928SBarry Smith } 5763e6d426SJed Brown if (mglevels[l-1]->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels[l-1]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 58f3fbd535SBarry Smith ierr = KSPSolve(mglevels[l-1]->smoothd,mglevels[l-1]->b,mglevels[l-1]->x);CHKERRQ(ierr); 59*c0decd05SBarry Smith ierr = KSPCheckSolve(mglevels[l-1]->smoothd,pc,mglevels[l-1]->x);CHKERRQ(ierr); 6063e6d426SJed Brown if (mglevels[l-1]->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels[l-1]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 614b9ad928SBarry Smith PetscFunctionReturn(0); 624b9ad928SBarry Smith } 634b9ad928SBarry Smith 644b9ad928SBarry Smith 65