EAMxx: Adds aerosols heterogeneous freezing calculations in P3 microphysics #6947

singhbalwinder · 2025-01-25T23:02:43Z

The heterogeneous freezing calculations from prognostics aerosols are
added to P3 microphysics. Setting use_hetfrz_classnuc to true
will turn on these calculations. Otherwise, P3 will use the default
prescribed aerosol calculations.

[BFB] for EAM and EAMxx

github-actions · 2025-01-25T23:04:53Z

PR Preview Action v1.6.0
🚀 View preview at https://E3SM-Project.github.io/E3SM/pr-preview/pr-6947/
Built to branch `gh-pages` at 2025-02-04 02:10 UTC. Preview will be ready when the GitHub Pages deployment is complete.

singhbalwinder · 2025-01-25T23:08:09Z

TODO:

Turn off this feature for default EAMxx
Revive commented-out P3 tests after adding missing arguments in various function signatures.

mahf708 · 2025-01-29T19:05:27Z

Qucik comments:

Please turn off the feature by default
Please follow for the do_ice_production procedure (as an example) for passing the flag inside
Please hide most (if not all) additions inside if-else guards (with the new flag), for example, the add_required calls and such
Please keep tests intact if you intend to integrate this
Also, ensure you don't break PAM/MMF2 (I am ~~100%~~ almost certain you're currently breaking it)

bartgol

I have a few comments. Mostly: why are lots of unit tests now commented?

bartgol · 2025-01-29T20:16:50Z

components/eamxx/cime_config/namelist_defaults_scream.xml

@@ -192,6 +192,7 @@ be lost if SCREAM_HACK_XML is not enabled.
    <!-- P3 microphysics -->
    <p3 inherit="atm_proc_base">
      <do_prescribed_ccn>true</do_prescribed_ccn>
+      <use_hetfrz_classnuc>true</use_hetfrz_classnuc>


Let's add metadata to new options. In particular, type and doc.

Sure. We will add the following info:

doc=“Switch to turn on heterogeneous freezing due to prognostic aerosols” type=”logical”

bartgol · 2025-01-29T20:23:37Z

components/eamxx/src/physics/p3/eamxx_p3_process_interface.cpp

@@ -44,6 +44,7 @@ void P3Microphysics::set_grids(const std::shared_ptr<const GridsManager> grids_m
  infrastructure.kte = m_num_levs-1;
  infrastructure.predictNc = m_params.get<bool>("do_predict_nc",true);
  infrastructure.prescribedCCN = m_params.get<bool>("do_prescribed_ccn",true);
+  infrastructure.use_hetfrz_classnuc = m_params.get<bool>("use_hetfrz_classnuc",true);


@jgfouca, @AaronDonahue, @tcclevenger, asking your opinion as well

I don't know how I feel about having defaults in the source code (regardless of whether their value). We should distinguish CIME runs from standalone.

CIME RUNS:

We already put defaults in namelist_defaults_scream.xml, so here we should not have to set defaults. Sure, since the value is in m_params (b/c of buildnml), the default is not used. But what if at some point a bug sneaks in, and the input yaml file does not contain the value. What if we start using the hard-coded default without knowing?

STANDALONE RUNS:

We prob don't want to specify ALL options in our input.yaml files, to keep them at a minimum. Perhaps a compromise would be to do this in the atm procs constructors:

P3Microphysics::P3Microphysics (const ekat::Comm& comm, const ekat::ParameterList& params) : AtmosphereProcess(comm, params) { #ifndef SCREAM_CIME_BUILD // set defaults for parameters, so unit tests can use small-ish input files m_params.get<bool>("do_predict_nc",true); // use get with default to "set if not set" ... #endif }

Then, we remove ALL hard-coded defaults in the rest of the parametrization.

@bartgol I think most (if not all) of these should be moved to the runtime struct and initialize stuff in that struct. I definitely will have a hard time approving this PR if use_hetfrz_classnuc is added here and not near the do_ice_production flag

We will revisit and mimic how the do_ice_production flag is used.

bartgol · 2025-01-29T20:25:57Z

components/eamxx/src/physics/p3/impl/p3_CNT_couple_impl.hpp

+      ninuc_cnt.set(mask, frzdep*1.0e6/rho, Zero);
+      qinuc_cnt.set(mask, ninuc_cnt*mi0, Zero);
+      break;
+    default:


Are you really ok doing nothing if Iflag is neither 1 nor 2? If that should be an error, please add EKAT_KERNEL_ERROR(...) in the default clause. If the default case should indeed do nothing, then maybe you should have an early return, avoiding pointless calculations before...

Thanks. We should add EKAT_KERNEL_ERROR(...) as in the current logic, its value can only be 1 or 2. If it has some other value, it should error out.

Since you're calling this function twice with 1 and and then with 2, why not just do the calculations in one function call and avoid switch case here?

if (use_hetfrz_classnuc){ CNT_couple(hetfrz_immersion_nucleation_tend(k), hetfrz_contact_nucleation_tend(k), hetfrz_deposition_nucleation_tend(k), rho(k), qc_incld(k), nc_incld(k), 1, ncheti_cnt, qcheti_cnt, nicnt, qicnt, ninuc_cnt, qinuc_cnt); CNT_couple(hetfrz_immersion_nucleation_tend(k), hetfrz_contact_nucleation_tend(k), hetfrz_deposition_nucleation_tend(k), rho(k), qc_incld(k), nc_incld(k), 2, ncheti_cnt, qcheti_cnt, nicnt, qicnt, ninuc_cnt, qinuc_cnt); ...

In other words, I think the added complexity and branching hinder readability and therefore maintainability, so I would recommend simplifying as much as possible

bartgol · 2025-01-29T20:26:47Z

components/eamxx/src/physics/p3/impl/p3_CNT_couple_impl.hpp

+{
+  constexpr Scalar pi       = C::Pi;
+  constexpr Scalar rho_h2o  = C::RHO_H2O;
+  constexpr Scalar qsmall   = 1.0e-18; // BAD_CONSTANT!


I don't think we should merge PRs that use constant marked as "BAD_CONSTANT!" ... If it's bad, let's change it before we merge,or else you know it's never going away...

There is Constants<S>::QSMALL. Does that fit your needs?

It is not "BAD" in the sense that it is wrong. We use the label BAD_CONSTANT as a reminder that it is not a good practice to have arbitrary values for these constants. We must revisit these and understand their impact.

For now, we intend to keep the same values as EAM for these small constants. In my previous experience, changing these constants can sometimes have a large impact on model results (even model climate!). It is better to evaluate the impact before making any change. I checked QSMALL is 1e-14. 1e-18 also has been used in P3, but it is used as a bare constant (numeric literal). Should we keep the BAD_CONSTANT label here or remove it?

Could you add instead a string "TODO: ..." ?

Yes, we can change it to a TODO.

bartgol · 2025-01-29T20:30:52Z

components/eamxx/src/physics/p3/impl/p3_CNT_couple_impl.hpp

+  const auto mask = qc_incld > qsmall;
+  switch (Iflag) {
+    case 1:  // cloud droplet immersion freezing
+      ncheti_cnt.set(mask, frzimm*1.0e6/rho /* frzimm input is in [#/cm3] */ , Zero);


In all these "set" calls, how often do you expect the mask to be true/false? If the mask could often be ALL false (not sometimes, often), then you may consider using if statements, to avoid computing the packs for the true case for nothing (e.g., in the 1st line we have to compute frzimm*1e6/rho regardless of whether we need it or not).

Note: this nano-opt makes sense only if you expect mask to be often false. I assume that's not the case, since qsmall is very small. But I don't know how qc_incld is computed, so maybe it's often 0?

That is a good question. I am not sure about that. @kaizhangpnl or @AaronDonahue might know if mask can often be false or not.

components/eamxx/src/physics/p3/tests/p3_main_unit_tests.cpp

bartgol · 2025-01-29T22:36:34Z

components/eamxx/src/physics/p3/tests/p3_nc_conservation_tests.cpp

@@ -57,7 +57,7 @@ struct UnitWrap::UnitTest<D>::TestNcConservation : public UnitWrap::UnitTest<D>:
        nc_selfcollect_tend[s] = cxx_device(vs).nc_selfcollect_tend;
      }

-      Functions::nc_conservation(nc, nc_selfcollect_tend, cxx_device(offset).dt, nc_collect_tend, nc2ni_immers_freeze_tend, nc_accret_tend, nc2nr_autoconv_tend);
+      //Functions::nc_conservation(nc, nc_selfcollect_tend, cxx_device(offset).dt, nc_collect_tend, nc2ni_immers_freeze_tend, nc_accret_tend, nc2nr_autoconv_tend);


Why is this commented?

Same question for all the other unit tests.

I can answer this for you since my recent experience with this code: tests are commented because they're all broken by the changes in the source tree under physics/p3

Thanks, @mahf708 and Luca for your reviews. I commented out these tests because the function signatures changed when I added this feature, and they failed. @overfelt is already working on reviving these tests (see my comment above about TODOs). He has revived some of the tests already. We will uncomment all these tests and revive them again before changing the status of this PR from "draft" to normal.

bartgol

I have a few comments. Mostly: why are lots of unit tests now commented?

mahf708 · 2025-02-02T05:48:07Z

components/eamxx/src/physics/p3/impl/p3_main_impl.hpp

@@ -243,6 +248,7 @@ ::p3_main_internal(

    p3_main_part2(
      team, nk_pack, runtime_options.max_total_ni, infrastructure.predictNc, infrastructure.prescribedCCN, infrastructure.dt, inv_dt,
+      runtime_options.use_hetfrz_classnuc, ohetfrz_immersion_nucleation_tend, ohetfrz_contact_nucleation_tend, ohetfrz_deposition_nucleation_tend,


Okay, better... but I would do this differently:

I wouldn't pass runtime_options.use_hetfrz_classnuc directly (note that runtime_options is passed in the call as the last arg).

Is there any chance we could make this name better? Imho, "hetfrz" is pretty hard to guess, even for me, so I recommend making the user-facing (namelist-facing) parameters as self-documenting as possible, and I wouldn't worry about the length of the name. Heterogenous, freezing, and nucleation should all appear in the name fully (I'd consider adding mam as well in the name if desired)

I would push the variables lower down as well (later args) but this is not as important

Good catch! Thanks. I didn't notice the runtime_options struct is sent as an arg. I will revisit this to fix it. I will change the variable names as well.

mahf708 · 2025-02-02T05:50:39Z

components/eamxx/src/physics/p3/impl/p3_CNT_couple_impl.hpp

Let's rename this to something better-documenting. If you look at all other names, you see none has ambiguous acronyms in the name. When I see "CNT" I actually first think about carbon nanotubes (my training is chemical engineering) ... I don't think this is about those CNTs, but maybe it is? 🤯

Consider "ice_classical_nucleation"

components/eamxx/src/physics/p3/impl/p3_main_impl_part2.hpp

mahf708 · 2025-02-02T06:01:19Z

components/eamxx/src/physics/p3/p3_functions.hpp

@@ -134,6 +134,7 @@ struct Functions
    bool set_cld_frac_l_to_one = false;
    bool set_cld_frac_i_to_one = false;
    bool set_cld_frac_r_to_one = false;
+    bool use_hetfrz_classnuc   = false; // Set to true to use het frz from MAM4xx-ACI


Comment doesn't belong here. I wouldn't want a reader to get confused and edit stuff here, especially that this is just setting a default that gets overwritten with a call to the following method (see two lines below). The comment belongs in the xml file, which you already have it written well:

<use_hetfrz_classnuc type="logical" doc="Switch to turn on heterogeneous freezing due to prognostic aerosols">false</use_hetfrz_classnuc>

This struct shouldn't have any comments in it, because it should be 1) pretty self-evident and 2) not really editable.

Btw, the first thing we do in initalize_impl is to call the method two lines below:

void P3Microphysics::initialize_impl (const RunType /* run_type */) { // Gather runtime options from file runtime_options.load_runtime_options_from_file(m_params); ...

mahf708 · 2025-02-02T18:17:11Z

Requesting reviews from @hassanbeydoun and @brhillman because I know they're very curious about and interested in this part of the p3 code

mahf708 · 2025-02-03T19:41:14Z

components/eamxx/src/physics/p3/eamxx_p3_process_interface.cpp

+  add_field<Required>("hetfrz_immersion_nucleation_tend", scalar3d_layout_mid,
+                      frz_unit, grid_name, ps);
+
+  // heterogeneous freezing by contact nucleation [cm^-3 s^-1]
+  add_field<Required>("hetfrz_contact_nucleation_tend", scalar3d_layout_mid, frz_unit,
+                      grid_name, ps);
+
+  // heterogeneous freezing by deposition nucleation [cm^-3 s^-1]
+  add_field<Required>("hetfrz_deposition_nucleation_tend", scalar3d_layout_mid,
+                      frz_unit, grid_name, ps);


@singhbalwinder this is the last potential important comment, otherwise, I think the PR is mergeable:

What happens if MAM is inactive? What happens to these variables? Maybe we should make these hidden behind if-else? IF so, make sure to do the same below (around 333):

// Inputs for the heteogeneous freezing diag_inputs.hetfrz_immersion_nucleation_tend = get_field_in("hetfrz_immersion_nucleation_tend").get_view<const Pack**>(); diag_inputs.hetfrz_contact_nucleation_tend = get_field_in("hetfrz_contact_nucleation_tend").get_view<const Pack**>(); diag_inputs.hetfrz_deposition_nucleation_tend = get_field_in("hetfrz_deposition_nucleation_tend").get_view<const Pack**>();

where in the unused/uneeded case, you can use the buffer.unused or something like that

When MAM is inactive, values for these variables will be picked up from the namelist and will stay the same for the entire simulation. The simulation will run as usual with no impact from these variables.

I can hide it behind an if/else if it is not the desired behavior. Is there a way for P3 to know if MAM4 is active or not? Previously, when I discussed this, the design principle did not allow this, as each parameterization should be able to run independently without the knowledge of other parameterizations/processes. Otherwise, it makes the logic complex (e.g., to decide if different combinations of processes are valid or not). If it has changed, please let me know, and I will add if/else blocks.

If they are defined, then there's no need to worry (I wasn't aware they would be defined; this is a little surprising to me tbh)

For hiding them, I would use the boolean you're adding here (hetfrz). It doesn't matter though, whether these statements are hidden or not, won't make a difference as far as I could tell. I was worried a CIME case would error out if these variables are here, but MAM is inactive. I guess we will pick this up in testing soon enough :)

singhbalwinder changed the title ~~Adds aerosols heterogeneous freezing calculations in P3 microphysics~~ EAMxx: Adds aerosols heterogeneous freezing calculations in P3 microphysics Jan 28, 2025

singhbalwinder requested review from mahf708, bartgol and kaizhangpnl January 29, 2025 19:01

mahf708 requested review from AaronDonahue and hassanbeydoun January 29, 2025 19:09

bartgol reviewed Jan 29, 2025

View reviewed changes

mahf708 reviewed Feb 2, 2025

View reviewed changes

components/eamxx/src/physics/p3/impl/p3_main_impl_part2.hpp Show resolved Hide resolved

mahf708 reviewed Feb 2, 2025

View reviewed changes

Adds heterogeneous freezing from aerosols in P3 microphysics

ede97b6

mahf708 requested a review from brhillman February 2, 2025 18:16

Passes all tests on Compy

9be599e

singhbalwinder force-pushed the jroverf/singhbalwinder/eamxx/add-het-frz-p3_rebase1_1 branch from 487f9b6 to 9be599e Compare February 3, 2025 19:25

mahf708 reviewed Feb 3, 2025

View reviewed changes

Turn of this feature by default, all P3 tests pass on PM-GPUs

3a76876

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

EAMxx: Adds aerosols heterogeneous freezing calculations in P3 microphysics #6947

EAMxx: Adds aerosols heterogeneous freezing calculations in P3 microphysics #6947

singhbalwinder commented Jan 25, 2025 •

edited

Loading

github-actions bot commented Jan 25, 2025 •

edited

Loading

Built to branch `gh-pages` at 2025-02-04 02:10 UTC.
Preview will be ready when the GitHub Pages deployment is complete.

singhbalwinder commented Jan 25, 2025

mahf708 commented Jan 29, 2025 •

edited

Loading

bartgol left a comment

bartgol Jan 29, 2025

singhbalwinder Jan 30, 2025

bartgol Jan 29, 2025

mahf708 Jan 29, 2025

singhbalwinder Jan 30, 2025

bartgol Jan 29, 2025

singhbalwinder Jan 30, 2025

mahf708 Feb 2, 2025

bartgol Jan 29, 2025

bartgol Jan 29, 2025

singhbalwinder Jan 30, 2025

AaronDonahue Jan 30, 2025

singhbalwinder Jan 30, 2025

bartgol Jan 29, 2025

singhbalwinder Jan 30, 2025

bartgol Jan 29, 2025

bartgol Jan 29, 2025

mahf708 Jan 29, 2025

singhbalwinder Jan 30, 2025

bartgol left a comment

mahf708 Feb 2, 2025

singhbalwinder Feb 2, 2025

mahf708 Feb 2, 2025

mahf708 Feb 2, 2025 •

edited

Loading

mahf708 commented Feb 2, 2025

mahf708 Feb 3, 2025

singhbalwinder Feb 3, 2025

mahf708 Feb 3, 2025 •

edited

Loading

EAMxx: Adds aerosols heterogeneous freezing calculations in P3 microphysics #6947

Are you sure you want to change the base?

EAMxx: Adds aerosols heterogeneous freezing calculations in P3 microphysics #6947

Conversation

singhbalwinder commented Jan 25, 2025 • edited Loading

github-actions bot commented Jan 25, 2025 • edited Loading

Built to branch gh-pages at 2025-02-04 02:10 UTC. Preview will be ready when the GitHub Pages deployment is complete.

singhbalwinder commented Jan 25, 2025

mahf708 commented Jan 29, 2025 • edited Loading

bartgol left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

bartgol left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mahf708 Feb 2, 2025 • edited Loading

Choose a reason for hiding this comment

mahf708 commented Feb 2, 2025

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mahf708 Feb 3, 2025 • edited Loading

Choose a reason for hiding this comment

singhbalwinder commented Jan 25, 2025 •

edited

Loading

github-actions bot commented Jan 25, 2025 •

edited

Loading

Built to branch `gh-pages` at 2025-02-04 02:10 UTC.
Preview will be ready when the GitHub Pages deployment is complete.

mahf708 commented Jan 29, 2025 •

edited

Loading

mahf708 Feb 2, 2025 •

edited

Loading

mahf708 Feb 3, 2025 •

edited

Loading