Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add NSSL 2/3-moment cloud microphysics #1113

Open
wants to merge 34 commits into
base: develop
Choose a base branch
from

Conversation

MicroTed
Copy link

@MicroTed MicroTed commented Aug 23, 2023

Add NSSL 2/3-moment cloud microphysics

Adds option to use the NSSL 2-moment cloud microphysics scheme with 6 hydrometeor species (droplets, rain, cloud ice, snow, graupel, and hail). It includes the option for 3-moment rain/graupel/hail (adds reflectivity moment) and an option to turn off hail (for non-severe convection). Details on using the scheme are in docs/README.NSSLmp. Also adds a 3D reflectivity field for optional output (refl10cm). The parameterization code is (or will be) at parity with the current WRF PR for the 3-moment scheme wrf-model/WRF#1876 (and eventual CCPP update).

Tested on daily 3-km regional CONUS forecasts (with v7.3), so similar performance is expected with v8.0

Outstanding issue: At least in the physics_mmm directory, the .mod files are not updated by recompiling newer .F files. They have to be deleted in order to be replaced. I don't know if this is general or a gfortran quirk.

…icrophysics (adds refl. moment variables for rain, graupel, hail: zrw, zgw, zhw).

 Also added code to compute the composite and 1-km reflectivity for NSSL2M (not actually tested yet)
  mpas_io_streams.F : pass ioformat to def_var
  mpas_io.F : add conditional to set deflate on variables if ioformat is netcdf4 (needs netcdf-4.7.4 or later)
…ophysics

ccpp_kinds.F : Add 'machine' module for CCPP compatibility
mpas_atmphys_driver_microphysics.F : Add support for NSSL microphysics (arrays, init/driver etc.)
mpas_atmphys_interface.F : Get/put data for NSSL microphysics
mpas_atmphys_packages.F : Set up packages for NSSL microphysics
physics_mmm/Makefile : add dependencies for NSSL-mp
physics_wrf/Makefile : remove NSSL-mp (moved to physics_mmm)
docs/README.NSSLmp : brief overview of options in the NSSL microphysics scheme
physics_mmm: module_mp_nssl_2mom.F, mp_nssl.F : new code file for NSSL-MP
@MicroTed
Copy link
Author

Hi, @ldfowler58 would you have time to look at this at some point? Let me know if this is the approach you want for physics_mmm or if there is more to shift from the model code to the interface. The physics module is an updated version from CCPP, so it has all the extra comments for documentation. The interface 'run' subroutine is from CCPP (for SCM/FV3), but isn't actually used since MPAS can call the microphysics directly with the "_p" arrays.
Thanks!

@mgduda
Copy link
Contributor

mgduda commented Aug 31, 2023

I think it would be good to separate the changes for NetCDF compression from the changes to add the NSSL microphysics scheme, as these (if I understand correctly) are unrelated changes.

@MicroTed
Copy link
Author

I agree that it would be cleaner to split it off the netcdf part, and it is only 2 files and I think no overlapping changes. It's a bare-bones implementation that assumes one would only choose the netcdf4 option because of compression being available. Otherwise the speed of pnetcdf for uncompressed writes really can't be beat.

@MicroTed MicroTed changed the title Add NSSL 2/3-moment cloud microphysics AND enable compression for netcdf-4.7.4 or greater Add NSSL 2/3-moment cloud microphysics Aug 31, 2023
@MicroTed
Copy link
Author

I moved the netcdf changes to a separate branch: https://github.com/MicroTed/MPAS-Model/tree/dev-nc4deflate
I could make a draft PR or just leave it on ice for now -- whatever you want to do given your priorities and limited time. Ideally I would add a netcdf version check and an extra flag for deflate level etc.

@ldfowler58
Copy link
Contributor

I started to look at the changes and add-ons made by Ted in order to include the NSSL cloud microphysics scheme. Here my first suggestions for the top files that I have been modified:

  1. In ./src/core_atmosphere/Registry.xml, I suggest to keep the "mp_" reference in front of the packages named nssl_hail_in, nssl_ccn_in, nssl3m_in, and nssl3m_hail_in. Users may appreciate that those packages point to cloud microphysics options. nssl3m_in can be changed to mp_nssl3m_in since you defined the NSSL 2-moment package as mp_nssl2m_in. Thanks.

  2. In ./src/core_atmosphere/physics/mpas_atmphys_packages.F, can you please keep everything "lower case" and remove extra blank spaces to ensure "continuity" in how the sourcecode has been originally written (if not much to ask). For instance, change "
    if(config_nssl_hail_on) then
    nssl_hail_in = .true.
    endif

  3. In ./src/core_atmosphere/physics/mpas_atmphys_control.F, change lines
    config_microp_scheme .eq. 'mp_thompson' .or. &
    config_microp_scheme .eq. 'mp_nssl2m' .or. &
    to
    config_microp_scheme .eq. 'mp_nssl2m' .or. &
    config_microp_scheme .eq. 'mp_thompson' .or. &

  4. In ./src/core_atmosphere/physics/mpas_atmphys_vars.F (lines 197-203), do you mind putting one variable per line as for lines 188-195 so that I can add the definition of the different variables later. Thanks.

  5. In ./src/core_atmosphere/physics/mpas_atmphys_manager.F, same comments as in 2) regarding lower case and space. Wouldn't it be simpler to use 75._RKIND instead of 1.25001*60._RKIND? Thanks.

Will look at the actual microphysics next. Thanks.

@MicroTed
Copy link
Author

MicroTed commented Jan 3, 2024

Thanks for the comments! My thinking for leaving off "mp_" from nssl3m_in and nssl_hail_in etc. is that they are package options rather than full packages. They are not seen by the user, but rather activated by config_ options. So maybe it doesn't really matter? The package is still the mp_nssl2m_in, and the others are options (subpackages?) to turn variables on or off.

Perhaps I could just rename mp_nssl2m_in to just mp_nssl_in? Or something like that? I also have a version with an option to run single-moment, but I'm not sure that's needed for a release code.

I can change the other formatting things to match the rest of the code. (Although I'm firmly in the camp of having a space after 'if' :) considering there is an unbalanced space in front of 'then', haha.)

  Registry.xml : fix ntiedtke vs. tiedtke
  physics/ccpp_kind_types.F : remove 'machine'
  physics/mpas_atmphys_interface.F : fix typo and error with refl10cm
  physics/physics_wrf/Makefile : add NSSL mp routines (for now)
  physics/physics_wrf/mp_nssl.F : switched module for kind_phys
@MicroTed
Copy link
Author

Updated (merged) to 8.2.1-hotfix. Also added function in init_atmosphere to generate LBCs for hydrometeor number concentrations when the LBC data only have mixing ratio. I still need to run some more tests and update the NSSL microphysics code. (Also some enhancements to the ideal initialization will be added.)

…l reordering of scalars

  mpas_atmphys_interface.F : For mp_nssl, changes to avoid problems with PG compilers and pointers
  module_mp_nssl_2mom.F : Updates to recent development version
     - Use more accurate form for water saturation
     - Add active breakup parameterization for 3-moment rain
     - Fixes problem with sedimentation when nz > 128
     - Cleans up some unused/unset rate arrays that could cause issues in debug compile
@MicroTed
Copy link
Author

I ran into an odd compiler issue with nvfortran/pgf90 at -O3 with mpas_atmphys_interface.F: If a pointer assignment is inside an IF structure, it may not work correctly. For example, the optional hail variables here:

             IF ( config_nssl_hail_on ) THEN
              call mpas_pool_get_dimension(state,'index_qh'  ,index_qh  )
              call mpas_pool_get_dimension(state,'index_nh'  ,index_nh  )
              call mpas_pool_get_dimension(state,'index_volh' ,index_volh  )
               qh   => scalars(index_qh,:,:)
               nh   => scalars(index_nh,:,:)
               volh => scalars(index_volh,:,:)
             ENDIF

At -O0 it works fine, but -O2 or 3, the pointer for qh seems to go randomly somewhere, maybe within the scalar structure or maybe not. I have only looked at the min/max values of qh. Obviously this is bad compiler behavior, but I found two work-arounds for it:

  1. Remove the IF clause and just always assign qh etc. even if the index is -1, but still have an IF clause around accessing the array. I don't know if this will cause debug problems, however (e.g., bounds check, although the array isn't being accessed)
  2. Just replace
    qh_p(i,k,j) = qh(k,i)
    with
    qh_p(i,k,j) = scalars(index_qh,k,i)

I currently have option 2 in the PR code, but can switch it to 1 or a better suggestion.

Also for this reason I have not yet included mp_nssl with the mass/number tendencies (e.g., rqcmpten). Those seem to be diagnostic only, but maybe there will be a need for those.

@ldfowler58
Copy link
Contributor

Hi Ted: Thanks for pointing this out. I do not use nvfortran or pgf90 at all but it is definitely something that we should look into. Is it a "feature" that only occurs with the NSSL cloud microphysics or a "feature" that could also occur with, for instance, the Thompson cloud microphysics scheme? I'll investigate with a different cloud microphysics scheme separately.
Laura

@MicroTed
Copy link
Author

Hi, Laura,
It seems to only be a compiler bug with the IF-ENDIF sections, not with the SELECT-CASE at all. So I think the mp-thompson sections are OK. I did test the regular mp_thompson and it ran fine. We have some AMD epyc computers locally, and the nvfortran/pgf runs pretty well for them. I haven't tried gfortran yet for comparison, but will try to do that.

That said, I don't yet fully understand the structure of the 'scalar' array and how it works to reference the first index. (Basically how it converts a s(numscalar)%vals2d(nlevs,ncells) kind of record into something that looks like a 3D array -- too fancy!)

@weiwangncar
Copy link
Contributor

@MicroTed Apparently not updating .mod file is a gfortran feature as I just learnt. See https://gcc.gnu.org/bugzilla/show_bug.cgi?id=31587

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants