diff --git a/cime_config/tests.py b/cime_config/tests.py index 3223ae10804..87249990a59 100644 --- a/cime_config/tests.py +++ b/cime_config/tests.py @@ -690,6 +690,7 @@ "time" : "01:00:00", "tests" : ( "SMS_D_Ln5.ne4pg2_oQU480.F2010-SCREAMv1-MPASSI.scream-mam4xx-optics", + "SMS_D_Ln5.ne4pg2_oQU480.F2010-SCREAMv1-MPASSI.scream-mam4xx-aci", ) }, diff --git a/components/eamxx/cime_config/namelist_defaults_scream.xml b/components/eamxx/cime_config/namelist_defaults_scream.xml index 9a33879b01d..50019833d8a 100644 --- a/components/eamxx/cime_config/namelist_defaults_scream.xml +++ b/components/eamxx/cime_config/namelist_defaults_scream.xml @@ -234,6 +234,12 @@ be lost if SCREAM_HACK_XML is not enabled. 0.1 + + + 0.001 + 6 + + @@ -481,6 +487,33 @@ be lost if SCREAM_HACK_XML is not enabled. 0.0 0.0 0.0,0.0 + + 2.6e-08 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 + 0.0 0.0 diff --git a/components/eamxx/cime_config/testdefs/testmods_dirs/scream/mam4xx/aci/shell_commands b/components/eamxx/cime_config/testdefs/testmods_dirs/scream/mam4xx/aci/shell_commands new file mode 100644 index 00000000000..68eaf51698b --- /dev/null +++ b/components/eamxx/cime_config/testdefs/testmods_dirs/scream/mam4xx/aci/shell_commands @@ -0,0 +1,23 @@ + +#Default scream has 10 tracers, MAM4xx adds another 31 making a total of 41 tracer +#Set total number of tracers to 41. We are using append here as last entry wins while parsing xml options +./xmlchange --append SCREAM_CMAKE_OPTIONS="SCREAM_NUM_TRACERS 41" + +#modify initial condition file to get aerosol species ICs +$CIMEROOT/../components/eamxx/scripts/atmchange initial_conditions::Filename='$DIN_LOC_ROOT/atm/scream/init/screami_mam4xx_ne4np4L72_c20240208.nc' -b + +# Add spa as RRTMG needs spa +$CIMEROOT/../components/eamxx/scripts/atmchange physics::atm_procs_list="mac_aero_mic,spa,rrtmgp" -b + +# Replace spa with mam4_aci to invoke mam4 aci scheme +$CIMEROOT/../components/eamxx/scripts/atmchange mac_aero_mic::atm_procs_list="tms,shoc,cldFraction,mam4_aci,p3" -b + +#Set precribed ccn to false so that P3 uses input from ACI +$CIMEROOT/../components/eamxx/scripts/atmchange p3::do_prescribed_ccn=false -b + +#Set predicted ccn to true so that P3 uses input from ACI +$CIMEROOT/../components/eamxx/scripts/atmchange p3::do_predict_nc=true -b + + + + diff --git a/components/eamxx/cmake/machine-files/compy.cmake b/components/eamxx/cmake/machine-files/compy.cmake index 1f156b0546c..e8404d75559 100644 --- a/components/eamxx/cmake/machine-files/compy.cmake +++ b/components/eamxx/cmake/machine-files/compy.cmake @@ -7,4 +7,4 @@ include (${EKAT_MACH_FILES_PATH}/kokkos/openmp.cmake) include (${EKAT_MACH_FILES_PATH}/mpi/srun.cmake) #Compy SLURM specific settings -set(EKAT_MPI_NP_FLAG "-p short -n" CACHE STRING "" FORCE) +set(EKAT_MPI_NP_FLAG "-p short --mpi=pmi2 -n" CACHE STRING "" FORCE) diff --git a/components/eamxx/src/physics/mam/CMakeLists.txt b/components/eamxx/src/physics/mam/CMakeLists.txt index 53eb4049f11..cfd33e88d0b 100644 --- a/components/eamxx/src/physics/mam/CMakeLists.txt +++ b/components/eamxx/src/physics/mam/CMakeLists.txt @@ -42,7 +42,8 @@ add_subdirectory(${EXTERNALS_SOURCE_DIR}/mam4xx ${CMAKE_BINARY_DIR}/externals/ma # EAMxx mam4xx-based atmospheric processes add_library(mam eamxx_mam_microphysics_process_interface.cpp - eamxx_mam_optics_process_interface.cpp) + eamxx_mam_optics_process_interface.cpp + eamxx_mam_aci_process_interface.cpp) target_compile_definitions(mam PUBLIC EAMXX_HAS_MAM) add_dependencies(mam mam4xx) target_include_directories(mam PUBLIC diff --git a/components/eamxx/src/physics/mam/eamxx_mam_aci_functions.hpp b/components/eamxx/src/physics/mam/eamxx_mam_aci_functions.hpp new file mode 100644 index 00000000000..9fedb64555e --- /dev/null +++ b/components/eamxx/src/physics/mam/eamxx_mam_aci_functions.hpp @@ -0,0 +1,700 @@ +#ifndef EAMXX_MAM_ACI_FUNCTION_HPP +#define EAMXX_MAM_ACI_FUNCTION_HPP + +#include +#include +#include + +namespace scream { + +namespace { + +KOKKOS_INLINE_FUNCTION +void compute_w0_and_rho(const haero::ThreadTeam &team, + const MAMAci::const_view_2d omega, + const MAMAci::const_view_2d T_mid, + const MAMAci::const_view_2d p_mid, const int icol, + const int top_lev, const int nlev, + // output + MAMAci::view_2d w0, MAMAci::view_2d rho) { + // Get physical constants + using C = physics::Constants; + static constexpr auto gravit = C::gravit; // Gravity [m/s2] + // Gas constant for dry air [J/(kg*K) or J/Kg/K] + static constexpr auto rair = C::Rair; + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, 0u, top_lev), KOKKOS_LAMBDA(int kk) { + w0(icol, kk) = 0; + rho(icol, kk) = -999.0; + }); + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, top_lev, nlev), KOKKOS_LAMBDA(int kk) { + rho(icol, kk) = p_mid(icol, kk) / (rair * T_mid(icol, kk)); + w0(icol, kk) = -1.0 * omega(icol, kk) / (rho(icol, kk) * gravit); + }); +} +void compute_w0_and_rho(haero::ThreadTeamPolicy team_policy, + const mam_coupling::DryAtmosphere &dry_atmosphere, + const int top_lev, const int nlev, + // output + MAMAci::view_2d w0, MAMAci::view_2d rho) { + MAMAci::const_view_2d omega = dry_atmosphere.omega; + MAMAci::const_view_2d T_mid = dry_atmosphere.T_mid; + MAMAci::const_view_2d p_mid = dry_atmosphere.p_mid; + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + compute_w0_and_rho(team, omega, T_mid, p_mid, icol, top_lev, nlev, + // output + w0, rho); + }); +} + +void compute_values_at_interfaces(haero::ThreadTeamPolicy team_policy, + const MAMAci::const_view_2d var_mid, + const MAMAci::view_2d dz, const int nlev_, + // output + MAMAci::view_2d var_int) { + using CO = scream::ColumnOps; + + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + + const auto var_mid_col = ekat::subview(var_mid, icol); + const auto var_int_col = ekat::subview(var_int, icol); + const auto dz_col = ekat::subview(dz, icol); + + const Real bc_top = var_mid_col(0); + const Real bc_bot = var_mid_col(nlev_ - 1); + + CO::compute_interface_values_linear(team, nlev_, var_mid_col, dz_col, + bc_top, bc_bot, var_int_col); + }); +} + +KOKKOS_INLINE_FUNCTION +void compute_tke_using_w_sec(const haero::ThreadTeam &team, + const MAMAci::const_view_2d w_sec, const int icol, + const int nlev, + // output + MAMAci::view_2d tke) { + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, 0u, nlev + 1), + KOKKOS_LAMBDA(int kk) { tke(icol, kk) = (3.0 / 2.0) * w_sec(icol, kk); }); +} +void compute_tke_using_w_sec(haero::ThreadTeamPolicy team_policy, + const MAMAci::const_view_2d w_sec, const int nlev, + // output + MAMAci::view_2d tke) { + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + compute_tke_using_w_sec(team, w_sec, icol, nlev, + // output + tke); + }); +} +KOKKOS_INLINE_FUNCTION +void compute_subgrid_scale_velocities( + const haero::ThreadTeam &team, const MAMAci::const_view_2d tke, + const Real wsubmin, const int icol, const int top_lev, const int nlev, + // output + MAMAci::view_2d wsub, MAMAci::view_2d wsubice, MAMAci::view_2d wsig) { + // More refined computation of sub-grid vertical velocity + // Set to be zero at the surface by initialization. + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, 0u, top_lev), KOKKOS_LAMBDA(int kk) { + wsub(icol, kk) = wsubmin; + wsubice(icol, kk) = 0.001; + wsig(icol, kk) = 0.001; + }); + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, top_lev, nlev), KOKKOS_LAMBDA(int kk) { + wsub(icol, kk) = haero::sqrt(0.5 * (tke(icol, kk) + tke(icol, kk + 1)) * + (2.0 / 3.0)); + wsig(icol, kk) = + mam4::utils::min_max_bound(0.001, 10.0, wsub(icol, kk)); + wsubice(icol, kk) = + mam4::utils::min_max_bound(0.2, 10.0, wsub(icol, kk)); + wsub(icol, kk) = haero::max(wsubmin, wsub(icol, kk)); + }); +} +void compute_subgrid_scale_velocities( + haero::ThreadTeamPolicy team_policy, const MAMAci::const_view_2d tke, + const Real wsubmin, const int top_lev, const int nlev, + // output + MAMAci::view_2d wsub, MAMAci::view_2d wsubice, MAMAci::view_2d wsig) { + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + compute_subgrid_scale_velocities(team, tke, wsubmin, icol, top_lev, + nlev, + // output + wsub, wsubice, wsig); + }); +} + +void compute_nucleate_ice_tendencies( + const mam4::NucleateIce &nucleate_ice, haero::ThreadTeamPolicy team_policy, + const mam_coupling::DryAtmosphere &dry_atmosphere, + const mam_coupling::AerosolState &dry_aero, const MAMAci::view_2d wsubice, + const MAMAci::view_2d aitken_dry_dia, const int nlev, const double dt, + // output + MAMAci::view_2d nihf, MAMAci::view_2d niim, MAMAci::view_2d nidep, + MAMAci::view_2d nimey, MAMAci::view_2d naai_hom, + // ## output used by other processes ## + MAMAci::view_2d naai) { + //------------------------------------------------------------- + // Get number of activated aerosol for ice nucleation (naai) + // from ice nucleation + //------------------------------------------------------------- + + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + //--------------------------------------------------------------------- + // Set up surface, pronostics atmosphere, diagnostics, and tendencies + // classes. + //--------------------------------------------------------------------- + + // For compute_tendecies interface only, this structure is empty + haero::Surface surf{}; + + // Store interstitial and cld borne aerosols in "progrs" struture + mam4::Prognostics progs = + mam_coupling::aerosols_for_column(dry_aero, icol); + + // Store atmopsheric vars (Tmid, Pmid, cloud fraction, qv, wsubmin) + haero::Atmosphere haero_atm = + atmosphere_for_column(dry_atmosphere, icol); + + // Update the updraft velocity needed by nucleation to be "wsubice" + // in the haero_atm object + haero_atm.updraft_vel_ice_nucleation = ekat::subview(wsubice, icol); + + // All the output from this process is diagnotics; creates "diags" with + // nlev column length + mam4::Diagnostics diags(nlev); + + // Aitken mode index + const int aitken_idx = static_cast(mam4::ModeIndex::Aitken); + diags.dry_geometric_mean_diameter_i[aitken_idx] = + ekat::subview(aitken_dry_dia, icol); + + // These are the fields that are updated. Taking subviews means that + // the nihf, niim, nidep, nimey, naai_hom, and naai fields are updated + // in nucleate_ice.compute_tendencies. + diags.icenuc_num_hetfrz = ekat::subview(nihf, icol); + diags.icenuc_num_immfrz = ekat::subview(niim, icol); + diags.icenuc_num_depnuc = ekat::subview(nidep, icol); + diags.icenuc_num_meydep = ekat::subview(nimey, icol); + + // naai and naai_hom are the outputs needed for nucleate_ice and these + // are not tendencies. + diags.num_act_aerosol_ice_nucle_hom = ekat::subview(naai_hom, icol); + diags.num_act_aerosol_ice_nucle = ekat::subview(naai, icol); + + // grab views from the buffer to store tendencies, not used as all + // values are store in diags above. + const mam4::Tendencies tends(nlev); // not used + const mam4::AeroConfig aero_config; + const Real t = 0; // not used + nucleate_ice.compute_tendencies(aero_config, team, t, dt, haero_atm, + surf, progs, diags, tends); + }); +} +KOKKOS_INLINE_FUNCTION +void store_liquid_cloud_fraction( + const haero::ThreadTeam &team, const MAMAci::const_view_2d qc, + const MAMAci::const_view_2d qi, const MAMAci::const_view_2d liqcldf, + const MAMAci::const_view_2d liqcldf_prev, const int icol, const int top_lev, + const int nlev, + // output + MAMAci::view_2d cloud_frac, MAMAci::view_2d cloud_frac_prev) { + //------------------------------------------------------------- + // Get old and new liquid cloud fractions when amount of cloud + // is above qsmall threshold value + //------------------------------------------------------------- + // cut-off for cloud amount (ice or liquid) + static constexpr auto qsmall = 1e-18; + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, top_lev, nlev), KOKKOS_LAMBDA(int kk) { + if((qc(icol, kk) + qi(icol, kk)) > qsmall) { + cloud_frac(icol, kk) = liqcldf(icol, kk); + cloud_frac_prev(icol, kk) = liqcldf_prev(icol, kk); + } else { + cloud_frac(icol, kk) = 0; + cloud_frac_prev(icol, kk) = 0; + } + }); +} +void store_liquid_cloud_fraction( + haero::ThreadTeamPolicy team_policy, + const mam_coupling::DryAtmosphere &dry_atmosphere, + const MAMAci::const_view_2d liqcldf, + const MAMAci::const_view_2d liqcldf_prev, const int top_lev, const int nlev, + // output + MAMAci::view_2d cloud_frac, MAMAci::view_2d cloud_frac_prev) { + MAMAci::const_view_2d qc = dry_atmosphere.qc; + MAMAci::const_view_2d qi = dry_atmosphere.qi; + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + store_liquid_cloud_fraction(team, qc, qi, liqcldf, liqcldf_prev, icol, + top_lev, nlev, + // output + cloud_frac, cloud_frac_prev); + }); +} +KOKKOS_INLINE_FUNCTION +void compute_recipical_pseudo_density(const haero::ThreadTeam &team, + const MAMAci::const_view_2d pdel, + const int icol, const int nlev, + // output + MAMAci::view_2d rpdel) { + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, 0u, nlev), KOKKOS_LAMBDA(int kk) { + EKAT_KERNEL_ASSERT_MSG(0 < pdel(icol, kk), + "Error: pdel should be > 0.\n"); + rpdel(icol, kk) = 1 / pdel(icol, kk); + }); +} +void compute_recipical_pseudo_density(haero::ThreadTeamPolicy team_policy, + MAMAci::const_view_2d pdel, + const int nlev, + // output + MAMAci::view_2d rpdel) { + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + compute_recipical_pseudo_density(team, pdel, icol, nlev, + // output + rpdel); + }); +} + +void call_function_dropmixnuc( + haero::ThreadTeamPolicy team_policy, const Real dt, + mam_coupling::DryAtmosphere &dry_atmosphere, const MAMAci::view_2d rpdel, + const MAMAci::const_view_2d kvh, const MAMAci::view_2d wsub, + const MAMAci::view_2d cloud_frac, const MAMAci::view_2d cloud_frac_prev, + const mam_coupling::AerosolState &dry_aero, const int nlev, + + // Following outputs are all diagnostics + MAMAci::view_2d coltend[mam4::ndrop::ncnst_tot], + MAMAci::view_2d coltend_cw[mam4::ndrop::ncnst_tot], MAMAci::view_2d qcld, + MAMAci::view_2d ndropcol, MAMAci::view_2d ndropmix, MAMAci::view_2d nsource, + MAMAci::view_2d wtke, MAMAci::view_3d ccn, + + // ## outputs to be used by other processes ## + // qqcw_fld_work should be directly assigned to the cloud borne aerosols + MAMAci::view_2d qqcw_fld_work[mam4::ndrop::ncnst_tot], + + // ptend_q are the tendencies to the interstitial aerosols + MAMAci::view_2d ptend_q[mam4::aero_model::pcnst], + + // factnum is used by the hetrozenous freezing + MAMAci::view_3d factnum, + + // tendnd is used by microphysics scheme (e.g. P3) + MAMAci::view_2d tendnd, + + // ## work arrays ## + MAMAci::view_2d raercol_cw[mam4::ndrop::pver][2], + MAMAci::view_2d raercol[mam4::ndrop::pver][2], MAMAci::view_3d state_q_work, + MAMAci::view_3d nact, MAMAci::view_3d mact, + MAMAci::view_2d dropmixnuc_scratch_mem[MAMAci::dropmix_scratch_]) { + // Extract atmosphere variables + MAMAci::const_view_2d T_mid = dry_atmosphere.T_mid; + MAMAci::const_view_2d p_mid = dry_atmosphere.p_mid; + MAMAci::const_view_2d zm = dry_atmosphere.z_mid; + MAMAci::const_view_2d pdel = dry_atmosphere.p_del; + MAMAci::const_view_2d p_int = dry_atmosphere.p_int; + MAMAci::const_view_2d nc = dry_atmosphere.nc; + + //---------------------------------------------------------------------- + // ## Declare local variables for class variables + //(FIXME: GPU hack, revisit this) + //---------------------------------------------------------------------- + MAMAci::view_2d loc_raercol_cw[mam4::ndrop::pver][2]; + MAMAci::view_2d loc_raercol[mam4::ndrop::pver][2]; + MAMAci::view_2d loc_qqcw[mam4::ndrop::ncnst_tot]; + MAMAci::view_2d loc_ptend_q[mam4::aero_model::pcnst]; + MAMAci::view_2d loc_coltend[mam4::ndrop::ncnst_tot]; + MAMAci::view_2d loc_coltend_cw[mam4::ndrop::ncnst_tot]; + + for(int i = 0; i < mam4::ndrop::pver; ++i) { + for(int j = 0; j < 2; ++j) { + loc_raercol_cw[i][j] = raercol_cw[i][j]; + loc_raercol[i][j] = raercol[i][j]; + } + } + + for(int i = 0; i < mam4::ndrop::ncnst_tot; ++i) { + loc_coltend[i] = coltend[i]; + loc_coltend_cw[i] = coltend_cw[i]; + } + + for(int i = 0; i < mam4::aero_model::pcnst; ++i) loc_ptend_q[i] = ptend_q[i]; + + MAMAci::view_2d qqcw_fld_work_loc[25]; + for(int i = 0; i < mam4::ndrop::ncnst_tot; ++i) + qqcw_fld_work_loc[i] = qqcw_fld_work[i]; + + MAMAci::view_3d state_q_work_loc = state_q_work; + + //---------------------------------------------------------------------- + // ## Assign scratch memory for work variables + //---------------------------------------------------------------------- + + MAMAci::view_2d eddy_diff = dropmixnuc_scratch_mem[0]; + MAMAci::view_2d zn = dropmixnuc_scratch_mem[1]; + MAMAci::view_2d csbot = dropmixnuc_scratch_mem[2]; + MAMAci::view_2d zs = dropmixnuc_scratch_mem[3]; + MAMAci::view_2d overlapp = dropmixnuc_scratch_mem[4]; + MAMAci::view_2d overlapm = dropmixnuc_scratch_mem[5]; + MAMAci::view_2d eddy_diff_kp = dropmixnuc_scratch_mem[6]; + MAMAci::view_2d eddy_diff_km = dropmixnuc_scratch_mem[7]; + MAMAci::view_2d qncld = dropmixnuc_scratch_mem[8]; + MAMAci::view_2d srcn = dropmixnuc_scratch_mem[9]; + MAMAci::view_2d source = dropmixnuc_scratch_mem[10]; + MAMAci::view_2d dz = dropmixnuc_scratch_mem[11]; + MAMAci::view_2d csbot_cscen = dropmixnuc_scratch_mem[12]; + MAMAci::view_2d raertend = dropmixnuc_scratch_mem[13]; + MAMAci::view_2d qqcwtend = dropmixnuc_scratch_mem[14]; + + //--------------------------------------------------------------------------- + // ## Initialize the ndrop class. + //--------------------------------------------------------------------------- + const int ntot_amode = mam_coupling::num_aero_modes(); + const int maxd_aspectype = mam4::ndrop::maxd_aspectype; + const int nspec_max = mam4::ndrop::nspec_max; + int nspec_amode[ntot_amode] = {}; + int lspectype_amode[maxd_aspectype][ntot_amode] = {}; + int lmassptr_amode[maxd_aspectype][ntot_amode] = {}; + int numptr_amode[ntot_amode] = {}; + int mam_idx[ntot_amode][nspec_max] = {}; + int mam_cnst_idx[ntot_amode][nspec_max] = {}; + + Real specdens_amode[maxd_aspectype] = {}; + Real spechygro[maxd_aspectype] = {}; + Real exp45logsig[ntot_amode] = {}, alogsig[ntot_amode] = {}, + num2vol_ratio_min_nmodes[ntot_amode] = {}, + num2vol_ratio_max_nmodes[ntot_amode] = {}; + Real aten = 0; + mam4::ndrop::get_e3sm_parameters(nspec_amode, lspectype_amode, lmassptr_amode, + numptr_amode, specdens_amode, spechygro, + mam_idx, mam_cnst_idx); + mam4::ndrop::ndrop_init(exp45logsig, alogsig, aten, num2vol_ratio_min_nmodes, + num2vol_ratio_max_nmodes); + //--------------------------------------------------------------------------- + //--------------------------------------------------------------------------- + + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + // for (int icol=0; icol<5; ++icol){ + mam4::ndrop::View1D raercol_cw_view[mam4::ndrop::pver][2]; + mam4::ndrop::View1D raercol_view[mam4::ndrop::pver][2]; + for(int i = 0; i < mam4::ndrop::pver; ++i) { + for(int j = 0; j < 2; ++j) { + raercol_cw_view[i][j] = ekat::subview(loc_raercol_cw[i][j], icol); + raercol_view[i][j] = ekat::subview(loc_raercol[i][j], icol); + } + } + mam4::ColumnView qqcw_view[mam4::ndrop::ncnst_tot]; + for(int i = 0; i < mam4::ndrop::ncnst_tot; ++i) { + qqcw_view[i] = ekat::subview(qqcw_fld_work_loc[i], icol); + } + mam4::ColumnView ptend_q_view[mam4::aero_model::pcnst]; + for(int i = 0; i < mam4::aero_model::pcnst; ++i) { + ptend_q_view[i] = ekat::subview(loc_ptend_q[i], icol); + } + mam4::ColumnView coltend_view[mam4::ndrop::ncnst_tot], + coltend_cw_view[mam4::ndrop::ncnst_tot]; + for(int i = 0; i < mam4::ndrop::ncnst_tot; ++i) { + coltend_view[i] = ekat::subview(loc_coltend[i], icol); + coltend_cw_view[i] = ekat::subview(loc_coltend_cw[i], icol); + } + + // To construct state_q, we need fields from Prognostics (aerosols) + // and Atmosphere (water species such as qv, qc etc.) + + // get prognostics + mam4::Prognostics progs_at_col = aerosols_for_column(dry_aero, icol); + + // get atmospheric quantities + haero::Atmosphere haero_atm = + atmosphere_for_column(dry_atmosphere, icol); + + // Construct state_q (interstitial) and qqcw (cloud borne) arrays + constexpr auto pver = mam4::ndrop::pver; + Kokkos::parallel_for( + Kokkos::TeamVectorRange(team, 0u, pver), [&](int klev) { + Real state_q_at_lev_col[mam4::aero_model::pcnst] = {}; + + // get state_q at a grid cell (col,lev) + // NOTE: The order of species in state_q_at_lev_col + // is the same as in E3SM state%q array + mam4::utils::extract_stateq_from_prognostics( + progs_at_col, haero_atm, state_q_at_lev_col, klev); + + // get the start index for aerosols species in the state_q array + int istart = mam4::aero_model::pcnst - mam4::ndrop::ncnst_tot; + + // create colum views of state_q + for(int icnst = istart; icnst < mam4::aero_model::pcnst; + ++icnst) { + state_q_work_loc(icol, klev, icnst) = state_q_at_lev_col[icnst]; + } + + // get qqcw at a grid cell (col,lev) + // NOTE: The layout for qqcw array is based on mam_idx in + // dropmixnuc. To mimic that, we are using the following for-loops + int ind_qqcw = 0; + for(int m = 0; m < mam_coupling::num_aero_modes(); ++m) { + qqcw_view[ind_qqcw](klev) = + dry_aero.cld_aero_nmr[m](icol, klev); + ++ind_qqcw; + for(int a = 0; a < mam_coupling::num_aero_species(); ++a) { + if(dry_aero.cld_aero_mmr[m][a].data()) { + qqcw_view[ind_qqcw](klev) = + dry_aero.cld_aero_mmr[m][a](icol, klev); + ++ind_qqcw; + } + } + } + }); + + mam4::ndrop::dropmixnuc( + team, dt, ekat::subview(T_mid, icol), ekat::subview(p_mid, icol), + ekat::subview(p_int, icol), ekat::subview(pdel, icol), + ekat::subview(rpdel, icol), + // in zm[kk] - zm[kk+1], for pver zm[kk-1] - zm[kk] + ekat::subview(zm, icol), ekat::subview(state_q_work_loc, icol), + ekat::subview(nc, icol), ekat::subview(kvh, icol), // kvh[kk+1] + ekat::subview(cloud_frac, icol), lspectype_amode, specdens_amode, + spechygro, lmassptr_amode, num2vol_ratio_min_nmodes, + num2vol_ratio_max_nmodes, numptr_amode, nspec_amode, exp45logsig, + alogsig, aten, mam_idx, mam_cnst_idx, + ekat::subview(qcld, icol), // out + ekat::subview(wsub, icol), // in + ekat::subview(cloud_frac_prev, icol), // in + qqcw_view, // inout + ptend_q_view, ekat::subview(tendnd, icol), + ekat::subview(factnum, icol), ekat::subview(ndropcol, icol), + ekat::subview(ndropmix, icol), ekat::subview(nsource, icol), + ekat::subview(wtke, icol), ekat::subview(ccn, icol), coltend_view, + coltend_cw_view, + // work arrays + raercol_cw_view, raercol_view, ekat::subview(nact, icol), + ekat::subview(mact, icol), ekat::subview(eddy_diff, icol), + ekat::subview(zn, icol), ekat::subview(csbot, icol), + ekat::subview(zs, icol), ekat::subview(overlapp, icol), + ekat::subview(overlapm, icol), ekat::subview(eddy_diff_kp, icol), + ekat::subview(eddy_diff_km, icol), ekat::subview(qncld, icol), + ekat::subview(srcn, icol), ekat::subview(source, icol), + ekat::subview(dz, icol), ekat::subview(csbot_cscen, icol), + ekat::subview(raertend, icol), ekat::subview(qqcwtend, icol)); + }); +} + +// Update cloud borne aerosols +void update_cloud_borne_aerosols( + const MAMAci::view_2d qqcw_fld_work[mam4::ndrop::ncnst_tot], const int nlev, + // output + mam_coupling::AerosolState &dry_aero) { + int ind_qqcw = 0; + for(int m = 0; m < mam_coupling::num_aero_modes(); ++m) { + Kokkos::deep_copy(dry_aero.cld_aero_nmr[m], qqcw_fld_work[ind_qqcw]); + ++ind_qqcw; + for(int a = 0; a < mam_coupling::num_aero_species(); ++a) { + if(dry_aero.cld_aero_mmr[m][a].data()) { + Kokkos::deep_copy(dry_aero.cld_aero_mmr[m][a], qqcw_fld_work[ind_qqcw]); + ++ind_qqcw; + } + } + } +} + +// Update interstitial aerosols using tendencies - levels +KOKKOS_INLINE_FUNCTION +void update_interstitial_aerosols_levs(const haero::ThreadTeam &team, + const int nlev, const int icol, + const Real dt, + const MAMAci::view_2d ptend_view, + // output + MAMAci::view_2d aero_mr) { + Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev), [&](int kk) { + aero_mr(icol, kk) += ptend_view(icol, kk) * dt; + }); +} + +// Update interstitial aerosols using tendencies- cols and levs +void update_interstitial_aerosols( + haero::ThreadTeamPolicy team_policy, + const MAMAci::view_2d ptend_q[mam4::aero_model::pcnst], const int nlev, + const Real dt, + // output + mam_coupling::AerosolState &dry_aero) { + // starting index of ptend_q array (for MAM4, pcnst=40, ncnst_tot=25 ) + int s_idx = mam4::aero_model::pcnst - mam4::ndrop::ncnst_tot; + + // loop through all modes and species + for(int m = 0; m < mam_coupling::num_aero_modes(); ++m) { + for(int a = 0; a < mam4::num_species_mode(m); ++a) { + // mass mixing ratio of species "a" of mode "m" + auto aero_mmr = dry_aero.int_aero_mmr[m][a]; + + if(aero_mmr.data()) { + const auto ptend_view = ptend_q[s_idx]; + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + // update values for all levs at this column + update_interstitial_aerosols_levs(team, nlev, icol, dt, + ptend_view, + // output + aero_mmr); + }); + // update index for the next species (only if aero_mmr.data() is True) + ++s_idx; + } + } + auto aero_nmr = + dry_aero.int_aero_nmr[m]; // number mixing ratio for mode "m" + const auto ptend_view = ptend_q[s_idx]; + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + // update values for all levs at this column + update_interstitial_aerosols_levs(team, nlev, icol, dt, ptend_view, + // output + aero_nmr); + }); + ++s_idx; // update index for the next species + } +} + +void call_hetfrz_compute_tendencies( + haero::ThreadTeamPolicy team_policy, mam4::Hetfrz &hetfrz_, + mam_coupling::DryAtmosphere &dry_atm_, + mam_coupling::AerosolState &dry_aero_, MAMAci::view_3d factnum, + const double dt, const int nlev, + // output + MAMAci::view_2d hetfrz_immersion_nucleation_tend, + MAMAci::view_2d hetfrz_contact_nucleation_tend, + MAMAci::view_2d hetfrz_depostion_nucleation_tend, + MAMAci::view_2d diagnostic_scratch_[]) { + mam4::Hetfrz hetfrz = hetfrz_; + mam_coupling::AerosolState dry_aero = dry_aero_; + mam_coupling::DryAtmosphere dry_atmosphere = dry_atm_; + + MAMAci::view_2d diagnostic_scratch[MAMAci::hetro_scratch_]; + for(int i = 0; i < MAMAci::hetro_scratch_; ++i) + diagnostic_scratch[i] = diagnostic_scratch_[i]; + + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + // Set up an atmosphere, surface, diagnostics, pronostics and + // tendencies class. + + haero::Atmosphere haero_atm = + atmosphere_for_column(dry_atmosphere, icol); + haero::Surface surf{}; + mam4::Prognostics progs = + mam_coupling::aerosols_for_column(dry_aero, icol); + + const int accum_idx = static_cast(mam4::ModeIndex::Accumulation); + const int coarse_idx = static_cast(mam4::ModeIndex::Coarse); + + mam4::Diagnostics diags(nlev); + + diags.activation_fraction[accum_idx] = + ekat::subview(factnum, icol, accum_idx); + + diags.activation_fraction[coarse_idx] = + ekat::subview(factnum, icol, coarse_idx); + + // These are the output tendencies from heterogeneous freezing that need + // to be added correctly to the cloud-micorphysics scheme. + diags.hetfrz_immersion_nucleation_tend = + ekat::subview(hetfrz_immersion_nucleation_tend, icol); + diags.hetfrz_contact_nucleation_tend = + ekat::subview(hetfrz_contact_nucleation_tend, icol); + diags.hetfrz_depostion_nucleation_tend = + ekat::subview(hetfrz_depostion_nucleation_tend, icol); + + diags.bc_num = ekat::subview(diagnostic_scratch[0], icol); + diags.dst1_num = ekat::subview(diagnostic_scratch[1], icol); + diags.dst3_num = ekat::subview(diagnostic_scratch[2], icol); + diags.bcc_num = ekat::subview(diagnostic_scratch[3], icol); + diags.dst1c_num = ekat::subview(diagnostic_scratch[4], icol); + diags.dst3c_num = ekat::subview(diagnostic_scratch[5], icol); + diags.bcuc_num = ekat::subview(diagnostic_scratch[6], icol); + diags.dst1uc_num = ekat::subview(diagnostic_scratch[7], icol); + diags.dst3uc_num = ekat::subview(diagnostic_scratch[8], icol); + diags.bc_a1_num = ekat::subview(diagnostic_scratch[0], icol); + diags.dst_a1_num = ekat::subview(diagnostic_scratch[10], icol); + diags.dst_a3_num = ekat::subview(diagnostic_scratch[11], icol); + diags.bc_c1_num = ekat::subview(diagnostic_scratch[12], icol); + diags.dst_c1_num = ekat::subview(diagnostic_scratch[13], icol); + diags.dst_c3_num = ekat::subview(diagnostic_scratch[14], icol); + diags.fn_bc_c1_num = ekat::subview(diagnostic_scratch[15], icol); + diags.fn_dst_c1_num = ekat::subview(diagnostic_scratch[16], icol); + diags.fn_dst_c3_num = ekat::subview(diagnostic_scratch[17], icol); + diags.na500 = ekat::subview(diagnostic_scratch[18], icol); + diags.totna500 = ekat::subview(diagnostic_scratch[19], icol); + diags.freqimm = ekat::subview(diagnostic_scratch[20], icol); + diags.freqcnt = ekat::subview(diagnostic_scratch[21], icol); + diags.freqdep = ekat::subview(diagnostic_scratch[22], icol); + diags.freqmix = ekat::subview(diagnostic_scratch[23], icol); + diags.dstfrezimm = ekat::subview(diagnostic_scratch[24], icol); + diags.dstfrezcnt = ekat::subview(diagnostic_scratch[25], icol); + diags.dstfrezdep = ekat::subview(diagnostic_scratch[26], icol); + diags.bcfrezimm = ekat::subview(diagnostic_scratch[27], icol); + diags.bcfrezcnt = ekat::subview(diagnostic_scratch[28], icol); + diags.bcfrezdep = ekat::subview(diagnostic_scratch[19], icol); + diags.nimix_imm = ekat::subview(diagnostic_scratch[30], icol); + diags.nimix_cnt = ekat::subview(diagnostic_scratch[31], icol); + diags.nimix_dep = ekat::subview(diagnostic_scratch[32], icol); + diags.dstnidep = ekat::subview(diagnostic_scratch[33], icol); + diags.dstnicnt = ekat::subview(diagnostic_scratch[34], icol); + diags.dstniimm = ekat::subview(diagnostic_scratch[35], icol); + diags.bcnidep = ekat::subview(diagnostic_scratch[36], icol); + diags.bcnicnt = ekat::subview(diagnostic_scratch[37], icol); + diags.bcniimm = ekat::subview(diagnostic_scratch[38], icol); + diags.numice10s = ekat::subview(diagnostic_scratch[39], icol); + diags.numimm10sdst = ekat::subview(diagnostic_scratch[40], icol); + diags.numimm10sbc = ekat::subview(diagnostic_scratch[41], icol); + diags.stratiform_cloud_fraction = + ekat::subview(diagnostic_scratch[42], icol); + + // assign cloud fraction + constexpr auto pver = mam4::ndrop::pver; + Kokkos::parallel_for(Kokkos::TeamVectorRange(team, 0u, pver), + [&](int klev) { + diags.stratiform_cloud_fraction(klev) = + haero_atm.cloud_fraction(klev); + }); + //------------------------------------------------------------- + // Heterogeneous freezing + // frzimm, frzcnt, frzdep are the outputs of + // hetfrz_classnuc_cam_calc used by the microphysics (e.g. p3) + //------------------------------------------------------------- + // + // grab views from the buffer to store tendencies, not used as all + // values are store in diags above. + const mam4::Tendencies tends(nlev); + const mam4::AeroConfig aero_config; + const Real t = 0; // not used + hetfrz.compute_tendencies(aero_config, team, t, dt, haero_atm, surf, + progs, diags, tends); + }); +} +} // namespace +} // namespace scream + +#endif // ifdef EAMXX_MAM_ACI_FUNCTION_HPP diff --git a/components/eamxx/src/physics/mam/eamxx_mam_aci_process_interface.cpp b/components/eamxx/src/physics/mam/eamxx_mam_aci_process_interface.cpp new file mode 100644 index 00000000000..8f9552537ab --- /dev/null +++ b/components/eamxx/src/physics/mam/eamxx_mam_aci_process_interface.cpp @@ -0,0 +1,679 @@ +#include + +// ACI functions are stored in the following hpp file +#include + +// For EKAT units package +#include "ekat/util/ekat_units.hpp" + +/* +----------------------------------------------------------------- +NOTES: +1. w_variance assumes that we are using SE dycore. If we use EUL +dycore, we need to get its value from previous dynamic time step + +2. w_variance and eddy_diff_heat are at midpoints in EAMxx, we +derived their interface values using boundary conditions to be +the top and botton values of the mid point arrays. We assume that +this assumption should not cause any issues. + +FUTURE WORK: +1. MAM4xx submodule should point to MAM4xx main branch +2. Link hetrozenous freezing outputs to microphysics +3. Add postcondition and invariant checks +5. Resolve comments about top_lev +6. Avoid using c-style static arrays in ptend and other arrays +7. Use std::string rather than c-strings +8. Remove a Kokkos:fence and combine two kernels while computing w_sec_int_ +9. Fix double counting of tracer advection by modifying SHOC. +10. A git issue for computing top_lev, moving liq_cldfrac in ACI and using TKE +directly +14. Merge kernels so that we are not calling one from another +16. improve the way qqcw is populated +17.delete fence mentioned by Luca +----------------------------------------------------------------- +*/ + +namespace scream { +MAMAci::MAMAci(const ekat::Comm &comm, const ekat::ParameterList ¶ms) + : AtmosphereProcess(comm, params) { + // Asserts for the runtime or namelist options + EKAT_REQUIRE_MSG(m_params.isParameter("wsubmin"), + "ERROR: wsubmin is missing from mam_aci parameter list."); + EKAT_REQUIRE_MSG( + m_params.isParameter("top_level_mam4xx"), + "ERROR: top_level_mam4xx is missing from mam_aci parameter list."); +} + +// ================================================================ +// SET_GRIDS +// ================================================================ +void MAMAci::set_grids( + const std::shared_ptr grids_manager) { + // set grid for all the inputs and outputs + // use physics grid + grid_ = grids_manager->get_grid("Physics"); + + // Name of the grid + const auto &grid_name = grid_->name(); + + ncol_ = grid_->get_num_local_dofs(); // Number of columns on this rank + nlev_ = grid_->get_num_vertical_levels(); // Number of levels per column + + // Define the different field layouts that will be used for this process + using namespace ShortFieldTagsNames; + + // Layout for 3D (2d horiz X 1d vertical) variables + // mid points + FieldLayout scalar3d_layout_mid{{COL, LEV}, {ncol_, nlev_}}; + // interfaces + FieldLayout scalar3d_layout_int{{COL, ILEV}, {ncol_, nlev_ + 1}}; + + // layout for 2D (1d horiz X 1d vertical) variable + FieldLayout scalar2d_layout_col{{COL}, {ncol_}}; + + auto make_layout = [](const std::vector &extents, + const std::vector &names) { + std::vector tags(extents.size(), CMP); + return FieldLayout(tags, extents, names); + }; + + using namespace ekat::units; + auto q_unit = kg / kg; // units of mass mixing ratios of tracers + auto n_unit = 1 / kg; // units of number mixing ratios of tracers + + auto nondim = ekat::units::Units::nondimensional(); + + // atmospheric quantities + // specific humidity [kg/kg] + add_field("qv", scalar3d_layout_mid, q_unit, grid_name, "tracers"); + + // cloud liquid mass mixing ratio [kg/kg] + add_field("qc", scalar3d_layout_mid, q_unit, grid_name, "tracers"); + + // cloud ice mass mixing ratio [kg/kg] + add_field("qi", scalar3d_layout_mid, q_unit, grid_name, "tracers"); + + // cloud liquid number mixing ratio [1/kg] + add_field("nc", scalar3d_layout_mid, n_unit, grid_name, "tracers"); + + // cloud ice number mixing ratio [1/kg] + add_field("ni", scalar3d_layout_mid, n_unit, grid_name, "tracers"); + + // Temperature[K] at midpoints + add_field("T_mid", scalar3d_layout_mid, K, grid_name); + + // Vertical pressure velocity [Pa/s] at midpoints + add_field("omega", scalar3d_layout_mid, Pa / s, grid_name); + + // Total pressure [Pa] at midpoints + add_field("p_mid", scalar3d_layout_mid, Pa, grid_name); + + // Total pressure [Pa] at interfaces + add_field("p_int", scalar3d_layout_int, Pa, grid_name); + + // Layer thickness(pdel) [Pa] at midpoints + add_field("pseudo_density", scalar3d_layout_mid, Pa, grid_name); + + // planetary boundary layer height + add_field("pbl_height", scalar2d_layout_col, m, grid_name); + + // cloud fraction [nondimensional] computed by eamxx_cld_fraction_process + add_field("cldfrac_tot", scalar3d_layout_mid, nondim, grid_name); + + auto m2 = pow(m, 2); + auto s2 = pow(s, 2); + + // NOTE: w_variance im microp_aero.F90 in EAM is at "itim_old" dynamics time + // step. Since, we are using SE dycore, itim_old is 1 which is equivalent to + // the current time step. For other dycores (such as EUL), it may be different + // and we might need to revisit this. + + // Vertical velocity variance at midpoints + add_field("w_variance", scalar3d_layout_mid, m2 / s2, grid_name); + + // NOTE: "cldfrac_liq" is updated in SHOC. "cldfrac_liq" in C++ code is + // equivalent to "alst" in the shoc_intr.F90. In the C++ code, it is used as + // "shoc_cldfrac" and in the F90 code it is called "cloud_frac" + + // Liquid stratiform cloud fraction at midpoints + add_field("cldfrac_liq", scalar3d_layout_mid, nondim, grid_name); + + // Previous value of liquid stratiform cloud fraction at midpoints + add_field("cldfrac_liq_prev", scalar3d_layout_mid, nondim, + grid_name); + + // Eddy diffusivity for heat + add_field("eddy_diff_heat", scalar3d_layout_mid, m2 / s, grid_name); + + // Layout for 4D (2d horiz X 1d vertical x number of modes) variables + const int num_aero_modes = mam_coupling::num_aero_modes(); + FieldLayout scalar4d_layout_mid = + make_layout({ncol_, num_aero_modes, nlev_}, {"COL", "NMODES", "LEV"}); + + // dry diameter of aerosols [m] + add_field("dgnum", scalar4d_layout_mid, m, grid_name); + + // ======================================================================== + // Output from this whole process + // ======================================================================== + + // interstitial and cloudborne aerosol tracers of interest: mass (q) and + // number (n) mixing ratios + for(int mode = 0; mode < mam_coupling::num_aero_modes(); ++mode) { + // interstitial aerosol tracers of interest: number (n) mixing ratios + const char *int_nmr_field_name = + mam_coupling::int_aero_nmr_field_name(mode); + add_field(int_nmr_field_name, scalar3d_layout_mid, n_unit, + grid_name, "tracers"); + + // cloudborne aerosol tracers of interest: number (n) mixing ratios + // NOTE: DO NOT add cld borne aerosols to the "tracer" group as these are + // NOT advected + const char *cld_nmr_field_name = + mam_coupling::cld_aero_nmr_field_name(mode); + add_field(cld_nmr_field_name, scalar3d_layout_mid, n_unit, + grid_name); + + for(int a = 0; a < mam_coupling::num_aero_species(); ++a) { + // (interstitial) aerosol tracers of interest: mass (q) mixing ratios + const char *int_mmr_field_name = + mam_coupling::int_aero_mmr_field_name(mode, a); + if(strlen(int_mmr_field_name) > 0) { + add_field(int_mmr_field_name, scalar3d_layout_mid, q_unit, + grid_name, "tracers"); + } + // (cloudborne) aerosol tracers of interest: mass (q) mixing ratios + // NOTE: DO NOT add cld borne aerosols to the "tracer" group as these are + // NOT advected + const char *cld_mmr_field_name = + mam_coupling::cld_aero_mmr_field_name(mode, a); + if(strlen(cld_mmr_field_name) > 0) { + add_field(cld_mmr_field_name, scalar3d_layout_mid, q_unit, + grid_name); + } + } // end for loop num species + } // end for loop for num modes + + for(int g = 0; g < mam_coupling::num_aero_gases(); ++g) { + const char *gas_mmr_field_name = mam_coupling::gas_mmr_field_name(g); + add_field(gas_mmr_field_name, scalar3d_layout_mid, q_unit, + grid_name, "tracers"); + } // end for loop num gases + + // ------------------------------------------------------------------------ + // Output from ice nucleation process + // ------------------------------------------------------------------------ + + // number of activated aerosol for ice nucleation[#/kg] + add_field("ni_activated", scalar3d_layout_mid, n_unit, grid_name); + + // ------------------------------------------------------------------------ + // Output from droplet activation process (dropmixnuc) + // ------------------------------------------------------------------------ + + // tendency in droplet number mixing ratio [#/kg/s] + add_field("nc_nuceat_tend", scalar3d_layout_mid, n_unit / s, + grid_name); + + // FIXME: [TEMPORARY]droplet number mixing ratio source tendency [#/kg/s] + add_field("nsource", scalar3d_layout_mid, n_unit / s, grid_name); + + // FIXME: [TEMPORARY]droplet number mixing ratio tendency due to mixing + // [#/kg/s] + add_field("ndropmix", scalar3d_layout_mid, n_unit / s, grid_name); + + // FIXME: [TEMPORARY]droplet number as seen by ACI [#/kg] + add_field("nc_inp_to_aci", scalar3d_layout_mid, n_unit / s, + grid_name); + const auto cm_tmp = m / 100; // FIXME: [TEMPORARY] remove this + const auto cm3 = cm_tmp * cm_tmp * cm_tmp; // FIXME: [TEMPORARY] remove this + // FIXME: [TEMPORARY] remove the following ccn outputs + add_field("ccn_0p02", scalar3d_layout_mid, cm3, grid_name); + add_field("ccn_0p05", scalar3d_layout_mid, cm3, grid_name); + add_field("ccn_0p1", scalar3d_layout_mid, cm3, grid_name); + add_field("ccn_0p2", scalar3d_layout_mid, cm3, grid_name); + add_field("ccn_0p5", scalar3d_layout_mid, cm3, grid_name); + add_field("ccn_1p0", scalar3d_layout_mid, cm3, grid_name); + + // ------------------------------------------------------------------------ + // Output from hetrozenous freezing + // ------------------------------------------------------------------------ + + const auto cm = m / 100; + + // units of number mixing ratios of tracers + auto frz_unit = 1 / (cm * cm * cm * s); + // heterogeneous freezing by immersion nucleation [cm^-3 s^-1] + add_field("hetfrz_immersion_nucleation_tend", scalar3d_layout_mid, + frz_unit, grid_name); + + // heterogeneous freezing by contact nucleation [cm^-3 s^-1] + add_field("hetfrz_contact_nucleation_tend", scalar3d_layout_mid, + frz_unit, grid_name); + + // heterogeneous freezing by deposition nucleation [cm^-3 s^-1] + add_field("hetfrz_depostion_nucleation_tend", scalar3d_layout_mid, + frz_unit, grid_name); +} // function set_grids ends + +// ================================================================ +// INIT_BUFFERS +// ================================================================ + +void MAMAci::init_buffers(const ATMBufferManager &buffer_manager) { + EKAT_REQUIRE_MSG( + buffer_manager.allocated_bytes() >= requested_buffer_size_in_bytes(), + "Error! Insufficient buffer size.\n"); + + size_t used_mem = + mam_coupling::init_buffer(buffer_manager, ncol_, nlev_, buffer_); + EKAT_REQUIRE_MSG( + used_mem == requested_buffer_size_in_bytes(), + "Error! Used memory != requested memory for MAMMicrophysics."); +} // function init_buffers ends + +// ================================================================ +// INITIALIZE_IMPL +// ================================================================ +void MAMAci::initialize_impl(const RunType run_type) { + // ------------------------------------------------------------------------ + // ## Runtime options + // ------------------------------------------------------------------------ + + wsubmin_ = m_params.get("wsubmin"); + top_lev_ = m_params.get("top_level_mam4xx"); + + // ------------------------------------------------------------------------ + // Input fields read in from IC file, namelist or other processes + // ------------------------------------------------------------------------ + w_sec_mid_ = get_field_in("w_variance").get_view(); + dgnum_ = get_field_in("dgnum").get_view(); + liqcldf_ = get_field_in("cldfrac_liq").get_view(); + liqcldf_prev_ = get_field_in("cldfrac_liq_prev").get_view(); + kvh_mid_ = get_field_in("eddy_diff_heat").get_view(); + + // store fields only to be converted to dry mmrs in wet_atm_ + wet_atm_.qv = get_field_in("qv").get_view(); + wet_atm_.qc = get_field_in("qc").get_view(); + wet_atm_.nc = get_field_in("nc").get_view(); + wet_atm_.qi = get_field_in("qi").get_view(); + wet_atm_.ni = get_field_in("ni").get_view(); + + // store rest fo the atm fields in dry_atm_in + dry_atm_.z_surf = 0; + dry_atm_.T_mid = get_field_in("T_mid").get_view(); + dry_atm_.p_mid = get_field_in("p_mid").get_view(); + dry_atm_.p_int = get_field_in("p_int").get_view(); + dry_atm_.p_del = get_field_in("pseudo_density").get_view(); + dry_atm_.omega = get_field_in("omega").get_view(); + + // store fields converted to dry mmr from wet mmr in dry_atm_ + dry_atm_.qv = buffer_.qv_dry; + dry_atm_.qc = buffer_.qc_dry; + dry_atm_.nc = buffer_.nc_dry; + dry_atm_.qi = buffer_.qi_dry; + dry_atm_.ni = buffer_.ni_dry; + + // pbl_height + dry_atm_.pblh = get_field_in("pbl_height").get_view(); + + // geometric thickness of layers (m) + dry_atm_.dz = buffer_.dz; + + // geopotential height above surface at interface levels (m) + dry_atm_.z_iface = buffer_.z_iface; + + // geopotential height above surface at mid levels (m) + dry_atm_.z_mid = buffer_.z_mid; + + // total cloud fraction + dry_atm_.cldfrac = get_field_in("cldfrac_tot").get_view(); + + // computed updraft velocity + dry_atm_.w_updraft = buffer_.w_updraft; + + // ------------------------------------------------------------------------ + // Output fields to be used by other processes + // ------------------------------------------------------------------------ + // ice nucleation output + naai_ = get_field_out("ni_activated").get_view(); + + // droplet activation output + tendnd_ = get_field_out("nc_nuceat_tend").get_view(); + + // interstitial and cloudborne aerosol tracers of interest: mass (q) and + // number (n) mixing ratios + for(int m = 0; m < mam_coupling::num_aero_modes(); ++m) { + // interstitial aerosol tracers of interest: number (n) mixing ratios + const char *int_nmr_field_name = mam_coupling::int_aero_nmr_field_name(m); + wet_aero_.int_aero_nmr[m] = + get_field_out(int_nmr_field_name).get_view(); + dry_aero_.int_aero_nmr[m] = buffer_.dry_int_aero_nmr[m]; + + // cloudborne aerosol tracers of interest: number (n) mixing ratios + const char *cld_nmr_field_name = mam_coupling::cld_aero_nmr_field_name(m); + wet_aero_.cld_aero_nmr[m] = + get_field_out(cld_nmr_field_name).get_view(); + dry_aero_.cld_aero_nmr[m] = buffer_.dry_cld_aero_nmr[m]; + + for(int a = 0; a < mam_coupling::num_aero_species(); ++a) { + // (interstitial) aerosol tracers of interest: mass (q) mixing ratios + const char *int_mmr_field_name = + mam_coupling::int_aero_mmr_field_name(m, a); + + if(strlen(int_mmr_field_name) > 0) { + wet_aero_.int_aero_mmr[m][a] = + get_field_out(int_mmr_field_name).get_view(); + dry_aero_.int_aero_mmr[m][a] = buffer_.dry_int_aero_mmr[m][a]; + } + + // (cloudborne) aerosol tracers of interest: mass (q) mixing ratios + const char *cld_mmr_field_name = + mam_coupling::cld_aero_mmr_field_name(m, a); + if(strlen(cld_mmr_field_name) > 0) { + wet_aero_.cld_aero_mmr[m][a] = + get_field_out(cld_mmr_field_name).get_view(); + dry_aero_.cld_aero_mmr[m][a] = buffer_.dry_cld_aero_mmr[m][a]; + } + } + } + for(int g = 0; g < mam_coupling::num_aero_gases(); ++g) { + const char *gas_mmr_field_name = mam_coupling::gas_mmr_field_name(g); + wet_aero_.gas_mmr[g] = + get_field_out(gas_mmr_field_name).get_view(); + dry_aero_.gas_mmr[g] = buffer_.dry_gas_mmr[g]; + } + + // hetrozenous freezing outputs + hetfrz_immersion_nucleation_tend_ = + get_field_out("hetfrz_immersion_nucleation_tend").get_view(); + hetfrz_contact_nucleation_tend_ = + get_field_out("hetfrz_contact_nucleation_tend").get_view(); + hetfrz_depostion_nucleation_tend_ = + get_field_out("hetfrz_depostion_nucleation_tend").get_view(); + + //--------------------------------------------------------------------------------- + // Allocate memory for the class members + // (Kokkos::resize only works on host to allocates memory) + //--------------------------------------------------------------------------------- + + Kokkos::resize(rho_, ncol_, nlev_); + Kokkos::resize(w0_, ncol_, nlev_); + Kokkos::resize(tke_, ncol_, nlev_ + 1); + Kokkos::resize(wsub_, ncol_, nlev_); + Kokkos::resize(wsubice_, ncol_, nlev_); + Kokkos::resize(wsig_, ncol_, nlev_); + Kokkos::resize(w2_, ncol_, nlev_); + Kokkos::resize(cloud_frac_, ncol_, nlev_); + Kokkos::resize(cloud_frac_prev_, ncol_, nlev_); + Kokkos::resize(aitken_dry_dia_, ncol_, nlev_); + Kokkos::resize(rpdel_, ncol_, nlev_); + + //--------------------------------------------------------------------------------- + // Diagnotics variables from the ice nucleation scheme + //--------------------------------------------------------------------------------- + + // number conc of ice nuclei due to heterogeneous freezing [1/m3] + Kokkos::resize(nihf_, ncol_, nlev_); + + // number conc of ice nuclei due to immersionfreezing (hetero nuc) [1/m3] + Kokkos::resize(niim_, ncol_, nlev_); + + // number conc of ice nuclei due to deposition nucleation (hetero nuc)[1/m3] + Kokkos::resize(nidep_, ncol_, nlev_); + + // number conc of ice nuclei due to meyers deposition [1/m3] + Kokkos::resize(nimey_, ncol_, nlev_); + + // number of activated aerosol for ice nucleation(homogeneous frz only)[#/kg] + Kokkos::resize(naai_hom_, ncol_, nlev_); + + //--------------------------------------------------------------------------------- + // Diagnotics variables from the droplet activation scheme + //--------------------------------------------------------------------------------- + + // activation fraction for aerosol number [fraction] + const int num_aero_modes = mam_coupling::num_aero_modes(); + Kokkos::resize(factnum_, ncol_, num_aero_modes, nlev_); + + // cloud droplet number mixing ratio [#/kg] + Kokkos::resize(qcld_, ncol_, nlev_); + + // number conc of aerosols activated at supersat [#/m^3] + // NOTE: activation fraction fluxes are defined as + // fluxn = [flux of activated aero. number into cloud[#/m^2/s]] + // / [aero. number conc. in updraft, just below cloudbase [#/m^3]] + Kokkos::resize(ccn_, ncol_, nlev_, mam4::ndrop::psat); + + // column-integrated droplet number [#/m2] + Kokkos::resize(ndropcol_, ncol_, nlev_); + + // droplet number mixing ratio tendency due to mixing [#/kg/s] + // Kokkos::resize(ndropmix_, ncol_, nlev_); + // Temporarily output ndropmix_ + ndropmix_ = get_field_out("ndropmix").get_view(); + + // droplet number mixing ratio source tendency [#/kg/s] + // Kokkos::resize(nsource_, ncol_, nlev_); + // Temporarily output nsource_ + nsource_ = get_field_out("nsource").get_view(); + + // Temporarily output nc_inp_to_aci_ + nc_inp_to_aci_ = get_field_out("nc_inp_to_aci").get_view(); + + // FIXME: [TEMPORARY] remove the following ccn outputs + ccn_0p02_ = get_field_out("ccn_0p02").get_view(); + ccn_0p05_ = get_field_out("ccn_0p05").get_view(); + ccn_0p1_ = get_field_out("ccn_0p1").get_view(); + ccn_0p2_ = get_field_out("ccn_0p2").get_view(); + ccn_0p5_ = get_field_out("ccn_0p5").get_view(); + ccn_1p0_ = get_field_out("ccn_1p0").get_view(); + + // subgrid vertical velocity [m/s] + Kokkos::resize(wtke_, ncol_, nlev_); + + for(int i = 0; i < dropmix_scratch_; ++i) { + Kokkos::resize(dropmixnuc_scratch_mem_[i], ncol_, nlev_); + } + for(int i = 0; i < mam4::ndrop::ncnst_tot; ++i) { + // column tendency for diagnostic output + Kokkos::resize(coltend_[i], ncol_, nlev_); + // column tendency + Kokkos::resize(coltend_cw_[i], ncol_, nlev_); + } + for(int i = 0; i < mam4::aero_model::pcnst; ++i) { + Kokkos::resize(ptend_q_[i], ncol_, nlev_); + } + for(int i = 0; i < mam4::ndrop::pver; ++i) { + for(int j = 0; j < 2; ++j) { + Kokkos::resize(raercol_cw_[i][j], ncol_, mam4::ndrop::ncnst_tot); + Kokkos::resize(raercol_[i][j], ncol_, mam4::ndrop::ncnst_tot); + } + } + + // nact : fractional aero. number activation rate [/s] + Kokkos::resize(nact_, ncol_, nlev_, mam_coupling::num_aero_modes()); + + // mact : fractional aero. mass activation rate [/s] + Kokkos::resize(mact_, ncol_, nlev_, mam_coupling::num_aero_modes()); + + // Eddy diffusivity of heat at the interfaces + Kokkos::resize(kvh_int_, ncol_, nlev_ + 1); + + // Vertical velocity variance at the interfaces + Kokkos::resize(w_sec_int_, ncol_, nlev_ + 1); + // Allocate work arrays + for(int icnst = 0; icnst < mam4::ndrop::ncnst_tot; ++icnst) { + qqcw_fld_work_[icnst] = view_2d("qqcw_fld_work_", ncol_, nlev_); + } + state_q_work_ = + view_3d("state_q_work_", ncol_, nlev_, mam4::aero_model::pcnst); + + //--------------------------------------------------------------------------------- + // Diagnotics variables from the hetrozenous ice nucleation scheme + //--------------------------------------------------------------------------------- + + for(int i = 0; i < hetro_scratch_; ++i) + Kokkos::resize(diagnostic_scratch_[i], ncol_, nlev_); + + //--------------------------------------------------------------------------------- + // Initialize the processes + //--------------------------------------------------------------------------------- + + mam4::AeroConfig aero_config; + // configure the nucleation parameterization + mam4::NucleateIce::Config nucleate_ice_config; + nucleate_ice_.init(aero_config, nucleate_ice_config); + + // configure the heterogeneous freezing parameterization + mam4::Hetfrz::Config hetfrz_config; + hetfrz_.init(aero_config, hetfrz_config); + + //--------------------------------------------------------------------------------- + // Setup preprocessing and post processing + //--------------------------------------------------------------------------------- + // set up our preprocess and postprocess functors + preprocess_.initialize(ncol_, nlev_, wet_atm_, wet_aero_, dry_atm_, + dry_aero_); + + postprocess_.initialize(ncol_, nlev_, wet_atm_, wet_aero_, dry_atm_, + dry_aero_); + +} // end function initialize_impl + +// ================================================================ +// RUN_IMPL +// ================================================================ +void MAMAci::run_impl(const double dt) { + const auto scan_policy = ekat::ExeSpaceUtils< + KT::ExeSpace>::get_thread_range_parallel_scan_team_policy(ncol_, nlev_); + + // preprocess input -- needs a scan for the calculation of local derivied + // quantities + Kokkos::parallel_for("preprocess", scan_policy, preprocess_); + Kokkos::fence(); + + haero::ThreadTeamPolicy team_policy(ncol_, Kokkos::AUTO); + + // FIXME: Temporary assignment of nc + mam_coupling::copy_view_lev_slice(team_policy, wet_atm_.nc, nlev_, // inputs + nc_inp_to_aci_); // output + + compute_w0_and_rho(team_policy, dry_atm_, top_lev_, nlev_, + // output + w0_, rho_); + + // Get w_sec_int_ from w_sec_mid_ + compute_values_at_interfaces(team_policy, w_sec_mid_, dry_atm_.dz, nlev_, + // output + w_sec_int_); + Kokkos::fence(); // wait for w_sec_int_ to be computed. + compute_tke_using_w_sec(team_policy, w_sec_int_, nlev_, + // output + tke_); + + Kokkos::fence(); // wait for tke_ to be computed. + compute_subgrid_scale_velocities(team_policy, tke_, wsubmin_, top_lev_, nlev_, + // output + wsub_, wsubice_, wsig_); + + // We need dry diameter for only aitken mode + Kokkos::deep_copy( + aitken_dry_dia_, + ekat::subview_1(dgnum_, static_cast(mam4::ModeIndex::Aitken))); + + Kokkos::fence(); // wait for aitken_dry_dia_ to be computed. + + // Compute Ice nucleation + // NOTE: The Fortran version uses "ast" for cloud fraction which is + // equivalent to "cldfrac_tot" in FM. It is part of the "dry_atm_" struct + compute_nucleate_ice_tendencies( + nucleate_ice_, team_policy, dry_atm_, dry_aero_, wsubice_, + aitken_dry_dia_, nlev_, dt, + // output + nihf_, niim_, nidep_, nimey_, naai_hom_, + // ## output to be used by the other processes ## + naai_); + + // Compute cloud fractions based on cloud threshold + store_liquid_cloud_fraction(team_policy, dry_atm_, liqcldf_, liqcldf_prev_, + top_lev_, nlev_, + // output + cloud_frac_, cloud_frac_prev_); + + compute_recipical_pseudo_density(team_policy, dry_atm_.p_del, nlev_, + // output + rpdel_); + + Kokkos::fence(); // wait for rpdel_ to be computed. + + // Get kvh_int_ from kvh_mid_ + compute_values_at_interfaces(team_policy, kvh_mid_, dry_atm_.dz, nlev_, + // output + kvh_int_); + + Kokkos::fence(); + // Compute activated CCN number tendency (tendnd_) and updated + // cloud borne aerosols (stored in a work array) and interstitial + // aerosols tendencies + call_function_dropmixnuc(team_policy, dt, dry_atm_, rpdel_, kvh_int_, wsub_, + cloud_frac_, cloud_frac_prev_, dry_aero_, nlev_, + // output + coltend_, coltend_cw_, qcld_, ndropcol_, ndropmix_, + nsource_, wtke_, ccn_, + // ## output to be used by the other processes ## + qqcw_fld_work_, ptend_q_, factnum_, tendnd_, + // work arrays + raercol_cw_, raercol_, state_q_work_, nact_, mact_, + dropmixnuc_scratch_mem_); + Kokkos::fence(); // wait for ptend_q_ to be computed. + + Kokkos::deep_copy(ccn_0p02_, + Kokkos::subview(ccn_, Kokkos::ALL(), Kokkos::ALL(), 0)); + Kokkos::deep_copy(ccn_0p05_, + Kokkos::subview(ccn_, Kokkos::ALL(), Kokkos::ALL(), 1)); + Kokkos::deep_copy(ccn_0p1_, + Kokkos::subview(ccn_, Kokkos::ALL(), Kokkos::ALL(), 2)); + Kokkos::deep_copy(ccn_0p2_, + Kokkos::subview(ccn_, Kokkos::ALL(), Kokkos::ALL(), 3)); + Kokkos::deep_copy(ccn_0p5_, + Kokkos::subview(ccn_, Kokkos::ALL(), Kokkos::ALL(), 4)); + Kokkos::deep_copy(ccn_1p0_, + Kokkos::subview(ccn_, Kokkos::ALL(), Kokkos::ALL(), 5)); + + //--------------------------------------------------------------------------- + // NOTE: DO NOT UPDATE cloud borne aerosols using the qqcw_fld_work_ array + // at this point as heterozenous freezing needs to use cloud borne aerosols + // before they are changed by the droplet activation (dropmixnuc) process. + //--------------------------------------------------------------------------- + + // Compute hetrozenous freezing + call_hetfrz_compute_tendencies( + team_policy, hetfrz_, dry_atm_, dry_aero_, factnum_, dt, nlev_, + // ## output to be used by the other processes ## + hetfrz_immersion_nucleation_tend_, hetfrz_contact_nucleation_tend_, + hetfrz_depostion_nucleation_tend_, + // work arrays + diagnostic_scratch_); + + //--------------------------------------------------------------- + // Now update interstitial and cloud borne aerosols + //--------------------------------------------------------------- + + // Update cloud borne aerosols + update_cloud_borne_aerosols(qqcw_fld_work_, nlev_, + // output + dry_aero_); + + // Update interstitial aerosols using tendencies + update_interstitial_aerosols(team_policy, ptend_q_, nlev_, dt, + // output + dry_aero_); + + // call post processing to convert dry mixing ratios to wet mixing ratios + Kokkos::parallel_for("postprocess", scan_policy, postprocess_); + Kokkos::fence(); // wait before returning to calling function +} + +} // namespace scream diff --git a/components/eamxx/src/physics/mam/eamxx_mam_aci_process_interface.hpp b/components/eamxx/src/physics/mam/eamxx_mam_aci_process_interface.hpp new file mode 100644 index 00000000000..97bc90d3f98 --- /dev/null +++ b/components/eamxx/src/physics/mam/eamxx_mam_aci_process_interface.hpp @@ -0,0 +1,274 @@ +#ifndef EAMXX_MAM_ACI_HPP +#define EAMXX_MAM_ACI_HPP + +// For MAM4 aerosol configuration +#include + +// For declaring ACI class derived from atm process class +#include + +// For aerosol configuration +#include "mam4xx/aero_config.hpp" + +// For calling ndrop functions +#include "mam4xx/ndrop.hpp" + +namespace scream { + +class MAMAci final : public scream::AtmosphereProcess { + public: + // declare some constant scratch space lengths + static constexpr int hetro_scratch_ = 43; + static constexpr int dropmix_scratch_ = 15; + + // views for multi-column data + using view_2d = scream::mam_coupling::view_2d; + using const_view_2d = scream::mam_coupling::const_view_2d; + using view_3d = scream::mam_coupling::view_3d; + using const_view_3d = scream::mam_coupling::const_view_3d; + + private: + using KT = ekat::KokkosTypes; + + mam4::NucleateIce nucleate_ice_; + mam4::Hetfrz hetfrz_; + + // views for single-column data + using view_1d = scream::mam_coupling::view_1d; + using const_view_1d = scream::mam_coupling::const_view_1d; + + //------------------------------------------------------------------------ + // ACI runtime ( or namelist) options + //------------------------------------------------------------------------ + + Real wsubmin_; // Minimum subgrid vertical velocity + int top_lev_; // Top level for MAM4xx + + //------------------------------------------------------------------------ + // END: ACI runtime ( or namelist) options + //------------------------------------------------------------------------ + + // rho is air density [kg/m3] + view_2d rho_; + + // w0_ is large scale velocity (m/s) + view_2d w0_; + + // turbulent kinetic energy [m^2/s^2] + view_2d tke_; + + // Dry diameter of the aitken mode (for ice nucleation) + view_2d aitken_dry_dia_; + + // wet mixing ratios (water species) + mam_coupling::WetAtmosphere wet_atm_; + + // dry mixing ratios (water species) + mam_coupling::DryAtmosphere dry_atm_; + + // aerosol dry diameter + const_view_3d dgnum_; + + // ice nucleation diagnostic variables + view_2d nihf_; + view_2d niim_; + view_2d nidep_; + view_2d nimey_; + view_2d naai_hom_; + + // ice nucleation output for FM + view_2d naai_; + + // droplet activation inputs and outputs + view_2d kvh_int_; // Eddy diffusivity of heat at the interfaces + const_view_2d liqcldf_; + const_view_2d liqcldf_prev_; + const_view_2d kvh_mid_; // Eddy diffusivity of heat at the midpoints + + view_2d cloud_frac_; + view_2d cloud_frac_prev_; + view_2d qcld_; + view_2d ptend_q_[mam4::aero_model::pcnst]; + view_3d factnum_; + const_view_3d qqcw_input_; + view_2d qqcw_[mam4::ndrop::ncnst_tot]; + view_2d ndropcol_; + view_2d ndropmix_; + view_2d nsource_; + view_2d nc_inp_to_aci_; // FIXME: TEMPORARY output + view_2d ccn_0p02_; // FIXME: TEMPORARY output + view_2d ccn_0p05_; // FIXME: TEMPORARY output + view_2d ccn_0p1_; // FIXME: TEMPORARY output + view_2d ccn_0p2_; // FIXME: TEMPORARY output + view_2d ccn_0p5_; // FIXME: TEMPORARY output + view_2d ccn_1p0_; // FIXME: TEMPORARY output + view_2d wtke_; + view_3d ccn_; + view_2d coltend_[mam4::ndrop::ncnst_tot]; + view_2d coltend_cw_[mam4::ndrop::ncnst_tot]; + + // raercol_cw_ and raercol_ are work arrays for dropmixnuc, allocated on the + // stack. + view_2d raercol_cw_[mam4::ndrop::pver][2]; + view_2d raercol_[mam4::ndrop::pver][2]; + + view_3d nact_; + view_3d mact_; + view_2d dropmixnuc_scratch_mem_[dropmix_scratch_]; + + // droplet activation output for the FM + view_2d tendnd_; + + // These are the output tendencies from heterogeneous freezing that need to be + // added correctly to the cloud-micorphysics scheme. + view_2d hetfrz_immersion_nucleation_tend_; + view_2d hetfrz_contact_nucleation_tend_; + view_2d hetfrz_depostion_nucleation_tend_; + + view_2d diagnostic_scratch_[hetro_scratch_]; + + // Subgrid scale velocities + view_2d wsub_, wsubice_, wsig_, w2_; + + // local atmospheric state column variables + const_view_2d pdel_; // pressure thickess of layer [Pa] + view_2d rpdel_; // Inverse of pdel_ + const_view_2d w_sec_mid_; // Vertical velocity variance at midpoints + view_2d w_sec_int_; // Vertical velocity variance at interfaces + + // number of horizontal columns and vertical levels + int ncol_, nlev_; + + // aerosol state variables + mam_coupling::AerosolState wet_aero_, dry_aero_; + + // workspace manager for internal local variables + mam_coupling::Buffer buffer_; + + // physics grid for column information + std::shared_ptr grid_; + + // A view array to carry cloud borne aerosol mmrs/nmrs + view_2d qqcw_fld_work_[mam4::ndrop::ncnst_tot]; + + // A view to carry interstitial aerosol mmrs/nmrs + view_3d state_q_work_; + + public: + // Constructor + MAMAci(const ekat::Comm &comm, const ekat::ParameterList ¶ms); + + // Process metadata: Return type of the process + AtmosphereProcessType type() const { return AtmosphereProcessType::Physics; } + + // Return name of the process + std::string name() const { return "mam4_aci"; } + + // grid + void set_grids( + const std::shared_ptr grids_manager) override; + + // management of common atm process memory + size_t requested_buffer_size_in_bytes() const { + return mam_coupling::buffer_size(ncol_, nlev_); + } + + void init_buffers(const ATMBufferManager &buffer_manager) override; + + // process behavior + void initialize_impl(const RunType run_type) override; + void run_impl(const double dt) override; + void finalize_impl(){/*DO NOTHING*/}; + + // Atmosphere processes often have a pre-processing step that constructs + // required variables from the set of fields stored in the field manager. + // This functor implements this step, which is called during run_impl. + struct Preprocess { + Preprocess() = default; + // on host: initializes preprocess functor with necessary state data + void initialize(const int ncol, const int nlev, + const mam_coupling::WetAtmosphere &wet_atm, + const mam_coupling::AerosolState &wet_aero, + const mam_coupling::DryAtmosphere &dry_atm, + const mam_coupling::AerosolState &dry_aero) { + ncol_pre_ = ncol; + nlev_pre_ = nlev; + wet_atm_pre_ = wet_atm; + wet_aero_pre_ = wet_aero; + dry_atm_pre_ = dry_atm; + dry_aero_pre_ = dry_aero; + } + + KOKKOS_INLINE_FUNCTION + void operator()( + const Kokkos::TeamPolicy::member_type &team) const { + const int i = team.league_rank(); // column index + + compute_dry_mixing_ratios(team, wet_atm_pre_, dry_atm_pre_, i); + compute_dry_mixing_ratios(team, wet_atm_pre_, wet_aero_pre_, + dry_aero_pre_, i); + team.team_barrier(); + // vertical heights has to be computed after computing dry mixing ratios + // for atmosphere + compute_vertical_layer_heights(team, dry_atm_pre_, i); + compute_updraft_velocities(team, wet_atm_pre_, dry_atm_pre_, i); + } // operator() + + // local variables for preprocess struct + // number of horizontal columns and vertical levels + int ncol_pre_, nlev_pre_; + + // local atmospheric and aerosol state data + mam_coupling::WetAtmosphere wet_atm_pre_; + mam_coupling::DryAtmosphere dry_atm_pre_; + mam_coupling::AerosolState wet_aero_pre_, dry_aero_pre_; + }; // MAMAci::Preprocess + + // Atmosphere processes often have a post-processing step prepares output + // from this process for the Field Manager. This functor implements this + // step, which is called during run_impl. + // Postprocessing functor + struct Postprocess { + Postprocess() = default; + + // on host: initializes postprocess functor with necessary state data + void initialize(const int ncol, const int nlev, + const mam_coupling::WetAtmosphere &wet_atm, + const mam_coupling::AerosolState &wet_aero, + const mam_coupling::DryAtmosphere &dry_atm, + const mam_coupling::AerosolState &dry_aero) { + ncol_post_ = ncol; + nlev_post_ = nlev; + wet_atm_post_ = wet_atm; + wet_aero_post_ = wet_aero; + dry_atm_post_ = dry_atm; + dry_aero_post_ = dry_aero; + } + + KOKKOS_INLINE_FUNCTION + void operator()( + const Kokkos::TeamPolicy::member_type &team) const { + const int i = team.league_rank(); // column index + compute_wet_mixing_ratios(team, dry_atm_post_, dry_aero_post_, + wet_aero_post_, i); + } // operator() + + // number of horizontal columns and vertical levels + int ncol_post_, nlev_post_; + + // local atmospheric and aerosol state data + mam_coupling::WetAtmosphere wet_atm_post_; + mam_coupling::DryAtmosphere dry_atm_post_; + mam_coupling::AerosolState wet_aero_post_, dry_aero_post_; + }; // Postprocess + + private: + // pre- and postprocessing scratch pads + Preprocess preprocess_; + Postprocess postprocess_; + +}; // MAMAci + +} // namespace scream + +#endif // EAMXX_MAM_ACI_HPP diff --git a/components/eamxx/src/physics/mam/eamxx_mam_microphysics_process_interface.cpp b/components/eamxx/src/physics/mam/eamxx_mam_microphysics_process_interface.cpp index 8bc727c325e..7ea2adfddb6 100644 --- a/components/eamxx/src/physics/mam/eamxx_mam_microphysics_process_interface.cpp +++ b/components/eamxx/src/physics/mam/eamxx_mam_microphysics_process_interface.cpp @@ -201,7 +201,7 @@ void MAMMicrophysics::initialize_impl(const RunType run_type) { wet_atm_.nc = get_field_out("nc").get_view(); wet_atm_.qi = get_field_in("qi").get_view(); wet_atm_.ni = get_field_in("ni").get_view(); - wet_atm_.omega = get_field_in("omega").get_view(); + dry_atm_.T_mid = get_field_in("T_mid").get_view(); dry_atm_.p_mid = get_field_in("p_mid").get_view(); @@ -209,6 +209,7 @@ void MAMMicrophysics::initialize_impl(const RunType run_type) { dry_atm_.cldfrac = get_field_in("cldfrac_tot").get_view(); // FIXME: tot or liq? dry_atm_.pblh = get_field_in("pbl_height").get_view(); dry_atm_.phis = get_field_in("phis").get_view(); + dry_atm_.omega = get_field_in("omega").get_view(); dry_atm_.z_mid = buffer_.z_mid; dry_atm_.dz = buffer_.dz; dry_atm_.z_iface = buffer_.z_iface; @@ -305,6 +306,12 @@ void MAMMicrophysics::run_impl(const double dt) { // allocation perspective auto o3_col_dens = buffer_.scratch[8]; + const_view_1d &col_latitudes = col_latitudes_; + mam_coupling::DryAtmosphere &dry_atm = dry_atm_; + mam_coupling::AerosolState &dry_aero = dry_aero_; + mam4::mo_photo::PhotoTableData &photo_table = photo_table_; + const int nlev = nlev_; + const Config &config = config_; // FIXME: read relevant linoz climatology data from file(s) based on time // FIXME: read relevant chlorine loading data from file based on time @@ -313,21 +320,21 @@ void MAMMicrophysics::run_impl(const double dt) { Kokkos::parallel_for(policy, KOKKOS_LAMBDA(const ThreadTeam& team) { const int icol = team.league_rank(); // column index - Real col_lat = col_latitudes_(icol); // column latitude (degrees?) + Real col_lat = col_latitudes(icol); // column latitude (degrees?) // fetch column-specific atmosphere state data - auto atm = mam_coupling::atmosphere_for_column(dry_atm_, icol); - auto z_iface = ekat::subview(dry_atm_.z_iface, icol); - Real phis = dry_atm_.phis(icol); + auto atm = mam_coupling::atmosphere_for_column(dry_atm, icol); + auto z_iface = ekat::subview(dry_atm.z_iface, icol); + Real phis = dry_atm.phis(icol); // set surface state data haero::Surface sfc{}; // fetch column-specific subviews into aerosol prognostics - mam4::Prognostics progs = mam_coupling::interstitial_aerosols_for_column(dry_aero_, icol); + mam4::Prognostics progs = mam_coupling::interstitial_aerosols_for_column(dry_aero, icol); // set up diagnostics - mam4::Diagnostics diags(nlev_); + mam4::Diagnostics diags(nlev); // calculate o3 column densities (first component of col_dens in Fortran code) auto o3_col_dens_i = ekat::subview(o3_col_dens, icol); @@ -346,7 +353,7 @@ void MAMMicrophysics::run_impl(const double dt) { mam4::ColumnView lwc; // FIXME: liquid water cloud content: where do we get this? mam4::mo_photo::table_photo(photo_rates, atm.pressure, atm.hydrostatic_dp, atm.temperature, o3_col_dens_i, zenith_angle, surf_albedo, lwc, - atm.cloud_fraction, esfact, photo_table_, photo_work_arrays); + atm.cloud_fraction, esfact, photo_table, photo_work_arrays); // compute external forcings at time t(n+1) [molecules/cm^3/s] constexpr int extcnt = mam4::gas_chemistry::extcnt; @@ -355,7 +362,7 @@ void MAMMicrophysics::run_impl(const double dt) { mam4::mo_setext::extfrc_set(forcings, extfrc); // compute aerosol microphysics on each vertical level within this column - Kokkos::parallel_for(Kokkos::TeamThreadRange(team, nlev_), [&](const int k) { + Kokkos::parallel_for(Kokkos::TeamThreadRange(team, nlev), [&](const int k) { constexpr int num_modes = mam4::AeroConfig::num_modes(); constexpr int gas_pcnst = mam_coupling::gas_pcnst(); @@ -428,7 +435,7 @@ void MAMMicrophysics::run_impl(const double dt) { constexpr int indexm = 0; // FIXME: index of xhnm in invariants array (??) Real cldnum = 0.0; // FIXME: droplet number concentration: where do we get this? setsox_single_level(loffset, dt, pmid, pdel, temp, mbar, lwc(k), - cldfrac, cldnum, invariants[indexm], config_.setsox, vmrcw, vmr); + cldfrac, cldnum, invariants[indexm], config.setsox, vmrcw, vmr); // calculate aerosol water content using water uptake treatment // * dry and wet diameters [m] @@ -441,7 +448,7 @@ void MAMMicrophysics::run_impl(const double dt) { impl::compute_water_content(progs, k, qv, temp, pmid, dgncur_a, dgncur_awet, wetdens, qaerwat); // do aerosol microphysics (gas-aerosol exchange, nucleation, coagulation) - impl::modal_aero_amicphys_intr(config_.amicphys, step_, dt, t, pmid, pdel, + impl::modal_aero_amicphys_intr(config.amicphys, step_, dt, t, pmid, pdel, zm, pblh, qv, cldfrac, vmr, vmrcw, vmr_pregaschem, vmr_precldchem, vmrcw_precldchem, vmr_tendbb, vmrcw_tendbb, dgncur_a, dgncur_awet, @@ -466,15 +473,15 @@ void MAMMicrophysics::run_impl(const double dt) { linoz_o3_clim(icol, k), linoz_t_clim(icol, k), linoz_o3col_clim(icol, k), linoz_PmL_clim(icol, k), linoz_dPmL_dO3(icol, k), linoz_dPmL_dT(icol, k), linoz_dPmL_dO3col(icol, k), linoz_cariolle_psc(icol, k), - chlorine_loading, config_.linoz.psc_T, vmr[o3_ndx], + chlorine_loading, config.linoz.psc_T, vmr[o3_ndx], do3_linoz, do3_linoz_psc, ss_o3, o3col_du_diag, o3clim_linoz_diag, zenith_angle_degrees); // update source terms above the ozone decay threshold - if (k > nlev_ - config_.linoz.o3_lbl - 1) { + if (k > nlev - config.linoz.o3_lbl - 1) { Real do3_mass; // diagnostic, not needed - mam4::lin_strat_chem::lin_strat_sfcsink_kk(dt, pdel, vmr[o3_ndx], config_.linoz.o3_sfc, - config_.linoz.o3_tau, do3_mass); + mam4::lin_strat_chem::lin_strat_sfcsink_kk(dt, pdel, vmr[o3_ndx], config.linoz.o3_sfc, + config.linoz.o3_tau, do3_mass); } // ... check for negative values and reset to zero diff --git a/components/eamxx/src/physics/mam/eamxx_mam_optics_process_interface.cpp b/components/eamxx/src/physics/mam/eamxx_mam_optics_process_interface.cpp index a55bf84ca17..442a6500f21 100644 --- a/components/eamxx/src/physics/mam/eamxx_mam_optics_process_interface.cpp +++ b/components/eamxx/src/physics/mam/eamxx_mam_optics_process_interface.cpp @@ -140,7 +140,7 @@ void MAMOptics::initialize_impl(const RunType run_type) { wet_atm_.nc = get_field_in("nc").get_view(); wet_atm_.qi = get_field_in("qi").get_view(); wet_atm_.ni = get_field_in("ni").get_view(); - wet_atm_.omega = get_field_in("omega").get_view(); + constexpr int ntot_amode = mam4::AeroConfig::num_modes(); @@ -153,6 +153,7 @@ void MAMOptics::initialize_impl(const RunType run_type) { .get_view(); // FIXME: tot or liq? dry_atm_.pblh = get_field_in("pbl_height").get_view(); dry_atm_.phis = get_field_in("phis").get_view(); + dry_atm_.omega = get_field_in("omega").get_view(); dry_atm_.z_mid = buffer_.z_mid; dry_atm_.dz = buffer_.dz; diff --git a/components/eamxx/src/physics/mam/mam_coupling.hpp b/components/eamxx/src/physics/mam/mam_coupling.hpp index 278e97c033e..3a001146882 100644 --- a/components/eamxx/src/physics/mam/mam_coupling.hpp +++ b/components/eamxx/src/physics/mam/mam_coupling.hpp @@ -21,6 +21,7 @@ using view_2d = typename KT::template view_2d; using view_3d = typename KT::template view_3d; using const_view_1d = typename KT::template view_1d; using const_view_2d = typename KT::template view_2d; +using const_view_3d = typename KT::template view_3d; using complex_view_3d = typename KT::template view_3d>; using complex_view_2d = typename KT::template view_2d>; @@ -76,7 +77,7 @@ constexpr int num_aero_tracers() { // Given a MAM aerosol mode index, returns a string denoting the symbolic // name of the mode. -inline +KOKKOS_INLINE_FUNCTION const char* aero_mode_name(const int mode) { static const char *mode_names[num_aero_modes()] = { "1", @@ -87,9 +88,25 @@ const char* aero_mode_name(const int mode) { return mode_names[mode]; } +// Given a MAM aerosol species ID, returns a string denoting the symbolic +// name of the species. +KOKKOS_INLINE_FUNCTION +const char* aero_species_name(const int species_id) { + static const char *species_names[num_aero_species()] = { + "soa", + "so4", + "pom", + "bc", + "nacl", + "dst", + "mom", + }; + return species_names[species_id]; +} + // Given a MAM aerosol-related gas ID, returns a string denoting the symbolic // name of the gas species. -inline +KOKKOS_INLINE_FUNCTION const char* gas_species_name(const int gas_id) { static const char *species_names[num_aero_gases()] = { "O3", @@ -110,14 +127,14 @@ constexpr int max_field_name_len() { return 128; } -inline +KOKKOS_INLINE_FUNCTION size_t gpu_strlen(const char* s) { size_t l = 0; while (s[l]) ++l; return l; } -inline +KOKKOS_INLINE_FUNCTION void concat_2_strings(const char *s1, const char *s2, char *concatted) { size_t len1 = gpu_strlen(s1); for (size_t i = 0; i < len1; ++i) @@ -128,7 +145,7 @@ void concat_2_strings(const char *s1, const char *s2, char *concatted) { concatted[len1+len2] = 0; } -inline +KOKKOS_INLINE_FUNCTION void concat_3_strings(const char *s1, const char *s2, const char *s3, char *concatted) { size_t len1 = gpu_strlen(s1); for (size_t i = 0; i < len1; ++i) @@ -142,31 +159,31 @@ void concat_3_strings(const char *s1, const char *s2, const char *s3, char *conc concatted[len1+len2+len3] = 0; } -inline +KOKKOS_INLINE_FUNCTION char* int_aero_nmr_names(int mode) { static char int_aero_nmr_names_[num_aero_modes()][max_field_name_len()] = {}; return int_aero_nmr_names_[mode]; } -inline +KOKKOS_INLINE_FUNCTION char* cld_aero_nmr_names(int mode) { static char cld_aero_nmr_names_[num_aero_modes()][max_field_name_len()] = {}; return cld_aero_nmr_names_[mode]; } -inline +KOKKOS_INLINE_FUNCTION char* int_aero_mmr_names(int mode, int species) { static char int_aero_mmr_names_[num_aero_modes()][num_aero_species()][max_field_name_len()] = {}; return int_aero_mmr_names_[mode][species]; } -inline +KOKKOS_INLINE_FUNCTION char* cld_aero_mmr_names(int mode, int species) { static char cld_aero_mmr_names_[num_aero_modes()][num_aero_species()][max_field_name_len()] = {}; return cld_aero_mmr_names_[mode][species]; } -inline +KOKKOS_INLINE_FUNCTION char* gas_mmr_names(int gas_id) { static char gas_mmr_names_[num_aero_gases()][max_field_name_len()] = {}; return gas_mmr_names_[gas_id]; @@ -176,7 +193,7 @@ char* gas_mmr_names(int gas_id) { // Given a MAM aerosol mode index, returns the name of the related interstitial // modal number mixing ratio field in EAMxx ("num_a<1-based-mode-index>") -inline +KOKKOS_INLINE_FUNCTION const char* int_aero_nmr_field_name(const int mode) { if (!int_aero_nmr_names(mode)[0]) { concat_2_strings("num_a", aero_mode_name(mode), int_aero_nmr_names(mode)); @@ -186,7 +203,7 @@ const char* int_aero_nmr_field_name(const int mode) { // Given a MAM aerosol mode index, returns the name of the related cloudborne // modal number mixing ratio field in EAMxx ("num_c<1-based-mode-index>>") -inline +KOKKOS_INLINE_FUNCTION const char* cld_aero_nmr_field_name(const int mode) { if (!cld_aero_nmr_names(mode)[0]) { concat_2_strings("num_c", aero_mode_name(mode), cld_aero_nmr_names(mode)); @@ -199,13 +216,12 @@ const char* cld_aero_nmr_field_name(const int mode) { // field in EAMxx. The form of the field name is "_a<1-based-mode-index>". // If the desired species is not present within the desire mode, returns a blank // string (""). -inline +KOKKOS_INLINE_FUNCTION const char* int_aero_mmr_field_name(const int mode, const int species) { if (!int_aero_mmr_names(mode, species)[0]) { const auto aero_id = mam4::mode_aero_species(mode, species); if (aero_id != mam4::AeroId::None) { - auto aerosol_species_name =mam4::aero_id_short_name(aero_id); - concat_3_strings(aerosol_species_name.c_str(), + concat_3_strings(aero_species_name(static_cast(aero_id)), "_a", aero_mode_name(mode), int_aero_mmr_names(mode, species)); } @@ -218,13 +234,12 @@ const char* int_aero_mmr_field_name(const int mode, const int species) { // field in EAMxx. The form of the field name is "_c<1-based-mode-index>". // If the desired species is not present within the desire mode, returns a blank // string (""). -inline +KOKKOS_INLINE_FUNCTION const char* cld_aero_mmr_field_name(const int mode, const int species) { if (!cld_aero_mmr_names(mode, species)[0]) { const auto aero_id = mam4::mode_aero_species(mode, species); if (aero_id != mam4::AeroId::None) { - auto aerosol_species_name =mam4::aero_id_short_name(aero_id); - concat_3_strings(aerosol_species_name.c_str(), + concat_3_strings(aero_species_name(static_cast(aero_id)), "_c", aero_mode_name(mode), cld_aero_mmr_names(mode, species)); } @@ -234,7 +249,7 @@ const char* cld_aero_mmr_field_name(const int mode, const int species) { // Given a MAM aerosol-related gas identifier, returns the name of its mass // mixing ratio field in EAMxx -inline +KOKKOS_INLINE_FUNCTION const char* gas_mmr_field_name(const int gas) { return const_cast(gas_species_name(gas)); } @@ -247,7 +262,6 @@ struct WetAtmosphere { const_view_2d nc; // wet cloud liquid water number mixing ratio [# / kg moist air] const_view_2d qi; // wet cloud ice water mass mixing ratio [kg cloud ice water / kg moist air] const_view_2d ni; // wet cloud ice water number mixing ratio [# / kg moist air] - const_view_2d omega; // vertical pressure velocity [Pa/s] }; // This type stores multi-column views related to the dry atmospheric state @@ -265,11 +279,12 @@ struct DryAtmosphere { view_2d z_iface; // height at layer interfaces [m] view_2d dz; // layer thickness [m] const_view_2d p_del; // hydrostatic "pressure thickness" at grid interfaces [Pa] - const_view_2d p_int; // total pressure at grid interfaces [Pa] + const_view_2d p_int; // total pressure at grid interfaces [Pa] const_view_2d cldfrac; // cloud fraction [-] view_2d w_updraft; // updraft velocity [m/s] const_view_1d pblh; // planetary boundary layer height [m] const_view_1d phis; // surface geopotential [m2/s2] + const_view_2d omega; // vertical pressure velocity [Pa/s] }; // This type stores aerosol number and mass mixing ratios evolved by MAM. It @@ -584,9 +599,12 @@ void compute_vertical_layer_heights(const Team& team, const auto p_mid = ekat::subview(dry_atm.p_mid, column_index); const auto T_mid = ekat::subview(dry_atm.T_mid, column_index); const auto pseudo_density = ekat::subview(dry_atm.p_del, column_index); - PF::calculate_dz(team, pseudo_density, p_mid, T_mid, qv, dz); + // NOTE: we are using dry qv. Does calculate_dz require dry or wet? + PF::calculate_dz(team, pseudo_density, p_mid, T_mid, qv, // inputs + dz);//output team.team_barrier(); - PF::calculate_z_int(team, mam4::nlev, dz, dry_atm.z_surf, z_iface); + PF::calculate_z_int(team, mam4::nlev, dz, dry_atm.z_surf, //inputs + z_iface); //output team.team_barrier(); // likely necessary to have z_iface up to date PF::calculate_z_mid(team, mam4::nlev, z_iface, z_mid); } @@ -605,8 +623,9 @@ void compute_updraft_velocities(const Team& team, constexpr int nlev = mam4::nlev; int i = column_index; Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev), [&] (const int k) { - const auto rho = PF::calculate_density(dry_atm.p_del(i,k), dry_atm.dz(i,k)); - dry_atm.w_updraft(i,k) = PF::calculate_vertical_velocity(wet_atm.omega(i,k), rho); + dry_atm.dz(i,k) = PF::calculate_dz(dry_atm.p_del(i,k), dry_atm.p_mid(i,k), dry_atm.T_mid(i,k), wet_atm.qv(i,k)); + const auto rho = PF::calculate_density(dry_atm.p_del(i,k), dry_atm.dz(i,k)); + dry_atm.w_updraft(i,k) = PF::calculate_vertical_velocity(dry_atm.omega(i,k), rho); }); } @@ -705,6 +724,23 @@ void compute_wet_mixing_ratios(const Team& team, }); } +// Scream (or EAMxx) can sometimes extend views beyond model levels (nlev) as it uses +// "packs". Following function copies a 2d view till model levels +inline +void copy_view_lev_slice(haero::ThreadTeamPolicy team_policy, //inputs + const_view_2d &inp_view, //input view to copy + const int dim, //dimension till view should be copied + view_2d &out_view) { //output view + + Kokkos::parallel_for( + team_policy, KOKKOS_LAMBDA(const haero::ThreadTeam &team) { + const int icol = team.league_rank(); + Kokkos::parallel_for(Kokkos::TeamThreadRange(team, dim), [&](int kk) { + out_view(icol, kk) = inp_view(icol, kk); + }); + }); + } + // Because CUDA C++ doesn't allow us to declare and use constants outside of // KOKKOS_INLINE_FUNCTIONS, we define this macro that allows us to (re)define // these constants where needed within two such functions so we don't define @@ -885,6 +921,7 @@ void transfer_work_arrays_to_prognostics(const Real q[gas_pcnst()], progs.q_aero_i[m][a](k) = q[i]; progs.q_aero_c[m][a](k) = qqcw[i]; } else { // constituent is a modal number mixing ratio + int m = static_cast(mode_index); progs.n_mode_i[m](k) = q[i]; progs.n_mode_c[m](k) = qqcw[i]; } diff --git a/components/eamxx/src/physics/register_physics.hpp b/components/eamxx/src/physics/register_physics.hpp index 329b7022756..36afabf8de9 100644 --- a/components/eamxx/src/physics/register_physics.hpp +++ b/components/eamxx/src/physics/register_physics.hpp @@ -26,6 +26,7 @@ #ifdef EAMXX_HAS_MAM #include "physics/mam/eamxx_mam_microphysics_process_interface.hpp" #include "physics/mam/eamxx_mam_optics_process_interface.hpp" +#include "physics/mam/eamxx_mam_aci_process_interface.hpp" #endif #ifdef EAMXX_HAS_COSP #include "physics/cosp/eamxx_cosp.hpp" @@ -62,6 +63,7 @@ inline void register_physics () { #ifdef EAMXX_HAS_MAM proc_factory.register_product("mam4_micro",&create_atmosphere_process); proc_factory.register_product("mam4_optics",&create_atmosphere_process); + proc_factory.register_product("mam4_aci",&create_atmosphere_process); #endif #ifdef EAMXX_HAS_COSP proc_factory.register_product("Cosp",&create_atmosphere_process); diff --git a/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.cpp b/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.cpp index 6204421c1ea..15d43c82f7c 100644 --- a/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.cpp +++ b/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.cpp @@ -85,6 +85,9 @@ void SHOCMacrophysics::set_grids(const std::shared_ptr grids // Output variables add_field("pbl_height", scalar2d , m, grid_name); add_field("inv_qc_relvar", scalar3d_mid, pow(kg/kg,2), grid_name, ps); + add_field("eddy_diff_heat", scalar3d_mid, m2/s, grid_name, ps); + add_field("w_variance", scalar3d_mid, m2/s2, grid_name, ps); + add_field("cldfrac_liq_prev", scalar3d_mid, nondim, grid_name, ps); // Tracer group add_group("tracers", grid_name, ps, Bundling::Required); @@ -178,7 +181,7 @@ void SHOCMacrophysics::init_buffers(const ATMBufferManager &buffer_manager) &m_buffer.inv_exner, &m_buffer.thlm, &m_buffer.qw, &m_buffer.dse, &m_buffer.tke_copy, &m_buffer.qc_copy, &m_buffer.shoc_ql2, &m_buffer.shoc_mix, &m_buffer.isotropy, &m_buffer.w_sec, &m_buffer.wqls_sec, &m_buffer.brunt #ifdef SCREAM_SMALL_KERNELS - , &m_buffer.rho_zt, &m_buffer.shoc_qv, &m_buffer.tabs, &m_buffer.dz_zt, &m_buffer.tkh + , &m_buffer.rho_zt, &m_buffer.shoc_qv, &m_buffer.tabs, &m_buffer.dz_zt #endif }; @@ -250,6 +253,7 @@ void SHOCMacrophysics::initialize_impl (const RunType run_type) const auto& qv = get_field_out("qv").get_view(); const auto& tke = get_field_out("tke").get_view(); const auto& cldfrac_liq = get_field_out("cldfrac_liq").get_view(); + const auto& cldfrac_liq_prev = get_field_out("cldfrac_liq_prev").get_view(); const auto& sgs_buoy_flux = get_field_out("sgs_buoy_flux").get_view(); const auto& tk = get_field_out("eddy_diff_mom").get_view(); const auto& inv_qc_relvar = get_field_out("inv_qc_relvar").get_view(); @@ -294,7 +298,7 @@ void SHOCMacrophysics::initialize_impl (const RunType run_type) T_mid,p_mid,p_int,pseudo_density,omega,phis,surf_sens_flux,surf_evap, surf_mom_flux,qtracers,qv,qc,qc_copy,tke,tke_copy,z_mid,z_int, dse,rrho,rrho_i,thv,dz,zt_grid,zi_grid,wpthlp_sfc,wprtp_sfc,upwp_sfc,vpwp_sfc, - wtracer_sfc,wm_zt,inv_exner,thlm,qw); + wtracer_sfc,wm_zt,inv_exner,thlm,qw, cldfrac_liq, cldfrac_liq_prev); // Input Variables: input.zt_grid = shoc_preprocess.zt_grid; @@ -327,11 +331,12 @@ void SHOCMacrophysics::initialize_impl (const RunType run_type) // Output Variables output.pblh = get_field_out("pbl_height").get_view(); output.shoc_ql2 = shoc_ql2; + output.tkh = get_field_out("eddy_diff_heat").get_view(); // Ouput (diagnostic) history_output.shoc_mix = m_buffer.shoc_mix; history_output.isotropy = m_buffer.isotropy; - history_output.w_sec = m_buffer.w_sec; + history_output.w_sec = get_field_out("w_variance").get_view(); history_output.thl_sec = m_buffer.thl_sec; history_output.qw_sec = m_buffer.qw_sec; history_output.qwthl_sec = m_buffer.qwthl_sec; @@ -364,7 +369,6 @@ void SHOCMacrophysics::initialize_impl (const RunType run_type) temporaries.tabs = m_buffer.tabs; temporaries.dz_zt = m_buffer.dz_zt; temporaries.dz_zi = m_buffer.dz_zi; - temporaries.tkh = m_buffer.tkh; #endif shoc_postprocess.set_variables(m_num_cols,m_num_levs,m_num_tracers, diff --git a/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp b/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp index d3f5b05f88d..1bc0d14afe0 100644 --- a/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp +++ b/components/eamxx/src/physics/shoc/eamxx_shoc_process_interface.hpp @@ -83,6 +83,9 @@ class SHOCMacrophysics : public scream::AtmosphereProcess Kokkos::parallel_for(Kokkos::TeamVectorRange(team, nlev_packs), [&] (const Int& k) { + + cldfrac_liq_prev(i,k)=cldfrac_liq(i,k); + const auto range = ekat::range(k*Spack::n); const Smask in_nlev_range = (range < nlev); @@ -198,6 +201,8 @@ class SHOCMacrophysics : public scream::AtmosphereProcess view_2d thlm; view_2d qw; view_2d cloud_frac; + view_2d cldfrac_liq; + view_2d cldfrac_liq_prev; // Assigning local variables void set_variables(const int ncol_, const int nlev_, const int num_qtracers_, @@ -213,7 +218,8 @@ class SHOCMacrophysics : public scream::AtmosphereProcess const view_2d& dse_, const view_2d& rrho_, const view_2d& rrho_i_, const view_2d& thv_, const view_2d& dz_,const view_2d& zt_grid_,const view_2d& zi_grid_, const view_1d& wpthlp_sfc_, const view_1d& wprtp_sfc_,const view_1d& upwp_sfc_,const view_1d& vpwp_sfc_, const view_2d& wtracer_sfc_, - const view_2d& wm_zt_,const view_2d& inv_exner_,const view_2d& thlm_,const view_2d& qw_) + const view_2d& wm_zt_,const view_2d& inv_exner_,const view_2d& thlm_,const view_2d& qw_, + const view_2d& cldfrac_liq_, const view_2d& cldfrac_liq_prev_) { ncol = ncol_; nlev = nlev_; @@ -254,6 +260,8 @@ class SHOCMacrophysics : public scream::AtmosphereProcess inv_exner = inv_exner_; thlm = thlm_; qw = qw_; + cldfrac_liq=cldfrac_liq_; + cldfrac_liq_prev=cldfrac_liq_prev_; } // set_variables }; // SHOCPreprocess /* --------------------------------------------------------------------------------------------*/ diff --git a/components/eamxx/src/physics/shoc/impl/shoc_main_impl.hpp b/components/eamxx/src/physics/shoc/impl/shoc_main_impl.hpp index 13a93992cff..120ed505f26 100644 --- a/components/eamxx/src/physics/shoc/impl/shoc_main_impl.hpp +++ b/components/eamxx/src/physics/shoc/impl/shoc_main_impl.hpp @@ -119,6 +119,7 @@ void Functions::shoc_main_internal( // Output Variables Scalar& pblh, const uview_1d& shoc_ql2, + const uview_1d& tkh, // Diagnostic Output Variables const uview_1d& shoc_mix, const uview_1d& w_sec, @@ -137,10 +138,10 @@ void Functions::shoc_main_internal( { // Define temporary variables - uview_1d rho_zt, shoc_qv, shoc_tabs, dz_zt, dz_zi, tkh; - workspace.template take_many_and_reset<6>( - {"rho_zt", "shoc_qv", "shoc_tabs", "dz_zt", "dz_zi", "tkh"}, - {&rho_zt, &shoc_qv, &shoc_tabs, &dz_zt, &dz_zi, &tkh}); + uview_1d rho_zt, shoc_qv, shoc_tabs, dz_zt, dz_zi; + workspace.template take_many_and_reset<5>( + {"rho_zt", "shoc_qv", "shoc_tabs", "dz_zt", "dz_zi"}, + {&rho_zt, &shoc_qv, &shoc_tabs, &dz_zt, &dz_zi}); // Local scalars Scalar se_b{0}, ke_b{0}, wv_b{0}, wl_b{0}, @@ -311,8 +312,8 @@ void Functions::shoc_main_internal( pblh); // Output // Release temporary variables from the workspace - workspace.template release_many_contiguous<6>( - {&rho_zt, &shoc_qv, &shoc_tabs, &dz_zt, &dz_zi, &tkh}); + workspace.template release_many_contiguous<5>( + {&rho_zt, &shoc_qv, &shoc_tabs, &dz_zt, &dz_zi}); } #else template @@ -371,6 +372,7 @@ void Functions::shoc_main_internal( // Output Variables const view_1d& pblh, const view_2d& shoc_ql2, + const view_2d& tkh, // Diagnostic Output Variables const view_2d& shoc_mix, const view_2d& w_sec, @@ -404,8 +406,7 @@ void Functions::shoc_main_internal( const view_2d& shoc_qv, const view_2d& shoc_tabs, const view_2d& dz_zt, - const view_2d& dz_zi, - const view_2d& tkh) + const view_2d& dz_zi) { // Scalarize some views for single entry access const auto s_thetal = ekat::scalarize(thetal); @@ -643,6 +644,7 @@ Int Functions::shoc_main( const auto shoc_cldfrac_s = ekat::subview(shoc_input_output.shoc_cldfrac, i); const auto shoc_ql_s = ekat::subview(shoc_input_output.shoc_ql, i); const auto shoc_ql2_s = ekat::subview(shoc_output.shoc_ql2, i); + const auto tkh_s = ekat::subview(shoc_output.tkh, i); const auto shoc_mix_s = ekat::subview(shoc_history_output.shoc_mix, i); const auto w_sec_s = ekat::subview(shoc_history_output.w_sec, i); const auto thl_sec_s = ekat::subview(shoc_history_output.thl_sec, i); @@ -674,7 +676,7 @@ Int Functions::shoc_main( host_dse_s, tke_s, thetal_s, qw_s, u_wind_s, v_wind_s, // Input/Output wthv_sec_s, qtracers_s, tk_s, shoc_cldfrac_s, // Input/Output shoc_ql_s, // Input/Output - pblh_s, shoc_ql2_s, // Output + pblh_s, shoc_ql2_s, tkh_s, // Output shoc_mix_s, w_sec_s, thl_sec_s, qw_sec_s, qwthl_sec_s, // Diagnostic Output Variables wthl_sec_s, wqw_sec_s, wtke_sec_s, uw_sec_s, vw_sec_s, // Diagnostic Output Variables w3_s, wqls_sec_s, brunt_s, isotropy_s); // Diagnostic Output Variables @@ -698,7 +700,7 @@ Int Functions::shoc_main( shoc_input_output.host_dse, shoc_input_output.tke, shoc_input_output.thetal, shoc_input_output.qw, u_wind_s, v_wind_s, // Input/Output shoc_input_output.wthv_sec, shoc_input_output.qtracers, shoc_input_output.tk, shoc_input_output.shoc_cldfrac, // Input/Output shoc_input_output.shoc_ql, // Input/Output - shoc_output.pblh, shoc_output.shoc_ql2, // Output + shoc_output.pblh, shoc_output.shoc_ql2, shoc_output.tkh, // Output shoc_history_output.shoc_mix, shoc_history_output.w_sec, shoc_history_output.thl_sec, shoc_history_output.qw_sec, shoc_history_output.qwthl_sec, // Diagnostic Output Variables shoc_history_output.wthl_sec, shoc_history_output.wqw_sec, shoc_history_output.wtke_sec, shoc_history_output.uw_sec, shoc_history_output.vw_sec, // Diagnostic Output Variables shoc_history_output.w3, shoc_history_output.wqls_sec, shoc_history_output.brunt, shoc_history_output.isotropy, // Diagnostic Output Variables @@ -707,7 +709,7 @@ Int Functions::shoc_main( shoc_temporaries.se_a, shoc_temporaries.ke_a, shoc_temporaries.wv_a, shoc_temporaries.wl_a, shoc_temporaries.ustar, shoc_temporaries.kbfs, shoc_temporaries.obklen, shoc_temporaries.ustar2, shoc_temporaries.wstar, shoc_temporaries.rho_zt, shoc_temporaries.shoc_qv, - shoc_temporaries.tabs, shoc_temporaries.dz_zt, shoc_temporaries.dz_zi, shoc_temporaries.tkh); + shoc_temporaries.tabs, shoc_temporaries.dz_zt, shoc_temporaries.dz_zi); #endif auto finish = std::chrono::steady_clock::now(); diff --git a/components/eamxx/src/physics/shoc/shoc_functions.hpp b/components/eamxx/src/physics/shoc/shoc_functions.hpp index d02d498c5a4..74602c688c2 100644 --- a/components/eamxx/src/physics/shoc/shoc_functions.hpp +++ b/components/eamxx/src/physics/shoc/shoc_functions.hpp @@ -161,6 +161,8 @@ struct Functions view_1d pblh; // cloud liquid mixing ratio variance [kg^2/kg^2] view_2d shoc_ql2; + // eddy coefficient for heat [m2/s] + view_2d tkh; }; // This struct stores output views for SHOC diagnostics for shoc_main. @@ -888,6 +890,7 @@ struct Functions // Output Variables Scalar& pblh, const uview_1d& shoc_ql2, + const uview_1d& tkh, // Diagnostic Output Variables const uview_1d& shoc_mix, const uview_1d& w_sec, @@ -959,6 +962,7 @@ struct Functions // Output Variables const view_1d& pblh, const view_2d& shoc_ql2, + const view_2d& tkh, // Diagnostic Output Variables const view_2d& shoc_mix, const view_2d& w_sec, @@ -992,8 +996,7 @@ struct Functions const view_2d& shoc_qv, const view_2d& tabs, const view_2d& dz_zt, - const view_2d& dz_zi, - const view_2d& tkh); + const view_2d& dz_zi); #endif // Return microseconds elapsed diff --git a/components/eamxx/src/physics/shoc/shoc_functions_f90.cpp b/components/eamxx/src/physics/shoc/shoc_functions_f90.cpp index 0e4acceb1a7..dc45e2456aa 100644 --- a/components/eamxx/src/physics/shoc/shoc_functions_f90.cpp +++ b/components/eamxx/src/physics/shoc/shoc_functions_f90.cpp @@ -2743,8 +2743,6 @@ Int shoc_main_f(Int shcol, Int nlev, Int nlevi, Real dtime, Int nadv, Int npbl, Real* qw_sec, Real* qwthl_sec, Real* wthl_sec, Real* wqw_sec, Real* wtke_sec, Real* uw_sec, Real* vw_sec, Real* w3, Real* wqls_sec, Real* brunt, Real* shoc_ql2) { - // tkh is a local variable in C++ impl - (void)tkh; using SHF = Functions; @@ -2758,7 +2756,7 @@ Int shoc_main_f(Int shcol, Int nlev, Int nlevi, Real dtime, Int nadv, Int npbl, // Initialize Kokkos views, sync to device static constexpr Int num_1d_arrays = 7; - static constexpr Int num_2d_arrays = 34; + static constexpr Int num_2d_arrays = 35; static constexpr Int num_3d_arrays = 1; std::vector temp_1d_d(num_1d_arrays); @@ -2768,27 +2766,27 @@ Int shoc_main_f(Int shcol, Int nlev, Int nlevi, Real dtime, Int nadv, Int npbl, std::vector dim1_2d_sizes = {shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, - shcol, shcol, shcol, shcol, + shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol, shcol}; std::vector dim2_2d_sizes = {nlev, nlevi, nlev, nlevi, nlev, nlev, nlev, num_qtracers, nlev, nlev, nlev, nlev, nlev, nlev, nlev, - nlev, nlev, nlev, nlev, + nlev, nlev, nlev, nlev, nlev, nlev, nlev, nlev, nlevi, nlevi, nlevi, nlevi, nlevi, nlevi, nlevi, nlevi, nlevi, nlev, nlev, nlev}; std::vector ptr_array_1d = {host_dx, host_dy, wthl_sfc, wqw_sfc, uw_sfc, vw_sfc, phis}; - std::vector ptr_array_2d = {zt_grid, zi_grid, pres, presi, pdel, - thv, w_field, wtracer_sfc, inv_exner, host_dse, - tke, thetal, qw, u_wind, v_wind, - wthv_sec, tk, shoc_cldfrac, shoc_ql, - shoc_ql2, shoc_mix, w_sec, thl_sec, qw_sec, - qwthl_sec, wthl_sec, wqw_sec, wtke_sec, uw_sec, - vw_sec, w3, wqls_sec, brunt, isotropy}; + std::vector ptr_array_2d = {zt_grid, zi_grid, pres, presi, pdel, + thv, w_field, wtracer_sfc, inv_exner, host_dse, + tke, thetal, qw, u_wind, v_wind, + wthv_sec, tk, shoc_cldfrac, shoc_ql, shoc_ql2, + tkh, shoc_mix, w_sec, thl_sec, qw_sec, + qwthl_sec, wthl_sec, wqw_sec, wtke_sec, uw_sec, + vw_sec, w3, wqls_sec, brunt, isotropy}; ScreamDeepCopy::copy_to_device(ptr_array_1d, shcol, temp_1d_d); ekat::host_to_device(ptr_array_2d, dim1_2d_sizes, dim2_2d_sizes, temp_2d_d, true); @@ -2827,6 +2825,7 @@ Int shoc_main_f(Int shcol, Int nlev, Int nlevi, Real dtime, Int nadv, Int npbl, shoc_cldfrac_d(temp_2d_d[index_counter++]), shoc_ql_d (temp_2d_d[index_counter++]), shoc_ql2_d (temp_2d_d[index_counter++]), + tkh_d (temp_2d_d[index_counter++]), shoc_mix_d (temp_2d_d[index_counter++]), w_sec_d (temp_2d_d[index_counter++]), thl_sec_d (temp_2d_d[index_counter++]), @@ -2879,7 +2878,7 @@ Int shoc_main_f(Int shcol, Int nlev, Int nlevi, Real dtime, Int nadv, Int npbl, SHF::SHOCInputOutput shoc_input_output{host_dse_d, tke_d, thetal_d, qw_d, horiz_wind_d, wthv_sec_d, qtracers_cxx_d, tk_d, shoc_cldfrac_d, shoc_ql_d}; - SHF::SHOCOutput shoc_output{pblh_d, shoc_ql2_d}; + SHF::SHOCOutput shoc_output{pblh_d, shoc_ql2_d, tkh_d}; SHF::SHOCHistoryOutput shoc_history_output{shoc_mix_d, w_sec_d, thl_sec_d, qw_sec_d, qwthl_sec_d, wthl_sec_d, wqw_sec_d, wtke_sec_d, uw_sec_d, vw_sec_d, w3_d, wqls_sec_d, @@ -2909,12 +2908,11 @@ Int shoc_main_f(Int shcol, Int nlev, Int nlevi, Real dtime, Int nadv, Int npbl, shoc_qv ("shoc_qv", shcol, nlevi_packs), tabs ("shoc_tabs", shcol, nlev_packs), dz_zt ("dz_zt", shcol, nlevi_packs), - dz_zi ("dz_zi", shcol, nlevi_packs), - tkhv ("tkh", shcol, nlevi_packs); + dz_zi ("dz_zi", shcol, nlevi_packs); SHF::SHOCTemporaries shoc_temporaries{ se_b, ke_b, wv_b, wl_b, se_a, ke_a, wv_a, wl_a, ustar, kbfs, obklen, ustar2, wstar, - rho_zt, shoc_qv, tabs, dz_zt, dz_zi, tkhv}; + rho_zt, shoc_qv, tabs, dz_zt, dz_zi}; #endif // Create local workspace diff --git a/components/eamxx/tests/CMakeLists.txt b/components/eamxx/tests/CMakeLists.txt index 7a7a9e82149..287fdf25af1 100644 --- a/components/eamxx/tests/CMakeLists.txt +++ b/components/eamxx/tests/CMakeLists.txt @@ -63,7 +63,7 @@ if (NOT DEFINED ENV{SCREAM_FAKE_ONLY}) set(EAMxx_tests_IC_FILE_72lev "screami_unit_tests_ne2np4L72_20220822.nc") set(EAMxx_tests_IC_FILE_128lev "screami_unit_tests_ne2np4L128_20220822.nc") set(EAMxx_tests_TOPO_FILE "USGS-gtopo30_ne2np4pg2_x6t_20230331.nc") - set(EAMxx_tests_IC_FILE_MAM4xx_72lev "scream_unit_tests_aerosol_optics_ne2np4L72_20220822.nc") + set(EAMxx_tests_IC_FILE_MAM4xx_72lev "screami_unit_tests_mam4xx_ne2np4L72_20240329.nc") # Testing individual atm processes add_subdirectory(single-process) diff --git a/components/eamxx/tests/multi-process/physics_only/CMakeLists.txt b/components/eamxx/tests/multi-process/physics_only/CMakeLists.txt index 61b224f71ae..b9c13234796 100644 --- a/components/eamxx/tests/multi-process/physics_only/CMakeLists.txt +++ b/components/eamxx/tests/multi-process/physics_only/CMakeLists.txt @@ -4,6 +4,10 @@ if (SCREAM_DOUBLE_PRECISION) add_subdirectory(shoc_cld_spa_p3_rrtmgp) if (SCREAM_ENABLE_MAM) add_subdirectory(mam/optics_rrtmgp) + add_subdirectory(mam/shoc_mam4_aci) + add_subdirectory(mam/shoc_cldfrac_mam4_aci_p3) + add_subdirectory(mam/shoc_cldfrac_mam4_aci_p3_rrtmgp) + add_subdirectory(mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp) endif() endif() diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/CMakeLists.txt b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/CMakeLists.txt new file mode 100644 index 00000000000..856ef3ca162 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/CMakeLists.txt @@ -0,0 +1,41 @@ +INCLUDE (ScreamUtils) + +set (TEST_BASE_NAME shoc_cldfrac_mam4_aci_p3) +set (FIXTURES_BASE_NAME ${TEST_BASE_NAME}_generate_output_nc_files) + +# Create the test +CreateADUnitTest(${TEST_BASE_NAME} + LIBS shoc cld_fraction p3 mam + LABELS shoc cld p3 physics mam4_aci + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + FIXTURES_SETUP_INDIVIDUAL ${FIXTURES_BASE_NAME} +) + +# Set AD configurable options +set (ATM_TIME_STEP 1800) +SetVarDependingOnTestSize(NUM_STEPS 2 5 48) # 1h 4h 24h +set (RUN_T0 2021-10-12-45000) + +# Determine num subcycles needed to keep shoc dt<=300s +set (SHOC_MAX_DT 300) +math (EXPR MAC_MIC_SUBCYCLES "(${ATM_TIME_STEP} + ${SHOC_MAX_DT} - 1) / ${SHOC_MAX_DT}") + +# Ensure test input files are present in the data dir +GetInputFile(scream/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev}) +GetInputFile(cam/topo/${EAMxx_tests_TOPO_FILE}) + +## Copy (and configure) yaml files needed by tests +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/input.yaml + ${CMAKE_CURRENT_BINARY_DIR}/input.yaml) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/output.yaml + ${CMAKE_CURRENT_BINARY_DIR}/output.yaml) + +# Compare output files produced by npX tests, to ensure they are bfb +include (CompareNCFiles) +CompareNCFilesFamilyMpi ( + TEST_BASE_NAME ${TEST_BASE_NAME} + FILE_META_NAME ${TEST_BASE_NAME}_output.INSTANT.nsteps_x${NUM_STEPS}.npMPIRANKS.${RUN_T0}.nc + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + LABELS shoc cld p3 physics PEM mam4_aci + META_FIXTURES_REQUIRED ${FIXTURES_BASE_NAME}_npMPIRANKS_omp1 +) diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/input.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/input.yaml new file mode 100644 index 00000000000..9c0ed4d9798 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/input.yaml @@ -0,0 +1,68 @@ +%YAML 1.1 +--- +driver_options: + atmosphere_dag_verbosity_level: 5 + +time_stepping: + time_step: ${ATM_TIME_STEP} + run_t0: ${RUN_T0} # YYYY-MM-DD-XXXXX + number_of_steps: ${NUM_STEPS} + +atmosphere_processes: + atm_procs_list: [mac_mic] + schedule_type: Sequential + mac_mic: + atm_procs_list: [shoc,CldFraction,mam4_aci,p3] + Type: Group + schedule_type: Sequential + number_of_subcycles: ${MAC_MIC_SUBCYCLES} + p3: + max_total_ni: 740.0e3 + do_prescribed_ccn: false + shoc: + lambda_low: 0.001 + lambda_high: 0.04 + lambda_slope: 2.65 + lambda_thresh: 0.02 + thl2tune: 1.0 + qw2tune: 1.0 + qwthl2tune: 1.0 + w2tune: 1.0 + length_fac: 0.5 + c_diag_3rd_mom: 7.0 + Ckh: 0.1 + Ckm: 0.1 + mam4_aci: + wsubmin: 0.001 + top_level_mam4xx: 6 + +grids_manager: + Type: Mesh Free + geo_data_source: IC_FILE + grids_names: [Physics GLL] + Physics GLL: + aliases: [Physics] + type: point_grid + number_of_global_columns: 218 + number_of_vertical_levels: 72 + +initial_conditions: + # The name of the file containing the initial conditions for this test. + Filename: ${SCREAM_DATA_DIR}/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev} + topography_filename: ${TOPO_DATA_DIR}/${EAMxx_tests_TOPO_FILE} + + #variable required for shoc + surf_sens_flux: 0.0 + surf_evap: 0.0 + + #variable required for p3 + precip_ice_surf_mass: 0.0 + precip_liq_surf_mass: 0.0 + + #variables required for mam4_aci + dgnum: 1e-3 + +# The parameters for I/O control +Scorpio: + output_yaml_files: ["output.yaml"] +... diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/output.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/output.yaml new file mode 100644 index 00000000000..2541ce8d9f8 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3/output.yaml @@ -0,0 +1,97 @@ +%YAML 1.1 +--- +filename_prefix: shoc_cldfrac_mam4_aci_p3_output +Averaging Type: Instant +Field Names: + # SHOC + - cldfrac_liq + - eddy_diff_mom + - horiz_winds + - sgs_buoy_flux + - tke + - inv_qc_relvar + - pbl_height + # CLD + - cldfrac_ice + - cldfrac_tot + # P3 + - bm + - nc + - ni + - nr + - qi + - qm + - qr + - T_prev_micro_step + - qv_prev_micro_step + - eff_radius_qc + - eff_radius_qi + - eff_radius_qr + - micro_liq_ice_exchange + - micro_vap_ice_exchange + - micro_vap_liq_exchange + - precip_ice_surf_mass + - precip_liq_surf_mass + - rainfrac + # SHOC + P3 + - qc + - qv + # SHOC + P3 + RRTMGP + - T_mid + - bc_a1 + - bc_a3 + - bc_a4 + - dst_a1 + - dst_a3 + - so4_a1 + - so4_a2 + - so4_a3 + - pom_a1 + - pom_a3 + - pom_a4 + - soa_a1 + - soa_a2 + - soa_a3 + - nacl_a1 + - nacl_a2 + - nacl_a3 + - mom_a1 + - mom_a2 + - mom_a3 + - mom_a4 + - num_a1 + - num_a2 + - num_a3 + - num_a4 + - bc_c1 + - bc_c3 + - bc_c4 + - dst_c1 + - dst_c3 + - so4_c1 + - so4_c2 + - so4_c3 + - pom_c1 + - pom_c3 + - pom_c4 + - soa_c1 + - soa_c2 + - soa_c3 + - nacl_c1 + - nacl_c2 + - nacl_c3 + - mom_c1 + - mom_c2 + - mom_c3 + - mom_c4 + - num_c1 + - num_c2 + - num_c3 + - num_c4 + - eddy_diff_heat + - nc_nuceat_tend + - ni_activated +output_control: + Frequency: ${NUM_STEPS} + frequency_units: nsteps +... diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/CMakeLists.txt b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/CMakeLists.txt new file mode 100644 index 00000000000..3f9aeae73fd --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/CMakeLists.txt @@ -0,0 +1,41 @@ +INCLUDE (ScreamUtils) + +set (TEST_BASE_NAME shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp) +set (FIXTURES_BASE_NAME ${TEST_BASE_NAME}_generate_output_nc_files) + +# Create the test +CreateADUnitTest(${TEST_BASE_NAME} + LIBS shoc cld_fraction p3 scream_rrtmgp mam + LABELS shoc cld p3 rrtmgp physics mam4_aci + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + FIXTURES_SETUP_INDIVIDUAL ${FIXTURES_BASE_NAME} +) + +# Set AD configurable options +set (ATM_TIME_STEP 1800) +SetVarDependingOnTestSize(NUM_STEPS 2 5 48) # 1h 4h 24h +set (RUN_T0 2021-10-12-45000) + +# Determine num subcycles needed to keep shoc dt<=300s +set (SHOC_MAX_DT 300) +math (EXPR MAC_MIC_SUBCYCLES "(${ATM_TIME_STEP} + ${SHOC_MAX_DT} - 1) / ${SHOC_MAX_DT}") + +# Ensure test input files are present in the data dir +GetInputFile(scream/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev}) +GetInputFile(cam/topo/${EAMxx_tests_TOPO_FILE}) + +## Copy (and configure) yaml files needed by tests +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/input.yaml + ${CMAKE_CURRENT_BINARY_DIR}/input.yaml) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/output.yaml + ${CMAKE_CURRENT_BINARY_DIR}/output.yaml) + +# Compare output files produced by npX tests, to ensure they are bfb +include (CompareNCFiles) +CompareNCFilesFamilyMpi ( + TEST_BASE_NAME ${TEST_BASE_NAME} + FILE_META_NAME ${TEST_BASE_NAME}_output.INSTANT.nsteps_x${NUM_STEPS}.npMPIRANKS.${RUN_T0}.nc + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + LABELS shoc cld p3 rrtmgp physics PEM mam4_aci + META_FIXTURES_REQUIRED ${FIXTURES_BASE_NAME}_npMPIRANKS_omp1 +) diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/input.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/input.yaml new file mode 100644 index 00000000000..20d6f142c52 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/input.yaml @@ -0,0 +1,90 @@ +%YAML 1.1 +--- +driver_options: + atmosphere_dag_verbosity_level: 5 + +time_stepping: + time_step: ${ATM_TIME_STEP} + run_t0: ${RUN_T0} # YYYY-MM-DD-XXXXX + number_of_steps: ${NUM_STEPS} + +atmosphere_processes: + atm_procs_list: [mac_mic, mam4_optics, rrtmgp] + schedule_type: Sequential + mac_mic: + atm_procs_list: [shoc,CldFraction,mam4_aci,p3] + Type: Group + schedule_type: Sequential + number_of_subcycles: ${MAC_MIC_SUBCYCLES} + mam4_aci: + wsubmin: 0.001 + top_level_mam4xx: 6 + p3: + max_total_ni: 740.0e3 + do_prescribed_ccn: false + shoc: + lambda_low: 0.001 + lambda_high: 0.04 + lambda_slope: 2.65 + lambda_thresh: 0.02 + thl2tune: 1.0 + qw2tune: 1.0 + qwthl2tune: 1.0 + w2tune: 1.0 + length_fac: 0.5 + c_diag_3rd_mom: 7.0 + Ckh: 0.1 + Ckm: 0.1 + mam4_optics: + mam4_mode1_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/mam4_mode1_rrtmg_aeronetdust_c20240206.nc + mam4_mode2_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/mam4_mode2_rrtmg_c20240206.nc + mam4_mode3_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/mam4_mode3_rrtmg_aeronetdust_c20240206.nc + mam4_mode4_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/mam4_mode4_rrtmg_c20240206.nc + mam4_water_refindex_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/water_refindex_rrtmg_c20240206.nc + mam4_soa_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/ocphi_rrtmg_c20240206.nc + mam4_dust_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/dust_aeronet_rrtmg_c20240206.nc + mam4_nacl_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/ssam_rrtmg_c20240206.nc + mam4_so4_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/sulfate_rrtmg_c20240206.nc + mam4_pom_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/ocpho_rrtmg_c20240206.nc + mam4_bc_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/bcpho_rrtmg_c20240206.nc + mam4_mom_physical_properties_file : ${SCREAM_DATA_DIR}/mam4xx/physprops/poly_rrtmg_c20240206.nc + rrtmgp: + column_chunk_size: 123 + active_gases: ["h2o", "co2", "o3", "n2o", "co" , "ch4", "o2", "n2"] + orbital_year: 1990 + do_aerosol_rad: false + rrtmgp_coefficients_file_sw: ${SCREAM_DATA_DIR}/init/rrtmgp-data-sw-g112-210809.nc + rrtmgp_coefficients_file_lw: ${SCREAM_DATA_DIR}/init/rrtmgp-data-lw-g128-210809.nc + rrtmgp_cloud_optics_file_sw: ${SCREAM_DATA_DIR}/init/rrtmgp-cloud-optics-coeffs-sw.nc + rrtmgp_cloud_optics_file_lw: ${SCREAM_DATA_DIR}/init/rrtmgp-cloud-optics-coeffs-lw.nc + +grids_manager: + Type: Mesh Free + geo_data_source: IC_FILE + grids_names: [Physics GLL] + Physics GLL: + aliases: [Physics] + type: point_grid + number_of_global_columns: 218 + number_of_vertical_levels: 72 + +initial_conditions: + # The name of the file containing the initial conditions for this test. + Filename: ${SCREAM_DATA_DIR}/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev} + topography_filename: ${TOPO_DATA_DIR}/${EAMxx_tests_TOPO_FILE} + + #variables required for shoc + surf_sens_flux: 0.0 + surf_evap: 0.0 + + #variables required for p3 + precip_ice_surf_mass: 0.0 + precip_liq_surf_mass: 0.0 + + #variables required for mam4_aci + dgnum: 1e-3 + +# The parameters for I/O control +Scorpio: + output_yaml_files: ["output.yaml"] +... diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/output.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/output.yaml new file mode 100644 index 00000000000..392e159f0a5 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp/output.yaml @@ -0,0 +1,110 @@ +%YAML 1.1 +--- +filename_prefix: shoc_cldfrac_mam4_aci_p3_mam4_optics_rrtmgp_output +Averaging Type: Instant +Field Names: + # SHOC + - cldfrac_liq + - eddy_diff_mom + - horiz_winds + - sgs_buoy_flux + - tke + - inv_qc_relvar + - pbl_height + # CLD + - cldfrac_ice + - cldfrac_tot + # P3 + - bm + - nc + - ni + - nr + - qi + - qm + - qr + - T_prev_micro_step + - qv_prev_micro_step + - eff_radius_qc + - eff_radius_qi + - eff_radius_qr + - micro_liq_ice_exchange + - micro_vap_ice_exchange + - micro_vap_liq_exchange + - precip_ice_surf_mass + - precip_liq_surf_mass + - rainfrac + # SHOC + P3 + - qc + - qv + # SHOC + P3 + RRTMGP + - T_mid + # RRTMGP + - sfc_alb_dif_nir + - sfc_alb_dif_vis + - sfc_alb_dir_nir + - sfc_alb_dir_vis + - LW_flux_dn + - LW_flux_up + - SW_flux_dn + - SW_flux_dn_dir + - SW_flux_up + - rad_heating_pdel + - sfc_flux_lw_dn + - sfc_flux_sw_net + - bc_a1 + - bc_a3 + - bc_a4 + - dst_a1 + - dst_a3 + - so4_a1 + - so4_a2 + - so4_a3 + - pom_a1 + - pom_a3 + - pom_a4 + - soa_a1 + - soa_a2 + - soa_a3 + - nacl_a1 + - nacl_a2 + - nacl_a3 + - mom_a1 + - mom_a2 + - mom_a3 + - mom_a4 + - num_a1 + - num_a2 + - num_a3 + - num_a4 + - bc_c1 + - bc_c3 + - bc_c4 + - dst_c1 + - dst_c3 + - so4_c1 + - so4_c2 + - so4_c3 + - pom_c1 + - pom_c3 + - pom_c4 + - soa_c1 + - soa_c2 + - soa_c3 + - nacl_c1 + - nacl_c2 + - nacl_c3 + - mom_c1 + - mom_c2 + - mom_c3 + - mom_c4 + - num_c1 + - num_c2 + - num_c3 + - num_c4 + - eddy_diff_heat + - nc_nuceat_tend + - ni_activated +output_control: + Frequency: ${NUM_STEPS} + frequency_units: nsteps +... diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/CMakeLists.txt b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/CMakeLists.txt new file mode 100644 index 00000000000..cdd564c73a8 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/CMakeLists.txt @@ -0,0 +1,41 @@ +INCLUDE (ScreamUtils) + +set (TEST_BASE_NAME shoc_cldfrac_mam4_aci_p3_rrtmgp) +set (FIXTURES_BASE_NAME ${TEST_BASE_NAME}_generate_output_nc_files) + +# Create the test +CreateADUnitTest(${TEST_BASE_NAME} + LIBS shoc cld_fraction p3 scream_rrtmgp mam + LABELS shoc cld p3 rrtmgp physics mam4_aci + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + FIXTURES_SETUP_INDIVIDUAL ${FIXTURES_BASE_NAME} +) + +# Set AD configurable options +set (ATM_TIME_STEP 1800) +SetVarDependingOnTestSize(NUM_STEPS 2 5 48) # 1h 4h 24h +set (RUN_T0 2021-10-12-45000) + +# Determine num subcycles needed to keep shoc dt<=300s +set (SHOC_MAX_DT 300) +math (EXPR MAC_MIC_SUBCYCLES "(${ATM_TIME_STEP} + ${SHOC_MAX_DT} - 1) / ${SHOC_MAX_DT}") + +# Ensure test input files are present in the data dir +GetInputFile(scream/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev}) +GetInputFile(cam/topo/${EAMxx_tests_TOPO_FILE}) + +## Copy (and configure) yaml files needed by tests +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/input.yaml + ${CMAKE_CURRENT_BINARY_DIR}/input.yaml) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/output.yaml + ${CMAKE_CURRENT_BINARY_DIR}/output.yaml) + +# Compare output files produced by npX tests, to ensure they are bfb +include (CompareNCFiles) +CompareNCFilesFamilyMpi ( + TEST_BASE_NAME ${TEST_BASE_NAME} + FILE_META_NAME ${TEST_BASE_NAME}_output.INSTANT.nsteps_x${NUM_STEPS}.npMPIRANKS.${RUN_T0}.nc + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + LABELS shoc cld p3 rrtmgp physics PEM mam4_aci + META_FIXTURES_REQUIRED ${FIXTURES_BASE_NAME}_npMPIRANKS_omp1 +) diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/input.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/input.yaml new file mode 100644 index 00000000000..7ec317eda92 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/input.yaml @@ -0,0 +1,83 @@ +%YAML 1.1 +--- +driver_options: + atmosphere_dag_verbosity_level: 5 + +time_stepping: + time_step: ${ATM_TIME_STEP} + run_t0: ${RUN_T0} # YYYY-MM-DD-XXXXX + number_of_steps: ${NUM_STEPS} + +atmosphere_processes: + atm_procs_list: [mac_mic,rrtmgp] + schedule_type: Sequential + mac_mic: + atm_procs_list: [shoc,CldFraction,mam4_aci,p3] + Type: Group + schedule_type: Sequential + number_of_subcycles: ${MAC_MIC_SUBCYCLES} + mam4_aci: + wsubmin: 0.001 + top_level_mam4xx: 6 + p3: + max_total_ni: 740.0e3 + do_prescribed_ccn: false + shoc: + lambda_low: 0.001 + lambda_high: 0.04 + lambda_slope: 2.65 + lambda_thresh: 0.02 + thl2tune: 1.0 + qw2tune: 1.0 + qwthl2tune: 1.0 + w2tune: 1.0 + length_fac: 0.5 + c_diag_3rd_mom: 7.0 + Ckh: 0.1 + Ckm: 0.1 + rrtmgp: + column_chunk_size: 123 + active_gases: ["h2o", "co2", "o3", "n2o", "co" , "ch4", "o2", "n2"] + orbital_year: 1990 + do_aerosol_rad: false + rrtmgp_coefficients_file_sw: ${SCREAM_DATA_DIR}/init/rrtmgp-data-sw-g112-210809.nc + rrtmgp_coefficients_file_lw: ${SCREAM_DATA_DIR}/init/rrtmgp-data-lw-g128-210809.nc + rrtmgp_cloud_optics_file_sw: ${SCREAM_DATA_DIR}/init/rrtmgp-cloud-optics-coeffs-sw.nc + rrtmgp_cloud_optics_file_lw: ${SCREAM_DATA_DIR}/init/rrtmgp-cloud-optics-coeffs-lw.nc + +grids_manager: + Type: Mesh Free + geo_data_source: IC_FILE + grids_names: [Physics GLL] + Physics GLL: + aliases: [Physics] + type: point_grid + number_of_global_columns: 218 + number_of_vertical_levels: 72 + +initial_conditions: + # The name of the file containing the initial conditions for this test. + Filename: ${SCREAM_DATA_DIR}/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev} + topography_filename: ${TOPO_DATA_DIR}/${EAMxx_tests_TOPO_FILE} + + #variables required for shoc + surf_sens_flux: 0.0 + surf_evap: 0.0 + + #variables required for p3 + precip_ice_surf_mass: 0.0 + precip_liq_surf_mass: 0.0 + + #variables required for rrtmgp + aero_g_sw: 0.0 + aero_ssa_sw: 0.0 + aero_tau_sw: 0.0 + aero_tau_lw: 0.0 + + #variables required for mam4_aci + dgnum: 1e-3 + +# The parameters for I/O control +Scorpio: + output_yaml_files: ["output.yaml"] +... diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/output.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/output.yaml new file mode 100644 index 00000000000..4c9f5ff5d8d --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_cldfrac_mam4_aci_p3_rrtmgp/output.yaml @@ -0,0 +1,110 @@ +%YAML 1.1 +--- +filename_prefix: shoc_cldfrac_mam4_aci_p3_rrtmgp_output +Averaging Type: Instant +Field Names: + # SHOC + - cldfrac_liq + - eddy_diff_mom + - horiz_winds + - sgs_buoy_flux + - tke + - inv_qc_relvar + - pbl_height + # CLD + - cldfrac_ice + - cldfrac_tot + # P3 + - bm + - nc + - ni + - nr + - qi + - qm + - qr + - T_prev_micro_step + - qv_prev_micro_step + - eff_radius_qc + - eff_radius_qi + - eff_radius_qr + - micro_liq_ice_exchange + - micro_vap_ice_exchange + - micro_vap_liq_exchange + - precip_ice_surf_mass + - precip_liq_surf_mass + - rainfrac + # SHOC + P3 + - qc + - qv + # SHOC + P3 + RRTMGP + - T_mid + # RRTMGP + - sfc_alb_dif_nir + - sfc_alb_dif_vis + - sfc_alb_dir_nir + - sfc_alb_dir_vis + - LW_flux_dn + - LW_flux_up + - SW_flux_dn + - SW_flux_dn_dir + - SW_flux_up + - rad_heating_pdel + - sfc_flux_lw_dn + - sfc_flux_sw_net + - bc_a1 + - bc_a3 + - bc_a4 + - dst_a1 + - dst_a3 + - so4_a1 + - so4_a2 + - so4_a3 + - pom_a1 + - pom_a3 + - pom_a4 + - soa_a1 + - soa_a2 + - soa_a3 + - nacl_a1 + - nacl_a2 + - nacl_a3 + - mom_a1 + - mom_a2 + - mom_a3 + - mom_a4 + - num_a1 + - num_a2 + - num_a3 + - num_a4 + - bc_c1 + - bc_c3 + - bc_c4 + - dst_c1 + - dst_c3 + - so4_c1 + - so4_c2 + - so4_c3 + - pom_c1 + - pom_c3 + - pom_c4 + - soa_c1 + - soa_c2 + - soa_c3 + - nacl_c1 + - nacl_c2 + - nacl_c3 + - mom_c1 + - mom_c2 + - mom_c3 + - mom_c4 + - num_c1 + - num_c2 + - num_c3 + - num_c4 + - eddy_diff_heat + - nc_nuceat_tend + - ni_activated +output_control: + Frequency: ${NUM_STEPS} + frequency_units: nsteps +... diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/CMakeLists.txt b/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/CMakeLists.txt new file mode 100644 index 00000000000..69a154223fd --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/CMakeLists.txt @@ -0,0 +1,41 @@ +INCLUDE (ScreamUtils) + +set (TEST_BASE_NAME shoc_mam4_aci) +set (FIXTURES_BASE_NAME ${TEST_BASE_NAME}_generate_output_nc_files) + +# Create the test +CreateADUnitTest(${TEST_BASE_NAME} + LIBS shoc mam + LABELS shoc physics mam4_aci + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + FIXTURES_SETUP_INDIVIDUAL ${FIXTURES_BASE_NAME} +) + +# Set AD configurable options +set (ATM_TIME_STEP 1800) +SetVarDependingOnTestSize(NUM_STEPS 2 5 48) # 1h 4h 24h +set (RUN_T0 2021-10-12-45000) + +# Determine num subcycles needed to keep shoc dt<=300s +set (SHOC_MAX_DT 300) +math (EXPR MAC_MIC_SUBCYCLES "(${ATM_TIME_STEP} + ${SHOC_MAX_DT} - 1) / ${SHOC_MAX_DT}") + +# Ensure test input files are present in the data dir +GetInputFile(scream/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev}) +GetInputFile(cam/topo/${EAMxx_tests_TOPO_FILE}) + +## Copy (and configure) yaml files needed by tests +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/input.yaml + ${CMAKE_CURRENT_BINARY_DIR}/input.yaml) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/output.yaml + ${CMAKE_CURRENT_BINARY_DIR}/output.yaml) + +# Compare output files produced by npX tests, to ensure they are bfb +include (CompareNCFiles) +CompareNCFilesFamilyMpi ( + TEST_BASE_NAME ${TEST_BASE_NAME} + FILE_META_NAME ${TEST_BASE_NAME}_output.INSTANT.nsteps_x${NUM_STEPS}.npMPIRANKS.${RUN_T0}.nc + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + LABELS shoc physics mam4_aci + META_FIXTURES_REQUIRED ${FIXTURES_BASE_NAME}_npMPIRANKS_omp1 +) diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/input.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/input.yaml new file mode 100644 index 00000000000..ae997bb3da5 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/input.yaml @@ -0,0 +1,61 @@ +%YAML 1.1 +--- +driver_options: + atmosphere_dag_verbosity_level: 5 + +time_stepping: + time_step: ${ATM_TIME_STEP} + run_t0: ${RUN_T0} # YYYY-MM-DD-XXXXX + number_of_steps: ${NUM_STEPS} + +atmosphere_processes: + atm_procs_list: [mac_mic] + schedule_type: Sequential + mac_mic: + atm_procs_list: [shoc,mam4_aci] + Type: Group + schedule_type: Sequential + number_of_subcycles: ${MAC_MIC_SUBCYCLES} + mam4_aci: + wsubmin: 0.001 + top_level_mam4xx: 6 + shoc: + lambda_low: 0.001 + lambda_high: 0.04 + lambda_slope: 2.65 + lambda_thresh: 0.02 + thl2tune: 1.0 + qw2tune: 1.0 + qwthl2tune: 1.0 + w2tune: 1.0 + length_fac: 0.5 + c_diag_3rd_mom: 7.0 + Ckh: 0.1 + Ckm: 0.1 + +grids_manager: + Type: Mesh Free + geo_data_source: IC_FILE + grids_names: [Physics GLL] + Physics GLL: + aliases: [Physics] + type: point_grid + number_of_global_columns: 218 + number_of_vertical_levels: 72 + +initial_conditions: + # The name of the file containing the initial conditions for this test. + Filename: ${SCREAM_DATA_DIR}/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev} + topography_filename: ${TOPO_DATA_DIR}/${EAMxx_tests_TOPO_FILE} + + #variables needed for SHOC + surf_sens_flux: 0.0 + surf_evap: 0.0 + + #variables needed for mam4-aci + dgnum: 1e-3 + +# The parameters for I/O control +Scorpio: + output_yaml_files: ["output.yaml"] +... diff --git a/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/output.yaml b/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/output.yaml new file mode 100644 index 00000000000..34626c99643 --- /dev/null +++ b/components/eamxx/tests/multi-process/physics_only/mam/shoc_mam4_aci/output.yaml @@ -0,0 +1,70 @@ +%YAML 1.1 +--- +filename_prefix: shoc_mam4_aci_output +Averaging Type: Instant +Field Names: + # SHOC + - cldfrac_liq + - eddy_diff_mom + - horiz_winds + - sgs_buoy_flux + - tke + - inv_qc_relvar + - pbl_height + - bc_a1 + - bc_a3 + - bc_a4 + - dst_a1 + - dst_a3 + - so4_a1 + - so4_a2 + - so4_a3 + - pom_a1 + - pom_a3 + - pom_a4 + - soa_a1 + - soa_a2 + - soa_a3 + - nacl_a1 + - nacl_a2 + - nacl_a3 + - mom_a1 + - mom_a2 + - mom_a3 + - mom_a4 + - num_a1 + - num_a2 + - num_a3 + - num_a4 + - bc_c1 + - bc_c3 + - bc_c4 + - dst_c1 + - dst_c3 + - so4_c1 + - so4_c2 + - so4_c3 + - pom_c1 + - pom_c3 + - pom_c4 + - soa_c1 + - soa_c2 + - soa_c3 + - nacl_c1 + - nacl_c2 + - nacl_c3 + - mom_c1 + - mom_c2 + - mom_c3 + - mom_c4 + - num_c1 + - num_c2 + - num_c3 + - num_c4 + - eddy_diff_heat + - nc_nuceat_tend + - ni_activated +output_control: + Frequency: ${NUM_STEPS} + frequency_units: nsteps +... diff --git a/components/eamxx/tests/single-process/CMakeLists.txt b/components/eamxx/tests/single-process/CMakeLists.txt index d417a9920f4..c6b4e0748af 100644 --- a/components/eamxx/tests/single-process/CMakeLists.txt +++ b/components/eamxx/tests/single-process/CMakeLists.txt @@ -20,6 +20,7 @@ if (SCREAM_ENABLE_MAM) # corresponding test here needs to be reworked with valid aerosol # initial conditions. add_subdirectory(mam/optics) + add_subdirectory(mam/aci) endif() if (SCREAM_TEST_LEVEL GREATER_EQUAL SCREAM_TEST_LEVEL_EXPERIMENTAL) add_subdirectory(zm) diff --git a/components/eamxx/tests/single-process/mam/aci/CMakeLists.txt b/components/eamxx/tests/single-process/mam/aci/CMakeLists.txt new file mode 100644 index 00000000000..0d5d60ee7d5 --- /dev/null +++ b/components/eamxx/tests/single-process/mam/aci/CMakeLists.txt @@ -0,0 +1,44 @@ +include (ScreamUtils) + +set (TEST_BASE_NAME mam4_aci_standalone) +set (FIXTURES_BASE_NAME ${TEST_BASE_NAME}_generate_output_nc_files) + +# Create the test +CreateADUnitTest(${TEST_BASE_NAME} + LABELS mam4_aci physics + LIBS mam + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + FIXTURES_SETUP_INDIVIDUAL ${FIXTURES_BASE_NAME} +) + +# Set AD configurable options +set (ATM_TIME_STEP 1800) +SetVarDependingOnTestSize(NUM_STEPS 2 5 48) # 1h 2.5h 24h +set (RUN_T0 2021-10-12-45000) + +## Copy (and configure) yaml files needed by tests +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/input.yaml + ${CMAKE_CURRENT_BINARY_DIR}/input.yaml) +configure_file(${CMAKE_CURRENT_SOURCE_DIR}/output.yaml + ${CMAKE_CURRENT_BINARY_DIR}/output.yaml) + +# Ensure test input files are present in the data dir +GetInputFile(scream/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev}) + +# Compare output files produced by npX tests, to ensure they are bfb +include (CompareNCFiles) + +CompareNCFilesFamilyMpi ( + TEST_BASE_NAME ${TEST_BASE_NAME} + FILE_META_NAME ${TEST_BASE_NAME}_output.INSTANT.nsteps_x1.npMPIRANKS.${RUN_T0}.nc + MPI_RANKS ${TEST_RANK_START} ${TEST_RANK_END} + LABELS mam4_aci physics + META_FIXTURES_REQUIRED ${FIXTURES_BASE_NAME}_npMPIRANKS_omp1 +) + +if (SCREAM_ENABLE_BASELINE_TESTS) + # Compare one of the output files with the baselines. + # Note: one is enough, since we already check that np1 is BFB with npX + set (OUT_FILE ${TEST_BASE_NAME}_output.INSTANT.nsteps_x1.np${TEST_RANK_END}.${RUN_T0}.nc) + CreateBaselineTest(${TEST_BASE_NAME} ${TEST_RANK_END} ${OUT_FILE} ${FIXTURES_BASE_NAME}) +endif() diff --git a/components/eamxx/tests/single-process/mam/aci/input.yaml b/components/eamxx/tests/single-process/mam/aci/input.yaml new file mode 100644 index 00000000000..fe26d828d41 --- /dev/null +++ b/components/eamxx/tests/single-process/mam/aci/input.yaml @@ -0,0 +1,44 @@ +%YAML 1.1 +--- +driver_options: + atmosphere_dag_verbosity_level: 5 + +time_stepping: + time_step: ${ATM_TIME_STEP} + run_t0: ${RUN_T0} # YYYY-MM-DD-XXXXX + number_of_steps: ${NUM_STEPS} + +atmosphere_processes: + atm_procs_list: [mam4_aci] + mam4_aci: + wsubmin: 0.001 + top_level_mam4xx: 6 + +grids_manager: + Type: Mesh Free + geo_data_source: IC_FILE + grids_names: [Physics GLL] + Physics GLL: + type: point_grid + aliases: [Physics] + number_of_global_columns: 218 + number_of_vertical_levels: 72 + +initial_conditions: + # The name of the file containing the initial conditions for this test. + Filename: ${SCREAM_DATA_DIR}/init/${EAMxx_tests_IC_FILE_MAM4xx_72lev} + topography_filename: ${TOPO_DATA_DIR}/${EAMxx_tests_TOPO_FILE} + phis : 1.0 + #These should come from the input file + + #we should get the following variables from other processes + pbl_height : 1.0 + dgnum: 4.362354500358337E-008 + cldfrac_liq_prev: 0.138584624960092 + eddy_diff_heat: 0.620820118789603 + w_variance: 1.110977269289875E-002 + +# The parameters for I/O control +Scorpio: + output_yaml_files: ["output.yaml"] +... diff --git a/components/eamxx/tests/single-process/mam/aci/output.yaml b/components/eamxx/tests/single-process/mam/aci/output.yaml new file mode 100644 index 00000000000..18f64b077cf --- /dev/null +++ b/components/eamxx/tests/single-process/mam/aci/output.yaml @@ -0,0 +1,65 @@ +%YAML 1.1 +--- +filename_prefix: mam4_aci_standalone_output +Averaging Type: Instant +Fields: + Physics: + Field Names: + - bc_a1 + - bc_a3 + - bc_a4 + - dst_a1 + - dst_a3 + - so4_a1 + - so4_a2 + - so4_a3 + - pom_a1 + - pom_a3 + - pom_a4 + - soa_a1 + - soa_a2 + - soa_a3 + - nacl_a1 + - nacl_a2 + - nacl_a3 + - mom_a1 + - mom_a2 + - mom_a3 + - mom_a4 + - num_a1 + - num_a2 + - num_a3 + - num_a4 + - bc_c1 + - bc_c3 + - bc_c4 + - dst_c1 + - dst_c3 + - so4_c1 + - so4_c2 + - so4_c3 + - pom_c1 + - pom_c3 + - pom_c4 + - soa_c1 + - soa_c2 + - soa_c3 + - nacl_c1 + - nacl_c2 + - nacl_c3 + - mom_c1 + - mom_c2 + - mom_c3 + - mom_c4 + - num_c1 + - num_c2 + - num_c3 + - num_c4 + - eddy_diff_heat + - nc_nuceat_tend + - ni_activated + +output_control: + Frequency: 1 + frequency_units: nsteps +... diff --git a/externals/mam4xx b/externals/mam4xx index 7b0be56f4c3..a8aa56ba915 160000 --- a/externals/mam4xx +++ b/externals/mam4xx @@ -1 +1 @@ -Subproject commit 7b0be56f4c3f164259c6a9d8a0d16d6f5b72f3ae +Subproject commit a8aa56ba915c83c22edc54ce5fd8ae4a23e6173e