1af0996ceSBarry Smith #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 28d59e960SJed Brown 35e4ed32fSEmil Constantinescu static PetscErrorCode TSAdaptChoose_CFL(TSAdapt adapt,TS ts,PetscReal h,PetscInt *next_sc,PetscReal *next_h,PetscBool *accept,PetscReal *wlte,PetscReal *wltea,PetscReal *wlter) 48d59e960SJed Brown { 58d59e960SJed Brown PetscErrorCode ierr; 61566a47fSLisandro Dalcin PetscReal hcfl,cfltimestep,ccfl; 71566a47fSLisandro Dalcin PetscInt ncandidates; 81566a47fSLisandro Dalcin const PetscReal *ccflarray; 98d59e960SJed Brown 108d59e960SJed Brown PetscFunctionBegin; 111566a47fSLisandro Dalcin ierr = TSGetCFLTime(ts,&cfltimestep);CHKERRQ(ierr); 121566a47fSLisandro Dalcin ierr = TSAdaptCandidatesGet(adapt,&ncandidates,NULL,NULL,&ccflarray,NULL);CHKERRQ(ierr); 131566a47fSLisandro Dalcin ccfl = (ncandidates > 0) ? ccflarray[0] : 1.0; 148d59e960SJed Brown 15*1917a363SLisandro Dalcin if (!adapt->always_accept) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_SUP,"Step rejection not implemented. The CFL implementation is incomplete/unusable"); 16*1917a363SLisandro Dalcin 171566a47fSLisandro Dalcin /* Determine whether the step is accepted of rejected */ 181566a47fSLisandro Dalcin *accept = PETSC_TRUE; 191566a47fSLisandro Dalcin if (h > cfltimestep * ccfl) { 20bf997491SLisandro Dalcin if (adapt->always_accept) { 211566a47fSLisandro Dalcin 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); 228d59e960SJed Brown } else { 231566a47fSLisandro Dalcin 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); 248d59e960SJed Brown *accept = PETSC_FALSE; 258d59e960SJed Brown } 268d59e960SJed Brown } 278d59e960SJed Brown 281566a47fSLisandro Dalcin /* The optimal new step based purely on CFL constraint for this step. */ 29*1917a363SLisandro Dalcin hcfl = adapt->safety * cfltimestep * ccfl; 301566a47fSLisandro Dalcin if (hcfl < adapt->dt_min) { 31*1917a363SLisandro Dalcin 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)adapt->safety,(double)adapt->dt_min,(double)ccfl);CHKERRQ(ierr); 321566a47fSLisandro Dalcin } 331566a47fSLisandro Dalcin 348d59e960SJed Brown *next_sc = 0; 358d59e960SJed Brown *next_h = PetscClipInterval(hcfl,adapt->dt_min,adapt->dt_max); 360b99f514SJed Brown *wlte = -1; /* Weighted local truncation error was not evaluated */ 375e4ed32fSEmil Constantinescu *wltea = -1; /* Weighted absolute local truncation error was not evaluated */ 385e4ed32fSEmil Constantinescu *wlter = -1; /* Weighted relative local truncation error was not evaluated */ 398d59e960SJed Brown PetscFunctionReturn(0); 408d59e960SJed Brown } 418d59e960SJed Brown 428d59e960SJed Brown /*MC 438d59e960SJed Brown TSADAPTCFL - CFL adaptive controller for time stepping 448d59e960SJed Brown 458d59e960SJed Brown Level: intermediate 468d59e960SJed Brown 478d59e960SJed Brown .seealso: TS, TSAdapt, TSSetAdapt() 488d59e960SJed Brown M*/ 498cc058d9SJed Brown PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt adapt) 508d59e960SJed Brown { 518d59e960SJed Brown PetscFunctionBegin; 528d59e960SJed Brown adapt->ops->choose = TSAdaptChoose_CFL; 538d59e960SJed Brown PetscFunctionReturn(0); 548d59e960SJed Brown } 55