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