1*dba47a55SKris Buschelman #define PETSCKSP_DLL 2*dba47a55SKris Buschelman 34b9ad928SBarry Smith /* 44b9ad928SBarry Smith Additive Multigrid V Cycle routine 54b9ad928SBarry Smith */ 64b9ad928SBarry Smith #include "src/ksp/pc/impls/mg/mgimpl.h" 74b9ad928SBarry Smith 84b9ad928SBarry Smith /* 94b9ad928SBarry Smith MGACycle_Private - Given an MG structure created with MGCreate() runs 104b9ad928SBarry Smith one cycle down through the levels and back up. Applys 114b9ad928SBarry Smith the smoothers in an additive manner. 124b9ad928SBarry Smith 134b9ad928SBarry Smith Iput Parameters: 144b9ad928SBarry Smith . mg - structure created with MGCreate(). 154b9ad928SBarry Smith 164b9ad928SBarry Smith */ 174b9ad928SBarry Smith #undef __FUNCT__ 184b9ad928SBarry Smith #define __FUNCT__ "MGACycle_Private" 19dfbe8321SBarry Smith PetscErrorCode MGACycle_Private(MG *mg) 204b9ad928SBarry Smith { 216849ba73SBarry Smith PetscErrorCode ierr; 2279416396SBarry Smith PetscInt i,l = mg[0]->levels; 234b9ad928SBarry Smith PetscScalar zero = 0.0; 244b9ad928SBarry Smith 254b9ad928SBarry Smith PetscFunctionBegin; 264b9ad928SBarry Smith /* compute RHS on each level */ 274b9ad928SBarry Smith for (i=l-1; i>0; i--) { 284b9ad928SBarry Smith ierr = MatRestrict(mg[i]->restrct,mg[i]->b,mg[i-1]->b);CHKERRQ(ierr); 294b9ad928SBarry Smith } 304b9ad928SBarry Smith /* solve seperately on each level */ 314b9ad928SBarry Smith for (i=0; i<l; i++) { 324b9ad928SBarry Smith ierr = VecSet(&zero,mg[i]->x);CHKERRQ(ierr); 334b9ad928SBarry Smith if (mg[i]->eventsolve) {ierr = PetscLogEventBegin(mg[i]->eventsolve,0,0,0,0);CHKERRQ(ierr);} 3423ce1328SBarry Smith ierr = KSPSolve(mg[i]->smoothd,mg[i]->b,mg[i]->x);CHKERRQ(ierr); 354b9ad928SBarry Smith if (mg[i]->eventsolve) {ierr = PetscLogEventEnd(mg[i]->eventsolve,0,0,0,0);CHKERRQ(ierr);} 364b9ad928SBarry Smith } 374b9ad928SBarry Smith for (i=1; i<l; i++) { 384b9ad928SBarry Smith ierr = MatInterpolateAdd(mg[i]->interpolate,mg[i-1]->x,mg[i]->x,mg[i]->x);CHKERRQ(ierr); 394b9ad928SBarry Smith } 404b9ad928SBarry Smith PetscFunctionReturn(0); 414b9ad928SBarry Smith } 42