ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
Public Member Functions | Static Public Attributes | Private Types | Private Member Functions | Private Attributes | List of all members
WSM6 Class Reference

#include <ERF_WSM6.H>

Inheritance diagram for WSM6:
Inheritance graph
[legend]
Collaboration diagram for WSM6:
Collaboration graph
[legend]

Public Member Functions

 WSM6 ()
 
virtual ~WSM6 ()=default
 
void Define (SolverChoice &sc) override
 
void Init (const amrex::MultiFab &cons_in, const amrex::BoxArray &grids, const amrex::Geometry &geom, const amrex::Real &dt_advance, std::unique_ptr< amrex::MultiFab > &z_phys_nd, std::unique_ptr< amrex::MultiFab > &detJ_cc) override
 
void Set_dzmin (const amrex::Real dz_min) override
 
void Copy_State_to_Micro (const amrex::MultiFab &cons_in) override
 
void Copy_Micro_to_State (amrex::MultiFab &cons_in) override
 
void Update_Micro_Vars (amrex::MultiFab &cons_in) override
 
void Update_State_Vars (amrex::MultiFab &cons_in, const amrex::MultiFab &) override
 
void Advance (const amrex::Real &dt_advance, const SolverChoice &solverChoice) override
 
amrex::MultiFab * Qmoist_Ptr (const int &varIdx) override
 
int Qmoist_Size () override
 
int Qstate_Moist_Size () override
 
int Qstate_Moist_NumConc_Size () override
 
void Qmoist_Restart_Vars (const SolverChoice &, std::vector< int > &a_idx, std::vector< std::string > &a_names) const override
 
- Public Member Functions inherited from NullMoist
 NullMoist ()
 
virtual ~NullMoist ()=default
 
virtual int Qstate_NonMoist_Size ()
 
virtual void GetPlotVarNames (amrex::Vector< std::string > &a_vec) const
 
virtual void GetPlotVar (const std::string &, amrex::MultiFab &) const
 
virtual void GetPlotVar (const std::string &a_name, amrex::MultiFab &a_mf, const int) const
 
virtual void SetCurrentLevel (const int &)
 
virtual void InitLevel (const int, const amrex::MultiFab &)
 
virtual int getDiagnosticsInterval () const
 
virtual void Set_RealWidth (const int)
 

Static Public Attributes

static constexpr amrex::Real dtcldcr = amrex::Real(120.0)
 
static constexpr amrex::Real n0r = amrex::Real(8.0e6)
 
static constexpr amrex::Real avtr = amrex::Real(841.9)
 
static constexpr amrex::Real bvtr = amrex::Real(0.8)
 
static constexpr amrex::Real r0 = amrex::Real(0.8e-5)
 
static constexpr amrex::Real peaut = amrex::Real(0.55)
 
static constexpr amrex::Real xncr = amrex::Real(3.0e8)
 
static constexpr amrex::Real xmyu = amrex::Real(1.718e-5)
 
static constexpr amrex::Real avts = amrex::Real(11.72)
 
static constexpr amrex::Real bvts = amrex::Real(0.41)
 
static constexpr amrex::Real lamdarmax = amrex::Real(8.0e4)
 
static constexpr amrex::Real lamdasmax = amrex::Real(1.0e5)
 
static constexpr amrex::Real dicon = amrex::Real(11.9)
 
static constexpr amrex::Real dimax = amrex::Real(500.0e-6)
 
static constexpr amrex::Real pfrz1 = amrex::Real(100.0)
 
static constexpr amrex::Real pfrz2 = amrex::Real(0.66)
 
static constexpr amrex::Real qcrmin = amrex::Real(1.0e-9)
 
static constexpr amrex::Real eacrc = amrex::Real(1.0)
 
static constexpr amrex::Real dens_snow = amrex::Real(100.0)
 
static constexpr amrex::Real qs0 = amrex::Real(6.0e-4)
 
static constexpr amrex::Real n0smax = amrex::Real(1.0e11)
 
static constexpr amrex::Real n0s = amrex::Real(2.0e6)
 
static constexpr amrex::Real alpha_wsm6 = amrex::Real(0.12)
 

Private Types

using FabPtr = std::shared_ptr< amrex::MultiFab >
 

Private Member Functions

void initialize_coeffs ()
 

Private Attributes

int m_qmoist_size = 3
 
int n_qstate_moist_size = 6
 
int n_qstate_moist_numconc_size = 0
 
amrex::Vector< int > MicVarMap
 
amrex::Geometry m_geom
 
amrex::Real dt {0.0}
 
amrex::Real m_dzmin {0.0}
 
int nlev {0}
 
int zlo {0}
 
int zhi {0}
 
int m_axis {2}
 
bool m_do_cond {true}
 
amrex::MultiFab * m_z_phys_nd {nullptr}
 
amrex::MultiFab * m_detJ_cc {nullptr}
 
amrex::Array< FabPtr, MicVar_WSM6::NumVars > mic_fab_vars
 
bool m_hail_opt {false}
 
amrex::Real m_n0g {0}
 
amrex::Real m_deng {0}
 
amrex::Real m_avtg {0}
 
amrex::Real m_bvtg {0}
 
amrex::Real m_lamdagmax {0}
 
amrex::Real m_pi_wsm6 {0}
 
amrex::Real m_xlv1 {0}
 
amrex::Real m_qc0 {0}
 
amrex::Real m_qck1 {0}
 
amrex::Real m_pidnc {0}
 
amrex::Real m_bvtr1 {0}
 
amrex::Real m_bvtr2 {0}
 
amrex::Real m_bvtr3 {0}
 
amrex::Real m_bvtr4 {0}
 
amrex::Real m_bvtr6 {0}
 
amrex::Real m_g1pbr {0}
 
amrex::Real m_g3pbr {0}
 
amrex::Real m_g4pbr {0}
 
amrex::Real m_g5pbro2 {0}
 
amrex::Real m_g6pbr {0}
 
amrex::Real m_pvtr {0}
 
amrex::Real m_eacrr {0}
 
amrex::Real m_pacrr {0}
 
amrex::Real m_precr1 {0}
 
amrex::Real m_precr2 {0}
 
amrex::Real m_roqimax {0}
 
amrex::Real m_bvts1 {0}
 
amrex::Real m_bvts2 {0}
 
amrex::Real m_bvts3 {0}
 
amrex::Real m_bvts4 {0}
 
amrex::Real m_g1pbs {0}
 
amrex::Real m_g3pbs {0}
 
amrex::Real m_g4pbs {0}
 
amrex::Real m_g5pbso2 {0}
 
amrex::Real m_pvts {0}
 
amrex::Real m_pacrs {0}
 
amrex::Real m_precs1 {0}
 
amrex::Real m_precs2 {0}
 
amrex::Real m_pidn0r {0}
 
amrex::Real m_pidn0s {0}
 
amrex::Real m_pacrc {0}
 
amrex::Real m_bvtg1 {0}
 
amrex::Real m_bvtg2 {0}
 
amrex::Real m_bvtg3 {0}
 
amrex::Real m_bvtg4 {0}
 
amrex::Real m_g1pbg {0}
 
amrex::Real m_g3pbg {0}
 
amrex::Real m_g4pbg {0}
 
amrex::Real m_g5pbgo2 {0}
 
amrex::Real m_pvtg {0}
 
amrex::Real m_pacrg {0}
 
amrex::Real m_precg1 {0}
 
amrex::Real m_precg2 {0}
 
amrex::Real m_pidn0g {0}
 
amrex::Real m_rslopermax {0}
 
amrex::Real m_rslopesmax {0}
 
amrex::Real m_rslopegmax {0}
 
amrex::Real m_rsloperbmax {0}
 
amrex::Real m_rslopesbmax {0}
 
amrex::Real m_rslopegbmax {0}
 
amrex::Real m_rsloper2max {0}
 
amrex::Real m_rslopes2max {0}
 
amrex::Real m_rslopeg2max {0}
 
amrex::Real m_rsloper3max {0}
 
amrex::Real m_rslopes3max {0}
 
amrex::Real m_rslopeg3max {0}
 

Member Typedef Documentation

◆ FabPtr

using WSM6::FabPtr = std::shared_ptr<amrex::MultiFab>
private

Constructor & Destructor Documentation

◆ WSM6()

WSM6::WSM6 ( )
inline
41 {}

◆ ~WSM6()

virtual WSM6::~WSM6 ( )
virtualdefault

Member Function Documentation

◆ Advance()

void WSM6::Advance ( const amrex::Real &  dt_advance,
const SolverChoice &  solverChoice 
)
overridevirtual

Reimplemented from NullMoist.

844 {
845  dt = dt_advance;
846 
847  int microphysics_debug = 0;
848  std::string micro_diag_mode = "canonical";
849  std::vector<std::string> micro_diag_tags = {"standing"};
850  std::vector<std::string> micro_diag_expr = {"standing"};
851  std::vector<std::string> micro_diag_store = {"standing"};
852  std::vector<int> micro_diag_target_column;
853  {
854  amrex::ParmParse pp("erf");
855  pp.query("microphysics_debug", microphysics_debug);
856  pp.query("micro_diag_mode", micro_diag_mode);
857  pp.queryarr("micro_diag_tags", micro_diag_tags);
858  pp.queryarr("micro_diag_expr", micro_diag_expr);
859  pp.queryarr("micro_diag_store", micro_diag_store);
860  pp.queryarr("micro_diag_target_column", micro_diag_target_column);
861  }
862  microphysics_debug = std::max(0, std::min(2, microphysics_debug));
863 #ifdef ERF_USE_WSM6_FORT
864  const std::string micro_diag_mode_lower = [&]{ std::string m = micro_diag_mode; std::transform(m.begin(), m.end(), m.begin(), [](unsigned char c){ return static_cast<char>(std::tolower(c)); }); return m; }();
865  const bool micro_diag_forensic = (micro_diag_mode_lower == "forensic" || micro_diag_mode_lower == "both");
866  const int microphysics_debug_bridge = micro_diag_forensic
867  ? microphysics_debug
868  : std::min(microphysics_debug, 1);
869  bool use_wsm6_cpp_answer = false;
870  { amrex::ParmParse pp("erf");
871  pp.query("use_wsm6_cpp_answer", use_wsm6_cpp_answer); }
872  bool run_wsm6_fort = !use_wsm6_cpp_answer;
873 
874  static bool wsm6_inited = false;
875 
876  // Minimal phase-1 initialization for single-moment WSM6.
877  if (!wsm6_inited) {
878  constexpr double den0 = 1.28; // Standard dry-air density (kg/m^3)
879  constexpr double denr = static_cast<double>(rhoh2o);
880  constexpr double dens = static_cast<double>(rhos);
881  constexpr double cl = static_cast<double>(Cp_l);
882  constexpr double cpv = static_cast<double>(Cp_v);
883  constexpr int hail_opt = 0; // Graupel mode
884  mp_wsm6_init_c(den0, denr, dens, cl, cpv, hail_opt);
885  wsm6_inited = true;
886  }
887 #endif
888 
889  constexpr double g = static_cast<double>(CONST_GRAV);
890  constexpr double cpd = static_cast<double>(Cp_d);
891  constexpr double cpv = static_cast<double>(Cp_v);
892  constexpr double rd = static_cast<double>(R_d);
893  constexpr double rv = static_cast<double>(R_v);
894  constexpr double t0c = 273.15;
895  constexpr double ep1 = static_cast<double>(R_v / R_d - one);
896  amrex::ignore_unused(g, rd, ep1);
897  constexpr double ep2 = static_cast<double>(R_d / R_v);
898  constexpr double qmin = 1.0e-12;
899  constexpr double xls = static_cast<double>(lsub);
900  constexpr double xlv0 = static_cast<double>(lat_vap);
901  constexpr double xlf0 = static_cast<double>(lat_ice);
902  constexpr double den0 = 1.28;
903  constexpr double denr = static_cast<double>(rhoh2o);
904  constexpr double cliq = static_cast<double>(Cp_l);
905  constexpr double cice = 2106.0;
906  constexpr double psat = 610.78;
907  for (MFIter mfi(*mic_fab_vars[MicVar_WSM6::qv], TileNoZ()); mfi.isValid(); ++mfi) {
908  const Box box = mfi.tilebox();
909  const Box fab_box = mfi.fabbox();
910 
911  auto const& t_arr = mic_fab_vars[MicVar_WSM6::tabs]->array(mfi);
912  auto const& qv_arr = mic_fab_vars[MicVar_WSM6::qv]->array(mfi);
913  auto const& qc_arr = mic_fab_vars[MicVar_WSM6::qc]->array(mfi);
914  auto const& qi_arr = mic_fab_vars[MicVar_WSM6::qi]->array(mfi);
915  auto const& qr_arr = mic_fab_vars[MicVar_WSM6::qr]->array(mfi);
916  auto const& qs_arr = mic_fab_vars[MicVar_WSM6::qs]->array(mfi);
917  auto const& qg_arr = mic_fab_vars[MicVar_WSM6::qg]->array(mfi);
918  auto const& den_arr = mic_fab_vars[MicVar_WSM6::rho]->array(mfi);
919  auto const& p_arr = mic_fab_vars[MicVar_WSM6::pres]->array(mfi);
920  auto rain_arr = mic_fab_vars[MicVar_WSM6::rain_accum]->array(mfi);
921  auto snow_arr = mic_fab_vars[MicVar_WSM6::snow_accum]->array(mfi);
922  auto graup_arr = mic_fab_vars[MicVar_WSM6::graup_accum]->array(mfi);
923 
924  const int ilo = box.smallEnd(0);
925  const int ihi = box.bigEnd(0);
926  const int jlo = box.smallEnd(1);
927  const int jhi = box.bigEnd(1);
928  const int klo = box.smallEnd(2);
929  const int khi = box.bigEnd(2);
930  const bool has_target_override = (micro_diag_target_column.size() == 2);
931  const int diag_i = has_target_override ? micro_diag_target_column[0] : ilo;
932  const int diag_j = has_target_override ? micro_diag_target_column[1] : jlo;
933 
934  const int imlo = fab_box.smallEnd(0);
935  const int imhi = fab_box.bigEnd(0);
936  const int jmlo = fab_box.smallEnd(1);
937  const int jmhi = fab_box.bigEnd(1);
938  const int kmlo = fab_box.smallEnd(2);
939  const int kmhi = fab_box.bigEnd(2);
940  amrex::ignore_unused(ihi, jhi, diag_i, diag_j, imlo, imhi, jmlo, jmhi, kmlo, kmhi);
941 
942  const Real dz_val = m_geom.CellSize(2);
943  FArrayBox delz_fab(fab_box, 1);
944  auto const& delz_arr = delz_fab.array();
945  ParallelFor(fab_box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
946  delz_arr(i,j,k) = dz_val;
947  });
948 
949  const Array4<const Real> z_arr = (m_z_phys_nd) ? m_z_phys_nd->const_array(mfi) : Array4<const Real> {};
950  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
951  delz_arr(i,j,k) = (z_arr) ? Real(0.25) * ( (z_arr(i ,j ,k+1) - z_arr(i ,j ,k))
952  + (z_arr(i+1,j ,k+1) - z_arr(i+1,j ,k))
953  + (z_arr(i ,j+1,k+1) - z_arr(i ,j+1,k))
954  + (z_arr(i+1,j+1,k+1) - z_arr(i+1,j+1,k)) ) : dz_val;
955  });
956 
957  Box box2d(box);
958  box2d.makeSlab(2, 0);
959  Box fab_box2d(fab_box);
960  fab_box2d.makeSlab(2, 0);
961 
962  // Fortran bridge uses ims:ime, jms:jme storage bounds; these buffers must
963  // therefore be allocated on fab_box extents even if C++ kernels only
964  // update the valid tile slab (box2d).
965  FArrayBox rainncv_fab(fab_box2d, 1);
966  FArrayBox sr_fab(fab_box2d, 1);
967  FArrayBox snowncv_fab(fab_box2d, 1);
968  FArrayBox graupelncv_fab(fab_box2d, 1);
969  FArrayBox rainacc_fab(fab_box2d, 1);
970  FArrayBox snowacc_fab(fab_box2d, 1);
971  FArrayBox graupacc_fab(fab_box2d, 1);
972 
973  auto const& rainncv_arr = rainncv_fab.array();
974  auto const& sr_arr = sr_fab.array();
975  auto const& snowncv_arr = snowncv_fab.array();
976  auto const& graupelncv_arr = graupelncv_fab.array();
977  auto const& rainacc_arr = rainacc_fab.array();
978  auto const& snowacc_arr = snowacc_fab.array();
979  auto const& graupacc_arr = graupacc_fab.array();
980  ParallelFor(fab_box2d, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
981  rainncv_arr(i,j,k) = Real(0.0);
982  sr_arr(i,j,k) = Real(0.0);
983  snowncv_arr(i,j,k) = Real(0.0);
984  graupelncv_arr(i,j,k) = Real(0.0);
985  rainacc_arr(i,j,k) = Real(0.0);
986  snowacc_arr(i,j,k) = Real(0.0);
987  graupacc_arr(i,j,k) = Real(0.0);
988  });
989  ParallelFor(box2d, [=] AMREX_GPU_DEVICE (int i, int j, int) {
990  rainacc_arr(i,j,0) = rain_arr(i,j,klo);
991  snowacc_arr(i,j,0) = snow_arr(i,j,klo);
992  graupacc_arr(i,j,0) = graup_arr(i,j,klo);
993  rainncv_arr(i,j,0) = Real(0.0);
994  sr_arr(i,j,0) = Real(0.0);
995  snowncv_arr(i,j,0) = Real(0.0);
996  graupelncv_arr(i,j,0) = Real(0.0);
997  });
998 
999 #ifdef ERF_USE_WSM6_FORT
1000  if (run_wsm6_fort) {
1001  mp_wsm6_run_c(
1002  t_arr.dataPtr(),
1003  qv_arr.dataPtr(), qc_arr.dataPtr(), qi_arr.dataPtr(),
1004  qr_arr.dataPtr(), qs_arr.dataPtr(), qg_arr.dataPtr(),
1005  den_arr.dataPtr(), p_arr.dataPtr(), delz_arr.dataPtr(),
1006  static_cast<double>(dt), g, cpd, cpv, rd, rv, t0c, ep1, ep2, qmin,
1007  xls, xlv0, xlf0, den0, denr, cliq, cice, psat,
1008  rainacc_arr.dataPtr(), rainncv_arr.dataPtr(), sr_arr.dataPtr(),
1009  snowacc_arr.dataPtr(), snowncv_arr.dataPtr(),
1010  graupacc_arr.dataPtr(), graupelncv_arr.dataPtr(),
1011  imlo, imhi, jmlo, jmhi, kmlo, kmhi,
1012  ilo, ihi, jlo, jhi, klo, khi, microphysics_debug_bridge, diag_i, diag_j);
1013  } else {
1014 #endif
1015  // --- Phase 4 native C++ kernel ---
1016 
1017  // box2d for 1D per-column arrays (already defined above)
1018  // delqrs1/2/3, delqi: surface precipitation flux accumulators
1019  FArrayBox delqrs1_fab(box2d,1);
1020  FArrayBox delqrs2_fab(box2d,1);
1021  FArrayBox delqrs3_fab(box2d,1);
1022  FArrayBox delqi_fab(box2d,1);
1023  FArrayBox tstepsnow_fab(box2d,1);
1024  FArrayBox tstepgraup_fab(box2d,1);
1025  auto const& delqrs1_arr = delqrs1_fab.array();
1026  auto const& delqrs2_arr = delqrs2_fab.array();
1027  auto const& delqrs3_arr = delqrs3_fab.array();
1028  auto const& delqi_arr = delqi_fab.array();
1029  auto const& tstepsnow_arr = tstepsnow_fab.array();
1030  auto const& tstepgraup_arr = tstepgraup_fab.array();
1031 
1032  ParallelFor(box2d, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1033  delqrs1_arr(i,j,k) = Real(0.0);
1034  delqrs2_arr(i,j,k) = Real(0.0);
1035  delqrs3_arr(i,j,k) = Real(0.0);
1036  delqi_arr(i,j,k) = Real(0.0);
1037  tstepsnow_arr(i,j,k) = Real(0.0);
1038  tstepgraup_arr(i,j,k) = Real(0.0);
1039  });
1040  // 3D working FABs
1041  FArrayBox denfac_fab(fab_box,1); FArrayBox xni_fab(fab_box,1);
1042  FArrayBox cpm_fab(fab_box,1); FArrayBox xl_fab(fab_box,1);
1043  FArrayBox qsatw_fab(fab_box,1); FArrayBox qsati_fab(fab_box,1);
1044  FArrayBox rhw_fab(fab_box,1); FArrayBox rhi_fab(fab_box,1);
1045  FArrayBox den_tmp_fab(fab_box,1); FArrayBox delz_tmp_fab(fab_box,1);
1046  FArrayBox n0sfac_fab(fab_box,1);
1047  FArrayBox qrs_tmp_r_fab(fab_box,1); FArrayBox qrs_tmp_s_fab(fab_box,1);
1048  FArrayBox qrs_tmp_g_fab(fab_box,1);
1049  FArrayBox rslope_r_fab(fab_box,1); FArrayBox rslope_s_fab(fab_box,1);
1050  FArrayBox rslope_g_fab(fab_box,1);
1051  FArrayBox rslope2_r_fab(fab_box,1); FArrayBox rslope2_s_fab(fab_box,1);
1052  FArrayBox rslope2_g_fab(fab_box,1);
1053  FArrayBox rslope3_r_fab(fab_box,1); FArrayBox rslope3_s_fab(fab_box,1);
1054  FArrayBox rslope3_g_fab(fab_box,1);
1055  FArrayBox rslopeb_r_fab(fab_box,1); FArrayBox rslopeb_s_fab(fab_box,1);
1056  FArrayBox rslopeb_g_fab(fab_box,1);
1057  FArrayBox work1_r_fab(fab_box,1); FArrayBox work1_s_fab(fab_box,1);
1058  FArrayBox work1_g_fab(fab_box,1);
1059  FArrayBox work2_fab(fab_box,1); FArrayBox workdiffw_fab(fab_box,1);
1060  FArrayBox workdiffi_fab(fab_box,1);
1061  FArrayBox workr_fab(fab_box,1); FArrayBox worka_fab(fab_box,1);
1062  FArrayBox work1c_fab(fab_box,1);
1063  FArrayBox denqrs1_fab(fab_box,1); FArrayBox denqrs2_fab(fab_box,1);
1064  FArrayBox denqrs3_fab(fab_box,1); FArrayBox denqci_fab(fab_box,1);
1065  FArrayBox fall_r_fab(fab_box,1); FArrayBox fall_s_fab(fab_box,1);
1066  FArrayBox fall_g_fab(fab_box,1); FArrayBox fallc_fab(fab_box,1);
1067  FArrayBox qsum_fab(fab_box,1);
1068  FArrayBox nislfv_r_diag_fab(fab_box,6);
1069  FArrayBox nislfv_sg_diag_fab(fab_box,6);
1070  FArrayBox sed_cell_scratch_fab(fab_box, WSM6SedCellScratch::NumComps);
1071  Box sed_node_box = amrex::surroundingNodes(fab_box, 2);
1072  FArrayBox sed_node_scratch_fab(sed_node_box, WSM6SedNodeScratch::NumComps);
1073  // process rates
1074  FArrayBox praut_fab(fab_box,1); FArrayBox pracw_fab(fab_box,1);
1075  FArrayBox prevp_fab(fab_box,1); FArrayBox psdep_fab(fab_box,1);
1076  FArrayBox pgdep_fab(fab_box,1); FArrayBox psaut_fab(fab_box,1);
1077  FArrayBox pgaut_fab(fab_box,1); FArrayBox praci_fab(fab_box,1);
1078  FArrayBox piacr_fab(fab_box,1); FArrayBox psaci_fab(fab_box,1);
1079  FArrayBox psacw_fab(fab_box,1); FArrayBox pgacw_fab(fab_box,1);
1080  FArrayBox pgaci_fab(fab_box,1); FArrayBox paacw_fab(fab_box,1);
1081  FArrayBox pracs_fab(fab_box,1); FArrayBox psacr_fab(fab_box,1);
1082  FArrayBox pgacr_fab(fab_box,1); FArrayBox pgacs_fab(fab_box,1);
1083  FArrayBox pigen_fab(fab_box,1); FArrayBox pidep_fab(fab_box,1);
1084  FArrayBox pcond_fab(fab_box,1); FArrayBox psmlt_fab(fab_box,1);
1085  FArrayBox pgmlt_fab(fab_box,1); FArrayBox pseml_fab(fab_box,1);
1086  FArrayBox pgeml_fab(fab_box,1); FArrayBox psevp_fab(fab_box,1);
1087  FArrayBox pgevp_fab(fab_box,1);
1088  FArrayBox pimlt_fab(fab_box,1); FArrayBox pihmf_fab(fab_box,1);
1089  FArrayBox pihtf_fab(fab_box,1); FArrayBox pgfrz_fab(fab_box,1);
1090 
1091  auto const& denfac_arr = denfac_fab.array();
1092  auto const& xni_arr = xni_fab.array();
1093  auto const& cpm_arr = cpm_fab.array();
1094  auto const& xl_arr = xl_fab.array();
1095  auto const& qsatw_arr = qsatw_fab.array();
1096  auto const& qsati_arr = qsati_fab.array();
1097  auto const& rhw_arr = rhw_fab.array();
1098  auto const& rhi_arr = rhi_fab.array();
1099  auto const& den_tmp_arr = den_tmp_fab.array();
1100  auto const& delz_tmp_arr = delz_tmp_fab.array();
1101  auto const& n0sfac_arr = n0sfac_fab.array();
1102  auto const& qrs_tmp_r_arr = qrs_tmp_r_fab.array();
1103  auto const& qrs_tmp_s_arr = qrs_tmp_s_fab.array();
1104  auto const& qrs_tmp_g_arr = qrs_tmp_g_fab.array();
1105  auto const& rslope_r_arr = rslope_r_fab.array();
1106  auto const& rslope_s_arr = rslope_s_fab.array();
1107  auto const& rslope_g_arr = rslope_g_fab.array();
1108  auto const& rslope2_r_arr = rslope2_r_fab.array();
1109  auto const& rslope2_s_arr = rslope2_s_fab.array();
1110  auto const& rslope2_g_arr = rslope2_g_fab.array();
1111  auto const& rslope3_r_arr = rslope3_r_fab.array();
1112  auto const& rslope3_s_arr = rslope3_s_fab.array();
1113  auto const& rslope3_g_arr = rslope3_g_fab.array();
1114  auto const& rslopeb_r_arr = rslopeb_r_fab.array();
1115  auto const& rslopeb_s_arr = rslopeb_s_fab.array();
1116  auto const& rslopeb_g_arr = rslopeb_g_fab.array();
1117  auto const& work1_r_arr = work1_r_fab.array();
1118  auto const& work1_s_arr = work1_s_fab.array();
1119  auto const& work1_g_arr = work1_g_fab.array();
1120  auto const& work2_arr = work2_fab.array();
1121  auto const& workdiffw_arr = workdiffw_fab.array();
1122  auto const& workdiffi_arr = workdiffi_fab.array();
1123  auto const& workr_arr = workr_fab.array();
1124  auto const& worka_arr = worka_fab.array();
1125  auto const& work1c_arr = work1c_fab.array();
1126  auto const& denqrs1_arr = denqrs1_fab.array();
1127  auto const& denqrs2_arr = denqrs2_fab.array();
1128  auto const& denqrs3_arr = denqrs3_fab.array();
1129  auto const& denqci_arr = denqci_fab.array();
1130  auto const& fall_r_arr = fall_r_fab.array();
1131  auto const& fall_s_arr = fall_s_fab.array();
1132  auto const& fall_g_arr = fall_g_fab.array();
1133  auto const& fallc_arr = fallc_fab.array();
1134  auto const& qsum_arr = qsum_fab.array();
1135  auto const& nislfv_r_diag_arr = nislfv_r_diag_fab.array();
1136  auto const& nislfv_sg_diag_arr = nislfv_sg_diag_fab.array();
1137  auto const& sed_cell_scratch_arr = sed_cell_scratch_fab.array();
1138  auto const& sed_node_scratch_arr = sed_node_scratch_fab.array();
1139 
1140  ParallelFor(fab_box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1141  work1c_arr(i,j,k) = Real(0.0);
1142  });
1143  ParallelFor(fab_box, nislfv_r_diag_fab.nComp(), [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
1144  nislfv_r_diag_arr(i,j,k,n) = Real(0.0);
1145  nislfv_sg_diag_arr(i,j,k,n) = Real(0.0);
1146  });
1147  auto const& praut_arr = praut_fab.array();
1148  auto const& pracw_arr = pracw_fab.array();
1149  auto const& prevp_arr = prevp_fab.array();
1150  auto const& psdep_arr = psdep_fab.array();
1151  auto const& pgdep_arr = pgdep_fab.array();
1152  auto const& psaut_arr = psaut_fab.array();
1153  auto const& pgaut_arr = pgaut_fab.array();
1154  auto const& praci_arr = praci_fab.array();
1155  auto const& piacr_arr = piacr_fab.array();
1156  auto const& psaci_arr = psaci_fab.array();
1157  auto const& psacw_arr = psacw_fab.array();
1158  auto const& pgacw_arr = pgacw_fab.array();
1159  auto const& pgaci_arr = pgaci_fab.array();
1160  auto const& paacw_arr = paacw_fab.array();
1161  auto const& pracs_arr = pracs_fab.array();
1162  auto const& psacr_arr = psacr_fab.array();
1163  auto const& pgacr_arr = pgacr_fab.array();
1164  auto const& pgacs_arr = pgacs_fab.array();
1165  auto const& pigen_arr = pigen_fab.array();
1166  auto const& pidep_arr = pidep_fab.array();
1167  auto const& pcond_arr = pcond_fab.array();
1168  auto const& psmlt_arr = psmlt_fab.array();
1169  auto const& pgmlt_arr = pgmlt_fab.array();
1170  auto const& pseml_arr = pseml_fab.array();
1171  auto const& pgeml_arr = pgeml_fab.array();
1172  auto const& psevp_arr = psevp_fab.array();
1173  auto const& pgevp_arr = pgevp_fab.array();
1174  auto const& pimlt_arr = pimlt_fab.array();
1175  auto const& pihmf_arr = pihmf_fab.array();
1176  auto const& pihtf_arr = pihtf_fab.array();
1177  auto const& pgfrz_arr = pgfrz_fab.array();
1178 
1179  // Groups A-E: pre-loop setup
1180  // Clamp negative values (Group A)
1181  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1182  qc_arr(i,j,k) = amrex::max(qc_arr(i,j,k), Real(0.0));
1183  qr_arr(i,j,k) = amrex::max(qr_arr(i,j,k), Real(0.0));
1184  qi_arr(i,j,k) = amrex::max(qi_arr(i,j,k), Real(0.0));
1185  qs_arr(i,j,k) = amrex::max(qs_arr(i,j,k), Real(0.0));
1186  qg_arr(i,j,k) = amrex::max(qg_arr(i,j,k), Real(0.0));
1187  den_tmp_arr(i,j,k) = den_arr(i,j,k);
1188  delz_tmp_arr(i,j,k) = delz_arr(i,j,k);
1189  });
1190 
1191  // Group B: cpm, xl — computed once from initial state [lines 455-460]
1192  const Real xlv1_loc = m_xlv1;
1193  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1194  cpm_arr(i,j,k) = wsm6_cpmcal(qv_arr(i,j,k), Real(qmin), Real(cpd), Real(cpv));
1195  xl_arr(i,j,k) = wsm6_xlcal(t_arr(i,j,k), Real(xlv0), xlv1_loc, Real(t0c));
1196  });
1197 
1198  // Outer minor timestep loop (Rule 29)
1199  const int wsm6_loops = std::max(
1200  static_cast<int>(std::round(dt / dtcldcr)), 1);
1201  const Real dtcld = dt / static_cast<Real>(wsm6_loops);
1202  const Real qc0 = m_qc0;
1203  const Real qck1 = m_qck1;
1204  const Real pvtr = m_pvtr;
1205  const Real pacrr = m_pacrr;
1206  const Real precr1 = m_precr1;
1207  const Real precr2 = m_precr2;
1208  const Real roqimax= m_roqimax;
1209  const Real pvts = m_pvts;
1210  const Real pacrc = m_pacrc;
1211  const Real precs1 = m_precs1;
1212  const Real precs2 = m_precs2;
1213  const Real g6pbr = m_g6pbr;
1214  const Real pvtg = m_pvtg;
1215  const Real pacrg = m_pacrg;
1216  const Real precg1 = m_precg1;
1217  const Real precg2 = m_precg2;
1218 
1219  for (int loop = 0; loop < wsm6_loops; ++loop) {
1220  // G1b: denfac = sqrt(den0/den) [lines 503-515]
1221  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1222  const Real invden = Real(1.0) / den_arr(i,j,k);
1223  denfac_arr(i,j,k) = std::sqrt(invden * Real(den0));
1224  });
1225  // G1c: qsatw, qsati, rhw, rhi [lines 517-549]
1226  {
1227  const Real ttp = Real(t0c) + Real(0.01);
1228  const Real dldt = Real(cpv) - Real(cliq);
1229  const Real xa = -dldt / Real(rv);
1230  const Real xb = xa + Real(xlv0) / (Real(rv)*ttp);
1231  const Real dldti= Real(cpv) - Real(cice);
1232  const Real xai = -dldti / Real(rv);
1233  const Real xbi = xai + Real(xls) / (Real(rv)*ttp);
1234  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1235  const Real tr = ttp / t_arr(i,j,k);
1236  Real qsw = Real(psat)*std::exp(std::log(tr)*xa)*std::exp(xb*(Real(1.0)-tr));
1237  qsw = amrex::min(qsw, Real(0.99)*p_arr(i,j,k));
1238  qsw = Real(ep2)*qsw / (p_arr(i,j,k) - qsw);
1239  qsw = amrex::max(qsw, Real(qmin));
1240  qsatw_arr(i,j,k) = qsw;
1241  rhw_arr(i,j,k) = amrex::max(qv_arr(i,j,k)/qsw, Real(qmin));
1242  Real qsi = (t_arr(i,j,k) < ttp)
1243  ? Real(psat)*std::exp(std::log(tr)*xai)*std::exp(xbi*(Real(1.0)-tr))
1244  : Real(psat)*std::exp(std::log(tr)*xa )*std::exp(xb *(Real(1.0)-tr));
1245  qsi = amrex::min(qsi, Real(0.99)*p_arr(i,j,k));
1246  qsi = Real(ep2)*qsi / (p_arr(i,j,k) - qsi);
1247  qsi = amrex::max(qsi, Real(qmin));
1248  qsati_arr(i,j,k) = qsi;
1249  rhi_arr(i,j,k) = amrex::max(qv_arr(i,j,k)/qsi, Real(qmin));
1250  });
1251  }
1252  // G2: zero all process rates each sub-step [lines 555-594]
1253  // WSM6-CPP TAG: RATES_ZERO
1254  // legacy_group: G2
1255  // process: Zero all process rates each sub-step
1256  // compare_vars: prevp, psdep, pgdep, praut, psaut, pgaut, pracw, praci, psaci, pracs, pidep, pcond, psmlt, pgmlt, pseml, psevp, fall_r, fall_s, fall_g, fallc
1257  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1258  prevp_arr(i,j,k) = Real(0.0); psdep_arr(i,j,k) = Real(0.0);
1259  pgdep_arr(i,j,k) = Real(0.0); praut_arr(i,j,k) = Real(0.0);
1260  psaut_arr(i,j,k) = Real(0.0); pgaut_arr(i,j,k) = Real(0.0);
1261  pracw_arr(i,j,k) = Real(0.0); praci_arr(i,j,k) = Real(0.0);
1262  piacr_arr(i,j,k) = Real(0.0); psaci_arr(i,j,k) = Real(0.0);
1263  psacw_arr(i,j,k) = Real(0.0); pracs_arr(i,j,k) = Real(0.0);
1264  psacr_arr(i,j,k) = Real(0.0); pgacw_arr(i,j,k) = Real(0.0);
1265  paacw_arr(i,j,k) = Real(0.0); pgaci_arr(i,j,k) = Real(0.0);
1266  pgacr_arr(i,j,k) = Real(0.0); pgacs_arr(i,j,k) = Real(0.0);
1267  pigen_arr(i,j,k) = Real(0.0); pidep_arr(i,j,k) = Real(0.0);
1268  pcond_arr(i,j,k) = Real(0.0); psmlt_arr(i,j,k) = Real(0.0);
1269  pgmlt_arr(i,j,k) = Real(0.0); pseml_arr(i,j,k) = Real(0.0);
1270  pgeml_arr(i,j,k) = Real(0.0); psevp_arr(i,j,k) = Real(0.0);
1271  pgevp_arr(i,j,k) = Real(0.0);
1272  fall_r_arr(i,j,k) = Real(0.0); fall_s_arr(i,j,k) = Real(0.0);
1273  fall_g_arr(i,j,k) = Real(0.0); fallc_arr(i,j,k) = Real(0.0);
1274  });
1275 
1276  // G3: xni ice crystal number concentration [lines 598-604]
1277  // WSM6-CPP TAG: XNI
1278  // legacy_group: G3
1279  // process: Ice crystal number concentration
1280  // compare_vars: xni, qi, den
1281  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1282  const Real tmp = den_arr(i,j,k)*amrex::max(qi_arr(i,j,k), Real(qmin));
1283  xni_arr(i,j,k) = amrex::min(
1284  amrex::max(Real(5.38e7)*std::sqrt(std::sqrt(tmp*tmp*tmp)), Real(1.e3)),
1285  Real(1.e6));
1286  });
1287 
1288  // G4: pack qrs_tmp, first slope_wsm6 [lines 610-618]
1289  // WSM6-CPP TAG: SLOPE1
1290  // legacy_group: G4
1291  // process: First slope calculation
1292  // compare_vars: rslope, rslope2, rslope3, rslopeb, falk, fall, work1
1293  const Real pidn0r_loc = m_pidn0r;
1294  const Real rslopermax_loc = m_rslopermax;
1295  const Real rsloperbmax_loc = m_rsloperbmax;
1296  const Real rsloper2max_loc = m_rsloper2max;
1297  const Real rsloper3max_loc = m_rsloper3max;
1298  const Real pidn0s_loc = m_pidn0s;
1299  const Real rslopesmax_loc = m_rslopesmax;
1300  const Real rslopesbmax_loc = m_rslopesbmax;
1301  const Real rslopes2max_loc = m_rslopes2max;
1302  const Real rslopes3max_loc = m_rslopes3max;
1303  const Real pidn0g_loc = m_pidn0g;
1304  const Real rslopegmax_loc = m_rslopegmax;
1305  const Real rslopegbmax_loc = m_rslopegbmax;
1306  const Real rslopeg2max_loc = m_rslopeg2max;
1307  const Real rslopeg3max_loc = m_rslopeg3max;
1308  const Real bvtg_loc = m_bvtg;
1309  const Real pi_wsm6_loc = m_pi_wsm6;
1310  const Real n0g_loc = m_n0g;
1311 
1312  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1313  qrs_tmp_r_arr(i,j,k) = qr_arr(i,j,k);
1314  qrs_tmp_s_arr(i,j,k) = qs_arr(i,j,k);
1315  qrs_tmp_g_arr(i,j,k) = qg_arr(i,j,k);
1316  Real dummy_n0sfac;
1317  wsm6_slope_rain_cell(
1318  qrs_tmp_r_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1319  pidn0r_loc, Real(qcrmin), rslopermax_loc, rsloperbmax_loc,
1320  rsloper2max_loc, rsloper3max_loc, Real(bvtr), Real(pvtr),
1321  rslope_r_arr(i,j,k), rslopeb_r_arr(i,j,k),
1322  rslope2_r_arr(i,j,k), rslope3_r_arr(i,j,k),
1323  work1_r_arr(i,j,k));
1324  wsm6_slope_snow_cell(
1325  qrs_tmp_s_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1326  t_arr(i,j,k), pidn0s_loc, Real(alpha_wsm6),
1327  Real(n0smax), Real(n0s), Real(t0c), Real(qcrmin),
1328  rslopesmax_loc, rslopesbmax_loc,
1329  rslopes2max_loc, rslopes3max_loc,
1330  Real(bvts), Real(pvts),
1331  rslope_s_arr(i,j,k), rslopeb_s_arr(i,j,k),
1332  rslope2_s_arr(i,j,k), rslope3_s_arr(i,j,k),
1333  work1_s_arr(i,j,k), dummy_n0sfac);
1334  wsm6_slope_graup_cell(
1335  qrs_tmp_g_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1336  pidn0g_loc, Real(qcrmin),
1337  rslopegmax_loc, rslopegbmax_loc,
1338  rslopeg2max_loc, rslopeg3max_loc,
1339  bvtg_loc, Real(pvtg),
1340  rslope_g_arr(i,j,k), rslopeb_g_arr(i,j,k),
1341  rslope2_g_arr(i,j,k), rslope3_g_arr(i,j,k),
1342  work1_g_arr(i,j,k));
1343  n0sfac_arr(i,j,k) = dummy_n0sfac;
1344  });
1345  // G5a-G5e: sedimentation setup, nislfv calls, and flux updates
1346  ParallelFor(box2d, [=] AMREX_GPU_DEVICE (int i, int j, int) {
1347  const int km_local = khi - klo + 1;
1348 
1349  constexpr Real qsum_min = Real(1.0e-15);
1350  auto den_col = [&](int k) -> amrex::Real& {
1351  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::den);
1352  };
1353  auto denfac_col = [&](int k) -> amrex::Real& {
1354  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::denfac);
1355  };
1356  auto t_col = [&](int k) -> amrex::Real& {
1357  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::tk);
1358  };
1359  auto dz_col = [&](int k) -> amrex::Real& {
1360  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::dz);
1361  };
1362  auto workr_col = [&](int k) -> amrex::Real& {
1363  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::workr_col);
1364  };
1365  auto worka_col = [&](int k) -> amrex::Real& {
1366  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::worka_col);
1367  };
1368  auto denqrs1_col = [&](int k) -> amrex::Real& {
1369  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::denqrs1_col);
1370  };
1371  auto denqrs2_col = [&](int k) -> amrex::Real& {
1372  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::denqrs2_col);
1373  };
1374  auto denqrs3_col = [&](int k) -> amrex::Real& {
1375  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::denqrs3_col);
1376  };
1377  auto qsum_col = [&](int k) -> amrex::Real& {
1378  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::qsum_col);
1379  };
1380  Real delqrs1_col = Real(0.0);
1381  Real delqrs2_col = Real(0.0);
1382  Real delqrs3_col = Real(0.0);
1383 
1384  // G5a: pack sedimentation work arrays
1385  for (int k = klo; k <= khi; ++k) {
1386  const int kk = k - klo;
1387  den_col(kk) = den_arr(i,j,k);
1388  denfac_col(kk) = denfac_arr(i,j,k);
1389  t_col(kk) = t_arr(i,j,k);
1390  dz_col(kk) = delz_tmp_arr(i,j,k);
1391  workr_col(kk) = work1_r_arr(i,j,k);
1392  qsum_col(kk) = amrex::max(qs_arr(i,j,k) + qg_arr(i,j,k), qsum_min);
1393  if (qsum_col(kk) > qsum_min) {
1394  worka_col(kk) = (work1_s_arr(i,j,k) * qs_arr(i,j,k)
1395  + work1_g_arr(i,j,k) * qg_arr(i,j,k))
1396  / qsum_col(kk);
1397  } else {
1398  worka_col(kk) = Real(0.0);
1399  }
1400  denqrs1_col(kk) = den_col(kk) * qr_arr(i,j,k);
1401  denqrs2_col(kk) = den_col(kk) * qs_arr(i,j,k);
1402  denqrs3_col(kk) = den_col(kk) * qg_arr(i,j,k);
1403  if (qr_arr(i,j,k) <= Real(0.0)) {
1404  workr_col(kk) = Real(0.0);
1405  }
1406  }
1407 
1408  // G5b: rain sedimentation
1409  wsm6_nislfv_rain_plm_scratch(
1410  km_local,
1411  WSM6SedCellScratch::workr_col,
1412  WSM6SedCellScratch::denqrs1_col,
1413  &delqrs1_col, dtcld, 1,
1414  sed_cell_scratch_arr, sed_node_scratch_arr, i, j, klo);
1415  // Strict Rule 30 snapshot: immediately after G5b
1416  for (int k = klo; k <= khi; ++k) {
1417  const int kk = k - klo;
1418  nislfv_r_diag_arr(i,j,k,0) = amrex::max(denqrs1_col(kk) / den_col(kk), Real(0.0));
1419  nislfv_r_diag_arr(i,j,k,1) = denqrs1_col(kk) * workr_col(kk) / delz_arr(i,j,k);
1420  nislfv_r_diag_arr(i,j,k,2) = workr_col(kk);
1421  nislfv_r_diag_arr(i,j,k,3) = denqrs1_col(kk);
1422  nislfv_r_diag_arr(i,j,k,4) = den_col(kk);
1423  nislfv_r_diag_arr(i,j,k,5) = denfac_col(kk);
1424  }
1425 
1426  // G5c: snow + graupel sedimentation
1427  wsm6_nislfv_rain_plm6_scratch(
1428  km_local,
1429  WSM6SedCellScratch::worka_col,
1430  WSM6SedCellScratch::denqrs2_col,
1431  WSM6SedCellScratch::denqrs3_col,
1432  &delqrs2_col, &delqrs3_col, dtcld, 1,
1433  sed_cell_scratch_arr, sed_node_scratch_arr, i, j, klo);
1434  // Strict Rule 30 snapshot: immediately after G5c
1435  for (int k = klo; k <= khi; ++k) {
1436  const int kk = k - klo;
1437  nislfv_sg_diag_arr(i,j,k,0) = amrex::max(denqrs2_col(kk) / den_col(kk), Real(0.0));
1438  nislfv_sg_diag_arr(i,j,k,1) = amrex::max(denqrs3_col(kk) / den_col(kk), Real(0.0));
1439  nislfv_sg_diag_arr(i,j,k,2) = denqrs2_col(kk) * worka_col(kk) / delz_arr(i,j,k);
1440  nislfv_sg_diag_arr(i,j,k,3) = denqrs3_col(kk) * worka_col(kk) / delz_arr(i,j,k);
1441  nislfv_sg_diag_arr(i,j,k,4) = denqrs2_col(kk);
1442  nislfv_sg_diag_arr(i,j,k,5) = denqrs3_col(kk);
1443  }
1444 
1445  // G5d: update species and fall speeds
1446  for (int k = klo; k <= khi; ++k) {
1447  const int kk = k - klo;
1448  qsum_arr(i,j,k) = qsum_col(kk);
1449  workr_arr(i,j,k) = workr_col(kk);
1450  worka_arr(i,j,k) = worka_col(kk);
1451  denqrs1_arr(i,j,k) = denqrs1_col(kk);
1452  denqrs2_arr(i,j,k) = denqrs2_col(kk);
1453  denqrs3_arr(i,j,k) = denqrs3_col(kk);
1454  qr_arr(i,j,k) = amrex::max(denqrs1_col(kk) / den_col(kk), Real(0.0));
1455  qs_arr(i,j,k) = amrex::max(denqrs2_col(kk) / den_col(kk), Real(0.0));
1456  qg_arr(i,j,k) = amrex::max(denqrs3_col(kk) / den_col(kk), Real(0.0));
1457  fall_r_arr(i,j,k) = denqrs1_col(kk) * workr_col(kk) / delz_arr(i,j,k);
1458  fall_s_arr(i,j,k) = denqrs2_col(kk) * worka_col(kk) / delz_arr(i,j,k);
1459  fall_g_arr(i,j,k) = denqrs3_col(kk) * worka_col(kk) / delz_arr(i,j,k);
1460  }
1461 
1462  // G5e: slab fall fluxes at the lower boundary
1463  delqrs1_arr(i,j,0) = delqrs1_col / delz_arr(i,j,klo) / dtcld;
1464  delqrs2_arr(i,j,0) = delqrs2_col / delz_arr(i,j,klo) / dtcld;
1465  delqrs3_arr(i,j,0) = delqrs3_col / delz_arr(i,j,klo) / dtcld;
1466  fall_r_arr(i,j,klo) = delqrs1_arr(i,j,0);
1467  fall_s_arr(i,j,klo) = delqrs2_arr(i,j,0);
1468  fall_g_arr(i,j,klo) = delqrs3_arr(i,j,0);
1469  });
1470 
1471  // G6: repack qrs_tmp, second slope_wsm6 [lines 655-663]
1472  // slope params updated after sedimentation moved mass
1473  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1474  qrs_tmp_r_arr(i,j,k) = qr_arr(i,j,k);
1475  qrs_tmp_s_arr(i,j,k) = qs_arr(i,j,k);
1476  qrs_tmp_g_arr(i,j,k) = qg_arr(i,j,k);
1477  Real dummy_n0sfac;
1478  wsm6_slope_rain_cell(
1479  qrs_tmp_r_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1480  pidn0r_loc, Real(qcrmin), rslopermax_loc, rsloperbmax_loc,
1481  rsloper2max_loc, rsloper3max_loc, Real(bvtr), Real(pvtr),
1482  rslope_r_arr(i,j,k), rslopeb_r_arr(i,j,k),
1483  rslope2_r_arr(i,j,k), rslope3_r_arr(i,j,k),
1484  work1_r_arr(i,j,k));
1485  wsm6_slope_snow_cell(
1486  qrs_tmp_s_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1487  t_arr(i,j,k), pidn0s_loc, Real(alpha_wsm6),
1488  Real(n0smax), Real(n0s), Real(t0c), Real(qcrmin),
1489  rslopesmax_loc, rslopesbmax_loc,
1490  rslopes2max_loc, rslopes3max_loc,
1491  Real(bvts), Real(pvts),
1492  rslope_s_arr(i,j,k), rslopeb_s_arr(i,j,k),
1493  rslope2_s_arr(i,j,k), rslope3_s_arr(i,j,k),
1494  work1_s_arr(i,j,k), dummy_n0sfac);
1495  wsm6_slope_graup_cell(
1496  qrs_tmp_g_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1497  pidn0g_loc, Real(qcrmin),
1498  rslopegmax_loc, rslopegbmax_loc,
1499  rslopeg2max_loc, rslopeg3max_loc,
1500  bvtg_loc, Real(pvtg),
1501  rslope_g_arr(i,j,k), rslopeb_g_arr(i,j,k),
1502  rslope2_g_arr(i,j,k), rslope3_g_arr(i,j,k),
1503  work1_g_arr(i,j,k));
1504  });
1505 
1506  // G7: melting (T>T0 only) [lines 665-704]
1507  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1508  const Real supcol = Real(t0c) - t_arr(i,j,k);
1509  n0sfac_arr(i,j,k) = amrex::max(
1510  amrex::min(std::exp(Real(alpha_wsm6) * supcol),
1511  Real(n0smax) / Real(n0s)),
1512  Real(1.0));
1513 
1514  if (t_arr(i,j,k) > Real(t0c)) {
1515  const Real xlf = Real(xlf0);
1516  work2_arr(i,j,k) = wsm6_venfac(
1517  p_arr(i,j,k), t_arr(i,j,k), den_arr(i,j,k), Real(den0));
1518 
1519  if (qs_arr(i,j,k) > Real(0.0)) {
1520  const Real coeres =
1521  rslope2_s_arr(i,j,k) *
1522  std::sqrt(rslope_s_arr(i,j,k) * rslopeb_s_arr(i,j,k));
1523  psmlt_arr(i,j,k) =
1524  wsm6_xka(t_arr(i,j,k), den_arr(i,j,k)) / xlf *
1525  (Real(t0c) - t_arr(i,j,k)) * pi_wsm6_loc * Real(0.5) *
1526  n0sfac_arr(i,j,k) *
1527  (Real(precs1) * rslope2_s_arr(i,j,k) +
1528  Real(precs2) * work2_arr(i,j,k) * coeres) /
1529  den_arr(i,j,k);
1530  psmlt_arr(i,j,k) = amrex::min(
1531  amrex::max(psmlt_arr(i,j,k) * dtcld,
1532  -qs_arr(i,j,k)),
1533  Real(0.0));
1534  qs_arr(i,j,k) = qs_arr(i,j,k) + psmlt_arr(i,j,k);
1535  qr_arr(i,j,k) = qr_arr(i,j,k) - psmlt_arr(i,j,k);
1536  t_arr(i,j,k) = t_arr(i,j,k) + xlf / cpm_arr(i,j,k) * psmlt_arr(i,j,k);
1537  }
1538 
1539  if (qg_arr(i,j,k) > Real(0.0)) {
1540  const Real coeres =
1541  rslope2_g_arr(i,j,k) *
1542  std::sqrt(rslope_g_arr(i,j,k) * rslopeb_g_arr(i,j,k));
1543  pgmlt_arr(i,j,k) =
1544  wsm6_xka(t_arr(i,j,k), den_arr(i,j,k)) / xlf *
1545  (Real(t0c) - t_arr(i,j,k)) *
1546  (Real(precg1) * rslope2_g_arr(i,j,k) +
1547  Real(precg2) * work2_arr(i,j,k) * coeres) /
1548  den_arr(i,j,k);
1549  pgmlt_arr(i,j,k) = amrex::min(
1550  amrex::max(pgmlt_arr(i,j,k) * dtcld,
1551  -qg_arr(i,j,k)),
1552  Real(0.0));
1553  qg_arr(i,j,k) = qg_arr(i,j,k) + pgmlt_arr(i,j,k);
1554  qr_arr(i,j,k) = qr_arr(i,j,k) - pgmlt_arr(i,j,k);
1555  t_arr(i,j,k) = t_arr(i,j,k) + xlf / cpm_arr(i,j,k) * pgmlt_arr(i,j,k);
1556  }
1557  }
1558  });
1559 
1560  // G8: cloud ice sedimentation/fallout [lines 708-735]
1561  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1562  if (qi_arr(i,j,k) <= Real(0.0)) {
1563  work1c_arr(i,j,k) = Real(0.0);
1564  } else {
1565  const Real tmp = den_arr(i,j,k) *
1566  amrex::max(qi_arr(i,j,k), Real(qmin));
1567  const Real xni = amrex::min(
1568  amrex::max(Real(5.38e7) *
1569  std::sqrt(std::sqrt(tmp*tmp*tmp)),
1570  Real(1.0e3)),
1571  Real(1.0e6));
1572  const Real xmi = den_arr(i,j,k) * qi_arr(i,j,k) / xni;
1573  const Real diameter = amrex::max(
1574  amrex::min(Real(dicon) * std::sqrt(xmi), Real(dimax)),
1575  Real(1.0e-25));
1576  work1c_arr(i,j,k) = Real(1.49e4) *
1577  std::exp(std::log(diameter) * Real(1.31));
1578  }
1579  denqci_arr(i,j,k) = den_arr(i,j,k) * qi_arr(i,j,k);
1580  });
1581 
1582  ParallelFor(box2d, [=] AMREX_GPU_DEVICE (int i, int j, int) {
1583  const int km_local = khi - klo + 1;
1584  auto den_col = [&](int k) -> amrex::Real& {
1585  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::den);
1586  };
1587  auto denfac_col = [&](int k) -> amrex::Real& {
1588  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::denfac);
1589  };
1590  auto t_col = [&](int k) -> amrex::Real& {
1591  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::tk);
1592  };
1593  auto dz_col = [&](int k) -> amrex::Real& {
1594  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::dz);
1595  };
1596  auto work1c_col = [&](int k) -> amrex::Real& {
1597  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::work1c_col);
1598  };
1599  auto denqci_col = [&](int k) -> amrex::Real& {
1600  return sed_cell_scratch_arr(i, j, klo + k, WSM6SedCellScratch::denqci_col);
1601  };
1602  Real delqi_col = Real(0.0);
1603 
1604  for (int k = klo; k <= khi; ++k) {
1605  const int kk = k - klo;
1606  den_col(kk) = den_arr(i,j,k);
1607  denfac_col(kk) = denfac_arr(i,j,k);
1608  t_col(kk) = t_arr(i,j,k);
1609  dz_col(kk) = delz_tmp_arr(i,j,k);
1610  work1c_col(kk) = work1c_arr(i,j,k);
1611  denqci_col(kk) = denqci_arr(i,j,k);
1612  }
1613 
1614  wsm6_nislfv_rain_plm_scratch(
1615  km_local,
1616  WSM6SedCellScratch::work1c_col,
1617  WSM6SedCellScratch::denqci_col,
1618  &delqi_col, dtcld, 0,
1619  sed_cell_scratch_arr, sed_node_scratch_arr, i, j, klo);
1620 
1621  for (int k = klo; k <= khi; ++k) {
1622  const int kk = k - klo;
1623  work1c_arr(i,j,k) = work1c_col(kk);
1624  denqci_arr(i,j,k) = denqci_col(kk);
1625  qi_arr(i,j,k) = amrex::max(
1626  denqci_col(kk) / den_col(kk), Real(0.0));
1627  }
1628 
1629  delqi_arr(i,j,0) = delqi_col / delz_arr(i,j,klo) / dtcld;
1630  fallc_arr(i,j,klo) = delqi_arr(i,j,0);
1631  });
1632 
1633  // G9: surface precipitation accumulation [lines 741-770]
1634  ParallelFor(box2d, [=] AMREX_GPU_DEVICE (int i, int j, int) {
1635  const Real fallsum =
1636  fall_r_arr(i,j,klo) + fall_s_arr(i,j,klo) +
1637  fall_g_arr(i,j,klo) + fallc_arr(i,j,klo);
1638  const Real fallsum_qsi = fall_s_arr(i,j,klo) + fallc_arr(i,j,klo);
1639  const Real fallsum_qg = fall_g_arr(i,j,klo);
1640  const Real precip = delz_arr(i,j,klo) / denr * dtcld * Real(1000.0);
1641 
1642  if (fallsum > Real(0.0)) {
1643  rainncv_arr(i,j,0) += fallsum * precip;
1644  rainacc_arr(i,j,0) += fallsum * precip;
1645  }
1646  if (fallsum_qsi > Real(0.0)) {
1647  tstepsnow_arr(i,j,0) += fallsum_qsi * precip;
1648  snowncv_arr(i,j,0) += fallsum_qsi * precip;
1649  snowacc_arr(i,j,0) += fallsum_qsi * precip;
1650  }
1651  if (fallsum_qg > Real(0.0)) {
1652  tstepgraup_arr(i,j,0) += fallsum_qg * precip;
1653  graupelncv_arr(i,j,0) += fallsum_qg * precip;
1654  graupacc_arr(i,j,0) += fallsum_qg * precip;
1655  }
1656  if (fallsum > Real(0.0)) {
1657  sr_arr(i,j,0) =
1658  (snowncv_arr(i,j,0) + graupelncv_arr(i,j,0)) /
1659  (rainncv_arr(i,j,0) + Real(1.0e-12));
1660  }
1661  });
1662 
1663  constexpr Real pi = Real(3.141592653589793238462643383279502884);
1664 
1665  // G10: instantaneous phase changes [lines 774-830]
1666  // WSM6-CPP TAG: PHASE
1667  // legacy_group: G10
1668  // process: Instantaneous phase changes
1669  // compare_vars: t, q, qci, qrs
1670  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1671  const Real supcol = Real(t0c) - t_arr(i,j,k);
1672  const Real xlf = (supcol < Real(0.0))
1673  ? Real(xlf0)
1674  : Real(xls) - xl_arr(i,j,k);
1675  pimlt_arr(i,j,k) = Real(0.0);
1676  pihmf_arr(i,j,k) = Real(0.0);
1677  pihtf_arr(i,j,k) = Real(0.0);
1678  pgfrz_arr(i,j,k) = Real(0.0);
1679 
1680  if (supcol < Real(0.0) && qi_arr(i,j,k) > Real(0.0)) {
1681  pimlt_arr(i,j,k) = qi_arr(i,j,k);
1682  qc_arr(i,j,k) = qc_arr(i,j,k) + qi_arr(i,j,k);
1683  t_arr(i,j,k) = t_arr(i,j,k) - xlf / cpm_arr(i,j,k) * qi_arr(i,j,k);
1684  qi_arr(i,j,k) = Real(0.0);
1685  }
1686 
1687  if (supcol > Real(40.0) && qc_arr(i,j,k) > Real(0.0)) {
1688  pihmf_arr(i,j,k) = qc_arr(i,j,k);
1689  qi_arr(i,j,k) = qi_arr(i,j,k) + qc_arr(i,j,k);
1690  t_arr(i,j,k) = t_arr(i,j,k) + xlf / cpm_arr(i,j,k) * qc_arr(i,j,k);
1691  qc_arr(i,j,k) = Real(0.0);
1692  }
1693 
1694  if (supcol > Real(0.0) && qc_arr(i,j,k) > Real(qmin)) {
1695  const Real supcolt = amrex::min(supcol, Real(50.0));
1696  const Real pfrzdtc = amrex::min(
1697  Real(pfrz1) *
1698  (std::exp(Real(pfrz2) * supcolt) - Real(1.0)) *
1699  den_arr(i,j,k) / Real(denr) / Real(WSM6::xncr) *
1700  qc_arr(i,j,k) * qc_arr(i,j,k) * dtcld,
1701  qc_arr(i,j,k));
1702  pihtf_arr(i,j,k) = pfrzdtc;
1703  qi_arr(i,j,k) = qi_arr(i,j,k) + pfrzdtc;
1704  t_arr(i,j,k) = t_arr(i,j,k) + xlf / cpm_arr(i,j,k) * pfrzdtc;
1705  qc_arr(i,j,k) = qc_arr(i,j,k) - pfrzdtc;
1706  }
1707 
1708  if (supcol > Real(0.0) && qr_arr(i,j,k) > Real(0.0)) {
1709  Real temp = rslope3_r_arr(i,j,k);
1710  temp = temp * temp * rslope_r_arr(i,j,k);
1711  const Real supcolt = amrex::min(supcol, Real(50.0));
1712  const Real pfrzdtr = amrex::min(
1713  Real(20.0) * pi * pi * Real(pfrz1) * Real(WSM6::n0r) *
1714  Real(denr) / den_arr(i,j,k) *
1715  (std::exp(Real(pfrz2) * supcolt) - Real(1.0)) *
1716  temp * dtcld,
1717  qr_arr(i,j,k));
1718  pgfrz_arr(i,j,k) = pfrzdtr;
1719  qg_arr(i,j,k) = qg_arr(i,j,k) + pfrzdtr;
1720  t_arr(i,j,k) = t_arr(i,j,k) + xlf / cpm_arr(i,j,k) * pfrzdtr;
1721  qr_arr(i,j,k) = qr_arr(i,j,k) - pfrzdtr;
1722  }
1723  });
1724  // G11: third slope_wsm6 call [lines 836-844]
1725  // WSM6-CPP TAG: SLOPE3
1726  // legacy_group: G11
1727  // process: Third slope calculation
1728  // compare_vars: rslope, rslope2, rslope3, rslopeb, falk, fall, work1
1729  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1730  qrs_tmp_r_arr(i,j,k) = qr_arr(i,j,k);
1731  qrs_tmp_s_arr(i,j,k) = qs_arr(i,j,k);
1732  qrs_tmp_g_arr(i,j,k) = qg_arr(i,j,k);
1733  Real dummy_n0sfac;
1734  wsm6_slope_rain_cell(
1735  qrs_tmp_r_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1736  pidn0r_loc, Real(qcrmin), rslopermax_loc, rsloperbmax_loc,
1737  rsloper2max_loc, rsloper3max_loc, Real(bvtr), Real(pvtr),
1738  rslope_r_arr(i,j,k), rslopeb_r_arr(i,j,k),
1739  rslope2_r_arr(i,j,k), rslope3_r_arr(i,j,k),
1740  work1_r_arr(i,j,k));
1741  wsm6_slope_snow_cell(
1742  qrs_tmp_s_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1743  t_arr(i,j,k), pidn0s_loc, Real(alpha_wsm6),
1744  Real(n0smax), Real(n0s), Real(t0c), Real(qcrmin),
1745  rslopesmax_loc, rslopesbmax_loc,
1746  rslopes2max_loc, rslopes3max_loc,
1747  Real(bvts), Real(pvts),
1748  rslope_s_arr(i,j,k), rslopeb_s_arr(i,j,k),
1749  rslope2_s_arr(i,j,k), rslope3_s_arr(i,j,k),
1750  work1_s_arr(i,j,k), dummy_n0sfac);
1751  wsm6_slope_graup_cell(
1752  qrs_tmp_g_arr(i,j,k), den_arr(i,j,k), denfac_arr(i,j,k),
1753  pidn0g_loc, Real(qcrmin),
1754  rslopegmax_loc, rslopegbmax_loc,
1755  rslopeg2max_loc, rslopeg3max_loc,
1756  bvtg_loc, Real(pvtg),
1757  rslope_g_arr(i,j,k), rslopeb_g_arr(i,j,k),
1758  rslope2_g_arr(i,j,k), rslope3_g_arr(i,j,k),
1759  work1_g_arr(i,j,k));
1760  n0sfac_arr(i,j,k) = dummy_n0sfac;
1761  });
1762  // G12: workdiffw, workdiffi, work2 [lines 851-857]
1763  // WSM6-CPP TAG: DIFF_PREP
1764  // legacy_group: G12
1765  // process: Prepare diffusion/work terms
1766  // compare_vars: workdiffw, workdiffi, work2
1767  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1768  workdiffw_arr(i,j,k) = wsm6_diffac(
1769  xl_arr(i,j,k), p_arr(i,j,k), t_arr(i,j,k),
1770  den_arr(i,j,k), qsatw_arr(i,j,k), Real(rv));
1771  workdiffi_arr(i,j,k) = wsm6_diffac(
1772  Real(xls), p_arr(i,j,k), t_arr(i,j,k),
1773  den_arr(i,j,k), qsati_arr(i,j,k), Real(rv));
1774  work2_arr(i,j,k) = wsm6_venfac(
1775  p_arr(i,j,k), t_arr(i,j,k), den_arr(i,j,k), Real(den0));
1776  });
1777 
1778  // G13a: warm rain — praut, pracw, prevp [lines 867-903]
1779  {
1780  const Real dtcld_l = dtcld;
1781  const Real qmin_l = Real(qmin);
1782  const Real qcrmin_l= Real(qcrmin);
1783  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1784  Real supsat = amrex::max(qv_arr(i,j,k), qmin_l)
1785  - qsatw_arr(i,j,k);
1786  Real satdt = supsat / dtcld_l;
1787  // praut: C->R autoconversion
1788  if (qc_arr(i,j,k) > Real(qc0)) {
1789  praut_arr(i,j,k) = amrex::min(
1790  Real(qck1)*std::pow(qc_arr(i,j,k), Real(7.0/3.0)),
1791  qc_arr(i,j,k)/dtcld_l);
1792  }
1793  // pracw: C->R accretion by rain
1794  if (qr_arr(i,j,k) > qcrmin_l && qc_arr(i,j,k) > qmin_l) {
1795  pracw_arr(i,j,k) = amrex::min(
1796  Real(pacrr)*rslope3_r_arr(i,j,k)
1797  *rslopeb_r_arr(i,j,k)
1798  *qc_arr(i,j,k)*denfac_arr(i,j,k),
1799  qc_arr(i,j,k)/dtcld_l);
1800  }
1801  // prevp: R evaporation/condensation
1802  if (qr_arr(i,j,k) > Real(0.0)) {
1803  Real coeres = rslope2_r_arr(i,j,k)*std::sqrt(
1804  rslope_r_arr(i,j,k)*rslopeb_r_arr(i,j,k));
1805  Real rate = (rhw_arr(i,j,k)-Real(1.0))
1806  *(Real(precr1)*rslope2_r_arr(i,j,k)
1807  + Real(precr2)*work2_arr(i,j,k)*coeres)
1808  / workdiffw_arr(i,j,k);
1809  if (rate < Real(0.0)) {
1810  rate = amrex::max(rate, -qr_arr(i,j,k)/dtcld_l);
1811  rate = amrex::max(rate, satdt/Real(2.0));
1812  } else {
1813  rate = amrex::min(rate, satdt/Real(2.0));
1814  }
1815  prevp_arr(i,j,k) = rate;
1816  }
1817  });
1818  }
1819  // G13b: cold-rain, mixed-phase, and ice deposition/nucleation
1820  // [lines 904-1062, 1075-1192]
1821  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
1822  const Real supcol = Real(t0c) - t_arr(i,j,k);
1823  n0sfac_arr(i,j,k) = amrex::max(
1824  amrex::min(std::exp(Real(alpha_wsm6) * supcol),
1825  Real(n0smax) / Real(n0s)),
1826  Real(1.0));
1827 
1828  const Real supsat = amrex::max(qv_arr(i,j,k), Real(qmin))
1829  - qsati_arr(i,j,k);
1830  const Real satdt = supsat / dtcld;
1831  int ifsat = 0;
1832 
1833  const Real tmp = den_arr(i,j,k)
1834  * amrex::max(qi_arr(i,j,k), Real(qmin));
1835  const Real temp = std::sqrt(std::sqrt(tmp * tmp * tmp));
1836  xni_arr(i,j,k) = amrex::min(
1837  amrex::max(Real(5.38e7) * temp, Real(1.0e3)),
1838  Real(1.0e6));
1839 
1840  const Real xni = xni_arr(i,j,k);
1841  const Real eacrs = std::exp(Real(0.07) * (-supcol));
1842 
1843  const Real xmi = den_arr(i,j,k) * qi_arr(i,j,k) / xni;
1844  const Real diameter = amrex::min(
1845  Real(dicon) * std::sqrt(xmi), Real(dimax));
1846  const Real vt2i = Real(1.49e4) * std::pow(diameter, Real(1.31));
1847  const Real vt2r = Real(pvtr) * rslopeb_r_arr(i,j,k)
1848  * denfac_arr(i,j,k);
1849  const Real vt2s = Real(pvts) * rslopeb_s_arr(i,j,k)
1850  * denfac_arr(i,j,k);
1851  const Real vt2g = Real(pvtg) * rslopeb_g_arr(i,j,k)
1852  * denfac_arr(i,j,k);
1853 
1854  const Real qsum = amrex::max(qsum_arr(i,j,k), Real(1.0e-15));
1855  const Real vt2ave = (qsum > Real(1.0e-15))
1856  ? (vt2s * qs_arr(i,j,k) + vt2g * qg_arr(i,j,k)) / qsum
1857  : Real(0.0);
1858 
1859  if (supcol > Real(0.0) && qi_arr(i,j,k) > Real(qmin)) {
1860  if (qr_arr(i,j,k) > Real(qcrmin)) {
1861  const Real acrfac =
1862  Real(2.0) * rslope3_r_arr(i,j,k)
1863  + Real(2.0) * diameter * rslope2_r_arr(i,j,k)
1864  + diameter * diameter * rslope_r_arr(i,j,k);
1865 
1866  // WSM6-CPP TAG: PRACI
1867  // legacy_group: G13b
1868  // process: Accretion of cloud ice by rain
1869  // compare_vars: praci, qi, qr, den
1870  praci_arr(i,j,k) = Real(pi) * qi_arr(i,j,k) * Real(n0r)
1871  * std::abs(vt2r - vt2i) * acrfac / Real(4.0);
1872  praci_arr(i,j,k) *= std::pow(
1873  amrex::min(
1874  amrex::max(Real(0.0), qr_arr(i,j,k) / qi_arr(i,j,k)),
1875  Real(1.0)),
1876  Real(2.0));
1877  praci_arr(i,j,k) = amrex::min(
1878  praci_arr(i,j,k), qi_arr(i,j,k) / dtcld);
1879 
1880  piacr_arr(i,j,k) = Real(pi) * Real(pi) * Real(avtr)
1881  * Real(n0r) * Real(denr) * xni * denfac_arr(i,j,k)
1882  * Real(g6pbr) * rslope3_r_arr(i,j,k)
1883  * rslope3_r_arr(i,j,k) * rslopeb_r_arr(i,j,k)
1884  / Real(24.0) / den_arr(i,j,k);
1885  piacr_arr(i,j,k) *= std::pow(
1886  amrex::min(
1887  amrex::max(Real(0.0), qi_arr(i,j,k) / qr_arr(i,j,k)),
1888  Real(1.0)),
1889  Real(2.0));
1890  piacr_arr(i,j,k) = amrex::min(
1891  piacr_arr(i,j,k), qr_arr(i,j,k) / dtcld);
1892  }
1893 
1894  if (qs_arr(i,j,k) > Real(qcrmin)) {
1895  const Real acrfac =
1896  Real(2.0) * rslope3_s_arr(i,j,k)
1897  + Real(2.0) * diameter * rslope2_s_arr(i,j,k)
1898  + diameter * diameter * rslope_s_arr(i,j,k);
1899  // WSM6-CPP TAG: PSACI
1900  // legacy_group: G13e
1901  // process: Accretion of cloud ice by snow
1902  // compare_vars: psaci, qi, qs, den
1903  psaci_arr(i,j,k) = Real(pi) * qi_arr(i,j,k) * eacrs
1904  * Real(n0s) * n0sfac_arr(i,j,k)
1905  * std::abs(vt2ave - vt2i) * acrfac / Real(4.0);
1906  psaci_arr(i,j,k) = amrex::min(
1907  psaci_arr(i,j,k), qi_arr(i,j,k) / dtcld);
1908  }
1909 
1910  if (qg_arr(i,j,k) > Real(qcrmin)) {
1911  const Real egi = std::exp(Real(0.07) * (-supcol));
1912  const Real acrfac =
1913  Real(2.0) * rslope3_g_arr(i,j,k)
1914  + Real(2.0) * diameter * rslope2_g_arr(i,j,k)
1915  + diameter * diameter * rslope_g_arr(i,j,k);
1916  pgaci_arr(i,j,k) = Real(pi) * egi * qi_arr(i,j,k)
1917  * n0g_loc * std::abs(vt2ave - vt2i) * acrfac
1918  / Real(4.0);
1919  pgaci_arr(i,j,k) = amrex::min(
1920  pgaci_arr(i,j,k), qi_arr(i,j,k) / dtcld);
1921  }
1922  }
1923 
1924  if (qs_arr(i,j,k) > Real(qcrmin) && qc_arr(i,j,k) > Real(qmin)) {
1925  psacw_arr(i,j,k) = amrex::min(
1926  Real(pacrc) * n0sfac_arr(i,j,k) * rslope3_s_arr(i,j,k)
1927  * rslopeb_s_arr(i,j,k)
1928  * std::pow(
1929  amrex::min(
1930  amrex::max(Real(0.0), qs_arr(i,j,k) / qc_arr(i,j,k)),
1931  Real(1.0)),
1932  Real(2.0))
1933  * qc_arr(i,j,k) * denfac_arr(i,j,k),
1934  qc_arr(i,j,k) / dtcld);
1935  }
1936 
1937  if (qg_arr(i,j,k) > Real(qcrmin) && qc_arr(i,j,k) > Real(qmin)) {
1938  pgacw_arr(i,j,k) = amrex::min(
1939  Real(pacrg) * rslope3_g_arr(i,j,k) * rslopeb_g_arr(i,j,k)
1940  * std::pow(
1941  amrex::min(
1942  amrex::max(Real(0.0), qg_arr(i,j,k) / qc_arr(i,j,k)),
1943  Real(1.0)),
1944  Real(2.0))
1945  * qc_arr(i,j,k) * denfac_arr(i,j,k),
1946  qc_arr(i,j,k) / dtcld);
1947  }
1948 
1949  if (qsum > Real(1.0e-15)) {
1950  paacw_arr(i,j,k) = (qs_arr(i,j,k) * psacw_arr(i,j,k)
1951  + qg_arr(i,j,k) * pgacw_arr(i,j,k))
1952  / qsum;
1953  }
1954 
1955  if (qs_arr(i,j,k) > Real(qcrmin) && qr_arr(i,j,k) > Real(qcrmin)) {
1956  if (supcol > Real(0.0)) {
1957  const Real acrfac =
1958  Real(5.0) * rslope3_s_arr(i,j,k) * rslope3_s_arr(i,j,k)
1959  * rslope_r_arr(i,j,k)
1960  + Real(2.0) * rslope3_s_arr(i,j,k) * rslope2_s_arr(i,j,k)
1961  * rslope2_r_arr(i,j,k)
1962  + Real(0.5) * rslope2_s_arr(i,j,k) * rslope2_s_arr(i,j,k)
1963  * rslope3_r_arr(i,j,k);
1964  // WSM6-CPP TAG: PRACS
1965  // legacy_group: G13f
1966  // process: Accretion of snow by rain / rain-snow interaction
1967  // compare_vars: pracs, qr, qs, den
1968  pracs_arr(i,j,k) = Real(pi) * Real(pi) * Real(n0r)
1969  * Real(n0s) * n0sfac_arr(i,j,k)
1970  * std::abs(vt2r - vt2ave) * (Real(dens_snow) / den_arr(i,j,k))
1971  * acrfac;
1972  pracs_arr(i,j,k) *= std::pow(
1973  amrex::min(
1974  amrex::max(Real(0.0), qr_arr(i,j,k) / qs_arr(i,j,k)),
1975  Real(1.0)),
1976  Real(2.0));
1977  pracs_arr(i,j,k) = amrex::min(
1978  pracs_arr(i,j,k), qs_arr(i,j,k) / dtcld);
1979  }
1980 
1981  {
1982  const Real acrfac =
1983  Real(5.0) * rslope3_r_arr(i,j,k) * rslope3_r_arr(i,j,k)
1984  * rslope_s_arr(i,j,k)
1985  + Real(2.0) * rslope3_r_arr(i,j,k) * rslope2_r_arr(i,j,k)
1986  * rslope2_s_arr(i,j,k)
1987  + Real(0.5) * rslope2_r_arr(i,j,k) * rslope2_r_arr(i,j,k)
1988  * rslope3_s_arr(i,j,k);
1989  psacr_arr(i,j,k) = Real(pi) * Real(pi) * Real(n0r)
1990  * Real(n0s) * n0sfac_arr(i,j,k)
1991  * std::abs(vt2ave - vt2r) * (Real(denr) / den_arr(i,j,k))
1992  * acrfac;
1993  psacr_arr(i,j,k) *= std::pow(
1994  amrex::min(
1995  amrex::max(Real(0.0), qs_arr(i,j,k) / qr_arr(i,j,k)),
1996  Real(1.0)),
1997  Real(2.0));
1998  psacr_arr(i,j,k) = amrex::min(
1999  psacr_arr(i,j,k), qr_arr(i,j,k) / dtcld);
2000  }
2001 
2002  if (qg_arr(i,j,k) > Real(qcrmin)) {
2003  const Real acrfac =
2004  Real(5.0) * rslope3_r_arr(i,j,k) * rslope3_r_arr(i,j,k)
2005  * rslope_g_arr(i,j,k)
2006  + Real(2.0) * rslope3_r_arr(i,j,k) * rslope2_r_arr(i,j,k)
2007  * rslope2_g_arr(i,j,k)
2008  + Real(0.5) * rslope2_r_arr(i,j,k) * rslope2_r_arr(i,j,k)
2009  * rslope3_g_arr(i,j,k);
2010  pgacr_arr(i,j,k) = Real(pi) * Real(pi) * Real(n0r)
2011  * n0g_loc * std::abs(vt2ave - vt2r)
2012  * (Real(denr) / den_arr(i,j,k)) * acrfac;
2013  pgacr_arr(i,j,k) *= std::pow(
2014  amrex::min(
2015  amrex::max(Real(0.0), qg_arr(i,j,k) / qr_arr(i,j,k)),
2016  Real(1.0)),
2017  Real(2.0));
2018  pgacr_arr(i,j,k) = amrex::min(
2019  pgacr_arr(i,j,k), qr_arr(i,j,k) / dtcld);
2020  }
2021  }
2022 
2023  if (qg_arr(i,j,k) > Real(qcrmin) && qs_arr(i,j,k) > Real(qcrmin)) {
2024  pgacs_arr(i,j,k) = Real(0.0);
2025  }
2026 
2027  if (supcol <= Real(0.0)) {
2028  const Real xlf = Real(xlf0);
2029  if (qs_arr(i,j,k) > Real(0.0)) {
2030  // WSM6-CPP TAG: PSEML
2031  // legacy_group: G13g
2032  // process: Snow evaporation/sublimation
2033  // compare_vars: pseml, qs, qv, qsat, den
2034  pseml_arr(i,j,k) = amrex::min(
2035  amrex::max(
2036  Real(cliq) * supcol
2037  * (paacw_arr(i,j,k) + psacr_arr(i,j,k)) / xlf,
2038  -qs_arr(i,j,k) / dtcld),
2039  Real(0.0));
2040  }
2041  if (qg_arr(i,j,k) > Real(0.0)) {
2042  pgeml_arr(i,j,k) = amrex::min(
2043  amrex::max(
2044  Real(cliq) * supcol
2045  * (paacw_arr(i,j,k) + pgacr_arr(i,j,k)) / xlf,
2046  -qg_arr(i,j,k) / dtcld),
2047  Real(0.0));
2048  }
2049  }
2050 
2051  if (supcol > Real(0.0)) {
2052  if (qi_arr(i,j,k) > Real(0.0) && ifsat != 1) {
2053  // WSM6-CPP TAG: PIDEP
2054  // legacy_group: G13h
2055  // process: Ice deposition/sublimation
2056  // compare_vars: pidep, qi, qv, qsat, den, t
2057  pidep_arr(i,j,k) = Real(4.0) * diameter * xni
2058  * (rhi_arr(i,j,k) - Real(1.0))
2059  / workdiffi_arr(i,j,k);
2060  Real supice = satdt - prevp_arr(i,j,k);
2061  if (pidep_arr(i,j,k) < Real(0.0)) {
2062  pidep_arr(i,j,k) = amrex::max(
2063  amrex::max(pidep_arr(i,j,k), satdt / Real(2.0)),
2064  supice);
2065  pidep_arr(i,j,k) = amrex::max(
2066  pidep_arr(i,j,k), -qi_arr(i,j,k) / dtcld);
2067  } else {
2068  pidep_arr(i,j,k) = amrex::min(
2069  amrex::min(pidep_arr(i,j,k), satdt / Real(2.0)),
2070  supice);
2071  }
2072  if (std::abs(prevp_arr(i,j,k) + pidep_arr(i,j,k))
2073  >= std::abs(satdt)) {
2074  ifsat = 1;
2075  }
2076  }
2077 
2078  if (qs_arr(i,j,k) > Real(0.0) && ifsat != 1) {
2079  const Real coeres = rslope2_s_arr(i,j,k)
2080  * std::sqrt(rslope_s_arr(i,j,k)
2081  * rslopeb_s_arr(i,j,k));
2082  psdep_arr(i,j,k) = (rhi_arr(i,j,k) - Real(1.0))
2083  * n0sfac_arr(i,j,k)
2084  * (Real(precs1) * rslope2_s_arr(i,j,k)
2085  + Real(precs2) * work2_arr(i,j,k) * coeres)
2086  / workdiffi_arr(i,j,k);
2087  Real supice = satdt - prevp_arr(i,j,k) - pidep_arr(i,j,k);
2088  if (psdep_arr(i,j,k) < Real(0.0)) {
2089  psdep_arr(i,j,k) = amrex::max(
2090  psdep_arr(i,j,k), -qs_arr(i,j,k) / dtcld);
2091  psdep_arr(i,j,k) = amrex::max(
2092  amrex::max(psdep_arr(i,j,k), satdt / Real(2.0)),
2093  supice);
2094  } else {
2095  psdep_arr(i,j,k) = amrex::min(
2096  amrex::min(psdep_arr(i,j,k), satdt / Real(2.0)),
2097  supice);
2098  }
2099  if (std::abs(prevp_arr(i,j,k) + pidep_arr(i,j,k)
2100  + psdep_arr(i,j,k))
2101  >= std::abs(satdt)) {
2102  ifsat = 1;
2103  }
2104  }
2105 
2106  if (qg_arr(i,j,k) > Real(0.0) && ifsat != 1) {
2107  const Real coeres = rslope2_g_arr(i,j,k)
2108  * std::sqrt(rslope_g_arr(i,j,k)
2109  * rslopeb_g_arr(i,j,k));
2110  pgdep_arr(i,j,k) = (rhi_arr(i,j,k) - Real(1.0))
2111  * (Real(precg1) * rslope2_g_arr(i,j,k)
2112  + Real(precg2) * work2_arr(i,j,k) * coeres)
2113  / workdiffi_arr(i,j,k);
2114  Real supice = satdt - prevp_arr(i,j,k)
2115  - pidep_arr(i,j,k) - psdep_arr(i,j,k);
2116  if (pgdep_arr(i,j,k) < Real(0.0)) {
2117  pgdep_arr(i,j,k) = amrex::max(
2118  pgdep_arr(i,j,k), -qg_arr(i,j,k) / dtcld);
2119  pgdep_arr(i,j,k) = amrex::max(
2120  amrex::max(pgdep_arr(i,j,k), satdt / Real(2.0)),
2121  supice);
2122  } else {
2123  pgdep_arr(i,j,k) = amrex::min(
2124  amrex::min(pgdep_arr(i,j,k), satdt / Real(2.0)),
2125  supice);
2126  }
2127  if (std::abs(prevp_arr(i,j,k) + pidep_arr(i,j,k)
2128  + psdep_arr(i,j,k) + pgdep_arr(i,j,k))
2129  >= std::abs(satdt)) {
2130  ifsat = 1;
2131  }
2132  }
2133 
2134  if (supsat > Real(0.0) && ifsat != 1) {
2135  const Real supice = satdt - prevp_arr(i,j,k)
2136  - pidep_arr(i,j,k)
2137  - psdep_arr(i,j,k)
2138  - pgdep_arr(i,j,k);
2139  const Real xni0 = Real(1.0e3) * std::exp(Real(0.1) * supcol);
2140  const Real roqi0 = Real(4.92e-11) * std::pow(xni0, Real(1.33));
2141  pigen_arr(i,j,k) = amrex::max(
2142  Real(0.0),
2143  (roqi0 / den_arr(i,j,k)
2144  - amrex::max(qi_arr(i,j,k), Real(0.0))) / dtcld);
2145  pigen_arr(i,j,k) = amrex::min(
2146  amrex::min(pigen_arr(i,j,k), satdt), supice);
2147  }
2148 
2149  if (qi_arr(i,j,k) > Real(0.0)) {
2150  const Real qimax = Real(roqimax) / den_arr(i,j,k);
2151  // WSM6-CPP TAG: PSAUT
2152  // legacy_group: G13i
2153  // process: Autoconversion to snow
2154  // compare_vars: psaut, qi, qs, den
2155  psaut_arr(i,j,k) = amrex::max(
2156  Real(0.0), (qi_arr(i,j,k) - qimax) / dtcld);
2157  }
2158 
2159  if (qs_arr(i,j,k) > Real(0.0)) {
2160  const Real alpha2 = Real(1.0e-3)
2161  * std::exp(Real(0.09) * (-supcol));
2162  pgaut_arr(i,j,k) = amrex::min(
2163  amrex::max(
2164  Real(0.0), alpha2 * (qs_arr(i,j,k) - Real(qs0))),
2165  qs_arr(i,j,k) / dtcld);
2166  }
2167  }
2168 
2169  if (supcol < Real(0.0)) {
2170  if (qs_arr(i,j,k) > Real(0.0)
2171  && rhw_arr(i,j,k) < Real(1.0)) {
2172  const Real coeres = rslope2_s_arr(i,j,k)
2173  * std::sqrt(rslope_s_arr(i,j,k)
2174  * rslopeb_s_arr(i,j,k));
2175  // WSM6-CPP TAG: PSEVP
2176  // legacy_group: G13j
2177  // process: Graupel evaporation/sublimation
2178  // compare_vars: psevp, qg, qv, qsat, den
2179  psevp_arr(i,j,k) = (rhw_arr(i,j,k) - Real(1.0))
2180  * n0sfac_arr(i,j,k)
2181  * (Real(precs1) * rslope2_s_arr(i,j,k)
2182  + Real(precs2) * work2_arr(i,j,k) * coeres)
2183  / workdiffw_arr(i,j,k);
2184  psevp_arr(i,j,k) = amrex::min(
2185  amrex::max(psevp_arr(i,j,k),
2186  -qs_arr(i,j,k) / dtcld),
2187  Real(0.0));
2188  }
2189 
2190  if (qg_arr(i,j,k) > Real(0.0)
2191  && rhw_arr(i,j,k) < Real(1.0)) {
2192  const Real coeres = rslope2_g_arr(i,j,k)
2193  * std::sqrt(rslope_g_arr(i,j,k)
2194  * rslopeb_g_arr(i,j,k));
2195  pgevp_arr(i,j,k) = (rhw_arr(i,j,k) - Real(1.0))
2196  * (Real(precg1) * rslope2_g_arr(i,j,k)
2197  + Real(precg2) * work2_arr(i,j,k) * coeres)
2198  / workdiffw_arr(i,j,k);
2199  pgevp_arr(i,j,k) = amrex::min(
2200  amrex::max(pgevp_arr(i,j,k),
2201  -qg_arr(i,j,k) / dtcld),
2202  Real(0.0));
2203  }
2204  }
2205  });
2206  // G14: mass conservation check and state update [lines 1200-1388]
2207  // WSM6-CPP TAG: UPDATE
2208  // legacy_group: G14
2209  // process: Mass conservation and state update
2210  // compare_vars: t, q, qci, qrs, qv
2211  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
2212  const Real qmin_l = Real(qmin);
2213  const Real qcrmin_l= Real(qcrmin);
2214  const Real t0c_l = Real(t0c);
2215 
2216  const Real delta2 =
2217  (qr_arr(i,j,k) < Real(1.0e-4) && qs_arr(i,j,k) < Real(1.0e-4))
2218  ? Real(1.0) : Real(0.0);
2219  const Real delta3 =
2220  (qr_arr(i,j,k) < Real(1.0e-4)) ? Real(1.0) : Real(0.0);
2221 
2222  if (t_arr(i,j,k) <= t0c_l) {
2223  Real value, source, factor, xlf, xlwork2;
2224 
2225  value = amrex::max(qmin_l, qc_arr(i,j,k));
2226  source = (praut_arr(i,j,k) + pracw_arr(i,j,k)
2227  + paacw_arr(i,j,k) + paacw_arr(i,j,k)) * dtcld;
2228 // + paacw_arr(i,j,k)) * dtcld;
2229  if (source > value) {
2230  factor = value / source;
2231  praut_arr(i,j,k) = praut_arr(i,j,k) * factor;
2232  pracw_arr(i,j,k) = pracw_arr(i,j,k) * factor;
2233  paacw_arr(i,j,k) = paacw_arr(i,j,k) * factor;
2234  }
2235 
2236  value = amrex::max(qmin_l, qi_arr(i,j,k));
2237  source = (psaut_arr(i,j,k) - pigen_arr(i,j,k)
2238  - pidep_arr(i,j,k) + praci_arr(i,j,k)
2239  + psaci_arr(i,j,k) + pgaci_arr(i,j,k)) * dtcld;
2240  if (source > value) {
2241  factor = value / source;
2242  psaut_arr(i,j,k) = psaut_arr(i,j,k) * factor;
2243  pigen_arr(i,j,k) = pigen_arr(i,j,k) * factor;
2244  pidep_arr(i,j,k) = pidep_arr(i,j,k) * factor;
2245  praci_arr(i,j,k) = praci_arr(i,j,k) * factor;
2246  psaci_arr(i,j,k) = psaci_arr(i,j,k) * factor;
2247  pgaci_arr(i,j,k) = pgaci_arr(i,j,k) * factor;
2248  }
2249 
2250  value = amrex::max(qmin_l, qr_arr(i,j,k));
2251  source = (-praut_arr(i,j,k) - prevp_arr(i,j,k)
2252  - pracw_arr(i,j,k) + piacr_arr(i,j,k)
2253  + psacr_arr(i,j,k) + pgacr_arr(i,j,k)) * dtcld;
2254  if (source > value) {
2255  factor = value / source;
2256  praut_arr(i,j,k) = praut_arr(i,j,k) * factor;
2257  prevp_arr(i,j,k) = prevp_arr(i,j,k) * factor;
2258  pracw_arr(i,j,k) = pracw_arr(i,j,k) * factor;
2259  piacr_arr(i,j,k) = piacr_arr(i,j,k) * factor;
2260  psacr_arr(i,j,k) = psacr_arr(i,j,k) * factor;
2261  pgacr_arr(i,j,k) = pgacr_arr(i,j,k) * factor;
2262  }
2263 
2264  value = amrex::max(qmin_l, qs_arr(i,j,k));
2265  source = -(psdep_arr(i,j,k) + psaut_arr(i,j,k)
2266  - pgaut_arr(i,j,k) + paacw_arr(i,j,k)
2267  + piacr_arr(i,j,k) * delta3
2268  + praci_arr(i,j,k) * delta3
2269  - pracs_arr(i,j,k) * (Real(1.0) - delta2)
2270  + psacr_arr(i,j,k) * delta2
2271  + psaci_arr(i,j,k) - pgacs_arr(i,j,k)) * dtcld;
2272  if (source > value) {
2273  factor = value / source;
2274  psdep_arr(i,j,k) = psdep_arr(i,j,k) * factor;
2275  psaut_arr(i,j,k) = psaut_arr(i,j,k) * factor;
2276  pgaut_arr(i,j,k) = pgaut_arr(i,j,k) * factor;
2277  paacw_arr(i,j,k) = paacw_arr(i,j,k) * factor;
2278  piacr_arr(i,j,k) = piacr_arr(i,j,k) * factor;
2279  praci_arr(i,j,k) = praci_arr(i,j,k) * factor;
2280  psaci_arr(i,j,k) = psaci_arr(i,j,k) * factor;
2281  pracs_arr(i,j,k) = pracs_arr(i,j,k) * factor;
2282  psacr_arr(i,j,k) = psacr_arr(i,j,k) * factor;
2283  pgacs_arr(i,j,k) = pgacs_arr(i,j,k) * factor;
2284  }
2285 
2286  value = amrex::max(qmin_l, qg_arr(i,j,k));
2287  source = -(pgdep_arr(i,j,k) + pgaut_arr(i,j,k)
2288  + piacr_arr(i,j,k) * (Real(1.0) - delta3)
2289  + praci_arr(i,j,k) * (Real(1.0) - delta3)
2290  + psacr_arr(i,j,k) * (Real(1.0) - delta2)
2291  + pracs_arr(i,j,k) * (Real(1.0) - delta2)
2292  + pgaci_arr(i,j,k) + paacw_arr(i,j,k)
2293  + pgacr_arr(i,j,k) + pgacs_arr(i,j,k)) * dtcld;
2294  if (source > value) {
2295  factor = value / source;
2296  pgdep_arr(i,j,k) = pgdep_arr(i,j,k) * factor;
2297  pgaut_arr(i,j,k) = pgaut_arr(i,j,k) * factor;
2298  piacr_arr(i,j,k) = piacr_arr(i,j,k) * factor;
2299  praci_arr(i,j,k) = praci_arr(i,j,k) * factor;
2300  psacr_arr(i,j,k) = psacr_arr(i,j,k) * factor;
2301  pracs_arr(i,j,k) = pracs_arr(i,j,k) * factor;
2302  paacw_arr(i,j,k) = paacw_arr(i,j,k) * factor;
2303  pgaci_arr(i,j,k) = pgaci_arr(i,j,k) * factor;
2304  pgacr_arr(i,j,k) = pgacr_arr(i,j,k) * factor;
2305  pgacs_arr(i,j,k) = pgacs_arr(i,j,k) * factor;
2306  }
2307 
2308  work2_arr(i,j,k) = -(prevp_arr(i,j,k) + psdep_arr(i,j,k)
2309  + pgdep_arr(i,j,k) + pigen_arr(i,j,k)
2310  + pidep_arr(i,j,k));
2311  qv_arr(i,j,k) = qv_arr(i,j,k) + work2_arr(i,j,k) * dtcld;
2312 // + paacw_arr(i,j,k))
2313  qc_arr(i,j,k) = amrex::max(
2314  qc_arr(i,j,k) - (praut_arr(i,j,k) + pracw_arr(i,j,k)
2315  + paacw_arr(i,j,k) + paacw_arr(i,j,k))
2316  * dtcld,
2317  Real(0.0));
2318  qr_arr(i,j,k) = amrex::max(
2319  qr_arr(i,j,k) + (praut_arr(i,j,k) + pracw_arr(i,j,k)
2320  + prevp_arr(i,j,k) - piacr_arr(i,j,k)
2321  - pgacr_arr(i,j,k) - psacr_arr(i,j,k))
2322  * dtcld,
2323  Real(0.0));
2324  qi_arr(i,j,k) = amrex::max(
2325  qi_arr(i,j,k) - (psaut_arr(i,j,k) + praci_arr(i,j,k)
2326  + psaci_arr(i,j,k) + pgaci_arr(i,j,k)
2327  - pigen_arr(i,j,k) - pidep_arr(i,j,k))
2328  * dtcld,
2329  Real(0.0));
2330  qs_arr(i,j,k) = amrex::max(
2331  qs_arr(i,j,k) + (psdep_arr(i,j,k) + psaut_arr(i,j,k)
2332  + paacw_arr(i,j,k) - pgaut_arr(i,j,k)
2333  + piacr_arr(i,j,k) * delta3
2334  + praci_arr(i,j,k) * delta3
2335  + psaci_arr(i,j,k) - pgacs_arr(i,j,k)
2336  - pracs_arr(i,j,k) * (Real(1.0) - delta2)
2337  + psacr_arr(i,j,k) * delta2) * dtcld,
2338  Real(0.0));
2339  qg_arr(i,j,k) = amrex::max(
2340  qg_arr(i,j,k) + (pgdep_arr(i,j,k) + pgaut_arr(i,j,k)
2341  + piacr_arr(i,j,k) * (Real(1.0) - delta3)
2342  + praci_arr(i,j,k) * (Real(1.0) - delta3)
2343  + psacr_arr(i,j,k) * (Real(1.0) - delta2)
2344  + pracs_arr(i,j,k) * (Real(1.0) - delta2)
2345  + pgaci_arr(i,j,k) + paacw_arr(i,j,k)
2346  + pgacr_arr(i,j,k) + pgacs_arr(i,j,k))
2347  * dtcld,
2348  Real(0.0));
2349  xlf = Real(xls) - xl_arr(i,j,k);
2350  xlwork2 = -Real(xls) * (psdep_arr(i,j,k) + pgdep_arr(i,j,k)
2351  + pidep_arr(i,j,k) + pigen_arr(i,j,k))
2352  - xl_arr(i,j,k) * prevp_arr(i,j,k)
2353  - xlf * (piacr_arr(i,j,k) + paacw_arr(i,j,k)
2354  + paacw_arr(i,j,k) + pgacr_arr(i,j,k)
2355  + psacr_arr(i,j,k));
2356  t_arr(i,j,k) = t_arr(i,j,k) - xlwork2 / cpm_arr(i,j,k) * dtcld;
2357  } else {
2358  Real value, source, factor, xlf, xlwork2;
2359 
2360  value = amrex::max(qmin_l, qc_arr(i,j,k));
2361  source = (praut_arr(i,j,k) + pracw_arr(i,j,k)
2362  + paacw_arr(i,j,k) + paacw_arr(i,j,k)) * dtcld;
2363 // + paacw_arr(i,j,k)) * dtcld;
2364  if (source > value) {
2365  factor = value / source;
2366  praut_arr(i,j,k) = praut_arr(i,j,k) * factor;
2367  pracw_arr(i,j,k) = pracw_arr(i,j,k) * factor;
2368  paacw_arr(i,j,k) = paacw_arr(i,j,k) * factor;
2369  }
2370 
2371  value = amrex::max(qmin_l, qr_arr(i,j,k));
2372  source = (-paacw_arr(i,j,k) - praut_arr(i,j,k)
2373  + pseml_arr(i,j,k) + pgeml_arr(i,j,k)
2374  - pracw_arr(i,j,k) - paacw_arr(i,j,k)
2375  - prevp_arr(i,j,k)) * dtcld;
2376  if (source > value) {
2377  factor = value / source;
2378  praut_arr(i,j,k) = praut_arr(i,j,k) * factor;
2379  prevp_arr(i,j,k) = prevp_arr(i,j,k) * factor;
2380  pracw_arr(i,j,k) = pracw_arr(i,j,k) * factor;
2381  paacw_arr(i,j,k) = paacw_arr(i,j,k) * factor;
2382  pseml_arr(i,j,k) = pseml_arr(i,j,k) * factor;
2383  pgeml_arr(i,j,k) = pgeml_arr(i,j,k) * factor;
2384  }
2385 
2386  value = amrex::max(qcrmin_l, qs_arr(i,j,k));
2387  source = (pgacs_arr(i,j,k) - pseml_arr(i,j,k)
2388  - psevp_arr(i,j,k)) * dtcld;
2389  if (source > value) {
2390  factor = value / source;
2391  pgacs_arr(i,j,k) = pgacs_arr(i,j,k) * factor;
2392  psevp_arr(i,j,k) = psevp_arr(i,j,k) * factor;
2393  pseml_arr(i,j,k) = pseml_arr(i,j,k) * factor;
2394  }
2395 
2396  value = amrex::max(qcrmin_l, qg_arr(i,j,k));
2397  source = -(pgacs_arr(i,j,k) + pgevp_arr(i,j,k)
2398  + pgeml_arr(i,j,k)) * dtcld;
2399  if (source > value) {
2400  factor = value / source;
2401  pgacs_arr(i,j,k) = pgacs_arr(i,j,k) * factor;
2402  pgevp_arr(i,j,k) = pgevp_arr(i,j,k) * factor;
2403  pgeml_arr(i,j,k) = pgeml_arr(i,j,k) * factor;
2404  }
2405 
2406  work2_arr(i,j,k) = -(prevp_arr(i,j,k) + psevp_arr(i,j,k)
2407  + pgevp_arr(i,j,k));
2408  qv_arr(i,j,k) = qv_arr(i,j,k) + work2_arr(i,j,k) * dtcld;
2409 // + paacw_arr(i,j,k))
2410  qc_arr(i,j,k) = amrex::max(
2411  qc_arr(i,j,k) - (praut_arr(i,j,k) + pracw_arr(i,j,k)
2412  + paacw_arr(i,j,k) + paacw_arr(i,j,k))
2413  * dtcld,
2414  Real(0.0));
2415  qr_arr(i,j,k) = amrex::max(
2416  qr_arr(i,j,k) + (praut_arr(i,j,k) + pracw_arr(i,j,k)
2417  + prevp_arr(i,j,k) + paacw_arr(i,j,k)
2418  + paacw_arr(i,j,k) - pseml_arr(i,j,k)
2419  - pgeml_arr(i,j,k)) * dtcld,
2420  Real(0.0));
2421  qs_arr(i,j,k) = amrex::max(
2422  qs_arr(i,j,k) + (psevp_arr(i,j,k) - pgacs_arr(i,j,k)
2423  + pseml_arr(i,j,k)) * dtcld,
2424  Real(0.0));
2425  qg_arr(i,j,k) = amrex::max(
2426  qg_arr(i,j,k) + (pgacs_arr(i,j,k) + pgevp_arr(i,j,k)
2427  + pgeml_arr(i,j,k)) * dtcld,
2428  Real(0.0));
2429  xlf = Real(xls) - xl_arr(i,j,k);
2430  xlwork2 = -xl_arr(i,j,k) * (prevp_arr(i,j,k)
2431  + psevp_arr(i,j,k)
2432  + pgevp_arr(i,j,k))
2433  - xlf * (pseml_arr(i,j,k) + pgeml_arr(i,j,k));
2434  t_arr(i,j,k) = t_arr(i,j,k) - xlwork2 / cpm_arr(i,j,k) * dtcld;
2435  }
2436  });
2437  // G15: second qsat computation [lines 1390-1420]
2438  // WSM6-CPP TAG: QSAT2
2439  // legacy_group: G15
2440  // process: Second saturation mixing ratio computation
2441  // compare_vars: qs, qvs, den, denfac, t, p
2442  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
2443  const Real ttp = Real(t0c) + Real(0.01);
2444  const Real dldt = Real(cpv) - Real(cliq);
2445  const Real xa = -dldt / Real(rv);
2446  const Real xb = xa + Real(xlv0) / (Real(rv) * ttp);
2447 
2448  Real tr = ttp / t_arr(i,j,k);
2449  Real qsw = Real(psat) * std::exp(std::log(tr) * xa)
2450  * std::exp(xb * (Real(1.0) - tr));
2451  qsw = amrex::min(qsw, Real(0.99) * p_arr(i,j,k));
2452  qsatw_arr(i,j,k) = Real(ep2) * qsw / (p_arr(i,j,k) - qsw);
2453  qsatw_arr(i,j,k) = amrex::max(qsatw_arr(i,j,k), Real(qmin));
2454  });
2455  // G16: pcond condensational/evaporational update [lines 1427-1437]
2456  // WSM6-CPP TAG: PCOND
2457  // legacy_group: G16
2458  // process: Condensation/evaporation update
2459  // compare_vars: pcond, t, qv, qc, qsat
2460  if (m_do_cond) {
2461  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
2462  const Real workcond = wsm6_conden(
2463  t_arr(i,j,k), qv_arr(i,j,k), qsatw_arr(i,j,k),
2464  xl_arr(i,j,k), cpm_arr(i,j,k), Real(qmin), Real(rv));
2465  const Real work2loc = qc_arr(i,j,k) + workcond;
2466  static_cast<void>(work2loc);
2467  pcond_arr(i,j,k) = amrex::min(
2468  amrex::max(workcond / dtcld, Real(0.0)),
2469  amrex::max(qv_arr(i,j,k), Real(0.0)) / dtcld);
2470  if (qc_arr(i,j,k) > Real(0.0) && workcond < Real(0.0)) {
2471  pcond_arr(i,j,k) = amrex::max(workcond, -qc_arr(i,j,k)) / dtcld;
2472  }
2473  qv_arr(i,j,k) = qv_arr(i,j,k) - pcond_arr(i,j,k) * dtcld;
2474  qc_arr(i,j,k) = amrex::max(
2475  qc_arr(i,j,k) + pcond_arr(i,j,k) * dtcld,
2476  Real(0.0));
2477  t_arr(i,j,k) = t_arr(i,j,k)
2478  + pcond_arr(i,j,k) * xl_arr(i,j,k)
2479  / cpm_arr(i,j,k) * dtcld;
2480  });
2481  }
2482  // G17: padding for small values [lines 1444-1449]
2483  // WSM6-CPP TAG: CLIP
2484  // legacy_group: G17
2485  // process: Padding/clipping for small values
2486  // compare_vars: qv, qc, qr, qi, qs, qg
2487  ParallelFor(box, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
2488  if (qc_arr(i,j,k) <= Real(qmin)) qc_arr(i,j,k) = Real(0.0);
2489  if (qi_arr(i,j,k) <= Real(qmin)) qi_arr(i,j,k) = Real(0.0);
2490  });
2491 
2492  }
2493 #ifdef ERF_USE_WSM6_FORT
2494  }
2495 #endif
2496  ParallelFor(box2d, [=] AMREX_GPU_DEVICE (int i, int j, int) {
2497  rain_arr(i,j,klo) = rainacc_arr(i,j,0);
2498  snow_arr(i,j,klo) = snowacc_arr(i,j,0);
2499  graup_arr(i,j,klo) = graupacc_arr(i,j,0);
2500  });
2501  }
2502 }
wsm6_xka
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real wsm6_xka(Real x, Real y)
Definition: ERF_AdvanceWSM6.cpp:44
wsm6_conden
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real wsm6_conden(Real a, Real b, Real c, Real d, Real e, Real qmin_arg, Real rv_arg)
Definition: ERF_AdvanceWSM6.cpp:64
wsm6_slope_snow_cell
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void wsm6_slope_snow_cell(Real qs, Real den, Real denfac, Real t, Real pidn0s_arg, Real alpha_arg, Real n0smax_arg, Real n0s_arg, Real t0c_arg, Real qcrmin_arg, Real rslopesmax_arg, Real rslopesbmax_arg, Real rslopes2max_arg, Real rslopes3max_arg, Real bvts_arg, Real pvts_arg, Real &rslope, Real &rslopeb, Real &rslope2, Real &rslope3, Real &vt, Real &n0sfac)
Definition: ERF_AdvanceWSM6.cpp:170
wsm6_slope_graup_cell
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void wsm6_slope_graup_cell(Real qg, Real den, Real denfac, Real pidn0g_arg, Real qcrmin_arg, Real rslopegmax_arg, Real rslopegbmax_arg, Real rslopeg2max_arg, Real rslopeg3max_arg, Real bvtg_arg, Real pvtg_arg, Real &rslope, Real &rslopeb, Real &rslope2, Real &rslope3, Real &vt)
Definition: ERF_AdvanceWSM6.cpp:200
wsm6_slope_rain_cell
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void wsm6_slope_rain_cell(Real qr, Real den, Real denfac, Real pidn0r_arg, Real qcrmin_arg, Real rslopermax_arg, Real rsloperbmax_arg, Real rsloper2max_arg, Real rsloper3max_arg, Real bvtr_arg, Real pvtr_arg, Real &rslope, Real &rslopeb, Real &rslope2, Real &rslope3, Real &vt)
Definition: ERF_AdvanceWSM6.cpp:145
wsm6_nislfv_rain_plm_scratch
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void wsm6_nislfv_rain_plm_scratch(int km, int ww_comp, int rq_comp, Real *precip, Real dt, int iter, Array4< Real > const &sed_cell, Array4< Real > const &sed_node, int i_s, int j_s, int klo_s)
Definition: ERF_AdvanceWSM6.cpp:224
wsm6_xlcal
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real wsm6_xlcal(Real x, Real xlv0_arg, Real xlv1_arg, Real t0c_arg)
Definition: ERF_AdvanceWSM6.cpp:29
wsm6_venfac
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real wsm6_venfac(Real a, Real b, Real c, Real den0_arg)
Definition: ERF_AdvanceWSM6.cpp:56
wsm6_cpmcal
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real wsm6_cpmcal(Real x, Real qmin_arg, Real cpd_arg, Real cpv_arg)
Definition: ERF_AdvanceWSM6.cpp:23
wsm6_nislfv_rain_plm6_scratch
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void wsm6_nislfv_rain_plm6_scratch(int km, int ww_comp, int rq_comp, int rq2_comp, Real *precip1, Real *precip2, Real dt, int iter, Array4< Real > const &sed_cell, Array4< Real > const &sed_node, int i_s, int j_s, int klo_s)
Definition: ERF_AdvanceWSM6.cpp:498
wsm6_diffac
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real wsm6_diffac(Real a, Real b, Real c, Real d, Real e, Real rv_arg)
Definition: ERF_AdvanceWSM6.cpp:49
R_v
constexpr amrex::Real R_v
Definition: ERF_Constants.H:30
lat_vap
constexpr amrex::Real lat_vap
Definition: ERF_Constants.H:104
Cp_d
constexpr amrex::Real Cp_d
Definition: ERF_Constants.H:31
lsub
constexpr amrex::Real lsub
Definition: ERF_Constants.H:87
rhoh2o
constexpr amrex::Real rhoh2o
Definition: ERF_Constants.H:113
one
constexpr amrex::Real one
Definition: ERF_Constants.H:7
lat_ice
constexpr amrex::Real lat_ice
Definition: ERF_Constants.H:105
rhos
constexpr amrex::Real rhos
Definition: ERF_Constants.H:48
CONST_GRAV
constexpr amrex::Real CONST_GRAV
Definition: ERF_Constants.H:40
Cp_l
constexpr amrex::Real Cp_l
Definition: ERF_Constants.H:33
Cp_v
constexpr amrex::Real Cp_v
Definition: ERF_Constants.H:32
R_d
constexpr amrex::Real R_d
Definition: ERF_Constants.H:29
value
amrex::Real value
Definition: ERF_HurricaneDiagnostics.H:20
pp
ParmParse pp("prob")
khi
const int khi
Definition: ERF_InitCustomPert_Bubble.H:21
qv_arr
auto qv_arr
Definition: ERF_InitCustomPert_MultiSpeciesBubble.H:210
ParallelFor
ParallelFor(grown_box, [=] AMREX_GPU_DEVICE(int i, int j, int k) { qrcuten_arr(i, j, k)=Real(0);qscuten_arr(i, j, k)=Real(0);qicuten_arr(i, j, k)=Real(0);})
Real
amrex::Real Real
Definition: ERF_ShocInterface.H:19
TileNoZ
AMREX_FORCE_INLINE amrex::IntVect TileNoZ()
Definition: ERF_TileNoZ.H:11
mp_wsm6_init_c
void mp_wsm6_init_c(double den0, double denr, double dens, double cl, double cpv, int hail_opt)
mp_wsm6_run_c
void mp_wsm6_run_c(double *t, double *qv, double *qc, double *qi, double *qr, double *qs, double *qg, double *den, double *p, double *delz, double delt, double g, double cpd, double cpv, double rd, double rv, double t0c, double ep1, double ep2, double qmin, double xls, double xlv0, double xlf0, double den0, double denr, double cliq, double cice, double psat, double *rain, double *rainncv, double *sr, double *snow, double *snowncv, double *graupel, double *graupelncv, int ims, int ime, int jms, int jme, int kms, int kme, int its, int ite, int jts, int jte, int kts, int kte, int microphysics_debug, int diag_i_dbg, int diag_j_dbg)
WSM6::m_roqimax
amrex::Real m_roqimax
Definition: ERF_WSM6.H:151
WSM6::m_z_phys_nd
amrex::MultiFab * m_z_phys_nd
Definition: ERF_WSM6.H:133
WSM6::m_rslopesbmax
amrex::Real m_rslopesbmax
Definition: ERF_WSM6.H:160
WSM6::avtr
static constexpr amrex::Real avtr
Definition: ERF_WSM6.H:53
WSM6::dicon
static constexpr amrex::Real dicon
Definition: ERF_WSM6.H:63
WSM6::m_rslopes2max
amrex::Real m_rslopes2max
Definition: ERF_WSM6.H:161
WSM6::m_rsloper2max
amrex::Real m_rsloper2max
Definition: ERF_WSM6.H:161
WSM6::qcrmin
static constexpr amrex::Real qcrmin
Definition: ERF_WSM6.H:67
WSM6::m_rslopermax
amrex::Real m_rslopermax
Definition: ERF_WSM6.H:159
WSM6::m_pacrc
amrex::Real m_pacrc
Definition: ERF_WSM6.H:155
WSM6::pfrz1
static constexpr amrex::Real pfrz1
Definition: ERF_WSM6.H:65
WSM6::m_rsloperbmax
amrex::Real m_rsloperbmax
Definition: ERF_WSM6.H:160
WSM6::m_pi_wsm6
amrex::Real m_pi_wsm6
Definition: ERF_WSM6.H:146
WSM6::n0smax
static constexpr amrex::Real n0smax
Definition: ERF_WSM6.H:71
WSM6::m_do_cond
bool m_do_cond
Definition: ERF_WSM6.H:131
WSM6::m_pidn0s
amrex::Real m_pidn0s
Definition: ERF_WSM6.H:155
WSM6::dimax
static constexpr amrex::Real dimax
Definition: ERF_WSM6.H:64
WSM6::pfrz2
static constexpr amrex::Real pfrz2
Definition: ERF_WSM6.H:66
WSM6::alpha_wsm6
static constexpr amrex::Real alpha_wsm6
Definition: ERF_WSM6.H:73
WSM6::m_pvtg
amrex::Real m_pvtg
Definition: ERF_WSM6.H:158
WSM6::m_qck1
amrex::Real m_qck1
Definition: ERF_WSM6.H:147
WSM6::m_precg2
amrex::Real m_precg2
Definition: ERF_WSM6.H:158
WSM6::m_precg1
amrex::Real m_precg1
Definition: ERF_WSM6.H:158
WSM6::m_rsloper3max
amrex::Real m_rsloper3max
Definition: ERF_WSM6.H:162
WSM6::n0s
static constexpr amrex::Real n0s
Definition: ERF_WSM6.H:72
WSM6::m_precr2
amrex::Real m_precr2
Definition: ERF_WSM6.H:151
WSM6::n0r
static constexpr amrex::Real n0r
Definition: ERF_WSM6.H:52
WSM6::m_n0g
amrex::Real m_n0g
Definition: ERF_WSM6.H:141
WSM6::m_pidn0g
amrex::Real m_pidn0g
Definition: ERF_WSM6.H:158
WSM6::dens_snow
static constexpr amrex::Real dens_snow
Definition: ERF_WSM6.H:69
WSM6::dtcldcr
static constexpr amrex::Real dtcldcr
Definition: ERF_WSM6.H:51
WSM6::m_rslopegmax
amrex::Real m_rslopegmax
Definition: ERF_WSM6.H:159
WSM6::bvts
static constexpr amrex::Real bvts
Definition: ERF_WSM6.H:60
WSM6::m_g6pbr
amrex::Real m_g6pbr
Definition: ERF_WSM6.H:149
WSM6::m_pacrg
amrex::Real m_pacrg
Definition: ERF_WSM6.H:158
WSM6::m_precs2
amrex::Real m_precs2
Definition: ERF_WSM6.H:154
WSM6::m_rslopeg3max
amrex::Real m_rslopeg3max
Definition: ERF_WSM6.H:162
WSM6::m_pvts
amrex::Real m_pvts
Definition: ERF_WSM6.H:154
WSM6::m_rslopegbmax
amrex::Real m_rslopegbmax
Definition: ERF_WSM6.H:160
WSM6::qs0
static constexpr amrex::Real qs0
Definition: ERF_WSM6.H:70
WSM6::m_geom
amrex::Geometry m_geom
Definition: ERF_WSM6.H:126
WSM6::m_precs1
amrex::Real m_precs1
Definition: ERF_WSM6.H:154
WSM6::dt
amrex::Real dt
Definition: ERF_WSM6.H:127
WSM6::xncr
static constexpr amrex::Real xncr
Definition: ERF_WSM6.H:57
WSM6::bvtr
static constexpr amrex::Real bvtr
Definition: ERF_WSM6.H:54
WSM6::m_rslopesmax
amrex::Real m_rslopesmax
Definition: ERF_WSM6.H:159
WSM6::m_qc0
amrex::Real m_qc0
Definition: ERF_WSM6.H:147
WSM6::mic_fab_vars
amrex::Array< FabPtr, MicVar_WSM6::NumVars > mic_fab_vars
Definition: ERF_WSM6.H:136
WSM6::m_bvtg
amrex::Real m_bvtg
Definition: ERF_WSM6.H:141
WSM6::m_pvtr
amrex::Real m_pvtr
Definition: ERF_WSM6.H:150
WSM6::m_pacrr
amrex::Real m_pacrr
Definition: ERF_WSM6.H:150
WSM6::m_precr1
amrex::Real m_precr1
Definition: ERF_WSM6.H:151
WSM6::m_pidn0r
amrex::Real m_pidn0r
Definition: ERF_WSM6.H:155
WSM6::m_rslopes3max
amrex::Real m_rslopes3max
Definition: ERF_WSM6.H:162
WSM6::m_xlv1
amrex::Real m_xlv1
Definition: ERF_WSM6.H:146
WSM6::m_rslopeg2max
amrex::Real m_rslopeg2max
Definition: ERF_WSM6.H:161
MORRInd::xlf
@ xlf
Definition: ERF_AdvanceMorrison.cpp:157
MicVar_WSM6::qr
@ qr
Definition: ERF_WSM6.H:27
MicVar_WSM6::qi
@ qi
Definition: ERF_WSM6.H:26
MicVar_WSM6::qs
@ qs
Definition: ERF_WSM6.H:28
MicVar_WSM6::qv
@ qv
Definition: ERF_WSM6.H:24
MicVar_WSM6::qc
@ qc
Definition: ERF_WSM6.H:25
MicVar_WSM6::rain_accum
@ rain_accum
Definition: ERF_WSM6.H:30
MicVar_WSM6::snow_accum
@ snow_accum
Definition: ERF_WSM6.H:31
MicVar_WSM6::rho
@ rho
Definition: ERF_WSM6.H:20
MicVar_WSM6::graup_accum
@ graup_accum
Definition: ERF_WSM6.H:32
MicVar_WSM6::qg
@ qg
Definition: ERF_WSM6.H:29
MicVar_WSM6::pres
@ pres
Definition: ERF_WSM6.H:23
MicVar_WSM6::tabs
@ tabs
Definition: ERF_WSM6.H:22
WSM6Ind::qsum
@ qsum
Definition: ERF_WSM6.H:219
WSM6Ind::xni
@ xni
Definition: ERF_WSM6.H:224
WSM6SedCellScratch::NumComps
@ NumComps
Definition: ERF_AdvanceWSM6.cpp:126
WSM6SedCellScratch::workr_col
@ workr_col
Definition: ERF_AdvanceWSM6.cpp:118
WSM6SedCellScratch::tmp
@ tmp
Definition: ERF_AdvanceWSM6.cpp:114
WSM6SedCellScratch::denqrs2_col
@ denqrs2_col
Definition: ERF_AdvanceWSM6.cpp:121
WSM6SedCellScratch::qsum_col
@ qsum_col
Definition: ERF_AdvanceWSM6.cpp:123
WSM6SedCellScratch::denqci_col
@ denqci_col
Definition: ERF_AdvanceWSM6.cpp:125
WSM6SedCellScratch::denqrs1_col
@ denqrs1_col
Definition: ERF_AdvanceWSM6.cpp:120
WSM6SedCellScratch::worka_col
@ worka_col
Definition: ERF_AdvanceWSM6.cpp:119
WSM6SedCellScratch::denqrs3_col
@ denqrs3_col
Definition: ERF_AdvanceWSM6.cpp:122
WSM6SedCellScratch::den
@ den
Definition: ERF_AdvanceWSM6.cpp:109
WSM6SedCellScratch::dz
@ dz
Definition: ERF_AdvanceWSM6.cpp:104
WSM6SedCellScratch::tk
@ tk
Definition: ERF_AdvanceWSM6.cpp:111
WSM6SedCellScratch::denfac
@ denfac
Definition: ERF_AdvanceWSM6.cpp:110
WSM6SedCellScratch::work1c_col
@ work1c_col
Definition: ERF_AdvanceWSM6.cpp:124
WSM6SedNodeScratch::NumComps
@ NumComps
Definition: ERF_AdvanceWSM6.cpp:140
module_model_constants::cice
real(c_double), parameter cice
Definition: ERF_module_model_constants.F90:30
module_model_constants::g
real(c_double), parameter g
Definition: ERF_module_model_constants.F90:19
module_model_constants::cliq
real(c_double), parameter cliq
Definition: ERF_module_model_constants.F90:29
module_model_constants::xlv0
real(c_double), parameter xlv0
Definition: ERF_module_model_constants.F90:51
module_model_constants::cpv
real(c_double), parameter cpv
Definition: ERF_module_model_constants.F90:26
module_model_constants::xls
real(c_double), parameter xls
Definition: ERF_module_model_constants.F90:56
module_model_constants::psat
real(c_double), parameter psat
Definition: ERF_module_model_constants.F90:31
module_mp_morr_two_moment::pi
real(c_double), parameter, private pi
Definition: ERF_module_mp_morr_two_moment.F90:100
mp_wsm6::precg2
real(kind=kind_phys), save precg2
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::precg1
real(kind=kind_phys), save precg1
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::roqimax
real(kind=kind_phys), save roqimax
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::dens
real(kind=kind_phys), parameter, private dens
Definition: ERF_module_mp_wsm6.F90:39
mp_wsm6::g6pbr
real(kind=kind_phys), save g6pbr
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::precr2
real(kind=kind_phys), save precr2
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::pacrr
real(kind=kind_phys), save pacrr
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::qck1
real(kind=kind_phys), save qck1
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::precr1
real(kind=kind_phys), save precr1
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::pacrg
real(kind=kind_phys), save pacrg
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::pvtg
real(kind=kind_phys), save pvtg
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::pvts
real(kind=kind_phys), save pvts
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::precs1
real(kind=kind_phys), save precs1
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::pvtr
real(kind=kind_phys), save pvtr
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::precs2
real(kind=kind_phys), save precs2
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::pacrc
real(kind=kind_phys), save pacrc
Definition: ERF_module_mp_wsm6.F90:46
mp_wsm6::qc0
real(kind=kind_phys), save qc0
Definition: ERF_module_mp_wsm6.F90:46
Here is the call graph for this function:

◆ Copy_Micro_to_State()

void WSM6::Copy_Micro_to_State ( amrex::MultiFab &  cons_in)
overridevirtual

Reimplemented from NullMoist.

8 {
9  for (MFIter mfi(cons, TilingIfNotGPU()); mfi.isValid(); ++mfi) {
10  const auto& box3d = mfi.tilebox();
11  auto states = cons.array(mfi);
12 
13  auto rho = mic_fab_vars[MicVar_WSM6::rho]->array(mfi);
14  auto theta = mic_fab_vars[MicVar_WSM6::theta]->array(mfi);
15  auto tabs = mic_fab_vars[MicVar_WSM6::tabs]->array(mfi);
16 
17  auto qv = mic_fab_vars[MicVar_WSM6::qv]->array(mfi);
18  auto qc = mic_fab_vars[MicVar_WSM6::qc]->array(mfi);
19  auto qi = mic_fab_vars[MicVar_WSM6::qi]->array(mfi);
20  auto qr = mic_fab_vars[MicVar_WSM6::qr]->array(mfi);
21  auto qs = mic_fab_vars[MicVar_WSM6::qs]->array(mfi);
22  auto qg = mic_fab_vars[MicVar_WSM6::qg]->array(mfi);
23 
24  ParallelFor(box3d, [=] AMREX_GPU_DEVICE(int i, int j, int k) {
25  theta(i,j,k) = getThgivenRandT(rho(i,j,k), tabs(i,j,k), R_d / Cp_d, qv(i,j,k));
26  states(i,j,k,RhoTheta_comp) = rho(i,j,k) * theta(i,j,k);
27  states(i,j,k,RhoQ1_comp) = rho(i,j,k) * amrex::max(Real(0), qv(i,j,k));
28  states(i,j,k,RhoQ2_comp) = rho(i,j,k) * amrex::max(Real(0), qc(i,j,k));
29  states(i,j,k,RhoQ3_comp) = rho(i,j,k) * amrex::max(Real(0), qi(i,j,k));
30  states(i,j,k,RhoQ4_comp) = rho(i,j,k) * amrex::max(Real(0), qr(i,j,k));
31  states(i,j,k,RhoQ5_comp) = rho(i,j,k) * amrex::max(Real(0), qs(i,j,k));
32  states(i,j,k,RhoQ6_comp) = rho(i,j,k) * amrex::max(Real(0), qg(i,j,k));
33  });
34  }
35 
36  cons.FillBoundary(m_geom.periodicity());
37 }
getThgivenRandT
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getThgivenRandT(const amrex::Real rho, const amrex::Real T, const amrex::Real rdOcp, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:64
RhoQ4_comp
#define RhoQ4_comp
Definition: ERF_IndexDefines.H:45
RhoTheta_comp
#define RhoTheta_comp
Definition: ERF_IndexDefines.H:37
RhoQ2_comp
#define RhoQ2_comp
Definition: ERF_IndexDefines.H:43
RhoQ3_comp
#define RhoQ3_comp
Definition: ERF_IndexDefines.H:44
RhoQ1_comp
#define RhoQ1_comp
Definition: ERF_IndexDefines.H:42
RhoQ6_comp
#define RhoQ6_comp
Definition: ERF_IndexDefines.H:47
RhoQ5_comp
#define RhoQ5_comp
Definition: ERF_IndexDefines.H:46
rho
rho
Definition: ERF_InitCustomPert_Bubble.H:106
LsmVar_MM5::theta
@ theta
Definition: ERF_MM5.H:20
MicVar_Kess::tabs
@ tabs
Definition: ERF_Kessler.H:25
MicVar_Kess::qv
@ qv
Definition: ERF_Kessler.H:29
MicVar_SatAdj::qc
@ qc
Definition: ERF_SatAdj.H:40
MicVar_WSM6::theta
@ theta
Definition: ERF_WSM6.H:21
Vars::cons
@ cons
Definition: ERF_IndexDefines.H:174
WSM6SedCellScratch::qr
@ qr
Definition: ERF_AdvanceWSM6.cpp:112

Referenced by Update_State_Vars().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ Copy_State_to_Micro()

void WSM6::Copy_State_to_Micro ( const amrex::MultiFab &  cons_in)
overridevirtual

Reimplemented from NullMoist.

40 {
41  for (MFIter mfi(cons_in); mfi.isValid(); ++mfi) {
42  // Match Morrison behavior: refresh microphysics ghost zones from state.
43  // WSM6 Fortran reads the full (ims:ime, jms:jme, kms:kme) slab.
44  const auto& box3d = mfi.growntilebox();
45  auto states = cons_in.array(mfi);
46 
47  auto rho = mic_fab_vars[MicVar_WSM6::rho]->array(mfi);
48  auto theta = mic_fab_vars[MicVar_WSM6::theta]->array(mfi);
49  auto tabs = mic_fab_vars[MicVar_WSM6::tabs]->array(mfi);
50  auto pres = mic_fab_vars[MicVar_WSM6::pres]->array(mfi);
51 
52  auto qv = mic_fab_vars[MicVar_WSM6::qv]->array(mfi);
53  auto qc = mic_fab_vars[MicVar_WSM6::qc]->array(mfi);
54  auto qi = mic_fab_vars[MicVar_WSM6::qi]->array(mfi);
55  auto qr = mic_fab_vars[MicVar_WSM6::qr]->array(mfi);
56  auto qs = mic_fab_vars[MicVar_WSM6::qs]->array(mfi);
57  auto qg = mic_fab_vars[MicVar_WSM6::qg]->array(mfi);
58 
59  ParallelFor(box3d, [=] AMREX_GPU_DEVICE(int i, int j, int k) {
60  rho(i,j,k) = states(i,j,k,Rho_comp);
61  theta(i,j,k) = states(i,j,k,RhoTheta_comp) / states(i,j,k,Rho_comp);
62 
63  qv(i,j,k) = amrex::max(Real(0.0), states(i,j,k,RhoQ1_comp) / states(i,j,k,Rho_comp));
64  qc(i,j,k) = amrex::max(Real(0.0), states(i,j,k,RhoQ2_comp) / states(i,j,k,Rho_comp));
65  qi(i,j,k) = amrex::max(Real(0.0), states(i,j,k,RhoQ3_comp) / states(i,j,k,Rho_comp));
66  qr(i,j,k) = amrex::max(Real(0.0), states(i,j,k,RhoQ4_comp) / states(i,j,k,Rho_comp));
67  qs(i,j,k) = amrex::max(Real(0.0), states(i,j,k,RhoQ5_comp) / states(i,j,k,Rho_comp));
68  qg(i,j,k) = amrex::max(Real(0.0), states(i,j,k,RhoQ6_comp) / states(i,j,k,Rho_comp));
69 
70  tabs(i,j,k) = getTgivenRandRTh(states(i,j,k,Rho_comp),
71  states(i,j,k,RhoTheta_comp),
72  qv(i,j,k));
73  pres(i,j,k) = getPgivenRTh(states(i,j,k,RhoTheta_comp), qv(i,j,k));
74  });
75  }
76 }
getTgivenRandRTh
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getTgivenRandRTh(const amrex::Real rho, const amrex::Real rhotheta, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:46
getPgivenRTh
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real getPgivenRTh(const amrex::Real rhotheta, const amrex::Real qv=amrex::Real(0))
Definition: ERF_EOS.H:81
Rho_comp
#define Rho_comp
Definition: ERF_IndexDefines.H:36
MicVar_Kess::pres
@ pres
Definition: ERF_Kessler.H:26

Referenced by Update_Micro_Vars().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ Define()

void WSM6::Define ( SolverChoice &  sc)
inlineoverridevirtual

Reimplemented from NullMoist.

45  {
46  m_axis = sc.ave_plane;
47  m_do_cond = (!sc.use_shoc);
48  }
WSM6::m_axis
int m_axis
Definition: ERF_WSM6.H:130
SolverChoice::use_shoc
bool use_shoc
Definition: ERF_DataStruct.H:1253
SolverChoice::ave_plane
int ave_plane
Definition: ERF_DataStruct.H:1313

◆ Init()

void WSM6::Init ( const amrex::MultiFab &  cons_in,
const amrex::BoxArray &  grids,
const amrex::Geometry &  geom,
const amrex::Real &  dt_advance,
std::unique_ptr< amrex::MultiFab > &  z_phys_nd,
std::unique_ptr< amrex::MultiFab > &  detJ_cc 
)
overridevirtual

Reimplemented from NullMoist.

15 {
16  dt = dt_advance;
17  m_geom = geom;
18 
19  m_z_phys_nd = z_phys_nd.get();
20  m_detJ_cc = detJ_cc.get();
21 
22  MicVarMap.resize(m_qmoist_size);
23  MicVarMap = {MicVar_WSM6::rain_accum, MicVar_WSM6::snow_accum, MicVar_WSM6::graup_accum};
24 
25  for (int ivar = 0; ivar < MicVar_WSM6::NumVars; ++ivar) {
26  mic_fab_vars[ivar] = std::make_shared<MultiFab>(cons_in.boxArray(), cons_in.DistributionMap(),
27  1, cons_in.nGrowVect());
28  mic_fab_vars[ivar]->setVal(0.0);
29  }
30 
31  nlev = m_geom.Domain().length(2);
32  zlo = m_geom.Domain().smallEnd(2);
33  zhi = m_geom.Domain().bigEnd(2);
34 
35  initialize_coeffs();
36 }
WSM6::m_detJ_cc
amrex::MultiFab * m_detJ_cc
Definition: ERF_WSM6.H:134
WSM6::zlo
int zlo
Definition: ERF_WSM6.H:129
WSM6::nlev
int nlev
Definition: ERF_WSM6.H:129
WSM6::zhi
int zhi
Definition: ERF_WSM6.H:129
WSM6::m_qmoist_size
int m_qmoist_size
Definition: ERF_WSM6.H:120
WSM6::initialize_coeffs
void initialize_coeffs()
Definition: ERF_InitWSM6.cpp:79
WSM6::MicVarMap
amrex::Vector< int > MicVarMap
Definition: ERF_WSM6.H:124
MicVar_WSM6::NumVars
@ NumVars
Definition: ERF_WSM6.H:33

◆ initialize_coeffs()

void WSM6::initialize_coeffs ( )
private
80 {
81  using amrex::Real;
82 
83  // Exact port of Fortran rgmma() — ERF_module_mp_wsm6.F90 lines 1472-1495
84  // Weierstrass infinite product form for Gamma(x).
85  // Special case x==1 returns 0.0 matching Fortran exactly.
86  // Never triggered in practice (bvtr1=1.8, bvts1=1.41, etc.)
87  // but preserved for bit-compatible validation against Fortran reference.
88  // CPU-only: called from initialize_coeffs(), not from GPU kernels.
89  auto rgmma = [](Real x) -> Real {
90  if (x == Real(1.0)) return Real(0.0);
91  constexpr Real euler = Real(0.577215664901532);
92  Real result = x * std::exp(euler * x);
93  for (int i = 1; i <= 10000; ++i) {
94  Real y = Real(i);
95  result = result * (Real(1.0) + x/y) * std::exp(-x/y);
96  }
97  return Real(1.0) / result;
98  };
99 
100  // Physical constants matching mp_wsm6_init argument list
101  // den0: reference air density at 850mb (kg/m^3)
102  // denr: liquid water density — rhoh2o from ERF_Constants.H
103  // dens_arg: snow density — dens_snow constexpr = 100.0
104  // cl: specific heat liquid water — Cp_l from ERF_Constants.H
105  // cpv: specific heat water vapor — Cp_v from ERF_Constants.H
106  const Real den0 = Real(1.28);
107  const Real denr = Real(rhoh2o);
108  const Real dens_arg = dens_snow;
109  const Real cl = Real(Cp_l);
110  const Real cpv_loc = Real(Cp_v);
111 
112  // hail_opt branch: 5 regime-dependent coefficients
113  if (m_hail_opt) {
114  m_n0g = Real(4.0e4);
115  m_deng = Real(700.0);
116  m_avtg = Real(285.0);
117  m_bvtg = Real(0.8);
118  m_lamdagmax = Real(2.0e4);
119  } else {
120  m_n0g = Real(4.0e6);
121  m_deng = Real(500.0);
122  m_avtg = Real(330.0);
123  m_bvtg = Real(0.8);
124  m_lamdagmax = Real(6.0e4);
125  }
126 
127  m_pi_wsm6 = Real(4.0) * std::atan(Real(1.0));
128  m_xlv1 = cl - cpv_loc;
129 
130  m_qc0 = Real(4.0)/Real(3.0) * m_pi_wsm6 * denr
131  * std::pow(r0, Real(3.0)) * xncr / den0;
132  m_qck1 = Real(0.104) * Real(9.8) * peaut
133  / std::pow(xncr * denr, Real(1.0)/Real(3.0))
134  / xmyu * std::pow(den0, Real(4.0)/Real(3.0));
135  m_pidnc = m_pi_wsm6 * denr / Real(6.0);
136 
137  // Rain coefficients
138  m_bvtr1 = Real(1.0) + bvtr;
139  m_bvtr2 = Real(2.5) + Real(0.5) * bvtr;
140  m_bvtr3 = Real(3.0) + bvtr;
141  m_bvtr4 = Real(4.0) + bvtr;
142  m_bvtr6 = Real(6.0) + bvtr;
143  m_g1pbr = rgmma(m_bvtr1);
144  m_g3pbr = rgmma(m_bvtr3);
145  m_g4pbr = rgmma(m_bvtr4);
146  m_g6pbr = rgmma(m_bvtr6);
147  m_g5pbro2 = rgmma(m_bvtr2);
148  m_pvtr = avtr * m_g4pbr / Real(6.0);
149  m_eacrr = Real(1.0);
150  m_pacrr = m_pi_wsm6 * n0r * avtr * m_g3pbr * Real(0.25) * m_eacrr;
151  m_precr1 = Real(2.0) * m_pi_wsm6 * n0r * Real(0.78);
152  m_precr2 = Real(2.0) * m_pi_wsm6 * n0r * Real(0.31)
153  * std::pow(avtr, Real(0.5)) * m_g5pbro2;
154  m_roqimax = Real(2.08e22) * std::pow(dimax, Real(8.0));
155 
156  // Snow coefficients
157  m_bvts1 = Real(1.0) + bvts;
158  m_bvts2 = Real(2.5) + Real(0.5) * bvts;
159  m_bvts3 = Real(3.0) + bvts;
160  m_bvts4 = Real(4.0) + bvts;
161  m_g1pbs = rgmma(m_bvts1);
162  m_g3pbs = rgmma(m_bvts3);
163  m_g4pbs = rgmma(m_bvts4);
164  m_g5pbso2 = rgmma(m_bvts2);
165  m_pvts = avts * m_g4pbs / Real(6.0);
166  m_pacrs = m_pi_wsm6 * n0s * avts * m_g3pbs * Real(0.25);
167  m_precs1 = Real(4.0) * n0s * Real(0.65);
168  m_precs2 = Real(4.0) * n0s * Real(0.44)
169  * std::pow(avts, Real(0.5)) * m_g5pbso2;
170  m_pidn0r = m_pi_wsm6 * denr * n0r;
171  m_pidn0s = m_pi_wsm6 * dens_arg * n0s;
172  m_pacrc = m_pi_wsm6 * n0s * avts * m_g3pbs * Real(0.25) * eacrc;
173 
174  // Graupel/hail coefficients
175  m_bvtg1 = Real(1.0) + m_bvtg;
176  m_bvtg2 = Real(2.5) + Real(0.5) * m_bvtg;
177  m_bvtg3 = Real(3.0) + m_bvtg;
178  m_bvtg4 = Real(4.0) + m_bvtg;
179  m_g1pbg = rgmma(m_bvtg1);
180  m_g3pbg = rgmma(m_bvtg3);
181  m_g4pbg = rgmma(m_bvtg4);
182  m_pacrg = m_pi_wsm6 * m_n0g * m_avtg * m_g3pbg * Real(0.25);
183  m_g5pbgo2 = rgmma(m_bvtg2);
184  m_pvtg = m_avtg * m_g4pbg / Real(6.0);
185  m_precg1 = Real(2.0) * m_pi_wsm6 * m_n0g * Real(0.78);
186  m_precg2 = Real(2.0) * m_pi_wsm6 * m_n0g * Real(0.31)
187  * std::pow(m_avtg, Real(0.5)) * m_g5pbgo2;
188  m_pidn0g = m_pi_wsm6 * m_deng * m_n0g;
189 
190  // Slope parameter limits
191  m_rslopermax = Real(1.0) / lamdarmax;
192  m_rslopesmax = Real(1.0) / lamdasmax;
193  m_rslopegmax = Real(1.0) / m_lamdagmax;
194  m_rsloperbmax = std::pow(m_rslopermax, bvtr);
195  m_rslopesbmax = std::pow(m_rslopesmax, bvts);
196  m_rslopegbmax = std::pow(m_rslopegmax, m_bvtg);
197  m_rsloper2max = m_rslopermax * m_rslopermax;
198  m_rslopes2max = m_rslopesmax * m_rslopesmax;
199  m_rslopeg2max = m_rslopegmax * m_rslopegmax;
200  m_rsloper3max = m_rsloper2max * m_rslopermax;
201  m_rslopes3max = m_rslopes2max * m_rslopesmax;
202  m_rslopeg3max = m_rslopeg2max * m_rslopegmax;
203 }
Coord::y
@ y
Coord::x
@ x
WSM6::m_g5pbgo2
amrex::Real m_g5pbgo2
Definition: ERF_WSM6.H:157
WSM6::m_bvts3
amrex::Real m_bvts3
Definition: ERF_WSM6.H:152
WSM6::m_hail_opt
bool m_hail_opt
Definition: ERF_WSM6.H:140
WSM6::m_bvtr1
amrex::Real m_bvtr1
Definition: ERF_WSM6.H:148
WSM6::m_g4pbs
amrex::Real m_g4pbs
Definition: ERF_WSM6.H:153
WSM6::m_g1pbg
amrex::Real m_g1pbg
Definition: ERF_WSM6.H:157
WSM6::m_g4pbg
amrex::Real m_g4pbg
Definition: ERF_WSM6.H:157
WSM6::m_g5pbro2
amrex::Real m_g5pbro2
Definition: ERF_WSM6.H:149
WSM6::m_g3pbr
amrex::Real m_g3pbr
Definition: ERF_WSM6.H:149
WSM6::m_bvtg1
amrex::Real m_bvtg1
Definition: ERF_WSM6.H:156
WSM6::m_bvts4
amrex::Real m_bvts4
Definition: ERF_WSM6.H:152
WSM6::m_deng
amrex::Real m_deng
Definition: ERF_WSM6.H:141
WSM6::m_g3pbg
amrex::Real m_g3pbg
Definition: ERF_WSM6.H:157
WSM6::lamdasmax
static constexpr amrex::Real lamdasmax
Definition: ERF_WSM6.H:62
WSM6::m_bvtr6
amrex::Real m_bvtr6
Definition: ERF_WSM6.H:148
WSM6::m_bvtr2
amrex::Real m_bvtr2
Definition: ERF_WSM6.H:148
WSM6::m_bvts1
amrex::Real m_bvts1
Definition: ERF_WSM6.H:152
WSM6::m_pacrs
amrex::Real m_pacrs
Definition: ERF_WSM6.H:154
WSM6::m_eacrr
amrex::Real m_eacrr
Definition: ERF_WSM6.H:150
WSM6::m_g1pbs
amrex::Real m_g1pbs
Definition: ERF_WSM6.H:153
WSM6::m_lamdagmax
amrex::Real m_lamdagmax
Definition: ERF_WSM6.H:141
WSM6::m_bvtg2
amrex::Real m_bvtg2
Definition: ERF_WSM6.H:156
WSM6::r0
static constexpr amrex::Real r0
Definition: ERF_WSM6.H:55
WSM6::m_bvtg4
amrex::Real m_bvtg4
Definition: ERF_WSM6.H:156
WSM6::m_bvtr4
amrex::Real m_bvtr4
Definition: ERF_WSM6.H:148
WSM6::m_bvts2
amrex::Real m_bvts2
Definition: ERF_WSM6.H:152
WSM6::m_bvtg3
amrex::Real m_bvtg3
Definition: ERF_WSM6.H:156
WSM6::m_pidnc
amrex::Real m_pidnc
Definition: ERF_WSM6.H:147
WSM6::m_g3pbs
amrex::Real m_g3pbs
Definition: ERF_WSM6.H:153
WSM6::lamdarmax
static constexpr amrex::Real lamdarmax
Definition: ERF_WSM6.H:61
WSM6::m_avtg
amrex::Real m_avtg
Definition: ERF_WSM6.H:141
WSM6::m_g4pbr
amrex::Real m_g4pbr
Definition: ERF_WSM6.H:149
WSM6::m_g1pbr
amrex::Real m_g1pbr
Definition: ERF_WSM6.H:149
WSM6::peaut
static constexpr amrex::Real peaut
Definition: ERF_WSM6.H:56
WSM6::avts
static constexpr amrex::Real avts
Definition: ERF_WSM6.H:59
WSM6::eacrc
static constexpr amrex::Real eacrc
Definition: ERF_WSM6.H:68
WSM6::m_g5pbso2
amrex::Real m_g5pbso2
Definition: ERF_WSM6.H:153
WSM6::xmyu
static constexpr amrex::Real xmyu
Definition: ERF_WSM6.H:58
WSM6::m_bvtr3
amrex::Real m_bvtr3
Definition: ERF_WSM6.H:148
mp_wsm6::rgmma
real(kind=kind_phys) function rgmma(x)
Definition: ERF_module_mp_wsm6.F90:1584
Here is the call graph for this function:

◆ Qmoist_Ptr()

amrex::MultiFab* WSM6::Qmoist_Ptr ( const int &  varIdx)
inlineoverridevirtual

Reimplemented from NullMoist.

102  {
103  AMREX_ALWAYS_ASSERT(varIdx < m_qmoist_size);
104  return mic_fab_vars[MicVarMap[varIdx]].get();
105  }
AMREX_ALWAYS_ASSERT
AMREX_ALWAYS_ASSERT(bx.length()[2]==khi+1)
Here is the call graph for this function:

◆ Qmoist_Restart_Vars()

void WSM6::Qmoist_Restart_Vars ( const SolverChoice &  ,
std::vector< int > &  a_idx,
std::vector< std::string > &  a_names 
) const
inlineoverridevirtual

Reimplemented from NullMoist.

114  {
115  a_idx = {0, 1, 2};
116  a_names = {"RainAccum", "SnowAccum", "GraupAccum"};
117  }

◆ Qmoist_Size()

int WSM6::Qmoist_Size ( )
inlineoverridevirtual

Reimplemented from NullMoist.

107 { return m_qmoist_size; }

◆ Qstate_Moist_NumConc_Size()

int WSM6::Qstate_Moist_NumConc_Size ( )
inlineoverridevirtual

Reimplemented from NullMoist.

109 { return n_qstate_moist_numconc_size; }
WSM6::n_qstate_moist_numconc_size
int n_qstate_moist_numconc_size
Definition: ERF_WSM6.H:122

◆ Qstate_Moist_Size()

int WSM6::Qstate_Moist_Size ( )
inlineoverridevirtual

Reimplemented from NullMoist.

108 { return n_qstate_moist_size; }
WSM6::n_qstate_moist_size
int n_qstate_moist_size
Definition: ERF_WSM6.H:121

◆ Set_dzmin()

void WSM6::Set_dzmin ( const amrex::Real  dz_min)
inlineoverridevirtual

Reimplemented from NullMoist.

82 { m_dzmin = dz_min; }
WSM6::m_dzmin
amrex::Real m_dzmin
Definition: ERF_WSM6.H:128

◆ Update_Micro_Vars()

void WSM6::Update_Micro_Vars ( amrex::MultiFab &  cons_in)
inlineoverridevirtual

Reimplemented from NullMoist.

88  {
89  Copy_State_to_Micro(cons_in);
90  }
WSM6::Copy_State_to_Micro
void Copy_State_to_Micro(const amrex::MultiFab &cons_in) override
Definition: ERF_InitWSM6.cpp:39
Here is the call graph for this function:

◆ Update_State_Vars()

void WSM6::Update_State_Vars ( amrex::MultiFab &  cons_in,
const amrex::MultiFab &   
)
inlineoverridevirtual

Reimplemented from NullMoist.

94  {
95  Copy_Micro_to_State(cons_in);
96  }
WSM6::Copy_Micro_to_State
void Copy_Micro_to_State(amrex::MultiFab &cons_in) override
Definition: ERF_UpdateWSM6.cpp:7
Here is the call graph for this function:

Member Data Documentation

◆ alpha_wsm6

constexpr amrex::Real WSM6::alpha_wsm6 = amrex::Real(0.12)
staticconstexpr

◆ avtr

constexpr amrex::Real WSM6::avtr = amrex::Real(841.9)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm_scratch().

◆ avts

constexpr amrex::Real WSM6::avts = amrex::Real(11.72)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm6_scratch().

◆ bvtr

constexpr amrex::Real WSM6::bvtr = amrex::Real(0.8)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm_scratch().

◆ bvts

constexpr amrex::Real WSM6::bvts = amrex::Real(0.41)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm6_scratch().

◆ dens_snow

constexpr amrex::Real WSM6::dens_snow = amrex::Real(100.0)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm6_scratch().

◆ dicon

constexpr amrex::Real WSM6::dicon = amrex::Real(11.9)
staticconstexpr

◆ dimax

constexpr amrex::Real WSM6::dimax = amrex::Real(500.0e-6)
staticconstexpr

◆ dt

amrex::Real WSM6::dt {0.0}
private

◆ dtcldcr

constexpr amrex::Real WSM6::dtcldcr = amrex::Real(120.0)
staticconstexpr

◆ eacrc

constexpr amrex::Real WSM6::eacrc = amrex::Real(1.0)
staticconstexpr

◆ lamdarmax

constexpr amrex::Real WSM6::lamdarmax = amrex::Real(8.0e4)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm_scratch().

◆ lamdasmax

constexpr amrex::Real WSM6::lamdasmax = amrex::Real(1.0e5)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm6_scratch().

◆ m_avtg

amrex::Real WSM6::m_avtg {0}
private

◆ m_axis

int WSM6::m_axis {2}
private

Referenced by Define().

◆ m_bvtg

amrex::Real WSM6::m_bvtg {0}
private

◆ m_bvtg1

amrex::Real WSM6::m_bvtg1 {0}
private

◆ m_bvtg2

amrex::Real WSM6::m_bvtg2 {0}
private

◆ m_bvtg3

amrex::Real WSM6::m_bvtg3 {0}
private

◆ m_bvtg4

amrex::Real WSM6::m_bvtg4 {0}
private

◆ m_bvtr1

amrex::Real WSM6::m_bvtr1 {0}
private

◆ m_bvtr2

amrex::Real WSM6::m_bvtr2 {0}
private

◆ m_bvtr3

amrex::Real WSM6::m_bvtr3 {0}
private

◆ m_bvtr4

amrex::Real WSM6::m_bvtr4 {0}
private

◆ m_bvtr6

amrex::Real WSM6::m_bvtr6 {0}
private

◆ m_bvts1

amrex::Real WSM6::m_bvts1 {0}
private

◆ m_bvts2

amrex::Real WSM6::m_bvts2 {0}
private

◆ m_bvts3

amrex::Real WSM6::m_bvts3 {0}
private

◆ m_bvts4

amrex::Real WSM6::m_bvts4 {0}
private

◆ m_deng

amrex::Real WSM6::m_deng {0}
private

◆ m_detJ_cc

amrex::MultiFab* WSM6::m_detJ_cc {nullptr}
private

◆ m_do_cond

bool WSM6::m_do_cond {true}
private

Referenced by Define().

◆ m_dzmin

amrex::Real WSM6::m_dzmin {0.0}
private

Referenced by Set_dzmin().

◆ m_eacrr

amrex::Real WSM6::m_eacrr {0}
private

◆ m_g1pbg

amrex::Real WSM6::m_g1pbg {0}
private

◆ m_g1pbr

amrex::Real WSM6::m_g1pbr {0}
private

◆ m_g1pbs

amrex::Real WSM6::m_g1pbs {0}
private

◆ m_g3pbg

amrex::Real WSM6::m_g3pbg {0}
private

◆ m_g3pbr

amrex::Real WSM6::m_g3pbr {0}
private

◆ m_g3pbs

amrex::Real WSM6::m_g3pbs {0}
private

◆ m_g4pbg

amrex::Real WSM6::m_g4pbg {0}
private

◆ m_g4pbr

amrex::Real WSM6::m_g4pbr {0}
private

◆ m_g4pbs

amrex::Real WSM6::m_g4pbs {0}
private

◆ m_g5pbgo2

amrex::Real WSM6::m_g5pbgo2 {0}
private

◆ m_g5pbro2

amrex::Real WSM6::m_g5pbro2 {0}
private

◆ m_g5pbso2

amrex::Real WSM6::m_g5pbso2 {0}
private

◆ m_g6pbr

amrex::Real WSM6::m_g6pbr {0}
private

◆ m_geom

amrex::Geometry WSM6::m_geom
private

◆ m_hail_opt

bool WSM6::m_hail_opt {false}
private

◆ m_lamdagmax

amrex::Real WSM6::m_lamdagmax {0}
private

◆ m_n0g

amrex::Real WSM6::m_n0g {0}
private

◆ m_pacrc

amrex::Real WSM6::m_pacrc {0}
private

◆ m_pacrg

amrex::Real WSM6::m_pacrg {0}
private

◆ m_pacrr

amrex::Real WSM6::m_pacrr {0}
private

◆ m_pacrs

amrex::Real WSM6::m_pacrs {0}
private

◆ m_pi_wsm6

amrex::Real WSM6::m_pi_wsm6 {0}
private

◆ m_pidn0g

amrex::Real WSM6::m_pidn0g {0}
private

◆ m_pidn0r

amrex::Real WSM6::m_pidn0r {0}
private

◆ m_pidn0s

amrex::Real WSM6::m_pidn0s {0}
private

◆ m_pidnc

amrex::Real WSM6::m_pidnc {0}
private

◆ m_precg1

amrex::Real WSM6::m_precg1 {0}
private

◆ m_precg2

amrex::Real WSM6::m_precg2 {0}
private

◆ m_precr1

amrex::Real WSM6::m_precr1 {0}
private

◆ m_precr2

amrex::Real WSM6::m_precr2 {0}
private

◆ m_precs1

amrex::Real WSM6::m_precs1 {0}
private

◆ m_precs2

amrex::Real WSM6::m_precs2 {0}
private

◆ m_pvtg

amrex::Real WSM6::m_pvtg {0}
private

◆ m_pvtr

amrex::Real WSM6::m_pvtr {0}
private

◆ m_pvts

amrex::Real WSM6::m_pvts {0}
private

◆ m_qc0

amrex::Real WSM6::m_qc0 {0}
private

◆ m_qck1

amrex::Real WSM6::m_qck1 {0}
private

◆ m_qmoist_size

int WSM6::m_qmoist_size = 3
private

Referenced by Qmoist_Ptr(), and Qmoist_Size().

◆ m_roqimax

amrex::Real WSM6::m_roqimax {0}
private

◆ m_rslopeg2max

amrex::Real WSM6::m_rslopeg2max {0}
private

◆ m_rslopeg3max

amrex::Real WSM6::m_rslopeg3max {0}
private

◆ m_rslopegbmax

amrex::Real WSM6::m_rslopegbmax {0}
private

◆ m_rslopegmax

amrex::Real WSM6::m_rslopegmax {0}
private

◆ m_rsloper2max

amrex::Real WSM6::m_rsloper2max {0}
private

◆ m_rsloper3max

amrex::Real WSM6::m_rsloper3max {0}
private

◆ m_rsloperbmax

amrex::Real WSM6::m_rsloperbmax {0}
private

◆ m_rslopermax

amrex::Real WSM6::m_rslopermax {0}
private

◆ m_rslopes2max

amrex::Real WSM6::m_rslopes2max {0}
private

◆ m_rslopes3max

amrex::Real WSM6::m_rslopes3max {0}
private

◆ m_rslopesbmax

amrex::Real WSM6::m_rslopesbmax {0}
private

◆ m_rslopesmax

amrex::Real WSM6::m_rslopesmax {0}
private

◆ m_xlv1

amrex::Real WSM6::m_xlv1 {0}
private

◆ m_z_phys_nd

amrex::MultiFab* WSM6::m_z_phys_nd {nullptr}
private

◆ mic_fab_vars

amrex::Array<FabPtr, MicVar_WSM6::NumVars> WSM6::mic_fab_vars
private

Referenced by Qmoist_Ptr().

◆ MicVarMap

amrex::Vector<int> WSM6::MicVarMap
private

Referenced by Qmoist_Ptr().

◆ n0r

constexpr amrex::Real WSM6::n0r = amrex::Real(8.0e6)
staticconstexpr

Referenced by Advance(), and wsm6_nislfv_rain_plm_scratch().

◆ n0s

constexpr amrex::Real WSM6::n0s = amrex::Real(2.0e6)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm6_scratch().

◆ n0smax

constexpr amrex::Real WSM6::n0smax = amrex::Real(1.0e11)
staticconstexpr

◆ n_qstate_moist_numconc_size

int WSM6::n_qstate_moist_numconc_size = 0
private

Referenced by Qstate_Moist_NumConc_Size().

◆ n_qstate_moist_size

int WSM6::n_qstate_moist_size = 6
private

Referenced by Qstate_Moist_Size().

◆ nlev

int WSM6::nlev {0}
private

◆ peaut

constexpr amrex::Real WSM6::peaut = amrex::Real(0.55)
staticconstexpr

◆ pfrz1

constexpr amrex::Real WSM6::pfrz1 = amrex::Real(100.0)
staticconstexpr

◆ pfrz2

constexpr amrex::Real WSM6::pfrz2 = amrex::Real(0.66)
staticconstexpr

◆ qcrmin

constexpr amrex::Real WSM6::qcrmin = amrex::Real(1.0e-9)
staticconstexpr

Referenced by wsm6_nislfv_rain_plm6_scratch(), and wsm6_nislfv_rain_plm_scratch().

◆ qs0

constexpr amrex::Real WSM6::qs0 = amrex::Real(6.0e-4)
staticconstexpr

◆ r0

constexpr amrex::Real WSM6::r0 = amrex::Real(0.8e-5)
staticconstexpr

◆ xmyu

constexpr amrex::Real WSM6::xmyu = amrex::Real(1.718e-5)
staticconstexpr

◆ xncr

constexpr amrex::Real WSM6::xncr = amrex::Real(3.0e8)
staticconstexpr

Referenced by Advance().

◆ zhi

int WSM6::zhi {0}
private

◆ zlo

int WSM6::zlo {0}
private
The documentation for this class was generated from the following files: