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: }