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