Skip to content

[Draft] Adding new nml options to gradually increase timestep size #1295

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
wants to merge 4 commits into
base: develop
Choose a base branch
from
Draft

[Draft] Adding new nml options to gradually increase timestep size #1295

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
ca5ca12
working
abishekg7 Mar 7, 2025
ed68202
correct results with dt=95s init and changing to 60s
abishekg7 Mar 8, 2025
851d60c
adding new namelist options to slowly ramp up the timestep size from …
abishekg7 Mar 21, 2025
149a00d
some cleanup
abishekg7 Mar 25, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
10 changes: 10 additions & 0 deletions src/core_atmosphere/Registry.xml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -72,6 +72,11 @@
description="Model time step, seconds"
possible_values="Positive real values"/>

<nml_option name="config_dt_init" type="real" default_value="0.0"
units="s"
description="Model time step at the start of ramp up period (seconds)"
possible_values="Positive real values"/>

<nml_option name="config_calendar_type" type="character" default_value="gregorian" in_defaults="false"
units="-"
description="Simulation calendar type"
Expand All @@ -92,6 +97,11 @@
description="Length of model simulation"
possible_values="[DDD_]hh:mm:ss"/>

<nml_option name="config_dt_rampup_duration" type="character" default_value="none"
units="-"
description="Length of model ramp up period"
possible_values="[DDD_]hh:mm:ss"/>

<nml_option name="config_split_dynamics_transport" type="logical" default_value="true"
units="-"
description="Whether to super-cycle scalar transport"
Expand Down
11 changes: 7 additions & 4 deletions src/core_atmosphere/dynamics/mpas_atm_iau.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,7 +19,7 @@ module mpas_atm_iau
contains

!==================================================================================================
real (kind=RKIND) function atm_iau_coef(configs, itimestep, dt) result(wgt_iau)
real (kind=RKIND) function atm_iau_coef(configs, itimestep, dt, xtime_s) result(wgt_iau)
!==================================================================================================
! Compute the coefficient (or weight) for the IAU forcing at itimestep.

Expand All @@ -28,6 +28,7 @@ real (kind=RKIND) function atm_iau_coef(configs, itimestep, dt) result(wgt_iau)
type (mpas_pool_type), intent(in) :: configs
integer, intent(in) :: itimestep
real (kind=RKIND), intent(in) :: dt
real (kind=RKIND), intent(in) :: xtime_s

integer :: nsteps_iau ! Total number of time steps where the IAU forcing is applied.
logical, parameter :: debug = .false.
Expand Down Expand Up @@ -67,7 +68,8 @@ real (kind=RKIND) function atm_iau_coef(configs, itimestep, dt) result(wgt_iau)
nsteps_iau = nint(time_window_sec / dt)
!if(debug) call mpas_log_write('atm_compute_iau_coef: nsteps_iau = $i', intArgs=(/nsteps_iau/))

if(itimestep <= nsteps_iau) then
!if(itimestep <= nsteps_iau) then
if(xtime_s <= time_window_sec) then
!wgt_iau = 1./nsteps_iau
wgt_iau = 1.0_RKIND / time_window_sec
if(debug) call mpas_log_write('atm_compute_iau_coef: wgt_iau = $r', realArgs=(/wgt_iau/))
Expand All @@ -78,7 +80,7 @@ real (kind=RKIND) function atm_iau_coef(configs, itimestep, dt) result(wgt_iau)
end function atm_iau_coef

!==================================================================================================
subroutine atm_add_tend_anal_incr (configs, structs, itimestep, dt, tend_ru, tend_rtheta, tend_rho)
subroutine atm_add_tend_anal_incr (configs, structs, itimestep, dt, xtime_s, tend_ru, tend_rtheta, tend_rho)
!==================================================================================================

implicit none
Expand All @@ -87,6 +89,7 @@ subroutine atm_add_tend_anal_incr (configs, structs, itimestep, dt, tend_ru, ten
type (mpas_pool_type), intent(inout) :: structs
integer, intent(in) :: itimestep
real (kind=RKIND), intent(in) :: dt
real (kind=RKIND), intent(in) :: xtime_s
real (kind=RKIND), dimension(:,:), intent(inout) :: tend_ru
real (kind=RKIND), dimension(:,:), intent(inout) :: tend_rtheta
real (kind=RKIND), dimension(:,:), intent(inout) :: tend_rho
Expand Down Expand Up @@ -114,7 +117,7 @@ subroutine atm_add_tend_anal_incr (configs, structs, itimestep, dt, tend_ru, ten
! Compute weight for IAU forcing in this timestep, and return if weight
! is essentially zero.
!
wgt_iau = atm_iau_coef(configs, itimestep, dt)
wgt_iau = atm_iau_coef(configs, itimestep, dt, xtime_s)
if (wgt_iau <= 1.0e-12_RKIND) then
return
end if
Expand Down
19 changes: 12 additions & 7 deletions src/core_atmosphere/dynamics/mpas_atm_time_integration.F
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -631,6 +631,7 @@ subroutine mpas_atm_dynamics_finalize(domain)
end subroutine mpas_atm_dynamics_finalize


!subroutine atm_timestep(domain, dt, nowTime, itimestep, xtime_s, exchange_halo_group)
subroutine atm_timestep(domain, dt, nowTime, itimestep, exchange_halo_group)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Advance model state forward in time by the specified time step
Expand All @@ -647,6 +648,8 @@ subroutine atm_timestep(domain, dt, nowTime, itimestep, exchange_halo_group)
real (kind=RKIND), intent(in) :: dt
type (MPAS_Time_type), intent(in) :: nowTime
integer, intent(in) :: itimestep
!real (kind=RKIND), intent(in) :: xtime_s
real (kind=RKIND) :: xtime_s
procedure (halo_exchange_routine) :: exchange_halo_group


Expand All @@ -659,6 +662,7 @@ subroutine atm_timestep(domain, dt, nowTime, itimestep, exchange_halo_group)
type (mpas_pool_type), pointer :: state
character (len=StrKIND), pointer :: config_time_integration

xtime_s = 0.0_RKIND

clock => domain % clock
block => domain % blocklist
Expand All @@ -667,7 +671,7 @@ subroutine atm_timestep(domain, dt, nowTime, itimestep, exchange_halo_group)
call mpas_pool_get_config(block % configs, 'config_apply_lbcs', config_apply_lbcs)

if (trim(config_time_integration) == 'SRK3') then
call atm_srk3(domain, dt, itimestep, exchange_halo_group)
call atm_srk3(domain, dt, itimestep, xtime_s, exchange_halo_group)
else
call mpas_log_write('Unknown time integration option '//trim(config_time_integration), messageType=MPAS_LOG_ERR)
call mpas_log_write('Currently, only ''SRK3'' is supported.', messageType=MPAS_LOG_CRIT)
Expand All @@ -692,7 +696,7 @@ subroutine atm_timestep(domain, dt, nowTime, itimestep, exchange_halo_group)
end subroutine atm_timestep


subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
subroutine atm_srk3(domain, dt, itimestep, xtime_s, exchange_halo_group)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Advance model state forward in time by the specified time step using
! time-split RK3 scheme
Expand All @@ -711,6 +715,7 @@ subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
type (domain_type), intent(inout) :: domain
real (kind=RKIND), intent(in) :: dt
integer, intent(in) :: itimestep
real (kind=RKIND), intent(in) :: xtime_s
procedure (halo_exchange_routine) :: exchange_halo_group

integer :: thread
Expand Down Expand Up @@ -742,7 +747,6 @@ subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
logical, pointer :: config_scalar_advection
logical, pointer :: config_positive_definite
logical, pointer :: config_monotonic
real (kind=RKIND), pointer :: config_dt
character (len=StrKIND), pointer :: config_microp_scheme
character (len=StrKIND), pointer :: config_convection_scheme

Expand Down Expand Up @@ -781,8 +785,7 @@ subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
call mpas_pool_get_config(block % configs, 'config_time_integration_order', config_time_integration_order)
call mpas_pool_get_config(block % configs, 'config_scalar_advection', config_scalar_advection)
call mpas_pool_get_config(block % configs, 'config_positive_definite', config_positive_definite)
call mpas_pool_get_config(block % configs, 'config_monotonic', config_monotonic)
call mpas_pool_get_config(block % configs, 'config_dt', config_dt)
call mpas_pool_get_config(block % configs, 'config_monotonic', config_monotonic)
call mpas_pool_get_config(block % configs, 'config_IAU_option', config_IAU_option)
! config variables for dynamics-transport splitting, WCS 18 November 2014
call mpas_pool_get_config(block % configs, 'config_split_dynamics_transport', config_split_dynamics_transport)
Expand Down Expand Up @@ -870,6 +873,8 @@ subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
! Initialize RK weights
!

call mpas_log_write('--- atm_srk3 time-step DT = $r',realArgs=(/dt/))

dynamics_split = config_dynamics_split
if (config_split_dynamics_transport) then
dt_dynamics = dt/real(dynamics_split)
Expand Down Expand Up @@ -974,7 +979,7 @@ subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
! IAU - Incremental Analysis Update
!
if (trim(config_IAU_option) /= 'off') then
call atm_add_tend_anal_incr(block % configs, block % structs, itimestep, dt, &
call atm_add_tend_anal_incr(block % configs, block % structs, itimestep, dt, xtime_s, &
tend_ru_physics, tend_rtheta_physics, tend_rho_physics)
end if

Expand Down Expand Up @@ -1482,7 +1487,7 @@ subroutine atm_srk3(domain, dt, itimestep, exchange_halo_group)
!requires that the subroutine atm_advance_scalars_mono was called on the third Runge Kutta step, so that a halo
!update for the scalars at time_levs(1) is applied. A halo update for the scalars at time_levs(2) is done above.
if (config_monotonic) then
rqvdynten(:,:) = ( scalars_2(index_qv,:,:) - scalars_1(index_qv,:,:) ) / config_dt
rqvdynten(:,:) = ( scalars_2(index_qv,:,:) - scalars_1(index_qv,:,:) ) / dt
else
rqvdynten(:,:) = 0._RKIND
end if
Expand Down
Loading