MetOffice / MetOffice/lfric_apps
CoMorph further KGO-preserving refactoring
@MichaelWhitall is already working on this.
Since Aug 13, 2026.
- Dominant language
- Fortran
- Stars
- 31
- Forks
- 118
- Avg merge
- 2d 4h
- Merged PRs (30d)
- 8
Description
### Are there any linked Issues or Pull Requests?
_No response_
### Brief description
The bulk of the refactoring of CoMorph in the comorph_dev branch was lodged in the lfric_apps trunk in Issue #178 / PR #292. But some more KGO-preserving refactoring has been done in the comorph_dev branch after the branch lodged for #292 was split-off from it. We now need to lodge these extra refactoring changes. We also need to add more of CoMorph's tune-able constants to the LFRic namelist.
lfric_apps branch:
[vn3.2 -> vn3.2_comorph_refact2](https://github.com/MetOffice/lfric_apps/compare/main...MichaelWhitall:lfric_apps:vn3.2_comorph_refact2)
Equivalent changes in the UM:
[vn14.2_comorph_dev_lilf -> vn14.2_comorph_refact2](https://github.com/MichaelWhitall/um/compare/vn14.2_comorph_dev_lilf..MichaelWhitall:um:vn14.2_comorph_refact2?w=1)
Note the base-line code for this diff is the comorph_dev "lilf" (lodged in LFRic) branch, which packages together all comorph-related changes which have been lodged already in LFRic but are not on the UM trunk.
### Further details of the issue.
**Main design changes:**
1. We port the comorph namelist from the UM into lfric_apps. So-far the only comorph run-time settings that can be altered via the namelist in LFRic were the few that differ between CoMA9 versus CoMA9_TB ("par_gen_mass_fac", "par_gen_rhpert", "par_radius_ppn_max", "resdep_precip_ramp", "dx_ref"). These were placed in the "convection" namelist block (confusingly mixed in with the 6A convection scheme settings, but trigger-ignored if not using comorph).
We now create a new "comorph" namelist block, defined in `interfaces/physics_schemes_interface/rose-meta/um-convection/HEAD/rose-meta.conf`, containing all the same sections and run-time settings that currently exist in the UM trunk comorph namelist. As in the UM, we set the "comorph" namelist block to be trigger-ignored if not using comorph (much neater to trigger-ignore the whole section instead of having to add triggering for each namelist entry individually). The upgrade macro added in `interfaces/physics_schemes_interface/rose-meta/um-convection/versions.py` does 3 things:
* Add the new "comorph" namelist block to existing apps.
* Migrate the existing CoMA9_TB-related namelist settings from the "convection" namelist to the new "comorph" namelist.
* Add all the new "comorph" namelist entries, set to the values they were previously hardwired to in `interfaces/physics_schemes_interface/source/support/um_physics_init_mod.f90`.
We then modify `um_physics_init_mod` to set the values of the run-time settings in `comorph_um_namelist_mod` using the values from the LFRic namelist instead of using hardwired defaults.
Note that lfric_apps contains a number of "example configuration" namelist (.nml) files used in various test-runs; these have to be updated by-hand as the upgrade macro is not applied to them. I noticed that the settings "dx_ref" and "resdep_precip_ramp" had been added to all of these, but weren't actually used in any of them as none of the example configurations use comorph. I have simply deleted them from the example configs. This is consistent with the fact these settings now reside in the "comorph" namelist block which is wholly trigger-ignored (since comorph is off), which would suppress any part of the "comorph" namelist being added to the actual namelist files by rose. For the same reason, I do not need to add the new "comorph" namelist block to any of the example configuration files in this PR.
Note that all the namelist entries added in `interfaces/physics_schemes_interface/rose-meta/um-convection/HEAD/rose-meta.conf` have `!kind=default` except for `par_radius_knob` which has `!kind=double`. This tells PsyClone to define just that one variable in double-precision in `comorph_config_mod` (the rest default to single precision). I found that, just for this one namelist item, reading it from the namelist in single-precision and then converting to double when copying to `comorph_um_namelist_mod` in `um_physics_init` slightly changed its value compared to the original code (which hardwired it to 0.45 in native double precision). Setting `!kind=double` for `par_radius_knob` avoids the precision conversion and so allows the branch for this PR to preserve KGO. It is otherwise completely unnecessary (change in answers is only at rounding-error level) and can be changed to `!kind=default` in another PR.
1. CoMorph's optional run-time checking for bad values of fields (subroutine `check_bad_values`) is upgraded to check whether a field is within plausible bounds (set for each specific variable). Previously it only checked whether fields were not NaN / Inf, or negative for variables that shouldn't be. This change was needed during earlier debugging work, to track down bugs which create values which are physically implausible but still representable in a 32-bit floating-point real. This change entails:
* `check_bad_values` itself is modified to take optional inputs for the acceptable min and max bounds `field_min`, `field_max` for the field being checked (if the optional inputs are not present, the bounds are set to +/- `HUGE()`, which will only trap NaN / Inf as before). These inputs replace the input logical `l_positive`, which indicated whether negative values should be flagged (setting the input lower bound to zero will have the same effect as setting `l_positive` to true in the old version).
* Every call to `check_bad_values` is updated to pass in plausible bounds for the field being checked (and instances where `l_positive` was passed in set to `.true.` are replaced with the lower bound passed in set to zero).
* Various comorph modules which handle the properties of fields within comorph (`fields_type_mod`, `parcel_type_mod`, `convcloud_type_mod`) have been updated to set the allowed lower and upper bounds for each of their fields. The lists of bounds (e.g. `field_min`, `field_max` for the fields in `fields_type_mod`) are then imported from the module to be passed into each call to `check_bad_values`. In the case of `fields_type_mod`, we had to set separate lists of lower and upper bounds `field_min_cons`, `field_max_cons` to apply when the fields have been converted to "conserved variable form" (e.g. when "temperature" includes factor of heat-capacity cp). In the case of `parcel_type_mod` we have added a new list containing the name of each parcel field, and the module's contained subroutine `parcel_check_bad_values` now uses this to embed the calls to `check_bad_values` within a loop over the parcel fields instead of having a separate call hardwired for each of "massflux_d", "radius", "edge_virt_temp".
* The plausible bounds for primary fields in `field_min`, `field_max` are used in the bad value checking within comorph's microphysics, which creates a dependency of the microphysics compile on `fields_type_mod`. This required adding this module to the compile in some of the comorph unit tests.
* The max possible representable floating-point real `max_float` is now declared and stored in `comorph_constants_mod`, so that `fields_type_mod` can use it to set the allowed bounds on unknown tracer fields in the list of primary fields.
1. Some existing hardwired parameters from `comorph_constants_mod` are changed to run-time variables and set via the namelist (see changes in the rose meta-data / upgrade macro, `comorph_um_namelist_mod`, `set_constants_from_um` to add the new namelist variables). The following entries are added to the namelist:
* `cf_area_coef`
* `min_cmr`, `max_cmr`
* `nconc_cl`/`rain`/`cf`/`snow`/`graup`
(see the new `rose-meta.conf` entries for details of what these are).
1. In `comorph_constants_mod`, the tuning factor `fac_tdep_n` (which sets the slope of the T-dependent ice number concentration within the comorph microphysics) has been moved to inside the `cond_params` derived-type structure, so that separate values can be set for each hydrometeor species. This replaces the `l_tdep_n` switch within the `cond_params` structure (`fac_tdep_n` can just be set to zero to disable T-dependent n for a given species). Then in the UM namelist, we replace `r_fac_tdep_n` (which was only used for ice-cloud) with separate T-dependent n e-folding temperature scales for ice-cloud and liquid-cloud, `tdep_n_cf` and `tdep_n_cl`. The net result is to add the functionality to allow the user to set a T-dependent liquid-cloud number (the upgrade macro sets this off by default to preserve existing science). Then rose meta-data / upgrade macro, `comorph_um_namelist_mod`, `set_constants_from_um` have been updated to handle this change to the namelist, and `calc_cond_properties` has been updated to use the T-dependence factors from `cond_params`.
1. In `parcel_type_mod`, subroutine `parcel_combine` (which computes the updated combined parcel properties where initiating parcels from each model-level are merged into the bulk plume), we have moved the calculations for the parcel core and edge properties inside their own subroutine, `core_combine`. Here, we have left the calculations exactly the same in their new location. But this rearrangement will facilitate future changes to:
* Fix some bugs in the core and edge combine calculations (see Issue #252).
* Use the same core and edge combine calculations elsewhere (see Issue #252).
* Add new options for how to combine parcel core properties...
1. We have moved the setting of the entrained air properties "ent_fields" from `conv_level_step` to inside the entrainment subroutine `set_ent`. Also the copying of "core_ent_ratio" into a diagnostics array for output is moved from `conv_level_step` to `set_ent`. Later PRs will add new functionality in the entrainment routine with new diagnostics to be copied alongside "core_ent_ratio".
1. Subroutine `set_par_cloudfrac` (which sets the in-parcel cloud-fraction fields by checking whether the parcel contains liquid or ice cloud-condensate) has been refactored to account for the possibility of multiple condensate species contributing to the liquid or ice cloud. It was hardwired to just check the species "q_cl" (for the liquid-cloud-fraction) and "q_cf" (for the ice-cloud-fraction). However, under future changes for CASIM coupling we will have 2 ice species contributing to the ice-cloud-fraction ("q_cf" and the not-yet-used "q_snow"). The existing code would wrongly set the ice-cloud fraction of the parcel (and hence of the detrained air) to zero if it contained "snow" but not "ice". This could create spurious instances of snow present with no fraction in the environment once the air is detrained.
To fix this, we now loop over condensate species and account for all that have "i_sg" set to "i_sg_frac_liq" for the liquid-cloud, and all that have "i_sg" set to "i_sg_frac_ice" for the ice-cloud (in `comorph_constants_mod` these are already set to put only "q_cl" in the liquid-cloud-fraction and "q_cf" and "q_snow" in the ice-cloud-fraction). All calls to `set_par_cloudfrac` are then altered to pass in the full list of condensate masses instead of just "q_cl" and "q_cf". Also the code to calculate the convective cloud-water (which gets output to radiation in the "ccw" array) in `set_diag_conv_cloud_a` is consistently modified to loop over all condensate fields that contribute to the liquid or ice cloud-fraction.
1. A unit test is added for the in-plume microphysics implicit phase-change solver (new files `test_phase_change_solve.F90` and its compile script `build_test_phase_change_solve` added in the `unit_tests` directory). This unit test was added to the comorph_dev branch to help debug a failure in the phase-change solver; print statements were used to extract all the inputs to `phase_change_solve` from the failing model-run at the grid-point / time when the problem occurred; these input values were then set in the unit test and successfully reproduced the failure there. The unit test can then easily be run in an interactive debugger to track down and fix the problem :)
**Miscellaneous minor changes**
1. In `comorph_constants_mod`, we correct the default setting for the sub-grid fraction in-which graupel is assumed to reside by changing `params_graup % i_sg` from `i_sg_homog` to `i_sg_frac_prec` (previously assumed homogeneous over the grid-box, now confined within the sub-grid precip fraction). Note this only affects the comorph stand-alone test; in the UM / LFRic this setting is overridden in `set_constants_from_um` based on host-model switches so there is no impact on model-runs.
1. The calculation of the mass-flux weighting to apply to comorph diagnostics which are averaged over multiple convection types / layers has been tweaked in `diags_super_type_mod`. This change avoids an "equality comparison for real32" warning when compiling the code with rigorous settings. The original code tested whether the current grid-point / conv type / conv layer has zero mass-flux in types / layers added-on so-far by testing whether the updated mass-flux equals the current type / layer's mass-flux at L805; this generates a warning. This has been refactored so that instead we test whether the mass-flux added so-far is > zero, and then put the special treatment of when it is zero in an `ELSE` block. This avoids the compiler warning and arguably makes the code more intuitive.
1. In `comorph_constants_mod`, the number of updraft types (`n_updraft_types`) has the parameter attribute removed. This makes it consistent with `n_dndraft_types` which was already run-time variable on the trunk. This change has meant that `n_updraft_types` needed to be added to OMP declarations in some places.
1. In `comorph_diags_um_mod`, now initialising the diagnostic pointer for `pressure_incr_env` to `NULL()` unconditionally (existing code only did this if using PC2; new code should be safer and more succinct). **Only applicable in the UM branch, not LFRic.**
1. The initialisation of outputs from `calc_cond_properties_cmpr` to zero in `microphysics_1` has been moved inside the subroutine (in `calc_cond_properties`). This simplifies the code a bit, removing the need for an IF test on whether any points have nonzero condensate in `microphysics_1`, and allowing the outputs from `calc_cond_properties_cmpr` to actually be declared with `intent(out)`.
1. In `ice_nucleation` and `ice_rain_to_graupel`, corrected dodgy comparison of real values with integer 0 in IF tests (now uses real-valued `zero` instead of just `0`).
1. In `fall_in`, the added reset of the hydrometeor fall-fluxes to zero is actually a bug-fix (the surrounding code in `moist_proc` assumes that `fall_in` will do this). This change has no effect on existing configs because this only affects fall-in of hydrometeor from the environment into the parcel, and that is hardwired off in `parcel_dyn` (pending implementation of the fall-in flux calculation, `flux_cond` is just set to zero before input to `moist_proc`).
1. In `moist_proc`, removed confusing reuse of the array `rho_wet` to store total-water mixing-ratio (new local array "q_tot" added for that). This then allowed moving the calculation of dry-density "rho_dry" and wet-density "rho_wet" slightly earlier in this routine (straight after the calculation of "q_tot", to put it all in one place).
1. The utility routine which calculates total-water mixing-ratio on uncompressed / host-model precision data has been changed from 3D to 2D in `calc_q_tot` (i.e. it now operates on data for a single model-level at a time instead of the full 3D domain). The call to this routine from the 3D calcaultion in `calc_virt_temp` has been updated accordingly so it now calls `calc_q_tot_2d` inside an OMP parallel loop over model-levels. This change removes the need for a full 3D array to store `q_tot` in `calc_virt_temp_3d`. It also paves the way for a not-yet-included option to include the moisture-dependence in the dry adiabatic lapse rate in `dry_adiabat_2d`.
1. In `calc_env_region_tq_nb`, the condensed water excess of each sub-grid region is precalculated and stored in its own array, instead of repeatedly recalculating it in various calculations. This also paves the way for an option not introduced yet, to ignore condensed water loading in the calculation of neutrally-buoyant sub-grid region T, q, by just setting the new condensed water excess variable to zero.
1. In `set_par_fields`, we remove a never-used option to set the parcel radius to an ad-hoc fixed number `par_gen_radius` when `l_turb_par_gen` is false (it is currently hardwired to true anyway). `par_gen_radius` is then deleted from `comorph_constants_mod`.
1. In `collision_ctl`, we've just renamed the hydrometeor species index for the collected species from "i_liq" to "i_cond1", and the index for the collecint species from "i_cond" to "i_cond2". This hopefully makes it clearer what the loops over hydrometeor species are doing; now "i_cond1" is collected by "i_cond2". This also paves the way for later changes to allow snow to collect ice (when snow is switched on for coupling with CASIM); labelling the collected species "i_liq" would be misleading if one of them is actually ice!
1. Commenting and code-style improved in a few places.
Contributor guide
No contributing guide indexed for this repository
First steps
- Read the whole issue, then the project's contributing guide.
- Comment on the issue to say you are picking it up — it saves two people doing the same work.
- Fork the repository and make your change on a branch.
- Open a pull request that references the issue number.
Assessment
This issue has not been assessed yet.