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