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:    Level: advanced

 28:    Notes:
 29:    `TSAdaptRegister()` may be called multiple times to add several user-defined families.

 31:    Sample usage:
 32: .vb
 33:    TSAdaptRegister("my_scheme",MySchemeCreate);
 34: .ve

 36:    Then, your scheme can be chosen with the procedural interface via
 37: $     TSAdaptSetType(ts,"my_scheme")
 38:    or at runtime via the option
 39: $     -ts_adapt_type my_scheme

 41: .seealso: [](chapter_ts), `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: [](chapter_ts), `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: [](chapter_ts), `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: [](chapter_ts), `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

109:   Logicially Collective onadapt

111:   Input Parameters:
112: + adapt - the `TS` adapter, most likely obtained with `TSGetAdapt()`
113: - type - one of the `TSAdaptType`

115:   Options Database Key:
116: . -ts_adapt_type <basic or dsp or none> - to set the adapter type

118:   Level: intermediate

120: .seealso: [](chapter_ts), `TSGetAdapt()`, `TSAdaptDestroy()`, `TSAdaptType`, `TSAdaptGetType()`, `TSAdaptType`
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: `TSAdapt`, `TSAdaptType`, `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

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:   Note:
184:    The type is determined by the data in the file, any type set into the `TSAdapt` before this call is ignored.

186: .seealso: [](chapter_ts), `PetscViewerBinaryOpen()`, `TSAdaptView()`, `MatLoad()`, `VecLoad()`, `TSAdapt`
187: @*/
188: PetscErrorCode TSAdaptLoad(TSAdapt adapt, PetscViewer viewer)
189: {
190:   PetscBool isbinary;
191:   char      type[256];

195:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);

198:   PetscViewerBinaryRead(viewer, type, 256, NULL, PETSC_CHAR);
199:   TSAdaptSetType(adapt, type);
200:   PetscTryTypeMethod(adapt, load, viewer);
201:   return 0;
202: }

204: PetscErrorCode TSAdaptView(TSAdapt adapt, PetscViewer viewer)
205: {
206:   PetscBool iascii, isbinary, isnone, isglee;

209:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)adapt), &viewer);
212:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
213:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
214:   if (iascii) {
215:     PetscObjectPrintClassNamePrefixType((PetscObject)adapt, viewer);
216:     PetscObjectTypeCompare((PetscObject)adapt, TSADAPTNONE, &isnone);
217:     PetscObjectTypeCompare((PetscObject)adapt, TSADAPTGLEE, &isglee);
218:     if (!isnone) {
219:       if (adapt->always_accept) PetscViewerASCIIPrintf(viewer, "  always accepting steps\n");
220:       PetscViewerASCIIPrintf(viewer, "  safety factor %g\n", (double)adapt->safety);
221:       PetscViewerASCIIPrintf(viewer, "  extra safety factor after step rejection %g\n", (double)adapt->reject_safety);
222:       PetscViewerASCIIPrintf(viewer, "  clip fastest increase %g\n", (double)adapt->clip[1]);
223:       PetscViewerASCIIPrintf(viewer, "  clip fastest decrease %g\n", (double)adapt->clip[0]);
224:       PetscViewerASCIIPrintf(viewer, "  maximum allowed timestep %g\n", (double)adapt->dt_max);
225:       PetscViewerASCIIPrintf(viewer, "  minimum allowed timestep %g\n", (double)adapt->dt_min);
226:       PetscViewerASCIIPrintf(viewer, "  maximum solution absolute value to be ignored %g\n", (double)adapt->ignore_max);
227:     }
228:     if (isglee) {
229:       if (adapt->glee_use_local) {
230:         PetscViewerASCIIPrintf(viewer, "  GLEE uses local error control\n");
231:       } else {
232:         PetscViewerASCIIPrintf(viewer, "  GLEE uses global error control\n");
233:       }
234:     }
235:     PetscViewerASCIIPushTab(viewer);
236:     PetscTryTypeMethod(adapt, view, viewer);
237:     PetscViewerASCIIPopTab(viewer);
238:   } else if (isbinary) {
239:     char type[256];

241:     /* need to save FILE_CLASS_ID for adapt class */
242:     PetscStrncpy(type, ((PetscObject)adapt)->type_name, 256);
243:     PetscViewerBinaryWrite(viewer, type, 256, PETSC_CHAR);
244:   } else PetscTryTypeMethod(adapt, view, viewer);
245:   return 0;
246: }

248: /*@
249:    TSAdaptReset - Resets a `TSAdapt` context to its defaults

251:    Collective

253:    Input Parameter:
254: .  adapt - the `TSAdapt` context obtained from `TSGetAdapt()` or `TSAdaptCreate()`

256:    Level: developer

258: .seealso: [](chapter_ts), `TSGetAdapt()`, `TSAdapt`, `TSAdaptCreate()`, `TSAdaptDestroy()`
259: @*/
260: PetscErrorCode TSAdaptReset(TSAdapt adapt)
261: {
263:   PetscTryTypeMethod(adapt, reset);
264:   return 0;
265: }

267: PetscErrorCode TSAdaptDestroy(TSAdapt *adapt)
268: {
269:   if (!*adapt) return 0;
271:   if (--((PetscObject)(*adapt))->refct > 0) {
272:     *adapt = NULL;
273:     return 0;
274:   }

276:   TSAdaptReset(*adapt);

278:   PetscTryTypeMethod((*adapt), destroy);
279:   PetscViewerDestroy(&(*adapt)->monitor);
280:   PetscHeaderDestroy(adapt);
281:   return 0;
282: }

284: /*@
285:    TSAdaptSetMonitor - Monitor the choices made by the adaptive controller

287:    Collective

289:    Input Parameters:
290: +  adapt - adaptive controller context
291: -  flg - `PETSC_TRUE` to active a monitor, `PETSC_FALSE` to disable

293:    Options Database Key:
294: .  -ts_adapt_monitor - to turn on monitoring

296:    Level: intermediate

298: .seealso: [](chapter_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptChoose()`
299: @*/
300: PetscErrorCode TSAdaptSetMonitor(TSAdapt adapt, PetscBool flg)
301: {
304:   if (flg) {
305:     if (!adapt->monitor) PetscViewerASCIIOpen(PetscObjectComm((PetscObject)adapt), "stdout", &adapt->monitor);
306:   } else {
307:     PetscViewerDestroy(&adapt->monitor);
308:   }
309:   return 0;
310: }

312: /*@C
313:    TSAdaptSetCheckStage - Set a callback to check convergence for a stage

315:    Logically collective

317:    Input Parameters:
318: +  adapt - adaptive controller context
319: -  func - stage check function

321:    Arguments of func:
322: $  PetscErrorCode func(TSAdapt adapt,TS ts,PetscBool *accept)

324: +  adapt - adaptive controller context
325: .  ts - time stepping context
326: -  accept - pending choice of whether to accept, can be modified by this routine

328:    Level: advanced

330: .seealso: [](chapter_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptChoose()`
331: @*/
332: PetscErrorCode TSAdaptSetCheckStage(TSAdapt adapt, PetscErrorCode (*func)(TSAdapt, TS, PetscReal, Vec, PetscBool *))
333: {
335:   adapt->checkstage = func;
336:   return 0;
337: }

339: /*@
340:    TSAdaptSetAlwaysAccept - Set whether to always accept steps regardless of
341:    any error or stability condition not meeting the prescribed goal.

343:    Logically collective

345:    Input Parameters:
346: +  adapt - time step adaptivity context, usually gotten with `TSGetAdapt()`
347: -  flag - whether to always accept steps

349:    Options Database Key:
350: .  -ts_adapt_always_accept - to always accept steps

352:    Level: intermediate

354: .seealso: [](chapter_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptChoose()`
355: @*/
356: PetscErrorCode TSAdaptSetAlwaysAccept(TSAdapt adapt, PetscBool flag)
357: {
360:   adapt->always_accept = flag;
361:   return 0;
362: }

364: /*@
365:    TSAdaptSetSafety - Set safety factors for time step adaptor

367:    Logically collective

369:    Input Parameters:
370: +  adapt - adaptive controller context
371: .  safety - safety factor relative to target error/stability goal
372: -  reject_safety - extra safety factor to apply if the last step was rejected

374:    Options Database Keys:
375: +  -ts_adapt_safety <safety> - to set safety factor
376: -  -ts_adapt_reject_safety <reject_safety> - to set reject safety factor

378:    Level: intermediate

380: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptGetSafety()`, `TSAdaptChoose()`
381: @*/
382: PetscErrorCode TSAdaptSetSafety(TSAdapt adapt, PetscReal safety, PetscReal reject_safety)
383: {
391:   if (safety != PETSC_DEFAULT) adapt->safety = safety;
392:   if (reject_safety != PETSC_DEFAULT) adapt->reject_safety = reject_safety;
393:   return 0;
394: }

396: /*@
397:    TSAdaptGetSafety - Get safety factors for time step adapter

399:    Not Collective

401:    Input Parameter:
402: .  adapt - adaptive controller context

404:    Output Parameters:
405: .  safety - safety factor relative to target error/stability goal
406: +  reject_safety - extra safety factor to apply if the last step was rejected

408:    Level: intermediate

410: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptSetSafety()`, `TSAdaptChoose()`
411: @*/
412: PetscErrorCode TSAdaptGetSafety(TSAdapt adapt, PetscReal *safety, PetscReal *reject_safety)
413: {
417:   if (safety) *safety = adapt->safety;
418:   if (reject_safety) *reject_safety = adapt->reject_safety;
419:   return 0;
420: }

422: /*@
423:    TSAdaptSetMaxIgnore - Set error estimation threshold. Solution components below this threshold value will not be considered when computing error norms
424:    for time step adaptivity (in absolute value). A negative value (default) of the threshold leads to considering all solution components.

426:    Logically collective

428:    Input Parameters:
429: +  adapt - adaptive controller context
430: -  max_ignore - threshold for solution components that are ignored during error estimation

432:    Options Database Key:
433: .  -ts_adapt_max_ignore <max_ignore> - to set the threshold

435:    Level: intermediate

437: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptGetMaxIgnore()`, `TSAdaptChoose()`
438: @*/
439: PetscErrorCode TSAdaptSetMaxIgnore(TSAdapt adapt, PetscReal max_ignore)
440: {
443:   adapt->ignore_max = max_ignore;
444:   return 0;
445: }

447: /*@
448:    TSAdaptGetMaxIgnore - Get 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).

451:    Not Collective

453:    Input Parameter:
454: .  adapt - adaptive controller context

456:    Output Parameter:
457: .  max_ignore - threshold for solution components that are ignored during error estimation

459:    Level: intermediate

461: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptSetMaxIgnore()`, `TSAdaptChoose()`
462: @*/
463: PetscErrorCode TSAdaptGetMaxIgnore(TSAdapt adapt, PetscReal *max_ignore)
464: {
467:   *max_ignore = adapt->ignore_max;
468:   return 0;
469: }

471: /*@
472:    TSAdaptSetClip - Sets the admissible decrease/increase factor in step size in the time step adapter

474:    Logically collective

476:    Input Parameters:
477: +  adapt - adaptive controller context
478: .  low - admissible decrease factor
479: -  high - admissible increase factor

481:    Options Database Key:
482: .  -ts_adapt_clip <low>,<high> - to set admissible time step decrease and increase factors

484:    Level: intermediate

486: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptGetClip()`, `TSAdaptSetScaleSolveFailed()`
487: @*/
488: PetscErrorCode TSAdaptSetClip(TSAdapt adapt, PetscReal low, PetscReal high)
489: {
496:   if (low != PETSC_DEFAULT) adapt->clip[0] = low;
497:   if (high != PETSC_DEFAULT) adapt->clip[1] = high;
498:   return 0;
499: }

501: /*@
502:    TSAdaptGetClip - Gets the admissible decrease/increase factor in step size in the time step adapter

504:    Not Collective

506:    Input Parameter:
507: .  adapt - adaptive controller context

509:    Output Parameters:
510: +  low - optional, admissible decrease factor
511: -  high - optional, admissible increase factor

513:    Level: intermediate

515: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptSetClip()`, `TSAdaptSetScaleSolveFailed()`
516: @*/
517: PetscErrorCode TSAdaptGetClip(TSAdapt adapt, PetscReal *low, PetscReal *high)
518: {
522:   if (low) *low = adapt->clip[0];
523:   if (high) *high = adapt->clip[1];
524:   return 0;
525: }

527: /*@
528:    TSAdaptSetScaleSolveFailed - Scale step size by this factor if solve fails

530:    Logically collective

532:    Input Parameters:
533: +  adapt - adaptive controller context
534: -  scale - scale

536:    Options Database Key:
537: .  -ts_adapt_scale_solve_failed <scale> - to set scale step by this factor if solve fails

539:    Level: intermediate

541: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptGetScaleSolveFailed()`, `TSAdaptGetClip()`
542: @*/
543: PetscErrorCode TSAdaptSetScaleSolveFailed(TSAdapt adapt, PetscReal scale)
544: {
549:   if (scale != PETSC_DEFAULT) adapt->scale_solve_failed = scale;
550:   return 0;
551: }

553: /*@
554:    TSAdaptGetScaleSolveFailed - Gets the admissible decrease/increase factor in step size

556:    Not Collective

558:    Input Parameter:
559: .  adapt - adaptive controller context

561:    Output Parameter:
562: .  scale - scale factor

564:    Level: intermediate

566: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptChoose()`, `TSAdaptSetScaleSolveFailed()`, `TSAdaptSetClip()`
567: @*/
568: PetscErrorCode TSAdaptGetScaleSolveFailed(TSAdapt adapt, PetscReal *scale)
569: {
572:   if (scale) *scale = adapt->scale_solve_failed;
573:   return 0;
574: }

576: /*@
577:    TSAdaptSetStepLimits - Set the minimum and maximum step sizes to be considered by the time step controller

579:    Logically collective

581:    Input Parameters:
582: +  adapt - time step adaptivity context, usually gotten with `TSGetAdapt()`
583: .  hmin - minimum time step
584: -  hmax - maximum time step

586:    Options Database Keys:
587: +  -ts_adapt_dt_min <min> - to set minimum time step
588: -  -ts_adapt_dt_max <max> - to set maximum time step

590:    Level: intermediate

592: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptGetStepLimits()`, `TSAdaptChoose()`
593: @*/
594: PetscErrorCode TSAdaptSetStepLimits(TSAdapt adapt, PetscReal hmin, PetscReal hmax)
595: {
601:   if (hmin != PETSC_DEFAULT) adapt->dt_min = hmin;
602:   if (hmax != PETSC_DEFAULT) adapt->dt_max = hmax;
603:   hmin = adapt->dt_min;
604:   hmax = adapt->dt_max;
606:   return 0;
607: }

609: /*@
610:    TSAdaptGetStepLimits - Get the minimum and maximum step sizes to be considered by the time step controller

612:    Not Collective

614:    Input Parameter:
615: .  adapt - time step adaptivity context, usually gotten with `TSGetAdapt()`

617:    Output Parameters:
618: +  hmin - minimum time step
619: -  hmax - maximum time step

621:    Level: intermediate

623: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptSetStepLimits()`, `TSAdaptChoose()`
624: @*/
625: PetscErrorCode TSAdaptGetStepLimits(TSAdapt adapt, PetscReal *hmin, PetscReal *hmax)
626: {
630:   if (hmin) *hmin = adapt->dt_min;
631:   if (hmax) *hmax = adapt->dt_max;
632:   return 0;
633: }

635: /*
636:    TSAdaptSetFromOptions - Sets various `TSAdapt` parameters from user options.

638:    Collective

640:    Input Parameter:
641: .  adapt - the `TSAdapt` context

643:    Options Database Keys:
644: +  -ts_adapt_type <type> - algorithm to use for adaptivity
645: .  -ts_adapt_always_accept - always accept steps regardless of error/stability goals
646: .  -ts_adapt_safety <safety> - safety factor relative to target error/stability goal
647: .  -ts_adapt_reject_safety <safety> - extra safety factor to apply if the last step was rejected
648: .  -ts_adapt_clip <low,high> - admissible time step decrease and increase factors
649: .  -ts_adapt_dt_min <min> - minimum timestep to use
650: .  -ts_adapt_dt_max <max> - maximum timestep to use
651: .  -ts_adapt_scale_solve_failed <scale> - scale timestep by this factor if a solve fails
652: .  -ts_adapt_wnormtype <2 or infinity> - type of norm for computing error estimates
653: -  -ts_adapt_time_step_increase_delay - number of timesteps to delay increasing the time step after it has been decreased due to failed solver

655:    Level: advanced

657:    Note:
658:    This function is automatically called by `TSSetFromOptions()`

660: .seealso: [](chapter_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptSetAlwaysAccept()`, `TSAdaptSetSafety()`,
661:           `TSAdaptSetClip()`, `TSAdaptSetScaleSolveFailed()`, `TSAdaptSetStepLimits()`, `TSAdaptSetMonitor()`
662: */
663: PetscErrorCode TSAdaptSetFromOptions(TSAdapt adapt, PetscOptionItems *PetscOptionsObject)
664: {
665:   char      type[256] = TSADAPTBASIC;
666:   PetscReal safety, reject_safety, clip[2], scale, hmin, hmax;
667:   PetscBool set, flg;
668:   PetscInt  two;

671:   /* This should use PetscOptionsBegin() if/when this becomes an object used outside of TS, but currently this
672:    * function can only be called from inside TSSetFromOptions()  */
673:   PetscOptionsHeadBegin(PetscOptionsObject, "TS Adaptivity options");
674:   PetscOptionsFList("-ts_adapt_type", "Algorithm to use for adaptivity", "TSAdaptSetType", TSAdaptList, ((PetscObject)adapt)->type_name ? ((PetscObject)adapt)->type_name : type, type, sizeof(type), &flg);
675:   if (flg || !((PetscObject)adapt)->type_name) TSAdaptSetType(adapt, type);

677:   PetscOptionsBool("-ts_adapt_always_accept", "Always accept the step", "TSAdaptSetAlwaysAccept", adapt->always_accept, &flg, &set);
678:   if (set) TSAdaptSetAlwaysAccept(adapt, flg);

680:   safety        = adapt->safety;
681:   reject_safety = adapt->reject_safety;
682:   PetscOptionsReal("-ts_adapt_safety", "Safety factor relative to target error/stability goal", "TSAdaptSetSafety", safety, &safety, &set);
683:   PetscOptionsReal("-ts_adapt_reject_safety", "Extra safety factor to apply if the last step was rejected", "TSAdaptSetSafety", reject_safety, &reject_safety, &flg);
684:   if (set || flg) TSAdaptSetSafety(adapt, safety, reject_safety);

686:   two     = 2;
687:   clip[0] = adapt->clip[0];
688:   clip[1] = adapt->clip[1];
689:   PetscOptionsRealArray("-ts_adapt_clip", "Admissible decrease/increase factor in step size", "TSAdaptSetClip", clip, &two, &set);
691:   if (set) TSAdaptSetClip(adapt, clip[0], clip[1]);

693:   hmin = adapt->dt_min;
694:   hmax = adapt->dt_max;
695:   PetscOptionsReal("-ts_adapt_dt_min", "Minimum time step considered", "TSAdaptSetStepLimits", hmin, &hmin, &set);
696:   PetscOptionsReal("-ts_adapt_dt_max", "Maximum time step considered", "TSAdaptSetStepLimits", hmax, &hmax, &flg);
697:   if (set || flg) TSAdaptSetStepLimits(adapt, hmin, hmax);

699:   PetscOptionsReal("-ts_adapt_max_ignore", "Adaptor ignores (absolute) solution values smaller than this value", "", adapt->ignore_max, &adapt->ignore_max, &set);
700:   PetscOptionsBool("-ts_adapt_glee_use_local", "GLEE adaptor uses local error estimation for step control", "", adapt->glee_use_local, &adapt->glee_use_local, &set);

702:   PetscOptionsReal("-ts_adapt_scale_solve_failed", "Scale step by this factor if solve fails", "TSAdaptSetScaleSolveFailed", adapt->scale_solve_failed, &scale, &set);
703:   if (set) TSAdaptSetScaleSolveFailed(adapt, scale);

705:   PetscOptionsEnum("-ts_adapt_wnormtype", "Type of norm computed for error estimation", "", NormTypes, (PetscEnum)adapt->wnormtype, (PetscEnum *)&adapt->wnormtype, NULL);

708:   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);

710:   PetscOptionsBool("-ts_adapt_monitor", "Print choices made by adaptive controller", "TSAdaptSetMonitor", adapt->monitor ? PETSC_TRUE : PETSC_FALSE, &flg, &set);
711:   if (set) TSAdaptSetMonitor(adapt, flg);

713:   PetscTryTypeMethod(adapt, setfromoptions, PetscOptionsObject);
714:   PetscOptionsHeadEnd();
715:   return 0;
716: }

718: /*@
719:    TSAdaptCandidatesClear - clear any previously set candidate schemes

721:    Logically collective

723:    Input Parameter:
724: .  adapt - adaptive controller

726:    Level: developer

728: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptCreate()`, `TSAdaptCandidateAdd()`, `TSAdaptChoose()`
729: @*/
730: PetscErrorCode TSAdaptCandidatesClear(TSAdapt adapt)
731: {
733:   PetscMemzero(&adapt->candidates, sizeof(adapt->candidates));
734:   return 0;
735: }

737: /*@C
738:    TSAdaptCandidateAdd - add a candidate scheme for the adaptive controller to select from

740:    Logically collective

742:    Input Parameters:
743: +  adapt - time step adaptivity context, obtained with `TSGetAdapt()` or `TSAdaptCreate()`
744: .  name - name of the candidate scheme to add
745: .  order - order of the candidate scheme
746: .  stageorder - stage order of the candidate scheme
747: .  ccfl - stability coefficient relative to explicit Euler, used for CFL constraints
748: .  cost - relative measure of the amount of work required for the candidate scheme
749: -  inuse - indicates that this scheme is the one currently in use, this flag can only be set for one scheme

751:    Level: developer

753:    Fortran Note:
754:    This routine is not available in Fortran.

756: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptCandidatesClear()`, `TSAdaptChoose()`
757: @*/
758: PetscErrorCode TSAdaptCandidateAdd(TSAdapt adapt, const char name[], PetscInt order, PetscInt stageorder, PetscReal ccfl, PetscReal cost, PetscBool inuse)
759: {
760:   PetscInt c;

764:   if (inuse) {
766:     adapt->candidates.inuse_set = PETSC_TRUE;
767:   }
768:   /* first slot if this is the current scheme, otherwise the next available slot */
769:   c = inuse ? 0 : !adapt->candidates.inuse_set + adapt->candidates.n;

771:   adapt->candidates.name[c]       = name;
772:   adapt->candidates.order[c]      = order;
773:   adapt->candidates.stageorder[c] = stageorder;
774:   adapt->candidates.ccfl[c]       = ccfl;
775:   adapt->candidates.cost[c]       = cost;
776:   adapt->candidates.n++;
777:   return 0;
778: }

780: /*@C
781:    TSAdaptCandidatesGet - Get the list of candidate orders of accuracy and cost

783:    Not Collective

785:    Input Parameter:
786: .  adapt - time step adaptivity context

788:    Output Parameters:
789: +  n - number of candidate schemes, always at least 1
790: .  order - the order of each candidate scheme
791: .  stageorder - the stage order of each candidate scheme
792: .  ccfl - the CFL coefficient of each scheme
793: -  cost - the relative cost of each scheme

795:    Level: developer

797:    Note:
798:    The current scheme is always returned in the first slot

800: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptCandidatesClear()`, `TSAdaptCandidateAdd()`, `TSAdaptChoose()`
801: @*/
802: PetscErrorCode TSAdaptCandidatesGet(TSAdapt adapt, PetscInt *n, const PetscInt **order, const PetscInt **stageorder, const PetscReal **ccfl, const PetscReal **cost)
803: {
805:   if (n) *n = adapt->candidates.n;
806:   if (order) *order = adapt->candidates.order;
807:   if (stageorder) *stageorder = adapt->candidates.stageorder;
808:   if (ccfl) *ccfl = adapt->candidates.ccfl;
809:   if (cost) *cost = adapt->candidates.cost;
810:   return 0;
811: }

813: /*@C
814:    TSAdaptChoose - choose which method and step size to use for the next step

816:    Collective

818:    Input Parameters:
819: +  adapt - adaptive controller
820: .  ts - time stepper
821: -  h - current step size

823:    Output Parameters:
824: +  next_sc - optional, scheme to use for the next step
825: .  next_h - step size to use for the next step
826: -  accept - `PETSC_TRUE` to accept the current step, `PETSC_FALSE` to repeat the current step with the new step size

828:    Level: developer

830:    Note:
831:    The input value of parameter accept is retained from the last time step, so it will be `PETSC_FALSE` if the step is
832:    being retried after an initial rejection.

834: .seealso: [](chapter_ts), `TSAdapt`, `TSAdaptCandidatesClear()`, `TSAdaptCandidateAdd()`
835: @*/
836: PetscErrorCode TSAdaptChoose(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept)
837: {
838:   PetscInt  ncandidates = adapt->candidates.n;
839:   PetscInt  scheme      = 0;
840:   PetscReal wlte        = -1.0;
841:   PetscReal wltea       = -1.0;
842:   PetscReal wlter       = -1.0;

849:   if (next_sc) *next_sc = 0;

851:   /* Do not mess with adaptivity while handling events*/
852:   if (ts->event && ts->event->status != TSEVENT_NONE) {
853:     *next_h = h;
854:     *accept = PETSC_TRUE;
855:     return 0;
856:   }

858:   PetscUseTypeMethod(adapt, choose, ts, h, &scheme, next_h, accept, &wlte, &wltea, &wlter);
861:   if (next_sc) *next_sc = scheme;

863:   if (*accept && ts->exact_final_time == TS_EXACTFINALTIME_MATCHSTEP) {
864:     /* Increase/reduce step size if end time of next step is close to or overshoots max time */
865:     PetscReal t = ts->ptime + ts->time_step, h = *next_h;
866:     PetscReal tend = t + h, tmax, hmax;
867:     PetscReal a    = (PetscReal)(1.0 + adapt->matchstepfac[0]);
868:     PetscReal b    = adapt->matchstepfac[1];

870:     if (ts->tspan) {
871:       if (PetscIsCloseAtTol(t, ts->tspan->span_times[ts->tspan->spanctr], ts->tspan->reltol * h + ts->tspan->abstol, 0)) /* hit a span time point */
872:         if (ts->tspan->spanctr + 1 < ts->tspan->num_span_times) tmax = ts->tspan->span_times[ts->tspan->spanctr + 1];
873:         else tmax = ts->max_time; /* hit the last span time point */
874:       else tmax = ts->tspan->span_times[ts->tspan->spanctr];
875:     } else tmax = ts->max_time;
876:     hmax = tmax - t;
877:     if (t < tmax && tend > tmax) *next_h = hmax;
878:     if (t < tmax && tend < tmax && h * b > hmax) *next_h = hmax / 2;
879:     if (t < tmax && tend < tmax && h * a > hmax) *next_h = hmax;
880:     /* if step size is changed to match a span time point */
881:     if (ts->tspan && h != *next_h && !adapt->dt_span_cached) adapt->dt_span_cached = h;
882:     /* reset time step after a span time point */
883:     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)) {
884:       *next_h               = adapt->dt_span_cached;
885:       adapt->dt_span_cached = 0;
886:     }
887:   }
888:   if (adapt->monitor) {
889:     const char *sc_name = (scheme < ncandidates) ? adapt->candidates.name[scheme] : "";
890:     PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
891:     if (wlte < 0) {
892:       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",
893:                                        (double)ts->ptime, (double)h, (double)*next_h));
894:     } else {
895:       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",
896:                                        (double)ts->ptime, (double)h, (double)*next_h, (double)wlte, (double)wltea, (double)wlter));
897:     }
898:     PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
899:   }
900:   return 0;
901: }

903: /*@
904:    TSAdaptSetTimeStepIncreaseDelay - The number of timesteps to wait after a decrease in the timestep due to failed solver
905:                                      before increasing the time step.

907:    Logicially Collective

909:    Input Parameters:
910: +  adapt - adaptive controller context
911: -  cnt - the number of timesteps

913:    Options Database Key:
914: .  -ts_adapt_time_step_increase_delay cnt - number of steps to delay the increase

916:    Level: advanced

918:    Notes:
919:    This is to prevent an adaptor from bouncing back and forth between two nearby timesteps. The default is 0.

921:    The successful use of this option is problem dependent

923:    Developer Note:
924:    There is no theory to support this option

926: .seealso: [](chapter_ts), `TSAdapt`
927: @*/
928: PetscErrorCode TSAdaptSetTimeStepIncreaseDelay(TSAdapt adapt, PetscInt cnt)
929: {
930:   adapt->timestepjustdecreased_delay = cnt;
931:   return 0;
932: }

934: /*@
935:    TSAdaptCheckStage - checks whether to accept a stage, (e.g. reject and change time step size if nonlinear solve fails or solution vector is infeasible)

937:    Collective

939:    Input Parameters:
940: +  adapt - adaptive controller context
941: .  ts - time stepper
942: .  t - Current simulation time
943: -  Y - Current solution vector

945:    Output Parameter:
946: .  accept - `PETSC_TRUE` to accept the stage, `PETSC_FALSE` to reject

948:    Level: developer

950: .seealso: [](chapter_ts), `TSAdapt`
951: @*/
952: PetscErrorCode TSAdaptCheckStage(TSAdapt adapt, TS ts, PetscReal t, Vec Y, PetscBool *accept)
953: {
954:   SNESConvergedReason snesreason = SNES_CONVERGED_ITERATING;


960:   if (ts->snes) SNESGetConvergedReason(ts->snes, &snesreason);
961:   if (snesreason < 0) {
962:     *accept = PETSC_FALSE;
963:     if (++ts->num_snes_failures >= ts->max_snes_failures && ts->max_snes_failures > 0) {
964:       ts->reason = TS_DIVERGED_NONLINEAR_SOLVE;
965:       PetscInfo(ts, "Step=%" PetscInt_FMT ", nonlinear solve failures %" PetscInt_FMT " greater than current TS allowed, stopping solve\n", ts->steps, ts->num_snes_failures);
966:       if (adapt->monitor) {
967:         PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
968:         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,
969:                                          (double)ts->ptime, (double)ts->time_step, ts->num_snes_failures));
970:         PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
971:       }
972:     }
973:   } else {
974:     *accept = PETSC_TRUE;
975:     TSFunctionDomainError(ts, t, Y, accept);
976:     if (*accept && adapt->checkstage) {
977:       (*adapt->checkstage)(adapt, ts, t, Y, accept);
978:       if (!*accept) {
979:         PetscInfo(ts, "Step=%" PetscInt_FMT ", solution rejected by user function provided by TSSetFunctionDomainError()\n", ts->steps);
980:         if (adapt->monitor) {
981:           PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
982:           PetscViewerASCIIPrintf(adapt->monitor, "    TSAdapt %s step %3" PetscInt_FMT " stage rejected by user function provided by TSSetFunctionDomainError()\n", ((PetscObject)adapt)->type_name, ts->steps);
983:           PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
984:         }
985:       }
986:     }
987:   }

989:   if (!(*accept) && !ts->reason) {
990:     PetscReal dt, new_dt;
991:     TSGetTimeStep(ts, &dt);
992:     new_dt = dt * adapt->scale_solve_failed;
993:     TSSetTimeStep(ts, new_dt);
994:     adapt->timestepjustdecreased += adapt->timestepjustdecreased_delay;
995:     if (adapt->monitor) {
996:       PetscViewerASCIIAddTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
997:       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);
998:       PetscViewerASCIISubtractTab(adapt->monitor, ((PetscObject)adapt)->tablevel);
999:     }
1000:   }
1001:   return 0;
1002: }

1004: /*@
1005:   TSAdaptCreate - create an adaptive controller context for time stepping

1007:   Collective

1009:   Input Parameter:
1010: . comm - The communicator

1012:   Output Parameter:
1013: . adapt - new `TSAdapt` object

1015:   Level: developer

1017:   Note:
1018:   `TSAdapt` creation is handled by `TS`, so users should not need to call this function.

1020: .seealso: [](chapter_ts), `TSAdapt`, `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptDestroy()`
1021: @*/
1022: PetscErrorCode TSAdaptCreate(MPI_Comm comm, TSAdapt *inadapt)
1023: {
1024:   TSAdapt adapt;

1027:   *inadapt = NULL;
1028:   TSAdaptInitializePackage();

1030:   PetscHeaderCreate(adapt, TSADAPT_CLASSID, "TSAdapt", "Time stepping adaptivity", "TS", comm, TSAdaptDestroy, TSAdaptView);

1032:   adapt->always_accept      = PETSC_FALSE;
1033:   adapt->safety             = 0.9;
1034:   adapt->reject_safety      = 0.5;
1035:   adapt->clip[0]            = 0.1;
1036:   adapt->clip[1]            = 10.;
1037:   adapt->dt_min             = 1e-20;
1038:   adapt->dt_max             = 1e+20;
1039:   adapt->ignore_max         = -1.0;
1040:   adapt->glee_use_local     = PETSC_TRUE;
1041:   adapt->scale_solve_failed = 0.25;
1042:   /* these two safety factors are not public, and they are used only in the TS_EXACTFINALTIME_MATCHSTEP case
1043:      to prevent from situations were unreasonably small time steps are taken in order to match the final time */
1044:   adapt->matchstepfac[0]             = 0.01; /* allow 1% step size increase in the last step */
1045:   adapt->matchstepfac[1]             = 2.0;  /* halve last step if it is greater than what remains divided this factor */
1046:   adapt->wnormtype                   = NORM_2;
1047:   adapt->timestepjustdecreased_delay = 0;

1049:   *inadapt = adapt;
1050:   return 0;
1051: }