1dba47a55SKris Buschelman 24b9ad928SBarry Smith /* 34b9ad928SBarry Smith Additive Multigrid V Cycle routine 44b9ad928SBarry Smith */ 5af0996ceSBarry Smith #include <petsc/private/pcmgimpl.h> 64b9ad928SBarry Smith 7*fcb023d4SJed Brown PetscErrorCode PCMGACycle_Private(PC pc,PC_MG_Levels **mglevels,PetscBool transpose) 84b9ad928SBarry Smith { 96849ba73SBarry Smith PetscErrorCode ierr; 10f3fbd535SBarry Smith PetscInt i,l = mglevels[0]->levels; 114b9ad928SBarry Smith 124b9ad928SBarry Smith PetscFunctionBegin; 134b9ad928SBarry Smith /* compute RHS on each level */ 144b9ad928SBarry Smith for (i=l-1; i>0; i--) { 1563e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 16*fcb023d4SJed Brown if (!transpose) { 17f3fbd535SBarry Smith ierr = MatRestrict(mglevels[i]->restrct,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr); 18*fcb023d4SJed Brown } else { 19*fcb023d4SJed Brown ierr = MatRestrict(mglevels[i]->interpolate,mglevels[i]->b,mglevels[i-1]->b);CHKERRQ(ierr); 20*fcb023d4SJed Brown } 2163e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 224b9ad928SBarry Smith } 23a8c7a070SBarry Smith /* solve separately on each level */ 244b9ad928SBarry Smith for (i=0; i<l; i++) { 25f3fbd535SBarry Smith ierr = VecSet(mglevels[i]->x,0.0);CHKERRQ(ierr); 2663e6d426SJed Brown if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventBegin(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 27*fcb023d4SJed Brown if (!transpose) { 28f3fbd535SBarry Smith ierr = KSPSolve(mglevels[i]->smoothd,mglevels[i]->b,mglevels[i]->x);CHKERRQ(ierr); 29c0decd05SBarry Smith ierr = KSPCheckSolve(mglevels[i]->smoothd,pc,mglevels[i]->x);CHKERRQ(ierr); 30*fcb023d4SJed Brown } else { 31*fcb023d4SJed Brown ierr = KSPSolveTranspose(mglevels[i]->smoothu,mglevels[i]->b,mglevels[i]->x);CHKERRQ(ierr); 32*fcb023d4SJed Brown ierr = KSPCheckSolve(mglevels[i]->smoothu,pc,mglevels[i]->x);CHKERRQ(ierr); 33*fcb023d4SJed Brown } 3463e6d426SJed Brown if (mglevels[i]->eventsmoothsolve) {ierr = PetscLogEventEnd(mglevels[i]->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 354b9ad928SBarry Smith } 364b9ad928SBarry Smith for (i=1; i<l; i++) { 3763e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventBegin(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 38*fcb023d4SJed Brown if (!transpose) { 39f3fbd535SBarry Smith ierr = MatInterpolateAdd(mglevels[i]->interpolate,mglevels[i-1]->x,mglevels[i]->x,mglevels[i]->x);CHKERRQ(ierr); 40*fcb023d4SJed Brown } else { 41*fcb023d4SJed Brown ierr = MatInterpolateAdd(mglevels[i]->restrct,mglevels[i-1]->x,mglevels[i]->x,mglevels[i]->x);CHKERRQ(ierr); 42*fcb023d4SJed Brown } 4363e6d426SJed Brown if (mglevels[i]->eventinterprestrict) {ierr = PetscLogEventEnd(mglevels[i]->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 444b9ad928SBarry Smith } 454b9ad928SBarry Smith PetscFunctionReturn(0); 464b9ad928SBarry Smith } 47