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