From c359584686d9ca8107f96758c8c5268ee3d6c3e9 Mon Sep 17 00:00:00 2001 From: Hassan Beydoun Date: Thu, 11 Jul 2024 10:08:40 -0700 Subject: [PATCH 1/3] Lin Lin identified a bug in computing rain drop effective diameter as well as a bug in multiplying the immersion freezing number rate by cld_frac_l rather than the ooverlapped fraction between liquid and ice. Both bugs are fixed here --- components/eam/src/physics/p3/scream/micro_p3.F90 | 4 ++-- .../src/physics/p3/impl/p3_back_to_cell_average_impl.hpp | 2 +- components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp | 2 +- 3 files changed, 4 insertions(+), 4 deletions(-) diff --git a/components/eam/src/physics/p3/scream/micro_p3.F90 b/components/eam/src/physics/p3/scream/micro_p3.F90 index f345e5abcf7..17c33a9f4fc 100644 --- a/components/eam/src/physics/p3/scream/micro_p3.F90 +++ b/components/eam/src/physics/p3/scream/micro_p3.F90 @@ -1030,7 +1030,7 @@ subroutine p3_main_part3(kts, kte, kbot, ktop, kdir, & ze_rain(k) = nr(k)*(mu_r(k)+6._rtype)*(mu_r(k)+5._rtype)*(mu_r(k)+4._rtype)* & (mu_r(k)+3._rtype)*(mu_r(k)+2._rtype)*(mu_r(k)+1._rtype)/bfb_pow(lamr(k), 6._rtype) ze_rain(k) = max(ze_rain(k),1.e-22_rtype) - diag_eff_radius_qr(k) = 1.5_rtype/lamr(k) + diag_eff_radius_qr(k) = 0.5_rtype*(mu_r(k)+3._rtype)/lamr(k) else qv(k) = qv(k)+qr(k) th_atm(k) = th_atm(k)-inv_exner(k)*qr(k)*latent_heat_vapor(k)*inv_cp @@ -2839,7 +2839,7 @@ subroutine back_to_cell_average(cld_frac_l,cld_frac_r,cld_frac_i, ni2nr_melt_tend = ni2nr_melt_tend*cld_frac_i ! Change in number due to melting nc_collect_tend = nc_collect_tend*il_cldm ! Cloud # change due to collection of cld water by ice ncshdc = ncshdc*il_cldm ! Number change due to shedding, occurs when ice and liquid interact - nc2ni_immers_freeze_tend = nc2ni_immers_freeze_tend*cld_frac_l ! Number change associated with freexzing of cld drops + nc2ni_immers_freeze_tend = nc2ni_immers_freeze_tend*il_cldm ! Number change associated with freexzing of cld drops nr_collect_tend = nr_collect_tend*ir_cldm ! Rain number change due to collection from ice ni_selfcollect_tend = ni_selfcollect_tend*cld_frac_i ! Ice self collection qidep = qidep*cld_frac_i ! Vapor deposition to ice phase diff --git a/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp b/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp index 0ba2e8b412a..3f3f9e262a2 100644 --- a/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp +++ b/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp @@ -62,7 +62,7 @@ ::back_to_cell_average( ni2nr_melt_tend.set(context, ni2nr_melt_tend * cld_frac_i); // Change in number due to melting nc_collect_tend.set(context, nc_collect_tend * il_cldm); // Cloud # change due to collection of cld water by ice ncshdc.set(context, ncshdc * il_cldm); // Number change due to shedding, occurs when ice and liquid interact - nc2ni_immers_freeze_tend.set(context, nc2ni_immers_freeze_tend * cld_frac_l); // Number change associated with freexzing of cld drops + nc2ni_immers_freeze_tend.set(context, nc2ni_immers_freeze_tend * il_cldm); // Number change associated with freexzing of cld drops nr_collect_tend.set(context, nr_collect_tend * ir_cldm); // Rain number change due to collection from ice ni_selfcollect_tend.set(context, ni_selfcollect_tend * cld_frac_i); // Ice self collection qv2qi_vapdep_tend.set(context, qv2qi_vapdep_tend * cld_frac_i); // Vapor deposition to ice phase diff --git a/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp b/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp index febb00ef4d3..855c6c673ba 100644 --- a/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp +++ b/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp @@ -118,7 +118,7 @@ ::p3_main_part3( ze_rain(k).set(qr_gt_small, nr(k)*(mu_r(k)+6)*(mu_r(k)+5)*(mu_r(k)+4)* (mu_r(k)+3)*(mu_r(k)+2)*(mu_r(k)+1)/pow(lamr(k), sp(6.0))); // once f90 is gone, 6 can be int ze_rain(k).set(qr_gt_small, max(ze_rain(k), sp(1.e-22))); - diag_eff_radius_qr(k).set(qr_gt_small, sp(1.5) / lamr(k)); + diag_eff_radius_qc(k).set(qr_gt_small, sp(0.5) * (mu_r(k) + 3) / lamr(k)); } if (qr_small.any()) { From 3ef56352f09d3acf1801c456ec017e1e674f038d Mon Sep 17 00:00:00 2001 From: Hassan Beydoun Date: Fri, 12 Jul 2024 12:24:34 -0700 Subject: [PATCH 2/3] both qc2qi_immersion and nc2ni_immersion are now mutliplied by cld_frac_l. Also added comment on how to compute effective radius --- components/eam/src/physics/p3/scream/micro_p3.F90 | 4 ++-- .../src/physics/p3/impl/p3_back_to_cell_average_impl.hpp | 4 ++-- components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp | 1 + 3 files changed, 5 insertions(+), 4 deletions(-) diff --git a/components/eam/src/physics/p3/scream/micro_p3.F90 b/components/eam/src/physics/p3/scream/micro_p3.F90 index 17c33a9f4fc..6b6b6b49a68 100644 --- a/components/eam/src/physics/p3/scream/micro_p3.F90 +++ b/components/eam/src/physics/p3/scream/micro_p3.F90 @@ -2830,7 +2830,7 @@ subroutine back_to_cell_average(cld_frac_l,cld_frac_r,cld_frac_i, ! map ice-phase process rates to cell-avg qi2qv_sublim_tend = qi2qv_sublim_tend*cld_frac_i ! Sublimation of ice in ice cloud nr_ice_shed_tend = nr_ice_shed_tend*il_cldm ! Rain # increase due to shedding from rain-ice collisions, occurs when ice and liquid interact - qc2qi_hetero_freeze_tend = qc2qi_hetero_freeze_tend*il_cldm ! Immersion freezing of cloud drops + qc2qi_hetero_freeze_tend = qc2qi_hetero_freeze_tend*cld_frac_l ! Immersion freezing of cloud drops qrcol = qrcol*ir_cldm ! Collection of rain mass by ice qc2qr_ice_shed_tend = qc2qr_ice_shed_tend*il_cldm ! Rain mass growth due to shedding of fain drops after collisions with ice, occurs when ice and liquid interact qi2qr_melt_tend = qi2qr_melt_tend*cld_frac_i ! Melting of ice @@ -2839,7 +2839,7 @@ subroutine back_to_cell_average(cld_frac_l,cld_frac_r,cld_frac_i, ni2nr_melt_tend = ni2nr_melt_tend*cld_frac_i ! Change in number due to melting nc_collect_tend = nc_collect_tend*il_cldm ! Cloud # change due to collection of cld water by ice ncshdc = ncshdc*il_cldm ! Number change due to shedding, occurs when ice and liquid interact - nc2ni_immers_freeze_tend = nc2ni_immers_freeze_tend*il_cldm ! Number change associated with freexzing of cld drops + nc2ni_immers_freeze_tend = nc2ni_immers_freeze_tend*cld_frac_l ! Number change associated with freexzing of cld drops nr_collect_tend = nr_collect_tend*ir_cldm ! Rain number change due to collection from ice ni_selfcollect_tend = ni_selfcollect_tend*cld_frac_i ! Ice self collection qidep = qidep*cld_frac_i ! Vapor deposition to ice phase diff --git a/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp b/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp index 3f3f9e262a2..014793bd27a 100644 --- a/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp +++ b/components/eamxx/src/physics/p3/impl/p3_back_to_cell_average_impl.hpp @@ -53,7 +53,7 @@ ::back_to_cell_average( // map ice-phase process rates to cell-avg qi2qv_sublim_tend.set(context, qi2qv_sublim_tend * cld_frac_i); // Sublimation of ice in ice cloud nr_ice_shed_tend.set(context, nr_ice_shed_tend * il_cldm); // Rain # increase due to shedding from rain-ice collisions, occurs when ice and liquid interact - qc2qi_hetero_freeze_tend.set(context, qc2qi_hetero_freeze_tend * il_cldm); // Immersion freezing of cloud drops + qc2qi_hetero_freeze_tend.set(context, qc2qi_hetero_freeze_tend * cld_frac_l); // Immersion freezing of cloud drops qr2qi_collect_tend.set(context, qr2qi_collect_tend * ir_cldm); // Collection of rain mass by ice qc2qr_ice_shed_tend.set(context, qc2qr_ice_shed_tend * il_cldm); // Rain mass growth due to shedding of fain drops after collisions with ice, occurs when ice and liquid interact qi2qr_melt_tend.set(context, qi2qr_melt_tend * cld_frac_i); // Melting of ice @@ -62,7 +62,7 @@ ::back_to_cell_average( ni2nr_melt_tend.set(context, ni2nr_melt_tend * cld_frac_i); // Change in number due to melting nc_collect_tend.set(context, nc_collect_tend * il_cldm); // Cloud # change due to collection of cld water by ice ncshdc.set(context, ncshdc * il_cldm); // Number change due to shedding, occurs when ice and liquid interact - nc2ni_immers_freeze_tend.set(context, nc2ni_immers_freeze_tend * il_cldm); // Number change associated with freexzing of cld drops + nc2ni_immers_freeze_tend.set(context, nc2ni_immers_freeze_tend * cld_frac_l); // Number change associated with freexzing of cld drops nr_collect_tend.set(context, nr_collect_tend * ir_cldm); // Rain number change due to collection from ice ni_selfcollect_tend.set(context, ni_selfcollect_tend * cld_frac_i); // Ice self collection qv2qi_vapdep_tend.set(context, qv2qi_vapdep_tend * cld_frac_i); // Vapor deposition to ice phase diff --git a/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp b/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp index 855c6c673ba..31f6f0b7e39 100644 --- a/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp +++ b/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp @@ -89,6 +89,7 @@ ::p3_main_part3( if (qc_gt_small.any()) { nc(k).set(qc_gt_small,nc_incld*cld_frac_l(k)); //cld_dsd2 might have changed incld nc... need consistency. + //diag_eff_radius_qc is obtained by diving the 3rd and 2nd moments of the DSD e.g., eqn 5 of MG2008 diag_eff_radius_qc(k).set(qc_gt_small, sp(0.5) * (mu_c(k) + 3) / lamc(k)); } if (qc_small.any()) { From 19d31c7eb1320541fac4f9f3f71ae38435c81c5e Mon Sep 17 00:00:00 2001 From: Hassan Beydoun Date: Fri, 12 Jul 2024 12:32:34 -0700 Subject: [PATCH 3/3] Accidently assigned qr effective radius calculation to p3. fixed now --- components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp b/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp index 31f6f0b7e39..26a8399f71d 100644 --- a/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp +++ b/components/eamxx/src/physics/p3/impl/p3_main_impl_part3.hpp @@ -119,7 +119,7 @@ ::p3_main_part3( ze_rain(k).set(qr_gt_small, nr(k)*(mu_r(k)+6)*(mu_r(k)+5)*(mu_r(k)+4)* (mu_r(k)+3)*(mu_r(k)+2)*(mu_r(k)+1)/pow(lamr(k), sp(6.0))); // once f90 is gone, 6 can be int ze_rain(k).set(qr_gt_small, max(ze_rain(k), sp(1.e-22))); - diag_eff_radius_qc(k).set(qr_gt_small, sp(0.5) * (mu_r(k) + 3) / lamr(k)); + diag_eff_radius_qr(k).set(qr_gt_small, sp(0.5) * (mu_r(k) + 3) / lamr(k)); } if (qr_small.any()) {