1 2# ----------------------------------------------------------------------------- 3 4class TSType(object): 5 """The time stepping method.""" 6 # native 7 EULER = S_(TSEULER) 8 BEULER = S_(TSBEULER) 9 BASICSYMPLECTIC = S_(TSBASICSYMPLECTIC) 10 PSEUDO = S_(TSPSEUDO) 11 CN = S_(TSCN) 12 SUNDIALS = S_(TSSUNDIALS) 13 RK = S_(TSRK) 14 PYTHON = S_(TSPYTHON) 15 THETA = S_(TSTHETA) 16 ALPHA = S_(TSALPHA) 17 ALPHA2 = S_(TSALPHA2) 18 GLLE = S_(TSGLLE) 19 GLEE = S_(TSGLEE) 20 SSP = S_(TSSSP) 21 ARKIMEX = S_(TSARKIMEX) 22 DIRK = S_(TSDIRK) 23 ROSW = S_(TSROSW) 24 EIMEX = S_(TSEIMEX) 25 MIMEX = S_(TSMIMEX) 26 BDF = S_(TSBDF) 27 RADAU5 = S_(TSRADAU5) 28 MPRK = S_(TSMPRK) 29 DISCGRAD = S_(TSDISCGRAD) 30 # aliases 31 FE = EULER 32 BE = BEULER 33 TH = THETA 34 CRANK_NICOLSON = CN 35 RUNGE_KUTTA = RK 36 37 38class TSRKType(object): 39 """The *RK* subtype.""" 40 RK1FE = S_(TSRK1FE) 41 RK2A = S_(TSRK2A) 42 RK2B = S_(TSRK2B) 43 RK4 = S_(TSRK4) 44 RK3BS = S_(TSRK3BS) 45 RK3 = S_(TSRK3) 46 RK5F = S_(TSRK5F) 47 RK5DP = S_(TSRK5DP) 48 RK5BS = S_(TSRK5BS) 49 RK6VR = S_(TSRK6VR) 50 RK7VR = S_(TSRK7VR) 51 RK8VR = S_(TSRK8VR) 52 53 54class TSARKIMEXType(object): 55 """The *ARKIMEX* subtype.""" 56 ARKIMEX1BEE = S_(TSARKIMEX1BEE) 57 ARKIMEXA2 = S_(TSARKIMEXA2) 58 ARKIMEXL2 = S_(TSARKIMEXL2) 59 ARKIMEXARS122 = S_(TSARKIMEXARS122) 60 ARKIMEX2C = S_(TSARKIMEX2C) 61 ARKIMEX2D = S_(TSARKIMEX2D) 62 ARKIMEX2E = S_(TSARKIMEX2E) 63 ARKIMEXPRSSP2 = S_(TSARKIMEXPRSSP2) 64 ARKIMEX3 = S_(TSARKIMEX3) 65 ARKIMEXBPR3 = S_(TSARKIMEXBPR3) 66 ARKIMEXARS443 = S_(TSARKIMEXARS443) 67 ARKIMEX4 = S_(TSARKIMEX4) 68 ARKIMEX5 = S_(TSARKIMEX5) 69 70 71class TSDIRKType(object): 72 """The *DIRK* subtype.""" 73 DIRKS212 = S_(TSDIRKS212) 74 DIRKES122SAL = S_(TSDIRKES122SAL) 75 DIRKES213SAL = S_(TSDIRKES213SAL) 76 DIRKES324SAL = S_(TSDIRKES324SAL) 77 DIRKES325SAL = S_(TSDIRKES325SAL) 78 DIRK657A = S_(TSDIRK657A) 79 DIRKES648SA = S_(TSDIRKES648SA) 80 DIRK658A = S_(TSDIRK658A) 81 DIRKS659A = S_(TSDIRKS659A) 82 DIRK7510SAL = S_(TSDIRK7510SAL) 83 DIRKES7510SA = S_(TSDIRKES7510SA) 84 DIRK759A = S_(TSDIRK759A) 85 DIRKS7511SAL = S_(TSDIRKS7511SAL) 86 DIRK8614A = S_(TSDIRK8614A) 87 DIRK8616SAL = S_(TSDIRK8616SAL) 88 DIRKES8516SAL = S_(TSDIRKES8516SAL) 89 90 91class TSProblemType(object): 92 """Distinguishes linear and nonlinear problems.""" 93 LINEAR = TS_LINEAR 94 NONLINEAR = TS_NONLINEAR 95 96 97class TSEquationType(object): 98 """Distinguishes among types of explicit and implicit equations.""" 99 UNSPECIFIED = TS_EQ_UNSPECIFIED 100 EXPLICIT = TS_EQ_EXPLICIT 101 ODE_EXPLICIT = TS_EQ_ODE_EXPLICIT 102 DAE_SEMI_EXPLICIT_INDEX1 = TS_EQ_DAE_SEMI_EXPLICIT_INDEX1 103 DAE_SEMI_EXPLICIT_INDEX2 = TS_EQ_DAE_SEMI_EXPLICIT_INDEX2 104 DAE_SEMI_EXPLICIT_INDEX3 = TS_EQ_DAE_SEMI_EXPLICIT_INDEX3 105 DAE_SEMI_EXPLICIT_INDEXHI = TS_EQ_DAE_SEMI_EXPLICIT_INDEXHI 106 IMPLICIT = TS_EQ_IMPLICIT 107 ODE_IMPLICIT = TS_EQ_ODE_IMPLICIT 108 DAE_IMPLICIT_INDEX1 = TS_EQ_DAE_IMPLICIT_INDEX1 109 DAE_IMPLICIT_INDEX2 = TS_EQ_DAE_IMPLICIT_INDEX2 110 DAE_IMPLICIT_INDEX3 = TS_EQ_DAE_IMPLICIT_INDEX3 111 DAE_IMPLICIT_INDEXHI = TS_EQ_DAE_IMPLICIT_INDEXHI 112 113 114class TSExactFinalTime(object): 115 """The method for ending time stepping.""" 116 UNSPECIFIED = TS_EXACTFINALTIME_UNSPECIFIED 117 STEPOVER = TS_EXACTFINALTIME_STEPOVER 118 INTERPOLATE = TS_EXACTFINALTIME_INTERPOLATE 119 MATCHSTEP = TS_EXACTFINALTIME_MATCHSTEP 120 121 122class TSConvergedReason: 123 """The reason the time step is converging.""" 124 # iterating 125 CONVERGED_ITERATING = TS_CONVERGED_ITERATING 126 ITERATING = TS_CONVERGED_ITERATING 127 # converged 128 CONVERGED_TIME = TS_CONVERGED_TIME 129 CONVERGED_ITS = TS_CONVERGED_ITS 130 CONVERGED_USER = TS_CONVERGED_USER 131 CONVERGED_EVENT = TS_CONVERGED_EVENT 132 # diverged 133 DIVERGED_NONLINEAR_SOLVE = TS_DIVERGED_NONLINEAR_SOLVE 134 DIVERGED_STEP_REJECTED = TS_DIVERGED_STEP_REJECTED 135 136# ----------------------------------------------------------------------------- 137 138 139cdef class TS(Object): 140 """ODE integrator. 141 142 TS is described in the `PETSc manual <petsc:manual/ts>`. 143 144 See Also 145 -------- 146 petsc.TS 147 148 """ 149 Type = TSType 150 RKType = TSRKType 151 ARKIMEXType = TSARKIMEXType 152 DIRKType = TSDIRKType 153 ProblemType = TSProblemType 154 EquationType = TSEquationType 155 ExactFinalTime = TSExactFinalTime 156 ExactFinalTimeOption = TSExactFinalTime 157 ConvergedReason = TSConvergedReason 158 159 # --- xxx --- 160 161 def __cinit__(self): 162 self.obj = <PetscObject*> &self.ts 163 self.ts = NULL 164 165 # --- xxx --- 166 167 def view(self, Viewer viewer=None) -> None: 168 """Print the `TS` object. 169 170 Collective. 171 172 Parameters 173 ---------- 174 viewer 175 The visualization context. 176 177 Notes 178 ----- 179 ``-ts_view`` calls TSView at the end of TSStep 180 181 See Also 182 -------- 183 petsc.TSView 184 185 """ 186 cdef PetscViewer cviewer = NULL 187 if viewer is not None: cviewer = viewer.vwr 188 CHKERR(TSView(self.ts, cviewer)) 189 190 def load(self, Viewer viewer) -> None: 191 """Load a `TS` that has been stored in binary with `view`. 192 193 Collective. 194 195 Parameters 196 ---------- 197 viewer 198 The visualization context. 199 200 See Also 201 -------- 202 petsc.TSLoad 203 204 """ 205 CHKERR(TSLoad(self.ts, viewer.vwr)) 206 207 def destroy(self) -> Self: 208 """Destroy the `TS` that was created with `create`. 209 210 Collective. 211 212 See Also 213 -------- 214 petsc.TSDestroy 215 216 """ 217 CHKERR(TSDestroy(&self.ts)) 218 return self 219 220 def create(self, comm: Comm | None = None) -> Self: 221 """Create an empty `TS`. 222 223 Collective. 224 225 The problem type can then be set with `setProblemType` and the type of 226 solver can then be set with `setType`. 227 228 Parameters 229 ---------- 230 comm 231 MPI communicator, defaults to `Sys.getDefaultComm`. 232 233 See Also 234 -------- 235 petsc.TSCreate 236 237 """ 238 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 239 cdef PetscTS newts = NULL 240 CHKERR(TSCreate(ccomm, &newts)) 241 CHKERR(PetscCLEAR(self.obj)); self.ts = newts 242 return self 243 244 def clone(self) -> TS: 245 """Return a shallow clone of the `TS` object. 246 247 Collective. 248 249 See Also 250 -------- 251 petsc.TSClone 252 253 """ 254 cdef TS ts = TS() 255 CHKERR(TSClone(self.ts, &ts.ts)) 256 return ts 257 258 def setType(self, ts_type: Type | str) -> None: 259 """Set the method to be used as the `TS` solver. 260 261 Collective. 262 263 Parameters 264 ---------- 265 ts_type 266 The solver type. 267 268 Notes 269 ----- 270 ``-ts_type`` sets the method from the commandline 271 272 See Also 273 -------- 274 petsc.TSSetType 275 276 """ 277 cdef PetscTSType cval = NULL 278 ts_type = str2bytes(ts_type, &cval) 279 CHKERR(TSSetType(self.ts, cval)) 280 281 def setRKType(self, ts_type: RKType | str) -> None: 282 """Set the type of the *Runge-Kutta* scheme. 283 284 Logically collective. 285 286 Parameters 287 ---------- 288 ts_type 289 The type of scheme. 290 291 Notes 292 ----- 293 ``-ts_rk_type`` sets scheme type from the commandline. 294 295 See Also 296 -------- 297 petsc.TSRKSetType 298 299 """ 300 cdef PetscTSRKType cval = NULL 301 ts_type = str2bytes(ts_type, &cval) 302 CHKERR(TSRKSetType(self.ts, cval)) 303 304 def setARKIMEXType(self, ts_type: ARKIMEXType | str) -> None: 305 """Set the type of `Type.ARKIMEX` scheme. 306 307 Logically collective. 308 309 Parameters 310 ---------- 311 ts_type 312 The type of `Type.ARKIMEX` scheme. 313 314 Notes 315 ----- 316 ``-ts_arkimex_type`` sets scheme type from the commandline. 317 318 See Also 319 -------- 320 petsc.TSARKIMEXSetType 321 322 """ 323 cdef PetscTSARKIMEXType cval = NULL 324 ts_type = str2bytes(ts_type, &cval) 325 CHKERR(TSARKIMEXSetType(self.ts, cval)) 326 327 def setARKIMEXFullyImplicit(self, flag: bool) -> None: 328 """Solve both parts of the equation implicitly. 329 330 Logically collective. 331 332 Parameters 333 ---------- 334 flag 335 Set to True for fully implicit. 336 337 See Also 338 -------- 339 petsc.TSARKIMEXSetFullyImplicit 340 341 """ 342 cdef PetscBool bval = asBool(flag) 343 CHKERR(TSARKIMEXSetFullyImplicit(self.ts, bval)) 344 345 def setARKIMEXFastSlowSplit(self, flag: bool) -> None: 346 """Use ARKIMEX for solving a fast-slow system. 347 348 Logically collective. 349 350 Parameters 351 ---------- 352 flag 353 Set to True for fast-slow partitioned systems. 354 355 See Also 356 -------- 357 petsc.TSARKIMEXSetType 358 """ 359 cdef PetscBool bval = asBool(flag) 360 CHKERR(TSARKIMEXSetFastSlowSplit(self.ts, bval)) 361 362 def getType(self) -> str: 363 """Return the `TS` type. 364 365 Not collective. 366 367 See Also 368 -------- 369 petsc.TSGetType 370 371 """ 372 cdef PetscTSType cval = NULL 373 CHKERR(TSGetType(self.ts, &cval)) 374 return bytes2str(cval) 375 376 def getRKType(self) -> str: 377 """Return the `Type.RK` scheme. 378 379 Not collective. 380 381 See Also 382 -------- 383 petsc.TSRKGetType 384 385 """ 386 cdef PetscTSRKType cval = NULL 387 CHKERR(TSRKGetType(self.ts, &cval)) 388 return bytes2str(cval) 389 390 def getARKIMEXType(self) -> str: 391 """Return the `Type.ARKIMEX` scheme. 392 393 Not collective. 394 395 See Also 396 -------- 397 petsc.TSARKIMEXGetType 398 399 """ 400 cdef PetscTSARKIMEXType cval = NULL 401 CHKERR(TSARKIMEXGetType(self.ts, &cval)) 402 return bytes2str(cval) 403 404 def setDIRKType(self, ts_type: DIRKType | str) -> None: 405 """Set the type of `Type.DIRK` scheme. 406 407 Logically collective. 408 409 Parameters 410 ---------- 411 ts_type 412 The type of `Type.DIRK` scheme. 413 414 Notes 415 ----- 416 ``-ts_dirk_type`` sets scheme type from the commandline. 417 418 See Also 419 -------- 420 petsc.TSDIRKSetType 421 422 """ 423 cdef PetscTSDIRKType cval = NULL 424 ts_type = str2bytes(ts_type, &cval) 425 CHKERR(TSDIRKSetType(self.ts, cval)) 426 427 def getDIRKType(self) -> str: 428 """Return the `Type.DIRK` scheme. 429 430 Not collective. 431 432 See Also 433 -------- 434 setDIRKType, petsc.TSDIRKGetType 435 436 """ 437 cdef PetscTSDIRKType cval = NULL 438 CHKERR(TSDIRKGetType(self.ts, &cval)) 439 return bytes2str(cval) 440 441 def setProblemType(self, ptype: ProblemType) -> None: 442 """Set the type of problem to be solved. 443 444 Logically collective. 445 446 Parameters 447 ---------- 448 ptype 449 The type of problem of the forms. 450 451 See Also 452 -------- 453 petsc.TSSetProblemType 454 455 """ 456 CHKERR(TSSetProblemType(self.ts, ptype)) 457 458 def getProblemType(self) -> ProblemType: 459 """Return the type of problem to be solved. 460 461 Not collective. 462 463 See Also 464 -------- 465 petsc.TSGetProblemType 466 467 """ 468 cdef PetscTSProblemType ptype = TS_NONLINEAR 469 CHKERR(TSGetProblemType(self.ts, &ptype)) 470 return ptype 471 472 def setEquationType(self, eqtype: EquationType) -> None: 473 """Set the type of the equation that `TS` is solving. 474 475 Logically collective. 476 477 Parameters 478 ---------- 479 eqtype 480 The type of equation. 481 482 See Also 483 -------- 484 petsc.TSSetEquationType 485 486 """ 487 CHKERR(TSSetEquationType(self.ts, eqtype)) 488 489 def getEquationType(self) -> EquationType: 490 """Get the type of the equation that `TS` is solving. 491 492 Not collective. 493 494 See Also 495 -------- 496 petsc.TSGetEquationType 497 498 """ 499 cdef PetscTSEquationType eqtype = TS_EQ_UNSPECIFIED 500 CHKERR(TSGetEquationType(self.ts, &eqtype)) 501 return eqtype 502 503 def setOptionsPrefix(self, prefix : str | None) -> None: 504 """Set the prefix used for all the `TS` options. 505 506 Logically collective. 507 508 Parameters 509 ---------- 510 prefix 511 The prefix to prepend to all option names. 512 513 Notes 514 ----- 515 A hyphen must not be given at the beginning of the prefix name. 516 517 See Also 518 -------- 519 petsc_options, petsc.TSSetOptionsPrefix 520 521 """ 522 cdef const char *cval = NULL 523 prefix = str2bytes(prefix, &cval) 524 CHKERR(TSSetOptionsPrefix(self.ts, cval)) 525 526 def getOptionsPrefix(self) -> str: 527 """Return the prefix used for all the `TS` options. 528 529 Not collective. 530 531 See Also 532 -------- 533 petsc.TSGetOptionsPrefix 534 535 """ 536 cdef const char *cval = NULL 537 CHKERR(TSGetOptionsPrefix(self.ts, &cval)) 538 return bytes2str(cval) 539 540 def appendOptionsPrefix(self, prefix: str | None) -> None: 541 """Append to the prefix used for all the `TS` options. 542 543 Logically collective. 544 545 Parameters 546 ---------- 547 prefix 548 The prefix to append to the current prefix. 549 550 Notes 551 ----- 552 A hyphen must not be given at the beginning of the prefix name. 553 554 See Also 555 -------- 556 petsc_options, petsc.TSAppendOptionsPrefix 557 558 """ 559 cdef const char *cval = NULL 560 prefix = str2bytes(prefix, &cval) 561 CHKERR(TSAppendOptionsPrefix(self.ts, cval)) 562 563 def setFromOptions(self) -> None: 564 """Set various `TS` parameters from user options. 565 566 Collective. 567 568 See Also 569 -------- 570 petsc_options, petsc.TSSetFromOptions 571 572 """ 573 CHKERR(TSSetFromOptions(self.ts)) 574 575 # --- application context --- 576 577 def setAppCtx(self, appctx: Any) -> None: 578 """Set the application context. 579 580 Not collective. 581 582 Parameters 583 ---------- 584 appctx 585 The application context. 586 587 """ 588 self.set_attr('__appctx__', appctx) 589 590 def getAppCtx(self) -> Any: 591 """Return the application context.""" 592 return self.get_attr('__appctx__') 593 594 # --- user RHS Function/Jacobian routines --- 595 596 def setRHSFunction( 597 self, 598 function: TSRHSFunction | None, 599 Vec f=None, 600 args : tuple[Any, ...] | None = None, 601 kargs : dict[str, Any] | None = None) -> None: 602 """Set the routine for evaluating the function ``G`` in ``U_t = G(t, u)``. 603 604 Logically collective. 605 606 Parameters 607 ---------- 608 function 609 The right-hand side function. 610 f 611 The vector into which the right-hand side is computed. 612 args 613 Additional positional arguments for ``function``. 614 kargs 615 Additional keyword arguments for ``function``. 616 617 See Also 618 -------- 619 petsc.TSSetRHSFunction 620 621 """ 622 cdef PetscVec fvec=NULL 623 if f is not None: fvec = f.vec 624 if function is not None: 625 if args is None: args = () 626 if kargs is None: kargs = {} 627 context = (function, args, kargs) 628 self.set_attr('__rhsfunction__', context) 629 CHKERR(TSSetRHSFunction(self.ts, fvec, TS_RHSFunction, <void*>context)) 630 else: 631 CHKERR(TSSetRHSFunction(self.ts, fvec, NULL, NULL)) 632 633 def setRHSJacobian( 634 self, 635 jacobian: TSRHSJacobian | None, 636 Mat J=None, 637 Mat P=None, 638 args : tuple[Any, ...] | None = None, 639 kargs : dict[str, Any] | None = None) -> None: 640 """Set the function to compute the Jacobian of ``G`` in ``U_t = G(U, t)``. 641 642 Logically collective. 643 644 Parameters 645 ---------- 646 jacobian 647 The right-hand side function. 648 J 649 The matrix into which the jacobian is computed. 650 P 651 The matrix into which the preconditioner is computed. 652 args 653 Additional positional arguments for ``jacobian``. 654 kargs 655 Additional keyword arguments for ``jacobian``. 656 657 See Also 658 -------- 659 petsc.TSSetRHSJacobian 660 661 """ 662 cdef PetscMat Jmat=NULL 663 if J is not None: Jmat = J.mat 664 cdef PetscMat Pmat=Jmat 665 if P is not None: Pmat = P.mat 666 if jacobian is not None: 667 if args is None: args = () 668 if kargs is None: kargs = {} 669 context = (jacobian, args, kargs) 670 self.set_attr('__rhsjacobian__', context) 671 CHKERR(TSSetRHSJacobian(self.ts, Jmat, Pmat, TS_RHSJacobian, <void*>context)) 672 else: 673 CHKERR(TSSetRHSJacobian(self.ts, Jmat, Pmat, NULL, NULL)) 674 675 def computeRHSFunction(self, t: float, Vec x, Vec f) -> None: 676 """Evaluate the right-hand side function. 677 678 Collective. 679 680 Parameters 681 ---------- 682 t 683 The time at which to evaluate the RHS. 684 x 685 The state vector. 686 f 687 The Vec into which the RHS is computed. 688 689 See Also 690 -------- 691 petsc.TSComputeRHSFunction 692 693 """ 694 cdef PetscReal time = asReal(t) 695 CHKERR(TSComputeRHSFunction(self.ts, time, x.vec, f.vec)) 696 697 def computeRHSFunctionLinear(self, t: float, Vec x, Vec f) -> None: 698 """Evaluate the right-hand side via the user-provided Jacobian. 699 700 Collective. 701 702 Parameters 703 ---------- 704 t 705 The time at which to evaluate the RHS. 706 x 707 The state vector. 708 f 709 The Vec into which the RHS is computed. 710 711 See Also 712 -------- 713 petsc.TSComputeRHSFunctionLinear 714 715 """ 716 cdef PetscReal time = asReal(t) 717 CHKERR(TSComputeRHSFunctionLinear(self.ts, time, x.vec, f.vec, NULL)) 718 719 def computeRHSJacobian(self, t: float, Vec x, Mat J, Mat P=None) -> None: 720 """Compute the Jacobian matrix that has been set with `setRHSJacobian`. 721 722 Collective. 723 724 Parameters 725 ---------- 726 t 727 The time at which to evaluate the Jacobian. 728 x 729 The state vector. 730 J 731 The matrix into which the Jacobian is computed. 732 P 733 The optional matrix to use for building a preconditioner. 734 735 See Also 736 -------- 737 petsc.TSComputeRHSJacobian 738 739 """ 740 cdef PetscReal time = asReal(t) 741 cdef PetscMat jmat = J.mat, pmat = J.mat 742 if P is not None: pmat = P.mat 743 CHKERR(TSComputeRHSJacobian(self.ts, time, x.vec, jmat, pmat)) 744 745 def computeRHSJacobianConstant(self, t: float, Vec x, Mat J, Mat P=None) -> None: 746 """Reuse a Jacobian that is time-independent. 747 748 Collective. 749 750 Parameters 751 ---------- 752 t 753 The time at which to evaluate the Jacobian. 754 x 755 The state vector. 756 J 757 A pointer to the stored Jacobian. 758 P 759 An optional pointer to the matrix used to construct the preconditioner. 760 761 See Also 762 -------- 763 petsc.TSComputeRHSJacobianConstant 764 765 """ 766 cdef PetscReal time = asReal(t) 767 cdef PetscMat jmat = J.mat, pmat = J.mat 768 if P is not None: pmat = P.mat 769 CHKERR(TSComputeRHSJacobianConstant(self.ts, time, x.vec, jmat, pmat, NULL)) 770 771 def getRHSFunction(self) -> tuple[Vec, TSRHSFunction]: 772 """Return the vector where the rhs is stored and the function used to compute it. 773 774 Not collective. 775 776 See Also 777 -------- 778 petsc.TSGetRHSFunction 779 780 """ 781 cdef Vec f = Vec() 782 CHKERR(TSGetRHSFunction(self.ts, &f.vec, NULL, NULL)) 783 CHKERR(PetscINCREF(f.obj)) 784 cdef object function = self.get_attr('__rhsfunction__') 785 return (f, function) 786 787 def getRHSJacobian(self) -> tuple[Mat, Mat, TSRHSJacobian]: 788 """Return the Jacobian and the function used to compute them. 789 790 Not collective. 791 792 See Also 793 -------- 794 petsc.TSGetRHSJacobian 795 796 """ 797 cdef Mat J = Mat(), P = Mat() 798 CHKERR(TSGetRHSJacobian(self.ts, &J.mat, &P.mat, NULL, NULL)) 799 CHKERR(PetscINCREF(J.obj)); CHKERR(PetscINCREF(P.obj)) 800 cdef object jacobian = self.get_attr('__rhsjacobian__') 801 return (J, P, jacobian) 802 803 # --- user Implicit Function/Jacobian routines --- 804 805 def setIFunction( 806 self, 807 function: TSIFunction | None, 808 Vec f=None, 809 args : tuple[Any, ...] | None = None, 810 kargs : dict[str, Any] | None = None) -> None: 811 """Set the function representing the DAE to be solved. 812 813 Logically collective. 814 815 Parameters 816 ---------- 817 function 818 The right-hand side function. 819 f 820 The vector to store values or `None` to be created internally. 821 args 822 Additional positional arguments for ``function``. 823 kargs 824 Additional keyword arguments for ``function``. 825 826 See Also 827 -------- 828 petsc.TSSetIFunction 829 830 """ 831 cdef PetscVec fvec=NULL 832 if f is not None: fvec = f.vec 833 if function is not None: 834 if args is None: args = () 835 if kargs is None: kargs = {} 836 context = (function, args, kargs) 837 self.set_attr('__ifunction__', context) 838 CHKERR(TSSetIFunction(self.ts, fvec, TS_IFunction, <void*>context)) 839 else: 840 CHKERR(TSSetIFunction(self.ts, fvec, NULL, NULL)) 841 842 def setIJacobian( 843 self, 844 jacobian: TSIJacobian | None, 845 Mat J=None, 846 Mat P=None, 847 args : tuple[Any, ...] | None = None, 848 kargs : dict[str, Any] | None = None) -> None: 849 """Set the function to compute the Jacobian. 850 851 Logically collective. 852 853 Set the function to compute the matrix ``dF/dU + a*dF/dU_t`` where 854 ``F(t, U, U_t)`` is the function provided with `setIFunction`. 855 856 Parameters 857 ---------- 858 jacobian 859 The function which computes the Jacobian. 860 J 861 The matrix into which the Jacobian is computed. 862 P 863 The optional matrix to use for building a preconditioner matrix. 864 args 865 Additional positional arguments for ``jacobian``. 866 kargs 867 Additional keyword arguments for ``jacobian``. 868 869 See Also 870 -------- 871 petsc.TSSetIJacobian 872 873 """ 874 cdef PetscMat Jmat=NULL 875 if J is not None: Jmat = J.mat 876 cdef PetscMat Pmat=Jmat 877 if P is not None: Pmat = P.mat 878 if jacobian is not None: 879 if args is None: args = () 880 if kargs is None: kargs = {} 881 context = (jacobian, args, kargs) 882 self.set_attr('__ijacobian__', context) 883 CHKERR(TSSetIJacobian(self.ts, Jmat, Pmat, TS_IJacobian, <void*>context)) 884 else: 885 CHKERR(TSSetIJacobian(self.ts, Jmat, Pmat, NULL, NULL)) 886 887 def setIJacobianP( 888 self, 889 jacobian, 890 Mat J=None, 891 args : tuple[Any, ...] | None = None, 892 kargs : dict[str, Any] | None = None) -> None: 893 """Set the function that computes the Jacobian. 894 895 Logically collective. 896 897 Set the function that computes the Jacobian of ``F`` with respect to 898 the parameters ``P`` where ``F(Udot, U, t) = G(U, P, t)``, as well as the 899 location to store the matrix. 900 901 Parameters 902 ---------- 903 jacobian 904 The function which computes the Jacobian. 905 J 906 The matrix into which the Jacobian is computed. 907 args 908 Additional positional arguments for ``jacobian``. 909 kargs 910 Additional keyword arguments for ``jacobian``. 911 912 See Also 913 -------- 914 petsc.TSSetIJacobianP 915 916 """ 917 cdef PetscMat Jmat=NULL 918 if J is not None: Jmat = J.mat 919 if jacobian is not None: 920 if args is None: args = () 921 if kargs is None: kargs = {} 922 context = (jacobian, args, kargs) 923 self.set_attr('__ijacobianp__', context) 924 CHKERR(TSSetIJacobianP(self.ts, Jmat, TS_IJacobianP, <void*>context)) 925 else: 926 CHKERR(TSSetIJacobianP(self.ts, Jmat, NULL, NULL)) 927 928 def computeIFunction(self, 929 t: float, Vec x, Vec xdot, 930 Vec f, imex: bool = False) -> None: 931 """Evaluate the DAE residual written in implicit form. 932 933 Collective. 934 935 Parameters 936 ---------- 937 t 938 The current time. 939 x 940 The state vector. 941 xdot 942 The time derivative of the state vector. 943 f 944 The vector into which the residual is stored. 945 imex 946 A flag which indicates if the RHS should be kept separate. 947 948 See Also 949 -------- 950 petsc.TSComputeIFunction 951 952 """ 953 cdef PetscReal rval = asReal(t) 954 cdef PetscBool bval = imex 955 CHKERR(TSComputeIFunction(self.ts, rval, x.vec, xdot.vec, 956 f.vec, bval)) 957 958 def computeIJacobian(self, 959 t: float, Vec x, Vec xdot, a: float, 960 Mat J, Mat P = None, imex: bool = False) -> None: 961 """Evaluate the Jacobian of the DAE. 962 963 Collective. 964 965 If ``F(t, U, Udot)=0`` is the DAE, the required Jacobian is 966 ``dF/dU + shift*dF/dUdot`` 967 968 Parameters 969 ---------- 970 t 971 The current time. 972 x 973 The state vector. 974 xdot 975 The time derivative of the state vector. 976 a 977 The shift to apply 978 J 979 The matrix into which the Jacobian is computed. 980 P 981 The optional matrix to use for building a preconditioner. 982 imex 983 A flag which indicates if the RHS should be kept separate. 984 985 See Also 986 -------- 987 petsc.TSComputeIJacobian 988 989 """ 990 cdef PetscReal rval1 = asReal(t) 991 cdef PetscReal rval2 = asReal(a) 992 cdef PetscBool bval = imex 993 cdef PetscMat jmat = J.mat, pmat = J.mat 994 if P is not None: pmat = P.mat 995 CHKERR(TSComputeIJacobian(self.ts, rval1, x.vec, xdot.vec, rval2, 996 jmat, pmat, bval)) 997 998 def computeIJacobianP(self, 999 t: float, Vec x, Vec xdot, a: float, 1000 Mat J, imex: bool = False) -> None: 1001 """Evaluate the Jacobian with respect to parameters. 1002 1003 Collective. 1004 1005 Parameters 1006 ---------- 1007 t 1008 The current time. 1009 x 1010 The state vector. 1011 xdot 1012 The time derivative of the state vector. 1013 a 1014 The shift to apply 1015 J 1016 The matrix into which the Jacobian is computed. 1017 imex 1018 A flag which indicates if the RHS should be kept separate. 1019 1020 See Also 1021 -------- 1022 petsc.TSComputeIJacobianP 1023 1024 """ 1025 cdef PetscReal rval1 = asReal(t) 1026 cdef PetscReal rval2 = asReal(a) 1027 cdef PetscBool bval = asBool(imex) 1028 cdef PetscMat jmat = J.mat 1029 CHKERR(TSComputeIJacobianP(self.ts, rval1, x.vec, xdot.vec, rval2, 1030 jmat, bval)) 1031 1032 def getIFunction(self) -> tuple[Vec, TSIFunction]: 1033 """Return the vector and function which computes the implicit residual. 1034 1035 Not collective. 1036 1037 See Also 1038 -------- 1039 petsc.TSGetIFunction 1040 1041 """ 1042 cdef Vec f = Vec() 1043 CHKERR(TSGetIFunction(self.ts, &f.vec, NULL, NULL)) 1044 CHKERR(PetscINCREF(f.obj)) 1045 cdef object function = self.get_attr('__ifunction__') 1046 return (f, function) 1047 1048 def getIJacobian(self) -> tuple[Mat, Mat, TSIJacobian]: 1049 """Return the matrices and function which computes the implicit Jacobian. 1050 1051 Not collective. 1052 1053 See Also 1054 -------- 1055 petsc.TSGetIJacobian 1056 1057 """ 1058 cdef Mat J = Mat(), P = Mat() 1059 CHKERR(TSGetIJacobian(self.ts, &J.mat, &P.mat, NULL, NULL)) 1060 CHKERR(PetscINCREF(J.obj)); CHKERR(PetscINCREF(P.obj)) 1061 cdef object jacobian = self.get_attr('__ijacobian__') 1062 return (J, P, jacobian) 1063 1064 def setI2Function( 1065 self, 1066 function: TSI2Function | None, 1067 Vec f=None, 1068 args : tuple[Any, ...] | None = None, 1069 kargs : dict[str, Any] | None = None) -> None: 1070 """Set the function to compute the 2nd order DAE. 1071 1072 Logically collective. 1073 1074 Parameters 1075 ---------- 1076 function 1077 The right-hand side function. 1078 f 1079 The vector to store values or `None` to be created internally. 1080 args 1081 Additional positional arguments for ``function``. 1082 kargs 1083 Additional keyword arguments for ``function``. 1084 1085 See Also 1086 -------- 1087 petsc.TSSetI2Function 1088 1089 """ 1090 cdef PetscVec fvec=NULL 1091 if f is not None: fvec = f.vec 1092 if function is not None: 1093 if args is None: args = () 1094 if kargs is None: kargs = {} 1095 context = (function, args, kargs) 1096 self.set_attr('__i2function__', context) 1097 CHKERR(TSSetI2Function(self.ts, fvec, TS_I2Function, <void*>context)) 1098 else: 1099 CHKERR(TSSetI2Function(self.ts, fvec, NULL, NULL)) 1100 1101 def setI2Jacobian( 1102 self, 1103 jacobian: TSI2Jacobian | None, 1104 Mat J=None, 1105 Mat P=None, 1106 args=None, 1107 kargs=None) -> None: 1108 """Set the function to compute the Jacobian of the 2nd order DAE. 1109 1110 Logically collective. 1111 1112 Parameters 1113 ---------- 1114 jacobian 1115 The function which computes the Jacobian. 1116 J 1117 The matrix into which the Jacobian is computed. 1118 P 1119 The optional matrix to use for building a preconditioner. 1120 args 1121 Additional positional arguments for ``jacobian``. 1122 kargs 1123 Additional keyword arguments for ``jacobian``. 1124 1125 See Also 1126 -------- 1127 petsc.TSSetI2Jacobian 1128 1129 """ 1130 cdef PetscMat Jmat=NULL 1131 if J is not None: Jmat = J.mat 1132 cdef PetscMat Pmat=Jmat 1133 if P is not None: Pmat = P.mat 1134 if jacobian is not None: 1135 if args is None: args = () 1136 if kargs is None: kargs = {} 1137 context = (jacobian, args, kargs) 1138 self.set_attr('__i2jacobian__', context) 1139 CHKERR(TSSetI2Jacobian(self.ts, Jmat, Pmat, TS_I2Jacobian, <void*>context)) 1140 else: 1141 CHKERR(TSSetI2Jacobian(self.ts, Jmat, Pmat, NULL, NULL)) 1142 1143 def computeI2Function(self, t: float, Vec x, Vec xdot, Vec xdotdot, Vec f) -> None: 1144 """Evaluate the DAE residual in implicit form. 1145 1146 Collective. 1147 1148 Parameters 1149 ---------- 1150 t 1151 The current time. 1152 x 1153 The state vector. 1154 xdot 1155 The time derivative of the state vector. 1156 xdotdot 1157 The second time derivative of the state vector. 1158 f 1159 The vector into which the residual is stored. 1160 1161 See Also 1162 -------- 1163 petsc.TSComputeI2Function 1164 1165 """ 1166 cdef PetscReal rval = asReal(t) 1167 CHKERR(TSComputeI2Function(self.ts, rval, x.vec, xdot.vec, xdotdot.vec, 1168 f.vec)) 1169 1170 def computeI2Jacobian( 1171 self, 1172 t: float, 1173 Vec x, 1174 Vec xdot, 1175 Vec xdotdot, 1176 v: float, 1177 a: float, 1178 Mat J, 1179 Mat P=None) -> None: 1180 """Evaluate the Jacobian of the DAE. 1181 1182 Collective. 1183 1184 If ``F(t, U, V, A)=0`` is the DAE, 1185 the required Jacobian is ``dF/dU + v dF/dV + a dF/dA``. 1186 1187 Parameters 1188 ---------- 1189 t 1190 The current time. 1191 x 1192 The state vector. 1193 xdot 1194 The time derivative of the state vector. 1195 xdotdot 1196 The second time derivative of the state vector. 1197 v 1198 The shift to apply to the first derivative. 1199 a 1200 The shift to apply to the second derivative. 1201 J 1202 The matrix into which the Jacobian is computed. 1203 P 1204 The optional matrix to use for building a preconditioner. 1205 1206 See Also 1207 -------- 1208 petsc.TSComputeI2Jacobian 1209 1210 """ 1211 cdef PetscReal rval1 = asReal(t) 1212 cdef PetscReal rval2 = asReal(v) 1213 cdef PetscReal rval3 = asReal(a) 1214 cdef PetscMat jmat = J.mat, pmat = J.mat 1215 if P is not None: pmat = P.mat 1216 CHKERR(TSComputeI2Jacobian(self.ts, rval1, x.vec, xdot.vec, xdotdot.vec, rval2, rval3, 1217 jmat, pmat)) 1218 1219 def getI2Function(self) -> tuple[Vec, TSI2Function]: 1220 """Return the vector and function which computes the residual. 1221 1222 Not collective. 1223 1224 See Also 1225 -------- 1226 petsc.TSGetI2Function 1227 1228 """ 1229 cdef Vec f = Vec() 1230 CHKERR(TSGetI2Function(self.ts, &f.vec, NULL, NULL)) 1231 CHKERR(PetscINCREF(f.obj)) 1232 cdef object function = self.get_attr('__i2function__') 1233 return (f, function) 1234 1235 def getI2Jacobian(self) -> tuple[Mat, Mat, TSI2Jacobian]: 1236 """Return the matrices and function which computes the Jacobian. 1237 1238 Not collective. 1239 1240 See Also 1241 -------- 1242 petsc.TSGetI2Jacobian 1243 1244 """ 1245 cdef Mat J = Mat(), P = Mat() 1246 CHKERR(TSGetI2Jacobian(self.ts, &J.mat, &P.mat, NULL, NULL)) 1247 CHKERR(PetscINCREF(J.obj)); CHKERR(PetscINCREF(P.obj)) 1248 cdef object jacobian = self.get_attr('__i2jacobian__') 1249 return (J, P, jacobian) 1250 1251 # --- TSRHSSplit routines to support multirate and IMEX solvers --- 1252 def setRHSSplitIS(self, splitname: str, IS iss) -> None: 1253 """Set the index set for the specified split. 1254 1255 Logically collective. 1256 1257 Parameters 1258 ---------- 1259 splitname 1260 Name of this split, if `None` the number of the split is used. 1261 iss 1262 The index set for part of the solution vector. 1263 1264 See Also 1265 -------- 1266 petsc.TSRHSSplitSetIS 1267 1268 """ 1269 cdef const char *cname = NULL 1270 splitname = str2bytes(splitname, &cname) 1271 CHKERR(TSRHSSplitSetIS(self.ts, cname, iss.iset)) 1272 1273 def setRHSSplitRHSFunction( 1274 self, 1275 splitname: str, 1276 function: TSRHSFunction, 1277 Vec r=None, 1278 args : tuple[Any, ...] | None = None, 1279 kargs : dict[str, Any] | None = None) -> None: 1280 """Set the split right-hand-side functions. 1281 1282 Logically collective. 1283 1284 Parameters 1285 ---------- 1286 splitname 1287 Name of the split. 1288 function 1289 The RHS function evaluation routine. 1290 r 1291 Vector to hold the residual. 1292 args 1293 Additional positional arguments for ``function``. 1294 kargs 1295 Additional keyword arguments for ``function``. 1296 1297 See Also 1298 ------- 1299 petsc.TSRHSSplitSetRHSFunction 1300 1301 """ 1302 cdef const char *cname = NULL 1303 cdef char *aname = NULL 1304 splitname = str2bytes(splitname, <const char**>&cname) 1305 cdef PetscVec rvec=NULL 1306 if r is not None: rvec = r.vec 1307 if function is not None: 1308 if args is None: args = () 1309 if kargs is None: kargs = {} 1310 context = (function, args, kargs) 1311 str2bytes(function.__name__, <const char**>&aname) 1312 self.set_attr(aname, context) # to avoid being GCed 1313 CHKERR(TSRHSSplitSetRHSFunction(self.ts, cname, rvec, TS_RHSFunction, <void*>context)) 1314 else: 1315 CHKERR(TSRHSSplitSetRHSFunction(self.ts, cname, rvec, NULL, NULL)) 1316 1317 def setRHSSplitIFunction( 1318 self, 1319 splitname: str, 1320 function: TSIFunction, 1321 Vec r=None, 1322 args : tuple[Any, ...] | None = None, 1323 kargs : dict[str, Any] | None = None) -> None: 1324 """Set the split implicit functions. 1325 1326 Logically collective. 1327 1328 Parameters 1329 ---------- 1330 splitname 1331 Name of the split. 1332 function 1333 The implicit function evaluation routine. 1334 r 1335 Vector to hold the residual. 1336 args 1337 Additional positional arguments for ``function``. 1338 kargs 1339 Additional keyword arguments for ``function``. 1340 1341 See Also 1342 ------- 1343 petsc.TSRHSSplitSetIFunction 1344 1345 """ 1346 cdef const char *cname = NULL 1347 cdef char *aname = NULL 1348 splitname = str2bytes(splitname, &cname) 1349 cdef PetscVec rvec=NULL 1350 if r is not None: rvec = r.vec 1351 if function is not None: 1352 if args is None: args = () 1353 if kargs is None: kargs = {} 1354 context = (function, args, kargs) 1355 str2bytes(function.__name__, <const char**>&aname) 1356 self.set_attr(aname, context) # to avoid being GCed 1357 CHKERR(TSRHSSplitSetIFunction(self.ts, cname, rvec, TS_IFunction, <void*>context)) 1358 else: 1359 CHKERR(TSRHSSplitSetIFunction(self.ts, cname, rvec, NULL, NULL)) 1360 1361 def setRHSSplitIJacobian( 1362 self, 1363 splitname: str, 1364 jacobian: TSRHSJacobian, 1365 Mat J=None, 1366 Mat P=None, 1367 args : tuple[Any, ...] | None = None, 1368 kargs : dict[str, Any] | None = None) -> None: 1369 """Set the Jacobian for the split implicit function. 1370 1371 Logically collective. 1372 1373 Parameters 1374 ---------- 1375 splitname 1376 Name of the split. 1377 jacobian 1378 The Jacobian evaluation routine. 1379 J 1380 Matrix to store Jacobian entries computed by ``jacobian``. 1381 P 1382 Matrix used to compute preconditioner (usually the same as ``J``). 1383 args 1384 Additional positional arguments for ``jacobian``. 1385 kargs 1386 Additional keyword arguments for ``jacobian``. 1387 1388 See Also 1389 ------- 1390 petsc.TSRHSSplitSetIJacobian 1391 1392 """ 1393 cdef const char *cname = NULL 1394 cdef char *aname = NULL 1395 splitname = str2bytes(splitname, &cname) 1396 cdef PetscMat Jmat=NULL 1397 if J is not None: Jmat = J.mat 1398 cdef PetscMat Pmat=Jmat 1399 if P is not None: Pmat = P.mat 1400 if jacobian is not None: 1401 if args is None: args = () 1402 if kargs is None: kargs = {} 1403 context = (jacobian, args, kargs) 1404 str2bytes(jacobian.__name__, <const char**>&aname) 1405 self.set_attr(aname, context) # to avoid being GCed 1406 CHKERR(TSRHSSplitSetIJacobian(self.ts, cname, Jmat, Pmat, TS_IJacobian, <void*>context)) 1407 else: 1408 CHKERR(TSRHSSplitSetIJacobian(self.ts, cname, Jmat, Pmat, NULL, NULL)) 1409 1410 # --- solution vector --- 1411 1412 def setSolution(self, Vec u) -> None: 1413 """Set the initial solution vector. 1414 1415 Logically collective. 1416 1417 Parameters 1418 ---------- 1419 u 1420 The solution vector. 1421 1422 See Also 1423 -------- 1424 petsc.TSSetSolution 1425 1426 """ 1427 CHKERR(TSSetSolution(self.ts, u.vec)) 1428 1429 def getSolution(self) -> Vec: 1430 """Return the solution at the present timestep. 1431 1432 Not collective. 1433 1434 It is valid to call this routine inside the function that you are 1435 evaluating in order to move to the new timestep. This vector is not 1436 changed until the solution at the next timestep has been calculated. 1437 1438 See Also 1439 -------- 1440 petsc.TSGetSolution 1441 1442 """ 1443 cdef Vec u = Vec() 1444 CHKERR(TSGetSolution(self.ts, &u.vec)) 1445 CHKERR(PetscINCREF(u.obj)) 1446 return u 1447 1448 def setSolution2(self, Vec u, Vec v) -> None: 1449 """Set the initial solution and its time derivative. 1450 1451 Logically collective. 1452 1453 Parameters 1454 ---------- 1455 u 1456 The solution vector. 1457 v 1458 The time derivative vector. 1459 1460 See Also 1461 -------- 1462 petsc.TS2SetSolution 1463 1464 """ 1465 CHKERR(TS2SetSolution(self.ts, u.vec, v.vec)) 1466 1467 def getSolution2(self) -> tuple[Vec, Vec]: 1468 """Return the solution and time derivative at the present timestep. 1469 1470 Not collective. 1471 1472 It is valid to call this routine inside the function that you are 1473 evaluating in order to move to the new timestep. These vectors are not 1474 changed until the solution at the next timestep has been calculated. 1475 1476 See Also 1477 -------- 1478 petsc.TS2GetSolution 1479 1480 """ 1481 cdef Vec u = Vec() 1482 cdef Vec v = Vec() 1483 CHKERR(TS2GetSolution(self.ts, &u.vec, &v.vec)) 1484 CHKERR(PetscINCREF(u.obj)) 1485 CHKERR(PetscINCREF(v.obj)) 1486 return (u, v) 1487 1488 # --- evaluation times --- 1489 1490 def setEvaluationTimes(self, tspan: Sequence[float]) -> None: 1491 """Sets evaluation points where solution will be computed and stored. 1492 1493 Collective. 1494 1495 The solution will be computed and stored for each time 1496 requested. The times must be all increasing and correspond 1497 to the intermediate points for time integration. 1498 `ExactFinalTime.MATCHSTEP` must be used to make the last time step in 1499 each sub-interval match the intermediate points specified. The 1500 intermediate solutions are saved in a vector array that can be accessed 1501 with `getEvaluationSolutions`. 1502 1503 Parameters 1504 ---------- 1505 tspan 1506 The sequence of time points. The first element and the last element 1507 are the initial time and the final time respectively. 1508 1509 Notes 1510 ----- 1511 ``-ts_eval_times <t0, ..., tn>`` sets the time span from the commandline 1512 1513 See Also 1514 -------- 1515 getEvaluationTimes, petsc.TSGetEvaluationTimes 1516 1517 """ 1518 cdef PetscInt nt = 0 1519 cdef PetscReal *rtspan = NULL 1520 cdef unused = oarray_r(tspan, &nt, &rtspan) 1521 CHKERR(TSSetEvaluationTimes(self.ts, nt, rtspan)) 1522 1523 def getEvaluationTimes(self) -> ArrayReal: 1524 """Return the evaluation points. 1525 1526 Not collective. 1527 1528 See Also 1529 -------- 1530 setEvaluationTimes 1531 1532 """ 1533 cdef const PetscReal *rtspan = NULL 1534 cdef PetscInt nt = 0 1535 CHKERR(TSGetEvaluationTimes(self.ts, &nt, &rtspan)) 1536 cdef object tspan = array_r(nt, rtspan) 1537 return tspan 1538 1539 def getEvaluationSolutions(self) -> tuple[ArrayReal, list[Vec]]: 1540 """Return the solutions and the times they were recorded at. 1541 1542 Not collective. 1543 1544 See Also 1545 -------- 1546 setEvaluationTimes 1547 1548 """ 1549 cdef PetscInt nt = 0 1550 cdef PetscVec *sols = NULL 1551 cdef const PetscReal *rtspan = NULL 1552 CHKERR(TSGetEvaluationSolutions(self.ts, &nt, &rtspan, &sols)) 1553 cdef object sollist = None 1554 if sols != NULL: 1555 sollist = [ref_Vec(sols[i]) for i from 0 <= i < nt] 1556 cdef object tspan = array_r(nt, rtspan) 1557 return tspan, sollist 1558 1559 # --- time span --- 1560 1561 def setTimeSpan(self, tspan: Sequence[float]) -> None: 1562 """Set the time span and time points to evaluate solution at. 1563 1564 Collective. 1565 1566 The solution will be computed and stored for each time 1567 requested in the span. The times must be all increasing and correspond 1568 to the intermediate points for time integration. 1569 `ExactFinalTime.MATCHSTEP` must be used to make the last time step in 1570 each sub-interval match the intermediate points specified. The 1571 intermediate solutions are saved in a vector array that can be accessed 1572 with `getEvaluationSolutions`. 1573 1574 Parameters 1575 ---------- 1576 tspan 1577 The sequence of time points. The first element and the last element 1578 are the initial time and the final time respectively. 1579 1580 Notes 1581 ----- 1582 ``-ts_time_span <t0, ..., tf>`` sets the time span from the commandline 1583 1584 See Also 1585 -------- 1586 setEvaluationTimes, petsc.TSSetTimeSpan 1587 1588 """ 1589 cdef PetscInt nt = 0 1590 cdef PetscReal *rtspan = NULL 1591 cdef unused = oarray_r(tspan, &nt, &rtspan) 1592 CHKERR(TSSetTimeSpan(self.ts, nt, rtspan)) 1593 1594 getTimeSpan = getEvaluationTimes 1595 1596 def getTimeSpanSolutions(self) -> list[Vec]: 1597 """Return the solutions at the times in the time span. Deprecated. 1598 1599 Not collective. 1600 1601 See Also 1602 -------- 1603 setTimeSpan, setEvaluationTimes, getEvaluationSolutions 1604 1605 """ 1606 cdef PetscInt nt = 0 1607 cdef PetscVec *sols = NULL 1608 CHKERR(TSGetEvaluationSolutions(self.ts, &nt, NULL, &sols)) 1609 cdef object sollist = None 1610 if sols != NULL: 1611 sollist = [ref_Vec(sols[i]) for i from 0 <= i < nt] 1612 return sollist 1613 1614 # --- inner solver --- 1615 1616 def getSNES(self) -> SNES: 1617 """Return the `SNES` associated with the `TS`. 1618 1619 Not collective. 1620 1621 See Also 1622 -------- 1623 petsc.TSGetSNES 1624 1625 """ 1626 cdef SNES snes = SNES() 1627 CHKERR(TSGetSNES(self.ts, &snes.snes)) 1628 CHKERR(PetscINCREF(snes.obj)) 1629 return snes 1630 1631 def getKSP(self) -> KSP: 1632 """Return the `KSP` associated with the `TS`. 1633 1634 Not collective. 1635 1636 See Also 1637 -------- 1638 petsc.TSGetKSP 1639 1640 """ 1641 cdef KSP ksp = KSP() 1642 CHKERR(TSGetKSP(self.ts, &ksp.ksp)) 1643 CHKERR(PetscINCREF(ksp.obj)) 1644 return ksp 1645 1646 # --- discretization space --- 1647 1648 def getDM(self) -> DM: 1649 """Return the `DM` associated with the `TS`. 1650 1651 Not collective. 1652 1653 Only valid if nonlinear solvers or preconditioners are 1654 used which use the `DM`. 1655 1656 See Also 1657 -------- 1658 petsc.TSGetDM 1659 1660 """ 1661 cdef PetscDM newdm = NULL 1662 CHKERR(TSGetDM(self.ts, &newdm)) 1663 cdef DM dm = subtype_DM(newdm)() 1664 dm.dm = newdm 1665 CHKERR(PetscINCREF(dm.obj)) 1666 return dm 1667 1668 def setDM(self, DM dm) -> None: 1669 """Set the DM that may be used by some nonlinear solvers or preconditioners. 1670 1671 Logically collective. 1672 1673 Parameters 1674 ---------- 1675 dm 1676 The `DM` object. 1677 1678 See Also 1679 -------- 1680 petsc.TSSetDM 1681 1682 """ 1683 CHKERR(TSSetDM(self.ts, dm.dm)) 1684 1685 # --- customization --- 1686 1687 def setTime(self, t: float) -> None: 1688 """Set the time. 1689 1690 Logically collective. 1691 1692 Parameters 1693 ---------- 1694 t 1695 The time. 1696 1697 See Also 1698 -------- 1699 petsc.TSSetTime 1700 1701 """ 1702 cdef PetscReal rval = asReal(t) 1703 CHKERR(TSSetTime(self.ts, rval)) 1704 1705 def getTime(self) -> float: 1706 """Return the time of the most recently completed step. 1707 1708 Not collective. 1709 1710 When called during time step evaluation (e.g. during 1711 residual evaluation or via hooks set using `setPreStep` or 1712 `setPostStep`), the time returned is at the start of the step. 1713 1714 See Also 1715 -------- 1716 petsc.TSGetTime 1717 1718 """ 1719 cdef PetscReal rval = 0 1720 CHKERR(TSGetTime(self.ts, &rval)) 1721 return toReal(rval) 1722 1723 def getPrevTime(self) -> float: 1724 """Return the starting time of the previously completed step. 1725 1726 Not collective. 1727 1728 See Also 1729 -------- 1730 petsc.TSGetPrevTime 1731 1732 """ 1733 cdef PetscReal rval = 0 1734 CHKERR(TSGetPrevTime(self.ts, &rval)) 1735 return toReal(rval) 1736 1737 def getSolveTime(self) -> float: 1738 """Return the time after a call to `solve`. 1739 1740 Not collective. 1741 1742 This time corresponds to the final time set with 1743 `setMaxTime`. 1744 1745 See Also 1746 -------- 1747 petsc.TSGetSolveTime 1748 1749 """ 1750 cdef PetscReal rval = 0 1751 CHKERR(TSGetSolveTime(self.ts, &rval)) 1752 return toReal(rval) 1753 1754 def setTimeStep(self, time_step: float) -> None: 1755 """Set the duration of the timestep. 1756 1757 Logically collective. 1758 1759 Parameters 1760 ---------- 1761 time_step 1762 the duration of the timestep 1763 1764 See Also 1765 -------- 1766 petsc.TSSetTimeStep 1767 1768 """ 1769 cdef PetscReal rval = asReal(time_step) 1770 CHKERR(TSSetTimeStep(self.ts, rval)) 1771 1772 def getTimeStep(self) -> float: 1773 """Return the duration of the current timestep. 1774 1775 Not collective. 1776 1777 See Also 1778 -------- 1779 petsc.TSGetTimeStep 1780 1781 """ 1782 cdef PetscReal tstep = 0 1783 CHKERR(TSGetTimeStep(self.ts, &tstep)) 1784 return toReal(tstep) 1785 1786 def setStepNumber(self, step_number: int) -> None: 1787 """Set the number of steps completed. 1788 1789 Logically collective. 1790 1791 For most uses of the `TS` solvers the user need 1792 not explicitly call `setStepNumber`, as the step counter is 1793 appropriately updated in `solve`/`step`/`rollBack`. Power users may call 1794 this routine to reinitialize timestepping by setting the step counter to 1795 zero (and time to the initial time) to solve a similar problem with 1796 different initial conditions or parameters. It may also be used to 1797 continue timestepping from a previously interrupted run in such a way 1798 that `TS` monitors will be called with a initial nonzero step counter. 1799 1800 Parameters 1801 ---------- 1802 step_number 1803 the number of steps completed 1804 1805 See Also 1806 -------- 1807 petsc.TSSetStepNumber 1808 1809 """ 1810 cdef PetscInt ival = asInt(step_number) 1811 CHKERR(TSSetStepNumber(self.ts, ival)) 1812 1813 def getStepNumber(self) -> int: 1814 """Return the number of time steps completed. 1815 1816 Not collective. 1817 1818 See Also 1819 -------- 1820 petsc.TSGetStepNumber 1821 1822 """ 1823 cdef PetscInt ival = 0 1824 CHKERR(TSGetStepNumber(self.ts, &ival)) 1825 return toInt(ival) 1826 1827 def setMaxTime(self, max_time: float) -> None: 1828 """Set the maximum (final) time. 1829 1830 Logically collective. 1831 1832 Parameters 1833 ---------- 1834 max_time 1835 the final time 1836 1837 Notes 1838 ----- 1839 ``-ts_max_time`` sets the max time from the commandline 1840 1841 See Also 1842 -------- 1843 petsc.TSSetMaxTime 1844 1845 """ 1846 cdef PetscReal rval = asReal(max_time) 1847 CHKERR(TSSetMaxTime(self.ts, rval)) 1848 1849 def getMaxTime(self) -> float: 1850 """Return the maximum (final) time. 1851 1852 Not collective. 1853 1854 Defaults to ``5``. 1855 1856 See Also 1857 -------- 1858 petsc.TSGetMaxTime 1859 1860 """ 1861 cdef PetscReal rval = 0 1862 CHKERR(TSGetMaxTime(self.ts, &rval)) 1863 return toReal(rval) 1864 1865 def setMaxSteps(self, max_steps: int) -> None: 1866 """Set the maximum number of steps to use. 1867 1868 Logically collective. 1869 1870 Defaults to ``5000``. 1871 1872 Parameters 1873 ---------- 1874 max_steps 1875 The maximum number of steps to use. 1876 1877 See Also 1878 -------- 1879 petsc.TSSetMaxSteps 1880 1881 """ 1882 cdef PetscInt ival = asInt(max_steps) 1883 CHKERR(TSSetMaxSteps(self.ts, ival)) 1884 1885 def getMaxSteps(self) -> int: 1886 """Return the maximum number of steps to use. 1887 1888 Not collective. 1889 1890 See Also 1891 -------- 1892 petsc.TSGetMaxSteps 1893 1894 """ 1895 cdef PetscInt ival = 0 1896 CHKERR(TSGetMaxSteps(self.ts, &ival)) 1897 return toInt(ival) 1898 1899 def getSNESIterations(self) -> int: 1900 """Return the total number of nonlinear iterations used by the `TS`. 1901 1902 Not collective. 1903 1904 This counter is reset to zero for each successive call 1905 to `solve`. 1906 1907 See Also 1908 -------- 1909 petsc.TSGetSNESIterations 1910 1911 """ 1912 cdef PetscInt n = 0 1913 CHKERR(TSGetSNESIterations(self.ts, &n)) 1914 return toInt(n) 1915 1916 def getKSPIterations(self) -> int: 1917 """Return the total number of linear iterations used by the `TS`. 1918 1919 Not collective. 1920 1921 This counter is reset to zero for each successive call 1922 to `solve`. 1923 1924 See Also 1925 -------- 1926 petsc.TSGetKSPIterations 1927 1928 """ 1929 cdef PetscInt n = 0 1930 CHKERR(TSGetKSPIterations(self.ts, &n)) 1931 return toInt(n) 1932 1933 def setMaxStepRejections(self, n: int) -> None: 1934 """Set the maximum number of step rejections before a time step fails. 1935 1936 Not collective. 1937 1938 Parameters 1939 ---------- 1940 n 1941 The maximum number of rejected steps, use ``-1`` for unlimited. 1942 1943 Notes 1944 ----- 1945 ``-ts_max_step_rejections`` can be used to set this from the commandline 1946 1947 See Also 1948 -------- 1949 petsc.TSSetMaxStepRejections 1950 1951 """ 1952 cdef PetscInt rej = asInt(n) 1953 CHKERR(TSSetMaxStepRejections(self.ts, rej)) 1954 1955 def getStepRejections(self) -> int: 1956 """Return the total number of rejected steps. 1957 1958 Not collective. 1959 1960 This counter is reset to zero for each successive call 1961 to `solve`. 1962 1963 See Also 1964 -------- 1965 petsc.TSGetStepRejections 1966 1967 """ 1968 cdef PetscInt n = 0 1969 CHKERR(TSGetStepRejections(self.ts, &n)) 1970 return toInt(n) 1971 1972 def setMaxSNESFailures(self, n: int) -> None: 1973 """Set the maximum number of SNES solves failures allowed. 1974 1975 Not collective. 1976 1977 Parameters 1978 ---------- 1979 n 1980 The maximum number of failed nonlinear solver, use ``-1`` for unlimited. 1981 1982 See Also 1983 -------- 1984 petsc.TSSetMaxSNESFailures 1985 1986 """ 1987 cdef PetscInt fails = asInt(n) 1988 CHKERR(TSSetMaxSNESFailures(self.ts, fails)) 1989 1990 def getSNESFailures(self) -> int: 1991 """Return the total number of failed `SNES` solves in the `TS`. 1992 1993 Not collective. 1994 1995 This counter is reset to zero for each successive call 1996 to `solve`. 1997 1998 See Also 1999 -------- 2000 petsc.TSGetSNESFailures 2001 2002 """ 2003 cdef PetscInt n = 0 2004 CHKERR(TSGetSNESFailures(self.ts, &n)) 2005 return toInt(n) 2006 2007 def setErrorIfStepFails(self, flag: bool = True) -> None: 2008 """Immediately error is no step succeeds. 2009 2010 Not collective. 2011 2012 Parameters 2013 ---------- 2014 flag 2015 Enable to error if no step succeeds. 2016 2017 Notes 2018 ----- 2019 ``-ts_error_if_step_fails`` to enable from the commandline. 2020 2021 See Also 2022 -------- 2023 petsc.TSSetErrorIfStepFails 2024 2025 """ 2026 cdef PetscBool bval = flag 2027 CHKERR(TSSetErrorIfStepFails(self.ts, bval)) 2028 2029 def setTolerances(self, rtol: float | None = None, atol: float | None = None) -> None: 2030 """Set tolerances for local truncation error when using an adaptive controller. 2031 2032 Logically collective. 2033 2034 Parameters 2035 ---------- 2036 rtol 2037 The relative tolerance, `DETERMINE` to use the value 2038 when the object's type was set, or `None` to leave the 2039 current value. 2040 atol 2041 The absolute tolerance, `DETERMINE` to use the 2042 value when the object's type was set, or `None` to 2043 leave the current value. 2044 2045 Notes 2046 ----- 2047 ``-ts_rtol`` and ``-ts_atol`` may be used to set values from the commandline. 2048 2049 See Also 2050 -------- 2051 petsc.TSSetTolerances 2052 2053 """ 2054 cdef PetscReal rrtol = PETSC_CURRENT 2055 cdef PetscReal ratol = PETSC_CURRENT 2056 cdef PetscVec vrtol = NULL 2057 cdef PetscVec vatol = NULL 2058 if rtol is None: 2059 pass 2060 elif isinstance(rtol, Vec): 2061 vrtol = (<Vec>rtol).vec 2062 else: 2063 rrtol = asReal(rtol) 2064 if atol is None: 2065 pass 2066 elif isinstance(atol, Vec): 2067 vatol = (<Vec>atol).vec 2068 else: 2069 ratol = asReal(atol) 2070 CHKERR(TSSetTolerances(self.ts, ratol, vatol, rrtol, vrtol)) 2071 2072 def getTolerances(self) ->tuple[float, float]: 2073 """Return the tolerances for local truncation error. 2074 2075 Logically collective. 2076 2077 Returns 2078 ------- 2079 rtol : float 2080 the relative tolerance 2081 atol : float 2082 the absolute tolerance 2083 2084 See Also 2085 -------- 2086 petsc.TSGetTolerances 2087 2088 """ 2089 cdef PetscReal rrtol = PETSC_DETERMINE 2090 cdef PetscReal ratol = PETSC_DETERMINE 2091 cdef PetscVec vrtol = NULL 2092 cdef PetscVec vatol = NULL 2093 CHKERR(TSGetTolerances(self.ts, &ratol, &vatol, &rrtol, &vrtol)) 2094 cdef object rtol = None 2095 if vrtol != NULL: 2096 rtol = ref_Vec(vrtol) 2097 else: 2098 rtol = toReal(rrtol) 2099 cdef object atol = None 2100 if vatol != NULL: 2101 atol = ref_Vec(vatol) 2102 else: 2103 atol = toReal(ratol) 2104 return (rtol, atol) 2105 2106 def setExactFinalTime(self, option: ExactFinalTime) -> None: 2107 """Set method of computing the final time step. 2108 2109 Logically collective. 2110 2111 Parameters 2112 ---------- 2113 option 2114 The exact final time option 2115 2116 Notes 2117 ----- 2118 ``-ts_exact_final_time`` may be used to specify from the commandline. 2119 2120 See Also 2121 -------- 2122 petsc.TSSetExactFinalTime 2123 2124 """ 2125 cdef PetscTSExactFinalTimeOption oval = option 2126 CHKERR(TSSetExactFinalTime(self.ts, oval)) 2127 2128 def setConvergedReason(self, reason: ConvergedReason) -> None: 2129 """Set the reason for handling the convergence of `solve`. 2130 2131 Logically collective. 2132 2133 Can only be called when `solve` is active and 2134 ``reason`` must contain common value. 2135 2136 Parameters 2137 ---------- 2138 reason 2139 The reason for convergence. 2140 2141 See Also 2142 -------- 2143 petsc.TSSetConvergedReason 2144 2145 """ 2146 cdef PetscTSConvergedReason cval = reason 2147 CHKERR(TSSetConvergedReason(self.ts, cval)) 2148 2149 def getConvergedReason(self) -> ConvergedReason: 2150 """Return the reason the `TS` step was stopped. 2151 2152 Not collective. 2153 2154 Can only be called once `solve` is complete. 2155 2156 See Also 2157 -------- 2158 petsc.TSGetConvergedReason 2159 2160 """ 2161 cdef PetscTSConvergedReason reason = TS_CONVERGED_ITERATING 2162 CHKERR(TSGetConvergedReason(self.ts, &reason)) 2163 return reason 2164 2165 # --- monitoring --- 2166 2167 def setMonitor( 2168 self, 2169 monitor: TSMonitorFunction | None, 2170 args : tuple[Any, ...] | None = None, 2171 kargs : dict[str, Any] | None = None) -> None: 2172 """Set an additional monitor to the `TS`. 2173 2174 Logically collective. 2175 2176 Parameters 2177 ---------- 2178 monitor 2179 The custom monitor function. 2180 args 2181 Additional positional arguments for ``monitor``. 2182 kargs 2183 Additional keyword arguments for ``monitor``. 2184 2185 See Also 2186 -------- 2187 petsc.TSMonitorSet 2188 2189 """ 2190 if monitor is None: return 2191 cdef object monitorlist = self.get_attr('__monitor__') 2192 if monitorlist is None: 2193 monitorlist = [] 2194 self.set_attr('__monitor__', monitorlist) 2195 CHKERR(TSMonitorSet(self.ts, TS_Monitor, NULL, NULL)) 2196 if args is None: args = () 2197 if kargs is None: kargs = {} 2198 context = (monitor, args, kargs) 2199 monitorlist.append(context) 2200 2201 def getMonitor(self) -> list[tuple[TSMonitorFunction, tuple[Any, ...], dict[str, Any]]]: 2202 """Return the monitor. 2203 2204 Not collective. 2205 2206 See Also 2207 -------- 2208 setMonitor 2209 2210 """ 2211 return self.get_attr('__monitor__') 2212 2213 def monitorCancel(self) -> None: 2214 """Clear all the monitors that have been set. 2215 2216 Logically collective. 2217 2218 See Also 2219 -------- 2220 petsc.TSMonitorCancel 2221 2222 """ 2223 self.set_attr('__monitor__', None) 2224 CHKERR(TSMonitorCancel(self.ts)) 2225 2226 cancelMonitor = monitorCancel 2227 2228 def monitor(self, step: int, time: float, Vec u=None) -> None: 2229 """Monitor the solve. 2230 2231 Collective. 2232 2233 Parameters 2234 ---------- 2235 step 2236 The step number that has just completed. 2237 time 2238 The model time of the state. 2239 u 2240 The state at the current model time. 2241 2242 See Also 2243 -------- 2244 petsc.TSMonitor 2245 2246 """ 2247 cdef PetscInt ival = asInt(step) 2248 cdef PetscReal rval = asReal(time) 2249 cdef PetscVec uvec = NULL 2250 if u is not None: uvec = u.vec 2251 if uvec == NULL: 2252 CHKERR(TSGetSolution(self.ts, &uvec)) 2253 CHKERR(TSMonitor(self.ts, ival, rval, uvec)) 2254 2255 # --- event handling --- 2256 2257 def setEventHandler( 2258 self, 2259 direction: Sequence[int], 2260 terminate: Sequence[bool], 2261 indicator: TSIndicatorFunction | None, 2262 postevent: TSPostEventFunction | None = None, 2263 args: tuple[Any, ...] | None = None, 2264 kargs: dict[str, Any] | None = None) -> None: 2265 """Set a function used for detecting events. 2266 2267 Logically collective. 2268 2269 Parameters 2270 ---------- 2271 direction 2272 Direction of zero crossing to be detected {-1, 0, +1}. 2273 terminate 2274 Flags for each event to indicate stepping should be terminated. 2275 indicator 2276 Function for defining the indicator-functions marking the events 2277 postevent 2278 Function to execute after the event 2279 args 2280 Additional positional arguments for ``indicator``. 2281 kargs 2282 Additional keyword arguments for ``indicator``. 2283 2284 See Also 2285 -------- 2286 petsc.TSSetEventHandler 2287 2288 """ 2289 cdef PetscInt ndirs = 0 2290 cdef PetscInt *idirs = NULL 2291 direction = iarray_i(direction, &ndirs, &idirs) 2292 2293 cdef PetscInt nterm = 0 2294 cdef PetscBool *iterm = NULL 2295 terminate = iarray_b(terminate, &nterm, &iterm) 2296 assert nterm == ndirs 2297 2298 cdef PetscInt nevents = ndirs 2299 if indicator is not None: 2300 if args is None: args = () 2301 if kargs is None: kargs = {} 2302 self.set_attr('__indicator__', (indicator, args, kargs)) 2303 if postevent is not None: 2304 self.set_attr('__postevent__', (postevent, args, kargs)) 2305 CHKERR(TSSetEventHandler(self.ts, nevents, idirs, iterm, TS_Indicator, TS_PostEvent, <void*>NULL)) 2306 else: 2307 self.set_attr('__postevent__', None) 2308 CHKERR(TSSetEventHandler(self.ts, nevents, idirs, iterm, TS_Indicator, NULL, <void*>NULL)) 2309 else: 2310 CHKERR(TSSetEventHandler(self.ts, nevents, idirs, iterm, NULL, NULL, <void*>NULL)) 2311 2312 def setEventTolerances(self, tol: float | None = None, vtol: Sequence[float] | None = None) -> None: 2313 """Set tolerances for event zero crossings when using event handler. 2314 2315 Logically collective. 2316 2317 ``setEventHandler`` must have already been called. 2318 2319 Parameters 2320 ---------- 2321 tol 2322 The scalar tolerance or `None` to leave at the current value 2323 vtol 2324 A sequence of scalar tolerance for each event. Used in preference to 2325 ``tol`` if present. Set to `None` to leave at the current value. 2326 2327 Notes 2328 ----- 2329 ``-ts_event_tol`` can be used to set values from the commandline. 2330 2331 See Also 2332 -------- 2333 petsc.TSSetEventTolerances 2334 2335 """ 2336 cdef PetscInt nevents = 0 2337 cdef PetscReal tolr = PETSC_CURRENT 2338 cdef PetscInt ntolr = 0 2339 cdef PetscReal *vtolr = NULL 2340 if tol is not None: 2341 tolr = asReal(tol) 2342 if vtol is not None: 2343 CHKERR(TSGetNumEvents(self.ts, &nevents)) 2344 vtol = iarray_r(vtol, &ntolr, &vtolr) 2345 assert ntolr == nevents 2346 CHKERR(TSSetEventTolerances(self.ts, tolr, vtolr)) 2347 2348 def getNumEvents(self) -> int: 2349 """Return the number of events. 2350 2351 Logically collective. 2352 2353 See Also 2354 -------- 2355 petsc.TSGetNumEvents 2356 2357 """ 2358 cdef PetscInt nevents = 0 2359 CHKERR(TSGetNumEvents(self.ts, &nevents)) 2360 return toInt(nevents) 2361 2362 # --- solving --- 2363 2364 def setPreStep( 2365 self, 2366 prestep: TSPreStepFunction | None, 2367 args: tuple[Any, ...] | None = None, 2368 kargs: dict[str, Any] | None = None) -> None: 2369 """Set a function to be called at the beginning of each time step. 2370 2371 Logically collective. 2372 2373 Parameters 2374 ---------- 2375 prestep 2376 The function to be called at the beginning of each step. 2377 args 2378 Additional positional arguments for ``prestep``. 2379 kargs 2380 Additional keyword arguments for ``prestep``. 2381 2382 See Also 2383 -------- 2384 petsc.TSSetPreStep 2385 2386 """ 2387 if prestep is not None: 2388 if args is None: args = () 2389 if kargs is None: kargs = {} 2390 context = (prestep, args, kargs) 2391 self.set_attr('__prestep__', context) 2392 CHKERR(TSSetPreStep(self.ts, TS_PreStep)) 2393 else: 2394 self.set_attr('__prestep__', None) 2395 CHKERR(TSSetPreStep(self.ts, NULL)) 2396 2397 def getPreStep(self) -> tuple[TSPreStepFunction, tuple[Any, ...] | None, dict[str, Any] | None]: 2398 """Return the prestep function. 2399 2400 Not collective. 2401 2402 See Also 2403 -------- 2404 setPreStep 2405 2406 """ 2407 return self.get_attr('__prestep__') 2408 2409 def setPostStep(self, 2410 poststep: TSPostStepFunction | None, 2411 args: tuple[Any, ...] | None = None, 2412 kargs: dict[str, Any] | None = None) -> None: 2413 """Set a function to be called at the end of each time step. 2414 2415 Logically collective. 2416 2417 Parameters 2418 ---------- 2419 poststep 2420 The function to be called at the end of each step. 2421 args 2422 Additional positional arguments for ``poststep``. 2423 kargs 2424 Additional keyword arguments for ``poststep``. 2425 2426 See Also 2427 -------- 2428 petsc.TSSetPostStep 2429 2430 """ 2431 if poststep is not None: 2432 if args is None: args = () 2433 if kargs is None: kargs = {} 2434 context = (poststep, args, kargs) 2435 self.set_attr('__poststep__', context) 2436 CHKERR(TSSetPostStep(self.ts, TS_PostStep)) 2437 else: 2438 self.set_attr('__poststep__', None) 2439 CHKERR(TSSetPostStep(self.ts, NULL)) 2440 2441 def getPostStep(self) -> tuple[TSPostStepFunction, tuple[Any, ...] | None, dict[str, Any] | None]: 2442 """Return the poststep function.""" 2443 return self.get_attr('__poststep__') 2444 2445 def setUp(self) -> None: 2446 """Set up the internal data structures for the `TS`. 2447 2448 Collective. 2449 2450 See Also 2451 -------- 2452 petsc.TSSetUp 2453 2454 """ 2455 CHKERR(TSSetUp(self.ts)) 2456 2457 def reset(self) -> None: 2458 """Reset the `TS`, removing any allocated vectors and matrices. 2459 2460 Collective. 2461 2462 See Also 2463 -------- 2464 petsc.TSReset 2465 2466 """ 2467 CHKERR(TSReset(self.ts)) 2468 2469 def step(self) -> None: 2470 """Take one step. 2471 2472 Collective. 2473 2474 The preferred interface for the `TS` solvers is `solve`. If 2475 you need to execute code at the beginning or ending of each step, use 2476 `setPreStep` and `setPostStep` respectively. 2477 2478 See Also 2479 -------- 2480 petsc.TSStep 2481 2482 """ 2483 CHKERR(TSStep(self.ts)) 2484 2485 def restartStep(self) -> None: 2486 """Flag the solver to restart the next step. 2487 2488 Collective. 2489 2490 Multistep methods like TSBDF or Runge-Kutta methods with 2491 FSAL property require restarting the solver in the event of 2492 discontinuities. These discontinuities may be introduced as a 2493 consequence of explicitly modifications to the solution vector (which 2494 PETSc attempts to detect and handle) or problem coefficients (which 2495 PETSc is not able to detect). For the sake of correctness and maximum 2496 safety, users are expected to call TSRestart() whenever they introduce 2497 discontinuities in callback routines (e.g. prestep and poststep 2498 routines, or implicit/rhs function routines with discontinuous source 2499 terms). 2500 2501 See Also 2502 -------- 2503 petsc.TSRestartStep 2504 2505 """ 2506 CHKERR(TSRestartStep(self.ts)) 2507 2508 def rollBack(self) -> None: 2509 """Roll back one time step. 2510 2511 Collective. 2512 2513 See Also 2514 -------- 2515 petsc.TSRollBack 2516 2517 """ 2518 CHKERR(TSRollBack(self.ts)) 2519 2520 def solve(self, Vec u=None) -> None: 2521 """Step the requested number of timesteps. 2522 2523 Collective. 2524 2525 Parameters 2526 ---------- 2527 u 2528 The solution vector. Can be `None`. 2529 2530 See Also 2531 -------- 2532 petsc.TSSolve 2533 2534 """ 2535 cdef PetscVec uvec=NULL 2536 if u is not None: uvec = u.vec 2537 CHKERR(TSSolve(self.ts, uvec)) 2538 2539 def interpolate(self, t: float, Vec u) -> None: 2540 """Interpolate the solution to a given time. 2541 2542 Collective. 2543 2544 Parameters 2545 ---------- 2546 t 2547 The time to interpolate. 2548 u 2549 The state vector to interpolate. 2550 2551 See Also 2552 -------- 2553 petsc.TSInterpolate 2554 2555 """ 2556 cdef PetscReal rval = asReal(t) 2557 CHKERR(TSInterpolate(self.ts, rval, u.vec)) 2558 2559 def setStepLimits(self, hmin: float, hmax: float) -> None: 2560 """Set the minimum and maximum allowed step sizes. 2561 2562 Logically collective. 2563 2564 Parameters 2565 ---------- 2566 hmin 2567 the minimum step size 2568 hmax 2569 the maximum step size 2570 2571 See Also 2572 -------- 2573 petsc.TSAdaptSetStepLimits 2574 2575 """ 2576 cdef PetscTSAdapt tsadapt = NULL 2577 cdef PetscReal hminr = toReal(hmin) 2578 cdef PetscReal hmaxr = toReal(hmax) 2579 TSGetAdapt(self.ts, &tsadapt) 2580 CHKERR(TSAdaptSetStepLimits(tsadapt, hminr, hmaxr)) 2581 2582 def getStepLimits(self) -> tuple[float, float]: 2583 """Return the minimum and maximum allowed time step sizes. 2584 2585 Not collective. 2586 2587 See Also 2588 -------- 2589 petsc.TSAdaptGetStepLimits 2590 2591 """ 2592 cdef PetscTSAdapt tsadapt = NULL 2593 cdef PetscReal hminr = 0. 2594 cdef PetscReal hmaxr = 0. 2595 TSGetAdapt(self.ts, &tsadapt) 2596 CHKERR(TSAdaptGetStepLimits(tsadapt, &hminr, &hmaxr)) 2597 return (asReal(hminr), asReal(hmaxr)) 2598 2599 # --- Adjoint methods --- 2600 2601 def setSaveTrajectory(self) -> None: 2602 """Enable to save solutions as an internal `TS` trajectory. 2603 2604 Collective. 2605 2606 This routine should be called after all `TS` options have 2607 been set. 2608 2609 Notes 2610 ----- 2611 ``-ts_save_trajectory`` can be used to save a trajectory to a file. 2612 2613 See Also 2614 -------- 2615 petsc.TSSetSaveTrajectory 2616 2617 """ 2618 CHKERR(TSSetSaveTrajectory(self.ts)) 2619 2620 def removeTrajectory(self) -> None: 2621 """Remove the internal `TS` trajectory object. 2622 2623 Collective. 2624 2625 See Also 2626 -------- 2627 petsc.TSRemoveTrajectory 2628 2629 """ 2630 CHKERR(TSRemoveTrajectory(self.ts)) 2631 2632 def getCostIntegral(self) -> Vec: 2633 """Return a vector of values of the integral term in the cost functions. 2634 2635 Not collective. 2636 2637 See Also 2638 -------- 2639 petsc.TSGetCostIntegral 2640 2641 """ 2642 cdef Vec cost = Vec() 2643 CHKERR(TSGetCostIntegral(self.ts, &cost.vec)) 2644 CHKERR(PetscINCREF(cost.obj)) 2645 return cost 2646 2647 def setCostGradients( 2648 self, 2649 vl: Vec | Sequence[Vec] | None, 2650 vm: Vec | Sequence[Vec] | None = None) -> None: 2651 """Set the cost gradients. 2652 2653 Logically collective. 2654 2655 Parameters 2656 ---------- 2657 vl 2658 gradients with respect to the initial condition variables, the 2659 dimension and parallel layout of these vectors is the same as the 2660 ODE solution vector 2661 vm 2662 gradients with respect to the parameters, the number of entries in 2663 these vectors is the same as the number of parameters 2664 2665 See Also 2666 -------- 2667 petsc.TSSetCostGradients 2668 2669 """ 2670 cdef PetscInt n = 0 2671 cdef PetscVec *vecl = NULL 2672 cdef PetscVec *vecm = NULL 2673 cdef mem1 = None, mem2 = None 2674 if isinstance(vl, Vec): vl = [vl] 2675 if isinstance(vm, Vec): vm = [vm] 2676 if vl is not None: 2677 n = <PetscInt>len(vl) 2678 elif vm is not None: 2679 n = <PetscInt>len(vm) 2680 if vl is not None: 2681 assert len(vl) == <Py_ssize_t>n 2682 mem1 = oarray_p(empty_p(<PetscInt>n), NULL, <void**>&vecl) 2683 for i from 0 <= i < n: 2684 vecl[i] = (<Vec?>vl[i]).vec 2685 if vm is not None: 2686 assert len(vm) == <Py_ssize_t>n 2687 mem2 = oarray_p(empty_p(<PetscInt>n), NULL, <void**>&vecm) 2688 for i from 0 <= i < n: 2689 vecm[i] = (<Vec?>vm[i]).vec 2690 self.set_attr('__costgradients_memory', (mem1, mem2)) 2691 CHKERR(TSSetCostGradients(self.ts, n, vecl, vecm)) 2692 2693 def getCostGradients(self) -> tuple[list[Vec], list[Vec]]: 2694 """Return the cost gradients. 2695 2696 Not collective. 2697 2698 See Also 2699 -------- 2700 setCostGradients, petsc.TSGetCostGradients 2701 2702 """ 2703 cdef PetscInt n = 0 2704 cdef PetscVec *vecl = NULL 2705 cdef PetscVec *vecm = NULL 2706 CHKERR(TSGetCostGradients(self.ts, &n, &vecl, &vecm)) 2707 cdef object vl = None, vm = None 2708 if vecl != NULL: 2709 vl = [ref_Vec(vecl[i]) for i from 0 <= i < n] 2710 if vecm != NULL: 2711 vm = [ref_Vec(vecm[i]) for i from 0 <= i < n] 2712 return (vl, vm) 2713 2714 def setRHSJacobianP( 2715 self, 2716 jacobianp: TSRHSJacobianP | None, 2717 Mat A=None, 2718 args: tuple[Any, ...] | None = None, 2719 kargs: dict[str, Any] | None = None) -> None: 2720 """Set the function that computes the Jacobian with respect to the parameters. 2721 2722 Logically collective. 2723 2724 Parameters 2725 ---------- 2726 jacobianp 2727 The user-defined function. 2728 A 2729 The matrix into which the Jacobian will be computed. 2730 args 2731 Additional positional arguments for ``jacobianp``. 2732 kargs 2733 Additional keyword arguments for ``jacobianp``. 2734 2735 See Also 2736 -------- 2737 petsc.TSSetRHSJacobianP 2738 2739 """ 2740 cdef PetscMat Amat=NULL 2741 if A is not None: Amat = A.mat 2742 if jacobianp is not None: 2743 if args is None: args = () 2744 if kargs is None: kargs = {} 2745 context = (jacobianp, args, kargs) 2746 self.set_attr('__rhsjacobianp__', context) 2747 CHKERR(TSSetRHSJacobianP(self.ts, Amat, TS_RHSJacobianP, <void*>context)) 2748 else: 2749 CHKERR(TSSetRHSJacobianP(self.ts, Amat, NULL, NULL)) 2750 2751 def createQuadratureTS(self, forward: bool = True) -> TS: 2752 """Create a sub `TS` that evaluates integrals over time. 2753 2754 Collective. 2755 2756 Parameters 2757 ---------- 2758 forward 2759 Enable to evaluate forward in time. 2760 2761 See Also 2762 -------- 2763 petsc.TSCreateQuadratureTS 2764 2765 """ 2766 cdef TS qts = TS() 2767 cdef PetscBool fwd = forward 2768 CHKERR(TSCreateQuadratureTS(self.ts, fwd, &qts.ts)) 2769 CHKERR(PetscINCREF(qts.obj)) 2770 return qts 2771 2772 def getQuadratureTS(self) -> tuple[bool, TS]: 2773 """Return the sub `TS` that evaluates integrals over time. 2774 2775 Not collective. 2776 2777 Returns 2778 ------- 2779 forward : bool 2780 True if evaluating the integral forward in time 2781 qts : TS 2782 The sub `TS` 2783 2784 See Also 2785 -------- 2786 petsc.TSGetQuadratureTS 2787 2788 """ 2789 cdef TS qts = TS() 2790 cdef PetscBool fwd = PETSC_FALSE 2791 CHKERR(TSGetQuadratureTS(self.ts, &fwd, &qts.ts)) 2792 CHKERR(PetscINCREF(qts.obj)) 2793 return (toBool(fwd), qts) 2794 2795 def setRHSJacobianP( 2796 self, 2797 rhsjacobianp: TSRHSJacobianP | None, 2798 Mat A=None, 2799 args: tuple[Any, ...] | None = None, 2800 kargs: dict[str, Any] | None = None) -> None: 2801 """Set the function that computes the Jacobian with respect to the parameters. 2802 2803 Collective. 2804 2805 Parameters 2806 ---------- 2807 rhsjacobianp 2808 The function to compute the Jacobian 2809 A 2810 The JacobianP matrix 2811 args 2812 Additional positional arguments for ``rhsjacobianp``. 2813 kargs 2814 Additional keyword arguments for ``rhsjacobianp``. 2815 2816 See Also 2817 -------- 2818 petsc.TSSetRHSJacobianP 2819 2820 """ 2821 cdef PetscMat Amat=NULL 2822 if A is not None: Amat = A.mat 2823 if rhsjacobianp is not None: 2824 if args is None: args = () 2825 if kargs is None: kargs = {} 2826 context = (rhsjacobianp, args, kargs) 2827 self.set_attr('__rhsjacobianp__', context) 2828 CHKERR(TSSetRHSJacobianP(self.ts, Amat, TS_RHSJacobianP, <void*>context)) 2829 else: 2830 CHKERR(TSSetRHSJacobianP(self.ts, Amat, NULL, NULL)) 2831 2832 def computeRHSJacobianP(self, t: float, Vec x, Mat J) -> None: 2833 """Run the user-defined JacobianP function. 2834 2835 Collective. 2836 2837 Parameters 2838 ---------- 2839 t 2840 The time at which to compute the Jacobian. 2841 x 2842 The solution at which to compute the Jacobian. 2843 J 2844 The output Jacobian matrix. 2845 2846 See Also 2847 -------- 2848 petsc.TSComputeRHSJacobianP 2849 2850 """ 2851 cdef PetscReal rval = asReal(t) 2852 CHKERR(TSComputeRHSJacobianP(self.ts, rval, x.vec, J.mat)) 2853 2854 def adjointSetSteps(self, adjoint_steps: int) -> None: 2855 """Set the number of steps the adjoint solver should take backward in time. 2856 2857 Logically collective. 2858 2859 Parameters 2860 ---------- 2861 adjoint_steps 2862 The number of steps to take. 2863 2864 See Also 2865 -------- 2866 petsc.TSAdjointSetSteps 2867 2868 """ 2869 cdef PetscInt ival = asInt(adjoint_steps) 2870 CHKERR(TSAdjointSetSteps(self.ts, ival)) 2871 2872 def adjointSetUp(self) -> None: 2873 """Set up the internal data structures for the later use of an adjoint solver. 2874 2875 Collective. 2876 2877 See Also 2878 -------- 2879 petsc.TSAdjointSetUp 2880 2881 """ 2882 CHKERR(TSAdjointSetUp(self.ts)) 2883 2884 def adjointSolve(self) -> None: 2885 """Solve the discrete adjoint problem for an ODE/DAE. 2886 2887 Collective. 2888 2889 See Also 2890 -------- 2891 petsc.TSAdjointSolve 2892 2893 """ 2894 CHKERR(TSAdjointSolve(self.ts)) 2895 2896 def adjointStep(self) -> None: 2897 """Step one time step backward in the adjoint run. 2898 2899 Collective. 2900 2901 See Also 2902 -------- 2903 petsc.TSAdjointStep 2904 2905 """ 2906 CHKERR(TSAdjointStep(self.ts)) 2907 2908 def adjointReset(self) -> None: 2909 """Reset a `TS`, removing any allocated vectors and matrices. 2910 2911 Collective. 2912 2913 See Also 2914 -------- 2915 petsc.TSAdjointReset 2916 2917 """ 2918 CHKERR(TSAdjointReset(self.ts)) 2919 2920 # --- Python --- 2921 2922 def createPython(self, context: Any = None, comm: Comm | None = None) -> Self: 2923 """Create an integrator of Python type. 2924 2925 Collective. 2926 2927 Parameters 2928 ---------- 2929 context 2930 An instance of the Python class implementing the required methods. 2931 comm 2932 MPI communicator, defaults to `Sys.getDefaultComm`. 2933 2934 See Also 2935 -------- 2936 petsc_python_ts, setType, setPythonContext, Type.PYTHON 2937 2938 """ 2939 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 2940 cdef PetscTS newts = NULL 2941 CHKERR(TSCreate(ccomm, &newts)) 2942 CHKERR(PetscCLEAR(self.obj)); self.ts = newts 2943 CHKERR(TSSetType(self.ts, TSPYTHON)) 2944 CHKERR(TSPythonSetContext(self.ts, <void*>context)) 2945 return self 2946 2947 def setPythonContext(self, context: Any) -> None: 2948 """Set the instance of the class implementing the required Python methods. 2949 2950 Not collective. 2951 2952 See Also 2953 -------- 2954 petsc_python_ts, getPythonContext 2955 2956 """ 2957 CHKERR(TSPythonSetContext(self.ts, <void*>context)) 2958 2959 def getPythonContext(self) -> Any: 2960 """Return the instance of the class implementing the required Python methods. 2961 2962 Not collective. 2963 2964 See Also 2965 -------- 2966 petsc_python_ts, setPythonContext 2967 2968 """ 2969 cdef void *context = NULL 2970 CHKERR(TSPythonGetContext(self.ts, &context)) 2971 if context == NULL: return None 2972 else: return <object> context 2973 2974 def setPythonType(self, py_type: str) -> None: 2975 """Set the fully qualified Python name of the class to be used. 2976 2977 Collective. 2978 2979 See Also 2980 -------- 2981 petsc_python_ts, setPythonContext, getPythonType, petsc.TSPythonSetType 2982 2983 """ 2984 cdef const char *cval = NULL 2985 py_type = str2bytes(py_type, &cval) 2986 CHKERR(TSPythonSetType(self.ts, cval)) 2987 2988 def getPythonType(self) -> str: 2989 """Return the fully qualified Python name of the class used by the solver. 2990 2991 Not collective. 2992 2993 See Also 2994 -------- 2995 petsc_python_ts, setPythonContext, setPythonType, petsc.TSPythonGetType 2996 2997 """ 2998 cdef const char *cval = NULL 2999 CHKERR(TSPythonGetType(self.ts, &cval)) 3000 return bytes2str(cval) 3001 3002 # --- Theta --- 3003 3004 def setTheta(self, theta: float) -> None: 3005 """Set the abscissa of the stage in ``(0, 1]`` for `Type.THETA`. 3006 3007 Logically collective. 3008 3009 Parameters 3010 ---------- 3011 theta 3012 stage abscissa 3013 3014 Notes 3015 ----- 3016 ``-ts_theta_theta`` can be used to set a value from the commandline. 3017 3018 See Also 3019 -------- 3020 petsc.TSThetaSetTheta 3021 3022 """ 3023 cdef PetscReal rval = asReal(theta) 3024 CHKERR(TSThetaSetTheta(self.ts, rval)) 3025 3026 def getTheta(self) -> float: 3027 """Return the abscissa of the stage in ``(0, 1]`` for `Type.THETA`. 3028 3029 Not collective. 3030 3031 See Also 3032 -------- 3033 petsc.TSThetaGetTheta 3034 3035 """ 3036 cdef PetscReal rval = 0 3037 CHKERR(TSThetaGetTheta(self.ts, &rval)) 3038 return toReal(rval) 3039 3040 def setThetaEndpoint(self, flag=True) -> None: 3041 """Set to use the endpoint variant of `Type.THETA`. 3042 3043 Logically collective. 3044 3045 Parameters 3046 ---------- 3047 flag 3048 Enable to use the endpoint variant. 3049 3050 See Also 3051 -------- 3052 petsc.TSThetaSetEndpoint 3053 3054 """ 3055 cdef PetscBool bval = flag 3056 CHKERR(TSThetaSetEndpoint(self.ts, bval)) 3057 3058 def getThetaEndpoint(self) -> bool: 3059 """Return whether the endpoint variable of `Type.THETA` is used. 3060 3061 Not collective. 3062 3063 See Also 3064 -------- 3065 petsc.TSThetaGetEndpoint 3066 3067 """ 3068 cdef PetscBool flag = PETSC_FALSE 3069 CHKERR(TSThetaGetEndpoint(self.ts, &flag)) 3070 return toBool(flag) 3071 3072 # --- Alpha --- 3073 3074 def setAlphaRadius(self, radius: float) -> None: 3075 """Set the spectral radius for `Type.ALPHA`. 3076 3077 Logically collective. 3078 3079 Parameters 3080 ---------- 3081 radius 3082 the spectral radius 3083 3084 Notes 3085 ----- 3086 ``-ts_alpha_radius`` can be used to set this from the commandline. 3087 3088 See Also 3089 -------- 3090 petsc.TSAlphaSetRadius 3091 3092 """ 3093 cdef PetscReal rval = asReal(radius) 3094 CHKERR(TSAlphaSetRadius(self.ts, rval)) 3095 3096 def setAlphaParams( 3097 self, 3098 alpha_m: float | None = None, 3099 alpha_f: float | None = None, 3100 gamma: float | None = None) -> None: 3101 """Set the algorithmic parameters for `Type.ALPHA`. 3102 3103 Logically collective. 3104 3105 Users should call `setAlphaRadius`. 3106 3107 Parameters 3108 ---------- 3109 alpha_m 3110 Parameter, leave `None` to keep current value. 3111 alpha_f 3112 Parameter, leave `None` to keep current value. 3113 gamma 3114 Parameter, leave `None` to keep current value. 3115 3116 See Also 3117 -------- 3118 petsc.TSAlphaSetParams 3119 3120 """ 3121 cdef PetscReal rval1 = 0, rval2 = 0, rval3 = 0 3122 try: CHKERR(TSAlphaGetParams(self.ts, &rval1, &rval2, &rval3)) 3123 except PetscError: pass 3124 if alpha_m is not None: rval1 = asReal(alpha_m) 3125 if alpha_f is not None: rval2 = asReal(alpha_f) 3126 if gamma is not None: rval3 = asReal(gamma) 3127 CHKERR(TSAlphaSetParams(self.ts, rval1, rval2, rval3)) 3128 3129 def getAlphaParams(self) -> tuple[float, float, float]: 3130 """Return the algorithmic parameters for `Type.ALPHA`. 3131 3132 Not collective. 3133 3134 See Also 3135 -------- 3136 petsc.TSAlphaGetParams 3137 3138 """ 3139 cdef PetscReal rval1 = 0, rval2 = 0, rval3 = 0 3140 CHKERR(TSAlphaGetParams(self.ts, &rval1, &rval2, &rval3)) 3141 return (toReal(rval1), toReal(rval2), toReal(rval3)) 3142 3143 # --- application context --- 3144 3145 property appctx: 3146 """Application context.""" 3147 def __get__(self) -> Any: 3148 return self.getAppCtx() 3149 3150 def __set__(self, value) -> None: 3151 self.setAppCtx(value) 3152 3153 # --- discretization space --- 3154 3155 property dm: 3156 """The `DM`.""" 3157 def __get__(self) -> DM: 3158 return self.getDM() 3159 3160 def __set__(self, value) -> None: 3161 self.setDM(value) 3162 3163 # --- xxx --- 3164 3165 property problem_type: 3166 """The problem type.""" 3167 def __get__(self) -> ProblemType: 3168 return self.getProblemType() 3169 3170 def __set__(self, value) -> None: 3171 self.setProblemType(value) 3172 3173 property equation_type: 3174 """The equation type.""" 3175 def __get__(self) -> EquationType: 3176 return self.getEquationType() 3177 3178 def __set__(self, value) -> None: 3179 self.setEquationType(value) 3180 3181 property snes: 3182 """The `SNES`.""" 3183 def __get__(self) -> SNES: 3184 return self.getSNES() 3185 3186 property ksp: 3187 """The `KSP`.""" 3188 def __get__(self) -> KSP: 3189 return self.getKSP() 3190 3191 property vec_sol: 3192 """The solution vector.""" 3193 def __get__(self) -> Vec: 3194 return self.getSolution() 3195 3196 # --- xxx --- 3197 3198 property time: 3199 """The current time.""" 3200 def __get__(self) -> float: 3201 return self.getTime() 3202 3203 def __set__(self, value) -> None: 3204 self.setTime(value) 3205 3206 property time_step: 3207 """The current time step size.""" 3208 def __get__(self) -> None: 3209 return self.getTimeStep() 3210 3211 def __set__(self, value): 3212 self.setTimeStep(value) 3213 3214 property step_number: 3215 """The current step number.""" 3216 def __get__(self) -> int: 3217 return self.getStepNumber() 3218 3219 def __set__(self, value) -> None: 3220 self.setStepNumber(value) 3221 3222 property max_time: 3223 """The maximum time.""" 3224 def __get__(self) -> float: 3225 return self.getMaxTime() 3226 3227 def __set__(self, value) -> None: 3228 self.setMaxTime(value) 3229 3230 property max_steps: 3231 """The maximum number of steps.""" 3232 def __get__(self) -> int: 3233 return self.getMaxSteps() 3234 3235 def __set__(self, value) -> None: 3236 self.setMaxSteps(value) 3237 3238 # --- convergence --- 3239 3240 property rtol: 3241 """The relative tolerance.""" 3242 def __get__(self) -> float: 3243 return self.getTolerances()[0] 3244 3245 def __set__(self, value) -> None: 3246 self.setTolerances(rtol=value) 3247 3248 property atol: 3249 """The absolute tolerance.""" 3250 def __get__(self) -> float: 3251 return self.getTolerances()[1] 3252 3253 def __set__(self, value) -> None: 3254 self.setTolerances(atol=value) 3255 3256 property reason: 3257 """The converged reason.""" 3258 def __get__(self) -> ConvergedReason: 3259 return self.getConvergedReason() 3260 3261 def __set__(self, value) -> None: 3262 self.setConvergedReason(value) 3263 3264 property iterating: 3265 """Indicates the `TS` is still iterating.""" 3266 def __get__(self) -> bool: 3267 return self.reason == 0 3268 3269 property converged: 3270 """Indicates the `TS` has converged.""" 3271 def __get__(self) -> bool: 3272 return self.reason > 0 3273 3274 property diverged: 3275 """Indicates the `TS` has stopped.""" 3276 def __get__(self) -> bool: 3277 return self.reason < 0 3278 3279# ----------------------------------------------------------------------------- 3280 3281del TSType 3282del TSRKType 3283del TSARKIMEXType 3284del TSDIRKType 3285del TSProblemType 3286del TSEquationType 3287del TSExactFinalTime 3288del TSConvergedReason 3289 3290# ----------------------------------------------------------------------------- 3291