Actual source code: tsadapt.c
1: #include <petsc/private/tsimpl.h>
3: PetscClassId TSADAPT_CLASSID;
5: static PetscFunctionList TSAdaptList;
6: static PetscBool TSAdaptPackageInitialized;
7: static PetscBool TSAdaptRegisterAllCalled;
9: PETSC_EXTERN PetscErrorCode TSAdaptCreate_None(TSAdapt);
10: PETSC_EXTERN PetscErrorCode TSAdaptCreate_Basic(TSAdapt);
11: PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt);
12: PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt);
13: PETSC_EXTERN PetscErrorCode TSAdaptCreate_GLEE(TSAdapt);
14: PETSC_EXTERN PetscErrorCode TSAdaptCreate_History(TSAdapt);
16: /*@C
17: TSAdaptRegister - adds a TSAdapt implementation
19: Not Collective, No Fortran Support
21: Input Parameters:
22: + sname - name of user-defined adaptivity scheme
23: - function - routine to create method context
25: Level: advanced
27: Notes:
28: `TSAdaptRegister()` may be called multiple times to add several user-defined families.
30: Example Usage:
31: .vb
32: TSAdaptRegister("my_scheme", MySchemeCreate);
33: .ve
35: Then, your scheme can be chosen with the procedural interface via
36: .vb
37: TSAdaptSetType(ts, "my_scheme")
38: .ve
39: or at runtime via the option
40: .vb
41: -ts_adapt_type my_scheme
42: .ve
44: .seealso: [](ch_ts), `TSAdaptRegisterAll()`
45: @*/
46: PetscErrorCode TSAdaptRegister(const char sname[], PetscErrorCode (*function)(TSAdapt))
47: {
48: PetscFunctionBegin;
49: PetscCall(TSAdaptInitializePackage());
50: PetscCall(PetscFunctionListAdd(&TSAdaptList, sname, function));
51: PetscFunctionReturn(PETSC_SUCCESS);
52: }
54: /*@C
55: TSAdaptRegisterAll - Registers all of the adaptivity schemes in `TSAdapt`
57: Not Collective
59: Level: advanced
61: .seealso: [](ch_ts), `TSAdaptRegisterDestroy()`
62: @*/
63: PetscErrorCode TSAdaptRegisterAll(void)
64: {
65: PetscFunctionBegin;
66: if (TSAdaptRegisterAllCalled) PetscFunctionReturn(PETSC_SUCCESS);
67: TSAdaptRegisterAllCalled = PETSC_TRUE;
68: PetscCall(TSAdaptRegister(TSADAPTNONE, TSAdaptCreate_None));
69: PetscCall(TSAdaptRegister(TSADAPTBASIC, TSAdaptCreate_Basic));
70: PetscCall(TSAdaptRegister(TSADAPTDSP, TSAdaptCreate_DSP));
71: PetscCall(TSAdaptRegister(TSADAPTCFL, TSAdaptCreate_CFL));
72: PetscCall(TSAdaptRegister(TSADAPTGLEE, TSAdaptCreate_GLEE));
73: PetscCall(TSAdaptRegister(TSADAPTHISTORY, TSAdaptCreate_History));
74: PetscFunctionReturn(PETSC_SUCCESS);
75: }
77: /*@C
78: TSAdaptFinalizePackage - This function destroys everything in the `TS` package. It is
79: called from `PetscFinalize()`.
81: Level: developer
83: .seealso: [](ch_ts), `PetscFinalize()`
84: @*/
85: PetscErrorCode TSAdaptFinalizePackage(void)
86: {
87: PetscFunctionBegin;
88: PetscCall(PetscFunctionListDestroy(&TSAdaptList));
89: TSAdaptPackageInitialized = PETSC_FALSE;
90: TSAdaptRegisterAllCalled = PETSC_FALSE;
91: PetscFunctionReturn(PETSC_SUCCESS);
92: }
94: /*@C
95: TSAdaptInitializePackage - This function initializes everything in the `TSAdapt` package. It is
96: called from `TSInitializePackage()`.
98: Level: developer
100: .seealso: [](ch_ts), `PetscInitialize()`
101: @*/
102: PetscErrorCode TSAdaptInitializePackage(void)
103: {
104: PetscFunctionBegin;
105: if (TSAdaptPackageInitialized) PetscFunctionReturn(PETSC_SUCCESS);
106: TSAdaptPackageInitialized = PETSC_TRUE;
107: PetscCall(PetscClassIdRegister("TSAdapt", &TSADAPT_CLASSID));
108: PetscCall(TSAdaptRegisterAll());
109: PetscCall(PetscRegisterFinalize(TSAdaptFinalizePackage));
110: PetscFunctionReturn(PETSC_SUCCESS);
111: }
113: /*@
114: TSAdaptSetType - sets the approach used for the error adapter
116: Logicially Collective
118: Input Parameters:
119: + adapt - the `TS` adapter, most likely obtained with `TSGetAdapt()`
120: - type - one of the `TSAdaptType`
122: Options Database Key:
123: . -ts_adapt_type <basic or dsp or none> - to set the adapter type
125: Level: intermediate
127: .seealso: [](ch_ts), `TSGetAdapt()`, `TSAdaptDestroy()`, `TSAdaptType`, `TSAdaptGetType()`
128: @*/
129: PetscErrorCode TSAdaptSetType(TSAdapt adapt, TSAdaptType type)
130: {
131: PetscBool match;
132: PetscErrorCode (*r)(TSAdapt);
134: PetscFunctionBegin;
136: PetscAssertPointer(type, 2);
137: PetscCall(PetscObjectTypeCompare((PetscObject)adapt, type, &match));
138: if (match) PetscFunctionReturn(PETSC_SUCCESS);
139: PetscCall(PetscFunctionListFind(TSAdaptList, type, &r));
140: PetscCheck(r, PetscObjectComm((PetscObject)adapt), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown TSAdapt type \"%s\" given", type);
141: PetscTryTypeMethod(adapt, destroy);
142: PetscCall(PetscMemzero(adapt->ops, sizeof(struct _TSAdaptOps)));
143: PetscCall(PetscObjectChangeTypeName((PetscObject)adapt, type));
144: PetscCall((*r)(adapt));
145: PetscFunctionReturn(PETSC_SUCCESS);
146: }
148: /*@
149: TSAdaptGetType - gets the `TS` adapter method type (as a string).
151: Not Collective
153: Input Parameter:
154: . adapt - The `TS` adapter, most likely obtained with `TSGetAdapt()`
156: Output Parameter:
157: . type - The name of `TS` adapter method
159: Level: intermediate
161: .seealso: `TSAdapt`, `TSAdaptType`, `TSAdaptSetType()`
162: @*/
163: PetscErrorCode TSAdaptGetType(TSAdapt adapt, TSAdaptType *type)
164: {
165: PetscFunctionBegin;
167: PetscAssertPointer(type, 2);
168: *type = ((PetscObject)adapt)->type_name;
169: PetscFunctionReturn(PETSC_SUCCESS);
170: }
172: PetscErrorCode TSAdaptSetOptionsPrefix(TSAdapt adapt, const char prefix[])
173: {
174: PetscFunctionBegin;
176: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)adapt, prefix));
177: PetscFunctionReturn(PETSC_SUCCESS);
178: }
180: /*@
181: TSAdaptLoad - Loads a TSAdapt that has been stored in binary with `TSAdaptView()`.
183: Collective
185: Input Parameters:
186: + adapt - the newly loaded `TSAdapt`, this needs to have been created with `TSAdaptCreate()` or
187: some related function before a call to `TSAdaptLoad()`.
188: - viewer - binary file viewer, obtained from `PetscViewerBinaryOpen()` or
189: HDF5 file viewer, obtained from `PetscViewerHDF5Open()`
191: Level: intermediate
193: Note:
194: The type is determined by the data in the file, any type set into the `TSAdapt` before this call is ignored.
196: .seealso: [](ch_ts), `PetscViewerBinaryOpen()`, `TSAdaptView()`, `MatLoad()`, `VecLoad()`, `TSAdapt`
197: @*/
198: PetscErrorCode TSAdaptLoad(TSAdapt adapt, PetscViewer viewer)
199: {
200: PetscBool isbinary;
201: char type[256];
203: PetscFunctionBegin;
206: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
207: PetscCheck(isbinary, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid viewer; open viewer with PetscViewerBinaryOpen()");
209: PetscCall(PetscViewerBinaryRead(viewer, type, 256, NULL, PETSC_CHAR));
210: PetscCall(TSAdaptSetType(adapt, type));
211: PetscTryTypeMethod(adapt, load, viewer);
212: PetscFunctionReturn(PETSC_SUCCESS);
213: }
215: PetscErrorCode TSAdaptView(TSAdapt adapt, PetscViewer viewer)
216: {
217: PetscBool iascii, isbinary, isnone, isglee;
219: PetscFunctionBegin;
221: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)adapt), &viewer));
223: PetscCheckSameComm(adapt, 1, viewer, 2);
224: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
225: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
226: if (iascii) {
227: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)adapt, viewer));
228: PetscCall(PetscObjectTypeCompare((PetscObject)adapt, TSADAPTNONE, &isnone));
229: PetscCall(PetscObjectTypeCompare((PetscObject)adapt, TSADAPTGLEE, &isglee));
230: if (!isnone) {
231: if (adapt->always_accept) PetscCall(PetscViewerASCIIPrintf(viewer, " always accepting steps\n"));
232: PetscCall(PetscViewerASCIIPrintf(viewer, " safety factor %g\n", (double)adapt->safety));
233: PetscCall(PetscViewerASCIIPrintf(viewer, " extra safety factor after step rejection %g\n", (double)adapt->reject_safety));
234: PetscCall(PetscViewerASCIIPrintf(viewer, " clip fastest increase %g\n", (double)adapt->clip[1]));
235: PetscCall(PetscViewerASCIIPrintf(viewer, " clip fastest decrease %g\n", (double)adapt->clip[0]));
236: PetscCall(PetscViewerASCIIPrintf(viewer, " maximum allowed timestep %g\n", (double)adapt->dt_max));
237: PetscCall(PetscViewerASCIIPrintf(viewer, " minimum allowed timestep %g\n", (double)adapt->dt_min));
238: PetscCall(PetscViewerASCIIPrintf(viewer, " maximum solution absolute value to be ignored %g\n", (double)adapt->ignore_max));
239: }
240: if (isglee) {
241: if (adapt->glee_use_local) {
242: PetscCall(PetscViewerASCIIPrintf(viewer, " GLEE uses local error control\n"));
243: } else {
244: PetscCall(PetscViewerASCIIPrintf(viewer, " GLEE uses global error control\n"));
245: }
246: }
247: PetscCall(PetscViewerASCIIPushTab(viewer));
248: PetscTryTypeMethod(adapt, view, viewer);
249: PetscCall(PetscViewerASCIIPopTab(viewer));
250: } else if (isbinary) {
251: char type[256];
253: /* need to save FILE_CLASS_ID for adapt class */
254: PetscCall(PetscStrncpy(type, ((PetscObject)adapt)->type_name, 256));
255: PetscCall(PetscViewerBinaryWrite(viewer, type, 256, PETSC_CHAR));
256: } else PetscTryTypeMethod(adapt, view, viewer);
257: PetscFunctionReturn(PETSC_SUCCESS);
258: }
260: /*@
261: TSAdaptReset - Resets a `TSAdapt` context to its defaults
263: Collective
265: Input Parameter:
266: . adapt - the `TSAdapt` context obtained from `TSGetAdapt()` or `TSAdaptCreate()`
268: Level: developer
270: .seealso: [](ch_ts), `TSGetAdapt()`, `TSAdapt`, `TSAdaptCreate()`, `TSAdaptDestroy()`
271: @*/
272: PetscErrorCode TSAdaptReset(TSAdapt adapt)
273: {
274: PetscFunctionBegin;
276: PetscTryTypeMethod(adapt, reset);
277: PetscFunctionReturn(PETSC_SUCCESS);
278: }
280: PetscErrorCode TSAdaptDestroy(TSAdapt *adapt)
281: {
282: PetscFunctionBegin;
283: if (!*adapt) PetscFunctionReturn(PETSC_SUCCESS);
285: if (--((PetscObject)*adapt)->refct > 0) {
286: *adapt = NULL;
287: PetscFunctionReturn(PETSC_SUCCESS);
288: }
290: PetscCall(TSAdaptReset(*adapt));
292: PetscTryTypeMethod(*adapt, destroy);
293: PetscCall(PetscViewerDestroy(&(*adapt)->monitor));
294: PetscCall(PetscHeaderDestroy(adapt));
295: PetscFunctionReturn(PETSC_SUCCESS);
296: }
298: /*@
299: TSAdaptSetMonitor - Monitor the choices made by the adaptive controller
301: Collective
303: Input Parameters:
304: + adapt - adaptive controller context
305: - flg - `PETSC_TRUE` to active a monitor, `PETSC_FALSE` to disable
307: Options Database Key:
308: . -ts_adapt_monitor - to turn on monitoring
310: Level: intermediate
312: .seealso: [](ch_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptChoose()`
313: @*/
314: PetscErrorCode TSAdaptSetMonitor(TSAdapt adapt, PetscBool flg)
315: {
316: PetscFunctionBegin;
319: if (flg) {
320: if (!adapt->monitor) PetscCall(PetscViewerASCIIOpen(PetscObjectComm((PetscObject)adapt), "stdout", &adapt->monitor));
321: } else {
322: PetscCall(PetscViewerDestroy(&adapt->monitor));
323: }
324: PetscFunctionReturn(PETSC_SUCCESS);
325: }
327: /*@C
328: TSAdaptSetCheckStage - Set a callback to check convergence for a stage
330: Logically Collective
332: Input Parameters:
333: + adapt - adaptive controller context
334: - func - stage check function
336: Calling sequence:
337: + adapt - adaptive controller context
338: . ts - time stepping context
339: . t - current time
340: . Y - current solution vector
341: - accept - pending choice of whether to accept, can be modified by this routine
343: Level: advanced
345: .seealso: [](ch_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptChoose()`
346: @*/
347: PetscErrorCode TSAdaptSetCheckStage(TSAdapt adapt, PetscErrorCode (*func)(TSAdapt adapt, TS ts, PetscReal t, Vec Y, PetscBool *accept))
348: {
349: PetscFunctionBegin;
351: adapt->checkstage = func;
352: PetscFunctionReturn(PETSC_SUCCESS);
353: }
355: /*@
356: TSAdaptSetAlwaysAccept - Set whether to always accept steps regardless of
357: any error or stability condition not meeting the prescribed goal.
359: Logically Collective
361: Input Parameters:
362: + adapt - time step adaptivity context, usually gotten with `TSGetAdapt()`
363: - flag - whether to always accept steps
365: Options Database Key:
366: . -ts_adapt_always_accept - to always accept steps
368: Level: intermediate
370: .seealso: [](ch_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptChoose()`
371: @*/
372: PetscErrorCode TSAdaptSetAlwaysAccept(TSAdapt adapt, PetscBool flag)
373: {
374: PetscFunctionBegin;
377: adapt->always_accept = flag;
378: PetscFunctionReturn(PETSC_SUCCESS);
379: }
381: /*@
382: TSAdaptSetSafety - Set safety factors for time step adaptor
384: Logically Collective
386: Input Parameters:
387: + adapt - adaptive controller context
388: . safety - safety factor relative to target error/stability goal
389: - reject_safety - extra safety factor to apply if the last step was rejected
391: Options Database Keys:
392: + -ts_adapt_safety <safety> - to set safety factor
393: - -ts_adapt_reject_safety <reject_safety> - to set reject safety factor
395: Level: intermediate
397: Note:
398: Use `PETSC_CURRENT` to keep the current value for either parameter
400: Fortran Note:
401: Use `PETSC_CURRENT_REAL`
403: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptGetSafety()`, `TSAdaptChoose()`
404: @*/
405: PetscErrorCode TSAdaptSetSafety(TSAdapt adapt, PetscReal safety, PetscReal reject_safety)
406: {
407: PetscFunctionBegin;
411: PetscCheck(safety == (PetscReal)PETSC_CURRENT || safety >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Safety factor %g must be non negative", (double)safety);
412: PetscCheck(safety == (PetscReal)PETSC_CURRENT || safety <= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Safety factor %g must be less than one", (double)safety);
413: PetscCheck(reject_safety == (PetscReal)PETSC_CURRENT || reject_safety >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Reject safety factor %g must be non negative", (double)reject_safety);
414: PetscCheck(reject_safety == (PetscReal)PETSC_CURRENT || reject_safety <= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Reject safety factor %g must be less than one", (double)reject_safety);
415: if (safety != (PetscReal)PETSC_CURRENT) adapt->safety = safety;
416: if (reject_safety != (PetscReal)PETSC_CURRENT) adapt->reject_safety = reject_safety;
417: PetscFunctionReturn(PETSC_SUCCESS);
418: }
420: /*@
421: TSAdaptGetSafety - Get safety factors for time step adapter
423: Not Collective
425: Input Parameter:
426: . adapt - adaptive controller context
428: Output Parameters:
429: + safety - safety factor relative to target error/stability goal
430: - reject_safety - extra safety factor to apply if the last step was rejected
432: Level: intermediate
434: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptSetSafety()`, `TSAdaptChoose()`
435: @*/
436: PetscErrorCode TSAdaptGetSafety(TSAdapt adapt, PetscReal *safety, PetscReal *reject_safety)
437: {
438: PetscFunctionBegin;
440: if (safety) PetscAssertPointer(safety, 2);
441: if (reject_safety) PetscAssertPointer(reject_safety, 3);
442: if (safety) *safety = adapt->safety;
443: if (reject_safety) *reject_safety = adapt->reject_safety;
444: PetscFunctionReturn(PETSC_SUCCESS);
445: }
447: /*@
448: TSAdaptSetMaxIgnore - Set error estimation threshold. Solution components below this threshold value will not be considered when computing error norms
449: for time step adaptivity (in absolute value). A negative value (default) of the threshold leads to considering all solution components.
451: Logically Collective
453: Input Parameters:
454: + adapt - adaptive controller context
455: - max_ignore - threshold for solution components that are ignored during error estimation
457: Options Database Key:
458: . -ts_adapt_max_ignore <max_ignore> - to set the threshold
460: Level: intermediate
462: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptGetMaxIgnore()`, `TSAdaptChoose()`
463: @*/
464: PetscErrorCode TSAdaptSetMaxIgnore(TSAdapt adapt, PetscReal max_ignore)
465: {
466: PetscFunctionBegin;
469: adapt->ignore_max = max_ignore;
470: PetscFunctionReturn(PETSC_SUCCESS);
471: }
473: /*@
474: TSAdaptGetMaxIgnore - Get error estimation threshold. Solution components below this threshold value will not be considered when computing error norms
475: for time step adaptivity (in absolute value).
477: Not Collective
479: Input Parameter:
480: . adapt - adaptive controller context
482: Output Parameter:
483: . max_ignore - threshold for solution components that are ignored during error estimation
485: Level: intermediate
487: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptSetMaxIgnore()`, `TSAdaptChoose()`
488: @*/
489: PetscErrorCode TSAdaptGetMaxIgnore(TSAdapt adapt, PetscReal *max_ignore)
490: {
491: PetscFunctionBegin;
493: PetscAssertPointer(max_ignore, 2);
494: *max_ignore = adapt->ignore_max;
495: PetscFunctionReturn(PETSC_SUCCESS);
496: }
498: /*@
499: TSAdaptSetClip - Sets the admissible decrease/increase factor in step size in the time step adapter
501: Logically collective
503: Input Parameters:
504: + adapt - adaptive controller context
505: . low - admissible decrease factor
506: - high - admissible increase factor
508: Options Database Key:
509: . -ts_adapt_clip <low>,<high> - to set admissible time step decrease and increase factors
511: Level: intermediate
513: Note:
514: Use `PETSC_CURRENT` to keep the current value for either parameter
516: Fortran Note:
517: Use `PETSC_CURRENT_REAL`
519: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptGetClip()`, `TSAdaptSetScaleSolveFailed()`
520: @*/
521: PetscErrorCode TSAdaptSetClip(TSAdapt adapt, PetscReal low, PetscReal high)
522: {
523: PetscFunctionBegin;
527: PetscCheck(low == (PetscReal)PETSC_CURRENT || low >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Decrease factor %g must be non negative", (double)low);
528: PetscCheck(low == (PetscReal)PETSC_CURRENT || low <= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Decrease factor %g must be less than one", (double)low);
529: PetscCheck(high == (PetscReal)PETSC_CURRENT || high >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Increase factor %g must be greater than one", (double)high);
530: if (low != (PetscReal)PETSC_CURRENT) adapt->clip[0] = low;
531: if (high != (PetscReal)PETSC_CURRENT) adapt->clip[1] = high;
532: PetscFunctionReturn(PETSC_SUCCESS);
533: }
535: /*@
536: TSAdaptGetClip - Gets the admissible decrease/increase factor in step size in the time step adapter
538: Not Collective
540: Input Parameter:
541: . adapt - adaptive controller context
543: Output Parameters:
544: + low - optional, admissible decrease factor
545: - high - optional, admissible increase factor
547: Level: intermediate
549: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptSetClip()`, `TSAdaptSetScaleSolveFailed()`
550: @*/
551: PetscErrorCode TSAdaptGetClip(TSAdapt adapt, PetscReal *low, PetscReal *high)
552: {
553: PetscFunctionBegin;
555: if (low) PetscAssertPointer(low, 2);
556: if (high) PetscAssertPointer(high, 3);
557: if (low) *low = adapt->clip[0];
558: if (high) *high = adapt->clip[1];
559: PetscFunctionReturn(PETSC_SUCCESS);
560: }
562: /*@
563: TSAdaptSetScaleSolveFailed - Scale step size by this factor if solve fails
565: Logically Collective
567: Input Parameters:
568: + adapt - adaptive controller context
569: - scale - scale
571: Options Database Key:
572: . -ts_adapt_scale_solve_failed <scale> - to set scale step by this factor if solve fails
574: Level: intermediate
576: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptGetScaleSolveFailed()`, `TSAdaptGetClip()`
577: @*/
578: PetscErrorCode TSAdaptSetScaleSolveFailed(TSAdapt adapt, PetscReal scale)
579: {
580: PetscFunctionBegin;
583: PetscCheck(scale > 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Scale factor %g must be positive", (double)scale);
584: PetscCheck(scale <= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Scale factor %g must be less than one", (double)scale);
585: adapt->scale_solve_failed = scale;
586: PetscFunctionReturn(PETSC_SUCCESS);
587: }
589: /*@
590: TSAdaptGetScaleSolveFailed - Gets the admissible decrease/increase factor in step size
592: Not Collective
594: Input Parameter:
595: . adapt - adaptive controller context
597: Output Parameter:
598: . scale - scale factor
600: Level: intermediate
602: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptSetScaleSolveFailed()`, `TSAdaptSetClip()`
603: @*/
604: PetscErrorCode TSAdaptGetScaleSolveFailed(TSAdapt adapt, PetscReal *scale)
605: {
606: PetscFunctionBegin;
608: if (scale) PetscAssertPointer(scale, 2);
609: if (scale) *scale = adapt->scale_solve_failed;
610: PetscFunctionReturn(PETSC_SUCCESS);
611: }
613: /*@
614: TSAdaptSetStepLimits - Set the minimum and maximum step sizes to be considered by the time step controller
616: Logically Collective
618: Input Parameters:
619: + adapt - time step adaptivity context, usually gotten with `TSGetAdapt()`
620: . hmin - minimum time step
621: - hmax - maximum time step
623: Options Database Keys:
624: + -ts_adapt_dt_min <min> - to set minimum time step
625: - -ts_adapt_dt_max <max> - to set maximum time step
627: Level: intermediate
629: Note:
630: Use `PETSC_CURRENT` to keep the current value for either parameter
632: Fortran Note:
633: Use `PETSC_CURRENT_REAL`
635: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptGetStepLimits()`, `TSAdaptChoose()`
636: @*/
637: PetscErrorCode TSAdaptSetStepLimits(TSAdapt adapt, PetscReal hmin, PetscReal hmax)
638: {
639: PetscFunctionBegin;
643: PetscCheck(hmin == (PetscReal)PETSC_CURRENT || hmin >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Minimum time step %g must be non negative", (double)hmin);
644: PetscCheck(hmax == (PetscReal)PETSC_CURRENT || hmax >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Minimum time step %g must be non negative", (double)hmax);
645: if (hmin != (PetscReal)PETSC_CURRENT) adapt->dt_min = hmin;
646: if (hmax != (PetscReal)PETSC_CURRENT) adapt->dt_max = hmax;
647: hmin = adapt->dt_min;
648: hmax = adapt->dt_max;
649: PetscCheck(hmax > hmin, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Maximum time step %g must greater than minimum time step %g", (double)hmax, (double)hmin);
650: PetscFunctionReturn(PETSC_SUCCESS);
651: }
653: /*@
654: TSAdaptGetStepLimits - Get the minimum and maximum step sizes to be considered by the time step controller
656: Not Collective
658: Input Parameter:
659: . adapt - time step adaptivity context, usually gotten with `TSGetAdapt()`
661: Output Parameters:
662: + hmin - minimum time step
663: - hmax - maximum time step
665: Level: intermediate
667: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptSetStepLimits()`, `TSAdaptChoose()`
668: @*/
669: PetscErrorCode TSAdaptGetStepLimits(TSAdapt adapt, PetscReal *hmin, PetscReal *hmax)
670: {
671: PetscFunctionBegin;
673: if (hmin) PetscAssertPointer(hmin, 2);
674: if (hmax) PetscAssertPointer(hmax, 3);
675: if (hmin) *hmin = adapt->dt_min;
676: if (hmax) *hmax = adapt->dt_max;
677: PetscFunctionReturn(PETSC_SUCCESS);
678: }
680: /*@C
681: TSAdaptSetFromOptions - Sets various `TSAdapt` parameters from user options.
683: Collective
685: Input Parameters:
686: + adapt - the `TSAdapt` context
687: - PetscOptionsObject - object created by `PetscOptionsBegin()`
689: Options Database Keys:
690: + -ts_adapt_type <type> - algorithm to use for adaptivity
691: . -ts_adapt_always_accept - always accept steps regardless of error/stability goals
692: . -ts_adapt_safety <safety> - safety factor relative to target error/stability goal
693: . -ts_adapt_reject_safety <safety> - extra safety factor to apply if the last step was rejected
694: . -ts_adapt_clip <low,high> - admissible time step decrease and increase factors
695: . -ts_adapt_dt_min <min> - minimum timestep to use
696: . -ts_adapt_dt_max <max> - maximum timestep to use
697: . -ts_adapt_scale_solve_failed <scale> - scale timestep by this factor if a solve fails
698: . -ts_adapt_wnormtype <2 or infinity> - type of norm for computing error estimates
699: - -ts_adapt_time_step_increase_delay - number of timesteps to delay increasing the time step after it has been decreased due to failed solver
701: Level: advanced
703: Note:
704: This function is automatically called by `TSSetFromOptions()`
706: .seealso: [](ch_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptSetAlwaysAccept()`, `TSAdaptSetSafety()`,
707: `TSAdaptSetClip()`, `TSAdaptSetScaleSolveFailed()`, `TSAdaptSetStepLimits()`, `TSAdaptSetMonitor()`
708: @*/
709: PetscErrorCode TSAdaptSetFromOptions(TSAdapt adapt, PetscOptionItems PetscOptionsObject)
710: {
711: char type[256] = TSADAPTBASIC;
712: PetscReal safety, reject_safety, clip[2], scale, hmin, hmax;
713: PetscBool set, flg;
714: PetscInt two;
716: PetscFunctionBegin;
718: /* This should use PetscOptionsBegin() if/when this becomes an object used outside of TS, but currently this
719: * function can only be called from inside TSSetFromOptions() */
720: PetscOptionsHeadBegin(PetscOptionsObject, "TS Adaptivity options");
721: PetscCall(PetscOptionsFList("-ts_adapt_type", "Algorithm to use for adaptivity", "TSAdaptSetType", TSAdaptList, ((PetscObject)adapt)->type_name ? ((PetscObject)adapt)->type_name : type, type, sizeof(type), &flg));
722: if (flg || !((PetscObject)adapt)->type_name) PetscCall(TSAdaptSetType(adapt, type));
724: PetscCall(PetscOptionsBool("-ts_adapt_always_accept", "Always accept the step", "TSAdaptSetAlwaysAccept", adapt->always_accept, &flg, &set));
725: if (set) PetscCall(TSAdaptSetAlwaysAccept(adapt, flg));
727: safety = adapt->safety;
728: reject_safety = adapt->reject_safety;
729: PetscCall(PetscOptionsReal("-ts_adapt_safety", "Safety factor relative to target error/stability goal", "TSAdaptSetSafety", safety, &safety, &set));
730: PetscCall(PetscOptionsReal("-ts_adapt_reject_safety", "Extra safety factor to apply if the last step was rejected", "TSAdaptSetSafety", reject_safety, &reject_safety, &flg));
731: if (set || flg) PetscCall(TSAdaptSetSafety(adapt, safety, reject_safety));
733: two = 2;
734: clip[0] = adapt->clip[0];
735: clip[1] = adapt->clip[1];
736: PetscCall(PetscOptionsRealArray("-ts_adapt_clip", "Admissible decrease/increase factor in step size", "TSAdaptSetClip", clip, &two, &set));
737: PetscCheck(!set || (two == 2), PetscObjectComm((PetscObject)adapt), PETSC_ERR_ARG_OUTOFRANGE, "Must give exactly two values to -ts_adapt_clip");
738: if (set) PetscCall(TSAdaptSetClip(adapt, clip[0], clip[1]));
740: hmin = adapt->dt_min;
741: hmax = adapt->dt_max;
742: PetscCall(PetscOptionsReal("-ts_adapt_dt_min", "Minimum time step considered", "TSAdaptSetStepLimits", hmin, &hmin, &set));
743: PetscCall(PetscOptionsReal("-ts_adapt_dt_max", "Maximum time step considered", "TSAdaptSetStepLimits", hmax, &hmax, &flg));
744: if (set || flg) PetscCall(TSAdaptSetStepLimits(adapt, hmin, hmax));
746: PetscCall(PetscOptionsReal("-ts_adapt_max_ignore", "Adaptor ignores (absolute) solution values smaller than this value", "", adapt->ignore_max, &adapt->ignore_max, &set));
747: PetscCall(PetscOptionsBool("-ts_adapt_glee_use_local", "GLEE adaptor uses local error estimation for step control", "", adapt->glee_use_local, &adapt->glee_use_local, &set));
749: PetscCall(PetscOptionsReal("-ts_adapt_scale_solve_failed", "Scale step by this factor if solve fails", "TSAdaptSetScaleSolveFailed", adapt->scale_solve_failed, &scale, &set));
750: if (set) PetscCall(TSAdaptSetScaleSolveFailed(adapt, scale));
752: PetscCall(PetscOptionsEnum("-ts_adapt_wnormtype", "Type of norm computed for error estimation", "", NormTypes, (PetscEnum)adapt->wnormtype, (PetscEnum *)&adapt->wnormtype, NULL));
753: PetscCheck(adapt->wnormtype == NORM_2 || adapt->wnormtype == NORM_INFINITY, PetscObjectComm((PetscObject)adapt), PETSC_ERR_SUP, "Only 2-norm and infinite norm supported");
755: PetscCall(PetscOptionsInt("-ts_adapt_time_step_increase_delay", "Number of timesteps to delay increasing the time step after it has been decreased due to failed solver", "TSAdaptSetTimeStepIncreaseDelay", adapt->timestepjustdecreased_delay, &adapt->timestepjustdecreased_delay, NULL));
757: PetscCall(PetscOptionsBool("-ts_adapt_monitor", "Print choices made by adaptive controller", "TSAdaptSetMonitor", adapt->monitor ? PETSC_TRUE : PETSC_FALSE, &flg, &set));
758: if (set) PetscCall(TSAdaptSetMonitor(adapt, flg));
760: PetscTryTypeMethod(adapt, setfromoptions, PetscOptionsObject);
761: PetscOptionsHeadEnd();
762: PetscFunctionReturn(PETSC_SUCCESS);
763: }
765: /*@
766: TSAdaptCandidatesClear - clear any previously set candidate schemes
768: Logically Collective
770: Input Parameter:
771: . adapt - adaptive controller
773: Level: developer
775: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptCreate()`, `TSAdaptCandidateAdd()`, `TSAdaptChoose()`
776: @*/
777: PetscErrorCode TSAdaptCandidatesClear(TSAdapt adapt)
778: {
779: PetscFunctionBegin;
781: PetscCall(PetscMemzero(&adapt->candidates, sizeof(adapt->candidates)));
782: PetscFunctionReturn(PETSC_SUCCESS);
783: }
785: /*@C
786: TSAdaptCandidateAdd - add a candidate scheme for the adaptive controller to select from
788: Logically Collective; No Fortran Support
790: Input Parameters:
791: + adapt - time step adaptivity context, obtained with `TSGetAdapt()` or `TSAdaptCreate()`
792: . name - name of the candidate scheme to add
793: . order - order of the candidate scheme
794: . stageorder - stage order of the candidate scheme
795: . ccfl - stability coefficient relative to explicit Euler, used for CFL constraints
796: . cost - relative measure of the amount of work required for the candidate scheme
797: - inuse - indicates that this scheme is the one currently in use, this flag can only be set for one scheme
799: Level: developer
801: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptCandidatesClear()`, `TSAdaptChoose()`
802: @*/
803: PetscErrorCode TSAdaptCandidateAdd(TSAdapt adapt, const char name[], PetscInt order, PetscInt stageorder, PetscReal ccfl, PetscReal cost, PetscBool inuse)
804: {
805: PetscInt c;
807: PetscFunctionBegin;
809: PetscCheck(order >= 1, PetscObjectComm((PetscObject)adapt), PETSC_ERR_ARG_OUTOFRANGE, "Classical order %" PetscInt_FMT " must be a positive integer", order);
810: if (inuse) {
811: PetscCheck(!adapt->candidates.inuse_set, PetscObjectComm((PetscObject)adapt), PETSC_ERR_ARG_WRONGSTATE, "Cannot set the inuse method twice, maybe forgot to call TSAdaptCandidatesClear()");
812: adapt->candidates.inuse_set = PETSC_TRUE;
813: }
814: /* first slot if this is the current scheme, otherwise the next available slot */
815: c = inuse ? 0 : !adapt->candidates.inuse_set + adapt->candidates.n;
817: adapt->candidates.name[c] = name;
818: adapt->candidates.order[c] = order;
819: adapt->candidates.stageorder[c] = stageorder;
820: adapt->candidates.ccfl[c] = ccfl;
821: adapt->candidates.cost[c] = cost;
822: adapt->candidates.n++;
823: PetscFunctionReturn(PETSC_SUCCESS);
824: }
826: /*@C
827: TSAdaptCandidatesGet - Get the list of candidate orders of accuracy and cost
829: Not Collective
831: Input Parameter:
832: . adapt - time step adaptivity context
834: Output Parameters:
835: + n - number of candidate schemes, always at least 1
836: . order - the order of each candidate scheme
837: . stageorder - the stage order of each candidate scheme
838: . ccfl - the CFL coefficient of each scheme
839: - cost - the relative cost of each scheme
841: Level: developer
843: Note:
844: The current scheme is always returned in the first slot
846: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptCandidatesClear()`, `TSAdaptCandidateAdd()`, `TSAdaptChoose()`
847: @*/
848: PetscErrorCode TSAdaptCandidatesGet(TSAdapt adapt, PetscInt *n, const PetscInt **order, const PetscInt **stageorder, const PetscReal **ccfl, const PetscReal **cost)
849: {
850: PetscFunctionBegin;
852: if (n) *n = adapt->candidates.n;
853: if (order) *order = adapt->candidates.order;
854: if (stageorder) *stageorder = adapt->candidates.stageorder;
855: if (ccfl) *ccfl = adapt->candidates.ccfl;
856: if (cost) *cost = adapt->candidates.cost;
857: PetscFunctionReturn(PETSC_SUCCESS);
858: }
860: /*@C
861: TSAdaptChoose - choose which method and step size to use for the next step
863: Collective
865: Input Parameters:
866: + adapt - adaptive controller
867: . ts - time stepper
868: - h - current step size
870: Output Parameters:
871: + next_sc - optional, scheme to use for the next step
872: . next_h - step size to use for the next step
873: - accept - `PETSC_TRUE` to accept the current step, `PETSC_FALSE` to repeat the current step with the new step size
875: Level: developer
877: Note:
878: The input value of parameter accept is retained from the last time step, so it will be `PETSC_FALSE` if the step is
879: being retried after an initial rejection.
881: .seealso: [](ch_ts), `TSAdapt`, `TSAdaptCandidatesClear()`, `TSAdaptCandidateAdd()`
882: @*/
883: PetscErrorCode TSAdaptChoose(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept)
884: {
885: PetscInt ncandidates = adapt->candidates.n;
886: PetscInt scheme = 0;
887: PetscReal wlte = -1.0;
888: PetscReal wltea = -1.0;
889: PetscReal wlter = -1.0;
891: PetscFunctionBegin;
894: if (next_sc) PetscAssertPointer(next_sc, 4);
895: PetscAssertPointer(next_h, 5);
896: PetscAssertPointer(accept, 6);
897: if (next_sc) *next_sc = 0;
899: /* Do not mess with adaptivity while handling events */
900: if (ts->event && ts->event->processing) {
901: *next_h = h;
902: *accept = PETSC_TRUE;
903: if (adapt->monitor) {
904: PetscCall(PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
906: if (ts->event->iterctr == 0) {
907: /*
908: An event has been found, now finalising the event processing: performing the 1st and 2nd post-event steps.
909: Entering this if-branch means both these steps (set to either PETSC_DECIDE or numerical value) are managed
910: by the event handler. In this case the 1st post-event step is always accepted, without interference of TSAdapt.
911: Note: if the 2nd post-event step is not managed by the event handler (e.g. given 1st = numerical, 2nd = PETSC_DECIDE),
912: this if-branch is not entered, and TSAdapt may reject/adjust the proposed 1st post-event step.
913: */
914: PetscCall(PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt does not interfere, step %3" PetscInt_FMT " accepted. Processing post-event steps: 1-st accepted just now, 2-nd yet to come\n", ts->steps));
915: } else PetscCall(PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt does not interfere, step %3" PetscInt_FMT " accepted. Event handling in progress\n", ts->steps));
917: PetscCall(PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
918: }
919: PetscFunctionReturn(PETSC_SUCCESS);
920: }
922: PetscUseTypeMethod(adapt, choose, ts, h, &scheme, next_h, accept, &wlte, &wltea, &wlter);
923: PetscCheck(scheme >= 0 && (ncandidates <= 0 || scheme < ncandidates), PetscObjectComm((PetscObject)adapt), PETSC_ERR_ARG_OUTOFRANGE, "Chosen scheme %" PetscInt_FMT " not in valid range 0..%" PetscInt_FMT, scheme, ncandidates - 1);
924: PetscCheck(*next_h >= 0, PetscObjectComm((PetscObject)adapt), PETSC_ERR_ARG_OUTOFRANGE, "Computed step size %g must be positive", (double)*next_h);
925: if (next_sc) *next_sc = scheme;
927: if (*accept && ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
928: /* Increase/reduce step size if end time of next step is close to or overshoots max time */
929: PetscReal t = ts->ptime + ts->time_step, tend, tmax, h1, hmax;
930: PetscReal a = (PetscReal)(1.0 + adapt->matchstepfac[0]);
931: PetscReal b = adapt->matchstepfac[1];
932: PetscReal eps = 10 * PETSC_MACHINE_EPSILON;
934: /*
935: Logic in using 'dt_span_cached':
936: 1. It always overrides *next_h, except (any of):
937: a) the current step was rejected,
938: b) the adaptor proposed to decrease the next step,
939: c) the adaptor proposed *next_h > dt_span_cached.
940: 2. If *next_h was adjusted by eval_times points (or the final point):
941: -- when dt_span_cached is filled (>0), it keeps its value,
942: -- when dt_span_cached is clear (==0), it gets the unadjusted version of *next_h.
943: 3. If *next_h was not adjusted as in (2), dt_span_cached is cleared.
944: Note, if a combination (1.b || 1.c) && (3) takes place, this means that
945: dt_span_cached remains unused at the moment of clearing.
946: If (1.a) takes place, dt_span_cached keeps its value.
947: Also, dt_span_cached can be updated by the event handler, see tsevent.c.
948: */
949: if (h <= *next_h && *next_h <= adapt->dt_eval_times_cached) *next_h = adapt->dt_eval_times_cached; /* try employing the cache */
950: h1 = *next_h;
951: tend = t + h1;
953: if (ts->eval_times && ts->eval_times->time_point_idx < ts->eval_times->num_time_points) {
954: PetscCheck(ts->eval_times->worktol == 0, PetscObjectComm((PetscObject)adapt), PETSC_ERR_PLIB, "Unexpected state (tspan->worktol != 0) in TSAdaptChoose()");
955: ts->eval_times->worktol = ts->eval_times->reltol * h1 + ts->eval_times->abstol;
956: if (PetscIsCloseAtTol(t, ts->eval_times->time_points[ts->eval_times->time_point_idx], ts->eval_times->worktol, 0)) /* hit a span time point */
957: if (ts->eval_times->time_point_idx + 1 < ts->eval_times->num_time_points) tmax = ts->eval_times->time_points[ts->eval_times->time_point_idx + 1];
958: else tmax = ts->max_time; /* hit the last span time point */
959: else tmax = ts->eval_times->time_points[ts->eval_times->time_point_idx];
960: } else tmax = ts->max_time;
961: tmax = PetscMin(tmax, ts->max_time);
962: hmax = tmax - t;
963: PetscCheck((hmax > eps) || (PetscAbsReal(hmax) <= eps && PetscIsCloseAtTol(t, ts->max_time, eps, 0)), PetscObjectComm((PetscObject)adapt), PETSC_ERR_PLIB, "Unexpected state: bad hmax in TSAdaptChoose()");
965: if (t < tmax && tend > tmax) *next_h = hmax;
966: if (t < tmax && tend < tmax && h1 * b > hmax) *next_h = hmax / 2;
967: if (t < tmax && tend < tmax && h1 * a > hmax) *next_h = hmax;
968: if (ts->eval_times && h1 != *next_h && !adapt->dt_eval_times_cached) adapt->dt_eval_times_cached = h1; /* cache the step size if it is to be changed */
969: if (ts->eval_times && h1 == *next_h && adapt->dt_eval_times_cached) adapt->dt_eval_times_cached = 0; /* clear the cache if the step size is unchanged */
970: }
971: if (adapt->monitor) {
972: const char *sc_name = (scheme < ncandidates) ? adapt->candidates.name[scheme] : "";
973: PetscCall(PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
974: if (wlte < 0) {
975: PetscCall(PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt %s %s %" PetscInt_FMT ":%s step %3" PetscInt_FMT " %s t=%-11g+%10.3e dt=%-10.3e\n", ((PetscObject)adapt)->type_name, ((PetscObject)ts)->type_name, scheme, sc_name, ts->steps, *accept ? "accepted" : "rejected",
976: (double)ts->ptime, (double)h, (double)*next_h));
977: } else {
978: PetscCall(PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt %s %s %" PetscInt_FMT ":%s step %3" PetscInt_FMT " %s t=%-11g+%10.3e dt=%-10.3e wlte=%5.3g wltea=%5.3g wlter=%5.3g\n", ((PetscObject)adapt)->type_name, ((PetscObject)ts)->type_name, scheme, sc_name, ts->steps, *accept ? "accepted" : "rejected",
979: (double)ts->ptime, (double)h, (double)*next_h, (double)wlte, (double)wltea, (double)wlter));
980: }
981: PetscCall(PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
982: }
983: PetscFunctionReturn(PETSC_SUCCESS);
984: }
986: /*@
987: TSAdaptSetTimeStepIncreaseDelay - The number of timesteps to wait after a decrease in the timestep due to failed solver
988: before increasing the time step.
990: Logicially Collective
992: Input Parameters:
993: + adapt - adaptive controller context
994: - cnt - the number of timesteps
996: Options Database Key:
997: . -ts_adapt_time_step_increase_delay cnt - number of steps to delay the increase
999: Level: advanced
1001: Notes:
1002: This is to prevent an adaptor from bouncing back and forth between two nearby timesteps. The default is 0.
1004: The successful use of this option is problem dependent
1006: Developer Notes:
1007: There is no theory to support this option
1009: .seealso: [](ch_ts), `TSAdapt`
1010: @*/
1011: PetscErrorCode TSAdaptSetTimeStepIncreaseDelay(TSAdapt adapt, PetscInt cnt)
1012: {
1013: PetscFunctionBegin;
1014: adapt->timestepjustdecreased_delay = cnt;
1015: PetscFunctionReturn(PETSC_SUCCESS);
1016: }
1018: /*@
1019: TSAdaptCheckStage - checks whether to accept a stage, (e.g. reject and change time step size if nonlinear solve fails or solution vector is infeasible)
1021: Collective
1023: Input Parameters:
1024: + adapt - adaptive controller context
1025: . ts - time stepper
1026: . t - Current simulation time
1027: - Y - Current solution vector
1029: Output Parameter:
1030: . accept - `PETSC_TRUE` to accept the stage, `PETSC_FALSE` to reject
1032: Level: developer
1034: .seealso: [](ch_ts), `TSAdapt`
1035: @*/
1036: PetscErrorCode TSAdaptCheckStage(TSAdapt adapt, TS ts, PetscReal t, Vec Y, PetscBool *accept)
1037: {
1038: SNESConvergedReason snesreason = SNES_CONVERGED_ITERATING;
1039: PetscBool func_accept, snes_div_func;
1041: PetscFunctionBegin;
1044: PetscAssertPointer(accept, 5);
1046: PetscCall(TSFunctionDomainError(ts, t, Y, &func_accept));
1047: if (ts->snes) PetscCall(SNESGetConvergedReason(ts->snes, &snesreason));
1048: snes_div_func = (PetscBool)(snesreason == SNES_DIVERGED_FUNCTION_DOMAIN);
1049: if (func_accept && snesreason < 0 && !snes_div_func) {
1050: *accept = PETSC_FALSE;
1051: PetscCall(PetscInfo(ts, "Step=%" PetscInt_FMT ", nonlinear solve failure: %s\n", ts->steps, SNESConvergedReasons[snesreason]));
1052: if (++ts->num_snes_failures >= ts->max_snes_failures && ts->max_snes_failures != PETSC_UNLIMITED) {
1053: ts->reason = TS_DIVERGED_NONLINEAR_SOLVE;
1054: PetscCall(PetscInfo(ts, "Step=%" PetscInt_FMT ", nonlinear solve failures %" PetscInt_FMT " greater than current TS allowed, stopping solve\n", ts->steps, ts->num_snes_failures));
1055: if (adapt->monitor) {
1056: PetscCall(PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
1057: PetscCall(PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt %s step %3" PetscInt_FMT " stage rejected t=%-11g+%10.3e, nonlinear solve failures %" PetscInt_FMT " greater than current TS allowed\n", ((PetscObject)adapt)->type_name, ts->steps,
1058: (double)ts->ptime, (double)ts->time_step, ts->num_snes_failures));
1059: PetscCall(PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
1060: }
1061: }
1062: } else {
1063: *accept = (PetscBool)(func_accept && !snes_div_func);
1064: if (*accept && adapt->checkstage) PetscCall((*adapt->checkstage)(adapt, ts, t, Y, accept));
1065: if (!*accept) {
1066: const char *user_func = !func_accept ? "TSSetFunctionDomainError()" : "TSAdaptSetCheckStage";
1067: const char *snes_err = "SNES invalid function domain";
1068: const char *err_msg = snes_div_func && func_accept ? snes_err : user_func;
1069: PetscCall(PetscInfo(ts, "Step=%" PetscInt_FMT ", solution rejected by %s\n", ts->steps, err_msg));
1070: if (adapt->monitor) {
1071: PetscCall(PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
1072: PetscCall(PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt %s step %3" PetscInt_FMT " stage rejected by %s\n", ((PetscObject)adapt)->type_name, ts->steps, err_msg));
1073: PetscCall(PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
1074: }
1075: }
1076: }
1078: if (!*accept && !ts->reason) {
1079: PetscReal dt, new_dt;
1080: PetscCall(TSGetTimeStep(ts, &dt));
1081: new_dt = dt * adapt->scale_solve_failed;
1082: PetscCall(TSSetTimeStep(ts, new_dt));
1083: adapt->timestepjustdecreased += adapt->timestepjustdecreased_delay;
1084: if (adapt->monitor) {
1085: PetscCall(PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
1086: PetscCall(PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt %s step %3" PetscInt_FMT " stage rejected (SNES reason %s) t=%-11g+%10.3e retrying with dt=%-10.3e\n", ((PetscObject)adapt)->type_name, ts->steps, SNESConvergedReasons[snesreason],
1087: (double)ts->ptime, (double)dt, (double)new_dt));
1088: PetscCall(PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel));
1089: }
1090: }
1091: PetscFunctionReturn(PETSC_SUCCESS);
1092: }
1094: /*@
1095: TSAdaptCreate - create an adaptive controller context for time stepping
1097: Collective
1099: Input Parameter:
1100: . comm - The communicator
1102: Output Parameter:
1103: . inadapt - new `TSAdapt` object
1105: Level: developer
1107: Note:
1108: `TSAdapt` creation is handled by `TS`, so users should not need to call this function.
1110: .seealso: [](ch_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptDestroy()`
1111: @*/
1112: PetscErrorCode TSAdaptCreate(MPI_Comm comm, TSAdapt *inadapt)
1113: {
1114: TSAdapt adapt;
1116: PetscFunctionBegin;
1117: PetscAssertPointer(inadapt, 2);
1118: PetscCall(TSAdaptInitializePackage());
1120: PetscCall(PetscHeaderCreate(adapt, TSADAPT_CLASSID, "TSAdapt", "Time stepping adaptivity", "TS", comm, TSAdaptDestroy, TSAdaptView));
1121: adapt->always_accept = PETSC_FALSE;
1122: adapt->safety = 0.9;
1123: adapt->reject_safety = 0.5;
1124: adapt->clip[0] = 0.1;
1125: adapt->clip[1] = 10.;
1126: adapt->dt_min = 1e-20;
1127: adapt->dt_max = 1e+20;
1128: adapt->ignore_max = -1.0;
1129: adapt->glee_use_local = PETSC_TRUE;
1130: adapt->scale_solve_failed = 0.25;
1131: /* these two safety factors are not public, and they are used only in the TS_EXACTFINALTIME_MATCHSTEP case
1132: to prevent from situations were unreasonably small time steps are taken in order to match the final time */
1133: adapt->matchstepfac[0] = 0.01; /* allow 1% step size increase in the last step */
1134: adapt->matchstepfac[1] = 2.0; /* halve last step if it is greater than what remains divided this factor */
1135: adapt->wnormtype = NORM_2;
1136: adapt->timestepjustdecreased_delay = 0;
1137: *inadapt = adapt;
1138: PetscFunctionReturn(PETSC_SUCCESS);
1139: }