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