From 1ade5cab54362f94f64136ecc03553e46c731b04 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Tue, 21 Oct 2025 11:43:43 -0600 Subject: [PATCH 01/22] Do not hard code path to working directory in GitHub Actions The workflow will fail when a fork uses a different repository name. --- .github/workflows/unit-tests.yaml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/unit-tests.yaml b/.github/workflows/unit-tests.yaml index 45ef3a57..bb54666f 100644 --- a/.github/workflows/unit-tests.yaml +++ b/.github/workflows/unit-tests.yaml @@ -35,7 +35,7 @@ jobs: - name: Build atmospheric_physics run: | cmake \ - -DCMAKE_PREFIX_PATH=/home/runner/work/atmospheric_physics/atmospheric_physics/pFUnit/build/installed \ + -DCMAKE_PREFIX_PATH=$GITHUB_WORKSPACE/pFUnit/build/installed \ -DATMOSPHERIC_PHYSICS_ENABLE_CODE_COVERAGE=ON \ -B./build \ -S./test/unit-test From d289167636e78c0e05ce102904b89b76dad30df3 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Mon, 8 Dec 2025 13:41:34 -0700 Subject: [PATCH 02/22] Make GitHub Linguist recognize `*.pf` files as Fortran source code --- .gitattributes | 1 + 1 file changed, 1 insertion(+) create mode 100644 .gitattributes diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 00000000..846e6dbf --- /dev/null +++ b/.gitattributes @@ -0,0 +1 @@ +*.pf linguist-language=Fortran-Free-Form From 3c70286c69f67d42f3e0820cedaaa96a8c2d41a9 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Mon, 8 Dec 2025 13:59:24 -0700 Subject: [PATCH 03/22] Synchronize `CMakeLists.txt` of MMM physics with MPAS dycore in CAM-SIMA --- schemes/mmm/CMakeLists.txt | 14 ++++++++------ test/unit-test/tests/mmm/CMakeLists.txt | 6 ++++++ 2 files changed, 14 insertions(+), 6 deletions(-) diff --git a/schemes/mmm/CMakeLists.txt b/schemes/mmm/CMakeLists.txt index e52c1192..82f278de 100644 --- a/schemes/mmm/CMakeLists.txt +++ b/schemes/mmm/CMakeLists.txt @@ -1,8 +1,8 @@ -cmake_minimum_required(VERSION 3.20) +cmake_minimum_required(VERSION 3.21) # `mmm_physics_compat` has not been integrated into the CMake build of any top level projects yet, # and this CMakeLists.txt file is currently for unit testing purposes only. -# Making a change to this CMakeLists.txt file will not impact the build of a parent project at this time. +# Therefore, making a change to this CMakeLists.txt file will not impact the build of a parent project at this time. project(mmm_physics_compat VERSION 0.1.0 @@ -19,11 +19,13 @@ target_sources(mmm_physics_compat ccpp_kind_types.F90 mmm_physics_compat.F90 ) -target_compile_options(mmm_physics_compat - PRIVATE - $<$,$>:-fbacktrace -fcheck=all -std=f2018 -Wall -Wextra -Wpedantic> -) target_include_directories(mmm_physics_compat INTERFACE ${CMAKE_CURRENT_BINARY_DIR} ) +target_compile_options(mmm_physics_compat + PRIVATE + $<$,$>:-fbacktrace -fcheck=all -ffpe-trap=invalid,overflow,zero -fno-unsafe-math-optimizations -frounding-math -fsignaling-nans -std=f2018 -Wall -Wextra -Wpedantic> + $<$,$>:-check all -fp-model=precise -stand f18 -traceback -warn all> + $<$,$>:-check all -fp-model=precise -stand f18 -traceback -warn all> +) diff --git a/test/unit-test/tests/mmm/CMakeLists.txt b/test/unit-test/tests/mmm/CMakeLists.txt index 69528847..269c7b6e 100644 --- a/test/unit-test/tests/mmm/CMakeLists.txt +++ b/test/unit-test/tests/mmm/CMakeLists.txt @@ -6,3 +6,9 @@ add_pfunit_ctest(mmm_physics_compat_tests LINK_LIBRARIES mmm_physics_compat ) +target_compile_options(mmm_physics_compat_tests + PRIVATE + $<$,$>:-fbacktrace -fcheck=all -ffpe-trap=invalid,overflow,zero -fno-unsafe-math-optimizations -frounding-math -fsignaling-nans -std=f2018 -Wall -Wextra -Wpedantic> + $<$,$>:-check all -fp-model=precise -stand f18 -traceback -warn all> + $<$,$>:-check all -fp-model=precise -stand f18 -traceback -warn all> +) From 017eca22d4f71c7f48215b9443645dacae50e96e Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 5 Dec 2025 15:34:15 -0700 Subject: [PATCH 04/22] Fix incorrect standard name in new Tiedtke convection scheme `hfx` is simply the sensible heat flux, not the sum of latent heat flux and sensible heat flux. --- schemes/mmm/cu_ntiedtke_compat.F90 | 7 ++----- schemes/mmm/cu_ntiedtke_compat.meta | 20 +------------------- 2 files changed, 3 insertions(+), 24 deletions(-) diff --git a/schemes/mmm/cu_ntiedtke_compat.F90 b/schemes/mmm/cu_ntiedtke_compat.F90 index 19717226..bb18a5a1 100644 --- a/schemes/mmm/cu_ntiedtke_compat.F90 +++ b/schemes/mmm/cu_ntiedtke_compat.F90 @@ -14,21 +14,19 @@ module cu_ntiedtke_compat !! \htmlinclude cu_ntiedtke_compat_pre_run.html subroutine cu_ntiedtke_compat_pre_run( & cflx, exner, landfrac, & - lhf, shf, & rthdynten, rthblten, rthratenlw, rthratensw, & rqvdynten, rqvblten, & lndj, & - ptf, pqvf, hfx, evap, & + ptf, pqvf, evap, & errmsg, errflg) use ccpp_kinds, only: kind_phys use ccpp_scheme_utils, only: ccpp_constituent_index real(kind_phys), intent(in) :: cflx(:, :), exner(:, :), landfrac(:), & - lhf(:), shf(:), & rthdynten(:, :), rthblten(:, :), rthratenlw(:, :), rthratensw(:, :), & rqvdynten(:, :), rqvblten(:, :) integer, intent(out) :: lndj(:) - real(kind_phys), intent(out) :: ptf(:, :), pqvf(:, :), hfx(:), evap(:) + real(kind_phys), intent(out) :: ptf(:, :), pqvf(:, :), evap(:) character(*), intent(out) :: errmsg integer, intent(out) :: errflg @@ -42,7 +40,6 @@ subroutine cu_ntiedtke_compat_pre_run( & ptf(:, :) = (rthdynten(:, :) + rthblten(:, :) + rthratenlw(:, :) + rthratensw(:, :)) * exner(:, :) pqvf(:, :) = rqvdynten(:, :) + rqvblten(:, :) - hfx(:) = lhf(:) + shf(:) evap(:) = 0.0_kind_phys call ccpp_constituent_index( & diff --git a/schemes/mmm/cu_ntiedtke_compat.meta b/schemes/mmm/cu_ntiedtke_compat.meta index 6a09d240..5e7d5ed1 100644 --- a/schemes/mmm/cu_ntiedtke_compat.meta +++ b/schemes/mmm/cu_ntiedtke_compat.meta @@ -23,18 +23,6 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = in -[ lhf ] - standard_name = surface_upward_latent_heat_flux_from_coupler - units = W m-2 - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent) - intent = in -[ shf ] - standard_name = surface_upward_sensible_heat_flux_from_coupler - units = W m-2 - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent) - intent = in [ rthdynten ] standard_name = tendency_of_air_potential_temperature_due_to_dynamics units = K s-1 @@ -89,12 +77,6 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = out -[ hfx ] - standard_name = surface_upward_heat_flux - units = W m-2 - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent) - intent = out [ evap ] standard_name = surface_upward_water_vapor_flux units = kg m-2 s-1 @@ -271,7 +253,7 @@ dimensions = (horizontal_loop_extent) intent = in [ hfx ] - standard_name = surface_upward_heat_flux + standard_name = surface_upward_sensible_heat_flux_from_coupler units = W m-2 type = real | kind = kind_phys dimensions = (horizontal_loop_extent) From dd97a2a51aab18f82afa476c30b8072a5aeb8a71 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Mon, 8 Dec 2025 15:11:50 -0700 Subject: [PATCH 05/22] Update git submodule of MMM physics ... for fixing compilation errors with MYNN surface layer scheme. --- .gitmodules | 2 +- schemes/mmm/mmm_physics | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/.gitmodules b/.gitmodules index 6606b271..41de1b86 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,7 +1,7 @@ [submodule "mmm-physics"] path = schemes/mmm/mmm_physics url = https://github.com/NCAR/MMM-physics.git - fxtag = 20250616-MPASv8.3 + fxtag = 6cb29bc8f17bc3efdfaaa8a268ca53e2504aadbd fxrequired = AlwaysRequired fxDONOTUSEurl = https://github.com/NCAR/MMM-physics.git [submodule "pumas"] diff --git a/schemes/mmm/mmm_physics b/schemes/mmm/mmm_physics index a4baf7f3..6cb29bc8 160000 --- a/schemes/mmm/mmm_physics +++ b/schemes/mmm/mmm_physics @@ -1 +1 @@ -Subproject commit a4baf7f3243d1db0dbc5f63473f895bdbdc05c30 +Subproject commit 6cb29bc8f17bc3efdfaaa8a268ca53e2504aadbd From 5c866cefa2725383934da0a3ead510f05c4481f0 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 5 Dec 2025 15:45:12 -0700 Subject: [PATCH 06/22] Implement interstitial schemes for MYNN surface layer scheme --- schemes/mmm/mmm_physics_compat.F90 | 94 +++++++++++++++++- schemes/mmm/mmm_physics_compat.meta | 144 ++++++++++++++++++++++++++++ 2 files changed, 235 insertions(+), 3 deletions(-) diff --git a/schemes/mmm/mmm_physics_compat.F90 b/schemes/mmm/mmm_physics_compat.F90 index 1bd0dfcc..157af4e6 100644 --- a/schemes/mmm/mmm_physics_compat.F90 +++ b/schemes/mmm/mmm_physics_compat.F90 @@ -3,31 +3,61 @@ module mmm_physics_compat implicit none private + public :: mmm_physics_compat_init public :: mmm_physics_compat_run public :: mmm_physics_accumulate_tendencies_timestep_init public :: mmm_physics_accumulate_tendencies_run public :: mmm_physics_persist_states_init public :: mmm_physics_persist_states_timestep_final + public :: compute_air_density_run + public :: compute_atmosphere_layer_thickness_run public :: compute_characteristic_grid_length_scale_init public :: geopotential_height_wrt_sfc_at_if_to_msl_run public :: geopotential_height_wrt_sfc_to_msl_run contains + !> \section arg_table_mmm_physics_compat_init Argument Table + !! \htmlinclude mmm_physics_compat_init.html + pure subroutine mmm_physics_compat_init( & + isfflx, isftcflx, iz0tlnd, & + errmsg, errflg) + + integer, intent(out) :: isfflx, isftcflx, iz0tlnd + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + errmsg = '' + errflg = 0 + + ! Set options that are specific to MMM physics at model initialization. + isfflx = 1 + isftcflx = 0 + iz0tlnd = 0 + end subroutine mmm_physics_compat_init + !> \section arg_table_mmm_physics_compat_run Argument Table !! \htmlinclude mmm_physics_compat_run.html pure subroutine mmm_physics_compat_run( & nstep, & dt, & theta_curr, theta_prev, qv_curr, qv_prev, & + icefrac, landfrac, & scheme_name, & rthdynten, rqvdynten, & + xland, & errmsg, errflg) use ccpp_kinds, only: kind_phys + ! This threshold is hardcoded to the same value as in MMM physics. + ! It is named `xice_threshold` there. + real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys + integer, intent(in) :: nstep real(kind_phys), intent(in) :: dt, & - theta_curr(:, :), theta_prev(:, :), qv_curr(:, :), qv_prev(:, :) + theta_curr(:, :), theta_prev(:, :), qv_curr(:, :), qv_prev(:, :), & + icefrac(:), landfrac(:) character(256), intent(out) :: scheme_name - real(kind_phys), intent(out) :: rthdynten(:, :), rqvdynten(:, :) + real(kind_phys), intent(out) :: rthdynten(:, :), rqvdynten(:, :), & + xland(:) character(*), intent(out) :: errmsg integer, intent(out) :: errflg @@ -43,6 +73,16 @@ pure subroutine mmm_physics_compat_run( & rthdynten(:, :) = (theta_curr(:, :) - theta_prev(:, :)) / dt rqvdynten(:, :) = (qv_curr(:, :) - qv_prev(:, :)) / dt end if + + ! For MMM physics, land mask (`xland`) is defined as + ! * xland = 1.0 for land cells, including sea ice cells. + ! * xland = 2.0 for water cells. + where (landfrac >= 0.5_kind_phys .or. & + icefrac >= sea_ice_area_fraction_threshold) + xland = 1.0_kind_phys + elsewhere + xland = 2.0_kind_phys + end where end subroutine mmm_physics_compat_run !> \section arg_table_mmm_physics_accumulate_tendencies_timestep_init Argument Table @@ -193,10 +233,58 @@ pure subroutine mmm_physics_persist_states_timestep_final( & qv_prev(:, :) = qv_curr(:, :) end subroutine mmm_physics_persist_states_timestep_final + !> \section arg_table_compute_air_density_run Argument Table + !! \htmlinclude compute_air_density_run.html + pure subroutine compute_air_density_run( & + pdel, gravit, dz, & + rho, & + errmsg, errflg) + use ccpp_kinds, only: kind_phys + + real(kind_phys), intent(in) :: pdel(:, :), gravit, dz(:, :) + real(kind_phys), intent(out) :: rho(:, :) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + errmsg = '' + errflg = 0 + + ! Compute air density by hydrostatic equation. + rho(:, :) = pdel(:, :) / (gravit * dz(:, :)) + end subroutine compute_air_density_run + + !> \section arg_table_compute_atmosphere_layer_thickness_run Argument Table + !! \htmlinclude compute_atmosphere_layer_thickness_run.html + pure subroutine compute_atmosphere_layer_thickness_run( & + ncol, & + zisfc, & + dz, & + errmsg, errflg) + use ccpp_kinds, only: kind_phys + + integer, intent(in) :: ncol + real(kind_phys), intent(in) :: zisfc(:, :) + real(kind_phys), intent(out) :: dz(:, :) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + integer :: i + + errmsg = '' + errflg = 0 + + ! In CAM-SIMA, the first vertical index is at top of atmosphere. + ! The last one is at bottom of atmosphere. The resulting `dz` is positive. + do i = 1, ncol + dz(i, :) = zisfc(i, 1:size(zisfc, 2) - 1) - zisfc(i, 2:size(zisfc, 2)) + end do + end subroutine compute_atmosphere_layer_thickness_run + !> \section arg_table_compute_characteristic_grid_length_scale_init Argument Table !! \htmlinclude compute_characteristic_grid_length_scale_init.html pure subroutine compute_characteristic_grid_length_scale_init( & - omega, rearth, dx, & + omega, rearth, & + dx, & errmsg, errflg) use ccpp_kinds, only: kind_phys diff --git a/schemes/mmm/mmm_physics_compat.meta b/schemes/mmm/mmm_physics_compat.meta index 5cd2ae93..407e6ce5 100644 --- a/schemes/mmm/mmm_physics_compat.meta +++ b/schemes/mmm/mmm_physics_compat.meta @@ -2,6 +2,42 @@ name = mmm_physics_compat type = scheme +[ccpp-arg-table] + name = mmm_physics_compat_init + type = scheme +[ isfflx ] + standard_name = control_for_surface_flux_for_mmm_scheme + units = 1 + type = integer + dimensions = () + intent = out +[ isftcflx ] + standard_name = control_for_surface_roughness_length_over_water_for_mmm_scheme + units = 1 + type = integer + dimensions = () + intent = out +[ iz0tlnd ] + standard_name = control_for_surface_roughness_length_over_land_for_mmm_scheme + units = 1 + type = integer + dimensions = () + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + [ccpp-arg-table] name = mmm_physics_compat_run type = scheme @@ -41,6 +77,18 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = in +[ icefrac ] + standard_name = sea_ice_area_fraction_from_coupler + units = fraction + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ landfrac ] + standard_name = land_area_fraction_from_coupler + units = fraction + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in [ scheme_name ] standard_name = scheme_name units = none @@ -59,6 +107,12 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent, vertical_layer_dimension) intent = out +[ xland ] + standard_name = land_binary_mask_for_mmm_scheme + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out [ errmsg ] standard_name = ccpp_error_message long_name = Error message for error handling in CCPP @@ -412,6 +466,96 @@ # ---------- +[ccpp-table-properties] + name = compute_air_density + type = scheme + +[ccpp-arg-table] + name = compute_air_density_run + type = scheme +[ pdel ] + standard_name = air_pressure_thickness + units = Pa + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ gravit ] + standard_name = standard_gravitational_acceleration + units = m s-2 + type = real | kind = kind_phys + dimensions = () + intent = in +[ dz ] + standard_name = atmosphere_layer_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ rho ] + standard_name = air_density + units = kg m-3 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + +# ---------- + +[ccpp-table-properties] + name = compute_atmosphere_layer_thickness + type = scheme + +[ccpp-arg-table] + name = compute_atmosphere_layer_thickness_run + type = scheme +[ ncol ] + standard_name = horizontal_loop_extent + units = count + type = integer + dimensions = () + intent = in +[ zisfc ] + standard_name = geopotential_height_wrt_surface_at_interface + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_interface_dimension) + intent = in +[ dz ] + standard_name = atmosphere_layer_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + +# ---------- + [ccpp-table-properties] name = compute_characteristic_grid_length_scale type = scheme From 46aacf39d8d3e68796f6e52ee8040498fbf3b8b8 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Mon, 8 Dec 2025 13:33:57 -0700 Subject: [PATCH 07/22] Implement unit tests --- .../tests/mmm/mmm_physics_compat_tests.pf | 117 +++++++++++++++++- 1 file changed, 115 insertions(+), 2 deletions(-) diff --git a/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf b/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf index 0160fd60..69c88a5b 100644 --- a/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf +++ b/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf @@ -1,3 +1,28 @@ +@test +subroutine test_mmm_physics_compat_init() + use funit + use mmm_physics_compat, only: mmm_physics_compat_init + + integer :: isfflx, isftcflx, iz0tlnd + character(100) :: errmsg + integer :: errflg + + isfflx = huge(0) + isftcflx = huge(0) + iz0tlnd = huge(0) + + call mmm_physics_compat_init( & + isfflx, isftcflx, iz0tlnd, & + errmsg, errflg) + + ! Should set options that are specific to MMM physics. + @assertEqual(1, isfflx) + @assertEqual(0, isftcflx) + @assertEqual(0, iz0tlnd) + @assertEqual('', errmsg) + @assertEqual(0, errflg) +end subroutine test_mmm_physics_compat_init + @test subroutine test_mmm_physics_compat_run() use ccpp_kinds, only: kind_phys @@ -8,11 +33,15 @@ subroutine test_mmm_physics_compat_run() integer :: nstep real(kind_phys) :: dt real(kind_phys) :: theta_curr(ncol, pver), theta_prev(ncol, pver), qv_curr(ncol, pver), qv_prev(ncol, pver) + real(kind_phys) :: icefrac(ncol), landfrac(ncol) character(256) :: scheme_name real(kind_phys) :: rthdynten(ncol, pver), rqvdynten(ncol, pver) + real(kind_phys) :: xland(ncol) character(100) :: errmsg integer :: errflg + integer :: i + nstep = 0 dt = 20.0_kind_phys @@ -20,41 +49,65 @@ subroutine test_mmm_physics_compat_run() theta_prev(:, :) = 273.15_kind_phys qv_curr(:, :) = 0.03_kind_phys qv_prev(:, :) = 0.01_kind_phys + icefrac(1:50) = 0.0_kind_phys + icefrac(51:100) = [(0.01_kind_phys + real(i - 1, kind_phys) * 0.02_kind_phys, i = 1, 50)] + landfrac(1:50) = [(0.01_kind_phys + real(i - 1, kind_phys) * 0.02_kind_phys, i = 1, 50)] + landfrac(51:100) = 0.0_kind_phys + scheme_name = '' rthdynten(:, :) = huge(0.0_kind_phys) rqvdynten(:, :) = huge(0.0_kind_phys) + xland(:) = huge(0.0_kind_phys) call mmm_physics_compat_run( & nstep, & dt, & theta_curr, theta_prev, qv_curr, qv_prev, & + icefrac, landfrac, & scheme_name, & rthdynten, rqvdynten, & + xland, & errmsg, errflg) - ! Tendencies should be zero at start (nstep = 0). + ! Should set scheme name. @assertEqual('mmm_physics_compat_run', scheme_name) + ! Tendencies should be zero at start (nstep = 0). @assertEqual(0.0_kind_phys, rthdynten) @assertEqual(0.0_kind_phys, rqvdynten) + ! Should set land mask according to ice and land fractions. + @assertEqual(2.0_kind_phys, xland(1:25)) + @assertEqual(1.0_kind_phys, xland(26:50)) + @assertEqual(2.0_kind_phys, xland(51)) + @assertEqual(1.0_kind_phys, xland(52:100)) @assertEqual('', errmsg) @assertEqual(0, errflg) nstep = 1 + scheme_name = '' rthdynten(:, :) = huge(0.0_kind_phys) rqvdynten(:, :) = huge(0.0_kind_phys) + xland(:) = huge(0.0_kind_phys) call mmm_physics_compat_run( & nstep, & dt, & theta_curr, theta_prev, qv_curr, qv_prev, & + icefrac, landfrac, & scheme_name, & rthdynten, rqvdynten, & + xland, & errmsg, errflg) - ! Should compute tendencies correctly afterwards (nstep > 0). + ! Should set scheme name. @assertEqual('mmm_physics_compat_run', scheme_name) + ! Should compute tendencies correctly afterwards (nstep > 0). @assertEqual(1.0_kind_phys, rthdynten, epsilon(0.0_kind_phys)) @assertEqual(0.001_kind_phys, rqvdynten, epsilon(0.0_kind_phys)) + ! Should set land mask according to ice and land fractions. + @assertEqual(2.0_kind_phys, xland(1:25)) + @assertEqual(1.0_kind_phys, xland(26:50)) + @assertEqual(2.0_kind_phys, xland(51)) + @assertEqual(1.0_kind_phys, xland(52:100)) @assertEqual('', errmsg) @assertEqual(0, errflg) end subroutine test_mmm_physics_compat_run @@ -389,6 +442,66 @@ subroutine test_mmm_physics_persist_states_timestep_final() @assertEqual(0, errflg) end subroutine test_mmm_physics_persist_states_timestep_final +@test +subroutine test_compute_air_density_run() + use ccpp_kinds, only: kind_phys + use funit + use mmm_physics_compat, only: compute_air_density_run + + integer, parameter :: ncol = 100, pver = 10 + real(kind_phys) :: pdel(ncol, pver), gravit, dz(ncol, pver) + real(kind_phys) :: rho(ncol, pver) + character(100) :: errmsg + integer :: errflg + + pdel(:, :) = 98.0_kind_phys + gravit = 9.8_kind_phys + dz(:, :) = 10.0_kind_phys + rho(:, :) = huge(0.0_kind_phys) + + call compute_air_density_run( & + pdel, gravit, dz, & + rho, & + errmsg, errflg) + + ! Should compute air density correctly. + @assertEqual(1.0_kind_phys, rho, epsilon(0.0_kind_phys)) + @assertEqual('', errmsg) + @assertEqual(0, errflg) +end subroutine test_compute_air_density_run + +@test +subroutine test_compute_atmosphere_layer_thickness_run() + use ccpp_kinds, only: kind_phys + use funit + use mmm_physics_compat, only: compute_atmosphere_layer_thickness_run + + integer, parameter :: ncol = 100, pver = 10, pverp = 11 + real(kind_phys) :: zisfc(ncol, pverp) + real(kind_phys) :: dz(ncol, pver) + character(100) :: errmsg + integer :: errflg + + integer :: i + + do i = 1, pverp + zisfc(:, i) = 10000.0_kind_phys - real(i - 1, kind_phys) * 1000.0_kind_phys + end do + + dz(:, :) = huge(0.0_kind_phys) + + call compute_atmosphere_layer_thickness_run( & + ncol, & + zisfc, & + dz, & + errmsg, errflg) + + ! Should compute atmosphere layer thickness correctly. + @assertEqual(1000.0_kind_phys, dz, epsilon(0.0_kind_phys)) + @assertEqual('', errmsg) + @assertEqual(0, errflg) +end subroutine test_compute_atmosphere_layer_thickness_run + @test subroutine test_compute_characteristic_grid_length_scale_init() use ccpp_kinds, only: kind_phys From 5b9ebe4c695f1279e7c83b3c60fb06f3a23672b5 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 5 Dec 2025 15:47:00 -0700 Subject: [PATCH 08/22] Implement MYNN surface layer scheme --- schemes/mmm/sf_mynn_compat.F90 | 513 +++++++++++++++++ schemes/mmm/sf_mynn_compat.meta | 970 ++++++++++++++++++++++++++++++++ 2 files changed, 1483 insertions(+) create mode 100644 schemes/mmm/sf_mynn_compat.F90 create mode 100644 schemes/mmm/sf_mynn_compat.meta diff --git a/schemes/mmm/sf_mynn_compat.F90 b/schemes/mmm/sf_mynn_compat.F90 new file mode 100644 index 00000000..a1c60208 --- /dev/null +++ b/schemes/mmm/sf_mynn_compat.F90 @@ -0,0 +1,513 @@ +!> This module contains interstitial schemes that are specific to MYNN surface layer scheme, +!> which is part of MMM physics. +module sf_mynn_compat + implicit none + + private + public :: sf_mynn_compat_pre_run + public :: sf_mynn_compat_init + public :: sf_mynn_compat_run + public :: sf_mynn_diagnostics_init + public :: sf_mynn_diagnostics_run +contains + !> \section arg_table_sf_mynn_compat_pre_run Argument Table + !! \htmlinclude sf_mynn_compat_pre_run.html + subroutine sf_mynn_compat_pre_run( & + itimestep, & + u, v, t, qv, p, dz, rho, & + icefrac, landfrac, snowhice, snowhland, & + u1d, v1d, t1d, qv1d, p1d, dz8w1d, rho1d, & + u1d2, v1d2, dz2w1d, & + chs, chs2, cqs2, cpm, rmol, & + znt, ust, zol, mol, regime, psim, & + psih, qfx, & + flhc, flqc, snowh, qgh, qsfc, & + gz1oz0, wspd, br, svp1, svp2, & + svp3, svpt0, qcg, & + spp_pbl, rstoch1d, & + errmsg, errflg) + use ccpp_kinds, only: kind_phys + use ccpp_scheme_utils, only: ccpp_constituent_index + + ! This threshold is hardcoded to the same value as in MMM physics. + ! It is named `xice_threshold` there. + real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys + + integer, intent(in) :: itimestep + real(kind_phys), intent(in) :: u(:, :), v(:, :), t(:, :), qv(:, :), p(:, :), dz(:, :), rho(:, :), & + icefrac(:), landfrac(:), snowhice(:), snowhland(:) + logical, intent(out) :: spp_pbl + real(kind_phys), intent(out) :: u1d(:), v1d(:), t1d(:), qv1d(:), p1d(:), dz8w1d(:), rho1d(:), & + u1d2(:), v1d2(:), dz2w1d(:), & + chs(:), chs2(:), cqs2(:), cpm(:), rmol(:), & + znt(:), ust(:), zol(:), mol(:), regime(:), psim(:), & + psih(:), qfx(:), & + flhc(:), flqc(:), snowh(:), qgh(:), qsfc(:), & + gz1oz0(:), wspd(:), br(:), svp1, svp2, & + svp3, svpt0, qcg(:), & + rstoch1d(:) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + ! Provide first guess at the first time step. + if (itimestep == 1) then + ust(:) = max(0.04_kind_phys * sqrt(u(:, 1) ** 2 + v(:, 1) ** 2), 0.001_kind_phys) + mol(:) = 0.0_kind_phys + qsfc(:) = qv(:, 1) ! Note that this `qv` is wet. + end if + + ! In CAM-SIMA, the first vertical index is at top of atmosphere. + ! The last one is at bottom of atmosphere, which is what we want here. + + u1d(:) = u(:, size(u, 2)) + v1d(:) = v(:, size(v, 2)) + t1d(:) = t(:, size(t, 2)) + qv1d(:) = qv(:, size(qv, 2)) + p1d(:) = p(:, size(p, 2)) + dz8w1d(:) = dz(:, size(dz, 2)) + rho1d(:) = rho(:, size(rho, 2)) + + u1d2(:) = u(:, size(u, 2) - 1) + v1d2(:) = v(:, size(v, 2) - 1) + dz2w1d(:) = dz(:, size(dz, 2) - 1) + + chs(:) = 0.0_kind_phys + chs2(:) = 0.0_kind_phys + cqs2(:) = 0.0_kind_phys + cpm(:) = 0.0_kind_phys + rmol(:) = 0.0_kind_phys + + znt(:) = 0.0_kind_phys + zol(:) = 0.0_kind_phys + regime(:) = 0.0_kind_phys + psim(:) = 0.0_kind_phys + + psih(:) = 0.0_kind_phys + qfx(:) = 0.0_kind_phys + + flhc(:) = 0.0_kind_phys + flqc(:) = 0.0_kind_phys + snowh(:) = 0.0_kind_phys + + where (landfrac >= 0.5_kind_phys) + snowh = snowhland + end where + + where (icefrac >= sea_ice_area_fraction_threshold) + snowh = snowhice + end where + + qgh(:) = 0.0_kind_phys + + gz1oz0(:) = 0.0_kind_phys + wspd(:) = 0.0_kind_phys + br(:) = 0.0_kind_phys + + ! Constants in equation 10 from Bolton (1980). See + ! doi:10.1175/1520-0493(1980)108<1046:TCOEPT>2.0.CO;2. + svp1 = 6.112_kind_phys + svp2 = 17.67_kind_phys + svp3 = 29.65_kind_phys + svpt0 = 273.15_kind_phys + + qcg(:) = 0.0_kind_phys ! Not used but still appear in the argument list... + + spp_pbl = .false. + rstoch1d(:) = 0.0_kind_phys + + errmsg = '' + errflg = 0 + end subroutine sf_mynn_compat_pre_run + + !> \section arg_table_sf_mynn_compat_init Argument Table + !! \htmlinclude sf_mynn_compat_init.html + subroutine sf_mynn_compat_init( & + ust, mol, qsfc, & + errmsg, errflg) + use ccpp_kinds, only: kind_phys + use sf_mynn, only: sf_mynn_init + + real(kind_phys), intent(out) :: ust(:), mol(:), qsfc(:) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + ! Precompute lookup tables in MYNN surface layer scheme. + call sf_mynn_init(errmsg, errflg) + + ! MYNN surface layer scheme takes time averages of these variables internally. + ! As a result, they must be able to persist across time steps. + ust(:) = 0.0_kind_phys + mol(:) = 0.0_kind_phys + qsfc(:) = 0.0_kind_phys + + errmsg = '' + errflg = 0 + end subroutine sf_mynn_compat_init + + !> \section arg_table_sf_mynn_compat_run Argument Table + !! \htmlinclude sf_mynn_compat_run.html + subroutine sf_mynn_compat_run( & + ncol, cflx, icefrac, sst, & + u1d, v1d, t1d, qv1d, p1d, dz8w1d, rho1d, & + u1d2, v1d2, dz2w1d, cp, g, rovcp, r, xlv, & + psfcpa, chs, chs2, cqs2, cpm, pblh, rmol, & + znt, ust, mavail, zol, mol, regime, psim, & + psih, xland, hfx, qfx, tsk, u10, v10, th2, & + t2, q2, flhc, flqc, snowh, qgh, qsfc, lh, & + gz1oz0, wspd, br, isfflx, dx, svp1, svp2, & + svp3, svpt0, ep1, ep2, karman, qcg, & + itimestep, wstar, qstar, ustm, ck, cka, & + cd, cda, spp_pbl, rstoch1d, isftcflx, & + iz0tlnd, its, ite, & + errmsg, errflg) + use ccpp_kinds, only: kind_phys + use ccpp_scheme_utils, only: ccpp_constituent_index + use sf_mynn, only: sf_mynn_run + + ! This threshold is hardcoded to the same value as in MMM physics. + ! It is named `xice_threshold` there. + real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys + ! Typical threshold between sea surface temperature and ice surface temperature. See + ! Algorithm Theoretical Basis Document (ATBD) for the MODIS Snow and Sea Ice-Mapping Algorithms, + ! Section 4.4.3 Ice Surface Temperature (IST) Algorithm. + real(kind_phys), parameter :: sea_ice_temperature_threshold = 271.4_kind_phys + + integer, intent(in) :: ncol, & + isfflx, & + itimestep, & + isftcflx, & + iz0tlnd, its, ite + logical, intent(in) :: spp_pbl + real(kind_phys), intent(in) :: icefrac(:), sst(:), & + u1d(:), v1d(:), t1d(:), qv1d(:), p1d(:), dz8w1d(:), rho1d(:), & + u1d2(:), v1d2(:), dz2w1d(:), cp, g, rovcp, r, xlv, & + psfcpa(:), pblh(:), & + mavail(:), & + xland(:), tsk(:), & + snowh(:), & + dx(:), svp1, svp2, & + svp3, svpt0, ep1, ep2, karman, qcg(:), & + rstoch1d(:) + real(kind_phys), intent(inout) :: cflx(:, :), & + chs(:), chs2(:), cqs2(:), cpm(:), rmol(:), & + znt(:), ust(:), zol(:), mol(:), regime(:), psim(:), & + psih(:), hfx(:), qfx(:), & + flhc(:), flqc(:), qgh(:), qsfc(:), lh(:), & + gz1oz0(:), wspd(:), br(:) + real(kind_phys), intent(out) :: u10(:), v10(:), th2(:), & + t2(:), q2(:), & + wstar(:), qstar(:), ustm(:), ck(:), cka(:), & + cd(:), cda(:) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + integer :: water_vapor_mixing_ratio_index + real(kind_phys), allocatable :: ch(:) + + ! Special handling for sea ice cells like in MMM physics. + logical, allocatable :: mask_sea_ice_cell(:) + real(kind_phys), allocatable :: mavail_sea(:), & + xland_sea(:), tsk_sea(:) + real(kind_phys), allocatable :: chs_sea(:), chs2_sea(:), cqs2_sea(:), cpm_sea(:), rmol_sea(:), & + znt_sea(:), ust_sea(:), zol_sea(:), mol_sea(:), regime_sea(:), psim_sea(:), & + psih_sea(:), hfx_sea(:), qfx_sea(:), & + flhc_sea(:), flqc_sea(:), qgh_sea(:), qsfc_sea(:), lh_sea(:), & + gz1oz0_sea(:), wspd_sea(:), br_sea(:), & + ch_sea(:) + real(kind_phys), allocatable :: u10_sea(:), v10_sea(:), th2_sea(:), & + t2_sea(:), q2_sea(:), & + wstar_sea(:), qstar_sea(:), ustm_sea(:), ck_sea(:), cka_sea(:), & + cd_sea(:), cda_sea(:) + + ! `ch` is duplicate of `chs`, but for unknown reasons it is passed separately in the argument list + ! to MYNN surface layer scheme. + allocate(ch(ncol), errmsg=errmsg, stat=errflg) + + if (errflg /= 0) then + return + end if + + ch(:) = chs(:) + + qfx(:) = 0.0_kind_phys + + call ccpp_constituent_index( & + 'water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water', water_vapor_mixing_ratio_index, errflg, errmsg) + + if (errflg /= 0 .or. & + water_vapor_mixing_ratio_index < lbound(cflx, 2) .or. water_vapor_mixing_ratio_index > ubound(cflx, 2)) then + errmsg = 'Failed to find desired constituent flux from cflx' + errflg = 1 + + return + end if + + qfx(:) = cflx(:, water_vapor_mixing_ratio_index) + + allocate(mask_sea_ice_cell(ncol), errmsg=errmsg, stat=errflg) + + if (errflg /= 0) then + return + end if + + mask_sea_ice_cell(:) = (icefrac >= sea_ice_area_fraction_threshold) + + ! If there are sea ice cells, make local copies of the variables in preparation for the second pass. + if (any(mask_sea_ice_cell)) then + allocate(mavail_sea(ncol), xland_sea(ncol), tsk_sea(ncol), & + errmsg=errmsg, stat=errflg) + + if (errflg /= 0) then + return + end if + + mavail_sea(:) = mavail(:) + xland_sea(:) = xland(:) + tsk_sea(:) = tsk(:) + + allocate(chs_sea(ncol), chs2_sea(ncol), cqs2_sea(ncol), cpm_sea(ncol), rmol_sea(ncol), & + znt_sea(ncol), ust_sea(ncol), zol_sea(ncol), mol_sea(ncol), regime_sea(ncol), psim_sea(ncol), & + psih_sea(ncol), hfx_sea(ncol), qfx_sea(ncol), & + flhc_sea(ncol), flqc_sea(ncol), qgh_sea(ncol), qsfc_sea(ncol), lh_sea(ncol), & + gz1oz0_sea(ncol), wspd_sea(ncol), br_sea(ncol), & + ch_sea(ncol), & + errmsg=errmsg, stat=errflg) + + if (errflg /= 0) then + return + end if + + chs_sea(:) = chs(:) + chs2_sea(:) = chs2(:) + cqs2_sea(:) = cqs2(:) + cpm_sea(:) = cpm(:) + rmol_sea(:) = rmol(:) + znt_sea(:) = znt(:) + ust_sea(:) = ust(:) + zol_sea(:) = zol(:) + mol_sea(:) = mol(:) + regime_sea(:) = regime(:) + psim_sea(:) = psim(:) + psih_sea(:) = psih(:) + hfx_sea(:) = hfx(:) + qfx_sea(:) = qfx(:) + flhc_sea(:) = flhc(:) + flqc_sea(:) = flqc(:) + qgh_sea(:) = qgh(:) + qsfc_sea(:) = qsfc(:) + lh_sea(:) = lh(:) + gz1oz0_sea(:) = gz1oz0(:) + wspd_sea(:) = wspd(:) + br_sea(:) = br(:) + ch_sea(:) = ch(:) + + allocate(u10_sea(ncol), v10_sea(ncol), th2_sea(ncol), & + t2_sea(ncol), q2_sea(ncol), & + wstar_sea(ncol), qstar_sea(ncol), ustm_sea(ncol), ck_sea(ncol), cka_sea(ncol), & + cd_sea(ncol), cda_sea(ncol), & + errmsg=errmsg, stat=errflg) + + if (errflg /= 0) then + return + end if + + u10_sea(:) = u10(:) + v10_sea(:) = v10(:) + th2_sea(:) = th2(:) + t2_sea(:) = t2(:) + q2_sea(:) = q2(:) + wstar_sea(:) = wstar(:) + qstar_sea(:) = qstar(:) + ustm_sea(:) = ustm(:) + ck_sea(:) = ck(:) + cka_sea(:) = cka(:) + cd_sea(:) = cd(:) + cda_sea(:) = cda(:) + + where (mask_sea_ice_cell) + ! Set surface moisture availability to maximum. + mavail_sea = 1.0_kind_phys + + ! Impose minimum skin temperature. + tsk_sea = max(sea_ice_temperature_threshold, sst) + + ! Treat as water cells. + xland_sea = 2.0_kind_phys + + ! Set surface roughness length to calm water. + znt_sea = 0.0001_kind_phys + end where + end if + + ! First pass for all cells. + call sf_mynn_run( & + u1d, v1d, t1d, qv1d, p1d, dz8w1d, rho1d, & + u1d2, v1d2, dz2w1d, cp, g, rovcp, r, xlv, & + psfcpa, chs, chs2, cqs2, cpm, pblh, rmol, & + znt, ust, mavail, zol, mol, regime, psim, & + psih, xland, hfx, qfx, tsk, u10, v10, th2, & + t2, q2, flhc, flqc, snowh, qgh, qsfc, lh, & + gz1oz0, wspd, br, isfflx, dx, svp1, svp2, & + svp3, svpt0, ep1, ep2, karman, ch, qcg, & + itimestep, wstar, qstar, ustm, ck, cka, & + cd, cda, spp_pbl, rstoch1d, isftcflx, & + iz0tlnd, its, ite, & + errmsg, errflg) + + ! The first pass treated sea ice cells as land cells due to how `xland` is defined. + ! The second pass treats sea ice cells as water cells instead. + ! Then, for sea ice cells only, the final results are weighted averages according to their area fractions. + if (any(mask_sea_ice_cell)) then + call sf_mynn_run( & + u1d, v1d, t1d, qv1d, p1d, dz8w1d, rho1d, & + u1d2, v1d2, dz2w1d, cp, g, rovcp, r, xlv, & + psfcpa, chs_sea, chs2_sea, cqs2_sea, cpm_sea, pblh, rmol_sea, & + znt_sea, ust_sea, mavail_sea, zol_sea, mol_sea, regime_sea, psim_sea, & + psih_sea, xland_sea, hfx_sea, qfx_sea, tsk_sea, u10_sea, v10_sea, th2_sea, & + t2_sea, q2_sea, flhc_sea, flqc_sea, snowh, qgh_sea, qsfc_sea, lh_sea, & + gz1oz0_sea, wspd_sea, br_sea, isfflx, dx, svp1, svp2, & + svp3, svpt0, ep1, ep2, karman, ch_sea, qcg, & + itimestep, wstar_sea, qstar_sea, ustm_sea, ck_sea, cka_sea, & + cd_sea, cda_sea, spp_pbl, rstoch1d, isftcflx, & + iz0tlnd, its, ite, & + errmsg, errflg) + + ! Blend the final results between sea ice and sea water. + where (mask_sea_ice_cell) + chs = chs * icefrac + (1.0_kind_phys - icefrac) * chs_sea + chs2 = chs2 * icefrac + (1.0_kind_phys - icefrac) * chs2_sea + cqs2 = cqs2 * icefrac + (1.0_kind_phys - icefrac) * cqs2_sea + cpm = cpm * icefrac + (1.0_kind_phys - icefrac) * cpm_sea + rmol = rmol * icefrac + (1.0_kind_phys - icefrac) * rmol_sea + znt = znt * icefrac + (1.0_kind_phys - icefrac) * znt_sea + ust = ust * icefrac + (1.0_kind_phys - icefrac) * ust_sea + zol = zol * icefrac + (1.0_kind_phys - icefrac) * zol_sea + mol = mol * icefrac + (1.0_kind_phys - icefrac) * mol_sea + psim = psim * icefrac + (1.0_kind_phys - icefrac) * psim_sea + psih = psih * icefrac + (1.0_kind_phys - icefrac) * psih_sea + hfx = hfx * icefrac + (1.0_kind_phys - icefrac) * hfx_sea + qfx = qfx * icefrac + (1.0_kind_phys - icefrac) * qfx_sea + flhc = flhc * icefrac + (1.0_kind_phys - icefrac) * flhc_sea + flqc = flqc * icefrac + (1.0_kind_phys - icefrac) * flqc_sea + qgh = qgh * icefrac + (1.0_kind_phys - icefrac) * qgh_sea + qsfc = qsfc * icefrac + (1.0_kind_phys - icefrac) * qsfc_sea + lh = lh * icefrac + (1.0_kind_phys - icefrac) * lh_sea + gz1oz0 = gz1oz0 * icefrac + (1.0_kind_phys - icefrac) * gz1oz0_sea + wspd = wspd * icefrac + (1.0_kind_phys - icefrac) * wspd_sea + br = br * icefrac + (1.0_kind_phys - icefrac) * br_sea + + u10 = u10 * icefrac + (1.0_kind_phys - icefrac) * u10_sea + v10 = v10 * icefrac + (1.0_kind_phys - icefrac) * v10_sea + th2 = th2 * icefrac + (1.0_kind_phys - icefrac) * th2_sea + t2 = t2 * icefrac + (1.0_kind_phys - icefrac) * t2_sea + q2 = q2 * icefrac + (1.0_kind_phys - icefrac) * q2_sea + wstar = wstar * icefrac + (1.0_kind_phys - icefrac) * wstar_sea + qstar = qstar * icefrac + (1.0_kind_phys - icefrac) * qstar_sea + ustm = ustm * icefrac + (1.0_kind_phys - icefrac) * ustm_sea + ck = ck * icefrac + (1.0_kind_phys - icefrac) * ck_sea + cka = cka * icefrac + (1.0_kind_phys - icefrac) * cka_sea + cd = cd * icefrac + (1.0_kind_phys - icefrac) * cd_sea + cda = cda * icefrac + (1.0_kind_phys - icefrac) * cda_sea + end where + + ! `regime` is a categorical variable. Assign according to whether sea ice or sea water is dominant. + where (mask_sea_ice_cell .and. icefrac < 0.5_kind_phys) + regime = regime_sea + end where + end if + + cflx(:, water_vapor_mixing_ratio_index) = qfx(:) + + errmsg = '' + errflg = 0 + end subroutine sf_mynn_compat_run + + !> \section arg_table_sf_mynn_diagnostics_init Argument Table + !! \htmlinclude sf_mynn_diagnostics_init.html + subroutine sf_mynn_diagnostics_init( & + errmsg, errflg) + use cam_history, only: history_add_field + use cam_history_support, only: horiz_only + + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + call history_add_field('sf_mynn_lh', & + 'surface_upward_latent_heat_flux_from_coupler', horiz_only, 'avg', 'W m-2') + call history_add_field('sf_mynn_hfx', & + 'surface_upward_sensible_heat_flux_from_coupler', horiz_only, 'avg', 'W m-2') + call history_add_field('sf_mynn_qfx', & + 'surface_upward_water_vapor_flux', horiz_only, 'avg', 'kg m-2 s-1') + + call history_add_field('sf_mynn_cd', & + 'drag_coefficient_for_momentum_at_10m', horiz_only, 'avg', '1') + call history_add_field('sf_mynn_cda', & + 'drag_coefficient_for_momentum_at_surface_adjacent_layer', horiz_only, 'avg', '1') + call history_add_field('sf_mynn_ck', & + 'bulk_exchange_coefficient_for_enthalpy_at_2m', horiz_only, 'avg', '1') + call history_add_field('sf_mynn_cka', & + 'bulk_exchange_coefficient_for_enthalpy_at_surface_adjacent_layer', horiz_only, 'avg', '1') + + call history_add_field('sf_mynn_q2', & + 'water_vapor_mixing_ratio_wrt_dry_air_at_2m', horiz_only, 'avg', 'kg kg-1') + call history_add_field('sf_mynn_t2', & + 'air_temperature_at_2m', horiz_only, 'avg', 'K') + call history_add_field('sf_mynn_th2', & + 'air_potential_temperature_at_2m', horiz_only, 'avg', 'K') + call history_add_field('sf_mynn_u10', & + 'eastward_wind_at_10m', horiz_only, 'avg', 'm s-1') + call history_add_field('sf_mynn_v10', & + 'northward_wind_at_10m', horiz_only, 'avg', 'm s-1') + + call history_add_field('sf_mynn_ustm', & + 'surface_friction_velocity_assuming_no_correction_for_surface_convective_velocity_scale', horiz_only, 'avg', 'm s-1') + call history_add_field('sf_mynn_wstar', & + 'surface_convective_velocity_scale', horiz_only, 'avg', 'm s-1') + call history_add_field('sf_mynn_qstar', & + 'surface_moisture_scale', horiz_only, 'avg', 'g kg-1') + + errmsg = '' + errflg = 0 + end subroutine sf_mynn_diagnostics_init + + !> \section arg_table_sf_mynn_diagnostics_run Argument Table + !! \htmlinclude sf_mynn_diagnostics_run.html + subroutine sf_mynn_diagnostics_run( & + lh, hfx, qfx, & + cd, cda, ck, cka, & + q2, t2, th2, u10, v10, & + ustm, wstar, qstar, & + errmsg, errflg) + use cam_history, only: history_out_field + use ccpp_kinds, only: kind_phys + + real(kind_phys), intent(in) :: lh(:), hfx(:), qfx(:), & + cd(:), cda(:), ck(:), cka(:), & + q2(:), t2(:), th2(:), u10(:), v10(:), & + ustm(:), wstar(:), qstar(:) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + call history_out_field('sf_mynn_lh', lh) + call history_out_field('sf_mynn_hfx', hfx) + call history_out_field('sf_mynn_qfx', qfx) + + call history_out_field('sf_mynn_cd', cd) + call history_out_field('sf_mynn_cda', cda) + call history_out_field('sf_mynn_ck', ck) + call history_out_field('sf_mynn_cka', cka) + + call history_out_field('sf_mynn_q2', q2) + call history_out_field('sf_mynn_t2', t2) + call history_out_field('sf_mynn_th2', th2) + call history_out_field('sf_mynn_u10', u10) + call history_out_field('sf_mynn_v10', v10) + + call history_out_field('sf_mynn_ustm', ustm) + call history_out_field('sf_mynn_wstar', wstar) + call history_out_field('sf_mynn_qstar', qstar) + + errmsg = '' + errflg = 0 + end subroutine sf_mynn_diagnostics_run +end module sf_mynn_compat diff --git a/schemes/mmm/sf_mynn_compat.meta b/schemes/mmm/sf_mynn_compat.meta new file mode 100644 index 00000000..09c0620f --- /dev/null +++ b/schemes/mmm/sf_mynn_compat.meta @@ -0,0 +1,970 @@ +[ccpp-table-properties] + name = sf_mynn_compat_pre + type = scheme + +[ccpp-arg-table] + name = sf_mynn_compat_pre_run + type = scheme +[ itimestep ] + standard_name = current_timestep_number + units = count + type = integer + dimensions = () + intent = in +[ u ] + standard_name = eastward_wind + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ v ] + standard_name = northward_wind + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ t ] + standard_name = air_temperature + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ qv ] + standard_name = water_vapor_mixing_ratio_wrt_dry_air + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ p ] + standard_name = air_pressure + units = Pa + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ dz ] + standard_name = atmosphere_layer_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ rho ] + standard_name = air_density + units = kg m-3 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ icefrac ] + standard_name = sea_ice_area_fraction_from_coupler + units = fraction + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ landfrac ] + standard_name = land_area_fraction_from_coupler + units = fraction + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ snowhice ] + standard_name = lwe_surface_snow_depth_over_ice_from_coupler + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ snowhland ] + standard_name = lwe_surface_snow_depth_over_land_from_coupler + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ u1d ] + standard_name = eastward_wind_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ v1d ] + standard_name = northward_wind_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ t1d ] + standard_name = air_temperature_at_surface_adjacent_layer + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ qv1d ] + standard_name = water_vapor_mixing_ratio_wrt_dry_air_at_surface_adjacent_layer + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ p1d ] + standard_name = air_pressure_at_surface_adjacent_layer + units = Pa + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ dz8w1d ] + standard_name = atmosphere_layer_thickness_at_surface_adjacent_layer + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ rho1d ] + standard_name = air_density_at_surface_adjacent_layer + units = kg m-3 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ u1d2 ] + standard_name = eastward_wind_at_second_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ v1d2 ] + standard_name = northward_wind_at_second_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ dz2w1d ] + standard_name = atmosphere_layer_thickness_at_second_surface_adjacent_layer + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ chs ] + standard_name = exchange_coefficient_for_heat_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ chs2 ] + standard_name = exchange_coefficient_for_heat_at_2m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ cqs2 ] + standard_name = exchange_coefficient_for_moisture_at_2m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ cpm ] + standard_name = composition_dependent_specific_heat_of_dry_air_at_constant_pressure_at_surface_adjacent_layer + units = J kg-1 K-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ rmol ] + standard_name = reciprocal_of_monin_obukhov_length + units = m-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ znt ] + standard_name = surface_roughness_length + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ ust ] + standard_name = surface_friction_velocity + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ zol ] + standard_name = ratio_of_height_to_monin_obukhov_length + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ mol ] + standard_name = surface_temperature_scale + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ regime ] + standard_name = control_for_pbl_stability_regime + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ psim ] + standard_name = monin_obukhov_similarity_function_for_momentum + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ psih ] + standard_name = monin_obukhov_similarity_function_for_heat + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ qfx ] + standard_name = surface_upward_water_vapor_flux + units = kg m-2 s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ flhc ] + standard_name = surface_kinematic_exchange_coefficient_for_heat + units = W m-2 K-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ flqc ] + standard_name = surface_kinematic_exchange_coefficient_for_moisture + units = kg m-2 s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ snowh ] + standard_name = surface_snow_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ qgh ] + standard_name = water_vapor_mixing_ratio_wrt_dry_air_at_surface_adjacent_layer_assuming_saturation + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ qsfc ] + standard_name = water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water_at_surface + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ gz1oz0 ] + standard_name = ln_ratio_of_height_at_surface_adjacent_layer_to_surface_roughness_length + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ wspd ] + standard_name = wind_speed_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ br ] + standard_name = bulk_richardson_number_at_surface_adjacent_layer + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ svp1 ] + standard_name = constant_1_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = out +[ svp2 ] + standard_name = constant_2_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = out +[ svp3 ] + standard_name = constant_3_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = out +[ svpt0 ] + standard_name = constant_4_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = out +[ qcg ] + standard_name = cloud_liquid_water_mixing_ratio_wrt_dry_air_at_surface + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ spp_pbl ] + standard_name = flag_for_stochastically_perturbed_parameterization_for_mynn_scheme + units = flag + type = logical + dimensions = () + intent = out +[ rstoch1d ] + standard_name = stochastically_perturbed_weights_for_mynn_scheme + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + +# ---------- + +[ccpp-table-properties] + name = sf_mynn_compat + type = scheme + dependencies = ccpp_kind_types.F90, mmm_physics/mynn_shared.F90, mmm_physics/sf_mynn.F90 + +[ccpp-arg-table] + name = sf_mynn_compat_init + type = scheme +[ ust ] + standard_name = surface_friction_velocity + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_dimension) + intent = out +[ mol ] + standard_name = surface_temperature_scale + units = K + type = real | kind = kind_phys + dimensions = (horizontal_dimension) + intent = out +[ qsfc ] + standard_name = water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water_at_surface + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_dimension) + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + +[ccpp-arg-table] + name = sf_mynn_compat_run + type = scheme +[ ncol ] + standard_name = horizontal_loop_extent + units = count + type = integer + dimensions = () + intent = in +[ cflx ] + standard_name = surface_upward_ccpp_constituent_fluxes_from_coupler + units = kg m-2 s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, number_of_ccpp_constituents) + intent = inout +[ icefrac ] + standard_name = sea_ice_area_fraction_from_coupler + units = fraction + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ sst ] + standard_name = sea_surface_temperature_from_coupler + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ u1d ] + standard_name = eastward_wind_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ v1d ] + standard_name = northward_wind_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ t1d ] + standard_name = air_temperature_at_surface_adjacent_layer + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ qv1d ] + standard_name = water_vapor_mixing_ratio_wrt_dry_air_at_surface_adjacent_layer + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ p1d ] + standard_name = air_pressure_at_surface_adjacent_layer + units = Pa + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ dz8w1d ] + standard_name = atmosphere_layer_thickness_at_surface_adjacent_layer + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ rho1d ] + standard_name = air_density_at_surface_adjacent_layer + units = kg m-3 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ u1d2 ] + standard_name = eastward_wind_at_second_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ v1d2 ] + standard_name = northward_wind_at_second_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ dz2w1d ] + standard_name = atmosphere_layer_thickness_at_second_surface_adjacent_layer + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ cp ] + standard_name = specific_heat_of_dry_air_at_constant_pressure + units = J kg-1 K-1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ g ] + standard_name = standard_gravitational_acceleration + units = m s-2 + type = real | kind = kind_phys + dimensions = () + intent = in +[ rovcp ] + standard_name = ratio_of_dry_air_gas_constant_to_specific_heat_of_dry_air_at_constant_pressure + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ r ] + standard_name = gas_constant_of_dry_air + units = J kg-1 K-1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ xlv ] + standard_name = latent_heat_of_vaporization_of_water_at_0c + units = J kg-1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ psfcpa ] + standard_name = surface_air_pressure + units = Pa + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ chs ] + standard_name = exchange_coefficient_for_heat_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ chs2 ] + standard_name = exchange_coefficient_for_heat_at_2m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ cqs2 ] + standard_name = exchange_coefficient_for_moisture_at_2m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ cpm ] + standard_name = composition_dependent_specific_heat_of_dry_air_at_constant_pressure_at_surface_adjacent_layer + units = J kg-1 K-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ pblh ] + standard_name = atmosphere_boundary_layer_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ rmol ] + standard_name = reciprocal_of_monin_obukhov_length + units = m-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ znt ] + standard_name = surface_roughness_length + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ ust ] + standard_name = surface_friction_velocity + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ mavail ] + standard_name = surface_moisture_availability + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ zol ] + standard_name = ratio_of_height_to_monin_obukhov_length + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ mol ] + standard_name = surface_temperature_scale + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ regime ] + standard_name = control_for_pbl_stability_regime + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ psim ] + standard_name = monin_obukhov_similarity_function_for_momentum + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ psih ] + standard_name = monin_obukhov_similarity_function_for_heat + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ xland ] + standard_name = land_binary_mask_for_mmm_scheme + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ hfx ] + standard_name = surface_upward_sensible_heat_flux_from_coupler + units = W m-2 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ qfx ] + standard_name = surface_upward_water_vapor_flux + units = kg m-2 s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ tsk ] + standard_name = skin_temperature_at_surface + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ u10 ] + standard_name = eastward_wind_at_10m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ v10 ] + standard_name = northward_wind_at_10m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ th2 ] + standard_name = air_potential_temperature_at_2m + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ t2 ] + standard_name = air_temperature_at_2m + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ q2 ] + standard_name = water_vapor_mixing_ratio_wrt_dry_air_at_2m + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ flhc ] + standard_name = surface_kinematic_exchange_coefficient_for_heat + units = W m-2 K-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ flqc ] + standard_name = surface_kinematic_exchange_coefficient_for_moisture + units = kg m-2 s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ snowh ] + standard_name = surface_snow_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ qgh ] + standard_name = water_vapor_mixing_ratio_wrt_dry_air_at_surface_adjacent_layer_assuming_saturation + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ qsfc ] + standard_name = water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water_at_surface + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ lh ] + standard_name = surface_upward_latent_heat_flux_from_coupler + units = W m-2 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ gz1oz0 ] + standard_name = ln_ratio_of_height_at_surface_adjacent_layer_to_surface_roughness_length + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ wspd ] + standard_name = wind_speed_at_surface_adjacent_layer + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ br ] + standard_name = bulk_richardson_number_at_surface_adjacent_layer + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = inout +[ isfflx ] + standard_name = control_for_surface_flux_for_mmm_scheme + units = 1 + type = integer + dimensions = () + intent = in +[ dx ] + standard_name = characteristic_grid_lengthscale + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ svp1 ] + standard_name = constant_1_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ svp2 ] + standard_name = constant_2_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ svp3 ] + standard_name = constant_3_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ svpt0 ] + standard_name = constant_4_in_formula_for_water_vapor_partial_pressure_assuming_saturation + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ ep1 ] + standard_name = ratio_of_water_vapor_to_dry_air_gas_constants_minus_one + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ ep2 ] + standard_name = ratio_of_water_vapor_to_dry_air_molecular_weights + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ karman ] + standard_name = von_karman_constant + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in +[ qcg ] + standard_name = cloud_liquid_water_mixing_ratio_wrt_dry_air_at_surface + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ itimestep ] + standard_name = current_timestep_number + units = count + type = integer + dimensions = () + intent = in +[ wstar ] + standard_name = surface_convective_velocity_scale + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ qstar ] + standard_name = surface_moisture_scale + units = g kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ ustm ] + standard_name = surface_friction_velocity_assuming_no_correction_for_surface_convective_velocity_scale + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ ck ] + standard_name = bulk_exchange_coefficient_for_enthalpy_at_2m + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ cka ] + standard_name = bulk_exchange_coefficient_for_enthalpy_at_surface_adjacent_layer + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ cd ] + standard_name = drag_coefficient_for_momentum_at_10m + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ cda ] + standard_name = drag_coefficient_for_momentum_at_surface_adjacent_layer + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = out +[ spp_pbl ] + standard_name = flag_for_stochastically_perturbed_parameterization_for_mynn_scheme + units = flag + type = logical + dimensions = () + intent = in +[ rstoch1d ] + standard_name = stochastically_perturbed_weights_for_mynn_scheme + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ isftcflx ] + standard_name = control_for_surface_roughness_length_over_water_for_mmm_scheme + units = 1 + type = integer + dimensions = () + intent = in +[ iz0tlnd ] + standard_name = control_for_surface_roughness_length_over_land_for_mmm_scheme + units = 1 + type = integer + dimensions = () + intent = in +[ its ] + standard_name = horizontal_loop_begin + units = count + type = integer + dimensions = () + intent = in +[ ite ] + standard_name = horizontal_loop_end + units = count + type = integer + dimensions = () + intent = in +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + +# ---------- + +[ccpp-table-properties] + name = sf_mynn_diagnostics + type = scheme + +[ccpp-arg-table] + name = sf_mynn_diagnostics_init + type = scheme +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + +[ccpp-arg-table] + name = sf_mynn_diagnostics_run + type = scheme +[ lh ] + standard_name = surface_upward_latent_heat_flux_from_coupler + units = W m-2 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ hfx ] + standard_name = surface_upward_sensible_heat_flux_from_coupler + units = W m-2 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ qfx ] + standard_name = surface_upward_water_vapor_flux + units = kg m-2 s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ cd ] + standard_name = drag_coefficient_for_momentum_at_10m + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ cda ] + standard_name = drag_coefficient_for_momentum_at_surface_adjacent_layer + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ ck ] + standard_name = bulk_exchange_coefficient_for_enthalpy_at_2m + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ cka ] + standard_name = bulk_exchange_coefficient_for_enthalpy_at_surface_adjacent_layer + units = 1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ q2 ] + standard_name = water_vapor_mixing_ratio_wrt_dry_air_at_2m + units = kg kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ t2 ] + standard_name = air_temperature_at_2m + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ th2 ] + standard_name = air_potential_temperature_at_2m + units = K + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ u10 ] + standard_name = eastward_wind_at_10m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ v10 ] + standard_name = northward_wind_at_10m + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ ustm ] + standard_name = surface_friction_velocity_assuming_no_correction_for_surface_convective_velocity_scale + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ wstar ] + standard_name = surface_convective_velocity_scale + units = m s-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ qstar ] + standard_name = surface_moisture_scale + units = g kg-1 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent) + intent = in +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out From dbb780df7863d149077bd3a63586a6c9b6b1cdce Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 5 Dec 2025 15:50:44 -0700 Subject: [PATCH 09/22] Add MYNN surface layer scheme to convection-permitting suite --- test/test_suites/suite_convection_permitting.xml | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/test/test_suites/suite_convection_permitting.xml b/test/test_suites/suite_convection_permitting.xml index 4c91fda1..3b7adb3e 100644 --- a/test/test_suites/suite_convection_permitting.xml +++ b/test/test_suites/suite_convection_permitting.xml @@ -4,6 +4,8 @@ calc_exner compute_characteristic_grid_length_scale + compute_atmosphere_layer_thickness + compute_air_density geopotential_height_wrt_sfc_at_if_to_msl geopotential_height_wrt_sfc_to_msl temp_to_potential_temp @@ -13,6 +15,10 @@ mmm_physics_compat + sf_mynn_compat_pre + sf_mynn_compat + sf_mynn_diagnostics + bl_gwdo_compat_pre bl_gwdo_compat bl_gwdo_diagnostics From 0e2e7d091cb5c784851045c3ab13e56ab0ed2fce Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Thu, 22 Jan 2026 11:26:26 -0700 Subject: [PATCH 10/22] Remove unused variable --- schemes/mmm/sf_mynn_compat.F90 | 1 - 1 file changed, 1 deletion(-) diff --git a/schemes/mmm/sf_mynn_compat.F90 b/schemes/mmm/sf_mynn_compat.F90 index a1c60208..e365051d 100644 --- a/schemes/mmm/sf_mynn_compat.F90 +++ b/schemes/mmm/sf_mynn_compat.F90 @@ -27,7 +27,6 @@ subroutine sf_mynn_compat_pre_run( & spp_pbl, rstoch1d, & errmsg, errflg) use ccpp_kinds, only: kind_phys - use ccpp_scheme_utils, only: ccpp_constituent_index ! This threshold is hardcoded to the same value as in MMM physics. ! It is named `xice_threshold` there. From 3b0af3dcec7181edb22520bee3ca826131ad5e1c Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Thu, 22 Jan 2026 11:29:30 -0700 Subject: [PATCH 11/22] Move parameter to module level to avoid multiple declarations --- schemes/mmm/sf_mynn_compat.F90 | 13 ++++++------- 1 file changed, 6 insertions(+), 7 deletions(-) diff --git a/schemes/mmm/sf_mynn_compat.F90 b/schemes/mmm/sf_mynn_compat.F90 index e365051d..4a2cea40 100644 --- a/schemes/mmm/sf_mynn_compat.F90 +++ b/schemes/mmm/sf_mynn_compat.F90 @@ -1,6 +1,8 @@ !> This module contains interstitial schemes that are specific to MYNN surface layer scheme, !> which is part of MMM physics. module sf_mynn_compat + use ccpp_kinds, only: kind_phys + implicit none private @@ -9,6 +11,10 @@ module sf_mynn_compat public :: sf_mynn_compat_run public :: sf_mynn_diagnostics_init public :: sf_mynn_diagnostics_run + + ! This threshold is hardcoded to the same value as in MMM physics. + ! It is named `xice_threshold` there. + real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys contains !> \section arg_table_sf_mynn_compat_pre_run Argument Table !! \htmlinclude sf_mynn_compat_pre_run.html @@ -28,10 +34,6 @@ subroutine sf_mynn_compat_pre_run( & errmsg, errflg) use ccpp_kinds, only: kind_phys - ! This threshold is hardcoded to the same value as in MMM physics. - ! It is named `xice_threshold` there. - real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys - integer, intent(in) :: itimestep real(kind_phys), intent(in) :: u(:, :), v(:, :), t(:, :), qv(:, :), p(:, :), dz(:, :), rho(:, :), & icefrac(:), landfrac(:), snowhice(:), snowhland(:) @@ -163,9 +165,6 @@ subroutine sf_mynn_compat_run( & use ccpp_scheme_utils, only: ccpp_constituent_index use sf_mynn, only: sf_mynn_run - ! This threshold is hardcoded to the same value as in MMM physics. - ! It is named `xice_threshold` there. - real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys ! Typical threshold between sea surface temperature and ice surface temperature. See ! Algorithm Theoretical Basis Document (ATBD) for the MODIS Snow and Sea Ice-Mapping Algorithms, ! Section 4.4.3 Ice Surface Temperature (IST) Algorithm. From c2959dfcd0197ba7803c95878cba1c6867d5da96 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Thu, 22 Jan 2026 11:35:05 -0700 Subject: [PATCH 12/22] Add more checks for error flags --- schemes/mmm/sf_mynn_compat.F90 | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/schemes/mmm/sf_mynn_compat.F90 b/schemes/mmm/sf_mynn_compat.F90 index 4a2cea40..b11eaf63 100644 --- a/schemes/mmm/sf_mynn_compat.F90 +++ b/schemes/mmm/sf_mynn_compat.F90 @@ -135,6 +135,10 @@ subroutine sf_mynn_compat_init( & ! Precompute lookup tables in MYNN surface layer scheme. call sf_mynn_init(errmsg, errflg) + if (errflg /= 0) then + return + end if + ! MYNN surface layer scheme takes time averages of these variables internally. ! As a result, they must be able to persist across time steps. ust(:) = 0.0_kind_phys @@ -352,6 +356,10 @@ subroutine sf_mynn_compat_run( & iz0tlnd, its, ite, & errmsg, errflg) + if (errflg /= 0) then + return + end if + ! The first pass treated sea ice cells as land cells due to how `xland` is defined. ! The second pass treats sea ice cells as water cells instead. ! Then, for sea ice cells only, the final results are weighted averages according to their area fractions. @@ -370,6 +378,10 @@ subroutine sf_mynn_compat_run( & iz0tlnd, its, ite, & errmsg, errflg) + if (errflg /= 0) then + return + end if + ! Blend the final results between sea ice and sea water. where (mask_sea_ice_cell) chs = chs * icefrac + (1.0_kind_phys - icefrac) * chs_sea From a32c862031f5d5a2646da4420c3a05bd61c274a6 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Thu, 22 Jan 2026 11:46:10 -0700 Subject: [PATCH 13/22] Update standard name to be more clear However, this standard name change deviates from the ESM Standard Name Library. In the ESM Standard Name Library, "ratio_of_height_to_monin_obukhov_length" is defined. However, here "ratio_of_height_at_surface_adjacent_layer_to_monin_obukhov_length" is used. --- schemes/mmm/sf_mynn_compat.meta | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/schemes/mmm/sf_mynn_compat.meta b/schemes/mmm/sf_mynn_compat.meta index 09c0620f..408514f7 100644 --- a/schemes/mmm/sf_mynn_compat.meta +++ b/schemes/mmm/sf_mynn_compat.meta @@ -180,7 +180,7 @@ dimensions = (horizontal_loop_extent) intent = out [ zol ] - standard_name = ratio_of_height_to_monin_obukhov_length + standard_name = ratio_of_height_at_surface_adjacent_layer_to_monin_obukhov_length units = 1 type = real | kind = kind_phys dimensions = (horizontal_loop_extent) @@ -541,7 +541,7 @@ dimensions = (horizontal_loop_extent) intent = in [ zol ] - standard_name = ratio_of_height_to_monin_obukhov_length + standard_name = ratio_of_height_at_surface_adjacent_layer_to_monin_obukhov_length units = 1 type = real | kind = kind_phys dimensions = (horizontal_loop_extent) From 9e86e5133c0ce797378dc5cbdf892b78a0fd6d22 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Thu, 22 Jan 2026 12:07:44 -0700 Subject: [PATCH 14/22] Add comments to suite definition file --- test/test_suites/suite_convection_permitting.xml | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/test/test_suites/suite_convection_permitting.xml b/test/test_suites/suite_convection_permitting.xml index 3b7adb3e..6a3dac97 100644 --- a/test/test_suites/suite_convection_permitting.xml +++ b/test/test_suites/suite_convection_permitting.xml @@ -2,6 +2,7 @@ + calc_exner compute_characteristic_grid_length_scale compute_atmosphere_layer_thickness @@ -15,18 +16,22 @@ mmm_physics_compat + sf_mynn_compat_pre sf_mynn_compat sf_mynn_diagnostics + bl_gwdo_compat_pre bl_gwdo_compat bl_gwdo_diagnostics + cu_ntiedtke_compat_pre cu_ntiedtke_compat cu_ntiedtke_diagnostics + mmm_physics_accumulate_tendencies apply_tendency_of_eastward_wind apply_tendency_of_northward_wind From 1baf0397f9e13ae06ad3b90a898bf15ee66853c1 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Thu, 22 Jan 2026 12:53:56 -0700 Subject: [PATCH 15/22] Move computational procedures that are general enough to state converters * To conform to the existing naming convention in state converters, `compute_ ...` has been renamed to `calc_ ...`. * Fix incorrect units in a code comment in state converters. --- schemes/mmm/mmm_physics_compat.F90 | 49 ---------- schemes/mmm/mmm_physics_compat.meta | 90 ------------------- schemes/utilities/state_converters.F90 | 57 +++++++++++- schemes/utilities/state_converters.meta | 88 ++++++++++++++++++ .../suite_convection_permitting.xml | 4 +- .../tests/mmm/mmm_physics_compat_tests.pf | 60 ------------- .../tests/utilities/test_state_converters.pf | 60 +++++++++++++ 7 files changed, 205 insertions(+), 203 deletions(-) diff --git a/schemes/mmm/mmm_physics_compat.F90 b/schemes/mmm/mmm_physics_compat.F90 index 157af4e6..0a8d838d 100644 --- a/schemes/mmm/mmm_physics_compat.F90 +++ b/schemes/mmm/mmm_physics_compat.F90 @@ -9,8 +9,6 @@ module mmm_physics_compat public :: mmm_physics_accumulate_tendencies_run public :: mmm_physics_persist_states_init public :: mmm_physics_persist_states_timestep_final - public :: compute_air_density_run - public :: compute_atmosphere_layer_thickness_run public :: compute_characteristic_grid_length_scale_init public :: geopotential_height_wrt_sfc_at_if_to_msl_run public :: geopotential_height_wrt_sfc_to_msl_run @@ -233,53 +231,6 @@ pure subroutine mmm_physics_persist_states_timestep_final( & qv_prev(:, :) = qv_curr(:, :) end subroutine mmm_physics_persist_states_timestep_final - !> \section arg_table_compute_air_density_run Argument Table - !! \htmlinclude compute_air_density_run.html - pure subroutine compute_air_density_run( & - pdel, gravit, dz, & - rho, & - errmsg, errflg) - use ccpp_kinds, only: kind_phys - - real(kind_phys), intent(in) :: pdel(:, :), gravit, dz(:, :) - real(kind_phys), intent(out) :: rho(:, :) - character(*), intent(out) :: errmsg - integer, intent(out) :: errflg - - errmsg = '' - errflg = 0 - - ! Compute air density by hydrostatic equation. - rho(:, :) = pdel(:, :) / (gravit * dz(:, :)) - end subroutine compute_air_density_run - - !> \section arg_table_compute_atmosphere_layer_thickness_run Argument Table - !! \htmlinclude compute_atmosphere_layer_thickness_run.html - pure subroutine compute_atmosphere_layer_thickness_run( & - ncol, & - zisfc, & - dz, & - errmsg, errflg) - use ccpp_kinds, only: kind_phys - - integer, intent(in) :: ncol - real(kind_phys), intent(in) :: zisfc(:, :) - real(kind_phys), intent(out) :: dz(:, :) - character(*), intent(out) :: errmsg - integer, intent(out) :: errflg - - integer :: i - - errmsg = '' - errflg = 0 - - ! In CAM-SIMA, the first vertical index is at top of atmosphere. - ! The last one is at bottom of atmosphere. The resulting `dz` is positive. - do i = 1, ncol - dz(i, :) = zisfc(i, 1:size(zisfc, 2) - 1) - zisfc(i, 2:size(zisfc, 2)) - end do - end subroutine compute_atmosphere_layer_thickness_run - !> \section arg_table_compute_characteristic_grid_length_scale_init Argument Table !! \htmlinclude compute_characteristic_grid_length_scale_init.html pure subroutine compute_characteristic_grid_length_scale_init( & diff --git a/schemes/mmm/mmm_physics_compat.meta b/schemes/mmm/mmm_physics_compat.meta index 407e6ce5..1cb1adcf 100644 --- a/schemes/mmm/mmm_physics_compat.meta +++ b/schemes/mmm/mmm_physics_compat.meta @@ -466,96 +466,6 @@ # ---------- -[ccpp-table-properties] - name = compute_air_density - type = scheme - -[ccpp-arg-table] - name = compute_air_density_run - type = scheme -[ pdel ] - standard_name = air_pressure_thickness - units = Pa - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = in -[ gravit ] - standard_name = standard_gravitational_acceleration - units = m s-2 - type = real | kind = kind_phys - dimensions = () - intent = in -[ dz ] - standard_name = atmosphere_layer_thickness - units = m - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = in -[ rho ] - standard_name = air_density - units = kg m-3 - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = out -[ errmsg ] - standard_name = ccpp_error_message - long_name = Error message for error handling in CCPP - units = none - type = character | kind = len=* - dimensions = () - intent = out -[ errflg ] - standard_name = ccpp_error_code - long_name = Error flag for error handling in CCPP - units = 1 - type = integer - dimensions = () - intent = out - -# ---------- - -[ccpp-table-properties] - name = compute_atmosphere_layer_thickness - type = scheme - -[ccpp-arg-table] - name = compute_atmosphere_layer_thickness_run - type = scheme -[ ncol ] - standard_name = horizontal_loop_extent - units = count - type = integer - dimensions = () - intent = in -[ zisfc ] - standard_name = geopotential_height_wrt_surface_at_interface - units = m - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent, vertical_interface_dimension) - intent = in -[ dz ] - standard_name = atmosphere_layer_thickness - units = m - type = real | kind = kind_phys - dimensions = (horizontal_loop_extent, vertical_layer_dimension) - intent = out -[ errmsg ] - standard_name = ccpp_error_message - long_name = Error message for error handling in CCPP - units = none - type = character | kind = len=* - dimensions = () - intent = out -[ errflg ] - standard_name = ccpp_error_code - long_name = Error flag for error handling in CCPP - units = 1 - type = integer - dimensions = () - intent = out - -# ---------- - [ccpp-table-properties] name = compute_characteristic_grid_length_scale type = scheme diff --git a/schemes/utilities/state_converters.F90 b/schemes/utilities/state_converters.F90 index 5d49ec97..a4b5c578 100644 --- a/schemes/utilities/state_converters.F90 +++ b/schemes/utilities/state_converters.F90 @@ -10,9 +10,15 @@ module state_converters public :: temp_to_potential_temp_run public :: potential_temp_to_temp_run - ! Calculate density from equation of state/ideal gas law + ! Calculate dry air density by equation of state/ideal gas law public :: calc_dry_air_ideal_gas_density_run + ! Calculate air density by hydrostatic equation + public :: calc_hydrostatic_air_density_run + + ! Calculate atmosphere layer thickness + public :: calc_atmosphere_layer_thickness_run + ! Calculate exner public :: calc_exner_run @@ -75,7 +81,7 @@ end subroutine potential_temp_to_temp_run subroutine calc_dry_air_ideal_gas_density_run(ncol, nz, rair, pmiddry, temp, rho, errmsg, errflg) integer, intent(in) :: ncol ! Number of columns integer, intent(in) :: nz ! Number of vertical levels - real(kind_phys), intent(in) :: rair(:,:) ! gas constant for dry air (J kg-1) + real(kind_phys), intent(in) :: rair(:,:) ! Gas constant of dry air (J kg-1 K-1) real(kind_phys), intent(in) :: pmiddry(:,:) ! Air pressure of dry air (Pa) real(kind_phys), intent(in) :: temp(:,:) ! Air temperature (K) real(kind_phys), intent(out) :: rho(:,:) ! Dry air density (kg m-3) @@ -93,6 +99,53 @@ subroutine calc_dry_air_ideal_gas_density_run(ncol, nz, rair, pmiddry, temp, rho end subroutine calc_dry_air_ideal_gas_density_run + !> \section arg_table_calc_hydrostatic_air_density_run Argument Table + !! \htmlinclude calc_hydrostatic_air_density_run.html + pure subroutine calc_hydrostatic_air_density_run( & + pdel, gravit, dz, & + rho, & + errmsg, errflg) + use ccpp_kinds, only: kind_phys + + real(kind_phys), intent(in) :: pdel(:, :), gravit, dz(:, :) + real(kind_phys), intent(out) :: rho(:, :) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + ! Calculate air density by hydrostatic equation. + rho(:, :) = pdel(:, :) / (gravit * dz(:, :)) + + errmsg = '' + errflg = 0 + end subroutine calc_hydrostatic_air_density_run + + !> \section arg_table_calc_atmosphere_layer_thickness_run Argument Table + !! \htmlinclude calc_atmosphere_layer_thickness_run.html + pure subroutine calc_atmosphere_layer_thickness_run( & + ncol, & + zisfc, & + dz, & + errmsg, errflg) + use ccpp_kinds, only: kind_phys + + integer, intent(in) :: ncol + real(kind_phys), intent(in) :: zisfc(:, :) + real(kind_phys), intent(out) :: dz(:, :) + character(*), intent(out) :: errmsg + integer, intent(out) :: errflg + + integer :: i + + ! In CAM-SIMA, the first vertical index is at top of atmosphere. + ! The last one is at bottom of atmosphere. The resulting `dz` is positive. + do i = 1, ncol + dz(i, :) = zisfc(i, 1:size(zisfc, 2) - 1) - zisfc(i, 2:size(zisfc, 2)) + end do + + errmsg = '' + errflg = 0 + end subroutine calc_atmosphere_layer_thickness_run + !> \section arg_table_calc_exner_run Argument Table !! \htmlinclude calc_exner_run.html subroutine calc_exner_run(ncol, nz, cpair, rair, ref_pres, pmid, exner, & diff --git a/schemes/utilities/state_converters.meta b/schemes/utilities/state_converters.meta index 1bbfc5bb..9389d8f5 100644 --- a/schemes/utilities/state_converters.meta +++ b/schemes/utilities/state_converters.meta @@ -173,6 +173,94 @@ type = integer intent = out +######################################################### +[ccpp-table-properties] + name = calc_hydrostatic_air_density + type = scheme + +[ccpp-arg-table] + name = calc_hydrostatic_air_density_run + type = scheme +[ pdel ] + standard_name = air_pressure_thickness + units = Pa + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ gravit ] + standard_name = standard_gravitational_acceleration + units = m s-2 + type = real | kind = kind_phys + dimensions = () + intent = in +[ dz ] + standard_name = atmosphere_layer_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = in +[ rho ] + standard_name = air_density + units = kg m-3 + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + +######################################################### +[ccpp-table-properties] + name = calc_atmosphere_layer_thickness + type = scheme + +[ccpp-arg-table] + name = calc_atmosphere_layer_thickness_run + type = scheme +[ ncol ] + standard_name = horizontal_loop_extent + units = count + type = integer + dimensions = () + intent = in +[ zisfc ] + standard_name = geopotential_height_wrt_surface_at_interface + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_interface_dimension) + intent = in +[ dz ] + standard_name = atmosphere_layer_thickness + units = m + type = real | kind = kind_phys + dimensions = (horizontal_loop_extent, vertical_layer_dimension) + intent = out +[ errmsg ] + standard_name = ccpp_error_message + long_name = Error message for error handling in CCPP + units = none + type = character | kind = len=* + dimensions = () + intent = out +[ errflg ] + standard_name = ccpp_error_code + long_name = Error flag for error handling in CCPP + units = 1 + type = integer + dimensions = () + intent = out + ######################################################### [ccpp-table-properties] name = calc_exner diff --git a/test/test_suites/suite_convection_permitting.xml b/test/test_suites/suite_convection_permitting.xml index 6a3dac97..884a0c1f 100644 --- a/test/test_suites/suite_convection_permitting.xml +++ b/test/test_suites/suite_convection_permitting.xml @@ -3,10 +3,10 @@ + calc_atmosphere_layer_thickness calc_exner + calc_hydrostatic_air_density compute_characteristic_grid_length_scale - compute_atmosphere_layer_thickness - compute_air_density geopotential_height_wrt_sfc_at_if_to_msl geopotential_height_wrt_sfc_to_msl temp_to_potential_temp diff --git a/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf b/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf index 69c88a5b..e79e61a4 100644 --- a/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf +++ b/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf @@ -442,66 +442,6 @@ subroutine test_mmm_physics_persist_states_timestep_final() @assertEqual(0, errflg) end subroutine test_mmm_physics_persist_states_timestep_final -@test -subroutine test_compute_air_density_run() - use ccpp_kinds, only: kind_phys - use funit - use mmm_physics_compat, only: compute_air_density_run - - integer, parameter :: ncol = 100, pver = 10 - real(kind_phys) :: pdel(ncol, pver), gravit, dz(ncol, pver) - real(kind_phys) :: rho(ncol, pver) - character(100) :: errmsg - integer :: errflg - - pdel(:, :) = 98.0_kind_phys - gravit = 9.8_kind_phys - dz(:, :) = 10.0_kind_phys - rho(:, :) = huge(0.0_kind_phys) - - call compute_air_density_run( & - pdel, gravit, dz, & - rho, & - errmsg, errflg) - - ! Should compute air density correctly. - @assertEqual(1.0_kind_phys, rho, epsilon(0.0_kind_phys)) - @assertEqual('', errmsg) - @assertEqual(0, errflg) -end subroutine test_compute_air_density_run - -@test -subroutine test_compute_atmosphere_layer_thickness_run() - use ccpp_kinds, only: kind_phys - use funit - use mmm_physics_compat, only: compute_atmosphere_layer_thickness_run - - integer, parameter :: ncol = 100, pver = 10, pverp = 11 - real(kind_phys) :: zisfc(ncol, pverp) - real(kind_phys) :: dz(ncol, pver) - character(100) :: errmsg - integer :: errflg - - integer :: i - - do i = 1, pverp - zisfc(:, i) = 10000.0_kind_phys - real(i - 1, kind_phys) * 1000.0_kind_phys - end do - - dz(:, :) = huge(0.0_kind_phys) - - call compute_atmosphere_layer_thickness_run( & - ncol, & - zisfc, & - dz, & - errmsg, errflg) - - ! Should compute atmosphere layer thickness correctly. - @assertEqual(1000.0_kind_phys, dz, epsilon(0.0_kind_phys)) - @assertEqual('', errmsg) - @assertEqual(0, errflg) -end subroutine test_compute_atmosphere_layer_thickness_run - @test subroutine test_compute_characteristic_grid_length_scale_init() use ccpp_kinds, only: kind_phys diff --git a/test/unit-test/tests/utilities/test_state_converters.pf b/test/unit-test/tests/utilities/test_state_converters.pf index bc56f038..dcd7bc8e 100644 --- a/test/unit-test/tests/utilities/test_state_converters.pf +++ b/test/unit-test/tests/utilities/test_state_converters.pf @@ -25,3 +25,63 @@ subroutine test_temp_to_potential_temp() @assertEqual(0, errflg) end subroutine test_temp_to_potential_temp + +@test +subroutine test_calc_hydrostatic_air_density_run() + use ccpp_kinds, only: kind_phys + use funit + use state_converters, only: calc_hydrostatic_air_density_run + + integer, parameter :: ncol = 100, pver = 10 + real(kind_phys) :: pdel(ncol, pver), gravit, dz(ncol, pver) + real(kind_phys) :: rho(ncol, pver) + character(100) :: errmsg + integer :: errflg + + pdel(:, :) = 98.0_kind_phys + gravit = 9.8_kind_phys + dz(:, :) = 10.0_kind_phys + rho(:, :) = huge(0.0_kind_phys) + + call calc_hydrostatic_air_density_run( & + pdel, gravit, dz, & + rho, & + errmsg, errflg) + + ! Should compute air density correctly. + @assertEqual(1.0_kind_phys, rho, spacing(1.0_kind_phys)) + @assertEqual('', errmsg) + @assertEqual(0, errflg) +end subroutine test_calc_hydrostatic_air_density_run + +@test +subroutine test_calc_atmosphere_layer_thickness_run() + use ccpp_kinds, only: kind_phys + use funit + use state_converters, only: calc_atmosphere_layer_thickness_run + + integer, parameter :: ncol = 100, pver = 10, pverp = 11 + real(kind_phys) :: zisfc(ncol, pverp) + real(kind_phys) :: dz(ncol, pver) + character(100) :: errmsg + integer :: errflg + + integer :: i + + do i = 1, pverp + zisfc(:, i) = 10000.0_kind_phys - real(i - 1, kind_phys) * 1000.0_kind_phys + end do + + dz(:, :) = huge(0.0_kind_phys) + + call calc_atmosphere_layer_thickness_run( & + ncol, & + zisfc, & + dz, & + errmsg, errflg) + + ! Should compute atmosphere layer thickness correctly. + @assertEqual(1000.0_kind_phys, dz, spacing(1000.0_kind_phys)) + @assertEqual('', errmsg) + @assertEqual(0, errflg) +end subroutine test_calc_atmosphere_layer_thickness_run From e340a7cd6ea336d67fb7cc7afd28e4a008f5a694 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Thu, 22 Jan 2026 13:51:55 -0700 Subject: [PATCH 16/22] Indicate the originating procedure name in error messages --- schemes/mmm/cu_ntiedtke_compat.F90 | 2 +- schemes/mmm/sf_mynn_compat.F90 | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/schemes/mmm/cu_ntiedtke_compat.F90 b/schemes/mmm/cu_ntiedtke_compat.F90 index bb18a5a1..8df9f01a 100644 --- a/schemes/mmm/cu_ntiedtke_compat.F90 +++ b/schemes/mmm/cu_ntiedtke_compat.F90 @@ -47,7 +47,7 @@ subroutine cu_ntiedtke_compat_pre_run( & if (errflg /= 0 .or. & water_vapor_mixing_ratio_index < lbound(cflx, 2) .or. water_vapor_mixing_ratio_index > ubound(cflx, 2)) then - errmsg = 'Failed to find desired constituent flux from cflx' + errmsg = 'cu_ntiedtke_compat_pre_run: Failed to find desired constituent flux from cflx' errflg = 1 return diff --git a/schemes/mmm/sf_mynn_compat.F90 b/schemes/mmm/sf_mynn_compat.F90 index b11eaf63..613ef09d 100644 --- a/schemes/mmm/sf_mynn_compat.F90 +++ b/schemes/mmm/sf_mynn_compat.F90 @@ -238,7 +238,7 @@ subroutine sf_mynn_compat_run( & if (errflg /= 0 .or. & water_vapor_mixing_ratio_index < lbound(cflx, 2) .or. water_vapor_mixing_ratio_index > ubound(cflx, 2)) then - errmsg = 'Failed to find desired constituent flux from cflx' + errmsg = 'sf_mynn_compat_run: Failed to find desired constituent flux from cflx' errflg = 1 return From d98721370165a89c9f378337a0142d17fac06a4f Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 6 Feb 2026 12:39:39 -0700 Subject: [PATCH 17/22] Pin mmm-physics repository by tag instead of hash --- .gitmodules | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.gitmodules b/.gitmodules index 41de1b86..10431cb9 100644 --- a/.gitmodules +++ b/.gitmodules @@ -1,7 +1,7 @@ [submodule "mmm-physics"] path = schemes/mmm/mmm_physics url = https://github.com/NCAR/MMM-physics.git - fxtag = 6cb29bc8f17bc3efdfaaa8a268ca53e2504aadbd + fxtag = 20251208-CAM-SIMA fxrequired = AlwaysRequired fxDONOTUSEurl = https://github.com/NCAR/MMM-physics.git [submodule "pumas"] From fe2728db747aa4ff0e7c52304286a9f5402b81e0 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 6 Feb 2026 13:10:50 -0700 Subject: [PATCH 18/22] Centralize option definitions into one place The definitions of `isfflx`, `isftcflx`, `iz0tlnd`, `spp_pbl`, and `xice_threshold` have been moved to the `mmm_physics_compat_init` subroutine. This avoids duplicate definitions in multiple places. Also adjust misleading code comments. --- schemes/mmm/mmm_physics_compat.F90 | 20 +++++++++++-------- schemes/mmm/mmm_physics_compat.meta | 18 +++++++++++++++++ schemes/mmm/sf_mynn_compat.F90 | 31 ++++++++++++++++------------- schemes/mmm/sf_mynn_compat.meta | 24 ++++++++++++++++------ 4 files changed, 65 insertions(+), 28 deletions(-) diff --git a/schemes/mmm/mmm_physics_compat.F90 b/schemes/mmm/mmm_physics_compat.F90 index 0a8d838d..8ab26573 100644 --- a/schemes/mmm/mmm_physics_compat.F90 +++ b/schemes/mmm/mmm_physics_compat.F90 @@ -17,19 +17,27 @@ module mmm_physics_compat !! \htmlinclude mmm_physics_compat_init.html pure subroutine mmm_physics_compat_init( & isfflx, isftcflx, iz0tlnd, & + spp_pbl, & + xice_threshold, & errmsg, errflg) + use ccpp_kinds, only: kind_phys integer, intent(out) :: isfflx, isftcflx, iz0tlnd + logical, intent(out) :: spp_pbl + real(kind_phys), intent(out) :: xice_threshold character(*), intent(out) :: errmsg integer, intent(out) :: errflg errmsg = '' errflg = 0 - ! Set options that are specific to MMM physics at model initialization. + ! These options are hardcoded to the same values as in the MPAS physics driver. + ! There are other possible values for them in WRF, but the following combination is the only supported one in MPAS. isfflx = 1 isftcflx = 0 iz0tlnd = 0 + spp_pbl = .false. + xice_threshold = 0.02_kind_phys end subroutine mmm_physics_compat_init !> \section arg_table_mmm_physics_compat_run Argument Table @@ -38,21 +46,17 @@ pure subroutine mmm_physics_compat_run( & nstep, & dt, & theta_curr, theta_prev, qv_curr, qv_prev, & - icefrac, landfrac, & + icefrac, xice_threshold, landfrac, & scheme_name, & rthdynten, rqvdynten, & xland, & errmsg, errflg) use ccpp_kinds, only: kind_phys - ! This threshold is hardcoded to the same value as in MMM physics. - ! It is named `xice_threshold` there. - real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys - integer, intent(in) :: nstep real(kind_phys), intent(in) :: dt, & theta_curr(:, :), theta_prev(:, :), qv_curr(:, :), qv_prev(:, :), & - icefrac(:), landfrac(:) + icefrac(:), xice_threshold, landfrac(:) character(256), intent(out) :: scheme_name real(kind_phys), intent(out) :: rthdynten(:, :), rqvdynten(:, :), & xland(:) @@ -76,7 +80,7 @@ pure subroutine mmm_physics_compat_run( & ! * xland = 1.0 for land cells, including sea ice cells. ! * xland = 2.0 for water cells. where (landfrac >= 0.5_kind_phys .or. & - icefrac >= sea_ice_area_fraction_threshold) + icefrac >= xice_threshold) xland = 1.0_kind_phys elsewhere xland = 2.0_kind_phys diff --git a/schemes/mmm/mmm_physics_compat.meta b/schemes/mmm/mmm_physics_compat.meta index 1cb1adcf..39195953 100644 --- a/schemes/mmm/mmm_physics_compat.meta +++ b/schemes/mmm/mmm_physics_compat.meta @@ -23,6 +23,18 @@ type = integer dimensions = () intent = out +[ spp_pbl ] + standard_name = flag_for_stochastically_perturbed_parameterization_for_mynn_scheme + units = flag + type = logical + dimensions = () + intent = out +[ xice_threshold ] + standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = out [ errmsg ] standard_name = ccpp_error_message long_name = Error message for error handling in CCPP @@ -83,6 +95,12 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = in +[ xice_threshold ] + standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in [ landfrac ] standard_name = land_area_fraction_from_coupler units = fraction diff --git a/schemes/mmm/sf_mynn_compat.F90 b/schemes/mmm/sf_mynn_compat.F90 index 613ef09d..3566405e 100644 --- a/schemes/mmm/sf_mynn_compat.F90 +++ b/schemes/mmm/sf_mynn_compat.F90 @@ -11,17 +11,14 @@ module sf_mynn_compat public :: sf_mynn_compat_run public :: sf_mynn_diagnostics_init public :: sf_mynn_diagnostics_run - - ! This threshold is hardcoded to the same value as in MMM physics. - ! It is named `xice_threshold` there. - real(kind_phys), parameter :: sea_ice_area_fraction_threshold = 0.02_kind_phys contains !> \section arg_table_sf_mynn_compat_pre_run Argument Table !! \htmlinclude sf_mynn_compat_pre_run.html subroutine sf_mynn_compat_pre_run( & itimestep, & + spp_pbl, & u, v, t, qv, p, dz, rho, & - icefrac, landfrac, snowhice, snowhland, & + icefrac, xice_threshold, landfrac, snowhice, snowhland, & u1d, v1d, t1d, qv1d, p1d, dz8w1d, rho1d, & u1d2, v1d2, dz2w1d, & chs, chs2, cqs2, cpm, rmol, & @@ -30,14 +27,14 @@ subroutine sf_mynn_compat_pre_run( & flhc, flqc, snowh, qgh, qsfc, & gz1oz0, wspd, br, svp1, svp2, & svp3, svpt0, qcg, & - spp_pbl, rstoch1d, & + rstoch1d, & errmsg, errflg) use ccpp_kinds, only: kind_phys integer, intent(in) :: itimestep + logical, intent(in) :: spp_pbl real(kind_phys), intent(in) :: u(:, :), v(:, :), t(:, :), qv(:, :), p(:, :), dz(:, :), rho(:, :), & - icefrac(:), landfrac(:), snowhice(:), snowhland(:) - logical, intent(out) :: spp_pbl + icefrac(:), xice_threshold, landfrac(:), snowhice(:), snowhland(:) real(kind_phys), intent(out) :: u1d(:), v1d(:), t1d(:), qv1d(:), p1d(:), dz8w1d(:), rho1d(:), & u1d2(:), v1d2(:), dz2w1d(:), & chs(:), chs2(:), cqs2(:), cpm(:), rmol(:), & @@ -94,7 +91,7 @@ subroutine sf_mynn_compat_pre_run( & snowh = snowhland end where - where (icefrac >= sea_ice_area_fraction_threshold) + where (icefrac >= xice_threshold) snowh = snowhice end where @@ -113,8 +110,14 @@ subroutine sf_mynn_compat_pre_run( & qcg(:) = 0.0_kind_phys ! Not used but still appear in the argument list... - spp_pbl = .false. - rstoch1d(:) = 0.0_kind_phys + if (spp_pbl) then + errmsg = 'sf_mynn_compat_pre_run: Stochastically perturbed parameterization is not supported' + errflg = 1 + + return + else + rstoch1d(:) = 0.0_kind_phys + end if errmsg = '' errflg = 0 @@ -152,7 +155,7 @@ end subroutine sf_mynn_compat_init !> \section arg_table_sf_mynn_compat_run Argument Table !! \htmlinclude sf_mynn_compat_run.html subroutine sf_mynn_compat_run( & - ncol, cflx, icefrac, sst, & + ncol, cflx, icefrac, xice_threshold, sst, & u1d, v1d, t1d, qv1d, p1d, dz8w1d, rho1d, & u1d2, v1d2, dz2w1d, cp, g, rovcp, r, xlv, & psfcpa, chs, chs2, cqs2, cpm, pblh, rmol, & @@ -180,7 +183,7 @@ subroutine sf_mynn_compat_run( & isftcflx, & iz0tlnd, its, ite logical, intent(in) :: spp_pbl - real(kind_phys), intent(in) :: icefrac(:), sst(:), & + real(kind_phys), intent(in) :: icefrac(:), xice_threshold, sst(:), & u1d(:), v1d(:), t1d(:), qv1d(:), p1d(:), dz8w1d(:), rho1d(:), & u1d2(:), v1d2(:), dz2w1d(:), cp, g, rovcp, r, xlv, & psfcpa(:), pblh(:), & @@ -252,7 +255,7 @@ subroutine sf_mynn_compat_run( & return end if - mask_sea_ice_cell(:) = (icefrac >= sea_ice_area_fraction_threshold) + mask_sea_ice_cell(:) = (icefrac >= xice_threshold) ! If there are sea ice cells, make local copies of the variables in preparation for the second pass. if (any(mask_sea_ice_cell)) then diff --git a/schemes/mmm/sf_mynn_compat.meta b/schemes/mmm/sf_mynn_compat.meta index 408514f7..380716bb 100644 --- a/schemes/mmm/sf_mynn_compat.meta +++ b/schemes/mmm/sf_mynn_compat.meta @@ -11,6 +11,12 @@ type = integer dimensions = () intent = in +[ spp_pbl ] + standard_name = flag_for_stochastically_perturbed_parameterization_for_mynn_scheme + units = flag + type = logical + dimensions = () + intent = in [ u ] standard_name = eastward_wind units = m s-1 @@ -59,6 +65,12 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = in +[ xice_threshold ] + standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in [ landfrac ] standard_name = land_area_fraction_from_coupler units = fraction @@ -293,12 +305,6 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = out -[ spp_pbl ] - standard_name = flag_for_stochastically_perturbed_parameterization_for_mynn_scheme - units = flag - type = logical - dimensions = () - intent = out [ rstoch1d ] standard_name = stochastically_perturbed_weights_for_mynn_scheme units = 1 @@ -384,6 +390,12 @@ type = real | kind = kind_phys dimensions = (horizontal_loop_extent) intent = in +[ xice_threshold ] + standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme + units = 1 + type = real | kind = kind_phys + dimensions = () + intent = in [ sst ] standard_name = sea_surface_temperature_from_coupler units = K From 1eac9d77726b68ab9da476fdac40460a975d306b Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 6 Feb 2026 13:11:31 -0700 Subject: [PATCH 19/22] Adjust unit tests for changes in arguments --- .../tests/mmm/mmm_physics_compat_tests.pf | 18 ++++++++++++++---- 1 file changed, 14 insertions(+), 4 deletions(-) diff --git a/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf b/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf index e79e61a4..d1e3e0af 100644 --- a/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf +++ b/test/unit-test/tests/mmm/mmm_physics_compat_tests.pf @@ -1,24 +1,33 @@ @test subroutine test_mmm_physics_compat_init() + use ccpp_kinds, only: kind_phys use funit use mmm_physics_compat, only: mmm_physics_compat_init integer :: isfflx, isftcflx, iz0tlnd + logical :: spp_pbl + real(kind_phys) :: xice_threshold character(100) :: errmsg integer :: errflg isfflx = huge(0) isftcflx = huge(0) iz0tlnd = huge(0) + spp_pbl = .false. + xice_threshold = huge(0.0_kind_phys) call mmm_physics_compat_init( & isfflx, isftcflx, iz0tlnd, & + spp_pbl, & + xice_threshold, & errmsg, errflg) - ! Should set options that are specific to MMM physics. + ! Should set options to the same values as in the MPAS physics driver. @assertEqual(1, isfflx) @assertEqual(0, isftcflx) @assertEqual(0, iz0tlnd) + @assertEqual(.false., spp_pbl) + @assertEqual(0.02_kind_phys, xice_threshold) @assertEqual('', errmsg) @assertEqual(0, errflg) end subroutine test_mmm_physics_compat_init @@ -33,7 +42,7 @@ subroutine test_mmm_physics_compat_run() integer :: nstep real(kind_phys) :: dt real(kind_phys) :: theta_curr(ncol, pver), theta_prev(ncol, pver), qv_curr(ncol, pver), qv_prev(ncol, pver) - real(kind_phys) :: icefrac(ncol), landfrac(ncol) + real(kind_phys) :: icefrac(ncol), xice_threshold, landfrac(ncol) character(256) :: scheme_name real(kind_phys) :: rthdynten(ncol, pver), rqvdynten(ncol, pver) real(kind_phys) :: xland(ncol) @@ -51,6 +60,7 @@ subroutine test_mmm_physics_compat_run() qv_prev(:, :) = 0.01_kind_phys icefrac(1:50) = 0.0_kind_phys icefrac(51:100) = [(0.01_kind_phys + real(i - 1, kind_phys) * 0.02_kind_phys, i = 1, 50)] + xice_threshold = 0.02_kind_phys landfrac(1:50) = [(0.01_kind_phys + real(i - 1, kind_phys) * 0.02_kind_phys, i = 1, 50)] landfrac(51:100) = 0.0_kind_phys scheme_name = '' @@ -62,7 +72,7 @@ subroutine test_mmm_physics_compat_run() nstep, & dt, & theta_curr, theta_prev, qv_curr, qv_prev, & - icefrac, landfrac, & + icefrac, xice_threshold, landfrac, & scheme_name, & rthdynten, rqvdynten, & xland, & @@ -92,7 +102,7 @@ subroutine test_mmm_physics_compat_run() nstep, & dt, & theta_curr, theta_prev, qv_curr, qv_prev, & - icefrac, landfrac, & + icefrac, xice_threshold, landfrac, & scheme_name, & rthdynten, rqvdynten, & xland, & From deeff2b7963bcfcc24a62a74208f8148581e5e91 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 6 Feb 2026 13:13:02 -0700 Subject: [PATCH 20/22] Prefix subroutine name in allocation error messages --- schemes/mmm/sf_mynn_compat.F90 | 15 +++++++++++++++ 1 file changed, 15 insertions(+) diff --git a/schemes/mmm/sf_mynn_compat.F90 b/schemes/mmm/sf_mynn_compat.F90 index 3566405e..80440a42 100644 --- a/schemes/mmm/sf_mynn_compat.F90 +++ b/schemes/mmm/sf_mynn_compat.F90 @@ -229,6 +229,9 @@ subroutine sf_mynn_compat_run( & allocate(ch(ncol), errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'sf_mynn_compat_run: Failed to allocate "ch"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if @@ -252,6 +255,9 @@ subroutine sf_mynn_compat_run( & allocate(mask_sea_ice_cell(ncol), errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'sf_mynn_compat_run: Failed to allocate "mask_sea_ice_cell"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if @@ -263,6 +269,9 @@ subroutine sf_mynn_compat_run( & errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'sf_mynn_compat_run: Failed to allocate "mavail_sea", "xland_sea", "tsk_sea"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if @@ -279,6 +288,9 @@ subroutine sf_mynn_compat_run( & errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'sf_mynn_compat_run: Failed to allocate "chs_sea", "chs2_sea", "cqs2_sea", ...' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if @@ -313,6 +325,9 @@ subroutine sf_mynn_compat_run( & errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'sf_mynn_compat_run: Failed to allocate "u10_sea", "v10_sea", "th2_sea", ...' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if From 9ed16fb8acd7ba30fff8dcb10d7e95079bd24468 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Fri, 6 Feb 2026 13:18:06 -0700 Subject: [PATCH 21/22] Prefix subroutine name in allocation error messages retrospectively --- schemes/mmm/cu_ntiedtke_compat.F90 | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/schemes/mmm/cu_ntiedtke_compat.F90 b/schemes/mmm/cu_ntiedtke_compat.F90 index 8df9f01a..e98a8e25 100644 --- a/schemes/mmm/cu_ntiedtke_compat.F90 +++ b/schemes/mmm/cu_ntiedtke_compat.F90 @@ -123,36 +123,54 @@ subroutine cu_ntiedtke_compat_run( & allocate(pu_local, source=pu, errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'cu_ntiedtke_compat_run: Failed to allocate "pu_local"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if allocate(pv_local, source=pv, errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'cu_ntiedtke_compat_run: Failed to allocate "pv_local"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if allocate(pt_local, source=pt, errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'cu_ntiedtke_compat_run: Failed to allocate "pt_local"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if allocate(pqv_local, source=pqv, errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'cu_ntiedtke_compat_run: Failed to allocate "pqv_local"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if allocate(pqc_local, source=pqc, errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'cu_ntiedtke_compat_run: Failed to allocate "pqc_local"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if allocate(pqi_local, source=pqi, errmsg=errmsg, stat=errflg) if (errflg /= 0) then + errmsg = 'cu_ntiedtke_compat_run: Failed to allocate "pqi_local"' // new_line('') // & + 'Allocation returned with error: ' // trim(adjustl(errmsg)) + return end if From 0872abb5a8d814684d435d7c6533769c24306811 Mon Sep 17 00:00:00 2001 From: Kuan-Chih Wang Date: Tue, 10 Feb 2026 12:31:12 -0700 Subject: [PATCH 22/22] The units for variables named "*_area_fraction" should be "fraction" --- schemes/mmm/mmm_physics_compat.meta | 4 ++-- schemes/mmm/sf_mynn_compat.meta | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/schemes/mmm/mmm_physics_compat.meta b/schemes/mmm/mmm_physics_compat.meta index 39195953..f6c74782 100644 --- a/schemes/mmm/mmm_physics_compat.meta +++ b/schemes/mmm/mmm_physics_compat.meta @@ -31,7 +31,7 @@ intent = out [ xice_threshold ] standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme - units = 1 + units = fraction type = real | kind = kind_phys dimensions = () intent = out @@ -97,7 +97,7 @@ intent = in [ xice_threshold ] standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme - units = 1 + units = fraction type = real | kind = kind_phys dimensions = () intent = in diff --git a/schemes/mmm/sf_mynn_compat.meta b/schemes/mmm/sf_mynn_compat.meta index 380716bb..53463d22 100644 --- a/schemes/mmm/sf_mynn_compat.meta +++ b/schemes/mmm/sf_mynn_compat.meta @@ -67,7 +67,7 @@ intent = in [ xice_threshold ] standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme - units = 1 + units = fraction type = real | kind = kind_phys dimensions = () intent = in @@ -392,7 +392,7 @@ intent = in [ xice_threshold ] standard_name = sea_ice_area_fraction_threshold_for_mmm_scheme - units = 1 + units = fraction type = real | kind = kind_phys dimensions = () intent = in