1 #include <petsc/private/tshistoryimpl.h> /*I "petscts.h" I*/ 2 3 typedef struct { 4 TSHistory hist; 5 PetscBool bw; 6 } TSAdapt_History; 7 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 PetscInt step; 10 TSAdapt_History *thadapt = (TSAdapt_History *)adapt->data; 11 12 PetscFunctionBegin; 13 PetscCheck(thadapt->hist, PetscObjectComm((PetscObject)adapt), PETSC_ERR_ORDER, "Need to call TSAdaptHistorySetHistory() first"); 14 PetscCall(TSGetStepNumber(ts, &step)); 15 PetscCall(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 PetscFunctionReturn(0); 22 } 23 24 static PetscErrorCode TSAdaptReset_History(TSAdapt adapt) { 25 TSAdapt_History *thadapt = (TSAdapt_History *)adapt->data; 26 27 PetscFunctionBegin; 28 PetscCall(TSHistoryDestroy(&thadapt->hist)); 29 PetscFunctionReturn(0); 30 } 31 32 static PetscErrorCode TSAdaptDestroy_History(TSAdapt adapt) { 33 PetscFunctionBegin; 34 PetscCall(TSAdaptReset_History(adapt)); 35 PetscCall(PetscFree(adapt->data)); 36 PetscFunctionReturn(0); 37 } 38 39 /* this is not public, as TSHistory is not a public object */ 40 PetscErrorCode TSAdaptHistorySetTSHistory(TSAdapt adapt, TSHistory hist, PetscBool backward) { 41 PetscReal *hist_t; 42 PetscInt n; 43 PetscBool flg; 44 45 PetscFunctionBegin; 46 PetscValidHeaderSpecific(adapt, TSADAPT_CLASSID, 1); 47 PetscValidLogicalCollectiveBool(adapt, backward, 3); 48 PetscCall(PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg)); 49 if (!flg) PetscFunctionReturn(0); 50 PetscCall(TSHistoryGetHistory(hist, &n, (const PetscReal **)&hist_t, NULL, NULL)); 51 PetscCall(TSAdaptHistorySetHistory(adapt, n, hist_t, backward)); 52 PetscFunctionReturn(0); 53 } 54 55 /*@ 56 TSAdaptHistoryGetStep - Gets time and time step for a given step number in the history 57 58 Logically Collective on TSAdapt 59 60 Input Parameters: 61 + adapt - the TSAdapt context 62 - step - the step number 63 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) 67 68 Notes: 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(). 69 70 Level: advanced 71 72 .seealso: `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptHistorySetTrajectory()`, `TSADAPTHISTORY` 73 @*/ 74 PetscErrorCode TSAdaptHistoryGetStep(TSAdapt adapt, PetscInt step, PetscReal *t, PetscReal *dt) { 75 TSAdapt_History *thadapt; 76 PetscBool flg; 77 78 PetscFunctionBegin; 79 PetscValidHeaderSpecific(adapt, TSADAPT_CLASSID, 1); 80 PetscValidLogicalCollectiveInt(adapt, step, 2); 81 PetscCall(PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg)); 82 PetscCheck(flg, PetscObjectComm((PetscObject)adapt), PETSC_ERR_SUP, "Not for type %s", ((PetscObject)adapt)->type_name); 83 thadapt = (TSAdapt_History *)adapt->data; 84 PetscCall(TSHistoryGetTimeStep(thadapt->hist, thadapt->bw, step, dt)); 85 PetscCall(TSHistoryGetTime(thadapt->hist, thadapt->bw, step, t)); 86 PetscFunctionReturn(0); 87 } 88 89 /*@ 90 TSAdaptHistorySetHistory - Sets a time history in the adaptor 91 92 Logically Collective on TSAdapt 93 94 Input Parameters: 95 + adapt - the TSAdapt context 96 . n - size of the time history 97 . hist - the time history 98 - backward - if the time history has to be followed backward 99 100 Notes: 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(). 101 102 Level: advanced 103 104 .seealso: `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptHistorySetTrajectory()`, `TSADAPTHISTORY` 105 @*/ 106 PetscErrorCode TSAdaptHistorySetHistory(TSAdapt adapt, PetscInt n, PetscReal hist[], PetscBool backward) { 107 TSAdapt_History *thadapt; 108 PetscBool flg; 109 110 PetscFunctionBegin; 111 PetscValidHeaderSpecific(adapt, TSADAPT_CLASSID, 1); 112 PetscValidLogicalCollectiveInt(adapt, n, 2); 113 PetscValidRealPointer(hist, 3); 114 PetscValidLogicalCollectiveBool(adapt, backward, 4); 115 PetscCall(PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg)); 116 if (!flg) PetscFunctionReturn(0); 117 thadapt = (TSAdapt_History *)adapt->data; 118 PetscCall(TSHistoryDestroy(&thadapt->hist)); 119 PetscCall(TSHistoryCreate(PetscObjectComm((PetscObject)adapt), &thadapt->hist)); 120 PetscCall(TSHistorySetHistory(thadapt->hist, n, hist, NULL, PETSC_FALSE)); 121 thadapt->bw = backward; 122 PetscFunctionReturn(0); 123 } 124 125 /*@ 126 TSAdaptHistorySetTrajectory - Sets a time history in the adaptor from a given TSTrajectory 127 128 Logically Collective on TSAdapt 129 130 Input Parameters: 131 + adapt - the TSAdapt context 132 . tj - the TSTrajectory context 133 - backward - if the time history has to be followed backward 134 135 Notes: The time history is internally copied, and the user can destroy the TSTrajectory if not needed. The user still needs to specify the termination of the solve via TSSetMaxSteps(). 136 137 Level: advanced 138 139 .seealso: `TSGetAdapt()`, `TSAdaptSetType()`, `TSAdaptHistorySetHistory()`, `TSADAPTHISTORY` 140 @*/ 141 PetscErrorCode TSAdaptHistorySetTrajectory(TSAdapt adapt, TSTrajectory tj, PetscBool backward) { 142 PetscBool flg; 143 144 PetscFunctionBegin; 145 PetscValidHeaderSpecific(adapt, TSADAPT_CLASSID, 1); 146 PetscValidHeaderSpecific(tj, TSTRAJECTORY_CLASSID, 2); 147 PetscValidLogicalCollectiveBool(adapt, backward, 3); 148 PetscCall(PetscObjectTypeCompare((PetscObject)adapt, TSADAPTHISTORY, &flg)); 149 if (!flg) PetscFunctionReturn(0); 150 PetscCall(TSAdaptHistorySetTSHistory(adapt, tj->tsh, backward)); 151 PetscFunctionReturn(0); 152 } 153 154 /*MC 155 TSADAPTHISTORY - Time stepping controller that follows a given time history, used for Tangent Linear Model simulations 156 157 Level: developer 158 159 .seealso: `TS`, `TSAdapt`, `TSGetAdapt()`, `TSAdaptHistorySetHistory()` 160 M*/ 161 PETSC_EXTERN PetscErrorCode TSAdaptCreate_History(TSAdapt adapt) { 162 TSAdapt_History *thadapt; 163 164 PetscFunctionBegin; 165 PetscCall(PetscNew(&thadapt)); 166 adapt->matchstepfac[0] = PETSC_SMALL; /* prevent from accumulation errors */ 167 adapt->matchstepfac[1] = 0.0; /* we will always match the final step, prevent TSAdaptChoose to mess with it */ 168 adapt->data = thadapt; 169 170 adapt->ops->choose = TSAdaptChoose_History; 171 adapt->ops->reset = TSAdaptReset_History; 172 adapt->ops->destroy = TSAdaptDestroy_History; 173 PetscFunctionReturn(0); 174 } 175