Skip to content

Commit

Permalink
This brings in several updates and fixes for the PAM CRM used in the …
Browse files Browse the repository at this point in the history 
…MMF2 compsets to improve stability and throughput. Some of these required updates to the PAM submodule, but the most important change was the addition of adaptive sub-cycling in the PAM driver. This change allowed the CRM time step to be increased to 10 seconds, which dramatically improved the throughput while also making the model much more resilient to numerically unstable conditions due to unrealistically strong velocities (still not sure why these occur in the first place).

List of notable changes for PAM and/or the driver:

adaptive subcycling in PAM driver
enable variance transport by default
disable mean-state acceleration for dry density
use crm_nx=45 with exact interpolation for PAM dycor (requires odd CRM columns)
increase PAM hyperdiffusion timescale to 60 seconds to reduce over-smoothing
Fixes E3SM-Project#6363

[non-BFB] only for MMF2

Merge branch 'whannah/mmf/pam-updates'
  • Loading branch information
@brhillman
brhillman committed Aug 22, 2024
2 parents c6ac698 + a39296c commit 69f24bf
Show file tree
Hide file tree
Showing 11 changed files with 460 additions and 223 deletions.
2 changes: 1 addition & 1 deletion cime_config/tests.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -370,7 +370,7 @@
"ERP_Ln9.ne4pg2_oQU480.WCYCL20TRNS-MMF1.allactive-mmf_fixed_subcycle",
"ERS_Ln9.ne4pg2_ne4pg2.FRCE-MMF1.eam-cosp_nhtfrq9",
"SMS_Ln5.ne4_ne4.FSCM-ARM97-MMF1",
# "SMS_Ln3.ne4pg2_ne4pg2.F2010-MMF2", - temporarily disabled due to ongoing dev issues
"SMS_Ln3.ne4pg2_oQU480.F2010-MMF2",
)
},

Expand Down
11 changes: 6 additions & 5 deletions components/eam/cime_config/config_component.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -71,14 +71,15 @@

<!-- Multiscale modeling framework (MMF) -->
<value compset="_EAM%.*MMF1" >-crm samxx -crm_dt 10 </value>
<value compset="_EAM%.*MMF2_" >-crm pam -pam_dycor spam -crm_dt 8 </value>
<value compset="_EAM%.*MMF2-AWFL">-crm pam -pam_dycor awfl -crm_dt 8 </value>
<value compset="_EAM%.*MMF2_" >-crm pam -pam_dycor spam -crm_dt 10 </value>
<value compset="_EAM%.*MMF2-AWFL">-crm pam -pam_dycor awfl -crm_dt 10 </value>
<value compset="_EAM%.*MMF" >-use_MMF -nlev 60 -crm_nz 50</value>
<value compset="_EAM%.*MMF" >-crm_dx 2000 -crm_nx 64 -crm_ny 1 </value>
<value compset="_EAM%.*MMF1" >-crm_dx 2000 -crm_nx 64 -crm_ny 1 -crm_nx_rad 4 -crm_ny_rad 1</value>
<value compset="_EAM%.*MMF2" >-crm_dx 3000 -crm_nx 45 -crm_ny 1 -crm_nx_rad 5 -crm_ny_rad 1</value>
<value compset="_EAM%.*MMF1" >-MMF_microphysics_scheme sam1mom -chem none</value>
<value compset="_EAM%.*MMF2" >-MMF_microphysics_scheme p3 -chem none</value>
<value compset="_EAM%.*MMF" >-crm_nx_rad 4 -crm_ny_rad 1 -rad rrtmgp -rrtmgpxx</value>
<value compset="_EAM%.*MMF1" >-use_MMF_VT</value> <!-- CRM variance transport -->
<value compset="_EAM%.*MMF" >-rad rrtmgp -rrtmgpxx</value>
<value compset="_EAM%.*MMF" >-use_MMF_VT</value> <!-- CRM variance transport -->
<value compset="_EAM%.*MMF1" >-use_MMF_ESMT</value> <!-- explicit scalar momentum transport -->
<value compset="_EAM%AQP-MMF.*" >-aquaplanet</value>
<value compset="_EAM%RCE-MMF.*" >-aquaplanet -rce</value>
Expand Down
4 changes: 2 additions & 2 deletions components/eam/src/physics/crm/crm_physics.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1379,8 +1379,8 @@ subroutine crm_physics_tend(ztodt, state, tend, ptend, pbuf2d, cam_in, cam_out,
call pam_mirror_array_readwrite( 'output_ni_mean', crm_output%ni_mean, '' )
call pam_mirror_array_readwrite( 'output_qm_mean', crm_output%qm_mean, '' )
call pam_mirror_array_readwrite( 'output_bm_mean', crm_output%bm_mean, '' )
call pam_mirror_array_readwrite( 'output_rho_d_mean', crm_output%rho_d_mean, '' )
call pam_mirror_array_readwrite( 'output_rho_v_mean', crm_output%rho_v_mean, '' )
call pam_mirror_array_readwrite( 'output_rho_d_mean', crm_output%rho_d_mean, '' )
call pam_mirror_array_readwrite( 'output_rho_v_mean', crm_output%rho_v_mean, '' )

call pam_mirror_array_readwrite( 'output_qt_ls', crm_output%qt_ls, '' )
call pam_mirror_array_readwrite( 'output_t_ls', crm_output%t_ls, '' )
Expand Down
6 changes: 4 additions & 2 deletions components/eam/src/physics/crm/pam/CMakeLists.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,9 @@ set(PAM_DRIVER_SRC
add_library(pam_driver
${PAM_DRIVER_SRC})

set(SCREAM_HOME ${CMAKE_CURRENT_SOURCE_DIR}/../../../../../..)
add_library(eamxx_physics INTERFACE ${SCREAM_HOME}/components/eamxx/src/physics/)

if (${USE_CUDA})
# PAM will be CUDA-linked with device symbols resolved at library creation
set_target_properties(pam_driver
Expand All @@ -37,7 +40,6 @@ set(PAM_MICRO p3)
set(PAM_SGS shoc)
set(PAM_RAD forced)
set(PAM_SCREAM_USE_CXX True)
set(SCREAM_HOME ${CMAKE_CURRENT_SOURCE_DIR}/../../../../../..)

# removing this once the cime cmake_macros aren't using obsolete
# Kokkos settings like KOKKOS_OPTIONS.
Expand Down Expand Up @@ -124,7 +126,7 @@ target_compile_options(pam_driver PUBLIC
)

if (PAM_SCREAM_USE_CXX)
target_link_libraries(pam_driver pam_core physics dynamics pam_scream_cxx_interfaces ekat p3 shoc physics_share scream_share)
target_link_libraries(pam_driver pam_core physics dynamics pam_scream_cxx_interfaces ekat p3 shoc physics_share scream_share eamxx_physics)
else()
target_link_libraries(pam_driver pam_core physics dynamics )
endif()
Expand Down
2 changes: 1 addition & 1 deletion components/eam/src/physics/crm/pam/external
Submodule external updated 10 files
+25 −5 dynamics/spam/Dycore.h
+34 −0 dynamics/spam/src/core/params.h
+71 −0 dynamics/spam/src/models/extrudedmodel.h
+6 −0 dynamics/spam/src/models/layermodel.h
+5 −3 pam_core/modules/gcm_forcing.h
+3 −4 physics/scream_cxx_interfaces/scream_cxx_interface_p3.cpp
+5 −1 physics/scream_cxx_interfaces/scream_cxx_interface_shoc.cpp
+2 −1 physics/scream_cxx_interfaces/scream_cxx_interface_shoc.h
+0 −24 physics/scream_cxx_p3_shoc/scream_share/CMakeLists.txt
+4 −1 physics/sgs/shoc/SGS.h
85 changes: 39 additions & 46 deletions components/eam/src/physics/crm/pam/pam_accelerate.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,16 +17,22 @@ inline void pam_accelerate_init( pam::PamCoupler &coupler ) {
using yakl::c::parallel_for;
using yakl::c::SimpleBounds;
using yakl::atomicAdd;
//------------------------------------------------------------------------------------------------
// The ability to accelerate dry density is implemented, but tests suggested
// that this was contributing to instabilities in PAM. Combined with the
// uncertainty of how to handle dry density in the anelastic case, we will
// disable dry density mean-state acceleration by default.
coupler.set_option<bool>("crm_accel_rd",false);
//------------------------------------------------------------------------------------------------
auto &dm_device = coupler.get_data_manager_device_readwrite();
auto nens = coupler.get_option<int>("ncrms");
auto nz = coupler.get_option<int>("crm_nz");
auto ny = coupler.get_option<int>("crm_ny");
auto nx = coupler.get_option<int>("crm_nx");
auto crm_accel_uv = coupler.get_option<bool>("crm_accel_uv");
//------------------------------------------------------------------------------------------------
dm_device.register_and_allocate<real>("accel_save_t", "saved temperature for MSA", {nz,nens}, {"z","nens"} );
dm_device.register_and_allocate<real>("accel_save_r", "saved dry density for MSA", {nz,nens}, {"z","nens"} );
dm_device.register_and_allocate<real>("accel_save_q", "saved total water for MSA", {nz,nens}, {"z","nens"} );
dm_device.register_and_allocate<real>("accel_save_r", "saved dry density for MSA", {nz,nens}, {"z","nens"} );
dm_device.register_and_allocate<real>("accel_save_u", "saved uvel for MSA", {nz,nens}, {"z","nens"} );
dm_device.register_and_allocate<real>("accel_save_v", "saved vvel for MSA", {nz,nens}, {"z","nens"} );
//------------------------------------------------------------------------------------------------
Expand All @@ -42,40 +48,37 @@ inline void pam_accelerate_diagnose( pam::PamCoupler &coupler ) {
auto nz = coupler.get_option<int>("crm_nz");
auto ny = coupler.get_option<int>("crm_ny");
auto nx = coupler.get_option<int>("crm_nx");
auto crm_accel_rd = coupler.get_option<bool>("crm_accel_rd");
auto crm_accel_uv = coupler.get_option<bool>("crm_accel_uv");
//------------------------------------------------------------------------------------------------
auto temp = dm_device.get<real,4>("temp" );
auto rhod = dm_device.get<real,4>("density_dry");
auto rhov = dm_device.get<real,4>("water_vapor");
auto rhol = dm_device.get<real,4>("cloud_water");
auto rhoi = dm_device.get<real,4>("ice" );
auto rhod = dm_device.get<real,4>("density_dry");
auto uvel = dm_device.get<real,4>("uvel" );
auto vvel = dm_device.get<real,4>("vvel" );
auto accel_save_t = dm_device.get<real,2>("accel_save_t");
auto accel_save_r = dm_device.get<real,2>("accel_save_r");
auto accel_save_q = dm_device.get<real,2>("accel_save_q");
auto accel_save_r = dm_device.get<real,2>("accel_save_r");
auto accel_save_u = dm_device.get<real,2>("accel_save_u");
auto accel_save_v = dm_device.get<real,2>("accel_save_v");
//------------------------------------------------------------------------------------------------
// compute horizontal means needed later for mean-state acceleration
parallel_for( SimpleBounds<2>(nz,nens) , YAKL_LAMBDA (int k, int n) {
accel_save_t(k,n) = 0.0;
accel_save_r(k,n) = 0.0;
accel_save_q(k,n) = 0.0;
if (crm_accel_uv) {
accel_save_u(k,n) = 0.0;
accel_save_v(k,n) = 0.0;
}
if (crm_accel_rd) { accel_save_r(k,n) = 0.0; }
if (crm_accel_uv) { accel_save_u(k,n) = 0.0; }
if (crm_accel_uv) { accel_save_v(k,n) = 0.0; }
});
real r_nx_ny = 1._fp/(nx*ny); // precompute reciprocal to avoid costly divisions
parallel_for( SimpleBounds<4>(nz,ny,nx,nens) , YAKL_LAMBDA (int k, int j, int i, int n) {
yakl::atomicAdd( accel_save_t(k,n), temp(k,j,i,n) * r_nx_ny );
yakl::atomicAdd( accel_save_r(k,n), rhod(k,j,i,n) * r_nx_ny );
yakl::atomicAdd( accel_save_q(k,n), ( rhov(k,j,i,n) + rhol(k,j,i,n) + rhoi(k,j,i,n) ) * r_nx_ny );
if (crm_accel_uv) {
yakl::atomicAdd( accel_save_u(k,n), uvel(k,j,i,n) * r_nx_ny );
yakl::atomicAdd( accel_save_v(k,n), vvel(k,j,i,n) * r_nx_ny );
}
if (crm_accel_rd) { yakl::atomicAdd( accel_save_r(k,n), rhod(k,j,i,n) * r_nx_ny ); }
if (crm_accel_uv) { yakl::atomicAdd( accel_save_u(k,n), uvel(k,j,i,n) * r_nx_ny ); }
if (crm_accel_uv) { yakl::atomicAdd( accel_save_v(k,n), vvel(k,j,i,n) * r_nx_ny ); }
});
//------------------------------------------------------------------------------------------------
}
Expand All @@ -86,28 +89,28 @@ inline void pam_accelerate( pam::PamCoupler &coupler, int nstep, int &nstop ) {
using yakl::c::SimpleBounds;
using yakl::atomicAdd;
using yakl::ScalarLiveOut;
auto &dm_device = coupler.get_data_manager_device_readwrite();
auto nens = coupler.get_option<int>("ncrms");
auto nz = coupler.get_option<int>("crm_nz");
auto ny = coupler.get_option<int>("crm_ny");
auto nx = coupler.get_option<int>("crm_nx");
auto &dm_device = coupler.get_data_manager_device_readwrite();
auto nens = coupler.get_option<int>("ncrms");
auto nz = coupler.get_option<int>("crm_nz");
auto ny = coupler.get_option<int>("crm_ny");
auto nx = coupler.get_option<int>("crm_nx");
auto crm_accel_rd = coupler.get_option<bool>("crm_accel_rd");
auto crm_accel_uv = coupler.get_option<bool>("crm_accel_uv");
real crm_accel_factor = coupler.get_option<real>("crm_accel_factor");
//------------------------------------------------------------------------------------------------
auto temp = dm_device.get<real,4>("temp" );
auto rhod = dm_device.get<real,4>("density_dry");
auto rhov = dm_device.get<real,4>("water_vapor");
auto rhol = dm_device.get<real,4>("cloud_water");
auto rhoi = dm_device.get<real,4>("ice" );
auto rhod = dm_device.get<real,4>("density_dry");
auto uvel = dm_device.get<real,4>("uvel" );
auto vvel = dm_device.get<real,4>("vvel" );
auto accel_save_t = dm_device.get<real,2>("accel_save_t");
auto accel_save_r = dm_device.get<real,2>("accel_save_r");
auto accel_save_q = dm_device.get<real,2>("accel_save_q");
auto accel_save_r = dm_device.get<real,2>("accel_save_r");
auto accel_save_u = dm_device.get<real,2>("accel_save_u");
auto accel_save_v = dm_device.get<real,2>("accel_save_v");
//------------------------------------------------------------------------------------------------
bool crm_accel_uv = coupler.get_option<bool>("crm_accel_uv");
real crm_accel_factor = coupler.get_option<real>("crm_accel_factor");
//------------------------------------------------------------------------------------------------
real2d hmean_t ("hmean_t", nz,nens);
real2d hmean_r ("hmean_r", nz,nens);
real2d hmean_q ("hmean_q", nz,nens);
Expand All @@ -127,35 +130,29 @@ inline void pam_accelerate( pam::PamCoupler &coupler, int nstep, int &nstop ) {
// Compute the horizontal mean for each variable
parallel_for( SimpleBounds<2>(nz,nens) , YAKL_LAMBDA (int k, int n) {
hmean_t(k,n) = 0.0;
hmean_r(k,n) = 0.0;
hmean_q(k,n) = 0.0;
if (crm_accel_uv) {
hmean_u(k,n) = 0.0;
hmean_v(k,n) = 0.0;
}
if (crm_accel_rd) { hmean_r(k,n) = 0.0; }
if (crm_accel_uv) { hmean_u(k,n) = 0.0; }
if (crm_accel_uv) { hmean_v(k,n) = 0.0; }
});
real r_nx_ny = 1._fp/(nx*ny); // precompute reciprocal to avoid costly divisions
parallel_for( SimpleBounds<4>(nz,ny,nx,nens) , YAKL_LAMBDA (int k, int j, int i, int n) {
yakl::atomicAdd( hmean_t(k,n), temp(k,j,i,n) * r_nx_ny );
yakl::atomicAdd( hmean_r(k,n), rhod(k,j,i,n) * r_nx_ny );
yakl::atomicAdd( hmean_q(k,n), ( rhov(k,j,i,n) + rhol(k,j,i,n) + rhoi(k,j,i,n) ) * r_nx_ny );
if (crm_accel_uv) {
yakl::atomicAdd( hmean_u(k,n), uvel(k,j,i,n) * r_nx_ny );
yakl::atomicAdd( hmean_v(k,n), vvel(k,j,i,n) * r_nx_ny );
}
if (crm_accel_rd) { yakl::atomicAdd( hmean_r(k,n), rhod(k,j,i,n) * r_nx_ny ); }
if (crm_accel_uv) { yakl::atomicAdd( hmean_u(k,n), uvel(k,j,i,n) * r_nx_ny ); }
if (crm_accel_uv) { yakl::atomicAdd( hmean_v(k,n), vvel(k,j,i,n) * r_nx_ny ); }
});
//------------------------------------------------------------------------------------------------
// Compute the accelerated tendencies
// NOTE - these are tendencies multiplied by the time step (i.e. d/dt * crm_dt)
ScalarLiveOut<bool> ceaseflag_liveout(false);
parallel_for( SimpleBounds<2>(nz,nens) , YAKL_LAMBDA (int k, int n) {
ttend_acc(k,n) = hmean_t(k,n) - accel_save_t(k,n);
rtend_acc(k,n) = hmean_r(k,n) - accel_save_r(k,n);
qtend_acc(k,n) = hmean_q(k,n) - accel_save_q(k,n);
if (crm_accel_uv) {
utend_acc(k,n) = hmean_u(k,n) - accel_save_u(k,n);
vtend_acc(k,n) = hmean_v(k,n) - accel_save_v(k,n);
}
if (crm_accel_rd) { rtend_acc(k,n) = hmean_r(k,n) - accel_save_r(k,n); }
if (crm_accel_uv) { utend_acc(k,n) = hmean_u(k,n) - accel_save_u(k,n); }
if (crm_accel_uv) { vtend_acc(k,n) = hmean_v(k,n) - accel_save_v(k,n); }
if (abs(ttend_acc(k,n)) > dtemp_max) {
ceaseflag_liveout = true;
}
Expand Down Expand Up @@ -191,14 +188,11 @@ inline void pam_accelerate( pam::PamCoupler &coupler, int nstep, int &nstop ) {
//------------------------------------------------------------------------------------------------
// Apply the accelerated tendencies
parallel_for( SimpleBounds<4>(nz,ny,nx,nens) , YAKL_LAMBDA (int k, int j, int i, int n) {

temp(k,j,i,n) = temp(k,j,i,n) + crm_accel_factor * ttend_acc(k,n);
rhod(k,j,i,n) = rhod(k,j,i,n) + crm_accel_factor * rtend_acc(k,n);
rhov(k,j,i,n) = rhov(k,j,i,n) + crm_accel_factor * qtend_acc(k,n);
if (crm_accel_uv) {
uvel(k,j,i,n) = uvel(k,j,i,n) + crm_accel_factor * utend_acc(k,n);
vvel(k,j,i,n) = vvel(k,j,i,n) + crm_accel_factor * vtend_acc(k,n);
}
if (crm_accel_rd) { rhod(k,j,i,n) = rhod(k,j,i,n) + crm_accel_factor * rtend_acc(k,n); }
if (crm_accel_uv) { uvel(k,j,i,n) = uvel(k,j,i,n) + crm_accel_factor * utend_acc(k,n); }
if (crm_accel_uv) { vvel(k,j,i,n) = vvel(k,j,i,n) + crm_accel_factor * vtend_acc(k,n); }
// apply limiter on temperature
temp(k,j,i,n) = std::max( temp_min, temp(k,j,i,n) );
});
Expand Down Expand Up @@ -231,4 +225,3 @@ inline void pam_accelerate( pam::PamCoupler &coupler, int nstep, int &nstop ) {
});
//------------------------------------------------------------------------------------------------
}

Loading

0 comments on commit 69f24bf

Please sign in to comment.