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