1 #include <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 PetscErrorCode ierr; 17 Vec X,Y; 18 PetscReal enorm,hfac_lte,h_lte,safety; 19 PetscInt order,stepno; 20 21 PetscFunctionBegin; 22 ierr = TSGetTimeStepNumber(ts,&stepno);CHKERRQ(ierr); 23 ierr = TSGetSolution(ts,&X);CHKERRQ(ierr); 24 if (!basic->Y) {ierr = VecDuplicate(X,&basic->Y);CHKERRQ(ierr);} 25 Y = basic->Y; 26 order = adapt->candidates.order[0]; 27 ierr = TSEvaluateStep(ts,order-1,Y,PETSC_NULL);CHKERRQ(ierr); 28 29 safety = basic->safety; 30 ierr = TSErrorNormWRMS(ts,Y,&enorm);CHKERRQ(ierr); 31 if (enorm > 1.) { 32 if (!*accept) safety *= basic->reject_safety; /* The last attempt also failed, shorten more aggressively */ 33 if (h < (1 + PETSC_SQRT_MACHINE_EPSILON)*adapt->dt_min) { 34 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %G, accepting because step size %G is at minimum\n",enorm,h);CHKERRQ(ierr); 35 *accept = PETSC_TRUE; 36 } else if (basic->always_accept) { 37 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %G, accepting step of size %G because always_accept is set\n",enorm,h);CHKERRQ(ierr); 38 *accept = PETSC_TRUE; 39 } else { 40 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %G, rejecting step of size %G\n",enorm,h);CHKERRQ(ierr); 41 *accept = PETSC_FALSE; 42 } 43 } else { 44 ierr = PetscInfo2(adapt,"Estimated scaled local truncation error %G, accepting step of size %G\n",enorm,h);CHKERRQ(ierr); 45 *accept = PETSC_TRUE; 46 } 47 48 /* The optimal new step based purely on local truncation error for this step. */ 49 hfac_lte = safety * PetscRealPart(PetscPowScalar((PetscScalar)enorm,(PetscReal)(-1./order))); 50 h_lte = h * PetscClipInterval(hfac_lte,basic->clip[0],basic->clip[1]); 51 52 *next_sc = 0; 53 *next_h = PetscClipInterval(h_lte,adapt->dt_min,adapt->dt_max); 54 *wlte = enorm; 55 PetscFunctionReturn(0); 56 } 57 58 #undef __FUNCT__ 59 #define __FUNCT__ "TSAdaptDestroy_Basic" 60 static PetscErrorCode TSAdaptDestroy_Basic(TSAdapt adapt) 61 { 62 TSAdapt_Basic *basic = (TSAdapt_Basic*)adapt->data; 63 PetscErrorCode ierr; 64 65 PetscFunctionBegin; 66 ierr = VecDestroy(&basic->Y);CHKERRQ(ierr); 67 ierr = PetscFree(adapt->data);CHKERRQ(ierr); 68 PetscFunctionReturn(0); 69 } 70 71 #undef __FUNCT__ 72 #define __FUNCT__ "TSAdaptSetFromOptions_Basic" 73 static PetscErrorCode TSAdaptSetFromOptions_Basic(TSAdapt adapt) 74 { 75 TSAdapt_Basic *basic = (TSAdapt_Basic*)adapt->data; 76 PetscErrorCode ierr; 77 PetscInt two; 78 PetscBool set; 79 80 PetscFunctionBegin; 81 ierr = PetscOptionsHead("Basic adaptive controller options");CHKERRQ(ierr); 82 two = 2; 83 ierr = PetscOptionsRealArray("-ts_adapt_basic_clip","Admissible decrease/increase in step size","",basic->clip,&two,&set);CHKERRQ(ierr); 84 if (set && (two != 2 || basic->clip[0] > basic->clip[1])) 85 SETERRQ(((PetscObject)adapt)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Must give exactly two values to -ts_adapt_basic_clip"); 86 ierr = PetscOptionsReal("-ts_adapt_basic_safety","Safety factor relative to target error","",basic->safety,&basic->safety,PETSC_NULL);CHKERRQ(ierr); 87 ierr = PetscOptionsReal("-ts_adapt_basic_reject_safety","Extra safety factor to apply if the last step was rejected","",basic->reject_safety,&basic->reject_safety,PETSC_NULL);CHKERRQ(ierr); 88 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,PETSC_NULL);CHKERRQ(ierr); 89 ierr = PetscOptionsTail();CHKERRQ(ierr); 90 PetscFunctionReturn(0); 91 } 92 93 EXTERN_C_BEGIN 94 #undef __FUNCT__ 95 #define __FUNCT__ "TSAdaptCreate_Basic" 96 /*MC 97 TSADAPTBASIC - Basic adaptive controller for time stepping 98 99 Level: intermediate 100 101 .seealso: TS, TSAdapt, TSSetAdapt() 102 M*/ 103 PetscErrorCode TSAdaptCreate_Basic(TSAdapt adapt) 104 { 105 PetscErrorCode ierr; 106 TSAdapt_Basic *a; 107 108 PetscFunctionBegin; 109 ierr = PetscNewLog(adapt,TSAdapt_Basic,&a);CHKERRQ(ierr); 110 adapt->data = (void*)a; 111 adapt->ops->choose = TSAdaptChoose_Basic; 112 adapt->ops->setfromoptions = TSAdaptSetFromOptions_Basic; 113 adapt->ops->destroy = TSAdaptDestroy_Basic; 114 115 a->clip[0] = 0.1; 116 a->clip[1] = 10.; 117 a->safety = 0.9; 118 a->reject_safety = 0.5; 119 a->always_accept = PETSC_FALSE; 120 PetscFunctionReturn(0); 121 } 122 EXTERN_C_END 123