Lines Matching refs:rk
21 TS_RK *rk = (TS_RK *)ts->data; in TSReset_RK_MultirateNonsplit() local
22 RKTableau tab = rk->tableau; in TSReset_RK_MultirateNonsplit()
25 PetscCall(VecDestroy(&rk->X0)); in TSReset_RK_MultirateNonsplit()
26 PetscCall(VecDestroyVecs(tab->s, &rk->YdotRHS_slow)); in TSReset_RK_MultirateNonsplit()
32 TS_RK *rk = (TS_RK *)ts->data; in TSInterpolate_RK_MultirateNonsplit() local
33 PetscInt s = rk->tableau->s, p = rk->tableau->p, i, j; in TSInterpolate_RK_MultirateNonsplit()
37 const PetscReal *B = rk->tableau->binterp; in TSInterpolate_RK_MultirateNonsplit()
40 …tscObject)ts), PETSC_ERR_SUP, "TSRK %s does not have an interpolation formula", rk->tableau->name); in TSInterpolate_RK_MultirateNonsplit()
41 t = (itime - rk->ptime) / h; in TSInterpolate_RK_MultirateNonsplit()
47 PetscCall(VecCopy(rk->X0, X)); in TSInterpolate_RK_MultirateNonsplit()
48 PetscCall(VecMAXPY(X, s, b, rk->YdotRHS_slow)); in TSInterpolate_RK_MultirateNonsplit()
56 TS_RK *rk = (TS_RK *)ts->data; in TSStepRefine_RK_MultirateNonsplit() local
57 RKTableau tab = rk->tableau; in TSStepRefine_RK_MultirateNonsplit()
58 Vec *Y = rk->Y, *YdotRHS = rk->YdotRHS; in TSStepRefine_RK_MultirateNonsplit()
62 PetscScalar *w = rk->work; in TSStepRefine_RK_MultirateNonsplit()
68 previousts = rk->subts_current; in TSStepRefine_RK_MultirateNonsplit()
69 PetscCall(TSRHSSplitGetSubTS(rk->subts_current, "fast", &subts)); in TSStepRefine_RK_MultirateNonsplit()
71 for (k = 0; k < rk->dtratio; k++) { in TSStepRefine_RK_MultirateNonsplit()
73 …PetscCall(TSInterpolate_RK_MultirateNonsplit(ts, t + k * h / rk->dtratio + h / rk->dtratio * c[i],… in TSStepRefine_RK_MultirateNonsplit()
74 for (j = 0; j < i; j++) w[j] = h / rk->dtratio * A[i * s + j]; in TSStepRefine_RK_MultirateNonsplit()
79 PetscCall(VecGetSubVector(vec_fast, rk->is_fast, &subvec_fast)); in TSStepRefine_RK_MultirateNonsplit()
80 PetscCall(VecISCopy(Y[i], rk->is_fast, SCATTER_FORWARD, subvec_fast)); in TSStepRefine_RK_MultirateNonsplit()
81 PetscCall(VecRestoreSubVector(vec_fast, rk->is_fast, &subvec_fast)); in TSStepRefine_RK_MultirateNonsplit()
83 …PetscCall(TSComputeRHSFunction(ts, t + k * h / rk->dtratio + h / rk->dtratio * c[i], Y[i], YdotRHS… in TSStepRefine_RK_MultirateNonsplit()
88 PetscCall(VecGetSubVector(vec_fast, rk->is_fast, &subvec_fast)); in TSStepRefine_RK_MultirateNonsplit()
89 PetscCall(VecISCopy(ts->vec_sol, rk->is_fast, SCATTER_FORWARD, subvec_fast)); in TSStepRefine_RK_MultirateNonsplit()
90 PetscCall(VecRestoreSubVector(vec_fast, rk->is_fast, &subvec_fast)); in TSStepRefine_RK_MultirateNonsplit()
95 rk->subts_current = rk->subts_fast; in TSStepRefine_RK_MultirateNonsplit()
96 ts->ptime = t + k * h / rk->dtratio; in TSStepRefine_RK_MultirateNonsplit()
97 ts->time_step = h / rk->dtratio; in TSStepRefine_RK_MultirateNonsplit()
98 PetscCall(TSRHSSplitGetIS(rk->subts_current, "fast", &rk->is_fast)); in TSStepRefine_RK_MultirateNonsplit()
100 PetscCall(VecCopy(rk->YdotRHS_slow[i], YdotRHS_copy[i])); in TSStepRefine_RK_MultirateNonsplit()
101 PetscCall(VecCopy(YdotRHS[i], rk->YdotRHS_slow[i])); in TSStepRefine_RK_MultirateNonsplit()
106 rk->subts_current = previousts; in TSStepRefine_RK_MultirateNonsplit()
109 PetscCall(TSRHSSplitGetIS(previousts, "fast", &rk->is_fast)); in TSStepRefine_RK_MultirateNonsplit()
110 for (i = 0; i < s; i++) PetscCall(VecCopy(YdotRHS_copy[i], rk->YdotRHS_slow[i])); in TSStepRefine_RK_MultirateNonsplit()
120 TS_RK *rk = (TS_RK *)ts->data; in TSStep_RK_MultirateNonsplit() local
121 RKTableau tab = rk->tableau; in TSStep_RK_MultirateNonsplit()
122 Vec *Y = rk->Y, *YdotRHS = rk->YdotRHS, *YdotRHS_slow = rk->YdotRHS_slow; in TSStep_RK_MultirateNonsplit()
125 IS is_slow = rk->is_slow; in TSStep_RK_MultirateNonsplit()
128 PetscScalar *w = rk->work; in TSStep_RK_MultirateNonsplit()
129 PetscInt i, j, dtratio = rk->dtratio; in TSStep_RK_MultirateNonsplit()
133 rk->status = TS_STEP_INCOMPLETE; in TSStep_RK_MultirateNonsplit()
136 PetscCall(VecCopy(ts->vec_sol, rk->X0)); in TSStep_RK_MultirateNonsplit()
138 rk->stage_time = t + h * c[i]; in TSStep_RK_MultirateNonsplit()
139 PetscCall(TSPreStage(ts, rk->stage_time)); in TSStep_RK_MultirateNonsplit()
143 PetscCall(TSPostStage(ts, rk->stage_time, i, Y)); in TSStep_RK_MultirateNonsplit()
148 rk->YdotRHS = YdotRHS_slow; in TSStep_RK_MultirateNonsplit()
149 rk->dtratio = 1; in TSStep_RK_MultirateNonsplit()
151 rk->dtratio = dtratio; in TSStep_RK_MultirateNonsplit()
152 rk->YdotRHS = YdotRHS; in TSStep_RK_MultirateNonsplit()
158 rk->subts_current = ts; in TSStep_RK_MultirateNonsplit()
159 rk->ptime = t; in TSStep_RK_MultirateNonsplit()
160 rk->time_step = h; in TSStep_RK_MultirateNonsplit()
165 rk->status = TS_STEP_COMPLETE; in TSStep_RK_MultirateNonsplit()
175 TS_RK *rk = (TS_RK *)ts->data; in TSSetUp_RK_MultirateNonsplit() local
176 RKTableau tab = rk->tableau; in TSSetUp_RK_MultirateNonsplit()
179 PetscCall(TSRHSSplitGetIS(ts, "slow", &rk->is_slow)); in TSSetUp_RK_MultirateNonsplit()
180 PetscCall(TSRHSSplitGetIS(ts, "fast", &rk->is_fast)); in TSSetUp_RK_MultirateNonsplit()
181 …PetscCheck(rk->is_slow && rk->is_fast, PetscObjectComm((PetscObject)ts), PETSC_ERR_USER, "Must set… in TSSetUp_RK_MultirateNonsplit()
182 PetscCall(TSRHSSplitGetSubTS(ts, "slow", &rk->subts_slow)); in TSSetUp_RK_MultirateNonsplit()
183 PetscCall(TSRHSSplitGetSubTS(ts, "fast", &rk->subts_fast)); in TSSetUp_RK_MultirateNonsplit()
184 …PetscCheck(rk->subts_slow && rk->subts_fast, PetscObjectComm((PetscObject)ts), PETSC_ERR_USER, "Mu… in TSSetUp_RK_MultirateNonsplit()
185 PetscCall(VecDuplicate(ts->vec_sol, &rk->X0)); in TSSetUp_RK_MultirateNonsplit()
186 PetscCall(VecDuplicateVecs(ts->vec_sol, tab->s, &rk->YdotRHS_slow)); in TSSetUp_RK_MultirateNonsplit()
187 rk->subts_current = rk->subts_fast; in TSSetUp_RK_MultirateNonsplit()
214 TS_RK *rk = (TS_RK *)ts->data; in TSReset_RK_MultirateSplit() local
217 if (rk->subts_slow) { in TSReset_RK_MultirateSplit()
218 PetscCall(PetscFree(rk->subts_slow->data)); in TSReset_RK_MultirateSplit()
219 rk->subts_slow = NULL; in TSReset_RK_MultirateSplit()
221 if (rk->subts_fast) { in TSReset_RK_MultirateSplit()
222 PetscCall(PetscFree(rk->YdotRHS_fast)); in TSReset_RK_MultirateSplit()
223 PetscCall(PetscFree(rk->YdotRHS_slow)); in TSReset_RK_MultirateSplit()
224 PetscCall(VecDestroy(&rk->X0)); in TSReset_RK_MultirateSplit()
225 PetscCall(TSReset_RK_MultirateSplit(rk->subts_fast)); in TSReset_RK_MultirateSplit()
226 PetscCall(PetscFree(rk->subts_fast->data)); in TSReset_RK_MultirateSplit()
227 rk->subts_fast = NULL; in TSReset_RK_MultirateSplit()
234 TS_RK *rk = (TS_RK *)ts->data; in TSInterpolate_RK_MultirateSplit() local
236 PetscInt s = rk->tableau->s, p = rk->tableau->p, i, j; in TSInterpolate_RK_MultirateSplit()
240 const PetscReal *B = rk->tableau->binterp; in TSInterpolate_RK_MultirateSplit()
243 …tscObject)ts), PETSC_ERR_SUP, "TSRK %s does not have an interpolation formula", rk->tableau->name); in TSInterpolate_RK_MultirateSplit()
245 switch (rk->status) { in TSInterpolate_RK_MultirateSplit()
263 …for (i = 0; i < s; i++) PetscCall(VecGetSubVector(rk->YdotRHS[i], rk->is_slow, &rk->YdotRHS_slow[i… in TSInterpolate_RK_MultirateSplit()
264 PetscCall(VecGetSubVector(X, rk->is_slow, &Xslow)); in TSInterpolate_RK_MultirateSplit()
265 PetscCall(VecISCopy(rk->X0, rk->is_slow, SCATTER_REVERSE, Xslow)); in TSInterpolate_RK_MultirateSplit()
266 PetscCall(VecMAXPY(Xslow, s, b, rk->YdotRHS_slow)); in TSInterpolate_RK_MultirateSplit()
267 PetscCall(VecRestoreSubVector(X, rk->is_slow, &Xslow)); in TSInterpolate_RK_MultirateSplit()
268 …for (i = 0; i < s; i++) PetscCall(VecRestoreSubVector(rk->YdotRHS[i], rk->is_slow, &rk->YdotRHS_sl… in TSInterpolate_RK_MultirateSplit()
282 TS_RK *rk = (TS_RK *)ts->data; in TSEvaluateStep_RK_MultirateSplit() local
283 RKTableau tab = rk->tableau; in TSEvaluateStep_RK_MultirateSplit()
285 PetscScalar *w = rk->work; in TSEvaluateStep_RK_MultirateSplit()
291 if (rk->slow) { in TSEvaluateStep_RK_MultirateSplit()
293 PetscCall(VecGetSubVector(ts->vec_sol, rk->is_slow, &Xslow)); in TSEvaluateStep_RK_MultirateSplit()
294 PetscCall(VecMAXPY(Xslow, s, w, rk->YdotRHS_slow)); in TSEvaluateStep_RK_MultirateSplit()
295 PetscCall(VecRestoreSubVector(ts->vec_sol, rk->is_slow, &Xslow)); in TSEvaluateStep_RK_MultirateSplit()
297 for (j = 0; j < s; j++) w[j] = h / rk->dtratio * tab->b[j]; in TSEvaluateStep_RK_MultirateSplit()
298 PetscCall(VecGetSubVector(X, rk->is_fast, &Xfast)); in TSEvaluateStep_RK_MultirateSplit()
299 PetscCall(VecMAXPY(Xfast, s, w, rk->YdotRHS_fast)); in TSEvaluateStep_RK_MultirateSplit()
300 PetscCall(VecRestoreSubVector(X, rk->is_fast, &Xfast)); in TSEvaluateStep_RK_MultirateSplit()
307 TS_RK *rk = (TS_RK *)ts->data; in TSStepRefine_RK_MultirateSplit() local
308 TS subts_fast = rk->subts_fast, currentlevelts; in TSStepRefine_RK_MultirateSplit()
310 RKTableau tab = rk->tableau; in TSStepRefine_RK_MultirateSplit()
311 Vec *Y = rk->Y; in TSStepRefine_RK_MultirateSplit()
312 Vec *YdotRHS = rk->YdotRHS, *YdotRHS_fast = rk->YdotRHS_fast; in TSStepRefine_RK_MultirateSplit()
316 PetscScalar *w = rk->work; in TSStepRefine_RK_MultirateSplit()
321 for (k = 0; k < rk->dtratio; k++) { in TSStepRefine_RK_MultirateSplit()
322 PetscCall(VecGetSubVector(ts->vec_sol, rk->is_fast, &Xfast)); in TSStepRefine_RK_MultirateSplit()
323 for (i = 0; i < s; i++) PetscCall(VecGetSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i])); in TSStepRefine_RK_MultirateSplit()
327 …PetscCall(TSInterpolate_RK_MultirateSplit(rk->ts_root, t + k * h / rk->dtratio + h / rk->dtratio *… in TSStepRefine_RK_MultirateSplit()
328 currentlevelts = rk->ts_root; in TSStepRefine_RK_MultirateSplit()
331 …PetscCall(TSInterpolate_RK_MultirateSplit(currentlevelts, t + k * h / rk->dtratio + h / rk->dtrati… in TSStepRefine_RK_MultirateSplit()
333 for (j = 0; j < i; j++) w[j] = h / rk->dtratio * A[i * s + j]; in TSStepRefine_RK_MultirateSplit()
334 subrk_fast->stage_time = t + h / rk->dtratio * c[i]; in TSStepRefine_RK_MultirateSplit()
337 PetscCall(VecGetSubVector(Y[i], rk->is_fast, &Yfast)); in TSStepRefine_RK_MultirateSplit()
340 PetscCall(VecRestoreSubVector(Y[i], rk->is_fast, &Yfast)); in TSStepRefine_RK_MultirateSplit()
343 …PetscCall(TSComputeRHSFunction(subts_fast, t + k * h * rk->dtratio + h / rk->dtratio * c[i], Y[i],… in TSStepRefine_RK_MultirateSplit()
345 PetscCall(VecRestoreSubVector(ts->vec_sol, rk->is_fast, &Xfast)); in TSStepRefine_RK_MultirateSplit()
347 rk->slow = PETSC_FALSE; in TSStepRefine_RK_MultirateSplit()
349 … for (i = 0; i < s; i++) PetscCall(VecRestoreSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i])); in TSStepRefine_RK_MultirateSplit()
352 subts_fast->ptime = t + k * h / rk->dtratio; in TSStepRefine_RK_MultirateSplit()
353 subts_fast->time_step = h / rk->dtratio; in TSStepRefine_RK_MultirateSplit()
357 PetscCall(VecGetSubVector(ts->vec_sol, rk->is_fast, &Xfast)); in TSStepRefine_RK_MultirateSplit()
358 PetscCall(VecISCopy(rk->X0, rk->is_fast, SCATTER_FORWARD, Xfast)); in TSStepRefine_RK_MultirateSplit()
359 PetscCall(VecRestoreSubVector(ts->vec_sol, rk->is_fast, &Xfast)); in TSStepRefine_RK_MultirateSplit()
366 TS_RK *rk = (TS_RK *)ts->data; in TSStep_RK_MultirateSplit() local
367 RKTableau tab = rk->tableau; in TSStep_RK_MultirateSplit()
368 Vec *Y = rk->Y, *YdotRHS = rk->YdotRHS; in TSStep_RK_MultirateSplit()
369 Vec *YdotRHS_fast = rk->YdotRHS_fast, *YdotRHS_slow = rk->YdotRHS_slow; in TSStep_RK_MultirateSplit()
373 PetscScalar *w = rk->work; in TSStep_RK_MultirateSplit()
378 rk->status = TS_STEP_INCOMPLETE; in TSStep_RK_MultirateSplit()
380 PetscCall(VecGetSubVector(YdotRHS[i], rk->is_slow, &YdotRHS_slow[i])); in TSStep_RK_MultirateSplit()
381 PetscCall(VecGetSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i])); in TSStep_RK_MultirateSplit()
383 PetscCall(VecCopy(ts->vec_sol, rk->X0)); in TSStep_RK_MultirateSplit()
386 rk->stage_time = t + h * c[i]; in TSStep_RK_MultirateSplit()
387 PetscCall(TSPreStage(ts, rk->stage_time)); in TSStep_RK_MultirateSplit()
389 PetscCall(VecGetSubVector(Y[i], rk->is_fast, &Yfast)); in TSStep_RK_MultirateSplit()
390 PetscCall(VecGetSubVector(Y[i], rk->is_slow, &Yslow)); in TSStep_RK_MultirateSplit()
394 PetscCall(VecRestoreSubVector(Y[i], rk->is_fast, &Yfast)); in TSStep_RK_MultirateSplit()
395 PetscCall(VecRestoreSubVector(Y[i], rk->is_slow, &Yslow)); in TSStep_RK_MultirateSplit()
396 PetscCall(TSPostStage(ts, rk->stage_time, i, Y)); in TSStep_RK_MultirateSplit()
397 PetscCall(TSComputeRHSFunction(rk->subts_slow, t + h * c[i], Y[i], YdotRHS_slow[i])); in TSStep_RK_MultirateSplit()
398 PetscCall(TSComputeRHSFunction(rk->subts_fast, t + h * c[i], Y[i], YdotRHS_fast[i])); in TSStep_RK_MultirateSplit()
400 rk->slow = PETSC_TRUE; in TSStep_RK_MultirateSplit()
405 PetscCall(VecRestoreSubVector(YdotRHS[i], rk->is_slow, &YdotRHS_slow[i])); in TSStep_RK_MultirateSplit()
406 PetscCall(VecRestoreSubVector(YdotRHS[i], rk->is_fast, &YdotRHS_fast[i])); in TSStep_RK_MultirateSplit()
413 rk->status = TS_STEP_COMPLETE; in TSStep_RK_MultirateSplit()
419 TS_RK *rk = (TS_RK *)ts->data, *nextlevelrk, *currentlevelrk; in TSSetUp_RK_MultirateSplit() local
424 PetscCall(TSRHSSplitGetIS(ts, "slow", &rk->is_slow)); in TSSetUp_RK_MultirateSplit()
425 PetscCall(TSRHSSplitGetIS(ts, "fast", &rk->is_fast)); in TSSetUp_RK_MultirateSplit()
426 …PetscCheck(rk->is_slow && rk->is_fast, PetscObjectComm((PetscObject)ts), PETSC_ERR_USER, "Must set… in TSSetUp_RK_MultirateSplit()
427 PetscCall(TSRHSSplitGetSubTS(ts, "slow", &rk->subts_slow)); in TSSetUp_RK_MultirateSplit()
428 PetscCall(TSRHSSplitGetSubTS(ts, "fast", &rk->subts_fast)); in TSSetUp_RK_MultirateSplit()
429 …PetscCheck(rk->subts_slow && rk->subts_fast, PetscObjectComm((PetscObject)ts), PETSC_ERR_USER, "Mu… in TSSetUp_RK_MultirateSplit()
433 currentlevelrk = rk; in TSSetUp_RK_MultirateSplit()
435 PetscCall(PetscMalloc1(rk->tableau->s, ¤tlevelrk->YdotRHS_fast)); in TSSetUp_RK_MultirateSplit()
436 PetscCall(PetscMalloc1(rk->tableau->s, ¤tlevelrk->YdotRHS_slow)); in TSSetUp_RK_MultirateSplit()
444 nextlevelrk->tableau = rk->tableau; in TSSetUp_RK_MultirateSplit()
445 nextlevelrk->work = rk->work; in TSSetUp_RK_MultirateSplit()
446 nextlevelrk->Y = rk->Y; in TSSetUp_RK_MultirateSplit()
447 nextlevelrk->YdotRHS = rk->YdotRHS; in TSSetUp_RK_MultirateSplit()
455 nextlevelrk->tableau = rk->tableau; in TSSetUp_RK_MultirateSplit()
456 nextlevelrk->work = rk->work; in TSSetUp_RK_MultirateSplit()
457 nextlevelrk->Y = rk->Y; in TSSetUp_RK_MultirateSplit()
458 nextlevelrk->YdotRHS = rk->YdotRHS; in TSSetUp_RK_MultirateSplit()
459 nextlevelrk->dtratio = rk->dtratio; in TSSetUp_RK_MultirateSplit()
480 TS_RK *rk = (TS_RK *)ts->data; in TSRKSetMultirate_RK() local
484 rk->use_multirate = use_multirate; in TSRKSetMultirate_RK()
486 rk->dtratio = 2; in TSRKSetMultirate_RK()
492 rk->dtratio = 0; in TSRKSetMultirate_RK()
503 TS_RK *rk = (TS_RK *)ts->data; in TSRKGetMultirate_RK() local
507 *use_multirate = rk->use_multirate; in TSRKGetMultirate_RK()