1 /* 2 Code for timestepping with BDF methods 3 */ 4 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 5 6 static PetscBool cited = PETSC_FALSE; 7 static const char citation[] = 8 "@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 PetscReal shift; 20 Vec vec_dot; 21 Vec vec_lte; 22 23 PetscInt order; 24 PetscBool adapt; 25 TSStepStatus status; 26 } TS_BDF; 27 28 29 PETSC_STATIC_INLINE void LagrangeBasisVals(PetscInt n,PetscReal t,const PetscReal T[],PetscScalar L[]) 30 { 31 PetscInt k,j; 32 for (k=0; k<n; k++) 33 for (L[k]=1, j=0; j<n; j++) 34 if (j != k) 35 L[k] *= (t - T[j])/(T[k] - T[j]); 36 } 37 38 PETSC_STATIC_INLINE void LagrangeBasisDers(PetscInt n,PetscReal t,const PetscReal T[],PetscScalar dL[]) 39 { 40 PetscInt k,j,i; 41 for (k=0; k<n; k++) 42 for (dL[k]=0, j=0; j<n; j++) 43 if (j != k) { 44 PetscReal L = 1/(T[k] - T[j]); 45 for (i=0; i<n; i++) 46 if (i != j && i != k) 47 L *= (t - T[i])/(T[k] - T[i]); 48 dL[k] += L; 49 } 50 } 51 52 #undef __FUNCT__ 53 #define __FUNCT__ "TSBDF_Advance" 54 static PetscErrorCode TSBDF_Advance(TS ts,PetscReal t,Vec X) 55 { 56 TS_BDF *bdf = (TS_BDF*)ts->data; 57 PetscInt i,n = (PetscInt)(sizeof(bdf->work)/sizeof(Vec)); 58 Vec tail = bdf->work[n-1]; 59 PetscErrorCode ierr; 60 61 PetscFunctionBegin; 62 for (i=n-1; i>=2; i--) { 63 bdf->time[i] = bdf->time[i-1]; 64 bdf->work[i] = bdf->work[i-1]; 65 } 66 bdf->n = PetscMin(bdf->n+1,n-1); 67 bdf->time[1] = t; 68 bdf->work[1] = tail; 69 ierr = VecCopy(X,tail);CHKERRQ(ierr); 70 PetscFunctionReturn(0); 71 } 72 73 #undef __FUNCT__ 74 #define __FUNCT__ "TSBDF_VecDot" 75 static PetscErrorCode TSBDF_VecDot(TS ts,PetscInt order,PetscReal t,Vec X,Vec Xdot,PetscReal *shift) 76 { 77 TS_BDF *bdf = (TS_BDF*)ts->data; 78 PetscInt i,n = order+1; 79 PetscReal time[7]; 80 Vec vecs[7]; 81 PetscScalar alpha[7]; 82 PetscErrorCode ierr; 83 84 PetscFunctionBegin; 85 n = PetscMax(n,2); 86 time[0] = t; for (i=1; i<n; i++) time[i] = bdf->time[i]; 87 vecs[0] = X; for (i=1; i<n; i++) vecs[i] = bdf->work[i]; 88 LagrangeBasisDers(n,t,time,alpha); 89 ierr = VecZeroEntries(Xdot);CHKERRQ(ierr); 90 ierr = VecMAXPY(Xdot,n,alpha,vecs);CHKERRQ(ierr); 91 if (shift) *shift = PetscRealPart(alpha[0]); 92 PetscFunctionReturn(0); 93 } 94 95 #undef __FUNCT__ 96 #define __FUNCT__ "TSBDF_VecLTE" 97 static PetscErrorCode TSBDF_VecLTE(TS ts,PetscInt order,Vec lte) 98 { 99 TS_BDF *bdf = (TS_BDF*)ts->data; 100 PetscInt i,n = order+1; 101 PetscReal *time = bdf->time; 102 Vec *vecs = bdf->work; 103 PetscScalar a[8],b[8],alpha[8]; 104 PetscErrorCode ierr; 105 106 PetscFunctionBegin; 107 LagrangeBasisDers(n+0,time[0],time,a); a[n] =0; 108 LagrangeBasisDers(n+1,time[0],time,b); 109 for (i=0; i<n+1; i++) alpha[i] = (a[i]-b[i])/a[0]; 110 ierr = VecZeroEntries(lte);CHKERRQ(ierr); 111 ierr = VecMAXPY(lte,n+1,alpha,vecs);CHKERRQ(ierr); 112 PetscFunctionReturn(0); 113 } 114 115 #undef __FUNCT__ 116 #define __FUNCT__ "TSBDF_Extrapolate" 117 static PetscErrorCode TSBDF_Extrapolate(TS ts,PetscInt order,PetscReal t,Vec X) 118 { 119 TS_BDF *bdf = (TS_BDF*)ts->data; 120 PetscInt n = order+1; 121 PetscReal *time = bdf->time+1; 122 Vec *vecs = bdf->work+1; 123 PetscScalar alpha[7]; 124 PetscErrorCode ierr; 125 126 PetscFunctionBegin; 127 n = PetscMin(n,bdf->n); 128 LagrangeBasisVals(n,t,time,alpha); 129 ierr = VecZeroEntries(X);CHKERRQ(ierr); 130 ierr = VecMAXPY(X,n,alpha,vecs);CHKERRQ(ierr); 131 PetscFunctionReturn(0); 132 } 133 134 #undef __FUNCT__ 135 #define __FUNCT__ "TSBDF_Interpolate" 136 static PetscErrorCode TSBDF_Interpolate(TS ts,PetscInt order,PetscReal t,Vec X) 137 { 138 TS_BDF *bdf = (TS_BDF*)ts->data; 139 PetscInt n = order+1; 140 PetscReal *time = bdf->time; 141 Vec *vecs = bdf->work; 142 PetscScalar alpha[7]; 143 PetscErrorCode ierr; 144 145 PetscFunctionBegin; 146 LagrangeBasisVals(n,t,time,alpha); 147 ierr = VecZeroEntries(X);CHKERRQ(ierr); 148 ierr = VecMAXPY(X,n,alpha,vecs);CHKERRQ(ierr); 149 PetscFunctionReturn(0); 150 } 151 152 #undef __FUNCT__ 153 #define __FUNCT__ "TS_SNESSolve" 154 static PetscErrorCode TS_SNESSolve(TS ts,Vec b,Vec x) 155 { 156 PetscInt nits,lits; 157 PetscErrorCode ierr; 158 159 PetscFunctionBegin; 160 ierr = SNESSolve(ts->snes,b,x);CHKERRQ(ierr); 161 ierr = SNESGetIterationNumber(ts->snes,&nits);CHKERRQ(ierr); 162 ierr = SNESGetLinearSolveIterations(ts->snes,&lits);CHKERRQ(ierr); 163 ts->snes_its += nits; ts->ksp_its += lits; 164 PetscFunctionReturn(0); 165 } 166 167 #undef __FUNCT__ 168 #define __FUNCT__ "TSBDF_Restart" 169 static PetscErrorCode TSBDF_Restart(TS ts,PetscBool *accept) 170 { 171 TS_BDF *bdf = (TS_BDF*)ts->data; 172 PetscErrorCode ierr; 173 174 PetscFunctionBegin; 175 bdf->k = 1; bdf->n = 0; 176 ierr = TSBDF_Advance(ts,ts->ptime,ts->vec_sol);CHKERRQ(ierr); 177 178 bdf->time[0] = ts->ptime + ts->time_step/2; 179 ierr = VecCopy(bdf->work[1],bdf->work[0]);CHKERRQ(ierr); 180 ierr = TSPreStage(ts,bdf->time[0]);CHKERRQ(ierr); 181 ierr = TS_SNESSolve(ts,NULL,bdf->work[0]);CHKERRQ(ierr); 182 ierr = TSPostStage(ts,bdf->time[0],0,&bdf->work[0]);CHKERRQ(ierr); 183 ierr = TSAdaptCheckStage(ts->adapt,ts,bdf->time[0],bdf->work[0],accept);CHKERRQ(ierr); 184 if (!*accept) PetscFunctionReturn(0); 185 186 bdf->k = PetscMin(2,bdf->order); bdf->n++; 187 ierr = VecCopy(bdf->work[0],bdf->work[2]);CHKERRQ(ierr); 188 bdf->time[2] = bdf->time[0]; 189 PetscFunctionReturn(0); 190 } 191 192 static const char *const BDF_SchemeName[] = {"", "1", "2", "3", "4", "5", "6"}; 193 194 #undef __FUNCT__ 195 #define __FUNCT__ "TSStep_BDF" 196 static PetscErrorCode TSStep_BDF(TS ts) 197 { 198 TS_BDF *bdf = (TS_BDF*)ts->data; 199 PetscInt rejections = 0; 200 PetscBool stageok,accept = PETSC_TRUE; 201 PetscReal next_time_step = ts->time_step; 202 PetscErrorCode ierr; 203 204 PetscFunctionBegin; 205 ierr = PetscCitationsRegister(citation,&cited);CHKERRQ(ierr); 206 207 if (!ts->steprollback && !ts->steprestart) { 208 bdf->k = PetscMin(bdf->k+1,bdf->order); 209 ierr = TSBDF_Advance(ts,ts->ptime,ts->vec_sol);CHKERRQ(ierr); 210 } 211 212 bdf->status = TS_STEP_INCOMPLETE; 213 while (!ts->reason && bdf->status != TS_STEP_COMPLETE) { 214 215 if (ts->steprestart) { 216 ierr = TSBDF_Restart(ts,&stageok);CHKERRQ(ierr); 217 if (!stageok) goto reject_step; 218 } 219 220 bdf->time[0] = ts->ptime + ts->time_step; 221 ierr = TSBDF_Extrapolate(ts,bdf->k-(accept?0:1),bdf->time[0],bdf->work[0]);CHKERRQ(ierr); 222 ierr = TSPreStage(ts,bdf->time[0]);CHKERRQ(ierr); 223 ierr = TS_SNESSolve(ts,NULL,bdf->work[0]);CHKERRQ(ierr); 224 ierr = TSPostStage(ts,bdf->time[0],0,&bdf->work[0]);CHKERRQ(ierr); 225 ierr = TSAdaptCheckStage(ts->adapt,ts,bdf->time[0],bdf->work[0],&stageok);CHKERRQ(ierr); 226 if (!stageok) goto reject_step; 227 228 bdf->status = TS_STEP_PENDING; 229 ierr = TSAdaptCandidatesClear(ts->adapt);CHKERRQ(ierr); 230 ierr = TSAdaptCandidateAdd(ts->adapt,BDF_SchemeName[bdf->k],bdf->k,1,1.0,1.0,PETSC_TRUE);CHKERRQ(ierr); 231 ierr = TSAdaptChoose(ts->adapt,ts,ts->time_step,NULL,&next_time_step,&accept);CHKERRQ(ierr); 232 bdf->status = accept ? TS_STEP_COMPLETE : TS_STEP_INCOMPLETE; 233 if (!accept) { ts->time_step = next_time_step; goto reject_step; } 234 235 ierr = VecCopy(bdf->work[0],ts->vec_sol);CHKERRQ(ierr); 236 ts->ptime += ts->time_step; 237 ts->time_step = next_time_step; 238 break; 239 240 reject_step: 241 ts->reject++; accept = PETSC_FALSE; 242 if (!ts->reason && ++rejections > ts->max_reject && ts->max_reject >= 0) { 243 ierr = PetscInfo2(ts,"Step=%D, step rejections %D greater than current TS allowed, stopping solve\n",ts->steps,rejections);CHKERRQ(ierr); 244 ts->reason = TS_DIVERGED_STEP_REJECTED; 245 } 246 } 247 PetscFunctionReturn(0); 248 } 249 250 #undef __FUNCT__ 251 #define __FUNCT__ "TSInterpolate_BDF" 252 static PetscErrorCode TSInterpolate_BDF(TS ts,PetscReal t,Vec X) 253 { 254 TS_BDF *bdf = (TS_BDF*)ts->data; 255 PetscErrorCode ierr; 256 257 PetscFunctionBegin; 258 ierr = TSBDF_Interpolate(ts,bdf->k,t,X);CHKERRQ(ierr); 259 PetscFunctionReturn(0); 260 } 261 262 #undef __FUNCT__ 263 #define __FUNCT__ "TSEvaluateWLTE_BDF" 264 static PetscErrorCode TSEvaluateWLTE_BDF(TS ts,NormType wnormtype,PetscInt *order,PetscReal *wlte) 265 { 266 TS_BDF *bdf = (TS_BDF*)ts->data; 267 PetscInt k = bdf->k; 268 Vec X = bdf->work[0], Y = bdf->vec_lte; 269 PetscErrorCode ierr; 270 271 PetscFunctionBegin; 272 k = PetscMin(k,bdf->n-1); 273 ierr = TSBDF_VecLTE(ts,k,Y);CHKERRQ(ierr); 274 ierr = VecAXPY(Y,1,X);CHKERRQ(ierr); 275 ierr = TSErrorWeightedNorm(ts,X,Y,wnormtype,wlte);CHKERRQ(ierr); 276 if (order) *order = k + 1; 277 PetscFunctionReturn(0); 278 } 279 280 #undef __FUNCT__ 281 #define __FUNCT__ "TSRollBack_BDF" 282 static PetscErrorCode TSRollBack_BDF(TS ts) 283 { 284 TS_BDF *bdf = (TS_BDF*)ts->data; 285 PetscErrorCode ierr; 286 287 PetscFunctionBegin; 288 ierr = VecCopy(bdf->work[1],ts->vec_sol);CHKERRQ(ierr); 289 PetscFunctionReturn(0); 290 } 291 292 #undef __FUNCT__ 293 #define __FUNCT__ "SNESTSFormFunction_BDF" 294 static PetscErrorCode SNESTSFormFunction_BDF(PETSC_UNUSED SNES snes,Vec X,Vec F,TS ts) 295 { 296 TS_BDF *bdf = (TS_BDF*)ts->data; 297 PetscInt k = bdf->k; 298 PetscReal t = bdf->time[0]; 299 Vec V = bdf->vec_dot; 300 PetscErrorCode ierr; 301 302 PetscFunctionBegin; 303 ierr = TSBDF_VecDot(ts,k,t,X,V,&bdf->shift);CHKERRQ(ierr); 304 /* F = Function(t,X,V) */ 305 ierr = TSComputeIFunction(ts,t,X,V,F,PETSC_FALSE);CHKERRQ(ierr); 306 PetscFunctionReturn(0); 307 } 308 309 #undef __FUNCT__ 310 #define __FUNCT__ "SNESTSFormJacobian_BDF" 311 static PetscErrorCode SNESTSFormJacobian_BDF(PETSC_UNUSED SNES snes, 312 PETSC_UNUSED Vec X, 313 Mat J,Mat P, 314 TS ts) 315 { 316 TS_BDF *bdf = (TS_BDF*)ts->data; 317 PetscReal t = bdf->time[0]; 318 Vec V = bdf->vec_dot; 319 PetscReal dVdX = bdf->shift; 320 PetscErrorCode ierr; 321 322 PetscFunctionBegin; 323 /* J,P = Jacobian(t,X,V) */ 324 ierr = TSComputeIJacobian(ts,t,X,V,dVdX,J,P,PETSC_FALSE);CHKERRQ(ierr); 325 PetscFunctionReturn(0); 326 } 327 328 #undef __FUNCT__ 329 #define __FUNCT__ "TSReset_BDF" 330 static PetscErrorCode TSReset_BDF(TS ts) 331 { 332 TS_BDF *bdf = (TS_BDF*)ts->data; 333 size_t i,n = sizeof(bdf->work)/sizeof(Vec); 334 PetscErrorCode ierr; 335 336 PetscFunctionBegin; 337 for (i=0; i<n; i++) {ierr = VecDestroy(&bdf->work[i]);CHKERRQ(ierr);} 338 ierr = VecDestroy(&bdf->vec_dot);CHKERRQ(ierr); 339 ierr = VecDestroy(&bdf->vec_lte);CHKERRQ(ierr); 340 PetscFunctionReturn(0); 341 } 342 343 #undef __FUNCT__ 344 #define __FUNCT__ "TSDestroy_BDF" 345 static PetscErrorCode TSDestroy_BDF(TS ts) 346 { 347 PetscErrorCode ierr; 348 349 PetscFunctionBegin; 350 ierr = TSReset_BDF(ts);CHKERRQ(ierr); 351 ierr = PetscFree(ts->data);CHKERRQ(ierr); 352 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSBDFSetOrder_C",NULL);CHKERRQ(ierr); 353 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSBDFGetOrder_C",NULL);CHKERRQ(ierr); 354 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSBDFUseAdapt_C",NULL);CHKERRQ(ierr); 355 PetscFunctionReturn(0); 356 } 357 358 #undef __FUNCT__ 359 #define __FUNCT__ "TSSetUp_BDF" 360 static PetscErrorCode TSSetUp_BDF(TS ts) 361 { 362 TS_BDF *bdf = (TS_BDF*)ts->data; 363 size_t i,n = sizeof(bdf->work)/sizeof(Vec); 364 PetscErrorCode ierr; 365 366 PetscFunctionBegin; 367 bdf->k = bdf->n = 0; 368 for (i=0; i<n; i++) {ierr = VecDuplicate(ts->vec_sol,&bdf->work[i]);CHKERRQ(ierr);} 369 ierr = VecDuplicate(ts->vec_sol,&bdf->vec_dot);CHKERRQ(ierr); 370 ierr = VecDuplicate(ts->vec_sol,&bdf->vec_lte);CHKERRQ(ierr); 371 372 ierr = TSGetAdapt(ts,&ts->adapt);CHKERRQ(ierr); 373 ierr = TSAdaptCandidatesClear(ts->adapt);CHKERRQ(ierr); 374 if (!bdf->adapt) { 375 ierr = TSAdaptSetType(ts->adapt,TSADAPTNONE);CHKERRQ(ierr); 376 } else { 377 PetscReal low,high; 378 ierr = TSAdaptBasicGetClip(ts->adapt,&low,&high);CHKERRQ(ierr); 379 high = PetscMin(high,2.0); 380 ierr = TSAdaptBasicSetClip(ts->adapt,low,high);CHKERRQ(ierr); 381 if (ts->exact_final_time == TS_EXACTFINALTIME_UNSPECIFIED) 382 ts->exact_final_time = TS_EXACTFINALTIME_MATCHSTEP; 383 } 384 385 ierr = TSGetSNES(ts,&ts->snes);CHKERRQ(ierr); 386 PetscFunctionReturn(0); 387 } 388 389 #undef __FUNCT__ 390 #define __FUNCT__ "TSSetFromOptions_BDF" 391 static PetscErrorCode TSSetFromOptions_BDF(PetscOptionItems *PetscOptionsObject,TS ts) 392 { 393 TS_BDF *bdf = (TS_BDF*)ts->data; 394 PetscErrorCode ierr; 395 396 PetscFunctionBegin; 397 ierr = PetscOptionsHead(PetscOptionsObject,"BDF ODE solver options");CHKERRQ(ierr); 398 { 399 PetscBool flg; 400 PetscInt order = bdf->order; 401 PetscBool adapt = bdf->adapt; 402 ierr = PetscOptionsInt("-ts_bdf_order","Order of the BDF method","TSBDFSetOrder",order,&order,&flg);CHKERRQ(ierr); 403 if (flg) {ierr = TSBDFSetOrder(ts,order);CHKERRQ(ierr);} 404 ierr = PetscOptionsBool("-ts_bdf_adapt","Use time-step adaptivity with the BDF method","TSBDFUseAdapt",adapt,&adapt,&flg);CHKERRQ(ierr); 405 if (flg) {ierr = TSBDFUseAdapt(ts,adapt);CHKERRQ(ierr);} 406 } 407 ierr = PetscOptionsTail();CHKERRQ(ierr); 408 PetscFunctionReturn(0); 409 } 410 411 #undef __FUNCT__ 412 #define __FUNCT__ "TSView_BDF" 413 static PetscErrorCode TSView_BDF(TS ts,PetscViewer viewer) 414 { 415 TS_BDF *bdf = (TS_BDF*)ts->data; 416 PetscBool iascii; 417 PetscErrorCode ierr; 418 419 PetscFunctionBegin; 420 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 421 if (iascii) { 422 ierr = PetscViewerASCIIPrintf(viewer," Order=%D\n",bdf->order);CHKERRQ(ierr); 423 } 424 if (ts->adapt) {ierr = TSAdaptView(ts->adapt,viewer);CHKERRQ(ierr);} 425 if (ts->snes) {ierr = SNESView(ts->snes,viewer);CHKERRQ(ierr);} 426 PetscFunctionReturn(0); 427 } 428 429 /* ------------------------------------------------------------ */ 430 431 #undef __FUNCT__ 432 #define __FUNCT__ "TSBDFSetOrder_BDF" 433 static PetscErrorCode TSBDFSetOrder_BDF(TS ts,PetscInt order) 434 { 435 TS_BDF *bdf = (TS_BDF*)ts->data; 436 437 PetscFunctionBegin; 438 if (order == bdf->order) PetscFunctionReturn(0); 439 if (order < 1 || order > 6) SETERRQ1(PetscObjectComm((PetscObject)ts),PETSC_ERR_ARG_OUTOFRANGE,"BDF Order %D not implemented",order); 440 bdf->order = order; 441 PetscFunctionReturn(0); 442 } 443 444 #undef __FUNCT__ 445 #define __FUNCT__ "TSBDFGetOrder_BDF" 446 static PetscErrorCode TSBDFGetOrder_BDF(TS ts,PetscInt *order) 447 { 448 TS_BDF *bdf = (TS_BDF*)ts->data; 449 450 PetscFunctionBegin; 451 *order = bdf->order; 452 PetscFunctionReturn(0); 453 } 454 455 #undef __FUNCT__ 456 #define __FUNCT__ "TSBDFUseAdapt_BDF" 457 static PetscErrorCode TSBDFUseAdapt_BDF(TS ts,PetscBool use) 458 { 459 TS_BDF *bdf = (TS_BDF*)ts->data; 460 461 PetscFunctionBegin; 462 if (use == bdf->adapt) PetscFunctionReturn(0); 463 if (ts->setupcalled) SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_ORDER,"Cannot change adaptivity after TSSetUp()"); 464 bdf->adapt = use; 465 PetscFunctionReturn(0); 466 } 467 468 /* ------------------------------------------------------------ */ 469 470 /*MC 471 TSBDF - DAE solver using BDF methods 472 473 Level: beginner 474 475 .seealso: TS, TSCreate(), TSSetType() 476 M*/ 477 #undef __FUNCT__ 478 #define __FUNCT__ "TSCreate_BDF" 479 PETSC_EXTERN PetscErrorCode TSCreate_BDF(TS ts) 480 { 481 TS_BDF *bdf; 482 PetscErrorCode ierr; 483 484 PetscFunctionBegin; 485 ts->ops->reset = TSReset_BDF; 486 ts->ops->destroy = TSDestroy_BDF; 487 ts->ops->view = TSView_BDF; 488 ts->ops->setup = TSSetUp_BDF; 489 ts->ops->setfromoptions = TSSetFromOptions_BDF; 490 ts->ops->step = TSStep_BDF; 491 ts->ops->evaluatewlte = TSEvaluateWLTE_BDF; 492 ts->ops->rollback = TSRollBack_BDF; 493 ts->ops->interpolate = TSInterpolate_BDF; 494 ts->ops->snesfunction = SNESTSFormFunction_BDF; 495 ts->ops->snesjacobian = SNESTSFormJacobian_BDF; 496 497 ierr = PetscNewLog(ts,&bdf);CHKERRQ(ierr); 498 ts->data = (void*)bdf; 499 500 bdf->order = 2; 501 bdf->adapt = PETSC_FALSE; 502 bdf->status = TS_STEP_COMPLETE; 503 504 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSBDFSetOrder_C",TSBDFSetOrder_BDF);CHKERRQ(ierr); 505 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSBDFGetOrder_C",TSBDFGetOrder_BDF);CHKERRQ(ierr); 506 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSBDFUseAdapt_C",TSBDFUseAdapt_BDF);CHKERRQ(ierr); 507 PetscFunctionReturn(0); 508 } 509 510 /* ------------------------------------------------------------ */ 511 512 #undef __FUNCT__ 513 #define __FUNCT__ "TSBDFSetOrder" 514 /*@ 515 TSBDFSetOrder - Set the order of the BDF method 516 517 Logically Collective on TS 518 519 Input Parameter: 520 + ts - timestepping context 521 - order - order of the method 522 523 Options Database: 524 . -ts_bdf_order <order> 525 526 Level: intermediate 527 528 @*/ 529 PetscErrorCode TSBDFSetOrder(TS ts,PetscInt order) 530 { 531 PetscErrorCode ierr; 532 533 PetscFunctionBegin; 534 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 535 PetscValidLogicalCollectiveInt(ts,order,2); 536 ierr = PetscTryMethod(ts,"TSBDFSetOrder_C",(TS,PetscInt),(ts,order));CHKERRQ(ierr); 537 PetscFunctionReturn(0); 538 } 539 540 #undef __FUNCT__ 541 #define __FUNCT__ "TSBDFGetOrder" 542 /*@ 543 TSBDFGetOrder - Get the order of the BDF method 544 545 Not Collective 546 547 Input Parameter: 548 . ts - timestepping context 549 550 Output Parameter: 551 . order - order of the method 552 553 Level: intermediate 554 555 @*/ 556 PetscErrorCode TSBDFGetOrder(TS ts,PetscInt *order) 557 { 558 PetscErrorCode ierr; 559 560 PetscFunctionBegin; 561 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 562 PetscValidIntPointer(order,2); 563 ierr = PetscUseMethod(ts,"TSBDFGetOrder_C",(TS,PetscInt*),(ts,order));CHKERRQ(ierr); 564 PetscFunctionReturn(0); 565 } 566 567 #undef __FUNCT__ 568 #define __FUNCT__ "TSBDFUseAdapt" 569 /*@ 570 TSBDFUseAdapt - Use time-step adaptivity with the BDF method 571 572 Logically Collective on TS 573 574 Input Parameter: 575 + ts - timestepping context 576 - use - flag to use adaptivity 577 578 Options Database: 579 . -ts_bdf_adapt 580 581 Level: intermediate 582 583 .seealso: TSAdapt 584 @*/ 585 PetscErrorCode TSBDFUseAdapt(TS ts,PetscBool use) 586 { 587 PetscErrorCode ierr; 588 589 PetscFunctionBegin; 590 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 591 PetscValidLogicalCollectiveBool(ts,use,2); 592 ierr = PetscTryMethod(ts,"TSBDFUseAdapt_C",(TS,PetscBool),(ts,use));CHKERRQ(ierr); 593 PetscFunctionReturn(0); 594 } 595