Actual source code: fas.c
1: #include <../src/snes/impls/fas/fasimpls.h>
3: const char *const SNESFASTypes[] = {"MULTIPLICATIVE", "ADDITIVE", "FULL", "KASKADE", "SNESFASType", "SNES_FAS", NULL};
5: static PetscErrorCode SNESReset_FAS(SNES snes)
6: {
7: SNES_FAS *fas = (SNES_FAS *)snes->data;
9: PetscFunctionBegin;
10: PetscCall(SNESDestroy(&fas->smoothu));
11: PetscCall(SNESDestroy(&fas->smoothd));
12: PetscCall(MatDestroy(&fas->inject));
13: PetscCall(MatDestroy(&fas->interpolate));
14: PetscCall(MatDestroy(&fas->restrct));
15: PetscCall(VecDestroy(&fas->rscale));
16: PetscCall(VecDestroy(&fas->Xg));
17: PetscCall(VecDestroy(&fas->Fg));
18: if (fas->coarseCorrectionLineSearch) PetscCall(SNESLineSearchReset(fas->coarseCorrectionLineSearch));
19: if (fas->next) PetscCall(SNESReset(fas->next));
20: PetscFunctionReturn(PETSC_SUCCESS);
21: }
23: static PetscErrorCode SNESDestroy_FAS(SNES snes)
24: {
25: SNES_FAS *fas = (SNES_FAS *)snes->data;
27: PetscFunctionBegin;
28: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESFASSetUseCoarseCorrectionLineSearch_C", NULL));
29: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESFASGetCoarseCorrectionLineSearch_C", NULL));
30: /* recursively resets and then destroys */
31: PetscCall(SNESReset_FAS(snes));
32: PetscCall(SNESLineSearchDestroy(&fas->coarseCorrectionLineSearch));
33: PetscCall(SNESDestroy(&fas->next));
34: PetscCall(PetscFree(fas));
35: PetscFunctionReturn(PETSC_SUCCESS);
36: }
38: static PetscErrorCode SNESFASSetUpLineSearch_Private(SNES snes, SNES smooth)
39: {
40: SNESLineSearch linesearch;
41: SNESLineSearch slinesearch;
42: void *lsprectx, *lspostctx;
43: PetscErrorCode (*precheck)(SNESLineSearch, Vec, Vec, PetscBool *, void *);
44: PetscErrorCode (*postcheck)(SNESLineSearch, Vec, Vec, Vec, PetscBool *, PetscBool *, void *);
46: PetscFunctionBegin;
47: if (!snes->linesearch) PetscFunctionReturn(PETSC_SUCCESS);
48: PetscCall(SNESGetLineSearch(snes, &linesearch));
49: PetscCall(SNESGetLineSearch(smooth, &slinesearch));
50: PetscCall(SNESLineSearchGetPreCheck(linesearch, &precheck, &lsprectx));
51: PetscCall(SNESLineSearchGetPostCheck(linesearch, &postcheck, &lspostctx));
52: PetscCall(SNESLineSearchSetPreCheck(slinesearch, precheck, lsprectx));
53: PetscCall(SNESLineSearchSetPostCheck(slinesearch, postcheck, lspostctx));
54: PetscCall(PetscObjectCopyFortranFunctionPointers((PetscObject)linesearch, (PetscObject)slinesearch));
55: PetscFunctionReturn(PETSC_SUCCESS);
56: }
58: static PetscErrorCode SNESFASCycleSetUpSmoother_Private(SNES snes, SNES smooth)
59: {
60: SNES_FAS *fas = (SNES_FAS *)snes->data;
62: PetscFunctionBegin;
63: PetscCall(PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)smooth));
64: PetscCall(SNESSetFromOptions(smooth));
65: PetscCall(SNESFASSetUpLineSearch_Private(snes, smooth));
67: PetscCall(PetscObjectReference((PetscObject)snes->vec_sol));
68: PetscCall(PetscObjectReference((PetscObject)snes->vec_sol_update));
69: PetscCall(PetscObjectReference((PetscObject)snes->vec_func));
70: smooth->vec_sol = snes->vec_sol;
71: smooth->vec_sol_update = snes->vec_sol_update;
72: smooth->vec_func = snes->vec_func;
74: if (fas->eventsmoothsetup) PetscCall(PetscLogEventBegin(fas->eventsmoothsetup, smooth, 0, 0, 0));
75: PetscCall(SNESSetUp(smooth));
76: if (fas->eventsmoothsetup) PetscCall(PetscLogEventEnd(fas->eventsmoothsetup, smooth, 0, 0, 0));
77: PetscFunctionReturn(PETSC_SUCCESS);
78: }
80: static PetscErrorCode SNESSetUp_FAS(SNES snes)
81: {
82: SNES_FAS *fas = (SNES_FAS *)snes->data;
83: PetscInt dm_levels;
84: SNES next;
85: PetscBool isFine, hasCreateRestriction, hasCreateInjection;
87: PetscFunctionBegin;
88: PetscCall(SNESFASCycleIsFine(snes, &isFine));
89: if (fas->usedmfornumberoflevels && isFine) {
90: PetscCall(DMGetRefineLevel(snes->dm, &dm_levels));
91: dm_levels++;
92: if (dm_levels > fas->levels) {
93: /* reset the number of levels */
94: PetscCall(SNESFASSetLevels(snes, dm_levels, NULL));
95: PetscCall(SNESSetFromOptions(snes));
96: }
97: }
98: PetscCall(SNESFASCycleGetCorrection(snes, &next));
99: if (!isFine) snes->gridsequence = 0; /* no grid sequencing inside the multigrid hierarchy! */
101: PetscCall(SNESSetWorkVecs(snes, 2)); /* work vectors used for intergrid transfers */
103: /* set up the smoothers if they haven't already been set up */
104: if (!fas->smoothd) PetscCall(SNESFASCycleCreateSmoother_Private(snes, &fas->smoothd));
106: if (snes->dm) {
107: /* set the smoother DMs properly */
108: if (fas->smoothu) PetscCall(SNESSetDM(fas->smoothu, snes->dm));
109: PetscCall(SNESSetDM(fas->smoothd, snes->dm));
110: /* construct EVERYTHING from the DM -- including the progressive set of smoothers */
111: if (next) {
112: /* for now -- assume the DM and the evaluation functions have been set externally */
113: if (!next->dm) {
114: PetscCall(DMCoarsen(snes->dm, PetscObjectComm((PetscObject)next), &next->dm));
115: PetscCall(SNESSetDM(next, next->dm));
116: }
117: /* set the interpolation and restriction from the DM */
118: if (!fas->interpolate) {
119: PetscCall(DMCreateInterpolation(next->dm, snes->dm, &fas->interpolate, &fas->rscale));
120: if (!fas->restrct) {
121: PetscCall(DMHasCreateRestriction(next->dm, &hasCreateRestriction));
122: /* DM can create restrictions, use that */
123: if (hasCreateRestriction) {
124: PetscCall(DMCreateRestriction(next->dm, snes->dm, &fas->restrct));
125: } else {
126: PetscCall(PetscObjectReference((PetscObject)fas->interpolate));
127: fas->restrct = fas->interpolate;
128: }
129: }
130: }
131: /* set the injection from the DM */
132: if (!fas->inject) {
133: PetscCall(DMHasCreateInjection(next->dm, &hasCreateInjection));
134: if (hasCreateInjection) PetscCall(DMCreateInjection(next->dm, snes->dm, &fas->inject));
135: }
136: }
137: }
139: /*pass the smoother, function, and jacobian up to the next level if it's not user set already */
140: if (fas->galerkin) {
141: if (next) PetscCall(SNESSetFunction(next, NULL, SNESFASGalerkinFunctionDefault, next));
142: if (fas->smoothd && fas->level != fas->levels - 1) PetscCall(SNESSetFunction(fas->smoothd, NULL, SNESFASGalerkinFunctionDefault, snes));
143: if (fas->smoothu && fas->level != fas->levels - 1) PetscCall(SNESSetFunction(fas->smoothu, NULL, SNESFASGalerkinFunctionDefault, snes));
144: }
146: /* sets the down (pre) smoother's default norm and sets it from options */
147: if (fas->smoothd) {
148: if (fas->level == 0) {
149: PetscCall(SNESSetNormSchedule(fas->smoothd, SNES_NORM_ALWAYS));
150: } else {
151: PetscCall(SNESSetNormSchedule(fas->smoothd, SNES_NORM_FINAL_ONLY));
152: }
153: PetscCall(SNESFASCycleSetUpSmoother_Private(snes, fas->smoothd));
154: }
156: /* sets the up (post) smoother's default norm and sets it from options */
157: if (fas->smoothu) {
158: if (fas->level != fas->levels - 1) {
159: PetscCall(SNESSetNormSchedule(fas->smoothu, SNES_NORM_NONE));
160: } else {
161: PetscCall(SNESSetNormSchedule(fas->smoothu, SNES_NORM_FINAL_ONLY));
162: }
163: PetscCall(SNESFASCycleSetUpSmoother_Private(snes, fas->smoothu));
164: }
166: if (next) {
167: /* gotta set up the solution vector for this to work */
168: if (!next->vec_sol) PetscCall(SNESFASCreateCoarseVec(snes, &next->vec_sol));
169: if (!next->vec_rhs) PetscCall(SNESFASCreateCoarseVec(snes, &next->vec_rhs));
170: PetscCall(PetscObjectCopyFortranFunctionPointers((PetscObject)snes, (PetscObject)next));
171: PetscCall(SNESFASSetUpLineSearch_Private(snes, next));
172: PetscCall(SNESSetUp(next));
173: }
175: /* setup FAS work vectors */
176: if (fas->galerkin) {
177: PetscCall(VecDuplicate(snes->vec_sol, &fas->Xg));
178: PetscCall(VecDuplicate(snes->vec_sol, &fas->Fg));
179: }
180: PetscFunctionReturn(PETSC_SUCCESS);
181: }
183: static PetscErrorCode SNESFASSetUpCoarseCorrectionLineSearch_Private(SNES snes)
184: {
185: SNES_FAS *fas = (SNES_FAS *)snes->data;
186: const char *optionsprefix;
188: PetscFunctionBegin;
189: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
190: for (PetscInt lvl = 1; lvl < fas->levels; lvl++) {
191: SNES levelsnes;
192: SNES_FAS *lfas;
193: SNESLineSearch ls;
194: char lsprefix[128];
196: PetscCall(SNESFASGetCycleSNES(snes, lvl, &levelsnes));
197: lfas = (SNES_FAS *)levelsnes->data;
198: if (!lfas->coarseCorrectionLineSearch) {
199: PetscCall(SNESLineSearchCreate(PetscObjectComm((PetscObject)levelsnes), &ls));
200: PetscCall(SNESLineSearchAppendOptionsPrefix(ls, optionsprefix));
201: PetscCall(PetscSNPrintf(lsprefix, sizeof(lsprefix), "fas_coarse_correction_%" PetscInt_FMT "_", lvl));
202: PetscCall(SNESLineSearchAppendOptionsPrefix(ls, lsprefix));
203: PetscCall(PetscObjectIncrementTabLevel((PetscObject)ls, (PetscObject)levelsnes, 1));
204: PetscCall(SNESLineSearchSetSNES(ls, levelsnes));
205: PetscCall(SNESLineSearchSetType(ls, SNESLINESEARCHSECANT));
206: lfas->coarseCorrectionLineSearch = ls;
207: }
208: }
209: PetscFunctionReturn(PETSC_SUCCESS);
210: }
212: static PetscErrorCode SNESSetFromOptions_FAS(SNES snes, PetscOptionItems PetscOptionsObject)
213: {
214: SNES_FAS *fas = (SNES_FAS *)snes->data;
215: PetscInt levels = 1;
216: PetscBool flg = PETSC_FALSE, upflg = PETSC_FALSE, downflg = PETSC_FALSE, monflg = PETSC_FALSE, galerkinflg = PETSC_FALSE, continuationflg = PETSC_FALSE;
217: SNESFASType fastype;
218: const char *optionsprefix;
219: SNESLineSearch linesearch = NULL;
220: PetscInt m, n_up, n_down;
221: SNES next;
222: PetscBool isFine;
224: PetscFunctionBegin;
225: PetscCall(SNESFASCycleIsFine(snes, &isFine));
226: PetscOptionsHeadBegin(PetscOptionsObject, "SNESFAS Options-----------------------------------");
228: /* number of levels -- only process most options on the finest level */
229: if (isFine) {
230: PetscCall(PetscOptionsInt("-snes_fas_levels", "Number of Levels", "SNESFASSetLevels", levels, &levels, &flg));
231: if (!flg && snes->dm) {
232: PetscCall(DMGetRefineLevel(snes->dm, &levels));
233: levels++;
234: fas->usedmfornumberoflevels = PETSC_TRUE;
235: }
236: PetscCall(SNESFASSetLevels(snes, levels, NULL));
237: fastype = fas->fastype;
238: PetscCall(PetscOptionsEnum("-snes_fas_type", "FAS correction type", "SNESFASSetType", SNESFASTypes, (PetscEnum)fastype, (PetscEnum *)&fastype, &flg));
239: if (flg) PetscCall(SNESFASSetType(snes, fastype));
241: PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
242: PetscCall(PetscOptionsInt("-snes_fas_cycles", "Number of cycles", "SNESFASSetCycles", fas->n_cycles, &m, &flg));
243: if (flg) PetscCall(SNESFASSetCycles(snes, m));
244: PetscCall(PetscOptionsBool("-snes_fas_continuation", "Corrected grid-sequence continuation", "SNESFASSetContinuation", fas->continuation, &continuationflg, &flg));
245: if (flg) PetscCall(SNESFASSetContinuation(snes, continuationflg));
247: PetscCall(PetscOptionsBool("-snes_fas_galerkin", "Form coarse problems with Galerkin", "SNESFASSetGalerkin", fas->galerkin, &galerkinflg, &flg));
248: if (flg) PetscCall(SNESFASSetGalerkin(snes, galerkinflg));
250: if (fas->fastype == SNES_FAS_FULL) {
251: PetscCall(PetscOptionsBool("-snes_fas_full_downsweep", "Smooth on the initial down sweep for full FAS cycles", "SNESFASFullSetDownSweep", fas->full_downsweep, &fas->full_downsweep, &flg));
252: if (flg) PetscCall(SNESFASFullSetDownSweep(snes, fas->full_downsweep));
253: PetscCall(PetscOptionsBool("-snes_fas_full_total", "Use total restriction and interpolaton on the indial down and up sweeps for the full FAS cycle", "SNESFASFullSetUseTotal", fas->full_total, &fas->full_total, &flg));
254: if (flg) PetscCall(SNESFASFullSetTotal(snes, fas->full_total));
255: }
257: PetscCall(PetscOptionsInt("-snes_fas_smoothup", "Number of post-smoothing steps", "SNESFASSetNumberSmoothUp", fas->max_up_it, &n_up, &upflg));
259: PetscCall(PetscOptionsInt("-snes_fas_smoothdown", "Number of pre-smoothing steps", "SNESFASSetNumberSmoothDown", fas->max_down_it, &n_down, &downflg));
261: {
262: PetscViewer viewer;
263: PetscViewerFormat format;
264: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)snes), ((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_fas_monitor", &viewer, &format, &monflg));
265: if (monflg) {
266: PetscViewerAndFormat *vf;
267: PetscCall(PetscViewerAndFormatCreate(viewer, format, &vf));
268: PetscCall(PetscViewerDestroy(&viewer));
269: PetscCall(SNESFASSetMonitor(snes, vf, PETSC_TRUE));
270: }
271: }
272: flg = PETSC_FALSE;
273: monflg = PETSC_TRUE;
274: PetscCall(PetscOptionsBool("-snes_fas_log", "Log times for each FAS level", "SNESFASSetLog", monflg, &monflg, &flg));
275: if (flg) PetscCall(SNESFASSetLog(snes, monflg));
276: }
277: PetscCall(PetscOptionsBool("-snes_fas_monitor_correction", "View the tau correction at each iteration", "SNESFASCoarseCorrection", fas->monitorCorrection, &fas->monitorCorrection, &flg));
279: /* setup from the determined types if there is no pointwise procedure or smoother defined */
280: if (upflg) PetscCall(SNESFASSetNumberSmoothUp(snes, n_up));
281: if (downflg) PetscCall(SNESFASSetNumberSmoothDown(snes, n_down));
283: /* set up the default line search for coarse grid corrections */
284: if (fas->fastype == SNES_FAS_ADDITIVE) {
285: PetscCall(SNESGetLineSearch(snes, &linesearch));
286: if (!((PetscObject)linesearch)->type_name) PetscCall(SNESLineSearchSetType(linesearch, SNESLINESEARCHSECANT));
287: }
289: if (isFine)
290: PetscCall(PetscOptionsBool("-snes_fas_use_coarse_correction_linesearch", "Use a line search for the multiplicative coarse correction", "SNESFASSetUseCoarseCorrectionLineSearch", fas->useCoarseCorrectionLineSearch, &fas->useCoarseCorrectionLineSearch, &flg));
291: PetscOptionsHeadEnd();
293: if (isFine && fas->fastype == SNES_FAS_MULTIPLICATIVE && fas->useCoarseCorrectionLineSearch) {
294: PetscCall(SNESFASSetUpCoarseCorrectionLineSearch_Private(snes));
296: for (PetscInt lvl = 1; lvl < fas->levels; lvl++) {
297: SNES levelsnes;
298: SNES_FAS *lfas;
300: PetscCall(SNESFASGetCycleSNES(snes, lvl, &levelsnes));
301: lfas = (SNES_FAS *)levelsnes->data;
302: PetscCall(SNESLineSearchSetFromOptions(lfas->coarseCorrectionLineSearch));
303: }
304: }
306: /* recursive option setting for the smoothers */
307: PetscCall(SNESFASCycleGetCorrection(snes, &next));
308: if (next) PetscCall(SNESSetFromOptions(next));
309: PetscFunctionReturn(PETSC_SUCCESS);
310: }
312: #include <petscdraw.h>
313: static PetscErrorCode SNESView_FAS(SNES snes, PetscViewer viewer)
314: {
315: SNES_FAS *fas = (SNES_FAS *)snes->data;
316: SNES_FAS *levelfas;
317: PetscBool isFine, isascii, isdraw;
318: SNES smoothu, smoothd, levelsnes;
320: PetscFunctionBegin;
321: PetscCall(SNESFASCycleIsFine(snes, &isFine));
322: if (isFine) {
323: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
324: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
325: if (isascii) {
326: PetscCall(PetscViewerASCIIPrintf(viewer, " type is %s, levels=%" PetscInt_FMT ", cycles=%" PetscInt_FMT "\n", SNESFASTypes[fas->fastype], fas->levels, fas->n_cycles));
327: if (fas->galerkin) {
328: PetscCall(PetscViewerASCIIPrintf(viewer, " Using Galerkin computed coarse grid function evaluation\n"));
329: } else {
330: PetscCall(PetscViewerASCIIPrintf(viewer, " Not using Galerkin computed coarse grid function evaluation\n"));
331: }
332: for (PetscInt i = 0; i < fas->levels; i++) {
333: PetscCall(SNESFASGetCycleSNES(snes, i, &levelsnes));
334: levelfas = (SNES_FAS *)levelsnes->data;
335: PetscCall(SNESFASCycleGetSmootherUp(levelsnes, &smoothu));
336: PetscCall(SNESFASCycleGetSmootherDown(levelsnes, &smoothd));
337: if (!i) {
338: PetscCall(PetscViewerASCIIPrintf(viewer, " Coarse grid solver -- level %" PetscInt_FMT " -------------------------------\n", i));
339: } else {
340: PetscCall(PetscViewerASCIIPrintf(viewer, " Down solver (pre-smoother) on level %" PetscInt_FMT " -------------------------------\n", i));
341: }
342: PetscCall(PetscViewerASCIIPushTab(viewer));
343: if (smoothd) {
344: PetscCall(SNESView(smoothd, viewer));
345: } else {
346: PetscCall(PetscViewerASCIIPrintf(viewer, "Not yet available\n"));
347: }
348: PetscCall(PetscViewerASCIIPopTab(viewer));
349: if (i && (smoothd == smoothu)) {
350: PetscCall(PetscViewerASCIIPrintf(viewer, " Up solver (post-smoother) same as down solver (pre-smoother)\n"));
351: } else if (i) {
352: PetscCall(PetscViewerASCIIPrintf(viewer, " Up solver (post-smoother) on level %" PetscInt_FMT " -------------------------------\n", i));
353: PetscCall(PetscViewerASCIIPushTab(viewer));
354: if (smoothu) {
355: PetscCall(SNESView(smoothu, viewer));
356: } else {
357: PetscCall(PetscViewerASCIIPrintf(viewer, "Not yet available\n"));
358: }
359: PetscCall(PetscViewerASCIIPopTab(viewer));
360: }
361: if (i && levelfas->coarseCorrectionLineSearch) {
362: PetscCall(PetscViewerASCIIPrintf(viewer, " Coarse correction line search on level %" PetscInt_FMT " -------------------------------\n", i));
363: PetscCall(PetscViewerASCIIPushTab(viewer));
364: PetscCall(SNESLineSearchView(levelfas->coarseCorrectionLineSearch, viewer));
365: PetscCall(PetscViewerASCIIPopTab(viewer));
366: }
367: }
368: } else if (isdraw) {
369: PetscDraw draw;
370: PetscReal x, w, y, bottom, th, wth;
371: SNES_FAS *curfas = fas;
372: PetscCall(PetscViewerDrawGetDraw(viewer, 0, &draw));
373: PetscCall(PetscDrawGetCurrentPoint(draw, &x, &y));
374: PetscCall(PetscDrawStringGetSize(draw, &wth, &th));
375: bottom = y - th;
376: while (curfas) {
377: if (!curfas->smoothu) {
378: PetscCall(PetscDrawPushCurrentPoint(draw, x, bottom));
379: if (curfas->smoothd) PetscCall(SNESView(curfas->smoothd, viewer));
380: PetscCall(PetscDrawPopCurrentPoint(draw));
381: } else {
382: w = 0.5 * PetscMin(1.0 - x, x);
383: PetscCall(PetscDrawPushCurrentPoint(draw, x - w, bottom));
384: if (curfas->smoothd) PetscCall(SNESView(curfas->smoothd, viewer));
385: PetscCall(PetscDrawPopCurrentPoint(draw));
386: PetscCall(PetscDrawPushCurrentPoint(draw, x + w, bottom));
387: if (curfas->smoothu) PetscCall(SNESView(curfas->smoothu, viewer));
388: PetscCall(PetscDrawPopCurrentPoint(draw));
389: }
390: /* this is totally bogus but we have no way of knowing how low the previous one was draw to */
391: bottom -= 5 * th;
392: if (curfas->next) curfas = (SNES_FAS *)curfas->next->data;
393: else curfas = NULL;
394: }
395: }
396: }
397: PetscFunctionReturn(PETSC_SUCCESS);
398: }
400: /*
401: Defines the action of the downsmoother
402: */
403: static PetscErrorCode SNESFASDownSmooth_Private(SNES snes, Vec B, Vec X, Vec F, PetscReal *fnorm)
404: {
405: SNESConvergedReason reason;
406: Vec FPC;
407: SNES smoothd;
408: PetscBool flg;
409: SNES_FAS *fas = (SNES_FAS *)snes->data;
411: PetscFunctionBegin;
412: PetscCall(SNESFASCycleGetSmootherDown(snes, &smoothd));
413: PetscCall(SNESSetInitialFunction(smoothd, F));
414: if (fas->eventsmoothsolve) PetscCall(PetscLogEventBegin(fas->eventsmoothsolve, smoothd, B, X, 0));
415: PetscCall(SNESSolve(smoothd, B, X));
416: if (fas->eventsmoothsolve) PetscCall(PetscLogEventEnd(fas->eventsmoothsolve, smoothd, B, X, 0));
417: /* check convergence reason for the smoother */
418: PetscCall(SNESGetConvergedReason(smoothd, &reason));
419: if (reason < 0 && !(reason == SNES_DIVERGED_MAX_IT || reason == SNES_DIVERGED_LOCAL_MIN || reason == SNES_DIVERGED_LINE_SEARCH)) {
420: snes->reason = SNES_DIVERGED_INNER;
421: PetscFunctionReturn(PETSC_SUCCESS);
422: }
424: PetscCall(SNESGetFunction(smoothd, &FPC, NULL, NULL));
425: PetscCall(SNESGetAlwaysComputesFinalResidual(smoothd, &flg));
426: if (!flg) PetscCall(SNESComputeFunction(smoothd, X, FPC));
427: PetscCall(VecCopy(FPC, F));
428: if (fnorm) {
429: PetscCall(VecNorm(F, NORM_2, fnorm));
430: SNESCheckFunctionDomainError(snes, *fnorm);
431: }
432: PetscFunctionReturn(PETSC_SUCCESS);
433: }
435: /*
436: Defines the action of the upsmoother
437: */
438: static PetscErrorCode SNESFASUpSmooth_Private(SNES snes, Vec B, Vec X, Vec F, PetscReal *fnorm)
439: {
440: SNESConvergedReason reason;
441: Vec FPC;
442: SNES smoothu;
443: PetscBool flg;
444: SNES_FAS *fas = (SNES_FAS *)snes->data;
446: PetscFunctionBegin;
447: PetscCall(SNESFASCycleGetSmootherUp(snes, &smoothu));
448: if (fas->eventsmoothsolve) PetscCall(PetscLogEventBegin(fas->eventsmoothsolve, smoothu, 0, 0, 0));
449: PetscCall(SNESSolve(smoothu, B, X));
450: if (fas->eventsmoothsolve) PetscCall(PetscLogEventEnd(fas->eventsmoothsolve, smoothu, 0, 0, 0));
451: /* check convergence reason for the smoother */
452: PetscCall(SNESGetConvergedReason(smoothu, &reason));
453: if (reason < 0 && !(reason == SNES_DIVERGED_MAX_IT || reason == SNES_DIVERGED_LOCAL_MIN || reason == SNES_DIVERGED_LINE_SEARCH)) {
454: snes->reason = SNES_DIVERGED_INNER;
455: PetscFunctionReturn(PETSC_SUCCESS);
456: }
457: PetscCall(SNESGetFunction(smoothu, &FPC, NULL, NULL));
458: PetscCall(SNESGetAlwaysComputesFinalResidual(smoothu, &flg));
459: if (!flg) PetscCall(SNESComputeFunction(smoothu, X, FPC));
460: PetscCall(VecCopy(FPC, F));
461: if (fnorm) {
462: PetscCall(VecNorm(F, NORM_2, fnorm));
463: SNESCheckFunctionDomainError(snes, *fnorm);
464: }
465: PetscFunctionReturn(PETSC_SUCCESS);
466: }
468: /*@
469: SNESFASCreateCoarseVec - create a `Vec` corresponding to a state vector on one level coarser than the current level
471: Collective
473: Input Parameter:
474: . snes - `SNESFAS` object
476: Output Parameter:
477: . Xcoarse - vector on level one coarser than the current level
479: Level: developer
481: .seealso: [](ch_snes), `SNESFASSetRestriction()`, `SNESFASRestrict()`, `SNESFAS`
482: @*/
483: PetscErrorCode SNESFASCreateCoarseVec(SNES snes, Vec *Xcoarse)
484: {
485: SNES_FAS *fas;
487: PetscFunctionBegin;
489: PetscAssertPointer(Xcoarse, 2);
490: fas = (SNES_FAS *)snes->data;
491: if (fas->rscale) {
492: PetscCall(VecDuplicate(fas->rscale, Xcoarse));
493: } else if (fas->interpolate) {
494: PetscCall(MatCreateVecs(fas->interpolate, Xcoarse, NULL));
495: } else SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Must set rscale or interpolation");
496: PetscFunctionReturn(PETSC_SUCCESS);
497: }
499: /*@
500: SNESFASRestrict - restrict a `Vec` to the next coarser level
502: Collective
504: Input Parameters:
505: + fine - `SNES` from which to restrict
506: - Xfine - vector to restrict
508: Output Parameter:
509: . Xcoarse - result of restriction
511: Level: developer
513: .seealso: [](ch_snes), `SNES`, `SNESFAS`, `SNESFASSetRestriction()`, `SNESFASSetInjection()`, `SNESFASCreateCoarseVec()`
514: @*/
515: PetscErrorCode SNESFASRestrict(SNES fine, Vec Xfine, Vec Xcoarse)
516: {
517: SNES_FAS *fas;
519: PetscFunctionBegin;
523: fas = (SNES_FAS *)fine->data;
524: if (fas->inject) {
525: PetscCall(MatRestrict(fas->inject, Xfine, Xcoarse));
526: } else {
527: PetscCall(MatRestrict(fas->restrct, Xfine, Xcoarse));
528: PetscCall(VecPointwiseMult(Xcoarse, fas->rscale, Xcoarse));
529: }
530: PetscFunctionReturn(PETSC_SUCCESS);
531: }
533: /*
534: Performs a variant of FAS using the interpolated total coarse solution
536: fine problem: F(x) = b
537: coarse problem: F^c(x^c) = Rb, Initial guess Rx
538: interpolated solution: x^f = I x^c (total solution interpolation
539: */
540: static PetscErrorCode SNESFASInterpolatedCoarseSolution(SNES snes, Vec X, Vec X_new)
541: {
542: Vec X_c, B_c;
543: SNESConvergedReason reason;
544: SNES next;
545: Mat restrct, interpolate;
546: SNES_FAS *fasc;
548: PetscFunctionBegin;
549: PetscCall(SNESFASCycleGetCorrection(snes, &next));
550: if (next) {
551: fasc = (SNES_FAS *)next->data;
553: PetscCall(SNESFASCycleGetRestriction(snes, &restrct));
554: PetscCall(SNESFASCycleGetInterpolation(snes, &interpolate));
556: X_c = next->vec_sol;
558: if (fasc->eventinterprestrict) PetscCall(PetscLogEventBegin(fasc->eventinterprestrict, snes, 0, 0, 0));
559: /* restrict the total solution: Rb */
560: PetscCall(SNESFASRestrict(snes, X, X_c));
561: B_c = next->vec_rhs;
562: if (snes->vec_rhs) {
563: /* restrict the total rhs defect: Rb */
564: PetscCall(MatRestrict(restrct, snes->vec_rhs, B_c));
565: } else {
566: PetscCall(VecSet(B_c, 0.));
567: }
568: if (fasc->eventinterprestrict) PetscCall(PetscLogEventEnd(fasc->eventinterprestrict, snes, 0, 0, 0));
570: PetscCall(SNESSolve(next, B_c, X_c));
571: PetscCall(SNESGetConvergedReason(next, &reason));
572: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
573: snes->reason = SNES_DIVERGED_INNER;
574: PetscFunctionReturn(PETSC_SUCCESS);
575: }
576: /* x^f <- Ix^c*/
577: DM dmc, dmf;
579: PetscCall(SNESGetDM(next, &dmc));
580: PetscCall(SNESGetDM(snes, &dmf));
581: if (fasc->eventinterprestrict) PetscCall(PetscLogEventBegin(fasc->eventinterprestrict, snes, 0, 0, 0));
582: PetscCall(DMInterpolateSolution(dmc, dmf, interpolate, X_c, X_new));
583: if (fasc->eventinterprestrict) PetscCall(PetscLogEventEnd(fasc->eventinterprestrict, snes, 0, 0, 0));
584: PetscCall(PetscObjectSetName((PetscObject)X_c, "Coarse solution"));
585: PetscCall(VecViewFromOptions(X_c, NULL, "-fas_coarse_solution_view"));
586: PetscCall(PetscObjectSetName((PetscObject)X_new, "Updated Fine solution"));
587: PetscCall(VecViewFromOptions(X_new, NULL, "-fas_levels_1_solution_view"));
588: }
589: PetscFunctionReturn(PETSC_SUCCESS);
590: }
592: /*
593: Performs the FAS coarse correction as:
595: fine problem: F(x) = b
596: coarse problem: F^c(x^c) = b^c
598: b^c = F^c(Rx) - R(F(x) - b)
599: = tau + R b
600: */
601: static PetscErrorCode SNESFASCoarseCorrection(SNES snes, Vec X, Vec F)
602: {
603: SNES_FAS *fas = (SNES_FAS *)snes->data;
604: SNES_FAS *fasc;
605: PetscBool monitorCorrection = fas->monitorCorrection;
606: PetscBool uselinesearch;
607: PetscReal xonorm = 0.0;
608: Vec X_c, Xo_c, F_c, B_c, Xhat;
609: SNESConvergedReason reason;
610: SNES next;
611: Mat restrct, interpolate;
613: PetscFunctionBegin;
614: PetscCall(SNESFASCycleGetCorrection(snes, &next));
615: if (next) {
616: fasc = (SNES_FAS *)next->data;
618: PetscCall(SNESFASCycleGetRestriction(snes, &restrct));
619: PetscCall(SNESFASCycleGetInterpolation(snes, &interpolate));
621: X_c = next->vec_sol;
622: Xo_c = next->work[0];
623: F_c = next->vec_func;
624: B_c = next->vec_rhs;
625: Xhat = snes->work[1];
627: if (fasc->eventinterprestrict) PetscCall(PetscLogEventBegin(fasc->eventinterprestrict, snes, 0, 0, 0));
628: PetscCall(SNESFASRestrict(snes, X, Xo_c));
629: // restrict the defect: R(F(x) - b)
630: PetscCall(MatRestrict(restrct, F, B_c));
631: if (monitorCorrection) {
632: PetscReal fnorm, rnorm;
634: PetscCall(VecNorm(F, NORM_2, &fnorm));
635: PetscCall(VecNorm(B_c, NORM_2, &rnorm));
636: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), "||F(X)|| %g\n||R (F(X) - b)|| %g\n", (double)fnorm, (double)rnorm));
637: }
638: if (fasc->eventinterprestrict) PetscCall(PetscLogEventEnd(fasc->eventinterprestrict, snes, 0, 0, 0));
640: if (fasc->eventresidual) PetscCall(PetscLogEventBegin(fasc->eventresidual, next, 0, 0, 0));
641: // F_c = F^c(Rx) - R(F(x) - b) since the second term was sitting in next->vec_rhs
642: PetscCall(SNESComputeFunction(next, Xo_c, F_c));
643: if (monitorCorrection) {
644: PetscReal xnorm, fnorm;
646: PetscCall(VecNorm(Xo_c, NORM_2, &xnorm));
647: PetscCall(VecNorm(F_c, NORM_2, &fnorm));
648: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), "||R(X)|| %g\n||tau + Rb = F_c (R X) - R(F(X) - b)|| %g\n", (double)xnorm, (double)fnorm));
649: PetscCall(PetscObjectSetName((PetscObject)F_c, "tau"));
650: PetscCall(VecViewFromOptions(F_c, (PetscObject)next, "-tau_view"));
651: }
652: if (fasc->eventresidual) PetscCall(PetscLogEventEnd(fasc->eventresidual, next, 0, 0, 0));
654: // solve the coarse problem corresponding to F^c(x^c) = b^c = F^c(Rx) - R(F(x) - b)
655: PetscCall(VecCopy(B_c, X_c)); // Now x^c holds R(F(x) - b)
656: PetscCall(VecCopy(F_c, B_c)); // Now b^c holds F^c(Rx) - R(F(x) - b)
657: PetscCall(VecCopy(X_c, F_c)); // Now F^c holds R(F(x) - b)
658: // set initial guess of the coarse problem to the projected fine solution
659: PetscCall(VecCopy(Xo_c, X_c));
661: /*
662: recurse to the next level
663: So x^c_0 = R x and F^c_0 = F^c(x^c_0) - b^c = F^c(x^c_0) - (F^c(x^c_0) - R(F(x) - b)) = R(F(x) - b)
664: */
665: PetscCall(SNESSetInitialFunction(next, F_c));
666: PetscCall(SNESSolve(next, B_c, X_c));
667: PetscCall(SNESGetConvergedReason(next, &reason));
668: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
669: snes->reason = SNES_DIVERGED_INNER;
670: PetscFunctionReturn(PETSC_SUCCESS);
671: }
672: // correct as x <- x + I(x^c - Rx), optionally via the coarse correction line search
673: uselinesearch = (PetscBool)(fas->fastype == SNES_FAS_MULTIPLICATIVE && fas->coarseCorrectionLineSearch != NULL);
674: if (monitorCorrection) PetscCall(VecNorm(X, NORM_2, &xonorm));
675: PetscCall(VecAXPY(X_c, -1.0, Xo_c));
676: if (uselinesearch) {
677: if (fasc->eventinterprestrict) PetscCall(PetscLogEventBegin(fasc->eventinterprestrict, snes, 0, 0, 0));
678: PetscCall(MatInterpolate(interpolate, X_c, Xhat));
679: if (fasc->eventinterprestrict) PetscCall(PetscLogEventEnd(fasc->eventinterprestrict, snes, 0, 0, 0));
680: PetscCall(VecScale(Xhat, -1.0)); // SNESLineSearch expects the descent direction as -Y, so negate the correction
681: PetscCall(SNESLineSearchApply(fas->coarseCorrectionLineSearch, X, F, &snes->norm, Xhat));
682: SNESCheckLineSearchFailure(snes, fas->coarseCorrectionLineSearch);
683: if (snes->reason < 0) PetscFunctionReturn(PETSC_SUCCESS);
684: } else {
685: if (fasc->eventinterprestrict) PetscCall(PetscLogEventBegin(fasc->eventinterprestrict, snes, 0, 0, 0));
686: PetscCall(MatInterpolateAdd(interpolate, X_c, X, X));
687: if (fasc->eventinterprestrict) PetscCall(PetscLogEventEnd(fasc->eventinterprestrict, snes, 0, 0, 0));
688: }
689: if (monitorCorrection) {
690: PetscReal xnorm, inorm;
692: PetscCall(VecNorm(X_c, NORM_2, &xnorm));
693: PetscCall(VecNorm(X, NORM_2, &inorm));
694: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)snes), "||X_c - Xo_c|| %g\n||X|| %g\n||X + lambda (X_c - X_co)|| %g\n", (double)xnorm, (double)xonorm, (double)inorm));
695: }
696: // TODO Check for snes->b, Technically we should strip out R b here if it is nonzero
697: PetscCall(PetscObjectSetName((PetscObject)B_c, "Tau correction"));
698: PetscCall(VecViewFromOptions(B_c, NULL, "-fas_tau_correction_view"));
699: PetscCall(PetscObjectSetName((PetscObject)X_c, "Coarse correction"));
700: PetscCall(VecViewFromOptions(X_c, NULL, "-fas_coarse_solution_view"));
701: PetscCall(PetscObjectSetName((PetscObject)X, "Updated Fine solution"));
702: PetscCall(VecViewFromOptions(X, NULL, "-fas_levels_1_solution_view"));
703: }
704: PetscFunctionReturn(PETSC_SUCCESS);
705: }
707: /*
708: The additive cycle is:
710: xhat = x
711: xhat = dS(x, b)
712: x = coarsecorrection(xhat, b_d)
713: x = x + nu*(xhat - x);
714: (optional) x = uS(x, b)
716: With the coarse RHS (defect correction) as below.
717: */
718: static PetscErrorCode SNESFASCycle_Additive(SNES snes, Vec X)
719: {
720: Vec F, B, Xhat;
721: Vec X_c, Xo_c, F_c, B_c;
722: SNESConvergedReason reason;
723: PetscReal xnorm, fnorm, ynorm;
724: SNES next;
725: Mat restrct, interpolate;
726: SNES_FAS *fas = (SNES_FAS *)snes->data, *fasc;
728: PetscFunctionBegin;
729: PetscCall(SNESFASCycleGetCorrection(snes, &next));
730: F = snes->vec_func;
731: B = snes->vec_rhs;
732: Xhat = snes->work[1];
733: PetscCall(VecCopy(X, Xhat));
734: /* recurse first */
735: if (next) {
736: fasc = (SNES_FAS *)next->data;
737: PetscCall(SNESFASCycleGetRestriction(snes, &restrct));
738: PetscCall(SNESFASCycleGetInterpolation(snes, &interpolate));
739: if (fas->eventresidual) PetscCall(PetscLogEventBegin(fas->eventresidual, snes, 0, 0, 0));
740: PetscCall(SNESComputeFunction(snes, Xhat, F));
741: if (fas->eventresidual) PetscCall(PetscLogEventEnd(fas->eventresidual, snes, 0, 0, 0));
742: PetscCall(VecNorm(F, NORM_2, &fnorm));
743: SNESCheckFunctionDomainError(snes, fnorm);
744: X_c = next->vec_sol;
745: Xo_c = next->work[0];
746: F_c = next->vec_func;
747: B_c = next->vec_rhs;
749: PetscCall(SNESFASRestrict(snes, Xhat, Xo_c));
750: /* restrict the defect */
751: PetscCall(MatRestrict(restrct, F, B_c));
753: /* solve the coarse problem corresponding to F^c(x^c) = b^c = Rb + F^c(Rx) - RF(x) */
754: if (fasc->eventresidual) PetscCall(PetscLogEventBegin(fasc->eventresidual, next, 0, 0, 0));
755: PetscCall(SNESComputeFunction(next, Xo_c, F_c));
756: if (fasc->eventresidual) PetscCall(PetscLogEventEnd(fasc->eventresidual, next, 0, 0, 0));
757: PetscCall(VecCopy(B_c, X_c));
758: PetscCall(VecCopy(F_c, B_c));
759: PetscCall(VecCopy(X_c, F_c));
760: /* set initial guess of the coarse problem to the projected fine solution */
761: PetscCall(VecCopy(Xo_c, X_c));
763: /* recurse */
764: PetscCall(SNESSetInitialFunction(next, F_c));
765: PetscCall(SNESSolve(next, B_c, X_c));
767: /* smooth on this level */
768: PetscCall(SNESFASDownSmooth_Private(snes, B, X, F, &fnorm));
770: PetscCall(SNESGetConvergedReason(next, &reason));
771: if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) {
772: snes->reason = SNES_DIVERGED_INNER;
773: PetscFunctionReturn(PETSC_SUCCESS);
774: }
776: /* correct as x <- x + I(x^c - Rx)*/
777: PetscCall(VecAYPX(X_c, -1.0, Xo_c));
778: PetscCall(MatInterpolate(interpolate, X_c, Xhat));
780: /* additive correction of the coarse direction*/
781: PetscCall(SNESLineSearchApply(snes->linesearch, X, F, &fnorm, Xhat));
782: PetscCall(SNESLineSearchGetNorms(snes->linesearch, &xnorm, &snes->norm, &ynorm));
783: SNESCheckLineSearchFailure(snes, snes->linesearch);
784: } else {
785: PetscCall(SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm));
786: }
787: PetscFunctionReturn(PETSC_SUCCESS);
788: }
790: /*
791: Defines the FAS cycle as:
793: fine problem: F(x) = b
794: coarse problem: F^c(x) = b^c
796: b^c = F^c(Rx) - R(F(x) - b)
797: = F^c(Rx) - R(F(x)) + R b
798: = tau + R b
800: correction:
802: x = x + I(x^c - Rx)
803: */
804: static PetscErrorCode SNESFASCycle_Multiplicative(SNES snes, Vec X)
805: {
806: Vec F, B;
807: SNES next;
809: PetscFunctionBegin;
810: F = snes->vec_func;
811: B = snes->vec_rhs;
812: /* pre-smooth -- just update using the pre-smoother */
813: PetscCall(SNESFASCycleGetCorrection(snes, &next));
814: PetscCall(SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm));
815: if (next) {
816: PetscCall(SNESFASCoarseCorrection(snes, X, F));
817: PetscCall(SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm));
818: }
819: PetscFunctionReturn(PETSC_SUCCESS);
820: }
822: static PetscErrorCode SNESFASCycleSetupPhase_Full(SNES snes)
823: {
824: SNES next;
825: SNES_FAS *fas = (SNES_FAS *)snes->data;
826: PetscBool isFine;
828: PetscFunctionBegin;
829: /* pre-smooth -- just update using the pre-smoother */
830: PetscCall(SNESFASCycleIsFine(snes, &isFine));
831: PetscCall(SNESFASCycleGetCorrection(snes, &next));
832: fas->full_stage = 0;
833: if (next) PetscCall(SNESFASCycleSetupPhase_Full(next));
834: PetscFunctionReturn(PETSC_SUCCESS);
835: }
837: static PetscErrorCode SNESFASCycle_Full(SNES snes, Vec X)
838: {
839: Vec F, B;
840: SNES_FAS *fas = (SNES_FAS *)snes->data;
841: PetscBool isFine;
842: SNES next;
844: PetscFunctionBegin;
845: F = snes->vec_func;
846: B = snes->vec_rhs;
847: PetscCall(SNESFASCycleIsFine(snes, &isFine));
848: PetscCall(SNESFASCycleGetCorrection(snes, &next));
850: if (isFine) PetscCall(SNESFASCycleSetupPhase_Full(snes));
852: if (fas->full_stage == 0) {
853: /* downsweep */
854: if (next) {
855: if (fas->level != 1) next->max_its += 1;
856: if (fas->full_downsweep) PetscCall(SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm));
857: fas->full_downsweep = PETSC_TRUE;
858: if (fas->full_total) PetscCall(SNESFASInterpolatedCoarseSolution(snes, X, X));
859: else PetscCall(SNESFASCoarseCorrection(snes, X, F));
860: fas->full_total = PETSC_FALSE;
861: PetscCall(SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm));
862: if (fas->level != 1) next->max_its -= 1;
863: } else {
864: /* The smoother on the coarse level is the coarse solver */
865: PetscCall(SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm));
866: }
867: fas->full_stage = 1;
868: } else if (fas->full_stage == 1) {
869: if (snes->iter == 0) PetscCall(SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm));
870: if (next) {
871: PetscCall(SNESFASCoarseCorrection(snes, X, F));
872: PetscCall(SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm));
873: }
874: }
875: /* final v-cycle */
876: if (isFine) {
877: if (next) {
878: PetscCall(SNESFASCoarseCorrection(snes, X, F));
879: PetscCall(SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm));
880: }
881: }
882: PetscFunctionReturn(PETSC_SUCCESS);
883: }
885: static PetscErrorCode SNESFASCycle_Kaskade(SNES snes, Vec X)
886: {
887: Vec F, B;
888: SNES next;
890: PetscFunctionBegin;
891: F = snes->vec_func;
892: B = snes->vec_rhs;
893: PetscCall(SNESFASCycleGetCorrection(snes, &next));
894: if (next) {
895: PetscCall(SNESFASCoarseCorrection(snes, X, F));
896: PetscCall(SNESFASUpSmooth_Private(snes, B, X, F, &snes->norm));
897: } else {
898: PetscCall(SNESFASDownSmooth_Private(snes, B, X, F, &snes->norm));
899: }
900: PetscFunctionReturn(PETSC_SUCCESS);
901: }
903: PetscBool SNEScite = PETSC_FALSE;
904: const char SNESCitation[] = "@Article{bruneknepleysmithtu15,"
905: " title = {Composing Scalable Nonlinear Algebraic Solvers},"
906: " author = {Peter R. Brune and Matthew G. Knepley and Barry F. Smith and Xuemin Tu},"
907: " journal = {SIAM Review},"
908: " volume = {57},"
909: " number = {4},"
910: " pages = {535--565},"
911: " doi = {10.1137/130936725},"
912: " year = {2015}\n";
914: static PetscErrorCode SNESSolve_FAS(SNES snes)
915: {
916: Vec X, F;
917: PetscReal fnorm;
918: SNES_FAS *fas = (SNES_FAS *)snes->data, *ffas;
919: DM dm;
920: PetscBool isFine;
922: PetscFunctionBegin;
923: PetscCheck(!snes->xl && !snes->xu && !snes->ops->computevariablebounds, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name);
925: PetscCall(PetscCitationsRegister(SNESCitation, &SNEScite));
926: snes->reason = SNES_CONVERGED_ITERATING;
927: X = snes->vec_sol;
928: F = snes->vec_func;
930: PetscCall(SNESFASCycleIsFine(snes, &isFine));
931: /* norm setup */
932: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
933: snes->iter = 0;
934: snes->norm = 0;
935: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
936: if (!snes->vec_func_init_set) {
937: if (fas->eventresidual) PetscCall(PetscLogEventBegin(fas->eventresidual, snes, 0, 0, 0));
938: PetscCall(SNESComputeFunction(snes, X, F));
939: if (fas->eventresidual) PetscCall(PetscLogEventEnd(fas->eventresidual, snes, 0, 0, 0));
940: } else snes->vec_func_init_set = PETSC_FALSE;
942: PetscCall(VecNorm(F, NORM_2, &fnorm)); /* fnorm <- ||F|| */
943: SNESCheckFunctionDomainError(snes, fnorm);
944: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
945: snes->norm = fnorm;
946: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
947: PetscCall(SNESLogConvergenceHistory(snes, fnorm, 0));
949: /* test convergence */
950: PetscCall(SNESConverged(snes, snes->iter, 0.0, 0.0, fnorm));
951: PetscCall(SNESMonitor(snes, snes->iter, fnorm));
952: if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
954: if (isFine) {
955: /* propagate scale-dependent data up the hierarchy */
956: PetscCall(SNESGetDM(snes, &dm));
957: for (ffas = fas; ffas->next; ffas = (SNES_FAS *)ffas->next->data) {
958: DM dmcoarse;
959: PetscCall(SNESGetDM(ffas->next, &dmcoarse));
960: PetscCall(DMRestrict(dm, ffas->restrct, ffas->rscale, ffas->inject, dmcoarse));
961: dm = dmcoarse;
962: }
963: }
965: for (PetscInt i = 0; i < snes->max_its; i++) {
966: /* Call general purpose update function */
967: PetscTryTypeMethod(snes, update, snes->iter);
969: if (fas->fastype == SNES_FAS_MULTIPLICATIVE) PetscCall(SNESFASCycle_Multiplicative(snes, X));
970: else if (fas->fastype == SNES_FAS_ADDITIVE) PetscCall(SNESFASCycle_Additive(snes, X));
971: else if (fas->fastype == SNES_FAS_FULL) PetscCall(SNESFASCycle_Full(snes, X));
972: else {
973: PetscCheck(fas->fastype == SNES_FAS_KASKADE, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Unsupported FAS type");
974: PetscCall(SNESFASCycle_Kaskade(snes, X));
975: }
977: /* check for FAS cycle divergence */
978: if (snes->reason != SNES_CONVERGED_ITERATING) PetscFunctionReturn(PETSC_SUCCESS);
980: /* Monitor convergence */
981: PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
982: snes->iter = i + 1;
983: PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
984: PetscCall(SNESLogConvergenceHistory(snes, snes->norm, 0));
985: PetscCall(SNESConverged(snes, snes->iter, 0.0, 0.0, snes->norm));
986: PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
987: if (snes->reason) break;
988: }
989: PetscFunctionReturn(PETSC_SUCCESS);
990: }
992: static PetscErrorCode SNESFASSetUseCoarseCorrectionLineSearch_FAS(SNES snes, PetscBool use)
993: {
994: SNES_FAS *fas = (SNES_FAS *)snes->data;
996: PetscFunctionBegin;
997: PetscCheck(snes->setfromoptionscalled == 0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Must be called before SNESSetFromOptions()");
998: fas->useCoarseCorrectionLineSearch = use;
999: PetscFunctionReturn(PETSC_SUCCESS);
1000: }
1002: static PetscErrorCode SNESFASGetCoarseCorrectionLineSearch_FAS(SNES snes, PetscInt level, SNESLineSearch *ls)
1003: {
1004: SNES_FAS *fas;
1005: SNES levelsnes;
1007: PetscFunctionBegin;
1008: PetscCheck(snes->setfromoptionscalled > 0, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Must be called after SNESSetFromOptions()");
1009: PetscCall(SNESFASGetCycleSNES(snes, level, &levelsnes));
1010: fas = (SNES_FAS *)levelsnes->data;
1011: *ls = fas->coarseCorrectionLineSearch;
1012: PetscFunctionReturn(PETSC_SUCCESS);
1013: }
1015: /*MC
1016: SNESFAS - An implementation of the Full Approximation Scheme nonlinear multigrid solver, FAS, or nonlinear multigrid {cite}`bruneknepleysmithtu15` for
1017: solving nonlinear systems of equations with `SNES`.
1019: The nonlinear problem is solved by correction using coarse versions
1020: of the nonlinear problem. This problem is perturbed so that a projected
1021: solution of the fine problem elicits no correction from the coarse problem.
1023: Options Database Keys and Prefixes:
1024: + -snes_fas_levels l - The number of levels
1025: . -snes_fas_cycles (1|2) - The number of cycles -- 1 for V, 2 for W
1026: . -snes_fas_type (additive|multiplicative|full|kaskade) - Additive or multiplicative cycle
1027: . -snes_fas_galerkin (true|false) - Form coarse problems by projection back upon the fine problem
1028: . -snes_fas_smoothup u - The number of iterations of the post-smoother
1029: . -snes_fas_smoothdown d - The number of iterations of the pre-smoother
1030: . -snes_fas_monitor - Monitor progress of all of the levels
1031: . -snes_fas_full_downsweep (true|false) - Call the downsmooth on the initial downsweep of full FAS
1032: . -snes_fas_use_coarse_correction_linesearch (true|false) - Use a line search for the multiplicative coarse correction update
1033: on every level except the coarsest
1034: . -fas_coarse_correction_snes_linesearch_type (none|bt|secant|cp|nleqerr|bisection|shell) - Type of the enabled coarse correction line search, on all levels
1035: . -fas_coarse_correction_i_snes_linesearch_type (none|bt|secant|cp|nleqerr|bisection|shell) - Override the coarse correction line search type on level i only
1036: . -fas_levels_snes_ - Prefix for `SNES` options for all smoothers
1037: . -fas_levels_cycle_snes_ - Prefix for `SNES` options for all cycles
1038: . -fas_levels_i_snes_ - Prefix `SNES` options for the smoothers on level i
1039: . -fas_levels_i_cycle_snes_ - Prefix for `SNES` options for the cycle on level i
1040: - -fas_coarse_snes_ - Prefix for `SNES` options for the coarsest smoother
1042: Level: beginner
1044: Notes:
1045: The organization of the `SNESFAS` solver is slightly different from the organization of `PCMG`
1046: As each level has smoother `SNES` instances(down and potentially up) and a cycle `SNES` instance.
1047: The cycle `SNES` instance may be used for monitoring convergence on a particular level.
1049: The coarse correction $X += I(x^c - Rx)$ is applied directly (a unit step, equivalent to the original FAS
1050: correction) unless `-snes_fas_use_coarse_correction_linesearch` or `SNESFASSetUseCoarseCorrectionLineSearch()` is set, in which case a dedicated
1051: `SNESLineSearch` is created on every level, except the coarsest, (one per recursive `SNESFAS` instance) and applied instead.
1052: MG-Opt {cite}`nash2000mgopt` generalizes the correction this way, which is useful when monotonic decrease
1053: in the residual norm or energy functional is desired (e.g., non-convex problems where the full coarse
1054: correction may overshoot without an explicit convergence control strategy). As with `PCMG`'s
1055: `-mg_levels_pc_type` vs. `-mg_levels_2_pc_type` pattern, `-fas_coarse_correction_snes_linesearch_type`
1056: sets the line search type for all levels while `-fas_coarse_correction_i_snes_linesearch_type` overrides
1057: only level i; the line search can be retrieved via `SNESFASGetCoarseCorrectionLineSearch()`.
1059: .seealso: [](ch_snes), `PCMG`, `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESType`, `SNESFASSetRestriction()`, `SNESFASSetInjection()`,
1060: `SNESFASFullGetTotal()`, `SNESFASSetType()`, `SNESFASGetType()`, `SNESFASSetLevels()`, `SNESFASGetLevels()`, `SNESFASGetCycleSNES()`,
1061: `SNESFASSetNumberSmoothUp()`, `SNESFASSetNumberSmoothDown()`, `SNESFASSetContinuation()`, `SNESFASSetCycles()`, `SNESFASSetMonitor()`,
1062: `SNESFASSetLog()`, `SNESFASCycleSetCycles()`, `SNESFASCycleGetSmoother()`, `SNESFASCycleGetSmootherUp()`, `SNESFASCycleGetSmootherDown()`,
1063: `SNESFASCycleGetCorrection()`, `SNESFASCycleGetInterpolation()`, `SNESFASCycleGetRestriction()`, `SNESFASCycleGetInjection()`,
1064: `SNESFASCycleGetRScale()`, `SNESFASCycleIsFine()`, `SNESFASSetInterpolation()`, `SNESFASGetInterpolation()`,
1065: `SNESFASGetRestriction()`, `SNESFASGetInjection()`, `SNESFASSetRScale()`, `SNESFASGetSmoother()`,
1066: `SNESFASGetSmootherDown()`, `SNESFASGetSmootherUp()`, `SNESFASGetCoarseSolve()`, `SNESFASFullSetDownSweep()`, `SNESFASFullSetTotal()`,
1067: `SNESFASSetUseCoarseCorrectionLineSearch()`, `SNESFASGetCoarseCorrectionLineSearch()`
1068: M*/
1070: PETSC_EXTERN PetscErrorCode SNESCreate_FAS(SNES snes)
1071: {
1072: SNES_FAS *fas;
1074: PetscFunctionBegin;
1075: snes->ops->destroy = SNESDestroy_FAS;
1076: snes->ops->setup = SNESSetUp_FAS;
1077: snes->ops->setfromoptions = SNESSetFromOptions_FAS;
1078: snes->ops->view = SNESView_FAS;
1079: snes->ops->solve = SNESSolve_FAS;
1080: snes->ops->reset = SNESReset_FAS;
1082: snes->usesksp = PETSC_FALSE;
1083: snes->usesnpc = PETSC_FALSE;
1085: PetscObjectParameterSetDefault(snes, max_funcs, 30000);
1086: PetscObjectParameterSetDefault(snes, max_its, 10000);
1088: snes->alwayscomputesfinalresidual = PETSC_TRUE;
1090: PetscCall(PetscNew(&fas));
1092: snes->data = (void *)fas;
1093: fas->level = 0;
1094: fas->levels = 1;
1095: fas->n_cycles = 1;
1096: fas->max_up_it = 1;
1097: fas->max_down_it = 1;
1098: fas->smoothu = NULL;
1099: fas->smoothd = NULL;
1100: fas->next = NULL;
1101: fas->previous = NULL;
1102: fas->fine = snes;
1103: fas->interpolate = NULL;
1104: fas->restrct = NULL;
1105: fas->inject = NULL;
1106: fas->useCoarseCorrectionLineSearch = PETSC_FALSE;
1107: fas->coarseCorrectionLineSearch = NULL;
1108: fas->usedmfornumberoflevels = PETSC_FALSE;
1109: fas->fastype = SNES_FAS_MULTIPLICATIVE;
1110: fas->full_downsweep = PETSC_FALSE;
1111: fas->full_total = PETSC_FALSE;
1113: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESFASSetUseCoarseCorrectionLineSearch_C", SNESFASSetUseCoarseCorrectionLineSearch_FAS));
1114: PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESFASGetCoarseCorrectionLineSearch_C", SNESFASGetCoarseCorrectionLineSearch_FAS));
1116: fas->eventsmoothsetup = 0;
1117: fas->eventsmoothsolve = 0;
1118: fas->eventresidual = 0;
1119: fas->eventinterprestrict = 0;
1120: PetscFunctionReturn(PETSC_SUCCESS);
1121: }