Actual source code: mrk.c
1: /*
2: Code for time stepping with the multi-rate Runge-Kutta method
4: Notes:
5: 1) The general system is written as
6: Udot = F(t,U) for the nonsplit version of multi-rate RK,
7: user should give the indexes for both slow and fast components;
8: 2) The general system is written as
9: Usdot = Fs(t,Us,Uf)
10: Ufdot = Ff(t,Us,Uf) for multi-rate RK with RHS splits,
11: user should partition RHS by themselves and also provide the indexes for both slow and fast components.
12: */
14: #include <petsc/private/tsimpl.h>
15: #include <petscdm.h>
16: #include <../src/ts/impls/explicit/rk/rk.h>
17: #include <../src/ts/impls/explicit/rk/mrk.h>
19: static PetscErrorCode TSReset_RK_MultirateNonsplit(TS ts)
20: {
21: TS_RK *rk = (TS_RK *)ts->data;
22: RKTableau tab = rk->tableau;
24: VecDestroy(&rk->X0);
25: VecDestroyVecs(tab->s, &rk->YdotRHS_slow);
26: return 0;
27: }
29: static PetscErrorCode TSInterpolate_RK_MultirateNonsplit(TS ts, PetscReal itime, Vec X)
30: {
31: TS_RK *rk = (TS_RK *)ts->data;
32: PetscInt s = rk->tableau->s, p = rk->tableau->p, i, j;
33: PetscReal h = ts->time_step;
34: PetscReal tt, t;
35: PetscScalar *b;
36: const PetscReal *B = rk->tableau->binterp;
39: t = (itime - rk->ptime) / h;
40: PetscMalloc1(s, &b);
41: for (i = 0; i < s; i++) b[i] = 0;
42: for (j = 0, tt = t; j < p; j++, tt *= t) {
43: for (i = 0; i < s; i++) b[i] += h * B[i * p + j] * tt;
44: }
45: VecCopy(rk->X0, X);
46: VecMAXPY(X, s, b, rk->YdotRHS_slow);
47: PetscFree(b);
48: return 0;
49: }
51: static PetscErrorCode TSStepRefine_RK_MultirateNonsplit(TS ts)
52: {
53: TS previousts, subts;
54: TS_RK *rk = (TS_RK *)ts->data;
55: RKTableau tab = rk->tableau;
56: Vec *Y = rk->Y, *YdotRHS = rk->YdotRHS;
57: Vec vec_fast, subvec_fast;
58: const PetscInt s = tab->s;
59: const PetscReal *A = tab->A, *c = tab->c;
60: PetscScalar *w = rk->work;
61: PetscInt i, j, k;
62: PetscReal t = ts->ptime, h = ts->time_step;
64: VecDuplicate(ts->vec_sol, &vec_fast);
65: previousts = rk->subts_current;
66: TSRHSSplitGetSubTS(rk->subts_current, "fast", &subts);
67: TSRHSSplitGetSubTS(subts, "fast", &subts);
68: for (k = 0; k < rk->dtratio; k++) {
69: for (i = 0; i < s; i++) {
70: TSInterpolate_RK_MultirateNonsplit(ts, t + k * h / rk->dtratio + h / rk->dtratio * c[i], Y[i]);
71: for (j = 0; j < i; j++) w[j] = h / rk->dtratio * A[i * s + j];
72: /* update the fast components in the stage value, the slow components will be ignored, so we do not care the slow part in vec_fast */
73: VecCopy(ts->vec_sol, vec_fast);
74: VecMAXPY(vec_fast, i, w, YdotRHS);
75: /* update the fast components in the stage value */
76: VecGetSubVector(vec_fast, rk->is_fast, &subvec_fast);
77: VecISCopy(Y[i], rk->is_fast, SCATTER_FORWARD, subvec_fast);
78: VecRestoreSubVector(vec_fast, rk->is_fast, &subvec_fast);
79: /* compute the stage RHS */
80: TSComputeRHSFunction(ts, t + k * h / rk->dtratio + h / rk->dtratio * c[i], Y[i], YdotRHS[i]);
81: }
82: VecCopy(ts->vec_sol, vec_fast);
83: TSEvaluateStep(ts, tab->order, vec_fast, NULL);
84: /* update the fast components in the solution */
85: VecGetSubVector(vec_fast, rk->is_fast, &subvec_fast);
86: VecISCopy(ts->vec_sol, rk->is_fast, SCATTER_FORWARD, subvec_fast);
87: VecRestoreSubVector(vec_fast, rk->is_fast, &subvec_fast);
89: if (subts) {
90: Vec *YdotRHS_copy;
91: VecDuplicateVecs(ts->vec_sol, s, &YdotRHS_copy);
92: rk->subts_current = rk->subts_fast;
93: ts->ptime = t + k * h / rk->dtratio;
94: ts->time_step = h / rk->dtratio;
95: TSRHSSplitGetIS(rk->subts_current, "fast", &rk->is_fast);
96: for (i = 0; i < s; i++) {
97: VecCopy(rk->YdotRHS_slow[i], YdotRHS_copy[i]);
98: VecCopy(YdotRHS[i], rk->YdotRHS_slow[i]);
99: }
101: TSStepRefine_RK_MultirateNonsplit(ts);
103: rk->subts_current = previousts;
104: ts->ptime = t;
105: ts->time_step = h;
106: TSRHSSplitGetIS(previousts, "fast", &rk->is_fast);
107: for (i = 0; i < s; i++) VecCopy(YdotRHS_copy[i], rk->YdotRHS_slow[i]);
108: VecDestroyVecs(s, &YdotRHS_copy);
109: }
110: }
111: VecDestroy(&vec_fast);
112: return 0;
113: }
115: static PetscErrorCode TSStep_RK_MultirateNonsplit(TS ts)
116: {
117: TS_RK *rk = (TS_RK *)ts->data;
118: RKTableau tab = rk->tableau;
119: Vec *Y = rk->Y, *YdotRHS = rk->YdotRHS, *YdotRHS_slow = rk->YdotRHS_slow;
120: Vec stage_slow, sol_slow; /* vectors store the slow components */
121: Vec subvec_slow; /* sub vector to store the slow components */
122: IS is_slow = rk->is_slow;
123: const PetscInt s = tab->s;
124: const PetscReal *A = tab->A, *c = tab->c;
125: PetscScalar *w = rk->work;
126: PetscInt i, j, dtratio = rk->dtratio;
127: PetscReal next_time_step = ts->time_step, t = ts->ptime, h = ts->time_step;
129: rk->status = TS_STEP_INCOMPLETE;
130: VecDuplicate(ts->vec_sol, &stage_slow);
131: VecDuplicate(ts->vec_sol, &sol_slow);
132: VecCopy(ts->vec_sol, rk->X0);
133: for (i = 0; i < s; i++) {
134: rk->stage_time = t + h * c[i];
135: TSPreStage(ts, rk->stage_time);
136: VecCopy(ts->vec_sol, Y[i]);
137: for (j = 0; j < i; j++) w[j] = h * A[i * s + j];
138: VecMAXPY(Y[i], i, w, YdotRHS_slow);
139: TSPostStage(ts, rk->stage_time, i, Y);
140: /* compute the stage RHS */
141: TSComputeRHSFunction(ts, t + h * c[i], Y[i], YdotRHS_slow[i]);
142: }
143: /* update the slow components in the solution */
144: rk->YdotRHS = YdotRHS_slow;
145: rk->dtratio = 1;
146: TSEvaluateStep(ts, tab->order, sol_slow, NULL);
147: rk->dtratio = dtratio;
148: rk->YdotRHS = YdotRHS;
149: /* update the slow components in the solution */
150: VecGetSubVector(sol_slow, is_slow, &subvec_slow);
151: VecISCopy(ts->vec_sol, is_slow, SCATTER_FORWARD, subvec_slow);
152: VecRestoreSubVector(sol_slow, is_slow, &subvec_slow);
154: rk->subts_current = ts;
155: rk->ptime = t;
156: rk->time_step = h;
157: TSStepRefine_RK_MultirateNonsplit(ts);
159: ts->ptime = t + ts->time_step;
160: ts->time_step = next_time_step;
161: rk->status = TS_STEP_COMPLETE;
163: /* free memory of work vectors */
164: VecDestroy(&stage_slow);
165: VecDestroy(&sol_slow);
166: return 0;
167: }
169: static PetscErrorCode TSSetUp_RK_MultirateNonsplit(TS ts)
170: {
171: TS_RK *rk = (TS_RK *)ts->data;
172: RKTableau tab = rk->tableau;
174: TSRHSSplitGetIS(ts, "slow", &rk->is_slow);
175: TSRHSSplitGetIS(ts, "fast", &rk->is_fast);
177: TSRHSSplitGetSubTS(ts, "slow", &rk->subts_slow);
178: TSRHSSplitGetSubTS(ts, "fast", &rk->subts_fast);
180: VecDuplicate(ts->vec_sol, &rk->X0);
181: VecDuplicateVecs(ts->vec_sol, tab->s, &rk->YdotRHS_slow);
182: rk->subts_current = rk->subts_fast;
184: ts->ops->step = TSStep_RK_MultirateNonsplit;
185: ts->ops->interpolate = TSInterpolate_RK_MultirateNonsplit;
186: return 0;
187: }
189: /*
190: Copy DM from tssrc to tsdest, while keeping the original DMTS and DMSNES in tsdest.
191: */
192: static PetscErrorCode TSCopyDM(TS tssrc, TS tsdest)
193: {
194: DM newdm, dmsrc, dmdest;
196: TSGetDM(tssrc, &dmsrc);
197: DMClone(dmsrc, &newdm);
198: TSGetDM(tsdest, &dmdest);
199: DMCopyDMTS(dmdest, newdm);
200: DMCopyDMSNES(dmdest, newdm);
201: TSSetDM(tsdest, newdm);
202: DMDestroy(&newdm);
203: return 0;
204: }
206: static PetscErrorCode TSReset_RK_MultirateSplit(TS ts)
207: {
208: TS_RK *rk = (TS_RK *)ts->data;
210: if (rk->subts_slow) {
211: PetscFree(rk->subts_slow->data);
212: rk->subts_slow = NULL;
213: }
214: if (rk->subts_fast) {
215: PetscFree(rk->YdotRHS_fast);
216: PetscFree(rk->YdotRHS_slow);
217: VecDestroy(&rk->X0);
218: TSReset_RK_MultirateSplit(rk->subts_fast);
219: PetscFree(rk->subts_fast->data);
220: rk->subts_fast = NULL;
221: }
222: return 0;
223: }
225: static PetscErrorCode TSInterpolate_RK_MultirateSplit(TS ts, PetscReal itime, Vec X)
226: {
227: TS_RK *rk = (TS_RK *)ts->data;
228: Vec Xslow;
229: PetscInt s = rk->tableau->s, p = rk->tableau->p, i, j;
230: PetscReal h;
231: PetscReal tt, t;
232: PetscScalar *b;
233: const PetscReal *B = rk->tableau->binterp;
237: switch (rk->status) {
238: case TS_STEP_INCOMPLETE:
239: case TS_STEP_PENDING:
240: h = ts->time_step;
241: t = (itime - ts->ptime) / h;
242: break;
243: case TS_STEP_COMPLETE:
244: h = ts->ptime - ts->ptime_prev;
245: t = (itime - ts->ptime) / h + 1; /* In the interval [0,1] */
246: break;
247: default:
248: SETERRQ(PetscObjectComm((PetscObject)ts), PETSC_ERR_PLIB, "Invalid TSStepStatus");
249: }
250: PetscMalloc1(s, &b);
251: for (i = 0; i < s; i++) b[i] = 0;
252: for (j = 0, tt = t; j < p; j++, tt *= t) {
253: for (i = 0; i < s; i++) b[i] += h * B[i * p + j] * tt;
254: }
255: for (i = 0; i < s; i++) VecGetSubVector(rk->YdotRHS[i], rk->is_slow, &rk->YdotRHS_slow[i]);
256: VecGetSubVector(X, rk->is_slow, &Xslow);
257: VecISCopy(rk->X0, rk->is_slow, SCATTER_REVERSE, Xslow);
258: VecMAXPY(Xslow, s, b, rk->YdotRHS_slow);
259: VecRestoreSubVector(X, rk->is_slow, &Xslow);
260: for (i = 0; i < s; i++) VecRestoreSubVector(rk->YdotRHS[i], rk->is_slow, &rk->YdotRHS_slow[i]);
261: PetscFree(b);
262: return 0;
263: }
265: /*
266: This is for partitioned RHS multirate RK method
267: The step completion formula is
269: x1 = x0 + h b^T YdotRHS
271: */
272: static PetscErrorCode TSEvaluateStep_RK_MultirateSplit(TS ts, PetscInt order, Vec X, PetscBool *done)
273: {
274: TS_RK *rk = (TS_RK *)ts->data;
275: RKTableau tab = rk->tableau;
276: Vec Xslow, Xfast; /* subvectors of X which store slow components and fast components respectively */
277: PetscScalar *w = rk->work;
278: PetscReal h = ts->time_step;
279: PetscInt s = tab->s, j;
281: VecCopy(ts->vec_sol, X);
282: if (rk->slow) {
283: for (j = 0; j < s; j++) w[j] = h * tab->b[j];
284: VecGetSubVector(ts->vec_sol, rk->is_slow, &Xslow);
285: VecMAXPY(Xslow, s, w, rk->YdotRHS_slow);
286: VecRestoreSubVector(ts->vec_sol, rk->is_slow, &Xslow);
287: } else {
288: for (j = 0; j < s; j++) w[j] = h / rk->dtratio * tab->b[j];
289: VecGetSubVector(X, rk->is_fast, &Xfast);
290: VecMAXPY(Xfast, s, w, rk->YdotRHS_fast);
291: VecRestoreSubVector(X, rk->is_fast, &Xfast);
292: }
293: return 0;
294: }
296: static PetscErrorCode TSStepRefine_RK_MultirateSplit(TS ts)
297: {
298: TS_RK *rk = (TS_RK *)ts->data;
299: TS subts_fast = rk->subts_fast, currentlevelts;
300: TS_RK *subrk_fast = (TS_RK *)subts_fast->data;
301: RKTableau tab = rk->tableau;
302: Vec *Y = rk->Y;
303: Vec *YdotRHS = rk->YdotRHS, *YdotRHS_fast = rk->YdotRHS_fast;
304: Vec Yfast, Xfast;
305: const PetscInt s = tab->s;
306: const PetscReal *A = tab->A, *c = tab->c;
307: PetscScalar *w = rk->work;
308: PetscInt i, j, k;
309: PetscReal t = ts->ptime, h = ts->time_step;
311: for (k = 0; k < rk->dtratio; k++) {
312: VecGetSubVector(ts->vec_sol, rk->is_fast, &Xfast);
313: for (i = 0; i < s; i++) VecGetSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i]);
314: /* propagate fast component using small time steps */
315: for (i = 0; i < s; i++) {
316: /* stage value for slow components */
317: TSInterpolate_RK_MultirateSplit(rk->ts_root, t + k * h / rk->dtratio + h / rk->dtratio * c[i], Y[i]);
318: currentlevelts = rk->ts_root;
319: while (currentlevelts != ts) { /* all the slow parts need to be interpolated separated */
320: currentlevelts = ((TS_RK *)currentlevelts->data)->subts_fast;
321: TSInterpolate_RK_MultirateSplit(currentlevelts, t + k * h / rk->dtratio + h / rk->dtratio * c[i], Y[i]);
322: }
323: for (j = 0; j < i; j++) w[j] = h / rk->dtratio * A[i * s + j];
324: subrk_fast->stage_time = t + h / rk->dtratio * c[i];
325: TSPreStage(subts_fast, subrk_fast->stage_time);
326: /* stage value for fast components */
327: VecGetSubVector(Y[i], rk->is_fast, &Yfast);
328: VecCopy(Xfast, Yfast);
329: VecMAXPY(Yfast, i, w, YdotRHS_fast);
330: VecRestoreSubVector(Y[i], rk->is_fast, &Yfast);
331: TSPostStage(subts_fast, subrk_fast->stage_time, i, Y);
332: /* compute the stage RHS for fast components */
333: TSComputeRHSFunction(subts_fast, t + k * h * rk->dtratio + h / rk->dtratio * c[i], Y[i], YdotRHS_fast[i]);
334: }
335: VecRestoreSubVector(ts->vec_sol, rk->is_fast, &Xfast);
336: /* update the value of fast components using fast time step */
337: rk->slow = PETSC_FALSE;
338: TSEvaluateStep_RK_MultirateSplit(ts, tab->order, ts->vec_sol, NULL);
339: for (i = 0; i < s; i++) VecRestoreSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i]);
341: if (subrk_fast->subts_fast) {
342: subts_fast->ptime = t + k * h / rk->dtratio;
343: subts_fast->time_step = h / rk->dtratio;
344: TSStepRefine_RK_MultirateSplit(subts_fast);
345: }
346: /* update the fast components of the solution */
347: VecGetSubVector(ts->vec_sol, rk->is_fast, &Xfast);
348: VecISCopy(rk->X0, rk->is_fast, SCATTER_FORWARD, Xfast);
349: VecRestoreSubVector(ts->vec_sol, rk->is_fast, &Xfast);
350: }
351: return 0;
352: }
354: static PetscErrorCode TSStep_RK_MultirateSplit(TS ts)
355: {
356: TS_RK *rk = (TS_RK *)ts->data;
357: RKTableau tab = rk->tableau;
358: Vec *Y = rk->Y, *YdotRHS = rk->YdotRHS;
359: Vec *YdotRHS_fast = rk->YdotRHS_fast, *YdotRHS_slow = rk->YdotRHS_slow;
360: Vec Yslow, Yfast; /* subvectors store the stges of slow components and fast components respectively */
361: const PetscInt s = tab->s;
362: const PetscReal *A = tab->A, *c = tab->c;
363: PetscScalar *w = rk->work;
364: PetscInt i, j;
365: PetscReal next_time_step = ts->time_step, t = ts->ptime, h = ts->time_step;
367: rk->status = TS_STEP_INCOMPLETE;
368: for (i = 0; i < s; i++) {
369: VecGetSubVector(YdotRHS[i], rk->is_slow, &YdotRHS_slow[i]);
370: VecGetSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i]);
371: }
372: VecCopy(ts->vec_sol, rk->X0);
373: /* propagate both slow and fast components using large time steps */
374: for (i = 0; i < s; i++) {
375: rk->stage_time = t + h * c[i];
376: TSPreStage(ts, rk->stage_time);
377: VecCopy(ts->vec_sol, Y[i]);
378: VecGetSubVector(Y[i], rk->is_fast, &Yfast);
379: VecGetSubVector(Y[i], rk->is_slow, &Yslow);
380: for (j = 0; j < i; j++) w[j] = h * A[i * s + j];
381: VecMAXPY(Yfast, i, w, YdotRHS_fast);
382: VecMAXPY(Yslow, i, w, YdotRHS_slow);
383: VecRestoreSubVector(Y[i], rk->is_fast, &Yfast);
384: VecRestoreSubVector(Y[i], rk->is_slow, &Yslow);
385: TSPostStage(ts, rk->stage_time, i, Y);
386: TSComputeRHSFunction(rk->subts_slow, t + h * c[i], Y[i], YdotRHS_slow[i]);
387: TSComputeRHSFunction(rk->subts_fast, t + h * c[i], Y[i], YdotRHS_fast[i]);
388: }
389: rk->slow = PETSC_TRUE;
390: /* update the slow components of the solution using slow time step */
391: TSEvaluateStep_RK_MultirateSplit(ts, tab->order, ts->vec_sol, NULL);
392: /* YdotRHS will be used for interpolation during refinement */
393: for (i = 0; i < s; i++) {
394: VecRestoreSubVector(YdotRHS[i], rk->is_slow, &YdotRHS_slow[i]);
395: VecRestoreSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i]);
396: }
398: TSStepRefine_RK_MultirateSplit(ts);
400: ts->ptime = t + ts->time_step;
401: ts->time_step = next_time_step;
402: rk->status = TS_STEP_COMPLETE;
403: return 0;
404: }
406: static PetscErrorCode TSSetUp_RK_MultirateSplit(TS ts)
407: {
408: TS_RK *rk = (TS_RK *)ts->data, *nextlevelrk, *currentlevelrk;
409: TS nextlevelts;
410: Vec X0;
412: TSRHSSplitGetIS(ts, "slow", &rk->is_slow);
413: TSRHSSplitGetIS(ts, "fast", &rk->is_fast);
415: TSRHSSplitGetSubTS(ts, "slow", &rk->subts_slow);
416: TSRHSSplitGetSubTS(ts, "fast", &rk->subts_fast);
419: VecDuplicate(ts->vec_sol, &X0);
420: /* The TS at each level share the same tableau, work array, solution vector, stage values and stage derivatives */
421: currentlevelrk = rk;
422: while (currentlevelrk->subts_fast) {
423: PetscMalloc1(rk->tableau->s, ¤tlevelrk->YdotRHS_fast);
424: PetscMalloc1(rk->tableau->s, ¤tlevelrk->YdotRHS_slow);
425: PetscObjectReference((PetscObject)X0);
426: currentlevelrk->X0 = X0;
427: currentlevelrk->ts_root = ts;
429: /* set up the ts for the slow part */
430: nextlevelts = currentlevelrk->subts_slow;
431: PetscNew(&nextlevelrk);
432: nextlevelrk->tableau = rk->tableau;
433: nextlevelrk->work = rk->work;
434: nextlevelrk->Y = rk->Y;
435: nextlevelrk->YdotRHS = rk->YdotRHS;
436: nextlevelts->data = (void *)nextlevelrk;
437: TSCopyDM(ts, nextlevelts);
438: TSSetSolution(nextlevelts, ts->vec_sol);
440: /* set up the ts for the fast part */
441: nextlevelts = currentlevelrk->subts_fast;
442: PetscNew(&nextlevelrk);
443: nextlevelrk->tableau = rk->tableau;
444: nextlevelrk->work = rk->work;
445: nextlevelrk->Y = rk->Y;
446: nextlevelrk->YdotRHS = rk->YdotRHS;
447: nextlevelrk->dtratio = rk->dtratio;
448: TSRHSSplitGetIS(nextlevelts, "slow", &nextlevelrk->is_slow);
449: TSRHSSplitGetSubTS(nextlevelts, "slow", &nextlevelrk->subts_slow);
450: TSRHSSplitGetIS(nextlevelts, "fast", &nextlevelrk->is_fast);
451: TSRHSSplitGetSubTS(nextlevelts, "fast", &nextlevelrk->subts_fast);
452: nextlevelts->data = (void *)nextlevelrk;
453: TSCopyDM(ts, nextlevelts);
454: TSSetSolution(nextlevelts, ts->vec_sol);
456: currentlevelrk = nextlevelrk;
457: }
458: VecDestroy(&X0);
460: ts->ops->step = TSStep_RK_MultirateSplit;
461: ts->ops->evaluatestep = TSEvaluateStep_RK_MultirateSplit;
462: ts->ops->interpolate = TSInterpolate_RK_MultirateSplit;
463: return 0;
464: }
466: PetscErrorCode TSRKSetMultirate_RK(TS ts, PetscBool use_multirate)
467: {
468: TS_RK *rk = (TS_RK *)ts->data;
471: rk->use_multirate = use_multirate;
472: if (use_multirate) {
473: rk->dtratio = 2;
474: PetscObjectComposeFunction((PetscObject)ts, "TSSetUp_RK_MultirateSplit_C", TSSetUp_RK_MultirateSplit);
475: PetscObjectComposeFunction((PetscObject)ts, "TSReset_RK_MultirateSplit_C", TSReset_RK_MultirateSplit);
476: PetscObjectComposeFunction((PetscObject)ts, "TSSetUp_RK_MultirateNonsplit_C", TSSetUp_RK_MultirateNonsplit);
477: PetscObjectComposeFunction((PetscObject)ts, "TSReset_RK_MultirateNonsplit_C", TSReset_RK_MultirateNonsplit);
478: } else {
479: rk->dtratio = 0;
480: PetscObjectComposeFunction((PetscObject)ts, "TSSetUp_RK_MultirateSplit_C", NULL);
481: PetscObjectComposeFunction((PetscObject)ts, "TSReset_RK_MultirateSplit_C", NULL);
482: PetscObjectComposeFunction((PetscObject)ts, "TSSetUp_RK_MultirateNonsplit_C", NULL);
483: PetscObjectComposeFunction((PetscObject)ts, "TSReset_RK_MultirateNonsplit_C", NULL);
484: }
485: return 0;
486: }
488: PetscErrorCode TSRKGetMultirate_RK(TS ts, PetscBool *use_multirate)
489: {
490: TS_RK *rk = (TS_RK *)ts->data;
493: *use_multirate = rk->use_multirate;
494: return 0;
495: }