1 /* 2 Code for Timestepping with explicit SSP. 3 */ 4 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 5 6 PetscFunctionList TSSSPList = 0; 7 static PetscBool TSSSPPackageInitialized; 8 9 typedef struct { 10 PetscErrorCode (*onestep)(TS,PetscReal,PetscReal,Vec); 11 char *type_name; 12 PetscInt nstages; 13 Vec *work; 14 PetscInt nwork; 15 PetscBool workout; 16 } TS_SSP; 17 18 19 static PetscErrorCode TSSSPGetWorkVectors(TS ts,PetscInt n,Vec **work) 20 { 21 TS_SSP *ssp = (TS_SSP*)ts->data; 22 PetscErrorCode ierr; 23 24 PetscFunctionBegin; 25 if (ssp->workout) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Work vectors already gotten"); 26 if (ssp->nwork < n) { 27 if (ssp->nwork > 0) { 28 ierr = VecDestroyVecs(ssp->nwork,&ssp->work);CHKERRQ(ierr); 29 } 30 ierr = VecDuplicateVecs(ts->vec_sol,n,&ssp->work);CHKERRQ(ierr); 31 ssp->nwork = n; 32 } 33 *work = ssp->work; 34 ssp->workout = PETSC_TRUE; 35 PetscFunctionReturn(0); 36 } 37 38 static PetscErrorCode TSSSPRestoreWorkVectors(TS ts,PetscInt n,Vec **work) 39 { 40 TS_SSP *ssp = (TS_SSP*)ts->data; 41 42 PetscFunctionBegin; 43 if (!ssp->workout) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"Work vectors have not been gotten"); 44 if (*work != ssp->work) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Wrong work vectors checked out"); 45 ssp->workout = PETSC_FALSE; 46 *work = NULL; 47 PetscFunctionReturn(0); 48 } 49 50 /*MC 51 TSSSPRKS2 - Optimal second order SSP Runge-Kutta method, low-storage, c_eff=(s-1)/s 52 53 Pseudocode 2 of Ketcheson 2008 54 55 Level: beginner 56 57 .seealso: TSSSP, TSSSPSetType(), TSSSPSetNumStages() 58 M*/ 59 static PetscErrorCode TSSSPStep_RK_2(TS ts,PetscReal t0,PetscReal dt,Vec sol) 60 { 61 TS_SSP *ssp = (TS_SSP*)ts->data; 62 Vec *work,F; 63 PetscInt i,s; 64 PetscErrorCode ierr; 65 66 PetscFunctionBegin; 67 s = ssp->nstages; 68 ierr = TSSSPGetWorkVectors(ts,2,&work);CHKERRQ(ierr); 69 F = work[1]; 70 ierr = VecCopy(sol,work[0]);CHKERRQ(ierr); 71 for (i=0; i<s-1; i++) { 72 PetscReal stage_time = t0+dt*(i/(s-1.)); 73 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 74 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 75 ierr = VecAXPY(work[0],dt/(s-1.),F);CHKERRQ(ierr); 76 } 77 ierr = TSComputeRHSFunction(ts,t0+dt,work[0],F);CHKERRQ(ierr); 78 ierr = VecAXPBYPCZ(sol,(s-1.)/s,dt/s,1./s,work[0],F);CHKERRQ(ierr); 79 ierr = TSSSPRestoreWorkVectors(ts,2,&work);CHKERRQ(ierr); 80 PetscFunctionReturn(0); 81 } 82 83 /*MC 84 TSSSPRKS3 - Optimal third order SSP Runge-Kutta, low-storage, c_eff=(PetscSqrtReal(s)-1)/PetscSqrtReal(s), where PetscSqrtReal(s) is an integer 85 86 Pseudocode 2 of Ketcheson 2008 87 88 Level: beginner 89 90 .seealso: TSSSP, TSSSPSetType(), TSSSPSetNumStages() 91 M*/ 92 static PetscErrorCode TSSSPStep_RK_3(TS ts,PetscReal t0,PetscReal dt,Vec sol) 93 { 94 TS_SSP *ssp = (TS_SSP*)ts->data; 95 Vec *work,F; 96 PetscInt i,s,n,r; 97 PetscReal c,stage_time; 98 PetscErrorCode ierr; 99 100 PetscFunctionBegin; 101 s = ssp->nstages; 102 n = (PetscInt)(PetscSqrtReal((PetscReal)s)+0.001); 103 r = s-n; 104 if (n*n != s) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for optimal third order schemes with %d stages, must be a square number at least 4",s); 105 ierr = TSSSPGetWorkVectors(ts,3,&work);CHKERRQ(ierr); 106 F = work[2]; 107 ierr = VecCopy(sol,work[0]);CHKERRQ(ierr); 108 for (i=0; i<(n-1)*(n-2)/2; i++) { 109 c = (i<n*(n+1)/2) ? 1.*i/(s-n) : (1.*i-n)/(s-n); 110 stage_time = t0+c*dt; 111 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 112 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 113 ierr = VecAXPY(work[0],dt/r,F);CHKERRQ(ierr); 114 } 115 ierr = VecCopy(work[0],work[1]);CHKERRQ(ierr); 116 for (; i<n*(n+1)/2-1; i++) { 117 c = (i<n*(n+1)/2) ? 1.*i/(s-n) : (1.*i-n)/(s-n); 118 stage_time = t0+c*dt; 119 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 120 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 121 ierr = VecAXPY(work[0],dt/r,F);CHKERRQ(ierr); 122 } 123 { 124 c = (i<n*(n+1)/2) ? 1.*i/(s-n) : (1.*i-n)/(s-n); 125 stage_time = t0+c*dt; 126 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 127 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 128 ierr = VecAXPBYPCZ(work[0],1.*n/(2*n-1.),(n-1.)*dt/(r*(2*n-1)),(n-1.)/(2*n-1.),work[1],F);CHKERRQ(ierr); 129 i++; 130 } 131 for (; i<s; i++) { 132 c = (i<n*(n+1)/2) ? 1.*i/(s-n) : (1.*i-n)/(s-n); 133 stage_time = t0+c*dt; 134 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 135 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 136 ierr = VecAXPY(work[0],dt/r,F);CHKERRQ(ierr); 137 } 138 ierr = VecCopy(work[0],sol);CHKERRQ(ierr); 139 ierr = TSSSPRestoreWorkVectors(ts,3,&work);CHKERRQ(ierr); 140 PetscFunctionReturn(0); 141 } 142 143 /*MC 144 TSSSPRKS104 - Optimal fourth order SSP Runge-Kutta, low-storage (2N), c_eff=0.6 145 146 SSPRK(10,4), Pseudocode 3 of Ketcheson 2008 147 148 Level: beginner 149 150 .seealso: TSSSP, TSSSPSetType() 151 M*/ 152 static PetscErrorCode TSSSPStep_RK_10_4(TS ts,PetscReal t0,PetscReal dt,Vec sol) 153 { 154 const PetscReal c[10] = {0, 1./6, 2./6, 3./6, 4./6, 2./6, 3./6, 4./6, 5./6, 1}; 155 Vec *work,F; 156 PetscInt i; 157 PetscReal stage_time; 158 PetscErrorCode ierr; 159 160 PetscFunctionBegin; 161 ierr = TSSSPGetWorkVectors(ts,3,&work);CHKERRQ(ierr); 162 F = work[2]; 163 ierr = VecCopy(sol,work[0]);CHKERRQ(ierr); 164 for (i=0; i<5; i++) { 165 stage_time = t0+c[i]*dt; 166 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 167 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 168 ierr = VecAXPY(work[0],dt/6,F);CHKERRQ(ierr); 169 } 170 ierr = VecAXPBYPCZ(work[1],1./25,9./25,0,sol,work[0]);CHKERRQ(ierr); 171 ierr = VecAXPBY(work[0],15,-5,work[1]);CHKERRQ(ierr); 172 for (; i<9; i++) { 173 stage_time = t0+c[i]*dt; 174 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 175 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 176 ierr = VecAXPY(work[0],dt/6,F);CHKERRQ(ierr); 177 } 178 stage_time = t0+dt; 179 ierr = TSPreStage(ts,stage_time);CHKERRQ(ierr); 180 ierr = TSComputeRHSFunction(ts,stage_time,work[0],F);CHKERRQ(ierr); 181 ierr = VecAXPBYPCZ(work[1],3./5,dt/10,1,work[0],F);CHKERRQ(ierr); 182 ierr = VecCopy(work[1],sol);CHKERRQ(ierr); 183 ierr = TSSSPRestoreWorkVectors(ts,3,&work);CHKERRQ(ierr); 184 PetscFunctionReturn(0); 185 } 186 187 188 static PetscErrorCode TSSetUp_SSP(TS ts) 189 { 190 PetscErrorCode ierr; 191 192 PetscFunctionBegin; 193 ierr = TSCheckImplicitTerm(ts);CHKERRQ(ierr); 194 ierr = TSGetAdapt(ts,&ts->adapt);CHKERRQ(ierr); 195 ierr = TSAdaptCandidatesClear(ts->adapt);CHKERRQ(ierr); 196 PetscFunctionReturn(0); 197 } 198 199 static PetscErrorCode TSStep_SSP(TS ts) 200 { 201 TS_SSP *ssp = (TS_SSP*)ts->data; 202 Vec sol = ts->vec_sol; 203 PetscBool stageok,accept = PETSC_TRUE; 204 PetscReal next_time_step = ts->time_step; 205 PetscErrorCode ierr; 206 207 PetscFunctionBegin; 208 ierr = (*ssp->onestep)(ts,ts->ptime,ts->time_step,sol);CHKERRQ(ierr); 209 ierr = TSPostStage(ts,ts->ptime,0,&sol);CHKERRQ(ierr); 210 ierr = TSAdaptCheckStage(ts->adapt,ts,ts->ptime+ts->time_step,sol,&stageok);CHKERRQ(ierr); 211 if(!stageok) {ts->reason = TS_DIVERGED_STEP_REJECTED; PetscFunctionReturn(0);} 212 213 ierr = TSAdaptChoose(ts->adapt,ts,ts->time_step,NULL,&next_time_step,&accept);CHKERRQ(ierr); 214 if (!accept) {ts->reason = TS_DIVERGED_STEP_REJECTED; PetscFunctionReturn(0);} 215 216 ts->ptime += ts->time_step; 217 ts->time_step = next_time_step; 218 PetscFunctionReturn(0); 219 } 220 /*------------------------------------------------------------*/ 221 222 static PetscErrorCode TSReset_SSP(TS ts) 223 { 224 TS_SSP *ssp = (TS_SSP*)ts->data; 225 PetscErrorCode ierr; 226 227 PetscFunctionBegin; 228 if (ssp->work) {ierr = VecDestroyVecs(ssp->nwork,&ssp->work);CHKERRQ(ierr);} 229 ssp->nwork = 0; 230 ssp->workout = PETSC_FALSE; 231 PetscFunctionReturn(0); 232 } 233 234 static PetscErrorCode TSDestroy_SSP(TS ts) 235 { 236 TS_SSP *ssp = (TS_SSP*)ts->data; 237 PetscErrorCode ierr; 238 239 PetscFunctionBegin; 240 ierr = TSReset_SSP(ts);CHKERRQ(ierr); 241 ierr = PetscFree(ssp->type_name);CHKERRQ(ierr); 242 ierr = PetscFree(ts->data);CHKERRQ(ierr); 243 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPGetType_C",NULL);CHKERRQ(ierr); 244 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPSetType_C",NULL);CHKERRQ(ierr); 245 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPGetNumStages_C",NULL);CHKERRQ(ierr); 246 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPSetNumStages_C",NULL);CHKERRQ(ierr); 247 PetscFunctionReturn(0); 248 } 249 /*------------------------------------------------------------*/ 250 251 /*@C 252 TSSSPSetType - set the SSP time integration scheme to use 253 254 Logically Collective 255 256 Input Arguments: 257 ts - time stepping object 258 type - type of scheme to use 259 260 Options Database Keys: 261 -ts_ssp_type <rks2>: Type of SSP method (one of) rks2 rks3 rk104 262 -ts_ssp_nstages <5>: Number of stages 263 264 Level: beginner 265 266 .seealso: TSSSP, TSSSPGetType(), TSSSPSetNumStages(), TSSSPRKS2, TSSSPRKS3, TSSSPRK104 267 @*/ 268 PetscErrorCode TSSSPSetType(TS ts,TSSSPType type) 269 { 270 PetscErrorCode ierr; 271 272 PetscFunctionBegin; 273 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 274 ierr = PetscTryMethod(ts,"TSSSPSetType_C",(TS,TSSSPType),(ts,type));CHKERRQ(ierr); 275 PetscFunctionReturn(0); 276 } 277 278 /*@C 279 TSSSPGetType - get the SSP time integration scheme 280 281 Logically Collective 282 283 Input Argument: 284 ts - time stepping object 285 286 Output Argument: 287 type - type of scheme being used 288 289 Level: beginner 290 291 .seealso: TSSSP, TSSSPSettype(), TSSSPSetNumStages(), TSSSPRKS2, TSSSPRKS3, TSSSPRK104 292 @*/ 293 PetscErrorCode TSSSPGetType(TS ts,TSSSPType *type) 294 { 295 PetscErrorCode ierr; 296 297 PetscFunctionBegin; 298 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 299 ierr = PetscUseMethod(ts,"TSSSPGetType_C",(TS,TSSSPType*),(ts,type));CHKERRQ(ierr); 300 PetscFunctionReturn(0); 301 } 302 303 /*@ 304 TSSSPSetNumStages - set the number of stages to use with the SSP method 305 306 Logically Collective 307 308 Input Arguments: 309 ts - time stepping object 310 nstages - number of stages 311 312 Options Database Keys: 313 -ts_ssp_type <rks2>: NumStages of SSP method (one of) rks2 rks3 rk104 314 -ts_ssp_nstages <5>: Number of stages 315 316 Level: beginner 317 318 .seealso: TSSSP, TSSSPGetNumStages(), TSSSPSetNumStages(), TSSSPRKS2, TSSSPRKS3, TSSSPRK104 319 @*/ 320 PetscErrorCode TSSSPSetNumStages(TS ts,PetscInt nstages) 321 { 322 PetscErrorCode ierr; 323 324 PetscFunctionBegin; 325 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 326 ierr = PetscTryMethod(ts,"TSSSPSetNumStages_C",(TS,PetscInt),(ts,nstages));CHKERRQ(ierr); 327 PetscFunctionReturn(0); 328 } 329 330 /*@ 331 TSSSPGetNumStages - get the number of stages in the SSP time integration scheme 332 333 Logically Collective 334 335 Input Argument: 336 ts - time stepping object 337 338 Output Argument: 339 nstages - number of stages 340 341 Level: beginner 342 343 .seealso: TSSSP, TSSSPGetType(), TSSSPSetNumStages(), TSSSPRKS2, TSSSPRKS3, TSSSPRK104 344 @*/ 345 PetscErrorCode TSSSPGetNumStages(TS ts,PetscInt *nstages) 346 { 347 PetscErrorCode ierr; 348 349 PetscFunctionBegin; 350 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 351 ierr = PetscUseMethod(ts,"TSSSPGetNumStages_C",(TS,PetscInt*),(ts,nstages));CHKERRQ(ierr); 352 PetscFunctionReturn(0); 353 } 354 355 static PetscErrorCode TSSSPSetType_SSP(TS ts,TSSSPType type) 356 { 357 PetscErrorCode ierr,(*r)(TS,PetscReal,PetscReal,Vec); 358 TS_SSP *ssp = (TS_SSP*)ts->data; 359 360 PetscFunctionBegin; 361 ierr = PetscFunctionListFind(TSSSPList,type,&r);CHKERRQ(ierr); 362 if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown TS_SSP type %s given",type); 363 ssp->onestep = r; 364 ierr = PetscFree(ssp->type_name);CHKERRQ(ierr); 365 ierr = PetscStrallocpy(type,&ssp->type_name);CHKERRQ(ierr); 366 PetscFunctionReturn(0); 367 } 368 static PetscErrorCode TSSSPGetType_SSP(TS ts,TSSSPType *type) 369 { 370 TS_SSP *ssp = (TS_SSP*)ts->data; 371 372 PetscFunctionBegin; 373 *type = ssp->type_name; 374 PetscFunctionReturn(0); 375 } 376 static PetscErrorCode TSSSPSetNumStages_SSP(TS ts,PetscInt nstages) 377 { 378 TS_SSP *ssp = (TS_SSP*)ts->data; 379 380 PetscFunctionBegin; 381 ssp->nstages = nstages; 382 PetscFunctionReturn(0); 383 } 384 static PetscErrorCode TSSSPGetNumStages_SSP(TS ts,PetscInt *nstages) 385 { 386 TS_SSP *ssp = (TS_SSP*)ts->data; 387 388 PetscFunctionBegin; 389 *nstages = ssp->nstages; 390 PetscFunctionReturn(0); 391 } 392 393 static PetscErrorCode TSSetFromOptions_SSP(PetscOptionItems *PetscOptionsObject,TS ts) 394 { 395 char tname[256] = TSSSPRKS2; 396 TS_SSP *ssp = (TS_SSP*)ts->data; 397 PetscErrorCode ierr; 398 PetscBool flg; 399 400 PetscFunctionBegin; 401 ierr = PetscOptionsHead(PetscOptionsObject,"SSP ODE solver options");CHKERRQ(ierr); 402 { 403 ierr = PetscOptionsFList("-ts_ssp_type","Type of SSP method","TSSSPSetType",TSSSPList,tname,tname,sizeof(tname),&flg);CHKERRQ(ierr); 404 if (flg) { 405 ierr = TSSSPSetType(ts,tname);CHKERRQ(ierr); 406 } 407 ierr = PetscOptionsInt("-ts_ssp_nstages","Number of stages","TSSSPSetNumStages",ssp->nstages,&ssp->nstages,NULL);CHKERRQ(ierr); 408 } 409 ierr = PetscOptionsTail();CHKERRQ(ierr); 410 PetscFunctionReturn(0); 411 } 412 413 static PetscErrorCode TSView_SSP(TS ts,PetscViewer viewer) 414 { 415 TS_SSP *ssp = (TS_SSP*)ts->data; 416 PetscBool ascii; 417 PetscErrorCode ierr; 418 419 PetscFunctionBegin; 420 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&ascii);CHKERRQ(ierr); 421 if (ascii) {ierr = PetscViewerASCIIPrintf(viewer," Scheme: %s\n",ssp->type_name);CHKERRQ(ierr);} 422 PetscFunctionReturn(0); 423 } 424 425 /* ------------------------------------------------------------ */ 426 427 /*MC 428 TSSSP - Explicit strong stability preserving ODE solver 429 430 Most hyperbolic conservation laws have exact solutions that are total variation diminishing (TVD) or total variation 431 bounded (TVB) although these solutions often contain discontinuities. Spatial discretizations such as Godunov's 432 scheme and high-resolution finite volume methods (TVD limiters, ENO/WENO) are designed to preserve these properties, 433 but they are usually formulated using a forward Euler time discretization or by coupling the space and time 434 discretization as in the classical Lax-Wendroff scheme. When the space and time discretization is coupled, it is very 435 difficult to produce schemes with high temporal accuracy while preserving TVD properties. An alternative is the 436 semidiscrete formulation where we choose a spatial discretization that is TVD with forward Euler and then choose a 437 time discretization that preserves the TVD property. Such integrators are called strong stability preserving (SSP). 438 439 Let c_eff be the minimum number of function evaluations required to step as far as one step of forward Euler while 440 still being SSP. Some theoretical bounds 441 442 1. There are no explicit methods with c_eff > 1. 443 444 2. There are no explicit methods beyond order 4 (for nonlinear problems) and c_eff > 0. 445 446 3. There are no implicit methods with order greater than 1 and c_eff > 2. 447 448 This integrator provides Runge-Kutta methods of order 2, 3, and 4 with maximal values of c_eff. More stages allows 449 for larger values of c_eff which improves efficiency. These implementations are low-memory and only use 2 or 3 work 450 vectors regardless of the total number of stages, so e.g. 25-stage 3rd order methods may be an excellent choice. 451 452 Methods can be chosen with -ts_ssp_type {rks2,rks3,rk104} 453 454 rks2: Second order methods with any number s>1 of stages. c_eff = (s-1)/s 455 456 rks3: Third order methods with s=n^2 stages, n>1. c_eff = (s-n)/s 457 458 rk104: A 10-stage fourth order method. c_eff = 0.6 459 460 Level: beginner 461 462 References: 463 + 1. - Ketcheson, Highly efficient strong stability preserving Runge Kutta methods with low storage implementations, SISC, 2008. 464 - 2. - Gottlieb, Ketcheson, and Shu, High order strong stability preserving time discretizations, J Scientific Computing, 2009. 465 466 .seealso: TSCreate(), TS, TSSetType() 467 468 M*/ 469 PETSC_EXTERN PetscErrorCode TSCreate_SSP(TS ts) 470 { 471 TS_SSP *ssp; 472 PetscErrorCode ierr; 473 474 PetscFunctionBegin; 475 ierr = TSSSPInitializePackage();CHKERRQ(ierr); 476 477 ts->ops->setup = TSSetUp_SSP; 478 ts->ops->step = TSStep_SSP; 479 ts->ops->reset = TSReset_SSP; 480 ts->ops->destroy = TSDestroy_SSP; 481 ts->ops->setfromoptions = TSSetFromOptions_SSP; 482 ts->ops->view = TSView_SSP; 483 484 ierr = PetscNewLog(ts,&ssp);CHKERRQ(ierr); 485 ts->data = (void*)ssp; 486 487 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPGetType_C",TSSSPGetType_SSP);CHKERRQ(ierr); 488 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPSetType_C",TSSSPSetType_SSP);CHKERRQ(ierr); 489 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPGetNumStages_C",TSSSPGetNumStages_SSP);CHKERRQ(ierr); 490 ierr = PetscObjectComposeFunction((PetscObject)ts,"TSSSPSetNumStages_C",TSSSPSetNumStages_SSP);CHKERRQ(ierr); 491 492 ierr = TSSSPSetType(ts,TSSSPRKS2);CHKERRQ(ierr); 493 ssp->nstages = 5; 494 PetscFunctionReturn(0); 495 } 496 497 /*@C 498 TSSSPInitializePackage - This function initializes everything in the TSSSP package. It is called 499 from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to TSCreate_SSP() 500 when using static libraries. 501 502 Level: developer 503 504 .keywords: TS, TSSSP, initialize, package 505 .seealso: PetscInitialize() 506 @*/ 507 PetscErrorCode TSSSPInitializePackage(void) 508 { 509 PetscErrorCode ierr; 510 511 PetscFunctionBegin; 512 if (TSSSPPackageInitialized) PetscFunctionReturn(0); 513 TSSSPPackageInitialized = PETSC_TRUE; 514 ierr = PetscFunctionListAdd(&TSSSPList,TSSSPRKS2, TSSSPStep_RK_2);CHKERRQ(ierr); 515 ierr = PetscFunctionListAdd(&TSSSPList,TSSSPRKS3, TSSSPStep_RK_3);CHKERRQ(ierr); 516 ierr = PetscFunctionListAdd(&TSSSPList,TSSSPRK104,TSSSPStep_RK_10_4);CHKERRQ(ierr); 517 ierr = PetscRegisterFinalize(TSSSPFinalizePackage);CHKERRQ(ierr); 518 PetscFunctionReturn(0); 519 } 520 521 /*@C 522 TSSSPFinalizePackage - This function destroys everything in the TSSSP package. It is 523 called from PetscFinalize(). 524 525 Level: developer 526 527 .keywords: Petsc, destroy, package 528 .seealso: PetscFinalize() 529 @*/ 530 PetscErrorCode TSSSPFinalizePackage(void) 531 { 532 PetscErrorCode ierr; 533 534 PetscFunctionBegin; 535 TSSSPPackageInitialized = PETSC_FALSE; 536 ierr = PetscFunctionListDestroy(&TSSSPList);CHKERRQ(ierr); 537 PetscFunctionReturn(0); 538 } 539