1 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 2 3 typedef struct { 4 PetscBool always_accept; 5 PetscReal clip[2]; /* admissible decrease/increase factors */ 6 PetscReal safety; /* safety factor relative to target error */ 7 PetscReal reject_safety; /* extra safety factor if the last step was rejected */ 8 Vec Y; 9 } TSAdapt_Basic; 10 11 #undef __FUNCT__ 12 #define __FUNCT__ "TSAdaptChoose_Basic" 13 static PetscErrorCode TSAdaptChoose_Basic(TSAdapt adapt,TS ts,PetscReal h,PetscInt *next_sc,PetscReal *next_h,PetscBool *accept,PetscReal *wlte) 14 { 15 TSAdapt_Basic *basic = (TSAdapt_Basic*)adapt->data; 16 PetscInt order; 17 PetscReal safety = basic->safety; 18 PetscReal enorm=-1,hfac_lte,h_lte; 19 PetscErrorCode ierr; 20 21 PetscFunctionBegin; 22 *next_sc = 0; /* Reuse the same order scheme */ 23 24 if (ts->ops->evaluatewlte) { 25 order = PETSC_DECIDE; 26 ierr = TSEvaluateWLTE(ts,adapt->wnormtype,&order,&enorm);CHKERRQ(ierr); 27 if (enorm >= 0 && order < 1) SETERRQ1(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Computed error order %D must be positive",order); 28 } else if (ts->ops->evaluatestep) { 29 if (adapt->candidates.n < 1) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_WRONGSTATE,"No candidate has been registered"); 30 if (!adapt->candidates.inuse_set) SETERRQ1(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_WRONGSTATE,"The current in-use scheme is not among the %D candidates",adapt->candidates.n); 31 if (!basic->Y) {ierr = VecDuplicate(ts->vec_sol,&basic->Y);CHKERRQ(ierr);} 32 order = adapt->candidates.order[0]; 33 ierr = TSEvaluateStep(ts,order-1,basic->Y,NULL);CHKERRQ(ierr); 34 ierr = TSErrorWeightedNorm(ts,ts->vec_sol,basic->Y,adapt->wnormtype,&enorm);CHKERRQ(ierr); 35 } 36 37 if (enorm < 0) { 38 *accept = PETSC_TRUE; 39 *next_h = h; /* Reuse the old step */ 40 *wlte = -1; /* Weighted local truncation error was not evaluated */ 41 PetscFunctionReturn(0); 42 } 43 44 /* Determine whether the step is accepted of rejected */ 45 if (enorm > 1) { 46 if (!*accept) safety *= basic->reject_safety; /* The last attempt also failed, shorten more aggressively */ 47 if (h < (1 + PETSC_SQRT_MACHINE_EPSILON)*adapt->dt_min) { 48 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %g, accepting because step size %g is at minimum\n",(double)enorm,(double)h);CHKERRQ(ierr); 49 *accept = PETSC_TRUE; 50 } else if (basic->always_accept) { 51 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %g, accepting step of size %g because always_accept is set\n",(double)enorm,(double)h);CHKERRQ(ierr); 52 *accept = PETSC_TRUE; 53 } else { 54 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %g, rejecting step of size %g\n",(double)enorm,(double)h);CHKERRQ(ierr); 55 *accept = PETSC_FALSE; 56 } 57 } else { 58 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %g, accepting step of size %g\n",(double)enorm,(double)h);CHKERRQ(ierr); 59 *accept = PETSC_TRUE; 60 } 61 62 /* The optimal new step based purely on local truncation error for this step. */ 63 if (enorm > 0) 64 hfac_lte = safety * PetscPowReal(enorm,((PetscReal)-1)/order); 65 else 66 hfac_lte = safety * PETSC_INFINITY; 67 h_lte = h * PetscClipInterval(hfac_lte,basic->clip[0],basic->clip[1]); 68 69 *next_h = PetscClipInterval(h_lte,adapt->dt_min,adapt->dt_max); 70 *wlte = enorm; 71 PetscFunctionReturn(0); 72 } 73 74 #undef __FUNCT__ 75 #define __FUNCT__ "TSAdaptReset_Basic" 76 static PetscErrorCode TSAdaptReset_Basic(TSAdapt adapt) 77 { 78 TSAdapt_Basic *basic = (TSAdapt_Basic*)adapt->data; 79 PetscErrorCode ierr; 80 81 PetscFunctionBegin; 82 ierr = VecDestroy(&basic->Y);CHKERRQ(ierr); 83 PetscFunctionReturn(0); 84 } 85 86 #undef __FUNCT__ 87 #define __FUNCT__ "TSAdaptDestroy_Basic" 88 static PetscErrorCode TSAdaptDestroy_Basic(TSAdapt adapt) 89 { 90 PetscErrorCode ierr; 91 92 PetscFunctionBegin; 93 ierr = TSAdaptReset_Basic(adapt);CHKERRQ(ierr); 94 ierr = PetscFree(adapt->data);CHKERRQ(ierr); 95 PetscFunctionReturn(0); 96 } 97 98 #undef __FUNCT__ 99 #define __FUNCT__ "TSAdaptSetFromOptions_Basic" 100 static PetscErrorCode TSAdaptSetFromOptions_Basic(PetscOptionItems *PetscOptionsObject,TSAdapt adapt) 101 { 102 TSAdapt_Basic *basic = (TSAdapt_Basic*)adapt->data; 103 PetscErrorCode ierr; 104 PetscInt two; 105 PetscBool set; 106 107 PetscFunctionBegin; 108 ierr = PetscOptionsHead(PetscOptionsObject,"Basic adaptive controller options");CHKERRQ(ierr); 109 two = 2; 110 ierr = PetscOptionsRealArray("-ts_adapt_basic_clip","Admissible decrease/increase in step size","",basic->clip,&two,&set);CHKERRQ(ierr); 111 if (set && (two != 2 || basic->clip[0] > basic->clip[1])) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_ARG_OUTOFRANGE,"Must give exactly two values to -ts_adapt_basic_clip"); 112 ierr = PetscOptionsReal("-ts_adapt_basic_safety","Safety factor relative to target error","",basic->safety,&basic->safety,NULL);CHKERRQ(ierr); 113 ierr = PetscOptionsReal("-ts_adapt_basic_reject_safety","Extra safety factor to apply if the last step was rejected","",basic->reject_safety,&basic->reject_safety,NULL);CHKERRQ(ierr); 114 ierr = PetscOptionsBool("-ts_adapt_basic_always_accept","Always accept the step regardless of whether local truncation error meets goal","",basic->always_accept,&basic->always_accept,NULL);CHKERRQ(ierr); 115 ierr = PetscOptionsTail();CHKERRQ(ierr); 116 PetscFunctionReturn(0); 117 } 118 119 #undef __FUNCT__ 120 #define __FUNCT__ "TSAdaptView_Basic" 121 static PetscErrorCode TSAdaptView_Basic(TSAdapt adapt,PetscViewer viewer) 122 { 123 TSAdapt_Basic *basic = (TSAdapt_Basic*)adapt->data; 124 PetscErrorCode ierr; 125 PetscBool iascii; 126 127 PetscFunctionBegin; 128 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 129 if (iascii) { 130 if (basic->always_accept) {ierr = PetscViewerASCIIPrintf(viewer," Basic: always accepting steps\n");CHKERRQ(ierr);} 131 ierr = PetscViewerASCIIPrintf(viewer," Basic: clip fastest decrease %g, fastest increase %g\n",(double)basic->clip[0],(double)basic->clip[1]);CHKERRQ(ierr); 132 ierr = PetscViewerASCIIPrintf(viewer," Basic: safety factor %g, extra factor after step rejection %g\n",(double)basic->safety,(double)basic->reject_safety);CHKERRQ(ierr); 133 } 134 PetscFunctionReturn(0); 135 } 136 137 #undef __FUNCT__ 138 #define __FUNCT__ "TSAdaptCreate_Basic" 139 /*MC 140 TSADAPTBASIC - Basic adaptive controller for time stepping 141 142 Level: intermediate 143 144 .seealso: TS, TSAdapt, TSSetAdapt() 145 M*/ 146 PETSC_EXTERN PetscErrorCode TSAdaptCreate_Basic(TSAdapt adapt) 147 { 148 PetscErrorCode ierr; 149 TSAdapt_Basic *a; 150 151 PetscFunctionBegin; 152 ierr = PetscNewLog(adapt,&a);CHKERRQ(ierr); 153 adapt->data = (void*)a; 154 adapt->ops->choose = TSAdaptChoose_Basic; 155 adapt->ops->setfromoptions = TSAdaptSetFromOptions_Basic; 156 adapt->ops->destroy = TSAdaptDestroy_Basic; 157 adapt->ops->view = TSAdaptView_Basic; 158 159 a->clip[0] = 0.1; 160 a->clip[1] = 10.; 161 a->safety = 0.9; 162 a->reject_safety = 0.5; 163 a->always_accept = PETSC_FALSE; 164 PetscFunctionReturn(0); 165 } 166