/* Performs the FAS coarse correction as: fine problem: F(x) = b coarse problem: F^c(x^c) = b^c b^c = F^c(Rx) - R(F(x) - b) */ PetscErrorCode SNESFASCoarseCorrection(SNES snes, Vec X, Vec F, Vec X_new) { PetscErrorCode ierr; Vec X_c, Xo_c, F_c, B_c; SNESConvergedReason reason; SNES next; Mat restrct, interpolate; SNES_FAS *fasc; PetscFunctionBegin; ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); if (next) { fasc = (SNES_FAS*)next->data; ierr = SNESFASCycleGetRestriction(snes, &restrct);CHKERRQ(ierr); ierr = SNESFASCycleGetInterpolation(snes, &interpolate);CHKERRQ(ierr); X_c = next->vec_sol; Xo_c = next->work[0]; F_c = next->vec_func; B_c = next->vec_rhs; if (fasc->eventinterprestrict) {ierr = PetscLogEventBegin(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} ierr = SNESFASRestrict(snes,X,Xo_c);CHKERRQ(ierr); /* restrict the defect: R(F(x) - b) */ ierr = MatRestrict(restrct, F, B_c);CHKERRQ(ierr); if (fasc->eventinterprestrict) {ierr = PetscLogEventEnd(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} /* F_c = F^c(Rx) - R(F(x) - b) since the second term was sitting in next->vec_rhs */ ierr = SNESComputeFunction(next, Xo_c, F_c);CHKERRQ(ierr); if (fasc->eventresidual) {ierr = PetscLogEventEnd(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} /* solve the coarse problem corresponding to F^c(x^c) = b^c = F^c(Rx) - R(F(x) - b) */ ierr = VecCopy(B_c, X_c);CHKERRQ(ierr); ierr = VecCopy(F_c, B_c);CHKERRQ(ierr); ierr = VecCopy(X_c, F_c);CHKERRQ(ierr); /* set initial guess of the coarse problem to the projected fine solution */ ierr = VecCopy(Xo_c, X_c);CHKERRQ(ierr); /* recurse to the next level */ ierr = SNESSetInitialFunction(next, F_c);CHKERRQ(ierr); ierr = SNESSolve(next, B_c, X_c);CHKERRQ(ierr); ierr = SNESGetConvergedReason(next,&reason);CHKERRQ(ierr); if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { snes->reason = SNES_DIVERGED_INNER; PetscFunctionReturn(0); } /* correct as x <- x + I(x^c - Rx)*/ ierr = VecAXPY(X_c, -1.0, Xo_c);CHKERRQ(ierr); if (fasc->eventinterprestrict) {ierr = PetscLogEventBegin(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} ierr = MatInterpolateAdd(interpolate, X_c, X, X_new);CHKERRQ(ierr); if (fasc->eventinterprestrict) {ierr = PetscLogEventEnd(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} } PetscFunctionReturn(0); }
PetscErrorCode SNESFASCycleSetupPhase_Full(SNES snes) { SNES next; SNES_FAS *fas = (SNES_FAS*)snes->data; PetscBool isFine; PetscErrorCode ierr; PetscFunctionBegin; /* pre-smooth -- just update using the pre-smoother */ ierr = SNESFASCycleIsFine(snes,&isFine);CHKERRQ(ierr); ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); fas->full_stage = 0; if (next) {ierr = SNESFASCycleSetupPhase_Full(next);CHKERRQ(ierr);} PetscFunctionReturn(0); }
PetscErrorCode SNESFASCycle_Full(SNES snes, Vec X) { PetscErrorCode ierr; Vec F,B; SNES_FAS *fas = (SNES_FAS*)snes->data; PetscBool isFine; SNES next; PetscFunctionBegin; F = snes->vec_func; B = snes->vec_rhs; ierr = SNESFASCycleIsFine(snes,&isFine);CHKERRQ(ierr); ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); if (isFine) { ierr = SNESFASCycleSetupPhase_Full(snes);CHKERRQ(ierr); } if (fas->full_stage == 0) { /* downsweep */ if (next) { if (fas->level != 1) next->max_its += 1; if (fas->full_downsweep||isFine) {ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr);} ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); if (fas->level != 1) next->max_its -= 1; } else { /* The smoother on the coarse level is the coarse solver */ ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); } fas->full_stage = 1; } else if (fas->full_stage == 1) { if (snes->iter == 0) {ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr);} if (next) { ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); } } /* final v-cycle */ if (isFine) { if (next) { ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); } } PetscFunctionReturn(0); }
PetscErrorCode SNESFASCycle_Kaskade(SNES snes, Vec X) { PetscErrorCode ierr; Vec F,B; SNES next; PetscFunctionBegin; F = snes->vec_func; B = snes->vec_rhs; ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); if (next) { ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); } else { ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); } PetscFunctionReturn(0); }
/* Defines the FAS cycle as: fine problem: F(x) = b coarse problem: F^c(x) = b^c b^c = F^c(Rx) - R(F(x) - b) correction: x = x + I(x^c - Rx) */ PetscErrorCode SNESFASCycle_Multiplicative(SNES snes, Vec X) { PetscErrorCode ierr; Vec F,B; SNES next; PetscFunctionBegin; F = snes->vec_func; B = snes->vec_rhs; /* pre-smooth -- just update using the pre-smoother */ ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); if (next) { ierr = SNESFASCoarseCorrection(snes, X, F, X);CHKERRQ(ierr); ierr = SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); } PetscFunctionReturn(0); }
/* The additive cycle looks like: xhat = x xhat = dS(x, b) x = coarsecorrection(xhat, b_d) x = x + nu*(xhat - x); (optional) x = uS(x, b) With the coarse RHS (defect correction) as below. */ PetscErrorCode SNESFASCycle_Additive(SNES snes, Vec X) { Vec F, B, Xhat; Vec X_c, Xo_c, F_c, B_c; PetscErrorCode ierr; SNESConvergedReason reason; PetscReal xnorm, fnorm, ynorm; PetscBool lssuccess; SNES next; Mat restrct, interpolate; SNES_FAS *fas = (SNES_FAS*)snes->data,*fasc; PetscFunctionBegin; ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); F = snes->vec_func; B = snes->vec_rhs; Xhat = snes->work[1]; ierr = VecCopy(X, Xhat);CHKERRQ(ierr); /* recurse first */ if (next) { fasc = (SNES_FAS*)next->data; ierr = SNESFASCycleGetRestriction(snes, &restrct);CHKERRQ(ierr); ierr = SNESFASCycleGetInterpolation(snes, &interpolate);CHKERRQ(ierr); if (fas->eventresidual) {ierr = PetscLogEventBegin(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} ierr = SNESComputeFunction(snes, Xhat, F);CHKERRQ(ierr); if (fas->eventresidual) {ierr = PetscLogEventEnd(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); X_c = next->vec_sol; Xo_c = next->work[0]; F_c = next->vec_func; B_c = next->vec_rhs; ierr = SNESFASRestrict(snes,Xhat,Xo_c);CHKERRQ(ierr); /* restrict the defect */ ierr = MatRestrict(restrct, F, B_c);CHKERRQ(ierr); /* solve the coarse problem corresponding to F^c(x^c) = b^c = Rb + F^c(Rx) - RF(x) */ if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} ierr = SNESComputeFunction(next, Xo_c, F_c);CHKERRQ(ierr); if (fasc->eventresidual) {ierr = PetscLogEventEnd(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} ierr = VecCopy(B_c, X_c);CHKERRQ(ierr); ierr = VecCopy(F_c, B_c);CHKERRQ(ierr); ierr = VecCopy(X_c, F_c);CHKERRQ(ierr); /* set initial guess of the coarse problem to the projected fine solution */ ierr = VecCopy(Xo_c, X_c);CHKERRQ(ierr); /* recurse */ ierr = SNESSetInitialFunction(next, F_c);CHKERRQ(ierr); ierr = SNESSolve(next, B_c, X_c);CHKERRQ(ierr); /* smooth on this level */ ierr = SNESFASDownSmooth_Private(snes, B, X, F, &fnorm);CHKERRQ(ierr); ierr = SNESGetConvergedReason(next,&reason);CHKERRQ(ierr); if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { snes->reason = SNES_DIVERGED_INNER; PetscFunctionReturn(0); } /* correct as x <- x + I(x^c - Rx)*/ ierr = VecAYPX(X_c, -1.0, Xo_c);CHKERRQ(ierr); ierr = MatInterpolate(interpolate, X_c, Xhat);CHKERRQ(ierr); /* additive correction of the coarse direction*/ ierr = SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Xhat);CHKERRQ(ierr); ierr = SNESLineSearchGetSuccess(snes->linesearch, &lssuccess);CHKERRQ(ierr); if (!lssuccess) { if (++snes->numFailures >= snes->maxFailures) { snes->reason = SNES_DIVERGED_LINE_SEARCH; PetscFunctionReturn(0); } } ierr = SNESLineSearchGetNorms(snes->linesearch, &xnorm, &snes->norm, &ynorm);CHKERRQ(ierr); } else { ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); } PetscFunctionReturn(0); }
PetscErrorCode SNESSetFromOptions_FAS(PetscOptions *PetscOptionsObject,SNES snes) { SNES_FAS *fas = (SNES_FAS*) snes->data; PetscInt levels = 1; PetscBool flg = PETSC_FALSE, upflg = PETSC_FALSE, downflg = PETSC_FALSE, monflg = PETSC_FALSE, galerkinflg = PETSC_FALSE,continuationflg = PETSC_FALSE; PetscErrorCode ierr; char monfilename[PETSC_MAX_PATH_LEN]; SNESFASType fastype; const char *optionsprefix; SNESLineSearch linesearch; PetscInt m, n_up, n_down; SNES next; PetscBool isFine; PetscFunctionBegin; ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); ierr = PetscOptionsHead(PetscOptionsObject,"SNESFAS Options-----------------------------------");CHKERRQ(ierr); /* number of levels -- only process most options on the finest level */ if (isFine) { ierr = PetscOptionsInt("-snes_fas_levels", "Number of Levels", "SNESFASSetLevels", levels, &levels, &flg);CHKERRQ(ierr); if (!flg && snes->dm) { ierr = DMGetRefineLevel(snes->dm,&levels);CHKERRQ(ierr); levels++; fas->usedmfornumberoflevels = PETSC_TRUE; } ierr = SNESFASSetLevels(snes, levels, NULL);CHKERRQ(ierr); fastype = fas->fastype; ierr = PetscOptionsEnum("-snes_fas_type","FAS correction type","SNESFASSetType",SNESFASTypes,(PetscEnum)fastype,(PetscEnum*)&fastype,&flg);CHKERRQ(ierr); if (flg) { ierr = SNESFASSetType(snes, fastype);CHKERRQ(ierr); } ierr = SNESGetOptionsPrefix(snes, &optionsprefix);CHKERRQ(ierr); ierr = PetscOptionsInt("-snes_fas_cycles","Number of cycles","SNESFASSetCycles",fas->n_cycles,&m,&flg);CHKERRQ(ierr); if (flg) { ierr = SNESFASSetCycles(snes, m);CHKERRQ(ierr); } ierr = PetscOptionsBool("-snes_fas_continuation","Corrected grid-sequence continuation","SNESFASSetContinuation",fas->continuation,&continuationflg,&flg);CHKERRQ(ierr); if (flg) { ierr = SNESFASSetContinuation(snes,continuationflg);CHKERRQ(ierr); } ierr = PetscOptionsBool("-snes_fas_galerkin", "Form coarse problems with Galerkin","SNESFASSetGalerkin",fas->galerkin,&galerkinflg,&flg);CHKERRQ(ierr); if (flg) { ierr = SNESFASSetGalerkin(snes, galerkinflg);CHKERRQ(ierr); } if (fas->fastype == SNES_FAS_FULL) { ierr = PetscOptionsBool("-snes_fas_full_downsweep","Smooth on the initial upsweep for full FAS cycles","SNESFASFullSetDownSweep",fas->full_downsweep,&fas->full_downsweep,&flg);CHKERRQ(ierr); if (flg) {SNESFASFullSetDownSweep(snes,fas->full_downsweep);CHKERRQ(ierr);} } ierr = PetscOptionsInt("-snes_fas_smoothup","Number of post-smoothing steps","SNESFASSetNumberSmoothUp",fas->max_up_it,&n_up,&upflg);CHKERRQ(ierr); ierr = PetscOptionsInt("-snes_fas_smoothdown","Number of pre-smoothing steps","SNESFASSetNumberSmoothDown",fas->max_down_it,&n_down,&downflg);CHKERRQ(ierr); ierr = PetscOptionsString("-snes_fas_monitor","Monitor FAS progress","SNESFASSetMonitor","stdout",monfilename,PETSC_MAX_PATH_LEN,&monflg);CHKERRQ(ierr); if (monflg) ierr = SNESFASSetMonitor(snes, PETSC_TRUE);CHKERRQ(ierr); flg = PETSC_FALSE; monflg = PETSC_TRUE; ierr = PetscOptionsBool("-snes_fas_log","Log times for each FAS level","SNESFASSetLog",monflg,&monflg,&flg);CHKERRQ(ierr); if (flg) {ierr = SNESFASSetLog(snes,monflg);CHKERRQ(ierr);} } ierr = PetscOptionsTail();CHKERRQ(ierr); /* setup from the determined types if there is no pointwise procedure or smoother defined */ if (upflg) { ierr = SNESFASSetNumberSmoothUp(snes,n_up);CHKERRQ(ierr); } if (downflg) { ierr = SNESFASSetNumberSmoothDown(snes,n_down);CHKERRQ(ierr); } /* set up the default line search for coarse grid corrections */ if (fas->fastype == SNES_FAS_ADDITIVE) { if (!snes->linesearch) { ierr = SNESGetLineSearch(snes, &linesearch);CHKERRQ(ierr); ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHL2);CHKERRQ(ierr); } } ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); /* recursive option setting for the smoothers */ if (next) {ierr = SNESSetFromOptions(next);CHKERRQ(ierr);} PetscFunctionReturn(0); }
PetscErrorCode SNESSetUp_FAS(SNES snes) { SNES_FAS *fas = (SNES_FAS*) snes->data; PetscErrorCode ierr; PetscInt dm_levels; Vec vec_sol, vec_func, vec_sol_update, vec_rhs; /* preserve these if they're set through the reset */ SNES next; PetscBool isFine; SNESLineSearch linesearch; SNESLineSearch slinesearch; void *lsprectx,*lspostctx; PetscErrorCode (*precheck)(SNESLineSearch,Vec,Vec,PetscBool*,void*); PetscErrorCode (*postcheck)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*); PetscFunctionBegin; ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); if (fas->usedmfornumberoflevels && isFine) { ierr = DMGetRefineLevel(snes->dm,&dm_levels);CHKERRQ(ierr); dm_levels++; if (dm_levels > fas->levels) { /* we don't want the solution and func vectors to be destroyed in the SNESReset when it's called in SNESFASSetLevels_FAS*/ vec_sol = snes->vec_sol; vec_func = snes->vec_func; vec_sol_update = snes->vec_sol_update; vec_rhs = snes->vec_rhs; snes->vec_sol = NULL; snes->vec_func = NULL; snes->vec_sol_update = NULL; snes->vec_rhs = NULL; /* reset the number of levels */ ierr = SNESFASSetLevels(snes,dm_levels,NULL);CHKERRQ(ierr); ierr = SNESSetFromOptions(snes);CHKERRQ(ierr); snes->vec_sol = vec_sol; snes->vec_func = vec_func; snes->vec_rhs = vec_rhs; snes->vec_sol_update = vec_sol_update; } } ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); if (!isFine) snes->gridsequence = 0; /* no grid sequencing inside the multigrid hierarchy! */ ierr = SNESSetWorkVecs(snes, 2);CHKERRQ(ierr); /* work vectors used for intergrid transfers */ /* set up the smoothers if they haven't already been set up */ if (!fas->smoothd) { ierr = SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);CHKERRQ(ierr); } if (snes->dm) { /* set the smoother DMs properly */ if (fas->smoothu) ierr = SNESSetDM(fas->smoothu, snes->dm);CHKERRQ(ierr); ierr = SNESSetDM(fas->smoothd, snes->dm);CHKERRQ(ierr); /* construct EVERYTHING from the DM -- including the progressive set of smoothers */ if (next) { /* for now -- assume the DM and the evaluation functions have been set externally */ if (!next->dm) { ierr = DMCoarsen(snes->dm, PetscObjectComm((PetscObject)next), &next->dm);CHKERRQ(ierr); ierr = SNESSetDM(next, next->dm);CHKERRQ(ierr); } /* set the interpolation and restriction from the DM */ if (!fas->interpolate) { ierr = DMCreateInterpolation(next->dm, snes->dm, &fas->interpolate, &fas->rscale);CHKERRQ(ierr); if (!fas->restrct) { ierr = PetscObjectReference((PetscObject)fas->interpolate);CHKERRQ(ierr); fas->restrct = fas->interpolate; } } /* set the injection from the DM */ if (!fas->inject) { ierr = DMCreateInjection(next->dm, snes->dm, &fas->inject);CHKERRQ(ierr); } } } /*pass the smoother, function, and jacobian up to the next level if it's not user set already */ if (fas->galerkin) { if (next) { ierr = SNESSetFunction(next, NULL, SNESFASGalerkinDefaultFunction, next);CHKERRQ(ierr); } if (fas->smoothd && fas->level != fas->levels - 1) { ierr = SNESSetFunction(fas->smoothd, NULL, SNESFASGalerkinDefaultFunction, snes);CHKERRQ(ierr); } if (fas->smoothu && fas->level != fas->levels - 1) { ierr = SNESSetFunction(fas->smoothu, NULL, SNESFASGalerkinDefaultFunction, snes);CHKERRQ(ierr); } } /* sets the down (pre) smoother's default norm and sets it from options */ if (fas->smoothd) { if (fas->level == 0 && fas->levels != 1) { ierr = SNESSetNormSchedule(fas->smoothd, SNES_NORM_NONE);CHKERRQ(ierr); } else { ierr = SNESSetNormSchedule(fas->smoothd, SNES_NORM_FINAL_ONLY);CHKERRQ(ierr); } ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)fas->smoothd);CHKERRQ(ierr); ierr = SNESSetFromOptions(fas->smoothd);CHKERRQ(ierr); ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); ierr = SNESGetLineSearch(fas->smoothd,&slinesearch);CHKERRQ(ierr); ierr = SNESLineSearchGetPreCheck(linesearch,&precheck,&lsprectx);CHKERRQ(ierr); ierr = SNESLineSearchGetPostCheck(linesearch,&postcheck,&lspostctx);CHKERRQ(ierr); ierr = SNESLineSearchSetPreCheck(slinesearch,precheck,lsprectx);CHKERRQ(ierr); ierr = SNESLineSearchSetPostCheck(slinesearch,postcheck,lspostctx);CHKERRQ(ierr); ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)slinesearch);CHKERRQ(ierr); fas->smoothd->vec_sol = snes->vec_sol; ierr = PetscObjectReference((PetscObject)snes->vec_sol);CHKERRQ(ierr); fas->smoothd->vec_sol_update = snes->vec_sol_update; ierr = PetscObjectReference((PetscObject)snes->vec_sol_update);CHKERRQ(ierr); fas->smoothd->vec_func = snes->vec_func; ierr = PetscObjectReference((PetscObject)snes->vec_func);CHKERRQ(ierr); if (fas->eventsmoothsetup) {ierr = PetscLogEventBegin(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} ierr = SNESSetUp(fas->smoothd);CHKERRQ(ierr); if (fas->eventsmoothsetup) {ierr = PetscLogEventEnd(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} } /* sets the up (post) smoother's default norm and sets it from options */ if (fas->smoothu) { if (fas->level != fas->levels - 1) { ierr = SNESSetNormSchedule(fas->smoothu, SNES_NORM_NONE);CHKERRQ(ierr); } else { ierr = SNESSetNormSchedule(fas->smoothu, SNES_NORM_FINAL_ONLY);CHKERRQ(ierr); } ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)fas->smoothu);CHKERRQ(ierr); ierr = SNESSetFromOptions(fas->smoothu);CHKERRQ(ierr); ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); ierr = SNESGetLineSearch(fas->smoothu,&slinesearch);CHKERRQ(ierr); ierr = SNESLineSearchGetPreCheck(linesearch,&precheck,&lsprectx);CHKERRQ(ierr); ierr = SNESLineSearchGetPostCheck(linesearch,&postcheck,&lspostctx);CHKERRQ(ierr); ierr = SNESLineSearchSetPreCheck(slinesearch,precheck,lsprectx);CHKERRQ(ierr); ierr = SNESLineSearchSetPostCheck(slinesearch,postcheck,lspostctx);CHKERRQ(ierr); ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)slinesearch);CHKERRQ(ierr); fas->smoothu->vec_sol = snes->vec_sol; ierr = PetscObjectReference((PetscObject)snes->vec_sol);CHKERRQ(ierr); fas->smoothu->vec_sol_update = snes->vec_sol_update; ierr = PetscObjectReference((PetscObject)snes->vec_sol_update);CHKERRQ(ierr); fas->smoothu->vec_func = snes->vec_func; ierr = PetscObjectReference((PetscObject)snes->vec_func);CHKERRQ(ierr); if (fas->eventsmoothsetup) {ierr = PetscLogEventBegin(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} ierr = SNESSetUp(fas->smoothu);CHKERRQ(ierr); if (fas->eventsmoothsetup) {ierr = PetscLogEventEnd(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} } if (next) { /* gotta set up the solution vector for this to work */ if (!next->vec_sol) {ierr = SNESFASCreateCoarseVec(snes,&next->vec_sol);CHKERRQ(ierr);} if (!next->vec_rhs) {ierr = SNESFASCreateCoarseVec(snes,&next->vec_rhs);CHKERRQ(ierr);} ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)next);CHKERRQ(ierr); ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); ierr = SNESGetLineSearch(fas->next,&slinesearch);CHKERRQ(ierr); ierr = SNESLineSearchGetPreCheck(linesearch,&precheck,&lsprectx);CHKERRQ(ierr); ierr = SNESLineSearchGetPostCheck(linesearch,&postcheck,&lspostctx);CHKERRQ(ierr); ierr = SNESLineSearchSetPreCheck(slinesearch,precheck,lsprectx);CHKERRQ(ierr); ierr = SNESLineSearchSetPostCheck(slinesearch,postcheck,lspostctx);CHKERRQ(ierr); ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)slinesearch);CHKERRQ(ierr); ierr = SNESSetUp(next);CHKERRQ(ierr); } /* setup FAS work vectors */ if (fas->galerkin) { ierr = VecDuplicate(snes->vec_sol, &fas->Xg);CHKERRQ(ierr); ierr = VecDuplicate(snes->vec_sol, &fas->Fg);CHKERRQ(ierr); } PetscFunctionReturn(0); }
PETSC_EXTERN void PETSC_STDCALL snesfascyclegetcorrection_(SNES snes,SNES *correction, int *__ierr ){ *__ierr = SNESFASCycleGetCorrection( (SNES)PetscToPointer((snes) ),correction); }