1 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 2 3 typedef struct { 4 PetscReal safety; /* safety factor relative to target CFL constraint */ 5 } TSAdapt_CFL; 6 7 static PetscErrorCode TSAdaptChoose_CFL(TSAdapt adapt,TS ts,PetscReal h,PetscInt *next_sc,PetscReal *next_h,PetscBool *accept,PetscReal *wlte,PetscReal *wltea,PetscReal *wlter) 8 { 9 TSAdapt_CFL *cfl = (TSAdapt_CFL*)adapt->data; 10 PetscErrorCode ierr; 11 PetscReal hcfl,cfltimestep,ccfl; 12 PetscInt ncandidates; 13 const PetscReal *ccflarray; 14 15 PetscFunctionBegin; 16 ierr = TSGetCFLTime(ts,&cfltimestep);CHKERRQ(ierr); 17 ierr = TSAdaptCandidatesGet(adapt,&ncandidates,NULL,NULL,&ccflarray,NULL);CHKERRQ(ierr); 18 ccfl = (ncandidates > 0) ? ccflarray[0] : 1.0; 19 20 /* Determine whether the step is accepted of rejected */ 21 *accept = PETSC_TRUE; 22 if (h > cfltimestep * ccfl) { 23 if (adapt->always_accept) { 24 ierr = PetscInfo3(adapt,"Step length %g with scheme of CFL coefficient %g did not satisfy user-provided CFL constraint %g, proceeding anyway\n",(double)h,(double)ccfl,(double)cfltimestep);CHKERRQ(ierr); 25 } else { 26 ierr = PetscInfo3(adapt,"Step length %g with scheme of CFL coefficient %g did not satisfy user-provided CFL constraint %g, step REJECTED\n",(double)h,(double)ccfl,(double)cfltimestep);CHKERRQ(ierr); 27 *accept = PETSC_FALSE; 28 } 29 } 30 31 /* The optimal new step based purely on CFL constraint for this step. */ 32 hcfl = cfl->safety * cfltimestep * ccfl; 33 if (hcfl < adapt->dt_min) { 34 ierr = PetscInfo4(adapt,"Cannot satisfy CFL constraint %g (with %g safety) at minimum time step %g with method coefficient %g, proceding anyway\n",(double)cfltimestep,(double)cfl->safety,(double)adapt->dt_min,(double)ccfl);CHKERRQ(ierr); 35 } 36 37 *next_sc = 0; 38 *next_h = PetscClipInterval(hcfl,adapt->dt_min,adapt->dt_max); 39 *wlte = -1; /* Weighted local truncation error was not evaluated */ 40 *wltea = -1; /* Weighted absolute local truncation error was not evaluated */ 41 *wlter = -1; /* Weighted relative local truncation error was not evaluated */ 42 PetscFunctionReturn(0); 43 } 44 45 static PetscErrorCode TSAdaptDestroy_CFL(TSAdapt adapt) 46 { 47 PetscErrorCode ierr; 48 49 PetscFunctionBegin; 50 ierr = PetscFree(adapt->data);CHKERRQ(ierr); 51 PetscFunctionReturn(0); 52 } 53 54 static PetscErrorCode TSAdaptSetFromOptions_CFL(PetscOptionItems *PetscOptionsObject,TSAdapt adapt) 55 { 56 TSAdapt_CFL *cfl = (TSAdapt_CFL*)adapt->data; 57 PetscErrorCode ierr; 58 59 PetscFunctionBegin; 60 ierr = PetscOptionsHead(PetscOptionsObject,"CFL adaptive controller options");CHKERRQ(ierr); 61 ierr = PetscOptionsReal("-ts_adapt_cfl_safety","Safety factor relative to target CFL constraint","",cfl->safety,&cfl->safety,NULL);CHKERRQ(ierr); 62 /* The TS implementations do not currently communicate CFL information to the controller. There is a placeholder, but 63 * we do not believe it to provide sufficiently rich information. That means the CFL adaptor is incomplete and 64 * unusable. Do not delete the guard below unless you have finished the implementation. */ 65 if (!adapt->always_accept) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_SUP,"Step rejection not implemented. The CFL implementation is incomplete/unusable"); 66 ierr = PetscOptionsTail();CHKERRQ(ierr); 67 PetscFunctionReturn(0); 68 } 69 70 /*MC 71 TSADAPTCFL - CFL adaptive controller for time stepping 72 73 Level: intermediate 74 75 .seealso: TS, TSAdapt, TSSetAdapt() 76 M*/ 77 PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt adapt) 78 { 79 PetscErrorCode ierr; 80 TSAdapt_CFL *a; 81 82 PetscFunctionBegin; 83 ierr = PetscNewLog(adapt,&a);CHKERRQ(ierr); 84 adapt->data = (void*)a; 85 86 adapt->ops->choose = TSAdaptChoose_CFL; 87 adapt->ops->setfromoptions = TSAdaptSetFromOptions_CFL; 88 adapt->ops->destroy = TSAdaptDestroy_CFL; 89 90 a->safety = 0.9; 91 PetscFunctionReturn(0); 92 } 93