Actual source code: tsadapt.c
2: #include <petsc/private/tsimpl.h>
4: PetscClassId TSADAPT_CLASSID;
6: static PetscFunctionList TSAdaptList;
7: static PetscBool TSAdaptPackageInitialized;
8: static PetscBool TSAdaptRegisterAllCalled;
10: PETSC_EXTERN PetscErrorCode TSAdaptCreate_None(TSAdapt);
11: PETSC_EXTERN PetscErrorCode TSAdaptCreate_Basic(TSAdapt);
12: PETSC_EXTERN PetscErrorCode TSAdaptCreate_DSP(TSAdapt);
13: PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt);
14: PETSC_EXTERN PetscErrorCode TSAdaptCreate_GLEE(TSAdapt);
15: PETSC_EXTERN PetscErrorCode TSAdaptCreate_History(TSAdapt);
17: /*@C
18: TSAdaptRegister - adds a TSAdapt implementation
20: Not Collective
22: Input Parameters:
23: + name_scheme - name of user-defined adaptivity scheme
24: - routine_create - routine to create method context
26: Notes:
27: TSAdaptRegister() may be called multiple times to add several user-defined families.
29: Sample usage:
30: .vb
31: TSAdaptRegister("my_scheme",MySchemeCreate);
32: .ve
34: Then, your scheme can be chosen with the procedural interface via
35: $ TSAdaptSetType(ts,"my_scheme")
36: or at runtime via the option
37: $ -ts_adapt_type my_scheme
39: Level: advanced
41: .seealso: `TSAdaptRegisterAll()`
42: @*/
43: PetscErrorCode TSAdaptRegister(const char sname[], PetscErrorCode (*function)(TSAdapt))
44: {
45: TSAdaptInitializePackage();
46: PetscFunctionListAdd(&TSAdaptList, sname, function);
47: return 0;
48: }
50: /*@C
51: TSAdaptRegisterAll - Registers all of the adaptivity schemes in TSAdapt
53: Not Collective
55: Level: advanced
57: .seealso: `TSAdaptRegisterDestroy()`
58: @*/
59: PetscErrorCode TSAdaptRegisterAll(void)
60: {
61: if (TSAdaptRegisterAllCalled) return 0;
62: TSAdaptRegisterAllCalled = PETSC_TRUE;
63: TSAdaptRegister(TSADAPTNONE, TSAdaptCreate_None);
64: TSAdaptRegister(TSADAPTBASIC, TSAdaptCreate_Basic);
65: TSAdaptRegister(TSADAPTDSP, TSAdaptCreate_DSP);
66: TSAdaptRegister(TSADAPTCFL, TSAdaptCreate_CFL);
67: TSAdaptRegister(TSADAPTGLEE, TSAdaptCreate_GLEE);
68: TSAdaptRegister(TSADAPTHISTORY, TSAdaptCreate_History);
69: return 0;
70: }
72: /*@C
73: TSAdaptFinalizePackage - This function destroys everything in the TS package. It is
74: called from PetscFinalize().
76: Level: developer
78: .seealso: `PetscFinalize()`
79: @*/
80: PetscErrorCode TSAdaptFinalizePackage(void)
81: {
82: PetscFunctionListDestroy(&TSAdaptList);
83: TSAdaptPackageInitialized = PETSC_FALSE;
84: TSAdaptRegisterAllCalled = PETSC_FALSE;
85: return 0;
86: }
88: /*@C
89: TSAdaptInitializePackage - This function initializes everything in the TSAdapt package. It is
90: called from TSInitializePackage().
92: Level: developer
94: .seealso: `PetscInitialize()`
95: @*/
96: PetscErrorCode TSAdaptInitializePackage(void)
97: {
98: if (TSAdaptPackageInitialized) return 0;
99: TSAdaptPackageInitialized = PETSC_TRUE;
100: PetscClassIdRegister("TSAdapt", &TSADAPT_CLASSID);
101: TSAdaptRegisterAll();
102: PetscRegisterFinalize(TSAdaptFinalizePackage);
103: return 0;
104: }
106: /*@C
107: TSAdaptSetType - sets the approach used for the error adapter, currently there is only TSADAPTBASIC and TSADAPTNONE
109: Logicially Collective on TSAdapt
111: Input Parameters:
112: + adapt - the TS adapter, most likely obtained with TSGetAdapt()
113: - type - either TSADAPTBASIC or TSADAPTNONE
115: Options Database:
116: . -ts_adapt_type <basic or dsp or none> - to set the adapter type
118: Level: intermediate
120: .seealso: `TSGetAdapt()`, `TSAdaptDestroy()`, `TSAdaptType`, `TSAdaptGetType()`
121: @*/
122: PetscErrorCode TSAdaptSetType(TSAdapt adapt, TSAdaptType type)
123: {
124: PetscBool match;
125: PetscErrorCode (*r)(TSAdapt);
129: PetscObjectTypeCompare((PetscObject)adapt, type, &match);
130: if (match) return 0;
131: PetscFunctionListFind(TSAdaptList, type, &r);
133: PetscTryTypeMethod(adapt, destroy);
134: PetscMemzero(adapt->ops, sizeof(struct _TSAdaptOps));
135: PetscObjectChangeTypeName((PetscObject)adapt, type);
136: (*r)(adapt);
137: return 0;
138: }
140: /*@C
141: TSAdaptGetType - gets the TS adapter method type (as a string).
143: Not Collective
145: Input Parameter:
146: . adapt - The TS adapter, most likely obtained with TSGetAdapt()
148: Output Parameter:
149: . type - The name of TS adapter method
151: Level: intermediate
153: .seealso `TSAdaptSetType()`
154: @*/
155: PetscErrorCode TSAdaptGetType(TSAdapt adapt, TSAdaptType *type)
156: {
159: *type = ((PetscObject)adapt)->type_name;
160: return 0;
161: }
163: PetscErrorCode TSAdaptSetOptionsPrefix(TSAdapt adapt, const char prefix[])
164: {
166: PetscObjectSetOptionsPrefix((PetscObject)adapt, prefix);
167: return 0;
168: }
170: /*@C
171: TSAdaptLoad - Loads a TSAdapt that has been stored in binary with TSAdaptView().
173: Collective on PetscViewer
175: Input Parameters:
176: + newdm - the newly loaded TSAdapt, this needs to have been created with TSAdaptCreate() or
177: some related function before a call to TSAdaptLoad().
178: - viewer - binary file viewer, obtained from PetscViewerBinaryOpen() or
179: HDF5 file viewer, obtained from PetscViewerHDF5Open()
181: Level: intermediate
183: Notes:
184: The type is determined by the data in the file, any type set into the TSAdapt before this call is ignored.
186: Notes for advanced users:
187: Most users should not need to know the details of the binary storage
188: format, since TSAdaptLoad() and TSAdaptView() completely hide these details.
189: But for anyone who's interested, the standard binary matrix storage
190: format is
191: .vb
192: has not yet been determined
193: .ve
195: .seealso: `PetscViewerBinaryOpen()`, `TSAdaptView()`, `MatLoad()`, `VecLoad()`
196: @*/
197: PetscErrorCode TSAdaptLoad(TSAdapt adapt, PetscViewer viewer)
198: {
199: PetscBool isbinary;
200: char type[256];
204: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
207: PetscViewerBinaryRead(viewer, type, 256, NULL, PETSC_CHAR);
208: TSAdaptSetType(adapt, type);
209: PetscTryTypeMethod(adapt, load, viewer);
210: return 0;
211: }
213: PetscErrorCode TSAdaptView(TSAdapt adapt, PetscViewer viewer)
214: {
215: PetscBool iascii, isbinary, isnone, isglee;
218: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)adapt), &viewer);
221: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
222: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
223: if (iascii) {
224: PetscObjectPrintClassNamePrefixType((PetscObject)adapt, viewer);
225: PetscObjectTypeCompare((PetscObject)adapt, TSADAPTNONE, &isnone);
226: PetscObjectTypeCompare((PetscObject)adapt, TSADAPTGLEE, &isglee);
227: if (!isnone) {
228: if (adapt->always_accept) PetscViewerASCIIPrintf(viewer, " always accepting steps\n");
229: PetscViewerASCIIPrintf(viewer, " safety factor %g\n", (double)adapt->safety);
230: PetscViewerASCIIPrintf(viewer, " extra safety factor after step rejection %g\n", (double)adapt->reject_safety);
231: PetscViewerASCIIPrintf(viewer, " clip fastest increase %g\n", (double)adapt->clip[1]);
232: PetscViewerASCIIPrintf(viewer, " clip fastest decrease %g\n", (double)adapt->clip[0]);
233: PetscViewerASCIIPrintf(viewer, " maximum allowed timestep %g\n", (double)adapt->dt_max);
234: PetscViewerASCIIPrintf(viewer, " minimum allowed timestep %g\n", (double)adapt->dt_min);
235: PetscViewerASCIIPrintf(viewer, " maximum solution absolute value to be ignored %g\n", (double)adapt->ignore_max);
236: }
237: if (isglee) {
238: if (adapt->glee_use_local) {
239: PetscViewerASCIIPrintf(viewer, " GLEE uses local error control\n");
240: } else {
241: PetscViewerASCIIPrintf(viewer, " GLEE uses global error control\n");
242: }
243: }
244: PetscViewerASCIIPushTab(viewer);
245: PetscTryTypeMethod(adapt, view, viewer);
246: PetscViewerASCIIPopTab(viewer);
247: } else if (isbinary) {
248: char type[256];
250: /* need to save FILE_CLASS_ID for adapt class */
251: PetscStrncpy(type, ((PetscObject)adapt)->type_name, 256);
252: PetscViewerBinaryWrite(viewer, type, 256, PETSC_CHAR);
253: } else PetscTryTypeMethod(adapt, view, viewer);
254: return 0;
255: }
257: /*@
258: TSAdaptReset - Resets a TSAdapt context.
260: Collective on TS
262: Input Parameter:
263: . adapt - the TSAdapt context obtained from TSAdaptCreate()
265: Level: developer
267: .seealso: `TSAdaptCreate()`, `TSAdaptDestroy()`
268: @*/
269: PetscErrorCode TSAdaptReset(TSAdapt adapt)
270: {
272: PetscTryTypeMethod(adapt, reset);
273: return 0;
274: }
276: PetscErrorCode TSAdaptDestroy(TSAdapt *adapt)
277: {
278: if (!*adapt) return 0;
280: if (--((PetscObject)(*adapt))->refct > 0) {
281: *adapt = NULL;
282: return 0;
283: }
285: TSAdaptReset(*adapt);
287: PetscTryTypeMethod((*adapt), destroy);
288: PetscViewerDestroy(&(*adapt)->monitor);
289: PetscHeaderDestroy(adapt);
290: return 0;
291: }
293: /*@
294: TSAdaptSetMonitor - Monitor the choices made by the adaptive controller
296: Collective on TSAdapt
298: Input Parameters:
299: + adapt - adaptive controller context
300: - flg - PETSC_TRUE to active a monitor, PETSC_FALSE to disable
302: Options Database Keys:
303: . -ts_adapt_monitor - to turn on monitoring
305: Level: intermediate
307: .seealso: `TSAdaptChoose()`
308: @*/
309: PetscErrorCode TSAdaptSetMonitor(TSAdapt adapt, PetscBool flg)
310: {
313: if (flg) {
314: if (!adapt->monitor) PetscViewerASCIIOpen(PetscObjectComm((PetscObject)adapt), "stdout", &adapt->monitor);
315: } else {
316: PetscViewerDestroy(&adapt->monitor);
317: }
318: return 0;
319: }
321: /*@C
322: TSAdaptSetCheckStage - Set a callback to check convergence for a stage
324: Logically collective on TSAdapt
326: Input Parameters:
327: + adapt - adaptive controller context
328: - func - stage check function
330: Arguments of func:
331: $ PetscErrorCode func(TSAdapt adapt,TS ts,PetscBool *accept)
333: + adapt - adaptive controller context
334: . ts - time stepping context
335: - accept - pending choice of whether to accept, can be modified by this routine
337: Level: advanced
339: .seealso: `TSAdaptChoose()`
340: @*/
341: PetscErrorCode TSAdaptSetCheckStage(TSAdapt adapt, PetscErrorCode (*func)(TSAdapt, TS, PetscReal, Vec, PetscBool *))
342: {
344: adapt->checkstage = func;
345: return 0;
346: }
348: /*@
349: TSAdaptSetAlwaysAccept - Set whether to always accept steps regardless of
350: any error or stability condition not meeting the prescribed goal.
352: Logically collective on TSAdapt
354: Input Parameters:
355: + adapt - time step adaptivity context, usually gotten with TSGetAdapt()
356: - flag - whether to always accept steps
358: Options Database Keys:
359: . -ts_adapt_always_accept - to always accept steps
361: Level: intermediate
363: .seealso: `TSAdapt`, `TSAdaptChoose()`
364: @*/
365: PetscErrorCode TSAdaptSetAlwaysAccept(TSAdapt adapt, PetscBool flag)
366: {
369: adapt->always_accept = flag;
370: return 0;
371: }
373: /*@
374: TSAdaptSetSafety - Set safety factors
376: Logically collective on TSAdapt
378: Input Parameters:
379: + adapt - adaptive controller context
380: . safety - safety factor relative to target error/stability goal
381: - reject_safety - extra safety factor to apply if the last step was rejected
383: Options Database Keys:
384: + -ts_adapt_safety <safety> - to set safety factor
385: - -ts_adapt_reject_safety <reject_safety> - to set reject safety factor
387: Level: intermediate
389: .seealso: `TSAdapt`, `TSAdaptGetSafety()`, `TSAdaptChoose()`
390: @*/
391: PetscErrorCode TSAdaptSetSafety(TSAdapt adapt, PetscReal safety, PetscReal reject_safety)
392: {
400: if (safety != PETSC_DEFAULT) adapt->safety = safety;
401: if (reject_safety != PETSC_DEFAULT) adapt->reject_safety = reject_safety;
402: return 0;
403: }
405: /*@
406: TSAdaptGetSafety - Get safety factors
408: Not Collective
410: Input Parameter:
411: . adapt - adaptive controller context
413: Output Parameters:
414: . safety - safety factor relative to target error/stability goal
415: + reject_safety - extra safety factor to apply if the last step was rejected
417: Level: intermediate
419: .seealso: `TSAdapt`, `TSAdaptSetSafety()`, `TSAdaptChoose()`
420: @*/
421: PetscErrorCode TSAdaptGetSafety(TSAdapt adapt, PetscReal *safety, PetscReal *reject_safety)
422: {
426: if (safety) *safety = adapt->safety;
427: if (reject_safety) *reject_safety = adapt->reject_safety;
428: return 0;
429: }
431: /*@
432: TSAdaptSetMaxIgnore - Set error estimation threshold. Solution components below this threshold value will not be considered when computing error norms for time step adaptivity (in absolute value). A negative value (default) of the threshold leads to considering all solution components.
434: Logically collective on TSAdapt
436: Input Parameters:
437: + adapt - adaptive controller context
438: - max_ignore - threshold for solution components that are ignored during error estimation
440: Options Database Keys:
441: . -ts_adapt_max_ignore <max_ignore> - to set the threshold
443: Level: intermediate
445: .seealso: `TSAdapt`, `TSAdaptGetMaxIgnore()`, `TSAdaptChoose()`
446: @*/
447: PetscErrorCode TSAdaptSetMaxIgnore(TSAdapt adapt, PetscReal max_ignore)
448: {
451: adapt->ignore_max = max_ignore;
452: return 0;
453: }
455: /*@
456: TSAdaptGetMaxIgnore - Get error estimation threshold. Solution components below this threshold value will not be considered when computing error norms for time step adaptivity (in absolute value).
458: Not Collective
460: Input Parameter:
461: . adapt - adaptive controller context
463: Output Parameter:
464: . max_ignore - threshold for solution components that are ignored during error estimation
466: Level: intermediate
468: .seealso: `TSAdapt`, `TSAdaptSetMaxIgnore()`, `TSAdaptChoose()`
469: @*/
470: PetscErrorCode TSAdaptGetMaxIgnore(TSAdapt adapt, PetscReal *max_ignore)
471: {
474: *max_ignore = adapt->ignore_max;
475: return 0;
476: }
478: /*@
479: TSAdaptSetClip - Sets the admissible decrease/increase factor in step size
481: Logically collective on TSAdapt
483: Input Parameters:
484: + adapt - adaptive controller context
485: . low - admissible decrease factor
486: - high - admissible increase factor
488: Options Database Keys:
489: . -ts_adapt_clip <low>,<high> - to set admissible time step decrease and increase factors
491: Level: intermediate
493: .seealso: `TSAdaptChoose()`, `TSAdaptGetClip()`, `TSAdaptSetScaleSolveFailed()`
494: @*/
495: PetscErrorCode TSAdaptSetClip(TSAdapt adapt, PetscReal low, PetscReal high)
496: {
503: if (low != PETSC_DEFAULT) adapt->clip[0] = low;
504: if (high != PETSC_DEFAULT) adapt->clip[1] = high;
505: return 0;
506: }
508: /*@
509: TSAdaptGetClip - Gets the admissible decrease/increase factor in step size
511: Not Collective
513: Input Parameter:
514: . adapt - adaptive controller context
516: Output Parameters:
517: + low - optional, admissible decrease factor
518: - high - optional, admissible increase factor
520: Level: intermediate
522: .seealso: `TSAdaptChoose()`, `TSAdaptSetClip()`, `TSAdaptSetScaleSolveFailed()`
523: @*/
524: PetscErrorCode TSAdaptGetClip(TSAdapt adapt, PetscReal *low, PetscReal *high)
525: {
529: if (low) *low = adapt->clip[0];
530: if (high) *high = adapt->clip[1];
531: return 0;
532: }
534: /*@
535: TSAdaptSetScaleSolveFailed - Scale step by this factor if solve fails
537: Logically collective on TSAdapt
539: Input Parameters:
540: + adapt - adaptive controller context
541: - scale - scale
543: Options Database Keys:
544: . -ts_adapt_scale_solve_failed <scale> - to set scale step by this factor if solve fails
546: Level: intermediate
548: .seealso: `TSAdaptChoose()`, `TSAdaptGetScaleSolveFailed()`, `TSAdaptGetClip()`
549: @*/
550: PetscErrorCode TSAdaptSetScaleSolveFailed(TSAdapt adapt, PetscReal scale)
551: {
556: if (scale != PETSC_DEFAULT) adapt->scale_solve_failed = scale;
557: return 0;
558: }
560: /*@
561: TSAdaptGetScaleSolveFailed - Gets the admissible decrease/increase factor in step size
563: Not Collective
565: Input Parameter:
566: . adapt - adaptive controller context
568: Output Parameter:
569: . scale - scale factor
571: Level: intermediate
573: .seealso: `TSAdaptChoose()`, `TSAdaptSetScaleSolveFailed()`, `TSAdaptSetClip()`
574: @*/
575: PetscErrorCode TSAdaptGetScaleSolveFailed(TSAdapt adapt, PetscReal *scale)
576: {
579: if (scale) *scale = adapt->scale_solve_failed;
580: return 0;
581: }
583: /*@
584: TSAdaptSetStepLimits - Set the minimum and maximum step sizes to be considered by the controller
586: Logically collective on TSAdapt
588: Input Parameters:
589: + adapt - time step adaptivity context, usually gotten with TSGetAdapt()
590: . hmin - minimum time step
591: - hmax - maximum time step
593: Options Database Keys:
594: + -ts_adapt_dt_min <min> - to set minimum time step
595: - -ts_adapt_dt_max <max> - to set maximum time step
597: Level: intermediate
599: .seealso: `TSAdapt`, `TSAdaptGetStepLimits()`, `TSAdaptChoose()`
600: @*/
601: PetscErrorCode TSAdaptSetStepLimits(TSAdapt adapt, PetscReal hmin, PetscReal hmax)
602: {
608: if (hmin != PETSC_DEFAULT) adapt->dt_min = hmin;
609: if (hmax != PETSC_DEFAULT) adapt->dt_max = hmax;
610: hmin = adapt->dt_min;
611: hmax = adapt->dt_max;
613: return 0;
614: }
616: /*@
617: TSAdaptGetStepLimits - Get the minimum and maximum step sizes to be considered by the controller
619: Not Collective
621: Input Parameter:
622: . adapt - time step adaptivity context, usually gotten with TSGetAdapt()
624: Output Parameters:
625: + hmin - minimum time step
626: - hmax - maximum time step
628: Level: intermediate
630: .seealso: `TSAdapt`, `TSAdaptSetStepLimits()`, `TSAdaptChoose()`
631: @*/
632: PetscErrorCode TSAdaptGetStepLimits(TSAdapt adapt, PetscReal *hmin, PetscReal *hmax)
633: {
637: if (hmin) *hmin = adapt->dt_min;
638: if (hmax) *hmax = adapt->dt_max;
639: return 0;
640: }
642: /*
643: TSAdaptSetFromOptions - Sets various TSAdapt parameters from user options.
645: Collective on TSAdapt
647: Input Parameter:
648: . adapt - the TSAdapt context
650: Options Database Keys:
651: + -ts_adapt_type <type> - algorithm to use for adaptivity
652: . -ts_adapt_always_accept - always accept steps regardless of error/stability goals
653: . -ts_adapt_safety <safety> - safety factor relative to target error/stability goal
654: . -ts_adapt_reject_safety <safety> - extra safety factor to apply if the last step was rejected
655: . -ts_adapt_clip <low,high> - admissible time step decrease and increase factors
656: . -ts_adapt_dt_min <min> - minimum timestep to use
657: . -ts_adapt_dt_max <max> - maximum timestep to use
658: . -ts_adapt_scale_solve_failed <scale> - scale timestep by this factor if a solve fails
659: . -ts_adapt_wnormtype <2 or infinity> - type of norm for computing error estimates
660: - -ts_adapt_time_step_increase_delay - number of timesteps to delay increasing the time step after it has been decreased due to failed solver
662: Level: advanced
664: Notes:
665: This function is automatically called by TSSetFromOptions()
667: .seealso: `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptSetAlwaysAccept()`, `TSAdaptSetSafety()`,
668: `TSAdaptSetClip()`, `TSAdaptSetScaleSolveFailed()`, `TSAdaptSetStepLimits()`, `TSAdaptSetMonitor()`
669: */
670: PetscErrorCode TSAdaptSetFromOptions(TSAdapt adapt, PetscOptionItems *PetscOptionsObject)
671: {
672: char type[256] = TSADAPTBASIC;
673: PetscReal safety, reject_safety, clip[2], scale, hmin, hmax;
674: PetscBool set, flg;
675: PetscInt two;
678: /* This should use PetscOptionsBegin() if/when this becomes an object used outside of TS, but currently this
679: * function can only be called from inside TSSetFromOptions() */
680: PetscOptionsHeadBegin(PetscOptionsObject, "TS Adaptivity options");
681: PetscOptionsFList("-ts_adapt_type", "Algorithm to use for adaptivity", "TSAdaptSetType", TSAdaptList, ((PetscObject)adapt)->type_name ? ((PetscObject)adapt)->type_name : type, type, sizeof(type), &flg);
682: if (flg || !((PetscObject)adapt)->type_name) TSAdaptSetType(adapt, type);
684: PetscOptionsBool("-ts_adapt_always_accept", "Always accept the step", "TSAdaptSetAlwaysAccept", adapt->always_accept, &flg, &set);
685: if (set) TSAdaptSetAlwaysAccept(adapt, flg);
687: safety = adapt->safety;
688: reject_safety = adapt->reject_safety;
689: PetscOptionsReal("-ts_adapt_safety", "Safety factor relative to target error/stability goal", "TSAdaptSetSafety", safety, &safety, &set);
690: PetscOptionsReal("-ts_adapt_reject_safety", "Extra safety factor to apply if the last step was rejected", "TSAdaptSetSafety", reject_safety, &reject_safety, &flg);
691: if (set || flg) TSAdaptSetSafety(adapt, safety, reject_safety);
693: two = 2;
694: clip[0] = adapt->clip[0];
695: clip[1] = adapt->clip[1];
696: PetscOptionsRealArray("-ts_adapt_clip", "Admissible decrease/increase factor in step size", "TSAdaptSetClip", clip, &two, &set);
698: if (set) TSAdaptSetClip(adapt, clip[0], clip[1]);
700: hmin = adapt->dt_min;
701: hmax = adapt->dt_max;
702: PetscOptionsReal("-ts_adapt_dt_min", "Minimum time step considered", "TSAdaptSetStepLimits", hmin, &hmin, &set);
703: PetscOptionsReal("-ts_adapt_dt_max", "Maximum time step considered", "TSAdaptSetStepLimits", hmax, &hmax, &flg);
704: if (set || flg) TSAdaptSetStepLimits(adapt, hmin, hmax);
706: PetscOptionsReal("-ts_adapt_max_ignore", "Adaptor ignores (absolute) solution values smaller than this value", "", adapt->ignore_max, &adapt->ignore_max, &set);
707: PetscOptionsBool("-ts_adapt_glee_use_local", "GLEE adaptor uses local error estimation for step control", "", adapt->glee_use_local, &adapt->glee_use_local, &set);
709: PetscOptionsReal("-ts_adapt_scale_solve_failed", "Scale step by this factor if solve fails", "TSAdaptSetScaleSolveFailed", adapt->scale_solve_failed, &scale, &set);
710: if (set) TSAdaptSetScaleSolveFailed(adapt, scale);
712: PetscOptionsEnum("-ts_adapt_wnormtype", "Type of norm computed for error estimation", "", NormTypes, (PetscEnum)adapt->wnormtype, (PetscEnum *)&adapt->wnormtype, NULL);
715: 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);
717: PetscOptionsBool("-ts_adapt_monitor", "Print choices made by adaptive controller", "TSAdaptSetMonitor", adapt->monitor ? PETSC_TRUE : PETSC_FALSE, &flg, &set);
718: if (set) TSAdaptSetMonitor(adapt, flg);
720: PetscTryTypeMethod(adapt, setfromoptions, PetscOptionsObject);
721: PetscOptionsHeadEnd();
722: return 0;
723: }
725: /*@
726: TSAdaptCandidatesClear - clear any previously set candidate schemes
728: Logically collective on TSAdapt
730: Input Parameter:
731: . adapt - adaptive controller
733: Level: developer
735: .seealso: `TSAdapt`, `TSAdaptCreate()`, `TSAdaptCandidateAdd()`, `TSAdaptChoose()`
736: @*/
737: PetscErrorCode TSAdaptCandidatesClear(TSAdapt adapt)
738: {
740: PetscMemzero(&adapt->candidates, sizeof(adapt->candidates));
741: return 0;
742: }
744: /*@C
745: TSAdaptCandidateAdd - add a candidate scheme for the adaptive controller to select from
747: Logically collective on TSAdapt
749: Input Parameters:
750: + adapt - time step adaptivity context, obtained with TSGetAdapt() or TSAdaptCreate()
751: . name - name of the candidate scheme to add
752: . order - order of the candidate scheme
753: . stageorder - stage order of the candidate scheme
754: . ccfl - stability coefficient relative to explicit Euler, used for CFL constraints
755: . cost - relative measure of the amount of work required for the candidate scheme
756: - inuse - indicates that this scheme is the one currently in use, this flag can only be set for one scheme
758: Note:
759: This routine is not available in Fortran.
761: Level: developer
763: .seealso: `TSAdaptCandidatesClear()`, `TSAdaptChoose()`
764: @*/
765: PetscErrorCode TSAdaptCandidateAdd(TSAdapt adapt, const char name[], PetscInt order, PetscInt stageorder, PetscReal ccfl, PetscReal cost, PetscBool inuse)
766: {
767: PetscInt c;
771: if (inuse) {
773: adapt->candidates.inuse_set = PETSC_TRUE;
774: }
775: /* first slot if this is the current scheme, otherwise the next available slot */
776: c = inuse ? 0 : !adapt->candidates.inuse_set + adapt->candidates.n;
778: adapt->candidates.name[c] = name;
779: adapt->candidates.order[c] = order;
780: adapt->candidates.stageorder[c] = stageorder;
781: adapt->candidates.ccfl[c] = ccfl;
782: adapt->candidates.cost[c] = cost;
783: adapt->candidates.n++;
784: return 0;
785: }
787: /*@C
788: TSAdaptCandidatesGet - Get the list of candidate orders of accuracy and cost
790: Not Collective
792: Input Parameter:
793: . adapt - time step adaptivity context
795: Output Parameters:
796: + n - number of candidate schemes, always at least 1
797: . order - the order of each candidate scheme
798: . stageorder - the stage order of each candidate scheme
799: . ccfl - the CFL coefficient of each scheme
800: - cost - the relative cost of each scheme
802: Level: developer
804: Note:
805: The current scheme is always returned in the first slot
807: .seealso: `TSAdaptCandidatesClear()`, `TSAdaptCandidateAdd()`, `TSAdaptChoose()`
808: @*/
809: PetscErrorCode TSAdaptCandidatesGet(TSAdapt adapt, PetscInt *n, const PetscInt **order, const PetscInt **stageorder, const PetscReal **ccfl, const PetscReal **cost)
810: {
812: if (n) *n = adapt->candidates.n;
813: if (order) *order = adapt->candidates.order;
814: if (stageorder) *stageorder = adapt->candidates.stageorder;
815: if (ccfl) *ccfl = adapt->candidates.ccfl;
816: if (cost) *cost = adapt->candidates.cost;
817: return 0;
818: }
820: /*@C
821: TSAdaptChoose - choose which method and step size to use for the next step
823: Collective on TSAdapt
825: Input Parameters:
826: + adapt - adaptive contoller
827: . ts - time stepper
828: - h - current step size
830: Output Parameters:
831: + next_sc - optional, scheme to use for the next step
832: . next_h - step size to use for the next step
833: - accept - PETSC_TRUE to accept the current step, PETSC_FALSE to repeat the current step with the new step size
835: Note:
836: The input value of parameter accept is retained from the last time step, so it will be PETSC_FALSE if the step is
837: being retried after an initial rejection.
839: Level: developer
841: .seealso: `TSAdapt`, `TSAdaptCandidatesClear()`, `TSAdaptCandidateAdd()`
842: @*/
843: PetscErrorCode TSAdaptChoose(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept)
844: {
845: PetscInt ncandidates = adapt->candidates.n;
846: PetscInt scheme = 0;
847: PetscReal wlte = -1.0;
848: PetscReal wltea = -1.0;
849: PetscReal wlter = -1.0;
856: if (next_sc) *next_sc = 0;
858: /* Do not mess with adaptivity while handling events*/
859: if (ts->event && ts->event->status != TSEVENT_NONE) {
860: *next_h = h;
861: *accept = PETSC_TRUE;
862: return 0;
863: }
865: PetscUseTypeMethod(adapt, choose, ts, h, &scheme, next_h, accept, &wlte, &wltea, &wlter);
868: if (next_sc) *next_sc = scheme;
870: if (*accept && ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
871: /* Increase/reduce step size if end time of next step is close to or overshoots max time */
872: PetscReal t = ts->ptime + ts->time_step, h = *next_h;
873: PetscReal tend = t + h, tmax, hmax;
874: PetscReal a = (PetscReal)(1.0 + adapt->matchstepfac[0]);
875: PetscReal b = adapt->matchstepfac[1];
877: if (ts->tspan) {
878: if (PetscIsCloseAtTol(t, ts->tspan->span_times[ts->tspan->spanctr], ts->tspan->reltol * h + ts->tspan->abstol, 0)) /* hit a span time point */
879: if (ts->tspan->spanctr + 1 < ts->tspan->num_span_times) tmax = ts->tspan->span_times[ts->tspan->spanctr + 1];
880: else tmax = ts->max_time; /* hit the last span time point */
881: else tmax = ts->tspan->span_times[ts->tspan->spanctr];
882: } else tmax = ts->max_time;
883: hmax = tmax - t;
884: if (t < tmax && tend > tmax) *next_h = hmax;
885: if (t < tmax && tend < tmax && h * b > hmax) *next_h = hmax / 2;
886: if (t < tmax && tend < tmax && h * a > hmax) *next_h = hmax;
887: /* if step size is changed to match a span time point */
888: if (ts->tspan && h != *next_h && !adapt->dt_span_cached) adapt->dt_span_cached = h;
889: /* reset time step after a span time point */
890: if (ts->tspan && h == *next_h && adapt->dt_span_cached && PetscIsCloseAtTol(t, ts->tspan->span_times[ts->tspan->spanctr], ts->tspan->reltol * h + ts->tspan->abstol, 0)) {
891: *next_h = adapt->dt_span_cached;
892: adapt->dt_span_cached = 0;
893: }
894: }
895: if (adapt->monitor) {
896: const char *sc_name = (scheme < ncandidates) ? adapt->candidates.name[scheme] : "";
897: PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
898: if (wlte < 0) {
899: 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",
900: (double)ts->ptime, (double)h, (double)*next_h));
901: } else {
902: 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",
903: (double)ts->ptime, (double)h, (double)*next_h, (double)wlte, (double)wltea, (double)wlter));
904: }
905: PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
906: }
907: return 0;
908: }
910: /*@
911: TSAdaptSetTimeStepIncreaseDelay - The number of timesteps to wait after a decrease in the timestep due to failed solver
912: before increasing the time step.
914: Logicially Collective on TSAdapt
916: Input Parameters:
917: + adapt - adaptive controller context
918: - cnt - the number of timesteps
920: Options Database Key:
921: . -ts_adapt_time_step_increase_delay cnt - number of steps to delay the increase
923: Notes: This is to prevent an adaptor from bouncing back and forth between two nearby timesteps. The default is 0.
924: The successful use of this option is problem dependent
926: Developer Note: there is no theory to support this option
928: Level: advanced
930: .seealso:
931: @*/
932: PetscErrorCode TSAdaptSetTimeStepIncreaseDelay(TSAdapt adapt, PetscInt cnt)
933: {
934: adapt->timestepjustdecreased_delay = cnt;
935: return 0;
936: }
938: /*@
939: TSAdaptCheckStage - checks whether to accept a stage, (e.g. reject and change time step size if nonlinear solve fails or solution vector is infeasible)
941: Collective on TSAdapt
943: Input Parameters:
944: + adapt - adaptive controller context
945: . ts - time stepper
946: . t - Current simulation time
947: - Y - Current solution vector
949: Output Parameter:
950: . accept - PETSC_TRUE to accept the stage, PETSC_FALSE to reject
952: Level: developer
954: .seealso:
955: @*/
956: PetscErrorCode TSAdaptCheckStage(TSAdapt adapt, TS ts, PetscReal t, Vec Y, PetscBool *accept)
957: {
958: SNESConvergedReason snesreason = SNES_CONVERGED_ITERATING;
964: if (ts->snes) SNESGetConvergedReason(ts->snes, &snesreason);
965: if (snesreason < 0) {
966: *accept = PETSC_FALSE;
967: if (++ts->num_snes_failures >= ts->max_snes_failures && ts->max_snes_failures > 0) {
968: ts->reason = TS_DIVERGED_NONLINEAR_SOLVE;
969: PetscInfo(ts, "Step=%" PetscInt_FMT ", nonlinear solve failures %" PetscInt_FMT " greater than current TS allowed, stopping solve\n", ts->steps, ts->num_snes_failures);
970: if (adapt->monitor) {
971: PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
972: 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,
973: (double)ts->ptime, (double)ts->time_step, ts->num_snes_failures));
974: PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
975: }
976: }
977: } else {
978: *accept = PETSC_TRUE;
979: TSFunctionDomainError(ts, t, Y, accept);
980: if (*accept && adapt->checkstage) {
981: (*adapt->checkstage)(adapt, ts, t, Y, accept);
982: if (!*accept) {
983: PetscInfo(ts, "Step=%" PetscInt_FMT ", solution rejected by user function provided by TSSetFunctionDomainError()\n", ts->steps);
984: if (adapt->monitor) {
985: PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
986: PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt %s step %3" PetscInt_FMT " stage rejected by user function provided by TSSetFunctionDomainError()\n", ((PetscObject)adapt)->type_name, ts->steps);
987: PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
988: }
989: }
990: }
991: }
993: if (!(*accept) && !ts->reason) {
994: PetscReal dt, new_dt;
995: TSGetTimeStep(ts, &dt);
996: new_dt = dt * adapt->scale_solve_failed;
997: TSSetTimeStep(ts, new_dt);
998: adapt->timestepjustdecreased += adapt->timestepjustdecreased_delay;
999: if (adapt->monitor) {
1000: PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
1001: PetscViewerASCIIPrintf(adapt->monitor, " TSAdapt %s step %3" PetscInt_FMT " stage rejected (%s) t=%-11g+%10.3e retrying with dt=%-10.3e\n", ((PetscObject)adapt)->type_name, ts->steps, SNESConvergedReasons[snesreason], (double)ts->ptime, (double)dt, (double)new_dt);
1002: PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
1003: }
1004: }
1005: return 0;
1006: }
1008: /*@
1009: TSAdaptCreate - create an adaptive controller context for time stepping
1011: Collective
1013: Input Parameter:
1014: . comm - The communicator
1016: Output Parameter:
1017: . adapt - new TSAdapt object
1019: Level: developer
1021: Notes:
1022: TSAdapt creation is handled by TS, so users should not need to call this function.
1024: .seealso: `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptDestroy()`
1025: @*/
1026: PetscErrorCode TSAdaptCreate(MPI_Comm comm, TSAdapt *inadapt)
1027: {
1028: TSAdapt adapt;
1031: *inadapt = NULL;
1032: TSAdaptInitializePackage();
1034: PetscHeaderCreate(adapt, TSADAPT_CLASSID, "TSAdapt", "Time stepping adaptivity", "TS", comm, TSAdaptDestroy, TSAdaptView);
1036: adapt->always_accept = PETSC_FALSE;
1037: adapt->safety = 0.9;
1038: adapt->reject_safety = 0.5;
1039: adapt->clip[0] = 0.1;
1040: adapt->clip[1] = 10.;
1041: adapt->dt_min = 1e-20;
1042: adapt->dt_max = 1e+20;
1043: adapt->ignore_max = -1.0;
1044: adapt->glee_use_local = PETSC_TRUE;
1045: adapt->scale_solve_failed = 0.25;
1046: /* these two safety factors are not public, and they are used only in the TS_EXACTFINALTIME_MATCHSTEP case
1047: to prevent from situations were unreasonably small time steps are taken in order to match the final time */
1048: adapt->matchstepfac[0] = 0.01; /* allow 1% step size increase in the last step */
1049: adapt->matchstepfac[1] = 2.0; /* halve last step if it is greater than what remains divided this factor */
1050: adapt->wnormtype = NORM_2;
1051: adapt->timestepjustdecreased_delay = 0;
1053: *inadapt = adapt;
1054: return 0;
1055: }