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 TSAdaptHistorySetHistory - Sets a time history in the adaptor 66 67 Logically Collective on TSAdapt 68 69 Input Parameters: 70 + adapt - the TSAdapt context 71 . n - size of the time history 72 . hist - the time history 73 - backward - if the time history has to be followed backward 74 75 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(). 76 77 Level: advanced 78 79 .keywords: TSAdapt 80 .seealso: TSGetAdapt(), TSAdaptSetType(), TSADAPTHISTORY 81 @*/ 82 PetscErrorCode TSAdaptHistorySetHistory(TSAdapt adapt, PetscInt n, PetscReal hist[], PetscBool backward) 83 { 84 TSAdapt_History *thadapt; 85 PetscBool flg; 86 PetscErrorCode ierr; 87 88 PetscFunctionBegin; 89 PetscValidHeaderSpecific(adapt,TSADAPT_CLASSID,1); 90 PetscValidLogicalCollectiveInt(adapt,n,2); 91 PetscValidRealPointer(hist,3); 92 PetscValidLogicalCollectiveBool(adapt,backward,4); 93 ierr = PetscObjectTypeCompare((PetscObject)adapt,TSADAPTHISTORY,&flg);CHKERRQ(ierr); 94 if (!flg) PetscFunctionReturn(0); 95 thadapt = (TSAdapt_History*)adapt->data; 96 ierr = TSHistoryDestroy(&thadapt->hist);CHKERRQ(ierr); 97 ierr = TSHistoryCreate(PetscObjectComm((PetscObject)adapt),&thadapt->hist);CHKERRQ(ierr); 98 ierr = TSHistorySetHistory(thadapt->hist,n,hist,NULL,PETSC_FALSE);CHKERRQ(ierr); 99 thadapt->bw = backward; 100 PetscFunctionReturn(0); 101 } 102 103 /*MC 104 TSADAPTHISTORY - Time stepping controller that follows a given time history, used for Tangent Linear Model simulations 105 106 Level: developer 107 108 .seealso: TS, TSAdapt, TSGetAdapt(), TSAdaptHistorySetHistory() 109 M*/ 110 PETSC_EXTERN PetscErrorCode TSAdaptCreate_History(TSAdapt adapt) 111 { 112 PetscErrorCode ierr; 113 TSAdapt_History *thadapt; 114 115 PetscFunctionBegin; 116 ierr = PetscNew(&thadapt);CHKERRQ(ierr); 117 adapt->matchstepfac[0] = PETSC_SMALL; /* prevent from accumulation errors */ 118 adapt->matchstepfac[1] = 0.0; /* we will always match the final step, prevent TSAdaptChoose to mess with it */ 119 adapt->data = thadapt; 120 121 adapt->ops->choose = TSAdaptChoose_History; 122 adapt->ops->reset = TSAdaptReset_History; 123 adapt->ops->destroy = TSAdaptDestroy_History; 124 PetscFunctionReturn(0); 125 } 126