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,0,0,0,0);CHKERRQ(ierr);} 166 ierr = SNESSetUp(fas->smoothd);CHKERRQ(ierr); 167 if (fas->eventsmoothsetup) {ierr = PetscLogEventEnd(fas->eventsmoothsetup,0,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,0,0,0,0);CHKERRQ(ierr);} 195 ierr = SNESSetUp(fas->smoothu);CHKERRQ(ierr); 196 if (fas->eventsmoothsetup) {ierr = PetscLogEventEnd(fas->eventsmoothsetup,0,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,0,0,0,0);CHKERRQ(ierr);} 417 ierr = SNESSolve(smoothd, B, X);CHKERRQ(ierr); 418 if (fas->eventsmoothsolve) {ierr = PetscLogEventEnd(fas->eventsmoothsolve,0,0,0,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,0,0,0,0);CHKERRQ(ierr);} 452 ierr = SNESSolve(smoothu, B, X);CHKERRQ(ierr); 453 if (fas->eventsmoothsolve) {ierr = PetscLogEventEnd(fas->eventsmoothsolve,0,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,0,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,0,0,0,0);CHKERRQ(ierr);} 570 571 if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,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,0,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,0,0,0,0);CHKERRQ(ierr);} 595 ierr = MatInterpolateAdd(interpolate, X_c, X, X_new);CHKERRQ(ierr); 596 if (fasc->eventinterprestrict) {ierr = PetscLogEventEnd(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 597 } 598 PetscFunctionReturn(0); 599 } 600 601 /* 602 603 The additive cycle looks like: 604 605 xhat = x 606 xhat = dS(x, b) 607 x = coarsecorrection(xhat, b_d) 608 x = x + nu*(xhat - x); 609 (optional) x = uS(x, b) 610 611 With the coarse RHS (defect correction) as below. 612 613 */ 614 static PetscErrorCode SNESFASCycle_Additive(SNES snes, Vec X) 615 { 616 Vec F, B, Xhat; 617 Vec X_c, Xo_c, F_c, B_c; 618 PetscErrorCode ierr; 619 SNESConvergedReason reason; 620 PetscReal xnorm, fnorm, ynorm; 621 SNESLineSearchReason lsresult; 622 SNES next; 623 Mat restrct, interpolate; 624 SNES_FAS *fas = (SNES_FAS*)snes->data,*fasc; 625 626 PetscFunctionBegin; 627 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 628 F = snes->vec_func; 629 B = snes->vec_rhs; 630 Xhat = snes->work[1]; 631 ierr = VecCopy(X, Xhat);CHKERRQ(ierr); 632 /* recurse first */ 633 if (next) { 634 fasc = (SNES_FAS*)next->data; 635 ierr = SNESFASCycleGetRestriction(snes, &restrct);CHKERRQ(ierr); 636 ierr = SNESFASCycleGetInterpolation(snes, &interpolate);CHKERRQ(ierr); 637 if (fas->eventresidual) {ierr = PetscLogEventBegin(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 638 ierr = SNESComputeFunction(snes, Xhat, F);CHKERRQ(ierr); 639 if (fas->eventresidual) {ierr = PetscLogEventEnd(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 640 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); 641 X_c = next->vec_sol; 642 Xo_c = next->work[0]; 643 F_c = next->vec_func; 644 B_c = next->vec_rhs; 645 646 ierr = SNESFASRestrict(snes,Xhat,Xo_c);CHKERRQ(ierr); 647 /* restrict the defect */ 648 ierr = MatRestrict(restrct, F, B_c);CHKERRQ(ierr); 649 650 /* solve the coarse problem corresponding to F^c(x^c) = b^c = Rb + F^c(Rx) - RF(x) */ 651 if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 652 ierr = SNESComputeFunction(next, Xo_c, F_c);CHKERRQ(ierr); 653 if (fasc->eventresidual) {ierr = PetscLogEventEnd(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 654 ierr = VecCopy(B_c, X_c);CHKERRQ(ierr); 655 ierr = VecCopy(F_c, B_c);CHKERRQ(ierr); 656 ierr = VecCopy(X_c, F_c);CHKERRQ(ierr); 657 /* set initial guess of the coarse problem to the projected fine solution */ 658 ierr = VecCopy(Xo_c, X_c);CHKERRQ(ierr); 659 660 /* recurse */ 661 ierr = SNESSetInitialFunction(next, F_c);CHKERRQ(ierr); 662 ierr = SNESSolve(next, B_c, X_c);CHKERRQ(ierr); 663 664 /* smooth on this level */ 665 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &fnorm);CHKERRQ(ierr); 666 667 ierr = SNESGetConvergedReason(next,&reason);CHKERRQ(ierr); 668 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 669 snes->reason = SNES_DIVERGED_INNER; 670 PetscFunctionReturn(0); 671 } 672 673 /* correct as x <- x + I(x^c - Rx)*/ 674 ierr = VecAYPX(X_c, -1.0, Xo_c);CHKERRQ(ierr); 675 ierr = MatInterpolate(interpolate, X_c, Xhat);CHKERRQ(ierr); 676 677 /* additive correction of the coarse direction*/ 678 ierr = SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Xhat);CHKERRQ(ierr); 679 ierr = SNESLineSearchGetReason(snes->linesearch, &lsresult);CHKERRQ(ierr); 680 ierr = SNESLineSearchGetNorms(snes->linesearch, &xnorm, &snes->norm, &ynorm);CHKERRQ(ierr); 681 if (lsresult) { 682 if (++snes->numFailures >= snes->maxFailures) { 683 snes->reason = SNES_DIVERGED_LINE_SEARCH; 684 PetscFunctionReturn(0); 685 } 686 } 687 } else { 688 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 689 } 690 PetscFunctionReturn(0); 691 } 692 693 /* 694 695 Defines the FAS cycle as: 696 697 fine problem: F(x) = b 698 coarse problem: F^c(x) = b^c 699 700 b^c = F^c(Rx) - R(F(x) - b) 701 702 correction: 703 704 x = x + I(x^c - Rx) 705 706 */ 707 static PetscErrorCode SNESFASCycle_Multiplicative(SNES snes, Vec X) 708 { 709 710 PetscErrorCode ierr; 711 Vec F,B; 712 SNES next; 713 714 PetscFunctionBegin; 715 F = snes->vec_func; 716 B = snes->vec_rhs; 717 /* pre-smooth -- just update using the pre-smoother */ 718 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 719 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 720 if (next) { 721 ierr = SNESFASCoarseCorrection(snes, X, F, X);CHKERRQ(ierr); 722 ierr = SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 723 } 724 PetscFunctionReturn(0); 725 } 726 727 static PetscErrorCode SNESFASCycleSetupPhase_Full(SNES snes) 728 { 729 SNES next; 730 SNES_FAS *fas = (SNES_FAS*)snes->data; 731 PetscBool isFine; 732 PetscErrorCode ierr; 733 734 PetscFunctionBegin; 735 /* pre-smooth -- just update using the pre-smoother */ 736 ierr = SNESFASCycleIsFine(snes,&isFine);CHKERRQ(ierr); 737 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 738 fas->full_stage = 0; 739 if (next) {ierr = SNESFASCycleSetupPhase_Full(next);CHKERRQ(ierr);} 740 PetscFunctionReturn(0); 741 } 742 743 static PetscErrorCode SNESFASCycle_Full(SNES snes, Vec X) 744 { 745 PetscErrorCode ierr; 746 Vec F,B; 747 SNES_FAS *fas = (SNES_FAS*)snes->data; 748 PetscBool isFine; 749 SNES next; 750 751 PetscFunctionBegin; 752 F = snes->vec_func; 753 B = snes->vec_rhs; 754 ierr = SNESFASCycleIsFine(snes,&isFine);CHKERRQ(ierr); 755 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 756 757 if (isFine) { 758 ierr = SNESFASCycleSetupPhase_Full(snes);CHKERRQ(ierr); 759 } 760 761 if (fas->full_stage == 0) { 762 /* downsweep */ 763 if (next) { 764 if (fas->level != 1) next->max_its += 1; 765 if (fas->full_downsweep||isFine) {ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr);} 766 ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); 767 ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 768 if (fas->level != 1) next->max_its -= 1; 769 } else { 770 /* The smoother on the coarse level is the coarse solver */ 771 ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 772 } 773 fas->full_stage = 1; 774 } else if (fas->full_stage == 1) { 775 if (snes->iter == 0) {ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr);} 776 if (next) { 777 ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); 778 ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 779 } 780 } 781 /* final v-cycle */ 782 if (isFine) { 783 if (next) { 784 ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); 785 ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 786 } 787 } 788 PetscFunctionReturn(0); 789 } 790 791 static PetscErrorCode SNESFASCycle_Kaskade(SNES snes, Vec X) 792 { 793 PetscErrorCode ierr; 794 Vec F,B; 795 SNES next; 796 797 PetscFunctionBegin; 798 F = snes->vec_func; 799 B = snes->vec_rhs; 800 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 801 if (next) { 802 ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); 803 ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 804 } else { 805 ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 806 } 807 PetscFunctionReturn(0); 808 } 809 810 PetscBool SNEScite = PETSC_FALSE; 811 const char SNESCitation[] = "@techreport{pbmkbsxt2012,\n" 812 " title = {Composing Scalable Nonlinear Algebraic Solvers},\n" 813 " author = {Peter Brune and Mathew Knepley and Barry Smith and Xuemin Tu},\n" 814 " year = 2013,\n" 815 " type = Preprint,\n" 816 " number = {ANL/MCS-P2010-0112},\n" 817 " institution = {Argonne National Laboratory}\n}\n"; 818 819 static PetscErrorCode SNESSolve_FAS(SNES snes) 820 { 821 PetscErrorCode ierr; 822 PetscInt i, maxits; 823 Vec X, F; 824 PetscReal fnorm; 825 SNES_FAS *fas = (SNES_FAS*)snes->data,*ffas; 826 DM dm; 827 PetscBool isFine; 828 829 PetscFunctionBegin; 830 831 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); 832 833 ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr); 834 maxits = snes->max_its; /* maximum number of iterations */ 835 snes->reason = SNES_CONVERGED_ITERATING; 836 X = snes->vec_sol; 837 F = snes->vec_func; 838 839 ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); 840 /*norm setup */ 841 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 842 snes->iter = 0; 843 snes->norm = 0.; 844 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 845 if (!snes->vec_func_init_set) { 846 if (fas->eventresidual) {ierr = PetscLogEventBegin(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 847 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 848 if (fas->eventresidual) {ierr = PetscLogEventEnd(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 849 } else snes->vec_func_init_set = PETSC_FALSE; 850 851 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 852 SNESCheckFunctionNorm(snes,fnorm); 853 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 854 snes->norm = fnorm; 855 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 856 ierr = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr); 857 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 858 859 /* test convergence */ 860 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 861 if (snes->reason) PetscFunctionReturn(0); 862 863 864 if (isFine) { 865 /* propagate scale-dependent data up the hierarchy */ 866 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 867 for (ffas=fas; ffas->next; ffas=(SNES_FAS*)ffas->next->data) { 868 DM dmcoarse; 869 ierr = SNESGetDM(ffas->next,&dmcoarse);CHKERRQ(ierr); 870 ierr = DMRestrict(dm,ffas->restrct,ffas->rscale,ffas->inject,dmcoarse);CHKERRQ(ierr); 871 dm = dmcoarse; 872 } 873 } 874 875 for (i = 0; i < maxits; i++) { 876 /* Call general purpose update function */ 877 878 if (snes->ops->update) { 879 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 880 } 881 if (fas->fastype == SNES_FAS_MULTIPLICATIVE) { 882 ierr = SNESFASCycle_Multiplicative(snes, X);CHKERRQ(ierr); 883 } else if (fas->fastype == SNES_FAS_ADDITIVE) { 884 ierr = SNESFASCycle_Additive(snes, X);CHKERRQ(ierr); 885 } else if (fas->fastype == SNES_FAS_FULL) { 886 ierr = SNESFASCycle_Full(snes, X);CHKERRQ(ierr); 887 } else if (fas->fastype ==SNES_FAS_KASKADE) { 888 ierr = SNESFASCycle_Kaskade(snes, X);CHKERRQ(ierr); 889 } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"Unsupported FAS type"); 890 891 /* check for FAS cycle divergence */ 892 if (snes->reason != SNES_CONVERGED_ITERATING) PetscFunctionReturn(0); 893 894 /* Monitor convergence */ 895 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 896 snes->iter = i+1; 897 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 898 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 899 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 900 /* Test for convergence */ 901 if (isFine) { 902 ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,snes->norm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 903 if (snes->reason) break; 904 } 905 } 906 if (i == maxits) { 907 ierr = PetscInfo1(snes, "Maximum number of iterations has been reached: %D\n", maxits);CHKERRQ(ierr); 908 if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT; 909 } 910 PetscFunctionReturn(0); 911 } 912 913 /*MC 914 915 SNESFAS - Full Approximation Scheme nonlinear multigrid solver. 916 917 The nonlinear problem is solved by correction using coarse versions 918 of the nonlinear problem. This problem is perturbed so that a projected 919 solution of the fine problem elicits no correction from the coarse problem. 920 921 Options Database: 922 + -snes_fas_levels - The number of levels 923 . -snes_fas_cycles<1> - The number of cycles -- 1 for V, 2 for W 924 . -snes_fas_type<additive,multiplicative,full,kaskade> - Additive or multiplicative cycle 925 . -snes_fas_galerkin<PETSC_FALSE> - Form coarse problems by projection back upon the fine problem 926 . -snes_fas_smoothup<1> - The number of iterations of the post-smoother 927 . -snes_fas_smoothdown<1> - The number of iterations of the pre-smoother 928 . -snes_fas_monitor - Monitor progress of all of the levels 929 . -snes_fas_full_downsweep<PETSC_FALSE> - call the downsmooth on the initial downsweep of full FAS 930 . -fas_levels_snes_ - SNES options for all smoothers 931 . -fas_levels_cycle_snes_ - SNES options for all cycles 932 . -fas_levels_i_snes_ - SNES options for the smoothers on level i 933 . -fas_levels_i_cycle_snes_ - SNES options for the cycle on level i 934 - -fas_coarse_snes_ - SNES options for the coarsest smoother 935 936 Notes: 937 The organization of the FAS solver is slightly different from the organization of PCMG 938 As each level has smoother SNES instances(down and potentially up) and a cycle SNES instance. 939 The cycle SNES instance may be used for monitoring convergence on a particular level. 940 941 Level: beginner 942 943 References: 944 . 1. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers", 945 SIAM Review, 57(4), 2015 946 947 .seealso: PCMG, SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 948 M*/ 949 950 PETSC_EXTERN PetscErrorCode SNESCreate_FAS(SNES snes) 951 { 952 SNES_FAS *fas; 953 PetscErrorCode ierr; 954 955 PetscFunctionBegin; 956 snes->ops->destroy = SNESDestroy_FAS; 957 snes->ops->setup = SNESSetUp_FAS; 958 snes->ops->setfromoptions = SNESSetFromOptions_FAS; 959 snes->ops->view = SNESView_FAS; 960 snes->ops->solve = SNESSolve_FAS; 961 snes->ops->reset = SNESReset_FAS; 962 963 snes->usesksp = PETSC_FALSE; 964 snes->usesnpc = PETSC_FALSE; 965 966 if (!snes->tolerancesset) { 967 snes->max_funcs = 30000; 968 snes->max_its = 10000; 969 } 970 971 snes->alwayscomputesfinalresidual = PETSC_TRUE; 972 973 ierr = PetscNewLog(snes,&fas);CHKERRQ(ierr); 974 975 snes->data = (void*) fas; 976 fas->level = 0; 977 fas->levels = 1; 978 fas->n_cycles = 1; 979 fas->max_up_it = 1; 980 fas->max_down_it = 1; 981 fas->smoothu = NULL; 982 fas->smoothd = NULL; 983 fas->next = NULL; 984 fas->previous = NULL; 985 fas->fine = snes; 986 fas->interpolate = NULL; 987 fas->restrct = NULL; 988 fas->inject = NULL; 989 fas->usedmfornumberoflevels = PETSC_FALSE; 990 fas->fastype = SNES_FAS_MULTIPLICATIVE; 991 fas->full_downsweep = PETSC_FALSE; 992 993 fas->eventsmoothsetup = 0; 994 fas->eventsmoothsolve = 0; 995 fas->eventresidual = 0; 996 fas->eventinterprestrict = 0; 997 PetscFunctionReturn(0); 998 } 999