diff --git a/bld/build-namelist b/bld/build-namelist
index 45a8e4f27b..e75ba11b86 100755
--- a/bld/build-namelist
+++ b/bld/build-namelist
@@ -4355,9 +4355,6 @@ if ($dyn =~ /se/) {
se_vert_remap_uvTq_alg
se_vert_remap_tracer_alg
se_fvm_supercycling
- se_fvm_supercycling_jet
- se_kmin_jet
- se_kmax_jet
se_molecular_diff
se_pgf_formulation
se_dribble_in_rsplit_loop
diff --git a/bld/namelist_files/namelist_defaults_cam.xml b/bld/namelist_files/namelist_defaults_cam.xml
index 38975992be..3d39f3c2b9 100644
--- a/bld/namelist_files/namelist_defaults_cam.xml
+++ b/bld/namelist_files/namelist_defaults_cam.xml
@@ -3240,17 +3240,13 @@
3.22D0
- 3
+-1
3
2
- 3
- 1
10
9
8
5
- 2
- 3
2
2
1
@@ -3262,16 +3258,15 @@
18
18
30
- 2
+ 3
+ 3
- 1
- 3
+-1
3
2
4
20
60
- 4
2
2
5
@@ -3282,12 +3277,13 @@
4
40
180
- 10
3
3
6
6
6
+ 1
+ 1
1
@@ -3372,7 +3368,7 @@
4
4
- 3
+-1
3
2
2
@@ -3380,17 +3376,20 @@
2
3
5
- 3
6
6
5
5
5
+ 3
+ 3
+ 3
+ 3
+ 3
+ 3
+ 3
-1
--1
--1
--1
3
diff --git a/bld/namelist_files/namelist_definition.xml b/bld/namelist_files/namelist_definition.xml
index 8d19021136..8c4e0e943c 100644
--- a/bld/namelist_files/namelist_definition.xml
+++ b/bld/namelist_files/namelist_definition.xml
@@ -8601,12 +8601,14 @@ Default: 0 (i.e., not used)
Number of hyperviscosity subcycles per dynamics timestep.
+Set to -1 to have the value selected automatically at runtime from the CFL stability limits.
Default: Set by build-namelist
Number of hyperviscosity subcycles per dynamics timestep in sponge del2 diffusion.
+Set to -1 to have the value selected automatically at runtime from the CFL stability limits.
Default: Set by build-namelist
@@ -8756,31 +8758,6 @@ Number of SE rsplit time-steps CSLAM supercycles
rsplit/se_fvm_supercycling must be an integer
-
-Number of SE rsplit time-steps CSLAM supercycles in the jet region
-(which is specified by kmin_jet and kmax_jet)
-rsplit/se_fvm_supercycling_jet must be an integer
-se_fvm_supercycling_jet/se_fvm_supercyling must be an integer
-
-
-
- Min level index where CSLAM runs with se_fvm_supercycling_jet
- (if se_fvm_supercycling_jet.ne.se_fvm_supercycling) or
- min index where are Courant number increment is active
- (se_large_Courant_incr=.true.)
-
-
-
- Max level index where CSLAM runs with se_fvm_supercycling_jet
- (if se_fvm_supercycling_jet.ne.se_fvm_supercycling) or
- max index where are Courant number increment is active
- (se_large_Courant_incr=.true.)
-
-
-
Tracer advection is done every qsplit dynamics timesteps.
@@ -8796,6 +8773,7 @@ Default: FALSE
Vertically lagrangian code vertically remaps every rsplit tracer timesteps.
+Set to -1 to have the value selected automatically at runtime from the CFL stability limits.
Default: Set by build-namelist.
@@ -8881,34 +8859,6 @@ differences into the simulation.
Default: .false.
-
-
-Nudging factor for prescribed winds in SE dycore
-Units: 1/sec
-Default: 2e-5
-
-
-
-Nudging factor for prescribed temperature in SE dycore
-Units: 1/sec
-Default: 2e-5
-
-
-
-Nudging factor for prescribed surface pressure in SE dycore
-Units: 1/sec
-Default: 0.0
-
-
-
-Switch to turn on/off time evolution of dynamics nudging
-Default: 0
-
-
diff --git a/doc/ChangeLog b/doc/ChangeLog
index 320c25ad23..7d57d56980 100644
--- a/doc/ChangeLog
+++ b/doc/ChangeLog
@@ -1,5 +1,119 @@
===============================================================
+Tag name:
+Originator(s): pel
+Date: August 22, 2026
+One-line Summary: Automatic dynamics-splitting defaults for the SE dycore; remove unused SE namelist options
+Github PR URL:
+
+Purpose of changes (include the issue number and title text for each relevant GitHub issue):
+
+ - #1655 Add automatic dynamics-splitting defaults to the SE dycore
+
+ - #969 Clean-up unused and/or unsupported namelist options in SE dycore
+
+Describe any changes made to build system:
+ - bld/build-namelist: removed se_fvm_supercycling_jet, se_kmin_jet, and
+ se_kmax_jet from the SE dycore namelist variable list.
+
+Describe any changes made to the namelist:
+
+ Changed defaults (automatic dynamics splitting):
+ - se_rsplit, se_hypervis_subcycle, se_hypervis_subcycle_sponge default
+ changed to -1 ("automatic"); the effective value is computed at runtime
+ from the CFL stability limits. WACCM, WACCM-X, variable-resolution, and
+ ne3 grids are pinned to their previous explicit values, so those
+ configurations are unchanged.
+
+ Removed namelist options (unused/untested):
+ - se_fvm_supercycling_jet, se_kmin_jet, se_kmax_jet
+ - se_met_nudge_u, se_met_nudge_t, se_met_nudge_p, se_met_tevolve
+
+List any changes to the defaults for the boundary datasets: N/A
+
+Describe any substantial timing or memory changes: N/A
+
+Code reviewed by:
+
+List all files eliminated:
+
+List all files added and what they do:
+
+List all existing files that have been modified, and describe the changes:
+
+M bld/build-namelist
+ - remove se_fvm_supercycling_jet/se_kmin_jet/se_kmax_jet from the
+ SE dycore variable list
+M bld/namelist_files/namelist_definition.xml
+ - document the -1 "automatic" sentinel for the split variables;
+ remove the jet-supercycling and se_met_nudge_* entries
+M bld/namelist_files/namelist_defaults_cam.xml
+ - set se_rsplit/se_hypervis_subcycle/se_hypervis_subcycle_sponge
+ defaults to -1; pin WACCM/WACCM-X/variable-resolution/ne3 to
+ their prior values; remove jet-supercycling defaults
+M src/dynamics/se/dyn_comp.F90
+ - accept the rsplit == -1 sentinel in dyn_readnl; remove the
+ se_fvm_supercycling_jet / se_kmin_jet / se_kmax_jet and
+ variable_nsplit plumbing
+M src/dynamics/se/dyn_grid.F90
+ - dyn_grid_init: resolve automatic (-1) rsplit and the fvm
+ supercycling default from the CFL stability limits
+M src/dynamics/se/dycore/global_norms_mod.F90
+ - compute the automatic hyperviscosity subcycle counts from the
+ grid stability constants in print_cfl
+M src/dynamics/se/dycore/prim_driver_mod.F90
+ - remove variable_nsplit machinery; the CFL dt computation is now
+ done in print_cfl
+M src/dynamics/se/dycore/prim_state_mod.F90
+ - remove variable_nsplit / nsplit_baseline / rsplit_baseline
+M src/dynamics/se/dycore/se_dyn_time_mod.F90
+ - remove variable_nsplit plumbing
+M src/dynamics/se/dycore/control_mod.F90
+ - remove 'variable_nsplit' and other unused control variables.
+M src/dynamics/se/dycore/dimensions_mod.F90
+ - remove kmin_jet, kmax_jet, and fvm_supercycling_jet declarations
+M src/dynamics/se/dycore/fvm_mod.F90
+ - remove jet-supercycling logic
+M src/dynamics/se/dycore/fvm_consistent_se_cslam.F90
+ - large_Courant_incr now operates over the full column (1..nlev)
+ instead of the removed [kmin_jet, kmax_jet] jet region
+
+If there were any failures reported from running test_driver.sh on any test
+platform, and checkin with these failures has been OK'd by the gatekeeper,
+then copy the lines from the td.*.status files for the failed tests to the
+appropriate machine below. All failed tests must be justified.
+
+derecho/intel/aux_cam:
+
+derecho/nvhpc/aux_cam:
+
+izumi/nag/aux_cam:
+
+izumi/gnu/aux_cam:
+
+CAM tag used for the baseline comparison tests if different than previous
+tag: cam6_4_196
+
+Summarize any changes to answers, i.e.,
+- what code configurations: SE dycore. Low/mid-top uniform-grid
+ configurations whose automatically selected se_hypervis_subcycle is smaller
+ than the previous hard-coded value change answers. Observed in the aux_cam
+ SE suite for:
+ ERI_D_Ln18.ne16pg3_ne16pg3_mt232.FHIST_C4 (hypervis_subcycle 3 -> 1)
+ ERP_Ln9.ne30pg3_ne30pg3_mg17.FCnudged (hypervis_subcycle 3 -> 2)
+ SMS_Ld1.ne30pg3_ne30pg3_mg17.FC2010climo (hypervis_subcycle 3 -> 2)
+ All other SE configurations (including all WACCM, WACCM-X,
+ variable-resolution, and ne3 grids, which are pinned) are bit-for-bit.
+ The namelist-option cleanup (#969) is bit-for-bit by construction (removal
+ of unused/dead code and orphaned definition entries).
+- what platforms/compilers: derecho/intel
+- nature of change: larger than roundoff (fewer hyperviscosity subcycles where
+ the CFL stability estimate permits); same climate .
+
+===============================================================
+===============================================================
+
Tag name: cam6_4_197
Originator(s): jimmielin
Date: August 22, 2026
diff --git a/src/dynamics/se/dycore/control_mod.F90 b/src/dynamics/se/dycore/control_mod.F90
index a5b2c95db0..d18c45781b 100644
--- a/src/dynamics/se/dycore/control_mod.F90
+++ b/src/dynamics/se/dycore/control_mod.F90
@@ -13,7 +13,6 @@ module control_mod
integer, public :: tstep_type= 0 ! 0 = leapfrog
! 1 = RK (foward-in-time)
- integer, public :: rk_stage_user = 0 ! number of RK stages to use
integer, public :: ftype = 2 ! Forcing Type
integer, public :: ftype_conserve = 1 !conserve momentum (dp*u)
integer, public :: dribble_in_rsplit_loop = 0
@@ -22,38 +21,26 @@ module control_mod
integer, public :: qsplit = 1 ! ratio of dynamics tsteps to tracer tsteps
integer, public :: rsplit =-1 ! for vertically lagrangian dynamics, apply remap
! every rsplit tracer timesteps
- logical, public :: variable_nsplit=.false.
logical, public :: refined_mesh
- integer, public :: vert_remap_q_alg = 10
-
-
integer, public :: cubed_sphere_map = -1 ! -1 = chosen at run time
! 0 = equi-angle Gnomonic (default)
! 1 = equi-spaced Gnomonic (not yet coded)
! 2 = element-local projection (for var-res)
! 3 = parametric (not yet coded)
-!tolerance to define smth small, was introduced for lim 8 in 2d and 3d
- real (kind=r8), public, parameter :: tol_limiter=1.0e-13_r8
-
integer , public :: limiter_option = 0
integer , public :: partmethod ! partition methods
character(len=MAX_STRING_LEN) , public :: topology ! options: "cube" is supported
- integer , public :: tasknum
- integer , public :: remapfreq ! remap frequency of synopsis of system state (steps)
- character(len=MAX_STRING_LEN) :: remap_type ! selected remapping option
integer , public :: statefreq ! output frequency of synopsis of system state (steps)
integer , public :: runtype
- integer , public :: timerdetail
integer , public :: numnodes
integer , public :: multilevel
character(len=MAX_STRING_LEN) , public :: columnpackage
- integer , public :: maxits ! max iterations of solver
real (kind=r8), public :: tol ! solver tolerance (convergence criteria)
integer , public :: fine_ne = -1 ! set for refined exodus meshes (variable viscosity)
@@ -83,8 +70,6 @@ module control_mod
integer, public :: hypervis_subcycle=1 ! number of subcycles for hyper viscsosity timestep
integer, public :: hypervis_subcycle_sponge=1 ! number of subcycles for hyper viscsosity timestep in sponge
integer, public :: hypervis_subcycle_q=1 ! number of subcycles for hyper viscsosity timestep on TRACERS
- integer, public :: psurf_vis = 0 ! 0 = use laplace on eta surfaces
- ! 1 = use (approx.) laplace on p surfaces
real (kind=r8), public :: hypervis_power=0 ! if not 0, use variable hyperviscosity based on element area
real (kind=r8), public :: hypervis_scaling=0 ! use tensor hyperviscosity
diff --git a/src/dynamics/se/dycore/dimensions_mod.F90 b/src/dynamics/se/dycore/dimensions_mod.F90
index 8da7f71e49..693dc790aa 100644
--- a/src/dynamics/se/dycore/dimensions_mod.F90
+++ b/src/dynamics/se/dycore/dimensions_mod.F90
@@ -55,9 +55,7 @@ module dimensions_mod
integer, parameter, public :: ubc = nc+nhc
logical, public :: large_Courant_incr
- integer, public :: kmin_jet,kmax_jet !min and max level index for the jet
- integer, public :: fvm_supercycling
- integer, public :: fvm_supercycling_jet
+ integer, public :: fvm_supercycling
integer, allocatable, public :: kord_tr(:), kord_tr_cslam(:)
diff --git a/src/dynamics/se/dycore/fvm_consistent_se_cslam.F90 b/src/dynamics/se/dycore/fvm_consistent_se_cslam.F90
index 3e56498b38..921b44d1cf 100644
--- a/src/dynamics/se/dycore/fvm_consistent_se_cslam.F90
+++ b/src/dynamics/se/dycore/fvm_consistent_se_cslam.F90
@@ -3,7 +3,6 @@ module fvm_consistent_se_cslam
use shr_kind_mod, only: r8=>shr_kind_r8
use dimensions_mod, only: nc, nhe, nlev, ntrac, np, nhr, nhc, ngpc, ns, nht
use dimensions_mod, only: irecons_tracer
- use dimensions_mod, only: kmin_jet,kmax_jet
use cam_abortutils, only: endrun
use cam_logfile, only: iulog
@@ -195,14 +194,11 @@ subroutine run_consistent_se_cslam(elem,fvm,hybrid,dt_fvm,tl,nets,nete,hvcoord,&
!
if (large_Courant_incr) then
if(FVM_TIMERS) call t_startf('fvm:fill_halo_fvm:large_Courant')
- !if (kmin_jetkmax) then
- ! call endrun('ERROR: kmax_jet must be .le. kmax passed to run_consistent_se_cslam')
- !end if
! Determine the extent of the JET that is owned by this thread
- ActiveJetThread = threadOwnsVertLevel(hybridnew,kmin_jet) .or. threadOwnsVertLevel(hybridnew,kmax_jet)
- kmin_jet_local = max(kmin_jet,kmin)
- kmax_jet_local = min(kmax_jet,kmax)
- klev = kmax_jet-kmin_jet+1
+ ActiveJetThread = threadOwnsVertLevel(hybridnew,1) .or. threadOwnsVertLevel(hybridnew,nlev)
+ kmin_jet_local = max(1,kmin)
+ kmax_jet_local = min(nlev,kmax)
+ klev = nlev
call fill_halo_fvm(ghostbufQ1,elem,fvm,hybridnew,nets,nete,1,kmin_jet_local,kmax_jet_local,klev,active=ActiveJetThread)
if(FVM_TIMERS) call t_stopf('fvm:fill_halo_fvm:large_Courant')
if(FVM_TIMERS) call t_startf('fvm:large_Courant_number_increment')
@@ -259,7 +255,6 @@ end subroutine run_consistent_se_cslam
subroutine swept_flux(elem,fvm,ilev,ctracer,irecons_tracer_actual,gsweights,gspts)
use fvm_analytic_mod , only: get_high_order_weights_over_areas
- use dimensions_mod, only : kmin_jet,kmax_jet
implicit none
type (element_t) , intent(in) :: elem
type (fvm_struct), intent(inout):: fvm
@@ -586,13 +581,6 @@ subroutine swept_flux(elem,fvm,ilev,ctracer,irecons_tracer_actual,gsweights,gspt
end if
end do
fvm%se_flux(i,j,iside,ilev) = mass_flux_se(i,j,iside)-flux
- if (fvm%se_flux(i,j,iside,ilev)>1.0E-13_r8.and.(ilevkmax_jet)) then
- write(iulog,*) "CN excess flux outside of pre-scribed jet region"
- write(iulog,*) "Increase jet region with kmin_jet and kmax_jet ",&
- ilev,fvm%se_flux(i,j,iside,ilev),mass_flux_se(i,j,iside),flux,flowcase,&
- kmin_jet,kmax_jet
- call endrun('ERROR in CSLAM: local Courant number is > 1; Increase kmin_jet/kmax_jet?')
- end if
fvm%dp_fvm(i ,j ,ilev ) = fvm%dp_fvm(i ,j ,ilev )-flux
fvm% c(i ,j ,ilev,1:ntrac) = fvm% c(i ,j ,ilev,1:ntrac)-flux_tracer(1:ntrac)
diff --git a/src/dynamics/se/dycore/fvm_mod.F90 b/src/dynamics/se/dycore/fvm_mod.F90
index e2f311ee81..81e2937b40 100644
--- a/src/dynamics/se/dycore/fvm_mod.F90
+++ b/src/dynamics/se/dycore/fvm_mod.F90
@@ -31,10 +31,7 @@ module fvm_mod
type (EdgeBuffer_t), public :: ghostBufQnhc_h
type (EdgeBuffer_t), public :: ghostBufQ1_h
type (EdgeBuffer_t), public :: ghostBufQ1_vh
-! type (EdgeBuffer_t), private :: ghostBufFlux_h
type (EdgeBuffer_t), public :: ghostBufFlux_vh
- type (EdgeBuffer_t), public :: ghostBufQnhcJet_h
- type (EdgeBuffer_t), public :: ghostBufFluxJet_h
type (EdgeBuffer_t), public :: ghostBufPG_s
type (EdgeBuffer_t), public :: ghostBuf_cslam2gll
@@ -291,10 +288,8 @@ subroutine fvm_init1(par,elem)
use cam_logfile, only: iulog
use control_mod, only: rsplit
use dimensions_mod, only: qsize, qsize_d
- use dimensions_mod, only: fvm_supercycling, fvm_supercycling_jet
+ use dimensions_mod, only: fvm_supercycling
use dimensions_mod, only: nc,nhe, nhc, nlev,ntrac, ntrac_d,ns, nhr, use_cslam
- use dimensions_mod, only: large_Courant_incr
- use dimensions_mod, only: kmin_jet,kmax_jet
type (parallel_t) :: par
type (element_t),intent(inout) :: elem(:)
@@ -318,8 +313,6 @@ subroutine fvm_init1(par,elem)
!
! PARAMETER ERROR CHECKING
!
- if (kmin_jet>kmax_jet) &
- call endrun("PARAMETER ERROR for fvm: kmin_jet must be < kmax_jet")
if (ntrac>ntrac_d) &
call endrun("PARAMETER ERROR for fvm: ntrac > ntrac_d")
@@ -332,23 +325,6 @@ subroutine fvm_init1(par,elem)
call endrun("PARAMETER ERROR for fvm: mod(rsplit,fvm_supercycling)<>0")
endif
- if (qsize>0.and.mod(rsplit,fvm_supercycling_jet).ne.0) then
- if (par%masterproc) then
- write(iulog,*)'cannot supercycle fvm tracers with respect to se tracers'
- write(iulog,*)'with this choice of rsplit =',rsplit
- write(iulog,*)'rsplit must be a multiple of fvm_supercycling_jet=',fvm_supercycling_jet
- end if
- call endrun("PARAMETER ERROR for fvm: mod(rsplit,fvm_supercycling_jet)<>0")
- endif
-
- if (large_Courant_incr.and.(fvm_supercycling.ne.fvm_supercycling_jet)) then
- if (par%masterproc) then
- write(iulog,*)'Large Courant number increment requires no level dependent supercycling'
- write(iulog,*)'i.e. fvm_supercycling must be equal to fvm_supercycling_jet'
- end if
- call endrun("PARAMETER ERROR for fvm: large_courant_incr requires fvm_supercycling=fvm_supercycling_jet")
- endif
-
if (par%masterproc) then
write(iulog,*) " "
write(iulog,*) "Done Tracer transport scheme information "
@@ -446,8 +422,6 @@ subroutine fvm_init2(elem,fvm,hybrid,nets,nete)
use bndry_mod, only: compute_ghost_corner_orientation
use dimensions_mod, only: nlev, nc, nhc, nhe, ntrac, ntrac_d, np
use dimensions_mod, only: nhc_phys, fv_nphys
- use dimensions_mod, only: fvm_supercycling, fvm_supercycling_jet
- use dimensions_mod, only: kmin_jet,kmax_jet
use hycoef, only: hyai, hybi, ps0
use derivative_mod, only: subcell_integration
use air_composition, only: thermodynamic_active_species_num
@@ -492,7 +466,7 @@ subroutine fvm_init2(elem,fvm,hybrid,nets,nete)
call initghostbuffer(hybrid%par,ghostBufQnhc_t1,elem,nlev, nhc,nc,nthreads=1)
call initghostbuffer(hybrid%par,ghostBufQnhc_h,elem,nlev*(ntrac+1),nhc,nc,nthreads=horz_num_threads)
call initghostbuffer(hybrid%par,ghostBufQnhc_vh,elem,nlev*(ntrac+1),nhc,nc,nthreads=vert_num_threads*horz_num_threads)
- klev = kmax_jet-kmin_jet+1
+ klev = nlev
call initghostbuffer(hybrid%par,ghostBufQ1_h,elem,klev*(ntrac+1),1,nc,nthreads=horz_num_threads)
call initghostbuffer(hybrid%par,ghostBufQ1_vh,elem,klev*(ntrac+1),1,nc,nthreads=vert_num_threads*horz_num_threads)
! call initghostbuffer(hybrid%par,ghostBufFlux_h,elem,4*nlev,nhe,nc,nthreads=horz_num_threads)
@@ -506,15 +480,6 @@ subroutine fvm_init2(elem,fvm,hybrid,nets,nete)
else
call initghostbuffer(hybrid%par,ghostBufPG_s,elem,nlev*3,nhc_phys,fv_nphys,nthreads=1)
end if
-
- if (fvm_supercycling.ne.fvm_supercycling_jet) then
- !
- ! buffers for running different fvm time-steps in the jet region
- !
- klev = kmax_jet-kmin_jet+1
- call initghostbuffer(hybrid%par,ghostBufQnhcJet_h,elem,klev*(ntrac+1),nhc,nc,nthreads=horz_num_threads)
- call initghostbuffer(hybrid%par,ghostBufFluxJet_h,elem,4*klev,nhe,nc,nthreads=horz_num_threads)
- end if
end subroutine fvm_init2
diff --git a/src/dynamics/se/dycore/global_norms_mod.F90 b/src/dynamics/se/dycore/global_norms_mod.F90
index 9699bfaa2f..80836e681a 100644
--- a/src/dynamics/se/dycore/global_norms_mod.F90
+++ b/src/dynamics/se/dycore/global_norms_mod.F90
@@ -17,6 +17,8 @@ module global_norms_mod
public :: linf_vnorm
public :: print_cfl
+ public :: auto_rsplit
+ public :: set_global_max_normDinv
public :: global_integral
public :: global_integrals_general
public :: wrap_repro_sum
@@ -24,6 +26,14 @@ module global_norms_mod
private :: global_maximum
type (EdgeBuffer_t), private :: edgebuf
+ real(kind=r8), parameter :: ugw = 342.0_r8 ! max gravity-wave speed [m/s]
+
+ ! Global max of the mesh metric elem(:)%normDinv. Set once by
+ ! set_global_max_normDinv (called from dyn_grid_init, before the dynamics
+ ! timestep is initialized) so that "automatic" (se_rsplit=-1) subcycling can
+ ! be resolved there.
+ real(kind=r8) :: global_max_normDinv = -1.0_r8
+
interface global_integral
module procedure global_integral_elem
module procedure global_integral_fvm
@@ -209,26 +219,7 @@ function global_integral_fvm(fvm,fld,hybrid,npts,nets,nete) result(I_sphere)
end function global_integral_fvm
-!------------------------------------------------------------------------------------
-
- ! ================================
- ! print_cfl:
- !
- ! Calculate / output CFL info
- ! (both advective and based on
- ! viscosity or hyperviscosity)
- !
- ! ================================
-
- subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
- dt_remap_actual,dt_tracer_fvm_actual,dt_tracer_se_actual,&
- dt_dyn_actual,dt_dyn_visco_actual,dt_dyn_del2_actual,dt_tracer_visco_actual,dt_phys)
- !
- ! estimate various CFL limits
- ! also, for variable resolution viscosity coefficient, make sure
- ! worse viscosity CFL (given by dtnu) is not violated by reducing
- ! viscosity coefficient in regions where CFL is violated
- !
+ subroutine print_cfl(elem,hybrid,nets,nete,ptop,pmid)
use hybrid_mod, only: hybrid_t
use element_mod, only: element_t
use dimensions_mod, only: np,ne,nelem,nc,nhe,use_cslam,nlev,large_Courant_incr
@@ -241,24 +232,28 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
use control_mod, only: nu, nu_div, nu_q, nu_p, nu_t, nu_top, fine_ne, max_hypervis_courant
use control_mod, only: tstep_type, hypervis_power, hypervis_scaling
use control_mod, only: sponge_del4_nu_div_fac, sponge_del4_nu_fac, sponge_del4_lev
+ use control_mod, only: rsplit, qsplit, hypervis_subcycle, hypervis_subcycle_sponge, hypervis_subcycle_q
+ use se_dyn_time_mod,only: tstep
+ use time_manager, only: get_step_size
use cam_abortutils, only: endrun
use parallel_mod, only: global_shared_buf, global_shared_sum
use edge_mod, only: initedgebuffer, FreeEdgeBuffer, edgeVpack, edgeVunpack
use bndry_mod, only: bndry_exchange
use mesh_mod, only: MeshUseMeshFile
- use dimensions_mod, only: ksponge_end, kmvis_ref, kmcnd_ref,rho_ref
+ use dimensions_mod, only: ksponge_end, kmvis_ref, kmcnd_ref,rho_ref, fvm_supercycling
use physconst, only: cpair
use std_atm_profile,only: std_atm_height
type(element_t) , intent(inout) :: elem(:)
integer , intent(in) :: nets,nete
type (hybrid_t) , intent(in) :: hybrid
- real (kind=r8), intent(in) :: dtnu, ptop, pmid(nlev)
+ real (kind=r8), intent(in) :: ptop, pmid(nlev)
!
- ! actual time-steps
+ ! actual time-steps (computed below, after "automatic" (se_*=-1) subcycling is resolved)
!
- real (kind=r8), intent(in) :: dt_remap_actual,dt_tracer_fvm_actual,dt_tracer_se_actual,&
- dt_dyn_actual,dt_dyn_visco_actual,dt_dyn_del2_actual, &
- dt_tracer_visco_actual, dt_phys
+ real (kind=r8) :: dtnu, dtime
+ real (kind=r8) :: dt_remap_actual,dt_tracer_fvm_actual,dt_tracer_se_actual,&
+ dt_dyn_actual,dt_dyn_visco_actual,dt_dyn_del2_actual, &
+ dt_tracer_visco_actual, dt_phys
! Element statisics
real (kind=r8) :: max_min_dx,min_min_dx,min_max_dx,max_unif_dx ! used for normalizing scalar HV
@@ -269,7 +264,7 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
real (kind=r8) :: normDinv_hypervis
real (kind=r8) :: x, y, noreast, nw, se, sw
real (kind=r8), dimension(np,np,nets:nete) :: zeta
- real (kind=r8) :: lambda_max, lambda_vis, min_gw, lambda,umax, ugw
+ real (kind=r8) :: lambda_max, lambda_vis, min_gw, lambda,umax
real (kind=r8) :: scale1, max_laplace,z(nlev)
integer :: ie, i, j, rowind, colind, k
type (quadrature_t) :: gp
@@ -285,33 +280,7 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
logical :: top_000_032km, top_032_042km, top_042_090km, top_090_140km, top_140_600km ! model top location ranges
logical :: nu_set,div_set,lev_set
- ! Eigenvalues calculated by folks at UMich (Paul U & Jared W)
- select case (np)
- case (2)
- lambda_max = 0.5_r8
- lambda_vis = 0.0_r8 ! need to compute this
- case (3)
- lambda_max = 1.5_r8
- lambda_vis = 12.0_r8
- case (4)
- lambda_max = 2.74_r8
- lambda_vis = 30.0_r8
- case (5)
- lambda_max = 4.18_r8
- lambda_vis = 91.6742_r8
- case (6)
- lambda_max = 5.86_r8
- lambda_vis = 190.1176_r8
- case (7)
- lambda_max = 7.79_r8
- lambda_vis = 374.7788_r8
- case (8)
- lambda_max = 10.0_r8
- lambda_vis = 652.3015_r8
- case DEFAULT
- lambda_max = 0.0_r8
- lambda_vis = 0.0_r8
- end select
+ call se_stability_constants(ptop, lambda_max, lambda_vis, s_rk, umax)
if ((lambda_max.eq.0_r8).and.(hybrid%masterthread)) then
print*, "lambda_max not calculated for NP = ",np
@@ -322,6 +291,10 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
print*, "Estimate of viscous CFLs will be incorrect"
end if
+ dtime = get_step_size()
+ ! timestep seen by the viscosity operators (ignores subcycling)
+ dtnu = max(tstep*max(nu,nu_div), tstep*qsplit*nu_q)
+
do ie=nets,nete
elem(ie)%variable_hyperviscosity = 1.0_r8
end do
@@ -370,7 +343,7 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
min_area = ParallelMin(min_area,hybrid)
max_area = ParallelMax(max_area,hybrid)
min_normDinv = ParallelMin(min_normDinv,hybrid)
- max_normDinv = ParallelMax(max_normDinv,hybrid)
+ max_normDinv = global_max_normDinv ! reduced once in dyn_grid_init (set_global_max_normDinv)
min_min_dx = ParallelMin(min_min_dx,hybrid)
max_min_dx = ParallelMax(max_min_dx,hybrid)
min_max_dx = ParallelMin(min_max_dx,hybrid)
@@ -652,15 +625,6 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
if (sponge_del4_nu_div_fac<0) sponge_del4_nu_div_fac = 1.0_r8
end if
- ! set max wind speed for diagnostics
- umax = 120.0_r8
- if (top_042_090km) then
- umax = 240._r8
- else if (top_090_140km) then
- umax = 300._r8
- else if (top_140_600km) then
- umax = 800._r8
- end if
!
! Log sponge layer configuration
!
@@ -720,16 +684,12 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
! S=time-step stability region (i.e. advection w/leapfrog: S=1, viscosity w/forward Euler: S=2)
!
if (tstep_type==1) then
- S_rk = 2.0_r8
rk_str = ' * RK2-SSP 3 stage (same as tracers)'
elseif (tstep_type==2) then
- S_rk = 2.0_r8
rk_str = ' * classic RK3'
elseif (tstep_type==3) then
- S_rk = 2.0_r8
rk_str = ' * Kinnmark&Gray RK4'
elseif (tstep_type==4) then
- S_rk = 3.0_r8
rk_str = ' * Kinnmark&Gray RK3 5 stage (3rd order)'
end if
if (hybrid%masterthread) then
@@ -745,8 +705,6 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
S_hypervis = 2.0_r8 !using forward Euler for hyperviscosity
S_rk_tracer = 2.0_r8
- ugw = 342.0_r8 !max gravity wave speed
-
dt_max_adv = S_rk/(umax*max_normDinv*lambda_max*ra)
dt_max_gw = S_rk/(ugw*max_normDinv*lambda_max*ra)
dt_max_tracer_se = S_rk_tracer*min_gw/(umax*max_normDinv*ra)
@@ -767,6 +725,24 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
max_laplace = MAX(max_laplace,MAXVAL(kmcnd_ref(:)/(cpair*rho_ref(:))))
dt_max_laplacian_top = 1.0_r8/(max_laplace*((ra*max_normDinv)**2)*lambda_vis)
+ ! Resolve "automatic" (se_*=-1) hyperviscosity subcycling from the stability
+ ! limits above. (rsplit=-1 is resolved earlier, in dyn_grid_init via
+ ! auto_rsplit)
+ !
+ if (hypervis_subcycle==-1) hypervis_subcycle = max(1, ceiling(tstep/(1.2_r8*dt_max_hypervis)))
+ if (hypervis_subcycle_sponge==-1) hypervis_subcycle_sponge = max(1, ceiling(tstep/dt_max_laplacian_top))
+ !
+ ! actual time-steps with final subcycling (printed below)
+ !
+ dt_remap_actual = tstep*qsplit*rsplit
+ dt_tracer_fvm_actual = tstep*qsplit*fvm_supercycling
+ dt_tracer_se_actual = tstep*qsplit
+ dt_dyn_actual = tstep
+ dt_dyn_visco_actual = tstep/hypervis_subcycle
+ dt_dyn_del2_actual = tstep/hypervis_subcycle_sponge
+ dt_tracer_visco_actual = tstep*qsplit/hypervis_subcycle_q
+ dt_phys = dtime
+
if (hybrid%masterthread) then
write(iulog,'(a,f10.2,a)') ' '
write(iulog,'(a,f10.2,a)') 'Estimates for maximum stable and actual time-steps for different aspects of algorithm:'
@@ -811,6 +787,8 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
end if
end do
write(iulog,*) ' '
+ write(iulog,'(a,3i6)') 'subcycling: rsplit, hypervis_subcycle, hypervis_subcycle_sponge =', &
+ rsplit,hypervis_subcycle,hypervis_subcycle_sponge
if (hypervis_power /= 0) then
write(iulog,'(a,3e11.4)')'Scalar hyperviscosity (dynamics): ave,min,max = ', &
nu*(/avg_hypervis**2,min_hypervis**2,max_hypervis**2/)
@@ -819,6 +797,86 @@ subroutine print_cfl(elem,hybrid,nets,nete,dtnu,ptop,pmid,&
end if
end subroutine print_cfl
+ !=============================================================================
+ ! se_stability_constants: single source for the stability constants used by
+ ! print_cfl and auto_rsplit.
+ !=============================================================================
+ pure subroutine se_stability_constants(ptop, lambda_max, lambda_vis, s_rk, umax)
+ use dimensions_mod, only: np
+ use control_mod, only: tstep_type
+
+ real(kind=r8), intent(in) :: ptop ! model top pressure [Pa]
+ real(kind=r8), intent(out) :: lambda_max ! max eigenvalue of gradient operator
+ real(kind=r8), intent(out) :: lambda_vis ! max eigenvalue of Laplace operator
+ real(kind=r8), intent(out) :: s_rk ! stability region of dynamics time-stepping
+ real(kind=r8), intent(out) :: umax ! max wind estimate [m/s] given model top
+
+ ! Eigenvalues calculated by folks at UMich (Paul U & Jared W)
+ select case (np)
+ case (2); lambda_max = 0.5_r8 ; lambda_vis = 0.0_r8 ! lambda_vis not computed
+ case (3); lambda_max = 1.5_r8 ; lambda_vis = 12.0_r8
+ case (4); lambda_max = 2.74_r8; lambda_vis = 30.0_r8
+ case (5); lambda_max = 4.18_r8; lambda_vis = 91.6742_r8
+ case (6); lambda_max = 5.86_r8; lambda_vis = 190.1176_r8
+ case (7); lambda_max = 7.79_r8; lambda_vis = 374.7788_r8
+ case (8); lambda_max = 10.0_r8; lambda_vis = 652.3015_r8
+ case DEFAULT; lambda_max = 0.0_r8 ; lambda_vis = 0.0_r8
+ end select
+ ! RK stability region (tstep_type 1-3: S=2; 4 = Kinnmark&Gray RK3 5-stage: S=3)
+ s_rk = 2.0_r8
+ if (tstep_type==4) s_rk = 3.0_r8
+ ! max wind estimate by model top location (same pressure ranges as the model
+ ! top damping classification in print_cfl)
+ if (ptop>100.0_r8) then
+ umax = 120._r8 ! model top below ~42km
+ else if (ptop>1e-1_r8) then
+ umax = 240._r8 ! CAM7 top (~42-90km)
+ else if (ptop>1e-4_r8) then
+ umax = 300._r8 ! WACCM (~90-140km)
+ else
+ umax = 800._r8 ! WACCM-x (>140km)
+ end if
+ end subroutine se_stability_constants
+
+ !=============================================================================
+ ! set_global_max_normDinv: form the global max of elem(:)%normDinv from the
+ ! task-local max (0 on tasks without elements). Called once from
+ ! dyn_grid_init, on all tasks, BEFORE the dynamics timestep is initialized.
+ !=============================================================================
+ subroutine set_global_max_normDinv(max_normDinv_loc, comm)
+ use spmd_utils, only: mpi_real8, mpi_max
+
+ real(kind=r8), intent(in) :: max_normDinv_loc ! task-local max of elem(:)%normDinv
+ integer, intent(in) :: comm ! communicator spanning all tasks
+ integer :: ierr
+
+ call MPI_Allreduce(max_normDinv_loc, global_max_normDinv, 1, mpi_real8, mpi_max, comm, ierr)
+ end subroutine set_global_max_normDinv
+
+ !=============================================================================
+ ! auto_rsplit: smallest rsplit keeping the dynamics timestep within the
+ ! advective/gravity-wave stability limit (as diagnosed in print_cfl). Called
+ ! from dyn_grid_init BEFORE tstep and the SE nstep counter are initialized, so
+ ! restart reading (TimeLevel_Qdp), the variable-resolution viscosity limiter
+ ! and all CFL diagnostics see the final timestep.
+ !=============================================================================
+ function auto_rsplit(dt_remap, ptop) result(rsplit_auto)
+ use physconst, only: ra
+ use cam_abortutils, only: endrun
+
+ real(kind=r8), intent(in) :: dt_remap ! vertical remap timestep (dtime/nsplit)
+ real(kind=r8), intent(in) :: ptop ! model top pressure [Pa]
+ integer :: rsplit_auto
+
+ real(kind=r8) :: lambda_max, lambda_vis, s_rk, umax
+
+ if (global_max_normDinv < 0.0_r8) call endrun('auto_rsplit: set_global_max_normDinv not called')
+ call se_stability_constants(ptop, lambda_max, lambda_vis, s_rk, umax)
+ if (lambda_max == 0.0_r8) call endrun('auto_rsplit (se_rsplit=-1): lambda_max not available for this np')
+ ! smallest rsplit with dt_remap/rsplit <= s_rk/(max(umax,ugw)*global_max_normDinv*lambda_max*ra)
+ rsplit_auto = max(1, ceiling(dt_remap*max(umax,ugw)*global_max_normDinv*lambda_max*ra/s_rk))
+ end function auto_rsplit
+
!
! ============================
! global_maximum:
diff --git a/src/dynamics/se/dycore/namelist_mod.F90 b/src/dynamics/se/dycore/namelist_mod.F90
index 7e375a097e..f28f27f6e1 100644
--- a/src/dynamics/se/dycore/namelist_mod.F90
+++ b/src/dynamics/se/dycore/namelist_mod.F90
@@ -11,8 +11,6 @@ module namelist_mod
partmethod, & ! Mesh partitioning method (METIS)
multilevel, &
numnodes, &
- tasknum, & ! used dg model in AIX machine
- remapfreq, & ! number of steps per remapping call
statefreq, & ! number of steps per printstate call
runtype, &
cubed_sphere_map, &
@@ -60,8 +58,6 @@ subroutine homme_set_defaults()
numnodes = -1
runtype = 0
statefreq = 1
- remapfreq = 240
- tasknum =-1
columnpackage = "none"
nu_top = 0
ne = 0
diff --git a/src/dynamics/se/dycore/prim_advance_mod.F90 b/src/dynamics/se/dycore/prim_advance_mod.F90
index 25323df956..f11f3333f9 100644
--- a/src/dynamics/se/dycore/prim_advance_mod.F90
+++ b/src/dynamics/se/dycore/prim_advance_mod.F90
@@ -108,7 +108,7 @@ subroutine prim_advance_exp(elem, fvm, deriv, hvcoord, hybrid,dt, tl, nets, net
call omp_set_nested(.true.)
! default weights for computing mean dynamics fluxes
- eta_ave_w = 1_r8/qsplit
+ eta_ave_w = 1.0_r8/qsplit
! ==================================
! Take timestep
diff --git a/src/dynamics/se/dycore/prim_driver_mod.F90 b/src/dynamics/se/dycore/prim_driver_mod.F90
index e2d470f616..521c63e747 100644
--- a/src/dynamics/se/dycore/prim_driver_mod.F90
+++ b/src/dynamics/se/dycore/prim_driver_mod.F90
@@ -28,11 +28,10 @@ subroutine prim_init2(elem, fvm, hybrid, nets, nete, tl, hvcoord)
use dimensions_mod, only: fv_nphys, nc
use parallel_mod, only: syncmp
use se_dyn_time_mod, only: timelevel_t, tstep, phys_tscale, nsplit, TimeLevel_Qdp
- use se_dyn_time_mod, only: nsplit_baseline,rsplit_baseline
use prim_state_mod, only: prim_printstate
- use control_mod, only: runtype, topology, rsplit, qsplit, rk_stage_user, &
+ use control_mod, only: runtype, topology, rsplit, qsplit, &
nu, nu_q, nu_div, hypervis_subcycle, hypervis_subcycle_q, &
- hypervis_subcycle_sponge, variable_nsplit
+ hypervis_subcycle_sponge
use fvm_mod, only: fill_halo_fvm,ghostBufQnhc_h
use thread_mod, only: omp_get_thread_num
use global_norms_mod, only: print_cfl
@@ -58,11 +57,6 @@ subroutine prim_init2(elem, fvm, hybrid, nets, nete, tl, hvcoord)
! Local variables
! ==================================
-! variables used to calculate CFL
- real (kind=r8) :: dtnu ! timestep*viscosity parameter
- real (kind=r8) :: dt_dyn_del2_sponge
- real (kind=r8) :: dt_tracer_vis ! viscosity timestep used in tracers
- real (kind=r8) :: dt_dyn_vis ! viscosity timestep
real (kind=r8) :: dt_remap ! remapping timestep
real (kind=r8) :: dp,dp0,T1,T0,pmid_ref(np,np)
@@ -83,23 +77,6 @@ subroutine prim_init2(elem, fvm, hybrid, nets, nete, tl, hvcoord)
! ==========================
!call prim_advance_init(hybrid%par,elem)
- ! compute most restrictive dt*nu for use by variable res viscosity:
- ! compute timestep seen by viscosity operator:
- dt_dyn_vis = tstep
- dt_dyn_del2_sponge = tstep
- dt_tracer_vis=tstep*qsplit
- dt_remap=dt_tracer_vis*rsplit
- ! compute most restrictive condition:
- ! note: dtnu ignores subcycling
- dtnu=max(dt_dyn_vis*max(nu,nu_div), dt_tracer_vis*nu_q)
- ! compute actual viscosity timesteps with subcycling
- dt_tracer_vis = dt_tracer_vis/hypervis_subcycle_q
- dt_dyn_vis = dt_dyn_vis/hypervis_subcycle
- dt_dyn_del2_sponge = dt_dyn_del2_sponge/hypervis_subcycle_sponge
- if (variable_nsplit) then
- nsplit_baseline=nsplit
- rsplit_baseline=rsplit
- end if
! ==================================
! Initialize derivative structure
! ==================================
@@ -122,13 +99,10 @@ subroutine prim_init2(elem, fvm, hybrid, nets, nete, tl, hvcoord)
! CAM has set tstep based on dtime before calling prim_init2(),
! so only now does HOMME learn the timstep. print them out:
- call print_cfl(elem,hybrid,nets,nete,dtnu,&
+ call print_cfl(elem,hybrid,nets,nete,&
!p top and p mid levels
- hvcoord%hyai(1)*hvcoord%ps0,hvcoord%hyam(:)*hvcoord%ps0+hvcoord%hybm(:)*pstd,&
- !dt_remap,dt_tracer_fvm,dt_tracer_se
- tstep*qsplit*rsplit,tstep*qsplit*fvm_supercycling,tstep*qsplit,&
- !dt_dyn,dt_dyn_visco,dt_tracer_visco, dt_phys
- tstep,dt_dyn_vis,dt_dyn_del2_sponge,dt_tracer_vis,tstep*nsplit*qsplit*rsplit)
+ hvcoord%hyai(1)*hvcoord%ps0,hvcoord%hyam(:)*hvcoord%ps0+hvcoord%hybm(:)*pstd)
+ dt_remap = tstep*qsplit*rsplit
if (hybrid%masterthread) then
if (phys_tscale/=0) then
@@ -219,10 +193,10 @@ subroutine prim_run_subcycle(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord,nsubst
!
use hybvcoord_mod, only : hvcoord_t
use se_dyn_time_mod, only: TimeLevel_t, timelevel_update, timelevel_qdp, nsplit
- use control_mod, only: statefreq,qsplit, rsplit, variable_nsplit, dribble_in_rsplit_loop
+ use control_mod, only: statefreq,qsplit, rsplit, dribble_in_rsplit_loop
use prim_advance_mod, only: applycamforcing
use prim_advance_mod, only: tot_energy_dyn,compute_omega
- use prim_state_mod, only: prim_printstate, adjust_nsplit
+ use prim_state_mod, only: prim_printstate
use prim_advection_mod, only: vertical_remap, deriv
use thread_mod, only: omp_get_thread_num
use perf_mod , only: t_startf, t_stopf
@@ -266,7 +240,7 @@ subroutine prim_run_subcycle(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord,nsubst
!
! initialize variables for computing vertical Courant number
!
- if (variable_nsplit.or.compute_diagnostics) then
+ if (compute_diagnostics) then
if (nsubstep==1) then
do ie=nets,nete
omega_cn(1,ie) = 0.0_r8
@@ -329,7 +303,7 @@ subroutine prim_run_subcycle(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord,nsubst
call tot_energy_dyn(elem,fvm,nets,nete,tl%np1,np1_qdp,'dAD')
- if (variable_nsplit.or.compute_diagnostics) then
+ if (compute_diagnostics) then
!
! initialize variables for computing vertical Courant number
!
@@ -374,7 +348,7 @@ subroutine prim_run_subcycle(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord,nsubst
!
! Compute vertical Courant numbers
!
- if (variable_nsplit.or.compute_diagnostics) then
+ if (compute_diagnostics) then
do ie=nets,nete
do k=1,nlev
do j=1,np
@@ -390,11 +364,6 @@ subroutine prim_run_subcycle(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord,nsubst
end do
end do
end do
- if (nsubstep==nsplit.and.variable_nsplit) then
- call t_startf('adjust_nsplit')
- call adjust_nsplit(elem, tl, hybrid,nets,nete, fvm, omega_cn)
- call t_stopf('adjust_nsplit')
- end if
end if
! ============================================================
@@ -440,10 +409,8 @@ subroutine prim_step(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord, rstep, last_s
use prim_advection_mod, only: prim_advec_tracers_remap, prim_advec_tracers_fvm, deriv
use derivative_mod, only: subcell_integration
use hybrid_mod, only: set_region_num_threads, config_thread_region, get_loop_ranges
- use dimensions_mod, only: use_cslam,fvm_supercycling,fvm_supercycling_jet
- use dimensions_mod, only: kmin_jet, kmax_jet
+ use dimensions_mod, only: use_cslam,fvm_supercycling
use fvm_mod, only: ghostBufQnhc_vh,ghostBufQ1_vh, ghostBufFlux_vh
- use fvm_mod, only: ghostBufQ1_h,ghostBufQnhcJet_h, ghostBufFluxJet_h
use se_dyn_time_mod, only: timelevel_qdp
use fvm_mapping, only: cslam2gll
#ifdef waccm_debug
@@ -582,7 +549,7 @@ subroutine prim_step(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord, rstep, last_s
!
! FVM transport
!
- if ((mod(rstep,fvm_supercycling) == 0).and.(mod(rstep,fvm_supercycling_jet) == 0)) then
+ if (mod(rstep,fvm_supercycling) == 0) then
! call omp_set_nested(.true.)
! !$OMP PARALLEL NUM_THREADS(vert_num_threads), DEFAULT(SHARED), PRIVATE(hybridnew2,kbeg,kend)
@@ -613,12 +580,6 @@ subroutine prim_step(elem, fvm, hybrid,nets,nete, dt, tl, hvcoord, rstep, last_s
end do
call TimeLevel_Qdp( tl, qsplit, n0_qdp, np1_qdp)
if (.not.last_step) call cslam2gll(elem, fvm, hybrid,nets,nete, tl%np1, np1_qdp)
- else if ((mod(rstep,fvm_supercycling_jet) == 0)) then
- !
- ! shorter fvm time-step in jet region
- !
- call Prim_Advec_Tracers_fvm(elem,fvm,hvcoord,hybrid,&
- dt_q,tl,nets,nete,ghostBufQnhcJet_h,ghostBufQ1_h, ghostBufFluxJet_h,kmin_jet,kmax_jet)
end if
#ifdef waccm_debug
diff --git a/src/dynamics/se/dycore/prim_state_mod.F90 b/src/dynamics/se/dycore/prim_state_mod.F90
index 0877834f97..0ad2bfa4ac 100644
--- a/src/dynamics/se/dycore/prim_state_mod.F90
+++ b/src/dynamics/se/dycore/prim_state_mod.F90
@@ -14,7 +14,7 @@ module prim_state_mod
implicit none
private
- public :: prim_printstate, adjust_nsplit
+ public :: prim_printstate
CONTAINS
@@ -345,91 +345,6 @@ subroutine prim_printstate_cslam_gamma(elem, tl,hybrid,nets,nete, fvm)
end subroutine prim_printstate_cslam_gamma
#endif
- subroutine adjust_nsplit(elem, tl,hybrid,nets,nete, fvm, omega_cn)
- use dimensions_mod, only: ksponge_end
- use dimensions_mod, only: fvm_supercycling, fvm_supercycling_jet
- use se_dyn_time_mod, only: tstep
- use control_mod, only: rsplit, qsplit
- use perf_mod, only: t_startf, t_stopf
- use se_dyn_time_mod, only: nsplit, nsplit_baseline,rsplit_baseline
- use control_mod, only: qsplit, rsplit
- use time_manager, only: get_step_size
- use cam_abortutils, only: endrun
- use control_mod, only: nu_top
- !
- type (element_t), intent(inout) :: elem(:)
- type (TimeLevel_t), target, intent(in) :: tl
- type (hybrid_t), intent(in) :: hybrid
- integer, intent(in) :: nets,nete
- type(fvm_struct), intent(inout) :: fvm(:)
- real (kind=r8), intent(in) :: omega_cn(2,nets:nete)
- ! Local variables...
- integer :: k,ie
- real (kind=r8), dimension(1) :: min_o
- real (kind=r8), dimension(1) :: max_o
- real (kind=r8) :: dtime
- character(len=128) :: errmsg
- real (kind=r8) :: threshold=0.90_r8
- real (kind=r8) :: max_abs_omega_cn(nets:nete)
- real (kind=r8) :: min_abs_omega_cn(nets:nete)
- !
- ! The threshold values for when to double nsplit are empirical.
- ! In FW2000climo runs the Courant numbers are large in the sponge
- !
- ! The model was found to be stable if regular del4 is increased
- ! in the sponge and nu_top is increased (when nsplit doubles)
- !
- !
- do ie=nets,nete
- max_abs_omega_cn(ie) = MAXVAL(ABS(omega_cn(:,ie)))
- end do
-
- !JMD This is a Thread Safe Reduction
- do k = 1,1
- max_o(k) = ParallelMax(max_abs_omega_cn(:),hybrid)
-! min_o(k) = ParallelMin(min_abs_omega_cn(:),hybrid)
- end do
- if (max_o(1)>threshold.and.nsplit==nsplit_baseline) then
- !
- ! change vertical remap time-step
- !
- nsplit=2*nsplit_baseline
- fvm_supercycling = rsplit
- fvm_supercycling_jet = rsplit
- nu_top=2.0_r8*nu_top
- !
- ! write diagnostics to log file
- !
- if(hybrid%masterthread) then
- !dynamics variables in n0 are at time = 'time': time=tl%nstep*tstep
- !dt=tstep*qsplit
- ! dt_remap = tstep*qsplit*rsplit ! vertical REMAP timestep
- !
- write(iulog,*) 'adj. nsplit: doubling nsplit; t=',Time_at(tl%nstep)/(24*3600)," [day]; max OMEGA",max_o(1)
- end if
- dtime = get_step_size()
- tstep = dtime / real(nsplit*qsplit*rsplit, r8)
-
- else if (nsplit.ne.nsplit_baseline.and.max_o(1)<0.4_r8*threshold) then
- !
- ! should nsplit be reduced again?
- !
- nsplit=nsplit_baseline
- rsplit=rsplit_baseline
- fvm_supercycling = rsplit
- fvm_supercycling_jet = rsplit
- nu_top=nu_top/2.0_r8
-
-! nu_div_scale_top(:) = 1.0_r8
-
- dtime = get_step_size()
- tstep = dtime / real(nsplit*qsplit*rsplit, r8)
- if(hybrid%masterthread) then
- write(iulog,*) 'adj. nsplit: reset nsplit ; t=',Time_at(tl%nstep)/(24*3600)," [day]; max OMEGA",max_o(1)
- end if
- end if
- end subroutine adjust_nsplit
-
subroutine prim_printstate_U(elem, tl,hybrid,nets,nete, fvm)
type (element_t), intent(inout) :: elem(:)
type(fvm_struct), intent(inout) :: fvm(:)
diff --git a/src/dynamics/se/dycore/se_dyn_time_mod.F90 b/src/dynamics/se/dycore/se_dyn_time_mod.F90
index cfe7ad2323..bc5091f3ec 100644
--- a/src/dynamics/se/dycore/se_dyn_time_mod.F90
+++ b/src/dynamics/se/dycore/se_dyn_time_mod.F90
@@ -4,8 +4,6 @@ module se_dyn_time_mod
!------------------
implicit none
integer,public :: nsplit=1
- integer,public :: nsplit_baseline=-1
- integer,public :: rsplit_baseline=-1
integer,public :: nmax ! Max number of timesteps
integer,public :: nEndStep ! Number of End Step
integer,public :: ndays ! Max number of days
diff --git a/src/dynamics/se/dyn_comp.F90 b/src/dynamics/se/dyn_comp.F90
index aef838adc2..7bb0e27d0e 100644
--- a/src/dynamics/se/dyn_comp.F90
+++ b/src/dynamics/se/dyn_comp.F90
@@ -110,9 +110,9 @@ subroutine dyn_readnl(NLFileName)
use control_mod, only: hypervis_subcycle_q, statefreq, runtype
use control_mod, only: nu, nu_div, nu_p, nu_q, nu_top, qsplit, rsplit
use control_mod, only: vert_remap_uvTq_alg, vert_remap_tracer_alg
- use control_mod, only: tstep_type, rk_stage_user
+ use control_mod, only: tstep_type
use control_mod, only: ftype, limiter_option, partmethod
- use control_mod, only: topology, variable_nsplit
+ use control_mod, only: topology
use control_mod, only: fine_ne, hypervis_power, hypervis_scaling
use control_mod, only: max_hypervis_courant, statediag_numtrac,refined_mesh
use control_mod, only: molecular_diff, pgf_formulation, dribble_in_rsplit_loop
@@ -120,8 +120,7 @@ subroutine dyn_readnl(NLFileName)
use control_mod, only: min_temperature
use dimensions_mod, only: ne, npart
use dimensions_mod, only: large_Courant_incr
- use dimensions_mod, only: fvm_supercycling, fvm_supercycling_jet
- use dimensions_mod, only: kmin_jet, kmax_jet
+ use dimensions_mod, only: fvm_supercycling
use params_mod, only: SFCURVE
use parallel_mod, only: initmpi
use thread_mod, only: initomp, max_num_threads
@@ -169,9 +168,6 @@ subroutine dyn_readnl(NLFileName)
logical :: se_write_restart_unstruct
logical :: se_large_Courant_incr
integer :: se_fvm_supercycling
- integer :: se_fvm_supercycling_jet
- integer :: se_kmin_jet
- integer :: se_kmax_jet
real(r8) :: se_molecular_diff
integer :: se_pgf_formulation
integer :: se_dribble_in_rsplit_loop
@@ -217,9 +213,6 @@ subroutine dyn_readnl(NLFileName)
se_write_restart_unstruct, &
se_large_Courant_incr, &
se_fvm_supercycling, &
- se_fvm_supercycling_jet, &
- se_kmin_jet, &
- se_kmax_jet, &
se_molecular_diff, &
se_pgf_formulation, &
se_dribble_in_rsplit_loop, &
@@ -292,9 +285,6 @@ subroutine dyn_readnl(NLFileName)
call MPI_bcast(se_write_restart_unstruct, 1, mpi_logical, masterprocid, mpicom, ierr)
call MPI_bcast(se_large_Courant_incr, 1, mpi_logical, masterprocid, mpicom, ierr)
call MPI_bcast(se_fvm_supercycling, 1, mpi_integer, masterprocid, mpicom, ierr)
- call MPI_bcast(se_fvm_supercycling_jet, 1, mpi_integer, masterprocid, mpicom, ierr)
- call MPI_bcast(se_kmin_jet, 1, mpi_integer, masterprocid, mpicom, ierr)
- call MPI_bcast(se_kmax_jet, 1, mpi_integer, masterprocid, mpicom, ierr)
call MPI_bcast(se_molecular_diff, 1, mpi_real8, masterprocid, mpicom, ierr)
call MPI_bcast(se_pgf_formulation, 1, mpi_integer, masterprocid, mpicom, ierr)
call MPI_bcast(se_dribble_in_rsplit_loop, 1, mpi_integer, masterprocid, mpicom, ierr)
@@ -309,9 +299,6 @@ subroutine dyn_readnl(NLFileName)
call initomp()
- if (se_fvm_supercycling < 0) se_fvm_supercycling = se_rsplit
- if (se_fvm_supercycling_jet < 0) se_fvm_supercycling_jet = se_rsplit
-
! Go ahead and enforce ne = 0 for refined mesh runs
if (se_refined_mesh) then
se_ne = 0
@@ -323,8 +310,6 @@ subroutine dyn_readnl(NLFileName)
partmethod = SFCURVE
npart = se_npes
! CAM requires forward-in-time, subcycled dynamics
- ! RK2 3 stage tracers, sign-preserving conservative
- rk_stage_user = 3
topology = "cube"
! Finally, set the HOMME variables which have different names
fine_ne = se_fine_ne
@@ -361,10 +346,6 @@ subroutine dyn_readnl(NLFileName)
fv_nphys = se_fv_nphys
large_Courant_incr = se_large_Courant_incr
fvm_supercycling = se_fvm_supercycling
- fvm_supercycling_jet = se_fvm_supercycling_jet
- kmin_jet = se_kmin_jet
- kmax_jet = se_kmax_jet
- variable_nsplit = .false.
molecular_diff = se_molecular_diff
pgf_formulation = se_pgf_formulation
dribble_in_rsplit_loop = se_dribble_in_rsplit_loop
@@ -384,8 +365,10 @@ subroutine dyn_readnl(NLFileName)
use_cslam = .false.
end if
- if (rsplit < 1) then
- call endrun('dyn_readnl: rsplit must be > 0')
+ ! rsplit == -1 is the "automatic" sentinel, resolved later in dyn_grid_init
+ ! (via auto_rsplit); only reject other non-positive values here.
+ if (rsplit < 1 .and. rsplit /= -1) then
+ call endrun('dyn_readnl: rsplit must be > 0 (or -1 for automatic subcycling)')
end if
! if restart or branch run
@@ -424,9 +407,6 @@ subroutine dyn_readnl(NLFileName)
write_restart_unstruct = se_write_restart_unstruct
- if (se_kmin_jet<0 ) kmin_jet = 1
- if (se_kmax_jet<0 ) kmax_jet = nlev
-
if (masterproc) then
write(iulog, '(a,i0)') 'dyn_readnl: se_ftype = ',ftype
write(iulog, '(a,i0)') 'dyn_readnl: se_statediag_numtrac = ',statediag_numtrac
@@ -460,9 +440,6 @@ subroutine dyn_readnl(NLFileName)
write(iulog, '(a,a)') 'dyn_readnl: se_vert_remap_uvTq_alg = ',trim(se_vert_remap_uvTq_alg)
write(iulog, '(a,a)') 'dyn_readnl: se_vert_remap_tracer_alg = ',trim(se_vert_remap_tracer_alg)
write(iulog, '(a,i0)') 'dyn_readnl: se_fvm_supercycling = ',fvm_supercycling
- write(iulog, '(a,i0)') 'dyn_readnl: se_fvm_supercycling_jet = ',fvm_supercycling_jet
- write(iulog, '(a,i0)') 'dyn_readnl: se_kmin_jet = ',kmin_jet
- write(iulog, '(a,i0)') 'dyn_readnl: se_kmax_jet = ',kmax_jet
write(iulog, *) 'dyn_readnl: se_sponge_del4_nu_fac = ',se_sponge_del4_nu_fac
if (se_sponge_del4_nu_fac < 0) write(iulog, '(a)') ' (automatically set based on model top location)'
diff --git a/src/dynamics/se/dyn_grid.F90 b/src/dynamics/se/dyn_grid.F90
index 69d9bbc520..e8ee0722ed 100644
--- a/src/dynamics/se/dyn_grid.F90
+++ b/src/dynamics/se/dyn_grid.F90
@@ -135,6 +135,8 @@ subroutine dyn_grid_init()
get_loop_ranges, config_thread_region
use control_mod, only: qsplit, rsplit
use se_dyn_time_mod, only: tstep, nsplit
+ use dimensions_mod, only: fvm_supercycling
+ use global_norms_mod, only: auto_rsplit, set_global_max_normDinv
use fvm_mod, only: fvm_init2, fvm_init3, fvm_pg_init
use dimensions_mod, only: irecons_tracer
use comp_gll_ctr_vol, only: gll_grid_write
@@ -150,6 +152,8 @@ subroutine dyn_grid_init()
type(hybrid_t) :: hybrid
integer :: ierr
integer :: dtime
+ integer :: ie
+ real(r8) :: max_normDinv_loc ! task-local max of the mesh metric
real(r8), allocatable ::clat(:), clon(:), areaa(:)
integer :: nets, nete
@@ -209,11 +213,30 @@ subroutine dyn_grid_init()
! nelemd (# of elements on this task) is set by prim_init1
call init_loop_ranges(nelemd)
+ ! Global max of the mesh metric normDinv, used for CFL-based "automatic"
+ ! (se_*=-1) subcycling below and reused by print_cfl (the reduction moved
+ ! here from print_cfl; none is added)
+ max_normDinv_loc = 0.0_r8
+ do ie = 1, nelemd
+ max_normDinv_loc = max(max_normDinv_loc, elem(ie)%normDinv)
+ end do
+ call set_global_max_normDinv(max_normDinv_loc, mpicom)
+
! Dynamics timestep
!
! Note: dtime = timestep for physics/dynamics coupling
! tstep = the dynamics timestep:
dtime = get_step_size()
+ if (rsplit == -1) then
+ ! "automatic" rsplit: smallest value keeping the dynamics timestep within
+ ! the advective/gravity-wave stability limit. Resolved here, before tstep
+ ! and TimeLevel%nstep are set, so restart reading and all downstream
+ ! initialization see the final timestep.
+ rsplit = auto_rsplit(dtime/real(nsplit,r8), hvcoord%hyai(1)*hvcoord%ps0)
+ if (masterproc) write(iulog,'(a,i0)') 'dyn_grid_init: automatic (se_rsplit=-1) rsplit = ', rsplit
+ end if
+ ! default (se_fvm_supercycling < 0): inherit rsplit
+ if (fvm_supercycling < 0) fvm_supercycling = rsplit
tstep = dtime / real(nsplit*qsplit*rsplit, r8)
TimeLevel%nstep = get_nstep()*nsplit*qsplit*rsplit