1 /* Defines the basic SNES object */ 2 #include <../src/snes/impls/fas/fasimpls.h> /*I "petscsnesfas.h" I*/ 3 4 const char *const SNESFASTypes[] = {"MULTIPLICATIVE","ADDITIVE","SNESFASType","SNES_FAS",0}; 5 6 extern PetscErrorCode SNESDestroy_FAS(SNES snes); 7 extern PetscErrorCode SNESSetUp_FAS(SNES snes); 8 extern PetscErrorCode SNESSetFromOptions_FAS(SNES snes); 9 extern PetscErrorCode SNESView_FAS(SNES snes, PetscViewer viewer); 10 extern PetscErrorCode SNESSolve_FAS(SNES snes); 11 extern PetscErrorCode SNESReset_FAS(SNES snes); 12 extern PetscErrorCode SNESFASGalerkinDefaultFunction(SNES, Vec, Vec, void*); 13 extern PetscErrorCode SNESFASCycleCreateSmoother_Private(SNES, SNES*); 14 15 /*MC 16 17 SNESFAS - Full Approximation Scheme nonlinear multigrid solver. 18 19 The nonlinear problem is solved by correction using coarse versions 20 of the nonlinear problem. This problem is perturbed so that a projected 21 solution of the fine problem elicits no correction from the coarse problem. 22 23 Options Database: 24 + -snes_fas_levels - The number of levels 25 . -snes_fas_cycles<1> - The number of cycles -- 1 for V, 2 for W 26 . -snes_fas_type<additive, multiplicative> - Additive or multiplicative cycle 27 . -snes_fas_galerkin<PETSC_FALSE> - Form coarse problems by projection back upon the fine problem 28 . -snes_fas_smoothup<1> - The number of iterations of the post-smoother 29 . -snes_fas_smoothdown<1> - The number of iterations of the pre-smoother 30 . -snes_fas_monitor - Monitor progress of all of the levels 31 . -fas_levels_snes_ - SNES options for all smoothers 32 . -fas_levels_cycle_snes_ - SNES options for all cycles 33 . -fas_levels_i_snes_ - SNES options for the smoothers on level i 34 . -fas_levels_i_cycle_snes_ - SNES options for the cycle on level i 35 - -fas_coarse_snes_ - SNES options for the coarsest smoother 36 37 Notes: 38 The organization of the FAS solver is slightly different from the organization of PCMG 39 As each level has smoother SNES instances(down and potentially up) and a cycle SNES instance. 40 The cycle SNES instance may be used for monitoring convergence on a particular level. 41 42 Level: beginner 43 44 .seealso: PCMG, SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 45 M*/ 46 47 #undef __FUNCT__ 48 #define __FUNCT__ "SNESCreate_FAS" 49 PETSC_EXTERN PetscErrorCode SNESCreate_FAS(SNES snes) 50 { 51 SNES_FAS *fas; 52 PetscErrorCode ierr; 53 54 PetscFunctionBegin; 55 snes->ops->destroy = SNESDestroy_FAS; 56 snes->ops->setup = SNESSetUp_FAS; 57 snes->ops->setfromoptions = SNESSetFromOptions_FAS; 58 snes->ops->view = SNESView_FAS; 59 snes->ops->solve = SNESSolve_FAS; 60 snes->ops->reset = SNESReset_FAS; 61 62 snes->usesksp = PETSC_FALSE; 63 snes->usespc = PETSC_FALSE; 64 65 if (!snes->tolerancesset) { 66 snes->max_funcs = 30000; 67 snes->max_its = 10000; 68 } 69 70 ierr = PetscNewLog(snes, SNES_FAS, &fas);CHKERRQ(ierr); 71 72 snes->data = (void*) fas; 73 fas->level = 0; 74 fas->levels = 1; 75 fas->n_cycles = 1; 76 fas->max_up_it = 1; 77 fas->max_down_it = 1; 78 fas->smoothu = NULL; 79 fas->smoothd = NULL; 80 fas->next = NULL; 81 fas->previous = NULL; 82 fas->fine = snes; 83 fas->interpolate = NULL; 84 fas->restrct = NULL; 85 fas->inject = NULL; 86 fas->monitor = NULL; 87 fas->usedmfornumberoflevels = PETSC_FALSE; 88 fas->fastype = SNES_FAS_MULTIPLICATIVE; 89 90 fas->eventsmoothsetup = 0; 91 fas->eventsmoothsolve = 0; 92 fas->eventresidual = 0; 93 fas->eventinterprestrict = 0; 94 PetscFunctionReturn(0); 95 } 96 97 #undef __FUNCT__ 98 #define __FUNCT__ "SNESReset_FAS" 99 PetscErrorCode SNESReset_FAS(SNES snes) 100 { 101 PetscErrorCode ierr = 0; 102 SNES_FAS * fas = (SNES_FAS*)snes->data; 103 104 PetscFunctionBegin; 105 ierr = SNESDestroy(&fas->smoothu);CHKERRQ(ierr); 106 ierr = SNESDestroy(&fas->smoothd);CHKERRQ(ierr); 107 ierr = MatDestroy(&fas->inject);CHKERRQ(ierr); 108 ierr = MatDestroy(&fas->interpolate);CHKERRQ(ierr); 109 ierr = MatDestroy(&fas->restrct);CHKERRQ(ierr); 110 ierr = VecDestroy(&fas->rscale);CHKERRQ(ierr); 111 if (fas->galerkin) { 112 ierr = VecDestroy(&fas->Xg);CHKERRQ(ierr); 113 ierr = VecDestroy(&fas->Fg);CHKERRQ(ierr); 114 } 115 if (fas->next) {ierr = SNESReset(fas->next);CHKERRQ(ierr);} 116 PetscFunctionReturn(0); 117 } 118 119 #undef __FUNCT__ 120 #define __FUNCT__ "SNESDestroy_FAS" 121 PetscErrorCode SNESDestroy_FAS(SNES snes) 122 { 123 SNES_FAS * fas = (SNES_FAS*)snes->data; 124 PetscErrorCode ierr = 0; 125 126 PetscFunctionBegin; 127 /* recursively resets and then destroys */ 128 ierr = SNESReset(snes);CHKERRQ(ierr); 129 if (fas->next) { 130 ierr = SNESDestroy(&fas->next);CHKERRQ(ierr); 131 } 132 ierr = PetscFree(fas);CHKERRQ(ierr); 133 PetscFunctionReturn(0); 134 } 135 136 #undef __FUNCT__ 137 #define __FUNCT__ "SNESSetUp_FAS" 138 PetscErrorCode SNESSetUp_FAS(SNES snes) 139 { 140 SNES_FAS *fas = (SNES_FAS*) snes->data; 141 PetscErrorCode ierr; 142 VecScatter injscatter; 143 PetscInt dm_levels; 144 Vec vec_sol, vec_func, vec_sol_update, vec_rhs; /* preserve these if they're set through the reset */ 145 SNES next; 146 PetscBool isFine; 147 SNESLineSearch linesearch; 148 SNESLineSearch slinesearch; 149 void *lsprectx,*lspostctx; 150 PetscErrorCode (*precheck)(SNESLineSearch,Vec,Vec,PetscBool*,void*); 151 PetscErrorCode (*postcheck)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*); 152 153 PetscFunctionBegin; 154 ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); 155 if (fas->usedmfornumberoflevels && isFine) { 156 ierr = DMGetRefineLevel(snes->dm,&dm_levels);CHKERRQ(ierr); 157 dm_levels++; 158 if (dm_levels > fas->levels) { 159 /* we don't want the solution and func vectors to be destroyed in the SNESReset when it's called in SNESFASSetLevels_FAS*/ 160 vec_sol = snes->vec_sol; 161 vec_func = snes->vec_func; 162 vec_sol_update = snes->vec_sol_update; 163 vec_rhs = snes->vec_rhs; 164 snes->vec_sol = NULL; 165 snes->vec_func = NULL; 166 snes->vec_sol_update = NULL; 167 snes->vec_rhs = NULL; 168 169 /* reset the number of levels */ 170 ierr = SNESFASSetLevels(snes,dm_levels,NULL);CHKERRQ(ierr); 171 ierr = SNESSetFromOptions(snes);CHKERRQ(ierr); 172 173 snes->vec_sol = vec_sol; 174 snes->vec_func = vec_func; 175 snes->vec_rhs = vec_rhs; 176 snes->vec_sol_update = vec_sol_update; 177 } 178 } 179 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 180 if (!isFine) snes->gridsequence = 0; /* no grid sequencing inside the multigrid hierarchy! */ 181 182 ierr = SNESSetWorkVecs(snes, 2);CHKERRQ(ierr); /* work vectors used for intergrid transfers */ 183 184 /* set up the smoothers if they haven't already been set up */ 185 if (!fas->smoothd) { 186 ierr = SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);CHKERRQ(ierr); 187 } 188 189 if (snes->dm) { 190 /* set the smoother DMs properly */ 191 if (fas->smoothu) ierr = SNESSetDM(fas->smoothu, snes->dm);CHKERRQ(ierr); 192 ierr = SNESSetDM(fas->smoothd, snes->dm);CHKERRQ(ierr); 193 /* construct EVERYTHING from the DM -- including the progressive set of smoothers */ 194 if (next) { 195 /* for now -- assume the DM and the evaluation functions have been set externally */ 196 if (!next->dm) { 197 ierr = DMCoarsen(snes->dm, PetscObjectComm((PetscObject)next), &next->dm);CHKERRQ(ierr); 198 ierr = SNESSetDM(next, next->dm);CHKERRQ(ierr); 199 } 200 /* set the interpolation and restriction from the DM */ 201 if (!fas->interpolate) { 202 ierr = DMCreateInterpolation(next->dm, snes->dm, &fas->interpolate, &fas->rscale);CHKERRQ(ierr); 203 if (!fas->restrct) { 204 ierr = PetscObjectReference((PetscObject)fas->interpolate);CHKERRQ(ierr); 205 fas->restrct = fas->interpolate; 206 } 207 } 208 /* set the injection from the DM */ 209 if (!fas->inject) { 210 ierr = DMCreateInjection(next->dm, snes->dm, &injscatter);CHKERRQ(ierr); 211 ierr = MatCreateScatter(PetscObjectComm((PetscObject)snes), injscatter, &fas->inject);CHKERRQ(ierr); 212 ierr = VecScatterDestroy(&injscatter);CHKERRQ(ierr); 213 } 214 } 215 } 216 /*pass the smoother, function, and jacobian up to the next level if it's not user set already */ 217 if (fas->galerkin) { 218 if (next) { 219 ierr = SNESSetFunction(next, NULL, SNESFASGalerkinDefaultFunction, next);CHKERRQ(ierr); 220 } 221 if (fas->smoothd && fas->level != fas->levels - 1) { 222 ierr = SNESSetFunction(fas->smoothd, NULL, SNESFASGalerkinDefaultFunction, snes);CHKERRQ(ierr); 223 } 224 if (fas->smoothu && fas->level != fas->levels - 1) { 225 ierr = SNESSetFunction(fas->smoothu, NULL, SNESFASGalerkinDefaultFunction, snes);CHKERRQ(ierr); 226 } 227 } 228 229 /* sets the down (pre) smoother's default norm and sets it from options */ 230 if (fas->smoothd) { 231 if (fas->level == 0 && fas->levels != 1) { 232 ierr = SNESSetNormSchedule(fas->smoothd, SNES_NORM_NONE);CHKERRQ(ierr); 233 } else { 234 ierr = SNESSetNormSchedule(fas->smoothd, SNES_NORM_FINAL_ONLY);CHKERRQ(ierr); 235 } 236 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)fas->smoothd);CHKERRQ(ierr); 237 ierr = SNESSetFromOptions(fas->smoothd);CHKERRQ(ierr); 238 ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); 239 ierr = SNESGetLineSearch(fas->smoothd,&slinesearch);CHKERRQ(ierr); 240 ierr = SNESLineSearchGetPreCheck(linesearch,&precheck,&lsprectx);CHKERRQ(ierr); 241 ierr = SNESLineSearchGetPostCheck(linesearch,&postcheck,&lspostctx);CHKERRQ(ierr); 242 ierr = SNESLineSearchSetPreCheck(slinesearch,precheck,lsprectx);CHKERRQ(ierr); 243 ierr = SNESLineSearchSetPostCheck(slinesearch,postcheck,lspostctx);CHKERRQ(ierr); 244 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)slinesearch);CHKERRQ(ierr); 245 246 fas->smoothd->vec_sol = snes->vec_sol; 247 ierr = PetscObjectReference((PetscObject)snes->vec_sol);CHKERRQ(ierr); 248 fas->smoothd->vec_sol_update = snes->vec_sol_update; 249 ierr = PetscObjectReference((PetscObject)snes->vec_sol_update);CHKERRQ(ierr); 250 fas->smoothd->vec_func = snes->vec_func; 251 ierr = PetscObjectReference((PetscObject)snes->vec_func);CHKERRQ(ierr); 252 253 if (fas->eventsmoothsetup) {ierr = PetscLogEventBegin(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} 254 ierr = SNESSetUp(fas->smoothd);CHKERRQ(ierr); 255 if (fas->eventsmoothsetup) {ierr = PetscLogEventEnd(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} 256 } 257 258 /* sets the up (post) smoother's default norm and sets it from options */ 259 if (fas->smoothu) { 260 if (fas->level != fas->levels - 1) { 261 ierr = SNESSetNormSchedule(fas->smoothu, SNES_NORM_NONE);CHKERRQ(ierr); 262 } else { 263 ierr = SNESSetNormSchedule(fas->smoothu, SNES_NORM_FINAL_ONLY);CHKERRQ(ierr); 264 } 265 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)fas->smoothu);CHKERRQ(ierr); 266 ierr = SNESSetFromOptions(fas->smoothu);CHKERRQ(ierr); 267 ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); 268 ierr = SNESGetLineSearch(fas->smoothu,&slinesearch);CHKERRQ(ierr); 269 ierr = SNESLineSearchGetPreCheck(linesearch,&precheck,&lsprectx);CHKERRQ(ierr); 270 ierr = SNESLineSearchGetPostCheck(linesearch,&postcheck,&lspostctx);CHKERRQ(ierr); 271 ierr = SNESLineSearchSetPreCheck(slinesearch,precheck,lsprectx);CHKERRQ(ierr); 272 ierr = SNESLineSearchSetPostCheck(slinesearch,postcheck,lspostctx);CHKERRQ(ierr); 273 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)slinesearch);CHKERRQ(ierr); 274 275 fas->smoothu->vec_sol = snes->vec_sol; 276 ierr = PetscObjectReference((PetscObject)snes->vec_sol);CHKERRQ(ierr); 277 fas->smoothu->vec_sol_update = snes->vec_sol_update; 278 ierr = PetscObjectReference((PetscObject)snes->vec_sol_update);CHKERRQ(ierr); 279 fas->smoothu->vec_func = snes->vec_func; 280 ierr = PetscObjectReference((PetscObject)snes->vec_func);CHKERRQ(ierr); 281 282 if (fas->eventsmoothsetup) {ierr = PetscLogEventBegin(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} 283 ierr = SNESSetUp(fas->smoothu);CHKERRQ(ierr); 284 if (fas->eventsmoothsetup) {ierr = PetscLogEventEnd(fas->eventsmoothsetup,0,0,0,0);CHKERRQ(ierr);} 285 286 } 287 288 if (next) { 289 /* gotta set up the solution vector for this to work */ 290 if (!next->vec_sol) {ierr = SNESFASCreateCoarseVec(snes,&next->vec_sol);CHKERRQ(ierr);} 291 if (!next->vec_rhs) {ierr = SNESFASCreateCoarseVec(snes,&next->vec_rhs);CHKERRQ(ierr);} 292 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)next);CHKERRQ(ierr); 293 ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); 294 ierr = SNESGetLineSearch(fas->next,&slinesearch);CHKERRQ(ierr); 295 ierr = SNESLineSearchGetPreCheck(linesearch,&precheck,&lsprectx);CHKERRQ(ierr); 296 ierr = SNESLineSearchGetPostCheck(linesearch,&postcheck,&lspostctx);CHKERRQ(ierr); 297 ierr = SNESLineSearchSetPreCheck(slinesearch,precheck,lsprectx);CHKERRQ(ierr); 298 ierr = SNESLineSearchSetPostCheck(slinesearch,postcheck,lspostctx);CHKERRQ(ierr); 299 ierr = PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)slinesearch);CHKERRQ(ierr); 300 ierr = SNESSetUp(next);CHKERRQ(ierr); 301 } 302 /* setup FAS work vectors */ 303 if (fas->galerkin) { 304 ierr = VecDuplicate(snes->vec_sol, &fas->Xg);CHKERRQ(ierr); 305 ierr = VecDuplicate(snes->vec_sol, &fas->Fg);CHKERRQ(ierr); 306 } 307 PetscFunctionReturn(0); 308 } 309 310 #undef __FUNCT__ 311 #define __FUNCT__ "SNESSetFromOptions_FAS" 312 PetscErrorCode SNESSetFromOptions_FAS(SNES snes) 313 { 314 SNES_FAS *fas = (SNES_FAS*) snes->data; 315 PetscInt levels = 1; 316 PetscBool flg = PETSC_FALSE, upflg = PETSC_FALSE, downflg = PETSC_FALSE, monflg = PETSC_FALSE, galerkinflg = PETSC_FALSE; 317 PetscErrorCode ierr; 318 char monfilename[PETSC_MAX_PATH_LEN]; 319 SNESFASType fastype; 320 const char *optionsprefix; 321 SNESLineSearch linesearch; 322 PetscInt m, n_up, n_down; 323 SNES next; 324 PetscBool isFine; 325 326 PetscFunctionBegin; 327 ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); 328 ierr = PetscOptionsHead("SNESFAS Options-----------------------------------");CHKERRQ(ierr); 329 330 /* number of levels -- only process most options on the finest level */ 331 if (isFine) { 332 ierr = PetscOptionsInt("-snes_fas_levels", "Number of Levels", "SNESFASSetLevels", levels, &levels, &flg);CHKERRQ(ierr); 333 if (!flg && snes->dm) { 334 ierr = DMGetRefineLevel(snes->dm,&levels);CHKERRQ(ierr); 335 levels++; 336 fas->usedmfornumberoflevels = PETSC_TRUE; 337 } 338 ierr = SNESFASSetLevels(snes, levels, NULL);CHKERRQ(ierr); 339 fastype = fas->fastype; 340 ierr = PetscOptionsEnum("-snes_fas_type","FAS correction type","SNESFASSetType",SNESFASTypes,(PetscEnum)fastype,(PetscEnum*)&fastype,&flg);CHKERRQ(ierr); 341 if (flg) { 342 ierr = SNESFASSetType(snes, fastype);CHKERRQ(ierr); 343 } 344 345 ierr = SNESGetOptionsPrefix(snes, &optionsprefix);CHKERRQ(ierr); 346 ierr = PetscOptionsInt("-snes_fas_cycles","Number of cycles","SNESFASSetCycles",fas->n_cycles,&m,&flg);CHKERRQ(ierr); 347 if (flg) { 348 ierr = SNESFASSetCycles(snes, m);CHKERRQ(ierr); 349 } 350 351 ierr = PetscOptionsBool("-snes_fas_galerkin", "Form coarse problems with Galerkin","SNESFASSetGalerkin",fas->galerkin,&galerkinflg,&flg);CHKERRQ(ierr); 352 if (flg) { 353 ierr = SNESFASSetGalerkin(snes, galerkinflg);CHKERRQ(ierr); 354 } 355 356 ierr = PetscOptionsInt("-snes_fas_smoothup","Number of post-smoothing steps","SNESFASSetNumberSmoothUp",fas->max_up_it,&n_up,&upflg);CHKERRQ(ierr); 357 358 ierr = PetscOptionsInt("-snes_fas_smoothdown","Number of pre-smoothing steps","SNESFASSetNumberSmoothDown",fas->max_down_it,&n_down,&downflg);CHKERRQ(ierr); 359 360 ierr = PetscOptionsString("-snes_fas_monitor","Monitor FAS progress","SNESFASSetMonitor","stdout",monfilename,PETSC_MAX_PATH_LEN,&monflg);CHKERRQ(ierr); 361 if (monflg) ierr = SNESFASSetMonitor(snes, PETSC_TRUE);CHKERRQ(ierr); 362 363 flg = PETSC_FALSE; 364 monflg = PETSC_TRUE; 365 ierr = PetscOptionsBool("-snes_fas_log","Log times for each FAS level","SNESFASSetLog",monflg,&monflg,&flg);CHKERRQ(ierr); 366 if (flg) {ierr = SNESFASSetLog(snes,monflg);CHKERRQ(ierr);} 367 } 368 369 ierr = PetscOptionsTail();CHKERRQ(ierr); 370 /* setup from the determined types if there is no pointwise procedure or smoother defined */ 371 if (upflg) { 372 ierr = SNESFASSetNumberSmoothUp(snes,n_up);CHKERRQ(ierr); 373 } 374 if (downflg) { 375 ierr = SNESFASSetNumberSmoothDown(snes,n_down);CHKERRQ(ierr); 376 } 377 378 /* set up the default line search for coarse grid corrections */ 379 if (fas->fastype == SNES_FAS_ADDITIVE) { 380 if (!snes->linesearch) { 381 ierr = SNESGetLineSearch(snes, &linesearch);CHKERRQ(ierr); 382 ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHL2);CHKERRQ(ierr); 383 } 384 } 385 386 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 387 /* recursive option setting for the smoothers */ 388 if (next) {ierr = SNESSetFromOptions(next);CHKERRQ(ierr);} 389 PetscFunctionReturn(0); 390 } 391 392 #include <petscdraw.h> 393 #undef __FUNCT__ 394 #define __FUNCT__ "SNESView_FAS" 395 PetscErrorCode SNESView_FAS(SNES snes, PetscViewer viewer) 396 { 397 SNES_FAS *fas = (SNES_FAS*) snes->data; 398 PetscBool isFine,iascii,isdraw; 399 PetscInt i; 400 PetscErrorCode ierr; 401 SNES smoothu, smoothd, levelsnes; 402 403 PetscFunctionBegin; 404 ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); 405 if (isFine) { 406 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 407 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);CHKERRQ(ierr); 408 if (iascii) { 409 ierr = PetscViewerASCIIPrintf(viewer, "FAS: type is %s, levels=%D, cycles=%D\n", SNESFASTypes[fas->fastype], fas->levels, fas->n_cycles);CHKERRQ(ierr); 410 if (fas->galerkin) { 411 ierr = PetscViewerASCIIPrintf(viewer," Using Galerkin computed coarse grid function evaluation\n");CHKERRQ(ierr); 412 } else { 413 ierr = PetscViewerASCIIPrintf(viewer," Not using Galerkin computed coarse grid function evaluation\n");CHKERRQ(ierr); 414 } 415 for (i=0; i<fas->levels; i++) { 416 ierr = SNESFASGetCycleSNES(snes, i, &levelsnes);CHKERRQ(ierr); 417 ierr = SNESFASCycleGetSmootherUp(levelsnes, &smoothu);CHKERRQ(ierr); 418 ierr = SNESFASCycleGetSmootherDown(levelsnes, &smoothd);CHKERRQ(ierr); 419 if (!i) { 420 ierr = PetscViewerASCIIPrintf(viewer,"Coarse grid solver -- level %D -------------------------------\n",i);CHKERRQ(ierr); 421 } else { 422 ierr = PetscViewerASCIIPrintf(viewer,"Down solver (pre-smoother) on level %D -------------------------------\n",i);CHKERRQ(ierr); 423 } 424 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 425 if (smoothd) { 426 ierr = SNESView(smoothd,viewer);CHKERRQ(ierr); 427 } else { 428 ierr = PetscViewerASCIIPrintf(viewer,"Not yet available\n");CHKERRQ(ierr); 429 } 430 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 431 if (i && (smoothd == smoothu)) { 432 ierr = PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) same as down solver (pre-smoother)\n");CHKERRQ(ierr); 433 } else if (i) { 434 ierr = PetscViewerASCIIPrintf(viewer,"Up solver (post-smoother) on level %D -------------------------------\n",i);CHKERRQ(ierr); 435 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 436 if (smoothu) { 437 ierr = SNESView(smoothu,viewer);CHKERRQ(ierr); 438 } else { 439 ierr = PetscViewerASCIIPrintf(viewer,"Not yet available\n");CHKERRQ(ierr); 440 } 441 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 442 } 443 } 444 } else if (isdraw) { 445 PetscDraw draw; 446 PetscReal x,w,y,bottom,th,wth; 447 SNES_FAS *curfas = fas; 448 ierr = PetscViewerDrawGetDraw(viewer,0,&draw);CHKERRQ(ierr); 449 ierr = PetscDrawGetCurrentPoint(draw,&x,&y);CHKERRQ(ierr); 450 ierr = PetscDrawStringGetSize(draw,&wth,&th);CHKERRQ(ierr); 451 bottom = y - th; 452 while (curfas) { 453 if (!curfas->smoothu) { 454 ierr = PetscDrawPushCurrentPoint(draw,x,bottom);CHKERRQ(ierr); 455 if (curfas->smoothd) ierr = SNESView(curfas->smoothd,viewer);CHKERRQ(ierr); 456 ierr = PetscDrawPopCurrentPoint(draw);CHKERRQ(ierr); 457 } else { 458 w = 0.5*PetscMin(1.0-x,x); 459 ierr = PetscDrawPushCurrentPoint(draw,x-w,bottom);CHKERRQ(ierr); 460 if (curfas->smoothd) ierr = SNESView(curfas->smoothd,viewer);CHKERRQ(ierr); 461 ierr = PetscDrawPopCurrentPoint(draw);CHKERRQ(ierr); 462 ierr = PetscDrawPushCurrentPoint(draw,x+w,bottom);CHKERRQ(ierr); 463 if (curfas->smoothu) ierr = SNESView(curfas->smoothu,viewer);CHKERRQ(ierr); 464 ierr = PetscDrawPopCurrentPoint(draw);CHKERRQ(ierr); 465 } 466 /* this is totally bogus but we have no way of knowing how low the previous one was draw to */ 467 bottom -= 5*th; 468 if (curfas->next) curfas = (SNES_FAS*)curfas->next->data; 469 else curfas = NULL; 470 } 471 } 472 } 473 PetscFunctionReturn(0); 474 } 475 476 #undef __FUNCT__ 477 #define __FUNCT__ "SNESFASDownSmooth_Private" 478 /* 479 Defines the action of the downsmoother 480 */ 481 PetscErrorCode SNESFASDownSmooth_Private(SNES snes, Vec B, Vec X, Vec F, PetscReal *fnorm) 482 { 483 PetscErrorCode ierr = 0; 484 SNESConvergedReason reason; 485 Vec FPC; 486 SNES smoothd; 487 SNES_FAS *fas = (SNES_FAS*) snes->data; 488 489 PetscFunctionBegin; 490 ierr = SNESFASCycleGetSmootherDown(snes, &smoothd);CHKERRQ(ierr); 491 ierr = SNESSetInitialFunction(smoothd, F);CHKERRQ(ierr); 492 if (fas->eventsmoothsolve) {ierr = PetscLogEventBegin(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 493 ierr = SNESSolve(smoothd, B, X);CHKERRQ(ierr); 494 if (fas->eventsmoothsolve) {ierr = PetscLogEventEnd(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 495 /* check convergence reason for the smoother */ 496 ierr = SNESGetConvergedReason(smoothd,&reason);CHKERRQ(ierr); 497 if (reason < 0 && !(reason == SNES_DIVERGED_MAX_IT || reason == SNES_DIVERGED_LOCAL_MIN)) { 498 snes->reason = SNES_DIVERGED_INNER; 499 PetscFunctionReturn(0); 500 } 501 ierr = SNESGetFunction(smoothd, &FPC, NULL, NULL);CHKERRQ(ierr); 502 ierr = VecCopy(FPC, F);CHKERRQ(ierr); 503 ierr = SNESGetFunctionNorm(smoothd, fnorm);CHKERRQ(ierr); 504 PetscFunctionReturn(0); 505 } 506 507 508 #undef __FUNCT__ 509 #define __FUNCT__ "SNESFASUpSmooth_Private" 510 /* 511 Defines the action of the upsmoother 512 */ 513 PetscErrorCode SNESFASUpSmooth_Private(SNES snes, Vec B, Vec X, Vec F, PetscReal *fnorm) 514 { 515 PetscErrorCode ierr = 0; 516 SNESConvergedReason reason; 517 Vec FPC; 518 SNES smoothu; 519 SNES_FAS *fas = (SNES_FAS*) snes->data; 520 521 PetscFunctionBegin; 522 ierr = SNESFASCycleGetSmootherUp(snes, &smoothu);CHKERRQ(ierr); 523 if (fas->eventsmoothsolve) {ierr = PetscLogEventBegin(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 524 ierr = SNESSolve(smoothu, B, X);CHKERRQ(ierr); 525 if (fas->eventsmoothsolve) {ierr = PetscLogEventEnd(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 526 /* check convergence reason for the smoother */ 527 ierr = SNESGetConvergedReason(smoothu,&reason);CHKERRQ(ierr); 528 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 529 snes->reason = SNES_DIVERGED_INNER; 530 PetscFunctionReturn(0); 531 } 532 ierr = SNESGetFunction(smoothu, &FPC, NULL, NULL);CHKERRQ(ierr); 533 ierr = VecCopy(FPC, F);CHKERRQ(ierr); 534 ierr = SNESGetFunctionNorm(smoothu, fnorm);CHKERRQ(ierr); 535 PetscFunctionReturn(0); 536 } 537 538 #undef __FUNCT__ 539 #define __FUNCT__ "SNESFASCreateCoarseVec" 540 /*@ 541 SNESFASCreateCoarseVec - create Vec corresponding to a state vector on one level coarser than current level 542 543 Collective 544 545 Input Arguments: 546 . snes - SNESFAS 547 548 Output Arguments: 549 . Xcoarse - vector on level one coarser than snes 550 551 Level: developer 552 553 .seealso: SNESFASSetRestriction(), SNESFASRestrict() 554 @*/ 555 PetscErrorCode SNESFASCreateCoarseVec(SNES snes,Vec *Xcoarse) 556 { 557 PetscErrorCode ierr; 558 SNES_FAS *fas = (SNES_FAS*)snes->data; 559 560 PetscFunctionBegin; 561 if (fas->rscale) { 562 ierr = VecDuplicate(fas->rscale,Xcoarse);CHKERRQ(ierr); 563 } else if (fas->inject) { 564 ierr = MatGetVecs(fas->inject,Xcoarse,NULL);CHKERRQ(ierr); 565 } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"Must set restriction or injection");CHKERRQ(ierr); 566 PetscFunctionReturn(0); 567 } 568 569 #undef __FUNCT__ 570 #define __FUNCT__ "SNESFASRestrict" 571 /*@ 572 SNESFASRestrict - restrict a Vec to the next coarser level 573 574 Collective 575 576 Input Arguments: 577 + fine - SNES from which to restrict 578 - Xfine - vector to restrict 579 580 Output Arguments: 581 . Xcoarse - result of restriction 582 583 Level: developer 584 585 .seealso: SNESFASSetRestriction(), SNESFASSetInjection() 586 @*/ 587 PetscErrorCode SNESFASRestrict(SNES fine,Vec Xfine,Vec Xcoarse) 588 { 589 PetscErrorCode ierr; 590 SNES_FAS *fas = (SNES_FAS*)fine->data; 591 592 PetscFunctionBegin; 593 PetscValidHeaderSpecific(fine,SNES_CLASSID,1); 594 PetscValidHeaderSpecific(Xfine,VEC_CLASSID,2); 595 PetscValidHeaderSpecific(Xcoarse,VEC_CLASSID,3); 596 if (fas->inject) { 597 ierr = MatRestrict(fas->inject,Xfine,Xcoarse);CHKERRQ(ierr); 598 } else { 599 ierr = MatRestrict(fas->restrct,Xfine,Xcoarse);CHKERRQ(ierr); 600 ierr = VecPointwiseMult(Xcoarse,fas->rscale,Xcoarse);CHKERRQ(ierr); 601 } 602 PetscFunctionReturn(0); 603 } 604 605 #undef __FUNCT__ 606 #define __FUNCT__ "SNESFASCoarseCorrection" 607 /* 608 609 Performs the FAS coarse correction as: 610 611 fine problem: F(x) = 0 612 coarse problem: F^c(x) = b^c 613 614 b^c = F^c(I^c_fx^f - I^c_fF(x)) 615 616 */ 617 PetscErrorCode SNESFASCoarseCorrection(SNES snes, Vec X, Vec F, Vec X_new) 618 { 619 PetscErrorCode ierr; 620 Vec X_c, Xo_c, F_c, B_c; 621 SNESConvergedReason reason; 622 SNES next; 623 Mat restrct, interpolate; 624 SNES_FAS *fasc; 625 626 PetscFunctionBegin; 627 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 628 if (next) { 629 fasc = (SNES_FAS*)next->data; 630 631 ierr = SNESFASCycleGetRestriction(snes, &restrct);CHKERRQ(ierr); 632 ierr = SNESFASCycleGetInterpolation(snes, &interpolate);CHKERRQ(ierr); 633 634 X_c = next->vec_sol; 635 Xo_c = next->work[0]; 636 F_c = next->vec_func; 637 B_c = next->vec_rhs; 638 639 if (fasc->eventinterprestrict) {ierr = PetscLogEventBegin(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 640 ierr = SNESFASRestrict(snes,X,Xo_c);CHKERRQ(ierr); 641 /* restrict the defect */ 642 ierr = MatRestrict(restrct, F, B_c);CHKERRQ(ierr); 643 if (fasc->eventinterprestrict) {ierr = PetscLogEventEnd(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 644 645 if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 646 ierr = SNESComputeFunction(next, Xo_c, F_c);CHKERRQ(ierr); 647 if (fasc->eventresidual) {ierr = PetscLogEventEnd(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 648 649 /* solve the coarse problem corresponding to F^c(x^c) = b^c = F^c(Rx) - R(F(x) - b) */ 650 ierr = VecCopy(B_c, X_c);CHKERRQ(ierr); 651 ierr = VecCopy(F_c, B_c);CHKERRQ(ierr); 652 ierr = VecCopy(X_c, F_c);CHKERRQ(ierr); 653 /* set initial guess of the coarse problem to the projected fine solution */ 654 ierr = VecCopy(Xo_c, X_c);CHKERRQ(ierr); 655 656 /* recurse to the next level */ 657 ierr = SNESSetInitialFunction(next, F_c);CHKERRQ(ierr); 658 ierr = SNESSolve(next, B_c, X_c);CHKERRQ(ierr); 659 ierr = SNESGetConvergedReason(next,&reason);CHKERRQ(ierr); 660 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 661 snes->reason = SNES_DIVERGED_INNER; 662 PetscFunctionReturn(0); 663 } 664 /* correct as x <- x + I(x^c - Rx)*/ 665 ierr = VecAXPY(X_c, -1.0, Xo_c);CHKERRQ(ierr); 666 667 if (fasc->eventinterprestrict) {ierr = PetscLogEventBegin(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 668 ierr = MatInterpolateAdd(interpolate, X_c, X, X_new);CHKERRQ(ierr); 669 if (fasc->eventinterprestrict) {ierr = PetscLogEventEnd(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 670 } 671 PetscFunctionReturn(0); 672 } 673 674 #undef __FUNCT__ 675 #define __FUNCT__ "SNESFASCycle_Additive" 676 /* 677 678 The additive cycle looks like: 679 680 xhat = x 681 xhat = dS(x, b) 682 x = coarsecorrection(xhat, b_d) 683 x = x + nu*(xhat - x); 684 (optional) x = uS(x, b) 685 686 With the coarse RHS (defect correction) as below. 687 688 */ 689 PetscErrorCode SNESFASCycle_Additive(SNES snes, Vec X) 690 { 691 Vec F, B, Xhat; 692 Vec X_c, Xo_c, F_c, B_c; 693 PetscErrorCode ierr; 694 SNESConvergedReason reason; 695 PetscReal xnorm, fnorm, ynorm; 696 PetscBool lssuccess; 697 SNES next; 698 Mat restrct, interpolate; 699 SNES_FAS *fas = (SNES_FAS*)snes->data,*fasc; 700 701 PetscFunctionBegin; 702 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 703 F = snes->vec_func; 704 B = snes->vec_rhs; 705 Xhat = snes->work[1]; 706 ierr = VecCopy(X, Xhat);CHKERRQ(ierr); 707 /* recurse first */ 708 if (next) { 709 fasc = (SNES_FAS*)next->data; 710 ierr = SNESFASCycleGetRestriction(snes, &restrct);CHKERRQ(ierr); 711 ierr = SNESFASCycleGetInterpolation(snes, &interpolate);CHKERRQ(ierr); 712 if (fas->eventresidual) {ierr = PetscLogEventBegin(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 713 ierr = SNESComputeFunction(snes, Xhat, F);CHKERRQ(ierr); 714 if (fas->eventresidual) {ierr = PetscLogEventEnd(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 715 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); 716 X_c = next->vec_sol; 717 Xo_c = next->work[0]; 718 F_c = next->vec_func; 719 B_c = next->vec_rhs; 720 721 ierr = SNESFASRestrict(snes,Xhat,Xo_c);CHKERRQ(ierr); 722 /* restrict the defect */ 723 ierr = MatRestrict(restrct, F, B_c);CHKERRQ(ierr); 724 725 /* solve the coarse problem corresponding to F^c(x^c) = b^c = Rb + F^c(Rx) - RF(x) */ 726 if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 727 ierr = SNESComputeFunction(next, Xo_c, F_c);CHKERRQ(ierr); 728 if (fasc->eventresidual) {ierr = PetscLogEventEnd(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 729 ierr = VecCopy(B_c, X_c);CHKERRQ(ierr); 730 ierr = VecCopy(F_c, B_c);CHKERRQ(ierr); 731 ierr = VecCopy(X_c, F_c);CHKERRQ(ierr); 732 /* set initial guess of the coarse problem to the projected fine solution */ 733 ierr = VecCopy(Xo_c, X_c);CHKERRQ(ierr); 734 735 /* recurse */ 736 ierr = SNESSetInitialFunction(next, F_c);CHKERRQ(ierr); 737 ierr = SNESSolve(next, B_c, X_c);CHKERRQ(ierr); 738 739 /* smooth on this level */ 740 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &fnorm);CHKERRQ(ierr); 741 742 ierr = SNESGetConvergedReason(next,&reason);CHKERRQ(ierr); 743 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 744 snes->reason = SNES_DIVERGED_INNER; 745 PetscFunctionReturn(0); 746 } 747 748 /* correct as x <- x + I(x^c - Rx)*/ 749 ierr = VecAYPX(X_c, -1.0, Xo_c);CHKERRQ(ierr); 750 ierr = MatInterpolate(interpolate, X_c, Xhat);CHKERRQ(ierr); 751 752 /* additive correction of the coarse direction*/ 753 ierr = SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Xhat);CHKERRQ(ierr); 754 ierr = SNESLineSearchGetSuccess(snes->linesearch, &lssuccess);CHKERRQ(ierr); 755 if (!lssuccess) { 756 if (++snes->numFailures >= snes->maxFailures) { 757 snes->reason = SNES_DIVERGED_LINE_SEARCH; 758 PetscFunctionReturn(0); 759 } 760 } 761 ierr = SNESLineSearchGetNorms(snes->linesearch, &xnorm, &snes->norm, &ynorm);CHKERRQ(ierr); 762 } else { 763 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 764 } 765 PetscFunctionReturn(0); 766 } 767 768 #undef __FUNCT__ 769 #define __FUNCT__ "SNESFASCycle_Multiplicative" 770 /* 771 772 Defines the FAS cycle as: 773 774 fine problem: F(x) = 0 775 coarse problem: F^c(x) = b^c 776 777 b^c = F^c(I^c_fx^f - I^c_fF(x)) 778 779 correction: 780 781 x = x + I(x^c - Rx) 782 783 */ 784 PetscErrorCode SNESFASCycle_Multiplicative(SNES snes, Vec X) 785 { 786 787 PetscErrorCode ierr; 788 Vec F,B; 789 SNES_FAS *fas = (SNES_FAS*)snes->data; 790 791 PetscFunctionBegin; 792 F = snes->vec_func; 793 B = snes->vec_rhs; 794 /* pre-smooth -- just update using the pre-smoother */ 795 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 796 if (fas->level != 0) { 797 ierr = SNESFASCoarseCorrection(snes, X, F, X);CHKERRQ(ierr); 798 ierr = SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 799 } 800 PetscFunctionReturn(0); 801 } 802 803 #undef __FUNCT__ 804 #define __FUNCT__ "SNESSolve_FAS" 805 806 PetscErrorCode SNESSolve_FAS(SNES snes) 807 { 808 PetscErrorCode ierr; 809 PetscInt i, maxits; 810 Vec X, F; 811 PetscReal fnorm; 812 SNES_FAS *fas = (SNES_FAS*)snes->data,*ffas; 813 DM dm; 814 PetscBool isFine; 815 816 PetscFunctionBegin; 817 maxits = snes->max_its; /* maximum number of iterations */ 818 snes->reason = SNES_CONVERGED_ITERATING; 819 X = snes->vec_sol; 820 F = snes->vec_func; 821 822 ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); 823 /*norm setup */ 824 ierr = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 825 snes->iter = 0; 826 snes->norm = 0.; 827 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 828 if (!snes->vec_func_init_set) { 829 if (fas->eventresidual) {ierr = PetscLogEventBegin(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 830 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 831 if (fas->eventresidual) {ierr = PetscLogEventEnd(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 832 if (snes->domainerror) { 833 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 834 PetscFunctionReturn(0); 835 } 836 } else snes->vec_func_init_set = PETSC_FALSE; 837 838 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 839 if (PetscIsInfOrNanReal(fnorm)) { 840 snes->reason = SNES_DIVERGED_FNORM_NAN; 841 PetscFunctionReturn(0); 842 } 843 844 ierr = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 845 snes->norm = fnorm; 846 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 847 ierr = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr); 848 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 849 850 /* test convergence */ 851 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 852 if (snes->reason) PetscFunctionReturn(0); 853 854 855 if (isFine) { 856 /* propagate scale-dependent data up the hierarchy */ 857 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 858 for (ffas=fas; ffas->next; ffas=(SNES_FAS*)ffas->next->data) { 859 DM dmcoarse; 860 ierr = SNESGetDM(ffas->next,&dmcoarse);CHKERRQ(ierr); 861 ierr = DMRestrict(dm,ffas->restrct,ffas->rscale,ffas->inject,dmcoarse);CHKERRQ(ierr); 862 dm = dmcoarse; 863 } 864 } 865 866 for (i = 0; i < maxits; i++) { 867 /* Call general purpose update function */ 868 869 if (snes->ops->update) { 870 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 871 } 872 if (fas->fastype == SNES_FAS_MULTIPLICATIVE) { 873 ierr = SNESFASCycle_Multiplicative(snes, X);CHKERRQ(ierr); 874 } else { 875 ierr = SNESFASCycle_Additive(snes, X);CHKERRQ(ierr); 876 } 877 878 /* check for FAS cycle divergence */ 879 if (snes->reason != SNES_CONVERGED_ITERATING) PetscFunctionReturn(0); 880 881 /* Monitor convergence */ 882 ierr = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 883 snes->iter = i+1; 884 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 885 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 886 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 887 /* Test for convergence */ 888 if (isFine) { 889 ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,snes->norm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 890 if (snes->reason) break; 891 } 892 } 893 if (i == maxits) { 894 ierr = PetscInfo1(snes, "Maximum number of iterations has been reached: %D\n", maxits);CHKERRQ(ierr); 895 if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT; 896 } 897 PetscFunctionReturn(0); 898 } 899