Extend diag_mediator to allow the piecemeal posting of diagnostics #809
Conversation
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Codecov ReportAttention: Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## dev/gfdl #809 +/- ##
===========================================
Coverage ? 21.81%
===========================================
Files ? 137
Lines ? 33019
Branches ? 5850
===========================================
Hits ? 7202
Misses ? 25260
Partials ? 557 ☔ View full report in Codecov by Sentry. |
Hallberg-NOAA
requested changes
Jan 23, 2025
There was a problem hiding this comment.
I think that this capability should be extended to include the ability to initialize a diagnostic buffer to land values to avoid the need to call the piecemeal posting for all columns, even over land.
Also, this PR should be revised to adhere to the guidance in the MOM6 style guide at https://github.com/NOAA-GFDL/MOM6/wiki/Code-style-guide, as described in specific comments.
|
@Hallberg-NOAA Thank you for the comments. I have added descriptions of the units (see the inline response) and made the if-statements conform the to the style guide. Additionally, the fill value for the arrays can now be specified; for the purposes of the piecemeal posting of the diagnostics these are set to the |
ashao
force-pushed
the
post_piecemeal
branch
2 times, most recently
from
5c5a238 to
a10cb9a
Compare
The MEKE GM computation of `src` and `src_GM`, a diagnostic array, were placed in a single loop. The similar RHS of each expression made it unfavorable to use FMAs on the `src` update. Older production runs depending on this FMA were seeing answer changes. This patch restores the FMA loop update of `src` by separating `src` and `src_GM` into separate loops.
This patch makes several adjustments to MOM_MEKE.F90 and MOM_hor_visc.F90 to ensure that the Laplacian and biharmonic friction coefficients are computed separately, and only if their respective terms are enabled. This resolves some subtle bugs where the default biharmonic value of -1 was applied to the Laplacian case, even when the biharmonic MEKE friction was disabled.
The if-test inside of the FrictWork loops are likely to impede performance. Even if the total work is reduced, they are likely to interrupt pipelines. When EY24_EBT_BS is disabled, they will clearly reduce performance. This patch moves those tests outside of the if-block and applies them separately. (Calculation would be slightly improved if the meaning of the flag were reversed, but I don't want to make additional changes.)
The damping MEKE loop also included updates to multiple diagnostics, even if they were not registered. This would presumably have a negative impact on performance. This patch moves each diagnostic into a separate loop. It also conditionally precomputes the damping and damp_rate parameters, which are now stored as 2d arrays rather than in-loop scalars. As before, the MEKE calculation is left unchanged in order to preserve bit reproducibility.
The redefining of int_tide_CS control struct in set_diffusivity_init caused errors in debug-mode for Intel compilers. The issue appears to be an internal function that expects a pointer rather than the type. This patch reverts this back to a pointer. We can revisit this if there is a need to reduce reliance on pointers.
This patch updates the expression for FrictWork_bh (biharmonic frictional work) when the FrictWork_bug flag is enabled. The new form is symmetric to rotations when FMA instructions are enabled.
Self-attraction and loading calculation in Boussinesq pressure gradient force is refactored to be consistent with the algorithm in non-Boussinesq version.
Using SSH to calculate self-attraction and loading (SAL) is only accurate for barotropic flows. Bottom pressure anomaly should really be used for general purposes. * New runtime parameter A runtime parameter SAL_USE_BPA is added to use bottom pressure anomaly. The option requires an input file for long term mean reference bottom pressure. The reference bottom pressure field is stored with SAL_CS. * Refactor SAL and tides in Boussinesq mode As the total bottom pressure is needed for bottom pressure anomaly, SAL calculation in Boussinesq mode needs to be refactored. In addition, there is a longstanding bug in Boussinesq mode that interface height is modified by SAL and tides, and the modified interface height is erroneously used for density and pressure calculation later on. Therefore the SAL and tides are refactored by moving their calculations after pressure is calculated. Tide answers before 20230630 is retained. Answers after 20230630 are changed for Boussinesq mode.
* Local variable SSH is renamed to pbot_anom. * Tides related local variable names are changed to be less confusing. Notably, `e_sal_tide` is renamed to `e_sal_and_tide` (the summation of sal and tide), not to be confused with `e_tide_sal`, which is renamed to `e_tidal_sal` (sal from tides). * Fix e_tidal diagnostics in Boussinesq mode. The term was not assigned if tide answer date is >20230630.
Add an alternative method for SAL and tides in Boussinesq mode. The current method adjusts interface heights with geopotential height anomaly for SAL and tides. For non-Boussinesq mode, the current method is algebraically the same as taking the gradient of SAL and tide geopotential (body forcing). For Boussinesq mode, there is a baroclinic component of tidal forcing and SAL. The alternative method is added to calculate the gradient of tidal forcing and SAL directly at the cost of additional multiplications. The new method is controlled by runtime parameter BOUSSINESQ_SAL_TIDES.
* Instead of rescaling bottom pressure to height unit, calc_Loving_scaling is modified to be conditionally dimensional. When calculating self-attraction and loading, Love numbers are now dimensional when bottom pressure anomaly is used as an input. This change eliminate Love numbers' dependence on mean seawater density. * A new coefficient called linear_scaling is added to SAL CS for the same purpose, although to use bottom pressure anomaly for scalar approximation is not quite justifiable. A WARNING is given when users try to do that. * Two separate field, SSH (sea surface height anomaly) and pbot (total bottom pressure), are used for calculating self attraction and loading when SAL_USE_BPA = False and SAL_USE_BPA = True, respectively. The change is to enhance code readability.
A new date for TIDES_ANSWER_DATE is added to recover answers for tides in Boussinesq mode before 20250201.
* Precalcualte a local field `coef_rhoE` to avoid in-loop division and if-blocks. The unit of coef_rhoE depends on use_bpa. * Fix a few unit description typos in SAL module and two other files.
Refactored MOM_opacity to replace hard-coded dimensional parameters for the Manizza and Morel opacity fits with run-time parameters, and also added the runtime parameter OPACITY_BAND_WAVELENGTHS to provide the ability to set the wavelengths of the bands, even though these are not actually used in MOM6. By default, these parameters are all set to the previous hard-coded values, using the recently added defaults argument to get_param_real_array(). The bounds of the frequency band label arrays with the MANIZZA_05 opacity scheme were also corrected when PEN_SW_NBANDS > 3, but it would not be typical to use so many bands for no purpose and these labeling arrays (optics%min_wavelength_band and optics%max_wavelength_band) do not appear to be used anywhere. In addition, the unused publicly visible routines opacity_manizza and opacity_morel were eliminated or made private. All answers are bitwise identical, but there are new entries in some MOM_parameter_doc files.
…om-ocean#794) * add option to specify horizontally varying decay * extend to vertical modes * fix comments * corrected description of decay_rate array --------- Co-authored-by: Raphael Dussin <Raphael.Dussin@noaa.gov>
* Add ice shelf pressure initialisation bug fix This commit adds a new parameter FRAC_DP_AT_POS_NEGATIVE_P_BUGFIX and associated bug fix. This bug occurs when TRIM_IC_FOR_P_SURF pressure initialisation is used for ice shelf, whilst in Boussinesq mode and a nonlinear equation of state. The subroutine trim_for_ice calls cut_off_column_top. This routine in Boussinesq mode calls find_depth_of_pressure_in_cell in MOM_density_integrals.F90. This subsequently calls the function frac_dp_at_pos which calls the density equation of state based on given salinity, temperature and depth, which is incorrectly converted into pressure by z position (which is negative) multiplied by g and rho0. The bug results in incorrect densities from the equation of state and therefore an imperfect initialisation of layer thicknesses and sea surface height due to the displacement of water by the ice. The bug is fixed by multiplying the pressure by -1. This commit introduces parameter FRAC_DP_AT_POS_NEGATIVE_P_BUGFIX with default value False to preserve answers. If true, the ice shelf initialisation is accurate to a high degree with a nonlinear equation of state. Answers should not change, except for the added parameter in MOM_parameter_doc. The same functions are called by a unit test in MOM_state_initialization; in this commit I set the MOM_state_initialization unit test to use the existing bug (with a false flag). * Resolve indenting and white space inconsitencies with MOM6 style for ice shelf pressure initialisation bug fix FRAC_DP_AT_POS_NEGATIVE_P_BUGFIX
Added the new optional argument old_name to the 8 get_param() routines. This new capability allows for an archaic parameter name to be specified and for appropriate warnings encouraging the user to migrate to using the new name while still setting the parameter as intended, or error messages in the case of inconsistent setting via the archaic name and the correct name. The logging inside of the MOM_parameter_doc files only uses the correct parameter name. Also added the new optional argument set to the 8 read_param routines, to indicate whether a parameter has been found and successfully set. The new set argument is now being used in read_param() calls in obsolete_int(), obsolete_real(), obsolete_char() and obsolete_logical(). Obsolete_logical() in particular was substantially simplified by the use of this new argument, and is now only about half as long as it was. The read_param() set argument is also used in all of the get_param() routines when they are given an old_name argument. The new old_name argument to get_param() is not yet being used in the version of the MOM6 code that is being checked in, but it has been tested extensively by adding or modifying get_param calls in a variant of the initialization code, and it will be used in an updated version of github.com/NOAA-GFDL/pull/725 to gracefully handle the deprecation of 4 parameter names. All answers are bitwise identical, but there are new optional arguments to two widely used interfaces.
Add the new element grid_unit_to_L to the ocean_grid_type and the dyn_horgrid_type, which can be used to convert the units of the geoLat and geoLon variables to rescaled horizontal distance units ([L ~> m]) when they are Cartesian coordinates. When Cartesian coordinates are not in use, G%grid_unit_to_L is set to 0. This new element of the grid type is used to test for inconsistent grids or to eliminate rescaling variables in set_rotation_beta_plane(), initialize_velocity_circular(), DOME_initialize_topography(), DOME_initialize_sponges(), DOME_set_OBC_data(), ISOMIP_initialize_topography(), idealized_hurricane_wind_forcing(), Kelvin_set_OBC_data(), Rossby_front_initialize_velocity(), soliton_initialize_thickness(), and soliton_initialize_velocity(). These are the instances where this new variable could be used and bitwise identical answers are recovered. There are a few other places where they should be used, but where answers would change, and these will be deferred to a subsequent commit. All answers are bitwise identical, but there are new elements in two transparent and widely used types.
…u, v or uh, vh The drifters package is able to run in 2 different modes: one where the particles are advected using u and v, and one where they are advected using uh and vh. When u and v are used (i.e. the drifters are advected using the resolved velocities only, with no sub-grd component), the particles package needs the layer thickness that is consistent with these velocities, so particles_run needs to be called before any subgrid scale parameterizations are applied. This was not implemented correctly in my last PR, and is corrected here. The particles package also needs the correct timestep for particle advection. This is added here.
Remove the unused element Rad_Earth from ocean_grid_type and dyn_horgrid_type. The dimensionally rescaled equivalent element G%Rad_Earth_L is extensively used, and it will continue to exist. G%Rad_Earth_L was introduced in November 27, 2021 as a dimensionally rescaled version of G%Rad_Earth, while the unscaled version was retained because at the time, the Rad_Earth element of the dyn_horgrid_type is also used in SIS2. However, SIS2 has not used G%Rad_Earth since December 23, 2021, so after 3 years we can now safely remove this unused element. Any cases on other branches that might be impacted by this change will not compile. All answers are bitwise identical, but there is one less element in two transparent types.
Some quantities in MOM6 are calculated in subroutines which expect slices of the model's 2d or 3d arrays. Diagnosing these quantities can be challenging because the usual post_data calls expect whole arrays. To solve this problem, the changes here introduce a dynamic buffer that can grow as needed over the course of an simulation. This buffer keeps track of the diagnostic IDs and the index in the buffer. When a slot in the buffer is no longer needed, for example if the whole array has been computed and posted, the slot is marked as "available" for overwriting. The buffer is only allowed to grow if all the currently allocated slots are in use. Any computational cost associated with growing this buffer will only happen in the first few steps of the model as post_data for requested diagnostics is called for the first time in that run.
The diag mediator has been extended to add a dynamic buffer to each axes group. Three new methods have also been added to enable the piecemeal posting (by column, by point) of a diagnostic and a 'final' method to allow the buffer to be reused later.
ePBL calculates the vertical diffusivity column by column. This provides a convenient sanity check of the new piecemeal posting of diagnostics. The original diagnostic Kd_ePBL is done by posting the full 3d prognostic array, whereas a new diagnostic Kd_ePBL_col_by_col posts the same array from within ePBL but does so column by column.
Fixes failures in the CI due to some procedures and type members not having docstrings.
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
Some quantities in MOM6 are calculated in subroutines which expect slices of the model's 2d or 3d arrays. Diagnosing these quantities can be challenging because the usual post_data calls expect whole arrays. To solve this problem, the changes here introduce a dynamic buffer that can grow as needed over the course of an simulation.
This buffer keeps track of the diagnostic IDs and the index in the buffer. When a slot in the buffer is no longer needed, for example if the whole array has been computed and posted, the slot is marked as "available" for overwriting. The buffer is only allowed to grow if all the currently allocated slots are in use.
Any computational cost associated with growing this buffer will only happen in the first few steps of the model as post_data for requested diagnostics is called for the first time in that run.
This new piecemeal posting of a diagnostic was implemented for
Kd_ePBL. The new diagnosticKd_ePBL_col_by_colis the same as the the original field in theBaltic_OM4_05case.@adcroft @marshallward: Two points in particular that I'd like your feedback on:
post_datacall to be done without an extra copy. I have included astoresubroutine that could be used if we made it private, but otherwise those should be deleted.