1b45d2f2cSJed Brown #include <petsc-private/tsimpl.h> /*I "petscts.h" I*/ 28d59e960SJed Brown 38d59e960SJed Brown typedef struct { 48d59e960SJed Brown PetscBool always_accept; 58d59e960SJed Brown PetscReal safety; /* safety factor relative to target error */ 68d59e960SJed Brown } TSAdapt_CFL; 78d59e960SJed Brown 88d59e960SJed Brown #undef __FUNCT__ 98d59e960SJed Brown #define __FUNCT__ "TSAdaptChoose_CFL" 100b99f514SJed Brown static PetscErrorCode TSAdaptChoose_CFL(TSAdapt adapt,TS ts,PetscReal h,PetscInt *next_sc,PetscReal *next_h,PetscBool *accept,PetscReal *wlte) 118d59e960SJed Brown { 128d59e960SJed Brown TSAdapt_CFL *cfl = (TSAdapt_CFL*)adapt->data; 138d59e960SJed Brown PetscErrorCode ierr; 148d59e960SJed Brown PetscReal hcfl,cfltime; 158d59e960SJed Brown PetscInt stepno,ncandidates; 168d59e960SJed Brown const PetscInt *order; 178d59e960SJed Brown const PetscReal *ccfl; 188d59e960SJed Brown 198d59e960SJed Brown PetscFunctionBegin; 208d59e960SJed Brown ierr = TSGetTimeStepNumber(ts,&stepno);CHKERRQ(ierr); 218d59e960SJed Brown ierr = TSGetCFLTime(ts,&cfltime);CHKERRQ(ierr); 220298fd71SBarry Smith ierr = TSAdaptCandidatesGet(adapt,&ncandidates,&order,NULL,&ccfl,NULL);CHKERRQ(ierr); 238d59e960SJed Brown 248d59e960SJed Brown hcfl = cfl->safety * cfltime * ccfl[0]; 258d59e960SJed Brown if (hcfl < adapt->dt_min) { 268d59e960SJed Brown ierr = PetscInfo4(adapt,"Cannot satisfy CFL constraint %G (with %G safety) at minimum time step %G with method coefficient %G, proceding anyway\n",cfltime,cfl->safety,adapt->dt_min,ccfl[0]);CHKERRQ(ierr); 278d59e960SJed Brown } 288d59e960SJed Brown 298d59e960SJed Brown if (h > cfltime * ccfl[0]) { 308d59e960SJed Brown if (cfl->always_accept) { 318d59e960SJed Brown ierr = PetscInfo3(adapt,"Step length %G with scheme of CFL coefficient %G did not satisfy user-provided CFL constraint %G, proceeding anyway\n",h,ccfl[0],cfltime);CHKERRQ(ierr); 328d59e960SJed Brown } else { 338d59e960SJed Brown ierr = PetscInfo3(adapt,"Step length %G with scheme of CFL coefficient %G did not satisfy user-provided CFL constraint %G, step REJECTED\n",h,ccfl[0],cfltime);CHKERRQ(ierr); 348d59e960SJed Brown *next_sc = 0; 358d59e960SJed Brown *next_h = PetscClipInterval(hcfl,adapt->dt_min,adapt->dt_max); 368d59e960SJed Brown *accept = PETSC_FALSE; 378d59e960SJed Brown } 388d59e960SJed Brown } 398d59e960SJed Brown 408d59e960SJed Brown *next_sc = 0; 418d59e960SJed Brown *next_h = PetscClipInterval(hcfl,adapt->dt_min,adapt->dt_max); 428d59e960SJed Brown *accept = PETSC_TRUE; 430b99f514SJed Brown *wlte = -1; /* Weighted local truncation error was not evaluated */ 448d59e960SJed Brown PetscFunctionReturn(0); 458d59e960SJed Brown } 468d59e960SJed Brown 478d59e960SJed Brown #undef __FUNCT__ 488d59e960SJed Brown #define __FUNCT__ "TSAdaptDestroy_CFL" 498d59e960SJed Brown static PetscErrorCode TSAdaptDestroy_CFL(TSAdapt adapt) 508d59e960SJed Brown { 518d59e960SJed Brown PetscErrorCode ierr; 528d59e960SJed Brown 538d59e960SJed Brown PetscFunctionBegin; 548d59e960SJed Brown ierr = PetscFree(adapt->data);CHKERRQ(ierr); 558d59e960SJed Brown PetscFunctionReturn(0); 568d59e960SJed Brown } 578d59e960SJed Brown 588d59e960SJed Brown #undef __FUNCT__ 598d59e960SJed Brown #define __FUNCT__ "TSAdaptSetFromOptions_CFL" 608d59e960SJed Brown static PetscErrorCode TSAdaptSetFromOptions_CFL(TSAdapt adapt) 618d59e960SJed Brown { 628d59e960SJed Brown TSAdapt_CFL *cfl = (TSAdapt_CFL*)adapt->data; 638d59e960SJed Brown PetscErrorCode ierr; 648d59e960SJed Brown 658d59e960SJed Brown PetscFunctionBegin; 668d59e960SJed Brown ierr = PetscOptionsHead("CFL adaptive controller options");CHKERRQ(ierr); 670298fd71SBarry Smith ierr = PetscOptionsReal("-ts_adapt_cfl_safety","Safety factor relative to target error","",cfl->safety,&cfl->safety,NULL);CHKERRQ(ierr); 680298fd71SBarry Smith ierr = PetscOptionsBool("-ts_adapt_cfl_always_accept","Always accept the step regardless of whether local truncation error meets goal","",cfl->always_accept,&cfl->always_accept,NULL);CHKERRQ(ierr); 69ce94432eSBarry Smith if (!cfl->always_accept) SETERRQ(PetscObjectComm((PetscObject)adapt),PETSC_ERR_SUP,"step rejection not implemented yet"); 708d59e960SJed Brown ierr = PetscOptionsTail();CHKERRQ(ierr); 718d59e960SJed Brown PetscFunctionReturn(0); 728d59e960SJed Brown } 738d59e960SJed Brown 748d59e960SJed Brown #undef __FUNCT__ 758d59e960SJed Brown #define __FUNCT__ "TSAdaptCreate_CFL" 768d59e960SJed Brown /*MC 778d59e960SJed Brown TSADAPTCFL - CFL adaptive controller for time stepping 788d59e960SJed Brown 798d59e960SJed Brown Level: intermediate 808d59e960SJed Brown 818d59e960SJed Brown .seealso: TS, TSAdapt, TSSetAdapt() 828d59e960SJed Brown M*/ 838cc058d9SJed Brown PETSC_EXTERN PetscErrorCode TSAdaptCreate_CFL(TSAdapt adapt) 848d59e960SJed Brown { 858d59e960SJed Brown PetscErrorCode ierr; 868d59e960SJed Brown TSAdapt_CFL *a; 878d59e960SJed Brown 888d59e960SJed Brown PetscFunctionBegin; 89*b00a9115SJed Brown ierr = PetscNewLog(adapt,&a);CHKERRQ(ierr); 908d59e960SJed Brown adapt->data = (void*)a; 918d59e960SJed Brown adapt->ops->choose = TSAdaptChoose_CFL; 928d59e960SJed Brown adapt->ops->setfromoptions = TSAdaptSetFromOptions_CFL; 938d59e960SJed Brown adapt->ops->destroy = TSAdaptDestroy_CFL; 948d59e960SJed Brown 958d59e960SJed Brown a->safety = 0.9; 968d59e960SJed Brown a->always_accept = PETSC_FALSE; 978d59e960SJed Brown PetscFunctionReturn(0); 988d59e960SJed Brown } 99