1 /* 2 Code for timestepping with BDF methods 3 */ 4 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 5 #include <petscdm.h> 6 7 static PetscBool cited = PETSC_FALSE; 8 static const char citation[] = "@book{Brenan1995,\n" 9 " title = {Numerical Solution of Initial-Value Problems in Differential-Algebraic Equations},\n" 10 " author = {Brenan, K. and Campbell, S. and Petzold, L.},\n" 11 " publisher = {Society for Industrial and Applied Mathematics},\n" 12 " year = {1995},\n" 13 " doi = {10.1137/1.9781611971224},\n}\n"; 14 15 typedef struct { 16 PetscInt k, n; 17 PetscReal time[6 + 2]; 18 Vec work[6 + 2]; 19 Vec tvwork[6 + 2]; 20 PetscReal shift; 21 Vec vec_dot; /* Xdot when !transientvar, else Cdot where C(X) is the transient variable. */ 22 Vec vec_wrk; 23 Vec vec_lte; 24 25 PetscBool transientvar; 26 PetscInt order; 27 TSStepStatus status; 28 } TS_BDF; 29 30 /* Compute Lagrange polynomials on T[:n] evaluated at t. 31 * If one has data (T[i], Y[i]), then the interpolation/extrapolation is f(t) = \sum_i L[i]*Y[i]. 32 */ 33 static inline void LagrangeBasisVals(PetscInt n, PetscReal t, const PetscReal T[], PetscScalar L[]) { 34 PetscInt k, j; 35 for (k = 0; k < n; k++) 36 for (L[k] = 1, j = 0; j < n; j++) 37 if (j != k) L[k] *= (t - T[j]) / (T[k] - T[j]); 38 } 39 40 static inline void LagrangeBasisDers(PetscInt n, PetscReal t, const PetscReal T[], PetscScalar dL[]) { 41 PetscInt k, j, i; 42 for (k = 0; k < n; k++) 43 for (dL[k] = 0, j = 0; j < n; j++) 44 if (j != k) { 45 PetscReal L = 1 / (T[k] - T[j]); 46 for (i = 0; i < n; i++) 47 if (i != j && i != k) L *= (t - T[i]) / (T[k] - T[i]); 48 dL[k] += L; 49 } 50 } 51 52 static PetscErrorCode TSBDF_GetVecs(TS ts, DM dm, Vec *Xdot, Vec *Ydot) { 53 TS_BDF *bdf = (TS_BDF *)ts->data; 54 55 PetscFunctionBegin; 56 if (dm && dm != ts->dm) { 57 PetscCall(DMGetNamedGlobalVector(dm, "TSBDF_Vec_Xdot", Xdot)); 58 PetscCall(DMGetNamedGlobalVector(dm, "TSBDF_Vec_Ydot", Ydot)); 59 } else { 60 *Xdot = bdf->vec_dot; 61 *Ydot = bdf->vec_wrk; 62 } 63 PetscFunctionReturn(0); 64 } 65 66 static PetscErrorCode TSBDF_RestoreVecs(TS ts, DM dm, Vec *Xdot, Vec *Ydot) { 67 TS_BDF *bdf = (TS_BDF *)ts->data; 68 69 PetscFunctionBegin; 70 if (dm && dm != ts->dm) { 71 PetscCall(DMRestoreNamedGlobalVector(dm, "TSBDF_Vec_Xdot", Xdot)); 72 PetscCall(DMRestoreNamedGlobalVector(dm, "TSBDF_Vec_Ydot", Ydot)); 73 } else { 74 PetscCheck(*Xdot == bdf->vec_dot, PetscObjectComm((PetscObject)ts), PETSC_ERR_ARG_INCOMP, "Vec does not match the cache"); 75 PetscCheck(*Ydot == bdf->vec_wrk, PetscObjectComm((PetscObject)ts), PETSC_ERR_ARG_INCOMP, "Vec does not match the cache"); 76 *Xdot = NULL; 77 *Ydot = NULL; 78 } 79 PetscFunctionReturn(0); 80 } 81 82 static PetscErrorCode DMCoarsenHook_TSBDF(DM fine, DM coarse, void *ctx) { 83 PetscFunctionBegin; 84 PetscFunctionReturn(0); 85 } 86 87 static PetscErrorCode DMRestrictHook_TSBDF(DM fine, Mat restrct, Vec rscale, Mat inject, DM coarse, void *ctx) { 88 TS ts = (TS)ctx; 89 Vec Ydot, Ydot_c; 90 Vec Xdot, Xdot_c; 91 92 PetscFunctionBegin; 93 PetscCall(TSBDF_GetVecs(ts, fine, &Xdot, &Ydot)); 94 PetscCall(TSBDF_GetVecs(ts, coarse, &Xdot_c, &Ydot_c)); 95 96 PetscCall(MatRestrict(restrct, Ydot, Ydot_c)); 97 PetscCall(VecPointwiseMult(Ydot_c, rscale, Ydot_c)); 98 99 PetscCall(TSBDF_RestoreVecs(ts, fine, &Xdot, &Ydot)); 100 PetscCall(TSBDF_RestoreVecs(ts, coarse, &Xdot_c, &Ydot_c)); 101 PetscFunctionReturn(0); 102 } 103 104 static PetscErrorCode TSBDF_Advance(TS ts, PetscReal t, Vec X) { 105 TS_BDF *bdf = (TS_BDF *)ts->data; 106 PetscInt i, n = (PetscInt)(sizeof(bdf->work) / sizeof(Vec)); 107 Vec tail = bdf->work[n - 1], tvtail = bdf->tvwork[n - 1]; 108 109 PetscFunctionBegin; 110 for (i = n - 1; i >= 2; i--) { 111 bdf->time[i] = bdf->time[i - 1]; 112 bdf->work[i] = bdf->work[i - 1]; 113 bdf->tvwork[i] = bdf->tvwork[i - 1]; 114 } 115 bdf->n = PetscMin(bdf->n + 1, n - 1); 116 bdf->time[1] = t; 117 bdf->work[1] = tail; 118 bdf->tvwork[1] = tvtail; 119 PetscCall(VecCopy(X, tail)); 120 PetscCall(TSComputeTransientVariable(ts, tail, tvtail)); 121 PetscFunctionReturn(0); 122 } 123 124 static PetscErrorCode TSBDF_VecLTE(TS ts, PetscInt order, Vec lte) { 125 TS_BDF *bdf = (TS_BDF *)ts->data; 126 PetscInt i, n = order + 1; 127 PetscReal *time = bdf->time; 128 Vec *vecs = bdf->work; 129 PetscScalar a[8], b[8], alpha[8]; 130 131 PetscFunctionBegin; 132 LagrangeBasisDers(n + 0, time[0], time, a); 133 a[n] = 0; 134 LagrangeBasisDers(n + 1, time[0], time, b); 135 for (i = 0; i < n + 1; i++) alpha[i] = (a[i] - b[i]) / a[0]; 136 PetscCall(VecZeroEntries(lte)); 137 PetscCall(VecMAXPY(lte, n + 1, alpha, vecs)); 138 PetscFunctionReturn(0); 139 } 140 141 static PetscErrorCode TSBDF_Extrapolate(TS ts, PetscInt order, PetscReal t, Vec X) { 142 TS_BDF *bdf = (TS_BDF *)ts->data; 143 PetscInt n = order + 1; 144 PetscReal *time = bdf->time + 1; 145 Vec *vecs = bdf->work + 1; 146 PetscScalar alpha[7]; 147 148 PetscFunctionBegin; 149 n = PetscMin(n, bdf->n); 150 LagrangeBasisVals(n, t, time, alpha); 151 PetscCall(VecZeroEntries(X)); 152 PetscCall(VecMAXPY(X, n, alpha, vecs)); 153 PetscFunctionReturn(0); 154 } 155 156 static PetscErrorCode TSBDF_Interpolate(TS ts, PetscInt order, PetscReal t, Vec X) { 157 TS_BDF *bdf = (TS_BDF *)ts->data; 158 PetscInt n = order + 1; 159 PetscReal *time = bdf->time; 160 Vec *vecs = bdf->work; 161 PetscScalar alpha[7]; 162 163 PetscFunctionBegin; 164 LagrangeBasisVals(n, t, time, alpha); 165 PetscCall(VecZeroEntries(X)); 166 PetscCall(VecMAXPY(X, n, alpha, vecs)); 167 PetscFunctionReturn(0); 168 } 169 170 /* Compute the affine term V0 such that Xdot = shift*X + V0. 171 * 172 * When using transient variables, we're computing Cdot = shift*C(X) + V0, and thus choose a linear combination of tvwork. 173 */ 174 static PetscErrorCode TSBDF_PreSolve(TS ts) { 175 TS_BDF *bdf = (TS_BDF *)ts->data; 176 PetscInt i, n = PetscMax(bdf->k, 1) + 1; 177 Vec V, V0; 178 Vec vecs[7]; 179 PetscScalar alpha[7]; 180 181 PetscFunctionBegin; 182 PetscCall(TSBDF_GetVecs(ts, NULL, &V, &V0)); 183 LagrangeBasisDers(n, bdf->time[0], bdf->time, alpha); 184 for (i = 1; i < n; i++) { vecs[i] = bdf->transientvar ? bdf->tvwork[i] : bdf->work[i]; } 185 PetscCall(VecZeroEntries(V0)); 186 PetscCall(VecMAXPY(V0, n - 1, alpha + 1, vecs + 1)); 187 bdf->shift = PetscRealPart(alpha[0]); 188 PetscCall(TSBDF_RestoreVecs(ts, NULL, &V, &V0)); 189 PetscFunctionReturn(0); 190 } 191 192 static PetscErrorCode TSBDF_SNESSolve(TS ts, Vec b, Vec x) { 193 PetscInt nits, lits; 194 195 PetscFunctionBegin; 196 PetscCall(TSBDF_PreSolve(ts)); 197 PetscCall(SNESSolve(ts->snes, b, x)); 198 PetscCall(SNESGetIterationNumber(ts->snes, &nits)); 199 PetscCall(SNESGetLinearSolveIterations(ts->snes, &lits)); 200 ts->snes_its += nits; 201 ts->ksp_its += lits; 202 PetscFunctionReturn(0); 203 } 204 205 static PetscErrorCode TSBDF_Restart(TS ts, PetscBool *accept) { 206 TS_BDF *bdf = (TS_BDF *)ts->data; 207 208 PetscFunctionBegin; 209 bdf->k = 1; 210 bdf->n = 0; 211 PetscCall(TSBDF_Advance(ts, ts->ptime, ts->vec_sol)); 212 213 bdf->time[0] = ts->ptime + ts->time_step / 2; 214 PetscCall(VecCopy(bdf->work[1], bdf->work[0])); 215 PetscCall(TSPreStage(ts, bdf->time[0])); 216 PetscCall(TSBDF_SNESSolve(ts, NULL, bdf->work[0])); 217 PetscCall(TSPostStage(ts, bdf->time[0], 0, &bdf->work[0])); 218 PetscCall(TSAdaptCheckStage(ts->adapt, ts, bdf->time[0], bdf->work[0], accept)); 219 if (!*accept) PetscFunctionReturn(0); 220 221 bdf->k = PetscMin(2, bdf->order); 222 bdf->n++; 223 PetscCall(VecCopy(bdf->work[0], bdf->work[2])); 224 bdf->time[2] = bdf->time[0]; 225 PetscCall(TSComputeTransientVariable(ts, bdf->work[2], bdf->tvwork[2])); 226 PetscFunctionReturn(0); 227 } 228 229 static const char *const BDF_SchemeName[] = {"", "1", "2", "3", "4", "5", "6"}; 230 231 static PetscErrorCode TSStep_BDF(TS ts) { 232 TS_BDF *bdf = (TS_BDF *)ts->data; 233 PetscInt rejections = 0; 234 PetscBool stageok, accept = PETSC_TRUE; 235 PetscReal next_time_step = ts->time_step; 236 237 PetscFunctionBegin; 238 PetscCall(PetscCitationsRegister(citation, &cited)); 239 240 if (!ts->steprollback && !ts->steprestart) { 241 bdf->k = PetscMin(bdf->k + 1, bdf->order); 242 PetscCall(TSBDF_Advance(ts, ts->ptime, ts->vec_sol)); 243 } 244 245 bdf->status = TS_STEP_INCOMPLETE; 246 while (!ts->reason && bdf->status != TS_STEP_COMPLETE) { 247 if (ts->steprestart) { 248 PetscCall(TSBDF_Restart(ts, &stageok)); 249 if (!stageok) goto reject_step; 250 } 251 252 bdf->time[0] = ts->ptime + ts->time_step; 253 PetscCall(TSBDF_Extrapolate(ts, bdf->k - (accept ? 0 : 1), bdf->time[0], bdf->work[0])); 254 PetscCall(TSPreStage(ts, bdf->time[0])); 255 PetscCall(TSBDF_SNESSolve(ts, NULL, bdf->work[0])); 256 PetscCall(TSPostStage(ts, bdf->time[0], 0, &bdf->work[0])); 257 PetscCall(TSAdaptCheckStage(ts->adapt, ts, bdf->time[0], bdf->work[0], &stageok)); 258 if (!stageok) goto reject_step; 259 260 bdf->status = TS_STEP_PENDING; 261 PetscCall(TSAdaptCandidatesClear(ts->adapt)); 262 PetscCall(TSAdaptCandidateAdd(ts->adapt, BDF_SchemeName[bdf->k], bdf->k, 1, 1.0, 1.0, PETSC_TRUE)); 263 PetscCall(TSAdaptChoose(ts->adapt, ts, ts->time_step, NULL, &next_time_step, &accept)); 264 bdf->status = accept ? TS_STEP_COMPLETE : TS_STEP_INCOMPLETE; 265 if (!accept) { 266 ts->time_step = next_time_step; 267 goto reject_step; 268 } 269 270 PetscCall(VecCopy(bdf->work[0], ts->vec_sol)); 271 ts->ptime += ts->time_step; 272 ts->time_step = next_time_step; 273 break; 274 275 reject_step: 276 ts->reject++; 277 accept = PETSC_FALSE; 278 if (!ts->reason && ++rejections > ts->max_reject && ts->max_reject >= 0) { 279 PetscCall(PetscInfo(ts, "Step=%" PetscInt_FMT ", step rejections %" PetscInt_FMT " greater than current TS allowed, stopping solve\n", ts->steps, rejections)); 280 ts->reason = TS_DIVERGED_STEP_REJECTED; 281 } 282 } 283 PetscFunctionReturn(0); 284 } 285 286 static PetscErrorCode TSInterpolate_BDF(TS ts, PetscReal t, Vec X) { 287 TS_BDF *bdf = (TS_BDF *)ts->data; 288 289 PetscFunctionBegin; 290 PetscCall(TSBDF_Interpolate(ts, bdf->k, t, X)); 291 PetscFunctionReturn(0); 292 } 293 294 static PetscErrorCode TSEvaluateWLTE_BDF(TS ts, NormType wnormtype, PetscInt *order, PetscReal *wlte) { 295 TS_BDF *bdf = (TS_BDF *)ts->data; 296 PetscInt k = bdf->k; 297 PetscReal wltea, wlter; 298 Vec X = bdf->work[0], Y = bdf->vec_lte; 299 300 PetscFunctionBegin; 301 k = PetscMin(k, bdf->n - 1); 302 PetscCall(TSBDF_VecLTE(ts, k, Y)); 303 PetscCall(VecAXPY(Y, 1, X)); 304 PetscCall(TSErrorWeightedNorm(ts, X, Y, wnormtype, wlte, &wltea, &wlter)); 305 if (order) *order = k + 1; 306 PetscFunctionReturn(0); 307 } 308 309 static PetscErrorCode TSRollBack_BDF(TS ts) { 310 TS_BDF *bdf = (TS_BDF *)ts->data; 311 312 PetscFunctionBegin; 313 PetscCall(VecCopy(bdf->work[1], ts->vec_sol)); 314 PetscFunctionReturn(0); 315 } 316 317 static PetscErrorCode SNESTSFormFunction_BDF(SNES snes, Vec X, Vec F, TS ts) { 318 TS_BDF *bdf = (TS_BDF *)ts->data; 319 DM dm, dmsave = ts->dm; 320 PetscReal t = bdf->time[0]; 321 PetscReal shift = bdf->shift; 322 Vec V, V0; 323 324 PetscFunctionBegin; 325 PetscCall(SNESGetDM(snes, &dm)); 326 PetscCall(TSBDF_GetVecs(ts, dm, &V, &V0)); 327 if (bdf->transientvar) { /* shift*C(X) + V0 */ 328 PetscCall(TSComputeTransientVariable(ts, X, V)); 329 PetscCall(VecAYPX(V, shift, V0)); 330 } else { /* shift*X + V0 */ 331 PetscCall(VecWAXPY(V, shift, X, V0)); 332 } 333 334 /* F = Function(t,X,V) */ 335 ts->dm = dm; 336 PetscCall(TSComputeIFunction(ts, t, X, V, F, PETSC_FALSE)); 337 ts->dm = dmsave; 338 339 PetscCall(TSBDF_RestoreVecs(ts, dm, &V, &V0)); 340 PetscFunctionReturn(0); 341 } 342 343 static PetscErrorCode SNESTSFormJacobian_BDF(SNES snes, Vec X, Mat J, Mat P, TS ts) { 344 TS_BDF *bdf = (TS_BDF *)ts->data; 345 DM dm, dmsave = ts->dm; 346 PetscReal t = bdf->time[0]; 347 PetscReal shift = bdf->shift; 348 Vec V, V0; 349 350 PetscFunctionBegin; 351 PetscCall(SNESGetDM(snes, &dm)); 352 PetscCall(TSBDF_GetVecs(ts, dm, &V, &V0)); 353 354 /* J,P = Jacobian(t,X,V) */ 355 ts->dm = dm; 356 PetscCall(TSComputeIJacobian(ts, t, X, V, shift, J, P, PETSC_FALSE)); 357 ts->dm = dmsave; 358 359 PetscCall(TSBDF_RestoreVecs(ts, dm, &V, &V0)); 360 PetscFunctionReturn(0); 361 } 362 363 static PetscErrorCode TSReset_BDF(TS ts) { 364 TS_BDF *bdf = (TS_BDF *)ts->data; 365 size_t i, n = sizeof(bdf->work) / sizeof(Vec); 366 367 PetscFunctionBegin; 368 bdf->k = bdf->n = 0; 369 for (i = 0; i < n; i++) { 370 PetscCall(VecDestroy(&bdf->work[i])); 371 PetscCall(VecDestroy(&bdf->tvwork[i])); 372 } 373 PetscCall(VecDestroy(&bdf->vec_dot)); 374 PetscCall(VecDestroy(&bdf->vec_wrk)); 375 PetscCall(VecDestroy(&bdf->vec_lte)); 376 if (ts->dm) PetscCall(DMCoarsenHookRemove(ts->dm, DMCoarsenHook_TSBDF, DMRestrictHook_TSBDF, ts)); 377 PetscFunctionReturn(0); 378 } 379 380 static PetscErrorCode TSDestroy_BDF(TS ts) { 381 PetscFunctionBegin; 382 PetscCall(TSReset_BDF(ts)); 383 PetscCall(PetscFree(ts->data)); 384 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSBDFSetOrder_C", NULL)); 385 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSBDFGetOrder_C", NULL)); 386 PetscFunctionReturn(0); 387 } 388 389 static PetscErrorCode TSSetUp_BDF(TS ts) { 390 TS_BDF *bdf = (TS_BDF *)ts->data; 391 size_t i, n = sizeof(bdf->work) / sizeof(Vec); 392 PetscReal low, high, two = 2; 393 394 PetscFunctionBegin; 395 PetscCall(TSHasTransientVariable(ts, &bdf->transientvar)); 396 bdf->k = bdf->n = 0; 397 for (i = 0; i < n; i++) { 398 PetscCall(VecDuplicate(ts->vec_sol, &bdf->work[i])); 399 if (i && bdf->transientvar) PetscCall(VecDuplicate(ts->vec_sol, &bdf->tvwork[i])); 400 } 401 PetscCall(VecDuplicate(ts->vec_sol, &bdf->vec_dot)); 402 PetscCall(VecDuplicate(ts->vec_sol, &bdf->vec_wrk)); 403 PetscCall(VecDuplicate(ts->vec_sol, &bdf->vec_lte)); 404 PetscCall(TSGetDM(ts, &ts->dm)); 405 PetscCall(DMCoarsenHookAdd(ts->dm, DMCoarsenHook_TSBDF, DMRestrictHook_TSBDF, ts)); 406 407 PetscCall(TSGetAdapt(ts, &ts->adapt)); 408 PetscCall(TSAdaptCandidatesClear(ts->adapt)); 409 PetscCall(TSAdaptGetClip(ts->adapt, &low, &high)); 410 PetscCall(TSAdaptSetClip(ts->adapt, low, PetscMin(high, two))); 411 412 PetscCall(TSGetSNES(ts, &ts->snes)); 413 PetscFunctionReturn(0); 414 } 415 416 static PetscErrorCode TSSetFromOptions_BDF(TS ts, PetscOptionItems *PetscOptionsObject) { 417 PetscFunctionBegin; 418 PetscOptionsHeadBegin(PetscOptionsObject, "BDF ODE solver options"); 419 { 420 PetscBool flg; 421 PetscInt order; 422 PetscCall(TSBDFGetOrder(ts, &order)); 423 PetscCall(PetscOptionsInt("-ts_bdf_order", "Order of the BDF method", "TSBDFSetOrder", order, &order, &flg)); 424 if (flg) PetscCall(TSBDFSetOrder(ts, order)); 425 } 426 PetscOptionsHeadEnd(); 427 PetscFunctionReturn(0); 428 } 429 430 static PetscErrorCode TSView_BDF(TS ts, PetscViewer viewer) { 431 TS_BDF *bdf = (TS_BDF *)ts->data; 432 PetscBool iascii; 433 434 PetscFunctionBegin; 435 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 436 if (iascii) PetscCall(PetscViewerASCIIPrintf(viewer, " Order=%" PetscInt_FMT "\n", bdf->order)); 437 PetscFunctionReturn(0); 438 } 439 440 /* ------------------------------------------------------------ */ 441 442 static PetscErrorCode TSBDFSetOrder_BDF(TS ts, PetscInt order) { 443 TS_BDF *bdf = (TS_BDF *)ts->data; 444 445 PetscFunctionBegin; 446 if (order == bdf->order) PetscFunctionReturn(0); 447 PetscCheck(order >= 1 && order <= 6, PetscObjectComm((PetscObject)ts), PETSC_ERR_ARG_OUTOFRANGE, "BDF Order %" PetscInt_FMT " not implemented", order); 448 bdf->order = order; 449 PetscFunctionReturn(0); 450 } 451 452 static PetscErrorCode TSBDFGetOrder_BDF(TS ts, PetscInt *order) { 453 TS_BDF *bdf = (TS_BDF *)ts->data; 454 455 PetscFunctionBegin; 456 *order = bdf->order; 457 PetscFunctionReturn(0); 458 } 459 460 /* ------------------------------------------------------------ */ 461 462 /*MC 463 TSBDF - DAE solver using BDF methods 464 465 Level: beginner 466 467 .seealso: `TS`, `TSCreate()`, `TSSetType()` 468 M*/ 469 PETSC_EXTERN PetscErrorCode TSCreate_BDF(TS ts) { 470 TS_BDF *bdf; 471 472 PetscFunctionBegin; 473 ts->ops->reset = TSReset_BDF; 474 ts->ops->destroy = TSDestroy_BDF; 475 ts->ops->view = TSView_BDF; 476 ts->ops->setup = TSSetUp_BDF; 477 ts->ops->setfromoptions = TSSetFromOptions_BDF; 478 ts->ops->step = TSStep_BDF; 479 ts->ops->evaluatewlte = TSEvaluateWLTE_BDF; 480 ts->ops->rollback = TSRollBack_BDF; 481 ts->ops->interpolate = TSInterpolate_BDF; 482 ts->ops->snesfunction = SNESTSFormFunction_BDF; 483 ts->ops->snesjacobian = SNESTSFormJacobian_BDF; 484 ts->default_adapt_type = TSADAPTBASIC; 485 486 ts->usessnes = PETSC_TRUE; 487 488 PetscCall(PetscNewLog(ts, &bdf)); 489 ts->data = (void *)bdf; 490 491 bdf->status = TS_STEP_COMPLETE; 492 493 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSBDFSetOrder_C", TSBDFSetOrder_BDF)); 494 PetscCall(PetscObjectComposeFunction((PetscObject)ts, "TSBDFGetOrder_C", TSBDFGetOrder_BDF)); 495 PetscCall(TSBDFSetOrder(ts, 2)); 496 PetscFunctionReturn(0); 497 } 498 499 /* ------------------------------------------------------------ */ 500 501 /*@ 502 TSBDFSetOrder - Set the order of the BDF method 503 504 Logically Collective on TS 505 506 Input Parameters: 507 + ts - timestepping context 508 - order - order of the method 509 510 Options Database: 511 . -ts_bdf_order <order> - select the order 512 513 Level: intermediate 514 515 @*/ 516 PetscErrorCode TSBDFSetOrder(TS ts, PetscInt order) { 517 PetscFunctionBegin; 518 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 519 PetscValidLogicalCollectiveInt(ts, order, 2); 520 PetscTryMethod(ts, "TSBDFSetOrder_C", (TS, PetscInt), (ts, order)); 521 PetscFunctionReturn(0); 522 } 523 524 /*@ 525 TSBDFGetOrder - Get the order of the BDF method 526 527 Not Collective 528 529 Input Parameter: 530 . ts - timestepping context 531 532 Output Parameter: 533 . order - order of the method 534 535 Level: intermediate 536 537 @*/ 538 PetscErrorCode TSBDFGetOrder(TS ts, PetscInt *order) { 539 PetscFunctionBegin; 540 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 541 PetscValidIntPointer(order, 2); 542 PetscUseMethod(ts, "TSBDFGetOrder_C", (TS, PetscInt *), (ts, order)); 543 PetscFunctionReturn(0); 544 } 545