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