Actual source code: alpha1.c

  1: /*
  2:   Code for timestepping with implicit generalized-\alpha method
  3:   for first order systems.
  4: */
  5: #include <petsc/private/tsimpl.h>

  7: static PetscBool  cited      = PETSC_FALSE;
  8: static const char citation[] = "@article{Jansen2000,\n"
  9:                                "  title   = {A generalized-$\\alpha$ method for integrating the filtered {N}avier--{S}tokes equations with a stabilized finite element method},\n"
 10:                                "  author  = {Kenneth E. Jansen and Christian H. Whiting and Gregory M. Hulbert},\n"
 11:                                "  journal = {Computer Methods in Applied Mechanics and Engineering},\n"
 12:                                "  volume  = {190},\n"
 13:                                "  number  = {3--4},\n"
 14:                                "  pages   = {305--319},\n"
 15:                                "  year    = {2000},\n"
 16:                                "  issn    = {0045-7825},\n"
 17:                                "  doi     = {http://dx.doi.org/10.1016/S0045-7825(00)00203-6}\n}\n";

 19: typedef struct {
 20:   PetscReal stage_time;
 21:   PetscReal shift_V;
 22:   PetscReal scale_F;
 23:   Vec       X0, Xa, X1;
 24:   Vec       V0, Va, V1;

 26:   PetscReal Alpha_m;
 27:   PetscReal Alpha_f;
 28:   PetscReal Gamma;
 29:   PetscInt  order;

 31:   Vec vec_sol_prev;
 32:   Vec vec_lte_work;

 34:   TSStepStatus status;
 35: } TS_Alpha;

 37: static PetscErrorCode TSAlpha_StageTime(TS ts)
 38: {
 39:   TS_Alpha *th      = (TS_Alpha *)ts->data;
 40:   PetscReal t       = ts->ptime;
 41:   PetscReal dt      = ts->time_step;
 42:   PetscReal Alpha_m = th->Alpha_m;
 43:   PetscReal Alpha_f = th->Alpha_f;
 44:   PetscReal Gamma   = th->Gamma;

 46:   th->stage_time = t + Alpha_f * dt;
 47:   th->shift_V    = Alpha_m / (Alpha_f * Gamma * dt);
 48:   th->scale_F    = 1 / Alpha_f;
 49:   return 0;
 50: }

 52: static PetscErrorCode TSAlpha_StageVecs(TS ts, Vec X)
 53: {
 54:   TS_Alpha *th = (TS_Alpha *)ts->data;
 55:   Vec       X1 = X, V1 = th->V1;
 56:   Vec       Xa = th->Xa, Va = th->Va;
 57:   Vec       X0 = th->X0, V0 = th->V0;
 58:   PetscReal dt      = ts->time_step;
 59:   PetscReal Alpha_m = th->Alpha_m;
 60:   PetscReal Alpha_f = th->Alpha_f;
 61:   PetscReal Gamma   = th->Gamma;

 63:   /* V1 = 1/(Gamma*dT)*(X1-X0) + (1-1/Gamma)*V0 */
 64:   VecWAXPY(V1, -1.0, X0, X1);
 65:   VecAXPBY(V1, 1 - 1 / Gamma, 1 / (Gamma * dt), V0);
 66:   /* Xa = X0 + Alpha_f*(X1-X0) */
 67:   VecWAXPY(Xa, -1.0, X0, X1);
 68:   VecAYPX(Xa, Alpha_f, X0);
 69:   /* Va = V0 + Alpha_m*(V1-V0) */
 70:   VecWAXPY(Va, -1.0, V0, V1);
 71:   VecAYPX(Va, Alpha_m, V0);
 72:   return 0;
 73: }

 75: static PetscErrorCode TSAlpha_SNESSolve(TS ts, Vec b, Vec x)
 76: {
 77:   PetscInt nits, lits;

 79:   SNESSolve(ts->snes, b, x);
 80:   SNESGetIterationNumber(ts->snes, &nits);
 81:   SNESGetLinearSolveIterations(ts->snes, &lits);
 82:   ts->snes_its += nits;
 83:   ts->ksp_its += lits;
 84:   return 0;
 85: }

 87: /*
 88:   Compute a consistent initial state for the generalized-alpha method.
 89:   - Solve two successive backward Euler steps with halved time step.
 90:   - Compute the initial time derivative using backward differences.
 91:   - If using adaptivity, estimate the LTE of the initial step.
 92: */
 93: static PetscErrorCode TSAlpha_Restart(TS ts, PetscBool *initok)
 94: {
 95:   TS_Alpha *th = (TS_Alpha *)ts->data;
 96:   PetscReal time_step;
 97:   PetscReal alpha_m, alpha_f, gamma;
 98:   Vec       X0 = ts->vec_sol, X1, X2 = th->X1;
 99:   PetscBool stageok;

101:   VecDuplicate(X0, &X1);

103:   /* Setup backward Euler with halved time step */
104:   TSAlphaGetParams(ts, &alpha_m, &alpha_f, &gamma);
105:   TSAlphaSetParams(ts, 1, 1, 1);
106:   TSGetTimeStep(ts, &time_step);
107:   ts->time_step = time_step / 2;
108:   TSAlpha_StageTime(ts);
109:   th->stage_time = ts->ptime;
110:   VecZeroEntries(th->V0);

112:   /* First BE step, (t0,X0) -> (t1,X1) */
113:   th->stage_time += ts->time_step;
114:   VecCopy(X0, th->X0);
115:   TSPreStage(ts, th->stage_time);
116:   VecCopy(th->X0, X1);
117:   TSAlpha_SNESSolve(ts, NULL, X1);
118:   TSPostStage(ts, th->stage_time, 0, &X1);
119:   TSAdaptCheckStage(ts->adapt, ts, th->stage_time, X1, &stageok);
120:   if (!stageok) goto finally;

122:   /* Second BE step, (t1,X1) -> (t2,X2) */
123:   th->stage_time += ts->time_step;
124:   VecCopy(X1, th->X0);
125:   TSPreStage(ts, th->stage_time);
126:   VecCopy(th->X0, X2);
127:   TSAlpha_SNESSolve(ts, NULL, X2);
128:   TSPostStage(ts, th->stage_time, 0, &X2);
129:   TSAdaptCheckStage(ts->adapt, ts, th->stage_time, X2, &stageok);
130:   if (!stageok) goto finally;

132:   /* Compute V0 ~ dX/dt at t0 with backward differences */
133:   VecZeroEntries(th->V0);
134:   VecAXPY(th->V0, -3 / ts->time_step, X0);
135:   VecAXPY(th->V0, +4 / ts->time_step, X1);
136:   VecAXPY(th->V0, -1 / ts->time_step, X2);

138:   /* Rough, lower-order estimate LTE of the initial step */
139:   if (th->vec_lte_work) {
140:     VecZeroEntries(th->vec_lte_work);
141:     VecAXPY(th->vec_lte_work, +2, X2);
142:     VecAXPY(th->vec_lte_work, -4, X1);
143:     VecAXPY(th->vec_lte_work, +2, X0);
144:   }

146: finally:
147:   /* Revert TSAlpha to the initial state (t0,X0) */
148:   if (initok) *initok = stageok;
149:   TSSetTimeStep(ts, time_step);
150:   TSAlphaSetParams(ts, alpha_m, alpha_f, gamma);
151:   VecCopy(ts->vec_sol, th->X0);

153:   VecDestroy(&X1);
154:   return 0;
155: }

157: static PetscErrorCode TSStep_Alpha(TS ts)
158: {
159:   TS_Alpha *th         = (TS_Alpha *)ts->data;
160:   PetscInt  rejections = 0;
161:   PetscBool stageok, accept = PETSC_TRUE;
162:   PetscReal next_time_step = ts->time_step;

164:   PetscCitationsRegister(citation, &cited);

166:   if (!ts->steprollback) {
167:     if (th->vec_sol_prev) VecCopy(th->X0, th->vec_sol_prev);
168:     VecCopy(ts->vec_sol, th->X0);
169:     VecCopy(th->V1, th->V0);
170:   }

172:   th->status = TS_STEP_INCOMPLETE;
173:   while (!ts->reason && th->status != TS_STEP_COMPLETE) {
174:     if (ts->steprestart) {
175:       TSAlpha_Restart(ts, &stageok);
176:       if (!stageok) goto reject_step;
177:     }

179:     TSAlpha_StageTime(ts);
180:     VecCopy(th->X0, th->X1);
181:     TSPreStage(ts, th->stage_time);
182:     TSAlpha_SNESSolve(ts, NULL, th->X1);
183:     TSPostStage(ts, th->stage_time, 0, &th->Xa);
184:     TSAdaptCheckStage(ts->adapt, ts, th->stage_time, th->Xa, &stageok);
185:     if (!stageok) goto reject_step;

187:     th->status = TS_STEP_PENDING;
188:     VecCopy(th->X1, ts->vec_sol);
189:     TSAdaptChoose(ts->adapt, ts, ts->time_step, NULL, &next_time_step, &accept);
190:     th->status = accept ? TS_STEP_COMPLETE : TS_STEP_INCOMPLETE;
191:     if (!accept) {
192:       VecCopy(th->X0, ts->vec_sol);
193:       ts->time_step = next_time_step;
194:       goto reject_step;
195:     }

197:     ts->ptime += ts->time_step;
198:     ts->time_step = next_time_step;
199:     break;

201:   reject_step:
202:     ts->reject++;
203:     accept = PETSC_FALSE;
204:     if (!ts->reason && ++rejections > ts->max_reject && ts->max_reject >= 0) {
205:       ts->reason = TS_DIVERGED_STEP_REJECTED;
206:       PetscInfo(ts, "Step=%" PetscInt_FMT ", step rejections %" PetscInt_FMT " greater than current TS allowed, stopping solve\n", ts->steps, rejections);
207:     }
208:   }
209:   return 0;
210: }

212: static PetscErrorCode TSEvaluateWLTE_Alpha(TS ts, NormType wnormtype, PetscInt *order, PetscReal *wlte)
213: {
214:   TS_Alpha *th = (TS_Alpha *)ts->data;
215:   Vec       X  = th->X1;           /* X = solution */
216:   Vec       Y  = th->vec_lte_work; /* Y = X + LTE  */
217:   PetscReal wltea, wlter;

219:   if (!th->vec_sol_prev) {
220:     *wlte = -1;
221:     return 0;
222:   }
223:   if (!th->vec_lte_work) {
224:     *wlte = -1;
225:     return 0;
226:   }
227:   if (ts->steprestart) {
228:     /* th->vec_lte_work is set to the LTE in TSAlpha_Restart() */
229:     VecAXPY(Y, 1, X);
230:   } else {
231:     /* Compute LTE using backward differences with non-constant time step */
232:     PetscReal   h = ts->time_step, h_prev = ts->ptime - ts->ptime_prev;
233:     PetscReal   a = 1 + h_prev / h;
234:     PetscScalar scal[3];
235:     Vec         vecs[3];
236:     scal[0] = +1 / a;
237:     scal[1] = -1 / (a - 1);
238:     scal[2] = +1 / (a * (a - 1));
239:     vecs[0] = th->X1;
240:     vecs[1] = th->X0;
241:     vecs[2] = th->vec_sol_prev;
242:     VecCopy(X, Y);
243:     VecMAXPY(Y, 3, scal, vecs);
244:   }
245:   TSErrorWeightedNorm(ts, X, Y, wnormtype, wlte, &wltea, &wlter);
246:   if (order) *order = 2;
247:   return 0;
248: }

250: static PetscErrorCode TSRollBack_Alpha(TS ts)
251: {
252:   TS_Alpha *th = (TS_Alpha *)ts->data;

254:   VecCopy(th->X0, ts->vec_sol);
255:   return 0;
256: }

258: static PetscErrorCode TSInterpolate_Alpha(TS ts, PetscReal t, Vec X)
259: {
260:   TS_Alpha *th = (TS_Alpha *)ts->data;
261:   PetscReal dt = t - ts->ptime;

263:   VecCopy(ts->vec_sol, X);
264:   VecAXPY(X, th->Gamma * dt, th->V1);
265:   VecAXPY(X, (1 - th->Gamma) * dt, th->V0);
266:   return 0;
267: }

269: static PetscErrorCode SNESTSFormFunction_Alpha(PETSC_UNUSED SNES snes, Vec X, Vec F, TS ts)
270: {
271:   TS_Alpha *th = (TS_Alpha *)ts->data;
272:   PetscReal ta = th->stage_time;
273:   Vec       Xa = th->Xa, Va = th->Va;

275:   TSAlpha_StageVecs(ts, X);
276:   /* F = Function(ta,Xa,Va) */
277:   TSComputeIFunction(ts, ta, Xa, Va, F, PETSC_FALSE);
278:   VecScale(F, th->scale_F);
279:   return 0;
280: }

282: static PetscErrorCode SNESTSFormJacobian_Alpha(PETSC_UNUSED SNES snes, PETSC_UNUSED Vec X, Mat J, Mat P, TS ts)
283: {
284:   TS_Alpha *th = (TS_Alpha *)ts->data;
285:   PetscReal ta = th->stage_time;
286:   Vec       Xa = th->Xa, Va = th->Va;
287:   PetscReal dVdX = th->shift_V;

289:   /* J,P = Jacobian(ta,Xa,Va) */
290:   TSComputeIJacobian(ts, ta, Xa, Va, dVdX, J, P, PETSC_FALSE);
291:   return 0;
292: }

294: static PetscErrorCode TSReset_Alpha(TS ts)
295: {
296:   TS_Alpha *th = (TS_Alpha *)ts->data;

298:   VecDestroy(&th->X0);
299:   VecDestroy(&th->Xa);
300:   VecDestroy(&th->X1);
301:   VecDestroy(&th->V0);
302:   VecDestroy(&th->Va);
303:   VecDestroy(&th->V1);
304:   VecDestroy(&th->vec_sol_prev);
305:   VecDestroy(&th->vec_lte_work);
306:   return 0;
307: }

309: static PetscErrorCode TSDestroy_Alpha(TS ts)
310: {
311:   TSReset_Alpha(ts);
312:   PetscFree(ts->data);

314:   PetscObjectComposeFunction((PetscObject)ts, "TSAlphaSetRadius_C", NULL);
315:   PetscObjectComposeFunction((PetscObject)ts, "TSAlphaSetParams_C", NULL);
316:   PetscObjectComposeFunction((PetscObject)ts, "TSAlphaGetParams_C", NULL);
317:   return 0;
318: }

320: static PetscErrorCode TSSetUp_Alpha(TS ts)
321: {
322:   TS_Alpha *th = (TS_Alpha *)ts->data;
323:   PetscBool match;

325:   VecDuplicate(ts->vec_sol, &th->X0);
326:   VecDuplicate(ts->vec_sol, &th->Xa);
327:   VecDuplicate(ts->vec_sol, &th->X1);
328:   VecDuplicate(ts->vec_sol, &th->V0);
329:   VecDuplicate(ts->vec_sol, &th->Va);
330:   VecDuplicate(ts->vec_sol, &th->V1);

332:   TSGetAdapt(ts, &ts->adapt);
333:   TSAdaptCandidatesClear(ts->adapt);
334:   PetscObjectTypeCompare((PetscObject)ts->adapt, TSADAPTNONE, &match);
335:   if (!match) {
336:     VecDuplicate(ts->vec_sol, &th->vec_sol_prev);
337:     VecDuplicate(ts->vec_sol, &th->vec_lte_work);
338:   }

340:   TSGetSNES(ts, &ts->snes);
341:   return 0;
342: }

344: static PetscErrorCode TSSetFromOptions_Alpha(TS ts, PetscOptionItems *PetscOptionsObject)
345: {
346:   TS_Alpha *th = (TS_Alpha *)ts->data;

348:   PetscOptionsHeadBegin(PetscOptionsObject, "Generalized-Alpha ODE solver options");
349:   {
350:     PetscBool flg;
351:     PetscReal radius = 1;
352:     PetscOptionsReal("-ts_alpha_radius", "Spectral radius (high-frequency dissipation)", "TSAlphaSetRadius", radius, &radius, &flg);
353:     if (flg) TSAlphaSetRadius(ts, radius);
354:     PetscOptionsReal("-ts_alpha_alpha_m", "Algorithmic parameter alpha_m", "TSAlphaSetParams", th->Alpha_m, &th->Alpha_m, NULL);
355:     PetscOptionsReal("-ts_alpha_alpha_f", "Algorithmic parameter alpha_f", "TSAlphaSetParams", th->Alpha_f, &th->Alpha_f, NULL);
356:     PetscOptionsReal("-ts_alpha_gamma", "Algorithmic parameter gamma", "TSAlphaSetParams", th->Gamma, &th->Gamma, NULL);
357:     TSAlphaSetParams(ts, th->Alpha_m, th->Alpha_f, th->Gamma);
358:   }
359:   PetscOptionsHeadEnd();
360:   return 0;
361: }

363: static PetscErrorCode TSView_Alpha(TS ts, PetscViewer viewer)
364: {
365:   TS_Alpha *th = (TS_Alpha *)ts->data;
366:   PetscBool iascii;

368:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
369:   if (iascii) PetscViewerASCIIPrintf(viewer, "  Alpha_m=%g, Alpha_f=%g, Gamma=%g\n", (double)th->Alpha_m, (double)th->Alpha_f, (double)th->Gamma);
370:   return 0;
371: }

373: static PetscErrorCode TSAlphaSetRadius_Alpha(TS ts, PetscReal radius)
374: {
375:   PetscReal alpha_m, alpha_f, gamma;

378:   alpha_m = (PetscReal)0.5 * (3 - radius) / (1 + radius);
379:   alpha_f = 1 / (1 + radius);
380:   gamma   = (PetscReal)0.5 + alpha_m - alpha_f;
381:   TSAlphaSetParams(ts, alpha_m, alpha_f, gamma);
382:   return 0;
383: }

385: static PetscErrorCode TSAlphaSetParams_Alpha(TS ts, PetscReal alpha_m, PetscReal alpha_f, PetscReal gamma)
386: {
387:   TS_Alpha *th  = (TS_Alpha *)ts->data;
388:   PetscReal tol = 100 * PETSC_MACHINE_EPSILON;
389:   PetscReal res = ((PetscReal)0.5 + alpha_m - alpha_f) - gamma;

391:   th->Alpha_m = alpha_m;
392:   th->Alpha_f = alpha_f;
393:   th->Gamma   = gamma;
394:   th->order   = (PetscAbsReal(res) < tol) ? 2 : 1;
395:   return 0;
396: }

398: static PetscErrorCode TSAlphaGetParams_Alpha(TS ts, PetscReal *alpha_m, PetscReal *alpha_f, PetscReal *gamma)
399: {
400:   TS_Alpha *th = (TS_Alpha *)ts->data;

402:   if (alpha_m) *alpha_m = th->Alpha_m;
403:   if (alpha_f) *alpha_f = th->Alpha_f;
404:   if (gamma) *gamma = th->Gamma;
405:   return 0;
406: }

408: /*MC
409:       TSALPHA - ODE/DAE solver using the implicit Generalized-Alpha method
410:                 for first-order systems

412:   Level: beginner

414:   References:
415: + * - K.E. Jansen, C.H. Whiting, G.M. Hulber, "A generalized-alpha
416:   method for integrating the filtered Navier-Stokes equations with a
417:   stabilized finite element method", Computer Methods in Applied
418:   Mechanics and Engineering, 190, 305-319, 2000.
419:   DOI: 10.1016/S0045-7825(00)00203-6.
420: - * -  J. Chung, G.M.Hubert. "A Time Integration Algorithm for Structural
421:   Dynamics with Improved Numerical Dissipation: The Generalized-alpha
422:   Method" ASME Journal of Applied Mechanics, 60, 371:375, 1993.

424: .seealso: `TS`, `TSCreate()`, `TSSetType()`, `TSAlphaSetRadius()`, `TSAlphaSetParams()`
425: M*/
426: PETSC_EXTERN PetscErrorCode TSCreate_Alpha(TS ts)
427: {
428:   TS_Alpha *th;

430:   ts->ops->reset          = TSReset_Alpha;
431:   ts->ops->destroy        = TSDestroy_Alpha;
432:   ts->ops->view           = TSView_Alpha;
433:   ts->ops->setup          = TSSetUp_Alpha;
434:   ts->ops->setfromoptions = TSSetFromOptions_Alpha;
435:   ts->ops->step           = TSStep_Alpha;
436:   ts->ops->evaluatewlte   = TSEvaluateWLTE_Alpha;
437:   ts->ops->rollback       = TSRollBack_Alpha;
438:   ts->ops->interpolate    = TSInterpolate_Alpha;
439:   ts->ops->snesfunction   = SNESTSFormFunction_Alpha;
440:   ts->ops->snesjacobian   = SNESTSFormJacobian_Alpha;
441:   ts->default_adapt_type  = TSADAPTNONE;

443:   ts->usessnes = PETSC_TRUE;

445:   PetscNew(&th);
446:   ts->data = (void *)th;

448:   th->Alpha_m = 0.5;
449:   th->Alpha_f = 0.5;
450:   th->Gamma   = 0.5;
451:   th->order   = 2;

453:   PetscObjectComposeFunction((PetscObject)ts, "TSAlphaSetRadius_C", TSAlphaSetRadius_Alpha);
454:   PetscObjectComposeFunction((PetscObject)ts, "TSAlphaSetParams_C", TSAlphaSetParams_Alpha);
455:   PetscObjectComposeFunction((PetscObject)ts, "TSAlphaGetParams_C", TSAlphaGetParams_Alpha);
456:   return 0;
457: }

459: /*@
460:   TSAlphaSetRadius - sets the desired spectral radius of the method
461:                      (i.e. high-frequency numerical damping)

463:   Logically Collective on TS

465:   The algorithmic parameters \alpha_m and \alpha_f of the
466:   generalized-\alpha method can be computed in terms of a specified
467:   spectral radius \rho in [0,1] for infinite time step in order to
468:   control high-frequency numerical damping:
469:     \alpha_m = 0.5*(3-\rho)/(1+\rho)
470:     \alpha_f = 1/(1+\rho)

472:   Input Parameters:
473: +  ts - timestepping context
474: -  radius - the desired spectral radius

476:   Options Database:
477: .  -ts_alpha_radius <radius> - set alpha radius

479:   Level: intermediate

481: .seealso: `TSAlphaSetParams()`, `TSAlphaGetParams()`
482: @*/
483: PetscErrorCode TSAlphaSetRadius(TS ts, PetscReal radius)
484: {
488:   PetscTryMethod(ts, "TSAlphaSetRadius_C", (TS, PetscReal), (ts, radius));
489:   return 0;
490: }

492: /*@
493:   TSAlphaSetParams - sets the algorithmic parameters for TSALPHA

495:   Logically Collective on TS

497:   Second-order accuracy can be obtained so long as:
498:     \gamma = 0.5 + alpha_m - alpha_f

500:   Unconditional stability requires:
501:     \alpha_m >= \alpha_f >= 0.5

503:   Backward Euler method is recovered with:
504:     \alpha_m = \alpha_f = gamma = 1

506:   Input Parameters:
507: +  ts - timestepping context
508: .  \alpha_m - algorithmic parameter
509: .  \alpha_f - algorithmic parameter
510: -  \gamma   - algorithmic parameter

512:    Options Database:
513: +  -ts_alpha_alpha_m <alpha_m> - set alpha_m
514: .  -ts_alpha_alpha_f <alpha_f> - set alpha_f
515: -  -ts_alpha_gamma   <gamma> - set gamma

517:   Note:
518:   Use of this function is normally only required to hack TSALPHA to
519:   use a modified integration scheme. Users should call
520:   TSAlphaSetRadius() to set the desired spectral radius of the methods
521:   (i.e. high-frequency damping) in order so select optimal values for
522:   these parameters.

524:   Level: advanced

526: .seealso: `TSAlphaSetRadius()`, `TSAlphaGetParams()`
527: @*/
528: PetscErrorCode TSAlphaSetParams(TS ts, PetscReal alpha_m, PetscReal alpha_f, PetscReal gamma)
529: {
534:   PetscTryMethod(ts, "TSAlphaSetParams_C", (TS, PetscReal, PetscReal, PetscReal), (ts, alpha_m, alpha_f, gamma));
535:   return 0;
536: }

538: /*@
539:   TSAlphaGetParams - gets the algorithmic parameters for TSALPHA

541:   Not Collective

543:   Input Parameter:
544: .  ts - timestepping context

546:   Output Parameters:
547: +  \alpha_m - algorithmic parameter
548: .  \alpha_f - algorithmic parameter
549: -  \gamma   - algorithmic parameter

551:   Note:
552:   Use of this function is normally only required to hack TSALPHA to
553:   use a modified integration scheme. Users should call
554:   TSAlphaSetRadius() to set the high-frequency damping (i.e. spectral
555:   radius of the method) in order so select optimal values for these
556:   parameters.

558:   Level: advanced

560: .seealso: `TSAlphaSetRadius()`, `TSAlphaSetParams()`
561: @*/
562: PetscErrorCode TSAlphaGetParams(TS ts, PetscReal *alpha_m, PetscReal *alpha_f, PetscReal *gamma)
563: {
568:   PetscUseMethod(ts, "TSAlphaGetParams_C", (TS, PetscReal *, PetscReal *, PetscReal *), (ts, alpha_m, alpha_f, gamma));
569:   return 0;
570: }