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