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 #undef __FUNCT__ 9 #define __FUNCT__ "SNESReset_FAS" 10 static PetscErrorCode SNESReset_FAS(SNES snes) 11 { 12 PetscErrorCode ierr = 0; 13 SNES_FAS * fas = (SNES_FAS*)snes->data; 14 15 PetscFunctionBegin; 16 ierr = SNESDestroy(&fas->smoothu);CHKERRQ(ierr); 17 ierr = SNESDestroy(&fas->smoothd);CHKERRQ(ierr); 18 ierr = MatDestroy(&fas->inject);CHKERRQ(ierr); 19 ierr = MatDestroy(&fas->interpolate);CHKERRQ(ierr); 20 ierr = MatDestroy(&fas->restrct);CHKERRQ(ierr); 21 ierr = VecDestroy(&fas->rscale);CHKERRQ(ierr); 22 if (fas->galerkin) { 23 ierr = VecDestroy(&fas->Xg);CHKERRQ(ierr); 24 ierr = VecDestroy(&fas->Fg);CHKERRQ(ierr); 25 } 26 if (fas->next) {ierr = SNESReset(fas->next);CHKERRQ(ierr);} 27 PetscFunctionReturn(0); 28 } 29 30 #undef __FUNCT__ 31 #define __FUNCT__ "SNESDestroy_FAS" 32 static PetscErrorCode SNESDestroy_FAS(SNES snes) 33 { 34 SNES_FAS * fas = (SNES_FAS*)snes->data; 35 PetscErrorCode ierr = 0; 36 37 PetscFunctionBegin; 38 /* recursively resets and then destroys */ 39 ierr = SNESReset(snes);CHKERRQ(ierr); 40 if (fas->next) { 41 ierr = SNESDestroy(&fas->next);CHKERRQ(ierr); 42 } 43 ierr = PetscFree(fas);CHKERRQ(ierr); 44 PetscFunctionReturn(0); 45 } 46 47 #undef __FUNCT__ 48 #define __FUNCT__ "SNESSetUp_FAS" 49 static PetscErrorCode SNESSetUp_FAS(SNES snes) 50 { 51 SNES_FAS *fas = (SNES_FAS*) snes->data; 52 PetscErrorCode ierr; 53 PetscInt dm_levels; 54 Vec vec_sol, vec_func, vec_sol_update, vec_rhs; /* preserve these if they're set through the reset */ 55 SNES next; 56 PetscBool isFine; 57 SNESLineSearch linesearch; 58 SNESLineSearch slinesearch; 59 void *lsprectx,*lspostctx; 60 PetscErrorCode (*precheck)(SNESLineSearch,Vec,Vec,PetscBool*,void*); 61 PetscErrorCode (*postcheck)(SNESLineSearch,Vec,Vec,Vec,PetscBool*,PetscBool*,void*); 62 63 PetscFunctionBegin; 64 ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); 65 if (fas->usedmfornumberoflevels && isFine) { 66 ierr = DMGetRefineLevel(snes->dm,&dm_levels);CHKERRQ(ierr); 67 dm_levels++; 68 if (dm_levels > fas->levels) { 69 /* we don't want the solution and func vectors to be destroyed in the SNESReset when it's called in SNESFASSetLevels_FAS*/ 70 vec_sol = snes->vec_sol; 71 vec_func = snes->vec_func; 72 vec_sol_update = snes->vec_sol_update; 73 vec_rhs = snes->vec_rhs; 74 snes->vec_sol = NULL; 75 snes->vec_func = NULL; 76 snes->vec_sol_update = NULL; 77 snes->vec_rhs = NULL; 78 79 /* reset the number of levels */ 80 ierr = SNESFASSetLevels(snes,dm_levels,NULL);CHKERRQ(ierr); 81 ierr = SNESSetFromOptions(snes);CHKERRQ(ierr); 82 83 snes->vec_sol = vec_sol; 84 snes->vec_func = vec_func; 85 snes->vec_rhs = vec_rhs; 86 snes->vec_sol_update = vec_sol_update; 87 } 88 } 89 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 90 if (!isFine) snes->gridsequence = 0; /* no grid sequencing inside the multigrid hierarchy! */ 91 92 ierr = SNESSetWorkVecs(snes, 2);CHKERRQ(ierr); /* work vectors used for intergrid transfers */ 93 94 /* set up the smoothers if they haven't already been set up */ 95 if (!fas->smoothd) { 96 ierr = SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd);CHKERRQ(ierr); 97 } 98 99 if (snes->dm) { 100 /* set the smoother DMs properly */ 101 if (fas->smoothu) ierr = SNESSetDM(fas->smoothu, snes->dm);CHKERRQ(ierr); 102 ierr = SNESSetDM(fas->smoothd, snes->dm);CHKERRQ(ierr); 103 /* construct EVERYTHING from the DM -- including the progressive set of smoothers */ 104 if (next) { 105 /* for now -- assume the DM and the evaluation functions have been set externally */ 106 if (!next->dm) { 107 ierr = DMCoarsen(snes->dm, PetscObjectComm((PetscObject)next), &next->dm);CHKERRQ(ierr); 108 ierr = SNESSetDM(next, next->dm);CHKERRQ(ierr); 109 } 110 /* set the interpolation and restriction from the DM */ 111 if (!fas->interpolate) { 112 ierr = DMCreateInterpolation(next->dm, snes->dm, &fas->interpolate, &fas->rscale);CHKERRQ(ierr); 113 if (!fas->restrct) { 114 ierr = PetscObjectReference((PetscObject)fas->interpolate);CHKERRQ(ierr); 115 fas->restrct = fas->interpolate; 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 #undef __FUNCT__ 219 #define __FUNCT__ "SNESSetFromOptions_FAS" 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 #undef __FUNCT__ 319 #define __FUNCT__ "SNESView_FAS" 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, "FAS: 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 #undef __FUNCT__ 402 #define __FUNCT__ "SNESFASDownSmooth_Private" 403 /* 404 Defines the action of the downsmoother 405 */ 406 static PetscErrorCode SNESFASDownSmooth_Private(SNES snes, Vec B, Vec X, Vec F, PetscReal *fnorm) 407 { 408 PetscErrorCode ierr = 0; 409 SNESConvergedReason reason; 410 Vec FPC; 411 SNES smoothd; 412 PetscBool flg; 413 SNES_FAS *fas = (SNES_FAS*) snes->data; 414 415 PetscFunctionBegin; 416 ierr = SNESFASCycleGetSmootherDown(snes, &smoothd);CHKERRQ(ierr); 417 ierr = SNESSetInitialFunction(smoothd, F);CHKERRQ(ierr); 418 if (fas->eventsmoothsolve) {ierr = PetscLogEventBegin(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 419 ierr = SNESSolve(smoothd, B, X);CHKERRQ(ierr); 420 if (fas->eventsmoothsolve) {ierr = PetscLogEventEnd(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 421 /* check convergence reason for the smoother */ 422 ierr = SNESGetConvergedReason(smoothd,&reason);CHKERRQ(ierr); 423 if (reason < 0 && !(reason == SNES_DIVERGED_MAX_IT || reason == SNES_DIVERGED_LOCAL_MIN || reason == SNES_DIVERGED_LINE_SEARCH)) { 424 snes->reason = SNES_DIVERGED_INNER; 425 PetscFunctionReturn(0); 426 } 427 428 ierr = SNESGetFunction(smoothd, &FPC, NULL, NULL);CHKERRQ(ierr); 429 ierr = SNESGetAlwaysComputesFinalResidual(smoothd, &flg);CHKERRQ(ierr); 430 if (!flg) { 431 ierr = SNESComputeFunction(smoothd, X, FPC);CHKERRQ(ierr); 432 } 433 ierr = VecCopy(FPC, F);CHKERRQ(ierr); 434 if (fnorm) {ierr = VecNorm(F,NORM_2,fnorm);CHKERRQ(ierr);} 435 PetscFunctionReturn(0); 436 } 437 438 439 #undef __FUNCT__ 440 #define __FUNCT__ "SNESFASUpSmooth_Private" 441 /* 442 Defines the action of the upsmoother 443 */ 444 static PetscErrorCode SNESFASUpSmooth_Private(SNES snes, Vec B, Vec X, Vec F, PetscReal *fnorm) 445 { 446 PetscErrorCode ierr = 0; 447 SNESConvergedReason reason; 448 Vec FPC; 449 SNES smoothu; 450 PetscBool flg; 451 SNES_FAS *fas = (SNES_FAS*) snes->data; 452 453 PetscFunctionBegin; 454 ierr = SNESFASCycleGetSmootherUp(snes, &smoothu);CHKERRQ(ierr); 455 if (fas->eventsmoothsolve) {ierr = PetscLogEventBegin(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 456 ierr = SNESSolve(smoothu, B, X);CHKERRQ(ierr); 457 if (fas->eventsmoothsolve) {ierr = PetscLogEventEnd(fas->eventsmoothsolve,0,0,0,0);CHKERRQ(ierr);} 458 /* check convergence reason for the smoother */ 459 ierr = SNESGetConvergedReason(smoothu,&reason);CHKERRQ(ierr); 460 if (reason < 0 && !(reason == SNES_DIVERGED_MAX_IT || reason == SNES_DIVERGED_LOCAL_MIN || reason == SNES_DIVERGED_LINE_SEARCH)) { 461 snes->reason = SNES_DIVERGED_INNER; 462 PetscFunctionReturn(0); 463 } 464 ierr = SNESGetFunction(smoothu, &FPC, NULL, NULL);CHKERRQ(ierr); 465 ierr = SNESGetAlwaysComputesFinalResidual(smoothu, &flg);CHKERRQ(ierr); 466 if (!flg) { 467 ierr = SNESComputeFunction(smoothu, X, FPC);CHKERRQ(ierr); 468 } 469 ierr = VecCopy(FPC, F);CHKERRQ(ierr); 470 if (fnorm) {ierr = VecNorm(F,NORM_2,fnorm);CHKERRQ(ierr);} 471 PetscFunctionReturn(0); 472 } 473 474 #undef __FUNCT__ 475 #define __FUNCT__ "SNESFASCreateCoarseVec" 476 /*@ 477 SNESFASCreateCoarseVec - create Vec corresponding to a state vector on one level coarser than current level 478 479 Collective 480 481 Input Arguments: 482 . snes - SNESFAS 483 484 Output Arguments: 485 . Xcoarse - vector on level one coarser than snes 486 487 Level: developer 488 489 .seealso: SNESFASSetRestriction(), SNESFASRestrict() 490 @*/ 491 PetscErrorCode SNESFASCreateCoarseVec(SNES snes,Vec *Xcoarse) 492 { 493 PetscErrorCode ierr; 494 SNES_FAS *fas = (SNES_FAS*)snes->data; 495 496 PetscFunctionBegin; 497 if (fas->rscale) { 498 ierr = VecDuplicate(fas->rscale,Xcoarse);CHKERRQ(ierr); 499 } else if (fas->inject) { 500 ierr = MatCreateVecs(fas->inject,Xcoarse,NULL);CHKERRQ(ierr); 501 } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"Must set restriction or injection");CHKERRQ(ierr); 502 PetscFunctionReturn(0); 503 } 504 505 #undef __FUNCT__ 506 #define __FUNCT__ "SNESFASRestrict" 507 /*@ 508 SNESFASRestrict - restrict a Vec to the next coarser level 509 510 Collective 511 512 Input Arguments: 513 + fine - SNES from which to restrict 514 - Xfine - vector to restrict 515 516 Output Arguments: 517 . Xcoarse - result of restriction 518 519 Level: developer 520 521 .seealso: SNESFASSetRestriction(), SNESFASSetInjection() 522 @*/ 523 PetscErrorCode SNESFASRestrict(SNES fine,Vec Xfine,Vec Xcoarse) 524 { 525 PetscErrorCode ierr; 526 SNES_FAS *fas = (SNES_FAS*)fine->data; 527 528 PetscFunctionBegin; 529 PetscValidHeaderSpecific(fine,SNES_CLASSID,1); 530 PetscValidHeaderSpecific(Xfine,VEC_CLASSID,2); 531 PetscValidHeaderSpecific(Xcoarse,VEC_CLASSID,3); 532 if (fas->inject) { 533 ierr = MatRestrict(fas->inject,Xfine,Xcoarse);CHKERRQ(ierr); 534 } else { 535 ierr = MatRestrict(fas->restrct,Xfine,Xcoarse);CHKERRQ(ierr); 536 ierr = VecPointwiseMult(Xcoarse,fas->rscale,Xcoarse);CHKERRQ(ierr); 537 } 538 PetscFunctionReturn(0); 539 } 540 541 #undef __FUNCT__ 542 #define __FUNCT__ "SNESFASCoarseCorrection" 543 /* 544 545 Performs the FAS coarse correction as: 546 547 fine problem: F(x) = b 548 coarse problem: F^c(x^c) = b^c 549 550 b^c = F^c(Rx) - R(F(x) - b) 551 552 */ 553 PetscErrorCode SNESFASCoarseCorrection(SNES snes, Vec X, Vec F, Vec X_new) 554 { 555 PetscErrorCode ierr; 556 Vec X_c, Xo_c, F_c, B_c; 557 SNESConvergedReason reason; 558 SNES next; 559 Mat restrct, interpolate; 560 SNES_FAS *fasc; 561 562 PetscFunctionBegin; 563 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 564 if (next) { 565 fasc = (SNES_FAS*)next->data; 566 567 ierr = SNESFASCycleGetRestriction(snes, &restrct);CHKERRQ(ierr); 568 ierr = SNESFASCycleGetInterpolation(snes, &interpolate);CHKERRQ(ierr); 569 570 X_c = next->vec_sol; 571 Xo_c = next->work[0]; 572 F_c = next->vec_func; 573 B_c = next->vec_rhs; 574 575 if (fasc->eventinterprestrict) {ierr = PetscLogEventBegin(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 576 ierr = SNESFASRestrict(snes,X,Xo_c);CHKERRQ(ierr); 577 /* restrict the defect: R(F(x) - b) */ 578 ierr = MatRestrict(restrct, F, B_c);CHKERRQ(ierr); 579 if (fasc->eventinterprestrict) {ierr = PetscLogEventEnd(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 580 581 if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 582 /* F_c = F^c(Rx) - R(F(x) - b) since the second term was sitting in next->vec_rhs */ 583 ierr = SNESComputeFunction(next, Xo_c, F_c);CHKERRQ(ierr); 584 if (fasc->eventresidual) {ierr = PetscLogEventEnd(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 585 586 /* solve the coarse problem corresponding to F^c(x^c) = b^c = F^c(Rx) - R(F(x) - b) */ 587 ierr = VecCopy(B_c, X_c);CHKERRQ(ierr); 588 ierr = VecCopy(F_c, B_c);CHKERRQ(ierr); 589 ierr = VecCopy(X_c, F_c);CHKERRQ(ierr); 590 /* set initial guess of the coarse problem to the projected fine solution */ 591 ierr = VecCopy(Xo_c, X_c);CHKERRQ(ierr); 592 593 /* recurse to the next level */ 594 ierr = SNESSetInitialFunction(next, F_c);CHKERRQ(ierr); 595 ierr = SNESSolve(next, B_c, X_c);CHKERRQ(ierr); 596 ierr = SNESGetConvergedReason(next,&reason);CHKERRQ(ierr); 597 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 598 snes->reason = SNES_DIVERGED_INNER; 599 PetscFunctionReturn(0); 600 } 601 /* correct as x <- x + I(x^c - Rx)*/ 602 ierr = VecAXPY(X_c, -1.0, Xo_c);CHKERRQ(ierr); 603 604 if (fasc->eventinterprestrict) {ierr = PetscLogEventBegin(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 605 ierr = MatInterpolateAdd(interpolate, X_c, X, X_new);CHKERRQ(ierr); 606 if (fasc->eventinterprestrict) {ierr = PetscLogEventEnd(fasc->eventinterprestrict,0,0,0,0);CHKERRQ(ierr);} 607 } 608 PetscFunctionReturn(0); 609 } 610 611 #undef __FUNCT__ 612 #define __FUNCT__ "SNESFASCycle_Additive" 613 /* 614 615 The additive cycle looks like: 616 617 xhat = x 618 xhat = dS(x, b) 619 x = coarsecorrection(xhat, b_d) 620 x = x + nu*(xhat - x); 621 (optional) x = uS(x, b) 622 623 With the coarse RHS (defect correction) as below. 624 625 */ 626 static PetscErrorCode SNESFASCycle_Additive(SNES snes, Vec X) 627 { 628 Vec F, B, Xhat; 629 Vec X_c, Xo_c, F_c, B_c; 630 PetscErrorCode ierr; 631 SNESConvergedReason reason; 632 PetscReal xnorm, fnorm, ynorm; 633 SNESLineSearchReason lsresult; 634 SNES next; 635 Mat restrct, interpolate; 636 SNES_FAS *fas = (SNES_FAS*)snes->data,*fasc; 637 638 PetscFunctionBegin; 639 ierr = SNESFASCycleGetCorrection(snes, &next);CHKERRQ(ierr); 640 F = snes->vec_func; 641 B = snes->vec_rhs; 642 Xhat = snes->work[1]; 643 ierr = VecCopy(X, Xhat);CHKERRQ(ierr); 644 /* recurse first */ 645 if (next) { 646 fasc = (SNES_FAS*)next->data; 647 ierr = SNESFASCycleGetRestriction(snes, &restrct);CHKERRQ(ierr); 648 ierr = SNESFASCycleGetInterpolation(snes, &interpolate);CHKERRQ(ierr); 649 if (fas->eventresidual) {ierr = PetscLogEventBegin(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 650 ierr = SNESComputeFunction(snes, Xhat, F);CHKERRQ(ierr); 651 if (fas->eventresidual) {ierr = PetscLogEventEnd(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 652 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); 653 X_c = next->vec_sol; 654 Xo_c = next->work[0]; 655 F_c = next->vec_func; 656 B_c = next->vec_rhs; 657 658 ierr = SNESFASRestrict(snes,Xhat,Xo_c);CHKERRQ(ierr); 659 /* restrict the defect */ 660 ierr = MatRestrict(restrct, F, B_c);CHKERRQ(ierr); 661 662 /* solve the coarse problem corresponding to F^c(x^c) = b^c = Rb + F^c(Rx) - RF(x) */ 663 if (fasc->eventresidual) {ierr = PetscLogEventBegin(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 664 ierr = SNESComputeFunction(next, Xo_c, F_c);CHKERRQ(ierr); 665 if (fasc->eventresidual) {ierr = PetscLogEventEnd(fasc->eventresidual,0,0,0,0);CHKERRQ(ierr);} 666 ierr = VecCopy(B_c, X_c);CHKERRQ(ierr); 667 ierr = VecCopy(F_c, B_c);CHKERRQ(ierr); 668 ierr = VecCopy(X_c, F_c);CHKERRQ(ierr); 669 /* set initial guess of the coarse problem to the projected fine solution */ 670 ierr = VecCopy(Xo_c, X_c);CHKERRQ(ierr); 671 672 /* recurse */ 673 ierr = SNESSetInitialFunction(next, F_c);CHKERRQ(ierr); 674 ierr = SNESSolve(next, B_c, X_c);CHKERRQ(ierr); 675 676 /* smooth on this level */ 677 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &fnorm);CHKERRQ(ierr); 678 679 ierr = SNESGetConvergedReason(next,&reason);CHKERRQ(ierr); 680 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 681 snes->reason = SNES_DIVERGED_INNER; 682 PetscFunctionReturn(0); 683 } 684 685 /* correct as x <- x + I(x^c - Rx)*/ 686 ierr = VecAYPX(X_c, -1.0, Xo_c);CHKERRQ(ierr); 687 ierr = MatInterpolate(interpolate, X_c, Xhat);CHKERRQ(ierr); 688 689 /* additive correction of the coarse direction*/ 690 ierr = SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Xhat);CHKERRQ(ierr); 691 ierr = SNESLineSearchGetReason(snes->linesearch, &lsresult);CHKERRQ(ierr); 692 ierr = SNESLineSearchGetNorms(snes->linesearch, &xnorm, &snes->norm, &ynorm);CHKERRQ(ierr); 693 if (lsresult) { 694 if (++snes->numFailures >= snes->maxFailures) { 695 snes->reason = SNES_DIVERGED_LINE_SEARCH; 696 PetscFunctionReturn(0); 697 } 698 } 699 } else { 700 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 701 } 702 PetscFunctionReturn(0); 703 } 704 705 #undef __FUNCT__ 706 #define __FUNCT__ "SNESFASCycle_Multiplicative" 707 /* 708 709 Defines the FAS cycle as: 710 711 fine problem: F(x) = b 712 coarse problem: F^c(x) = b^c 713 714 b^c = F^c(Rx) - R(F(x) - b) 715 716 correction: 717 718 x = x + I(x^c - Rx) 719 720 */ 721 static PetscErrorCode SNESFASCycle_Multiplicative(SNES snes, Vec X) 722 { 723 724 PetscErrorCode ierr; 725 Vec F,B; 726 SNES next; 727 728 PetscFunctionBegin; 729 F = snes->vec_func; 730 B = snes->vec_rhs; 731 /* pre-smooth -- just update using the pre-smoother */ 732 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 733 ierr = SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 734 if (next) { 735 ierr = SNESFASCoarseCorrection(snes, X, F, X);CHKERRQ(ierr); 736 ierr = SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm);CHKERRQ(ierr); 737 } 738 PetscFunctionReturn(0); 739 } 740 741 #undef __FUNCT__ 742 #define __FUNCT__ "SNESFASCycleSetupPhase_Full" 743 static PetscErrorCode SNESFASCycleSetupPhase_Full(SNES snes) 744 { 745 SNES next; 746 SNES_FAS *fas = (SNES_FAS*)snes->data; 747 PetscBool isFine; 748 PetscErrorCode ierr; 749 750 PetscFunctionBegin; 751 /* pre-smooth -- just update using the pre-smoother */ 752 ierr = SNESFASCycleIsFine(snes,&isFine);CHKERRQ(ierr); 753 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 754 fas->full_stage = 0; 755 if (next) {ierr = SNESFASCycleSetupPhase_Full(next);CHKERRQ(ierr);} 756 PetscFunctionReturn(0); 757 } 758 759 #undef __FUNCT__ 760 #define __FUNCT__ "SNESFASCycle_Full" 761 static PetscErrorCode SNESFASCycle_Full(SNES snes, Vec X) 762 { 763 PetscErrorCode ierr; 764 Vec F,B; 765 SNES_FAS *fas = (SNES_FAS*)snes->data; 766 PetscBool isFine; 767 SNES next; 768 769 PetscFunctionBegin; 770 F = snes->vec_func; 771 B = snes->vec_rhs; 772 ierr = SNESFASCycleIsFine(snes,&isFine);CHKERRQ(ierr); 773 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 774 775 if (isFine) { 776 ierr = SNESFASCycleSetupPhase_Full(snes);CHKERRQ(ierr); 777 } 778 779 if (fas->full_stage == 0) { 780 /* downsweep */ 781 if (next) { 782 if (fas->level != 1) next->max_its += 1; 783 if (fas->full_downsweep||isFine) {ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr);} 784 ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); 785 ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 786 if (fas->level != 1) next->max_its -= 1; 787 } else { 788 /* The smoother on the coarse level is the coarse solver */ 789 ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 790 } 791 fas->full_stage = 1; 792 } else if (fas->full_stage == 1) { 793 if (snes->iter == 0) {ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr);} 794 if (next) { 795 ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); 796 ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 797 } 798 } 799 /* final v-cycle */ 800 if (isFine) { 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 } 805 } 806 PetscFunctionReturn(0); 807 } 808 809 #undef __FUNCT__ 810 #define __FUNCT__ "SNESFASCycle_Kaskade" 811 static PetscErrorCode SNESFASCycle_Kaskade(SNES snes, Vec X) 812 { 813 PetscErrorCode ierr; 814 Vec F,B; 815 SNES next; 816 817 PetscFunctionBegin; 818 F = snes->vec_func; 819 B = snes->vec_rhs; 820 ierr = SNESFASCycleGetCorrection(snes,&next);CHKERRQ(ierr); 821 if (next) { 822 ierr = SNESFASCoarseCorrection(snes,X,F,X);CHKERRQ(ierr); 823 ierr = SNESFASUpSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 824 } else { 825 ierr = SNESFASDownSmooth_Private(snes,B,X,F,&snes->norm);CHKERRQ(ierr); 826 } 827 PetscFunctionReturn(0); 828 } 829 830 PetscBool SNEScite = PETSC_FALSE; 831 const char SNESCitation[] = "@techreport{pbmkbsxt2012,\n" 832 " title = {Composing Scalable Nonlinear Algebraic Solvers},\n" 833 " author = {Peter Brune and Mathew Knepley and Barry Smith and Xuemin Tu},\n" 834 " year = 2013,\n" 835 " type = Preprint,\n" 836 " number = {ANL/MCS-P2010-0112},\n" 837 " institution = {Argonne National Laboratory}\n}\n"; 838 839 #undef __FUNCT__ 840 #define __FUNCT__ "SNESSolve_FAS" 841 static PetscErrorCode SNESSolve_FAS(SNES snes) 842 { 843 PetscErrorCode ierr; 844 PetscInt i, maxits; 845 Vec X, F; 846 PetscReal fnorm; 847 SNES_FAS *fas = (SNES_FAS*)snes->data,*ffas; 848 DM dm; 849 PetscBool isFine; 850 851 PetscFunctionBegin; 852 853 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); 854 855 ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr); 856 maxits = snes->max_its; /* maximum number of iterations */ 857 snes->reason = SNES_CONVERGED_ITERATING; 858 X = snes->vec_sol; 859 F = snes->vec_func; 860 861 ierr = SNESFASCycleIsFine(snes, &isFine);CHKERRQ(ierr); 862 /*norm setup */ 863 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 864 snes->iter = 0; 865 snes->norm = 0.; 866 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 867 if (!snes->vec_func_init_set) { 868 if (fas->eventresidual) {ierr = PetscLogEventBegin(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 869 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 870 if (fas->eventresidual) {ierr = PetscLogEventEnd(fas->eventresidual,0,0,0,0);CHKERRQ(ierr);} 871 } else snes->vec_func_init_set = PETSC_FALSE; 872 873 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 874 SNESCheckFunctionNorm(snes,fnorm); 875 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 876 snes->norm = fnorm; 877 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 878 ierr = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr); 879 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 880 881 /* test convergence */ 882 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 883 if (snes->reason) PetscFunctionReturn(0); 884 885 886 if (isFine) { 887 /* propagate scale-dependent data up the hierarchy */ 888 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 889 for (ffas=fas; ffas->next; ffas=(SNES_FAS*)ffas->next->data) { 890 DM dmcoarse; 891 ierr = SNESGetDM(ffas->next,&dmcoarse);CHKERRQ(ierr); 892 ierr = DMRestrict(dm,ffas->restrct,ffas->rscale,ffas->inject,dmcoarse);CHKERRQ(ierr); 893 dm = dmcoarse; 894 } 895 } 896 897 for (i = 0; i < maxits; i++) { 898 /* Call general purpose update function */ 899 900 if (snes->ops->update) { 901 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 902 } 903 if (fas->fastype == SNES_FAS_MULTIPLICATIVE) { 904 ierr = SNESFASCycle_Multiplicative(snes, X);CHKERRQ(ierr); 905 } else if (fas->fastype == SNES_FAS_ADDITIVE) { 906 ierr = SNESFASCycle_Additive(snes, X);CHKERRQ(ierr); 907 } else if (fas->fastype == SNES_FAS_FULL) { 908 ierr = SNESFASCycle_Full(snes, X);CHKERRQ(ierr); 909 } else if (fas->fastype ==SNES_FAS_KASKADE) { 910 ierr = SNESFASCycle_Kaskade(snes, X);CHKERRQ(ierr); 911 } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"Unsupported FAS type"); 912 913 /* check for FAS cycle divergence */ 914 if (snes->reason != SNES_CONVERGED_ITERATING) PetscFunctionReturn(0); 915 916 /* Monitor convergence */ 917 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 918 snes->iter = i+1; 919 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 920 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 921 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 922 /* Test for convergence */ 923 if (isFine) { 924 ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,snes->norm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 925 if (snes->reason) break; 926 } 927 } 928 if (i == maxits) { 929 ierr = PetscInfo1(snes, "Maximum number of iterations has been reached: %D\n", maxits);CHKERRQ(ierr); 930 if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT; 931 } 932 PetscFunctionReturn(0); 933 } 934 935 /*MC 936 937 SNESFAS - Full Approximation Scheme nonlinear multigrid solver. 938 939 The nonlinear problem is solved by correction using coarse versions 940 of the nonlinear problem. This problem is perturbed so that a projected 941 solution of the fine problem elicits no correction from the coarse problem. 942 943 Options Database: 944 + -snes_fas_levels - The number of levels 945 . -snes_fas_cycles<1> - The number of cycles -- 1 for V, 2 for W 946 . -snes_fas_type<additive,multiplicative,full,kaskade> - Additive or multiplicative cycle 947 . -snes_fas_galerkin<PETSC_FALSE> - Form coarse problems by projection back upon the fine problem 948 . -snes_fas_smoothup<1> - The number of iterations of the post-smoother 949 . -snes_fas_smoothdown<1> - The number of iterations of the pre-smoother 950 . -snes_fas_monitor - Monitor progress of all of the levels 951 . -snes_fas_full_downsweep<PETSC_FALSE> - call the downsmooth on the initial downsweep of full FAS 952 . -fas_levels_snes_ - SNES options for all smoothers 953 . -fas_levels_cycle_snes_ - SNES options for all cycles 954 . -fas_levels_i_snes_ - SNES options for the smoothers on level i 955 . -fas_levels_i_cycle_snes_ - SNES options for the cycle on level i 956 - -fas_coarse_snes_ - SNES options for the coarsest smoother 957 958 Notes: 959 The organization of the FAS solver is slightly different from the organization of PCMG 960 As each level has smoother SNES instances(down and potentially up) and a cycle SNES instance. 961 The cycle SNES instance may be used for monitoring convergence on a particular level. 962 963 Level: beginner 964 965 References: 966 . 1. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers", 967 SIAM Review, 57(4), 2015 968 969 .seealso: PCMG, SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 970 M*/ 971 972 #undef __FUNCT__ 973 #define __FUNCT__ "SNESCreate_FAS" 974 PETSC_EXTERN PetscErrorCode SNESCreate_FAS(SNES snes) 975 { 976 SNES_FAS *fas; 977 PetscErrorCode ierr; 978 979 PetscFunctionBegin; 980 snes->ops->destroy = SNESDestroy_FAS; 981 snes->ops->setup = SNESSetUp_FAS; 982 snes->ops->setfromoptions = SNESSetFromOptions_FAS; 983 snes->ops->view = SNESView_FAS; 984 snes->ops->solve = SNESSolve_FAS; 985 snes->ops->reset = SNESReset_FAS; 986 987 snes->usesksp = PETSC_FALSE; 988 snes->usespc = PETSC_FALSE; 989 990 if (!snes->tolerancesset) { 991 snes->max_funcs = 30000; 992 snes->max_its = 10000; 993 } 994 995 snes->alwayscomputesfinalresidual = PETSC_TRUE; 996 997 ierr = PetscNewLog(snes,&fas);CHKERRQ(ierr); 998 999 snes->data = (void*) fas; 1000 fas->level = 0; 1001 fas->levels = 1; 1002 fas->n_cycles = 1; 1003 fas->max_up_it = 1; 1004 fas->max_down_it = 1; 1005 fas->smoothu = NULL; 1006 fas->smoothd = NULL; 1007 fas->next = NULL; 1008 fas->previous = NULL; 1009 fas->fine = snes; 1010 fas->interpolate = NULL; 1011 fas->restrct = NULL; 1012 fas->inject = NULL; 1013 fas->usedmfornumberoflevels = PETSC_FALSE; 1014 fas->fastype = SNES_FAS_MULTIPLICATIVE; 1015 fas->full_downsweep = PETSC_FALSE; 1016 1017 fas->eventsmoothsetup = 0; 1018 fas->eventsmoothsolve = 0; 1019 fas->eventresidual = 0; 1020 fas->eventinterprestrict = 0; 1021 PetscFunctionReturn(0); 1022 } 1023