Actual source code: adaptcfl.c

  1: #include <petsc/private/tsimpl.h>

  3: static PetscErrorCode TSAdaptChoose_CFL(TSAdapt adapt, TS ts, PetscReal h, PetscInt *next_sc, PetscReal *next_h, PetscBool *accept, PetscReal *wlte, PetscReal *wltea, PetscReal *wlter)
  4: {
  5:   PetscReal        hcfl, cfltimestep, ccfl;
  6:   PetscInt         ncandidates;
  7:   const PetscReal *ccflarray;

  9:   TSGetCFLTime(ts, &cfltimestep);
 10:   TSAdaptCandidatesGet(adapt, &ncandidates, NULL, NULL, &ccflarray, NULL);
 11:   ccfl = (ncandidates > 0) ? ccflarray[0] : 1.0;


 15:   /* Determine whether the step is accepted of rejected */
 16:   *accept = PETSC_TRUE;
 17:   if (h > cfltimestep * ccfl) {
 18:     if (adapt->always_accept) {
 19:       PetscInfo(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);
 20:     } else {
 21:       PetscInfo(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);
 22:       *accept = PETSC_FALSE;
 23:     }
 24:   }

 26:   /* The optimal new step based purely on CFL constraint for this step. */
 27:   hcfl = adapt->safety * cfltimestep * ccfl;
 28:   if (hcfl < adapt->dt_min) {
 29:     PetscInfo(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);
 30:   }

 32:   *next_sc = 0;
 33:   *next_h  = PetscClipInterval(hcfl, adapt->dt_min, adapt->dt_max);
 34:   *wlte    = -1; /* Weighted local truncation error was not evaluated */
 35:   *wltea   = -1; /* Weighted absolute local truncation error was not evaluated */
 36:   *wlter   = -1; /* Weighted relative local truncation error was not evaluated */
 37:   return 0;
 38: }

 40: /*MC
 41:    TSADAPTCFL - CFL adaptive controller for time stepping

 43:    Level: intermediate

 45: .seealso: [](chapter_ts), `TS`, `TSAdapt`, `TSGetAdapt()`, `TSAdaptType`
 46: M*/
 47: PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt adapt)
 48: {
 49:   adapt->ops->choose = TSAdaptChoose_CFL;
 50:   return 0;
 51: }