Actual source code: adapthist.c
1: #include <petsc/private/tshistoryimpl.h>
3: typedef struct {
4: TSHistory hist;
5: PetscBool bw;
6: } TSAdapt_History;
8: static PetscErrorCode TSAdaptChoose_History(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept, PetscReal *wlte, PetscReal *wltea, PetscReal *wlter)
9: {
10: PetscInt step;
11: TSAdapt_History *thadapt = (TSAdapt_History *)adapt->data;
14: TSGetStepNumber(ts, &step);
15: TSHistoryGetTimeStep(thadapt->hist, thadapt->bw, step + 1, next_h);
16: *accept = PETSC_TRUE;
17: *next_sc = 0;
18: *wlte = -1;
19: *wltea = -1;
20: *wlter = -1;
21: return 0;
22: }
24: static PetscErrorCode TSAdaptReset_History(TSAdapt adapt)
25: {
26: TSAdapt_History *thadapt = (TSAdapt_History *)adapt->data;
28: TSHistoryDestroy(&thadapt->hist);
29: return 0;
30: }
32: static PetscErrorCode TSAdaptDestroy_History(TSAdapt adapt)
33: {
34: TSAdaptReset_History(adapt);
35: PetscFree(adapt->data);
36: return 0;
37: }
39: /* this is not public, as TSHistory is not a public object */
40: PetscErrorCode TSAdaptHistorySetTSHistory(TSAdapt adapt, TSHistory hist, PetscBool backward)
41: {
42: PetscReal *hist_t;
43: PetscInt n;
44: PetscBool flg;
48: PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg);
49: if (!flg) return 0;
50: TSHistoryGetHistory(hist, &n, (const PetscReal **)&hist_t, NULL, NULL);
51: TSAdaptHistorySetHistory(adapt, n, hist_t, backward);
52: return 0;
53: }
55: /*@
56: TSAdaptHistoryGetStep - Gets time and time step for a given step number in the history
58: Logically Collective
60: Input Parameters:
61: + adapt - the TSAdapt context
62: - step - the step number
64: Output Parameters:
65: + t - the time corresponding to the requested step (can be NULL)
66: - dt - the time step to be taken at the requested step (can be NULL)
68: Level: advanced
70: Note:
71: The time history is internally copied, and the user can free the hist array. The user still needs to specify the termination of the solve via `TSSetMaxSteps()`.
73: .seealso: [](chapter_ts), `TS`, `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptHistorySetTrajectory()`, `TSADAPTHISTORY`
74: @*/
75: PetscErrorCode TSAdaptHistoryGetStep(TSAdapt adapt, PetscInt step, PetscReal *t, PetscReal *dt)
76: {
77: TSAdapt_History *thadapt;
78: PetscBool flg;
82: PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg);
84: thadapt = (TSAdapt_History *)adapt->data;
85: TSHistoryGetTimeStep(thadapt->hist, thadapt->bw, step, dt);
86: TSHistoryGetTime(thadapt->hist, thadapt->bw, step, t);
87: return 0;
88: }
90: /*@
91: TSAdaptHistorySetHistory - Sets the time history in the adaptor
93: Logically Collective
95: Input Parameters:
96: + adapt - the `TSAdapt` context
97: . n - size of the time history
98: . hist - the time history
99: - backward - if the time history has to be followed backward
101: Level: advanced
103: Note:
104: The time history is internally copied, and the user can free the hist array. The user still needs to specify the termination of the solve via `TSSetMaxSteps()`.
106: .seealso: [](chapter_ts), `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptHistorySetTrajectory()`, `TSADAPTHISTORY`
107: @*/
108: PetscErrorCode TSAdaptHistorySetHistory(TSAdapt adapt, PetscInt n, PetscReal hist[], PetscBool backward)
109: {
110: TSAdapt_History *thadapt;
111: PetscBool flg;
117: PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg);
118: if (!flg) return 0;
119: thadapt = (TSAdapt_History *)adapt->data;
120: TSHistoryDestroy(&thadapt->hist);
121: TSHistoryCreate(PetscObjectComm((PetscObject)adapt), &thadapt->hist);
122: TSHistorySetHistory(thadapt->hist, n, hist, NULL, PETSC_FALSE);
123: thadapt->bw = backward;
124: return 0;
125: }
127: /*@
128: TSAdaptHistorySetTrajectory - Sets a time history in the adaptor from a given `TSTrajectory`
130: Logically Collective
132: Input Parameters:
133: + adapt - the `TSAdapt` context
134: . tj - the `TSTrajectory` context
135: - backward - if the time history has to be followed backward
137: Level: advanced
139: Notes:
140: The time history is internally copied, and the user can destroy the `TSTrajectory` if not needed.
142: The user still needs to specify the termination of the solve via `TSSetMaxSteps()`.
144: .seealso: [](chapter_ts), `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptHistorySetHistory()`, `TSADAPTHISTORY`, `TSAdapt`
145: @*/
146: PetscErrorCode TSAdaptHistorySetTrajectory(TSAdapt adapt, TSTrajectory tj, PetscBool backward)
147: {
148: PetscBool flg;
153: PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg);
154: if (!flg) return 0;
155: TSAdaptHistorySetTSHistory(adapt, tj->tsh, backward);
156: return 0;
157: }
159: /*MC
160: TSADAPTHISTORY - Time stepping controller that follows a given time history, used for Tangent Linear Model simulations
162: Level: developer
164: .seealso: [](chapter_ts), `TS`, `TSAdapt`, `TSGetAdapt()`, `TSAdaptHistorySetHistory()`, `TSAdaptType`
165: M*/
166: PETSC_EXTERN PetscErrorCode TSAdaptCreate_History(TSAdapt adapt)
167: {
168: TSAdapt_History *thadapt;
170: PetscNew(&thadapt);
171: adapt->matchstepfac[0] = PETSC_SMALL; /* prevent from accumulation errors */
172: adapt->matchstepfac[1] = 0.0; /* we will always match the final step, prevent TSAdaptChoose to mess with it */
173: adapt->data = thadapt;
175: adapt->ops->choose = TSAdaptChoose_History;
176: adapt->ops->reset = TSAdaptReset_History;
177: adapt->ops->destroy = TSAdaptDestroy_History;
178: return 0;
179: }