From a3c7f4931dd41f904ac42b6e1ae33fb7b24c7294 Mon Sep 17 00:00:00 2001 From: Bill Skamarock Date: Mon, 20 Oct 2025 11:22:50 -0600 Subject: [PATCH] This commit provides access to the hybrid vertical coordinate and the transition height of the coordinate through namelist variables for the init_atmosphere core. Additionally, the smoothing coefficient formula for the hybrid coordinate now uses the transition height instead of the model top height. These changes affect only the real-data construction of the vertical mesh. The changes include two new namelist configurations variables in the init_atmosphere Registry file: a logical for activating the hybrid vertical coordinate and the transition height for the coordinate (where it transitions from terrain-following to constant height). These variables are config_hybrid_coordinate (logical) and config_hybrid_top_z (real, meters). The variables had been hardwired in the code up until now (true and 30 km). Also included are a few changes to provide further information concerning diagnostic output appearing in the init_atmosphere log files. --- src/core_init_atmosphere/Registry.xml | 10 +++++ .../mpas_init_atm_cases.F | 38 +++++++++++++------ 2 files changed, 36 insertions(+), 12 deletions(-) diff --git a/src/core_init_atmosphere/Registry.xml b/src/core_init_atmosphere/Registry.xml index cf4934a81b..fcb31aaab4 100644 --- a/src/core_init_atmosphere/Registry.xml +++ b/src/core_init_atmosphere/Registry.xml @@ -260,6 +260,16 @@ description="Whether to blend terrain along domain boundaries with first-guess terrain. Only useful for limited-area domains." possible_values="true or false"/> + + + + diff --git a/src/core_init_atmosphere/mpas_init_atm_cases.F b/src/core_init_atmosphere/mpas_init_atm_cases.F index 673ebfc525..aa400ee481 100644 --- a/src/core_init_atmosphere/mpas_init_atm_cases.F +++ b/src/core_init_atmosphere/mpas_init_atm_cases.F @@ -2761,6 +2761,8 @@ subroutine init_atm_case_gfs(block, mesh, nCells, nEdges, nVertLevels, fg, state integer, pointer :: config_nsm real (kind=RKIND), pointer :: config_dzmin real (kind=RKIND), pointer :: config_ztop + logical , pointer :: config_hybrid_coordinate + real (kind=RKIND), pointer :: config_hybrid_top_z logical, pointer :: config_tc_vertical_grid character (len=StrKIND), pointer :: config_specified_zeta_levels logical, pointer :: config_use_spechumd @@ -2837,6 +2839,8 @@ subroutine init_atm_case_gfs(block, mesh, nCells, nEdges, nVertLevels, fg, state call mpas_pool_get_config(configs, 'config_nsm', config_nsm) call mpas_pool_get_config(configs, 'config_dzmin', config_dzmin) call mpas_pool_get_config(configs, 'config_ztop', config_ztop) + call mpas_pool_get_config(configs, 'config_hybrid_coordinate' , config_hybrid_coordinate) + call mpas_pool_get_config(configs, 'config_hybrid_top_z' , config_hybrid_top_z) call mpas_pool_get_config(configs, 'config_tc_vertical_grid', config_tc_vertical_grid) call mpas_pool_get_config(configs, 'config_specified_zeta_levels', config_specified_zeta_levels) call mpas_pool_get_config(configs, 'config_use_spechumd', config_use_spechumd) @@ -3155,38 +3159,36 @@ subroutine init_atm_case_gfs(block, mesh, nCells, nEdges, nVertLevels, fg, state ! ah(k) = 1.-zw(k)/zt is the basic terrain-following coordinate ! ah(k) = 0 is a height coordinate - hybrid = .true. -! hybrid = .false. - kz = nz - if (hybrid) then + if (config_hybrid_coordinate) then - zh = 30000.0 -! zh = 0.5*zt - + zh = config_hybrid_top_z do k=1,nz if (zw(k) < zh) then ah(k) = cos(.5*pii*zw(k)/zh)**6 - -!!! ah(k) = ah(k)*(1.-zw(k)/zt) - else ah(k) = 0. kz = min(kz,k) end if end do + call mpas_log_write('Hybrid vertical coordinate is used. Transition height, ztop are $r $r ',realArgs=(/config_hybrid_top_z,zt/)) else + zh = zt do k=1,nz ah(k) = 1.-zw(k)/zt end do + call mpas_log_write('Hybrid vertical coordinate is not selected. ztop is ',realArgs=(/zt/)) end if + call mpas_log_write(' ') + call mpas_log_write(' interface k, height, ah value ') do k=1,nz call mpas_log_write('$i $r $r', intArgs=(/k/), realArgs=(/zw(k), ah(k)/)) end do + call mpas_log_write(' ') do k=1,nz1 dzw (k) = zw(k+1)-zw(k) @@ -3237,6 +3239,12 @@ subroutine init_atm_case_gfs(block, mesh, nCells, nEdges, nVertLevels, fg, state if (smooth) then + call mpas_log_write(' ') + call mpas_log_write(' surfaces will be smoothed, dzmin = $r ',realArgs=(/config_dzmin/)) + call mpas_log_write(' base number of smoothing passes = $i ',intArgs=(/config_nsm/)) + call mpas_log_write(' ') + call mpas_log_write(' level, smoothing steps, smoothing factor, smallest fractional dz ') + dzmin = config_dzmin allocate(sm0(nCells+1)) @@ -3257,7 +3265,7 @@ subroutine init_atm_case_gfs(block, mesh, nCells, nEdges, nVertLevels, fg, state do i=1,config_nsm + k do iCell=1,nCells - sm = sm0(iCell) * min((3.0*zw(k)/zt)**2.0, 1.0_RKIND) + sm = sm0(iCell) * min((3.0*zw(k)/zh)**2.0, 1.0_RKIND) hs1(iCell) = 0. do j = 1,nEdgesOnCell(iCell) @@ -3302,7 +3310,8 @@ subroutine init_atm_case_gfs(block, mesh, nCells, nEdges, nVertLevels, fg, state call mpas_dmpar_min_real(dminfo, dzminf, dzminf_global) call mpas_log_write('$i $i $r $r', intArgs=(/k,i/), realArgs=(/sm,dzminf_global/(zw(k)-zw(k-1))/)) end do - + call mpas_log_write(' ') + deallocate(sm0) do k=kz,nz @@ -3311,10 +3320,15 @@ subroutine init_atm_case_gfs(block, mesh, nCells, nEdges, nVertLevels, fg, state else + call mpas_log_write(' ') + call mpas_log_write(' surfaces are not smoothed ') + call mpas_log_write(' level, smallest fractional dz ') + do k=2,nz1 dzmina = minval(zw(k)+ah(k)*hx(k,:)-zw(k-1)-ah(k-1)*hx(k-1,:)) call mpas_log_write('$i $r', intArgs=(/k/), realArgs=(/dzmina/(zw(k)-zw(k-1))/)) end do + call mpas_log_write(' ') end if