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