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 partioned 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++) {
44: b[i] += h * B[i*p+j] * tt;
45: }
46: }
47: VecCopy(rk->X0,X);
48: VecMAXPY(X,s,b,rk->YdotRHS_slow);
49: PetscFree(b);
50: return 0;
51: }
53: static PetscErrorCode TSStepRefine_RK_MultirateNonsplit(TS ts)
54: {
55: TS previousts,subts;
56: TS_RK *rk = (TS_RK*)ts->data;
57: RKTableau tab = rk->tableau;
58: Vec *Y = rk->Y,*YdotRHS = rk->YdotRHS;
59: Vec vec_fast,subvec_fast;
60: const PetscInt s = tab->s;
61: const PetscReal *A = tab->A,*c = tab->c;
62: PetscScalar *w = rk->work;
63: PetscInt i,j,k;
64: PetscReal t = ts->ptime,h = ts->time_step;
66: VecDuplicate(ts->vec_sol,&vec_fast);
67: previousts = rk->subts_current;
68: TSRHSSplitGetSubTS(rk->subts_current,"fast",&subts);
69: TSRHSSplitGetSubTS(subts,"fast",&subts);
70: for (k=0; k<rk->dtratio; k++) {
71: for (i=0; i<s; i++) {
72: TSInterpolate_RK_MultirateNonsplit(ts,t+k*h/rk->dtratio+h/rk->dtratio*c[i],Y[i]);
73: for (j=0; j<i; j++) w[j] = h/rk->dtratio*A[i*s+j];
74: /* 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 */
75: VecCopy(ts->vec_sol,vec_fast);
76: VecMAXPY(vec_fast,i,w,YdotRHS);
77: /* update the fast components in the stage value */
78: VecGetSubVector(vec_fast,rk->is_fast,&subvec_fast);
79: VecISCopy(Y[i],rk->is_fast,SCATTER_FORWARD,subvec_fast);
80: VecRestoreSubVector(vec_fast,rk->is_fast,&subvec_fast);
81: /* compute the stage RHS */
82: TSComputeRHSFunction(ts,t+k*h/rk->dtratio+h/rk->dtratio*c[i],Y[i],YdotRHS[i]);
83: }
84: VecCopy(ts->vec_sol,vec_fast);
85: TSEvaluateStep(ts,tab->order,vec_fast,NULL);
86: /* update the fast components in the solution */
87: VecGetSubVector(vec_fast,rk->is_fast,&subvec_fast);
88: VecISCopy(ts->vec_sol,rk->is_fast,SCATTER_FORWARD,subvec_fast);
89: VecRestoreSubVector(vec_fast,rk->is_fast,&subvec_fast);
91: if (subts) {
92: Vec *YdotRHS_copy;
93: VecDuplicateVecs(ts->vec_sol,s,&YdotRHS_copy);
94: rk->subts_current = rk->subts_fast;
95: ts->ptime = t+k*h/rk->dtratio;
96: ts->time_step = h/rk->dtratio;
97: TSRHSSplitGetIS(rk->subts_current,"fast",&rk->is_fast);
98: for (i=0; i<s; i++) {
99: VecCopy(rk->YdotRHS_slow[i],YdotRHS_copy[i]);
100: VecCopy(YdotRHS[i],rk->YdotRHS_slow[i]);
101: }
103: TSStepRefine_RK_MultirateNonsplit(ts);
105: rk->subts_current = previousts;
106: ts->ptime = t;
107: ts->time_step = h;
108: TSRHSSplitGetIS(previousts,"fast",&rk->is_fast);
109: for (i=0; i<s; i++) {
110: VecCopy(YdotRHS_copy[i],rk->YdotRHS_slow[i]);
111: }
112: VecDestroyVecs(s,&YdotRHS_copy);
113: }
114: }
115: VecDestroy(&vec_fast);
116: return 0;
117: }
119: static PetscErrorCode TSStep_RK_MultirateNonsplit(TS ts)
120: {
121: TS_RK *rk = (TS_RK*)ts->data;
122: RKTableau tab = rk->tableau;
123: Vec *Y = rk->Y,*YdotRHS = rk->YdotRHS,*YdotRHS_slow = rk->YdotRHS_slow;
124: Vec stage_slow,sol_slow; /* vectors store the slow components */
125: Vec subvec_slow; /* sub vector to store the slow components */
126: IS is_slow = rk->is_slow;
127: const PetscInt s = tab->s;
128: const PetscReal *A = tab->A,*c = tab->c;
129: PetscScalar *w = rk->work;
130: PetscInt i,j,dtratio = rk->dtratio;
131: PetscReal next_time_step = ts->time_step,t = ts->ptime,h = ts->time_step;
133: rk->status = TS_STEP_INCOMPLETE;
134: VecDuplicate(ts->vec_sol,&stage_slow);
135: VecDuplicate(ts->vec_sol,&sol_slow);
136: VecCopy(ts->vec_sol,rk->X0);
137: for (i=0; i<s; i++) {
138: rk->stage_time = t + h*c[i];
139: TSPreStage(ts,rk->stage_time);
140: VecCopy(ts->vec_sol,Y[i]);
141: for (j=0; j<i; j++) w[j] = h*A[i*s+j];
142: VecMAXPY(Y[i],i,w,YdotRHS_slow);
143: TSPostStage(ts,rk->stage_time,i,Y);
144: /* compute the stage RHS */
145: TSComputeRHSFunction(ts,t+h*c[i],Y[i],YdotRHS_slow[i]);
146: }
147: /* update the slow components in the solution */
148: rk->YdotRHS = YdotRHS_slow;
149: rk->dtratio = 1;
150: TSEvaluateStep(ts,tab->order,sol_slow,NULL);
151: rk->dtratio = dtratio;
152: rk->YdotRHS = YdotRHS;
153: /* update the slow components in the solution */
154: VecGetSubVector(sol_slow,is_slow,&subvec_slow);
155: VecISCopy(ts->vec_sol,is_slow,SCATTER_FORWARD,subvec_slow);
156: VecRestoreSubVector(sol_slow,is_slow,&subvec_slow);
158: rk->subts_current = ts;
159: rk->ptime = t;
160: rk->time_step = h;
161: TSStepRefine_RK_MultirateNonsplit(ts);
163: ts->ptime = t + ts->time_step;
164: ts->time_step = next_time_step;
165: rk->status = TS_STEP_COMPLETE;
167: /* free memory of work vectors */
168: VecDestroy(&stage_slow);
169: VecDestroy(&sol_slow);
170: return 0;
171: }
173: static PetscErrorCode TSSetUp_RK_MultirateNonsplit(TS ts)
174: {
175: TS_RK *rk = (TS_RK*)ts->data;
176: RKTableau tab = rk->tableau;
178: TSRHSSplitGetIS(ts,"slow",&rk->is_slow);
179: TSRHSSplitGetIS(ts,"fast",&rk->is_fast);
181: TSRHSSplitGetSubTS(ts,"slow",&rk->subts_slow);
182: TSRHSSplitGetSubTS(ts,"fast",&rk->subts_fast);
184: VecDuplicate(ts->vec_sol,&rk->X0);
185: VecDuplicateVecs(ts->vec_sol,tab->s,&rk->YdotRHS_slow);
186: rk->subts_current = rk->subts_fast;
188: ts->ops->step = TSStep_RK_MultirateNonsplit;
189: ts->ops->interpolate = TSInterpolate_RK_MultirateNonsplit;
190: return 0;
191: }
193: /*
194: Copy DM from tssrc to tsdest, while keeping the original DMTS and DMSNES in tsdest.
195: */
196: static PetscErrorCode TSCopyDM(TS tssrc,TS tsdest)
197: {
198: DM newdm,dmsrc,dmdest;
200: TSGetDM(tssrc,&dmsrc);
201: DMClone(dmsrc,&newdm);
202: TSGetDM(tsdest,&dmdest);
203: DMCopyDMTS(dmdest,newdm);
204: DMCopyDMSNES(dmdest,newdm);
205: TSSetDM(tsdest,newdm);
206: DMDestroy(&newdm);
207: return 0;
208: }
210: static PetscErrorCode TSReset_RK_MultirateSplit(TS ts)
211: {
212: TS_RK *rk = (TS_RK*)ts->data;
214: if (rk->subts_slow) {
215: PetscFree(rk->subts_slow->data);
216: rk->subts_slow = NULL;
217: }
218: if (rk->subts_fast) {
219: PetscFree(rk->YdotRHS_fast);
220: PetscFree(rk->YdotRHS_slow);
221: VecDestroy(&rk->X0);
222: TSReset_RK_MultirateSplit(rk->subts_fast);
223: PetscFree(rk->subts_fast->data);
224: rk->subts_fast = NULL;
225: }
226: return 0;
227: }
229: static PetscErrorCode TSInterpolate_RK_MultirateSplit(TS ts,PetscReal itime,Vec X)
230: {
231: TS_RK *rk = (TS_RK*)ts->data;
232: Vec Xslow;
233: PetscInt s = rk->tableau->s,p = rk->tableau->p,i,j;
234: PetscReal h;
235: PetscReal tt,t;
236: PetscScalar *b;
237: const PetscReal *B = rk->tableau->binterp;
241: switch (rk->status) {
242: case TS_STEP_INCOMPLETE:
243: case TS_STEP_PENDING:
244: h = ts->time_step;
245: t = (itime - ts->ptime)/h;
246: break;
247: case TS_STEP_COMPLETE:
248: h = ts->ptime - ts->ptime_prev;
249: t = (itime - ts->ptime)/h + 1; /* In the interval [0,1] */
250: break;
251: default: SETERRQ(PetscObjectComm((PetscObject)ts),PETSC_ERR_PLIB,"Invalid TSStepStatus");
252: }
253: PetscMalloc1(s,&b);
254: for (i=0; i<s; i++) b[i] = 0;
255: for (j=0,tt=t; j<p; j++,tt*=t) {
256: for (i=0; i<s; i++) {
257: b[i] += h * B[i*p+j] * tt;
258: }
259: }
260: for (i=0; i<s; i++) {
261: VecGetSubVector(rk->YdotRHS[i],rk->is_slow,&rk->YdotRHS_slow[i]);
262: }
263: VecGetSubVector(X,rk->is_slow,&Xslow);
264: VecISCopy(rk->X0,rk->is_slow,SCATTER_REVERSE,Xslow);
265: VecMAXPY(Xslow,s,b,rk->YdotRHS_slow);
266: VecRestoreSubVector(X,rk->is_slow,&Xslow);
267: for (i=0; i<s; i++) {
268: VecRestoreSubVector(rk->YdotRHS[i],rk->is_slow,&rk->YdotRHS_slow[i]);
269: }
270: PetscFree(b);
271: return 0;
272: }
274: /*
275: This is for partitioned RHS multirate RK method
276: The step completion formula is
278: x1 = x0 + h b^T YdotRHS
280: */
281: static PetscErrorCode TSEvaluateStep_RK_MultirateSplit(TS ts,PetscInt order,Vec X,PetscBool *done)
282: {
283: TS_RK *rk = (TS_RK*)ts->data;
284: RKTableau tab = rk->tableau;
285: Vec Xslow,Xfast; /* subvectors of X which store slow components and fast components respectively */
286: PetscScalar *w = rk->work;
287: PetscReal h = ts->time_step;
288: PetscInt s = tab->s,j;
290: VecCopy(ts->vec_sol,X);
291: if (rk->slow) {
292: for (j=0; j<s; j++) w[j] = h*tab->b[j];
293: VecGetSubVector(ts->vec_sol,rk->is_slow,&Xslow);
294: VecMAXPY(Xslow,s,w,rk->YdotRHS_slow);
295: VecRestoreSubVector(ts->vec_sol,rk->is_slow,&Xslow);
296: } else {
297: for (j=0; j<s; j++) w[j] = h/rk->dtratio*tab->b[j];
298: VecGetSubVector(X,rk->is_fast,&Xfast);
299: VecMAXPY(Xfast,s,w,rk->YdotRHS_fast);
300: VecRestoreSubVector(X,rk->is_fast,&Xfast);
301: }
302: return 0;
303: }
305: static PetscErrorCode TSStepRefine_RK_MultirateSplit(TS ts)
306: {
307: TS_RK *rk = (TS_RK*)ts->data;
308: TS subts_fast = rk->subts_fast,currentlevelts;
309: TS_RK *subrk_fast = (TS_RK*)subts_fast->data;
310: RKTableau tab = rk->tableau;
311: Vec *Y = rk->Y;
312: Vec *YdotRHS = rk->YdotRHS,*YdotRHS_fast = rk->YdotRHS_fast;
313: Vec Yfast,Xfast;
314: const PetscInt s = tab->s;
315: const PetscReal *A = tab->A,*c = tab->c;
316: PetscScalar *w = rk->work;
317: PetscInt i,j,k;
318: PetscReal t = ts->ptime,h = ts->time_step;
320: for (k=0; k<rk->dtratio; k++) {
321: VecGetSubVector(ts->vec_sol,rk->is_fast,&Xfast);
322: for (i=0; i<s; i++) {
323: VecGetSubVector(YdotRHS[i],rk->is_fast,&YdotRHS_fast[i]);
324: }
325: /* propagate fast component using small time steps */
326: for (i=0; i<s; i++) {
327: /* stage value for slow components */
328: TSInterpolate_RK_MultirateSplit(rk->ts_root,t+k*h/rk->dtratio+h/rk->dtratio*c[i],Y[i]);
329: currentlevelts = rk->ts_root;
330: while (currentlevelts != ts) { /* all the slow parts need to be interpolated separated */
331: currentlevelts = ((TS_RK*)currentlevelts->data)->subts_fast;
332: TSInterpolate_RK_MultirateSplit(currentlevelts,t+k*h/rk->dtratio+h/rk->dtratio*c[i],Y[i]);
333: }
334: for (j=0; j<i; j++) w[j] = h/rk->dtratio*A[i*s+j];
335: subrk_fast->stage_time = t + h/rk->dtratio*c[i];
336: TSPreStage(subts_fast,subrk_fast->stage_time);
337: /* stage value for fast components */
338: VecGetSubVector(Y[i],rk->is_fast,&Yfast);
339: VecCopy(Xfast,Yfast);
340: VecMAXPY(Yfast,i,w,YdotRHS_fast);
341: VecRestoreSubVector(Y[i],rk->is_fast,&Yfast);
342: TSPostStage(subts_fast,subrk_fast->stage_time,i,Y);
343: /* compute the stage RHS for fast components */
344: TSComputeRHSFunction(subts_fast,t+k*h*rk->dtratio+h/rk->dtratio*c[i],Y[i],YdotRHS_fast[i]);
345: }
346: VecRestoreSubVector(ts->vec_sol,rk->is_fast,&Xfast);
347: /* update the value of fast components using fast time step */
348: rk->slow = PETSC_FALSE;
349: TSEvaluateStep_RK_MultirateSplit(ts,tab->order,ts->vec_sol,NULL);
350: for (i=0; i<s; i++) {
351: VecRestoreSubVector(YdotRHS[i],rk->is_fast,&YdotRHS_fast[i]);
352: }
354: if (subrk_fast->subts_fast) {
355: subts_fast->ptime = t+k*h/rk->dtratio;
356: subts_fast->time_step = h/rk->dtratio;
357: TSStepRefine_RK_MultirateSplit(subts_fast);
358: }
359: /* update the fast components of the solution */
360: VecGetSubVector(ts->vec_sol,rk->is_fast,&Xfast);
361: VecISCopy(rk->X0,rk->is_fast,SCATTER_FORWARD,Xfast);
362: VecRestoreSubVector(ts->vec_sol,rk->is_fast,&Xfast);
363: }
364: return 0;
365: }
367: static PetscErrorCode TSStep_RK_MultirateSplit(TS ts)
368: {
369: TS_RK *rk = (TS_RK*)ts->data;
370: RKTableau tab = rk->tableau;
371: Vec *Y = rk->Y,*YdotRHS = rk->YdotRHS;
372: Vec *YdotRHS_fast = rk->YdotRHS_fast,*YdotRHS_slow = rk->YdotRHS_slow;
373: Vec Yslow,Yfast; /* subvectors store the stges of slow components and fast components respectively */
374: const PetscInt s = tab->s;
375: const PetscReal *A = tab->A,*c = tab->c;
376: PetscScalar *w = rk->work;
377: PetscInt i,j;
378: PetscReal next_time_step = ts->time_step,t = ts->ptime,h = ts->time_step;
380: rk->status = TS_STEP_INCOMPLETE;
381: for (i=0; i<s; i++) {
382: VecGetSubVector(YdotRHS[i],rk->is_slow,&YdotRHS_slow[i]);
383: VecGetSubVector(YdotRHS[i],rk->is_fast,&YdotRHS_fast[i]);
384: }
385: VecCopy(ts->vec_sol,rk->X0);
386: /* propagate both slow and fast components using large time steps */
387: for (i=0; i<s; i++) {
388: rk->stage_time = t + h*c[i];
389: TSPreStage(ts,rk->stage_time);
390: VecCopy(ts->vec_sol,Y[i]);
391: VecGetSubVector(Y[i],rk->is_fast,&Yfast);
392: VecGetSubVector(Y[i],rk->is_slow,&Yslow);
393: for (j=0; j<i; j++) w[j] = h*A[i*s+j];
394: VecMAXPY(Yfast,i,w,YdotRHS_fast);
395: VecMAXPY(Yslow,i,w,YdotRHS_slow);
396: VecRestoreSubVector(Y[i],rk->is_fast,&Yfast);
397: VecRestoreSubVector(Y[i],rk->is_slow,&Yslow);
398: TSPostStage(ts,rk->stage_time,i,Y);
399: TSComputeRHSFunction(rk->subts_slow,t+h*c[i],Y[i],YdotRHS_slow[i]);
400: TSComputeRHSFunction(rk->subts_fast,t+h*c[i],Y[i],YdotRHS_fast[i]);
401: }
402: rk->slow = PETSC_TRUE;
403: /* update the slow components of the solution using slow time step */
404: TSEvaluateStep_RK_MultirateSplit(ts,tab->order,ts->vec_sol,NULL);
405: /* YdotRHS will be used for interpolation during refinement */
406: for (i=0; i<s; i++) {
407: VecRestoreSubVector(YdotRHS[i],rk->is_slow,&YdotRHS_slow[i]);
408: VecRestoreSubVector(YdotRHS[i],rk->is_fast,&YdotRHS_fast[i]);
409: }
411: TSStepRefine_RK_MultirateSplit(ts);
413: ts->ptime = t + ts->time_step;
414: ts->time_step = next_time_step;
415: rk->status = TS_STEP_COMPLETE;
416: return 0;
417: }
419: static PetscErrorCode TSSetUp_RK_MultirateSplit(TS ts)
420: {
421: TS_RK *rk = (TS_RK*)ts->data,*nextlevelrk,*currentlevelrk;
422: TS nextlevelts;
423: Vec X0;
425: TSRHSSplitGetIS(ts,"slow",&rk->is_slow);
426: TSRHSSplitGetIS(ts,"fast",&rk->is_fast);
428: TSRHSSplitGetSubTS(ts,"slow",&rk->subts_slow);
429: TSRHSSplitGetSubTS(ts,"fast",&rk->subts_fast);
432: VecDuplicate(ts->vec_sol,&X0);
433: /* The TS at each level share the same tableau, work array, solution vector, stage values and stage derivatives */
434: currentlevelrk = rk;
435: while (currentlevelrk->subts_fast) {
436: PetscMalloc1(rk->tableau->s,¤tlevelrk->YdotRHS_fast);
437: PetscMalloc1(rk->tableau->s,¤tlevelrk->YdotRHS_slow);
438: PetscObjectReference((PetscObject)X0);
439: currentlevelrk->X0 = X0;
440: currentlevelrk->ts_root = ts;
442: /* set up the ts for the slow part */
443: nextlevelts = currentlevelrk->subts_slow;
444: PetscNewLog(nextlevelts,&nextlevelrk);
445: nextlevelrk->tableau = rk->tableau;
446: nextlevelrk->work = rk->work;
447: nextlevelrk->Y = rk->Y;
448: nextlevelrk->YdotRHS = rk->YdotRHS;
449: nextlevelts->data = (void*)nextlevelrk;
450: TSCopyDM(ts,nextlevelts);
451: TSSetSolution(nextlevelts,ts->vec_sol);
453: /* set up the ts for the fast part */
454: nextlevelts = currentlevelrk->subts_fast;
455: PetscNewLog(nextlevelts,&nextlevelrk);
456: nextlevelrk->tableau = rk->tableau;
457: nextlevelrk->work = rk->work;
458: nextlevelrk->Y = rk->Y;
459: nextlevelrk->YdotRHS = rk->YdotRHS;
460: nextlevelrk->dtratio = rk->dtratio;
461: TSRHSSplitGetIS(nextlevelts,"slow",&nextlevelrk->is_slow);
462: TSRHSSplitGetSubTS(nextlevelts,"slow",&nextlevelrk->subts_slow);
463: TSRHSSplitGetIS(nextlevelts,"fast",&nextlevelrk->is_fast);
464: TSRHSSplitGetSubTS(nextlevelts,"fast",&nextlevelrk->subts_fast);
465: nextlevelts->data = (void*)nextlevelrk;
466: TSCopyDM(ts,nextlevelts);
467: TSSetSolution(nextlevelts,ts->vec_sol);
469: currentlevelrk = nextlevelrk;
470: }
471: VecDestroy(&X0);
473: ts->ops->step = TSStep_RK_MultirateSplit;
474: ts->ops->evaluatestep = TSEvaluateStep_RK_MultirateSplit;
475: ts->ops->interpolate = TSInterpolate_RK_MultirateSplit;
476: return 0;
477: }
479: PetscErrorCode TSRKSetMultirate_RK(TS ts,PetscBool use_multirate)
480: {
481: TS_RK *rk = (TS_RK*)ts->data;
484: rk->use_multirate = use_multirate;
485: if (use_multirate) {
486: rk->dtratio = 2;
487: PetscObjectComposeFunction((PetscObject)ts,"TSSetUp_RK_MultirateSplit_C",TSSetUp_RK_MultirateSplit);
488: PetscObjectComposeFunction((PetscObject)ts,"TSReset_RK_MultirateSplit_C",TSReset_RK_MultirateSplit);
489: PetscObjectComposeFunction((PetscObject)ts,"TSSetUp_RK_MultirateNonsplit_C",TSSetUp_RK_MultirateNonsplit);
490: PetscObjectComposeFunction((PetscObject)ts,"TSReset_RK_MultirateNonsplit_C",TSReset_RK_MultirateNonsplit);
491: } else {
492: rk->dtratio = 0;
493: PetscObjectComposeFunction((PetscObject)ts,"TSSetUp_RK_MultirateSplit_C",NULL);
494: PetscObjectComposeFunction((PetscObject)ts,"TSReset_RK_MultirateSplit_C",NULL);
495: PetscObjectComposeFunction((PetscObject)ts,"TSSetUp_RK_MultirateNonsplit_C",NULL);
496: PetscObjectComposeFunction((PetscObject)ts,"TSReset_RK_MultirateNonsplit_C",NULL);
497: }
498: return 0;
499: }
501: PetscErrorCode TSRKGetMultirate_RK(TS ts,PetscBool *use_multirate)
502: {
503: TS_RK *rk = (TS_RK*)ts->data;
506: *use_multirate = rk->use_multirate;
507: return 0;
508: }