Add GL90 parameterization for stacked shallow water #267
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
Created the new module remapping_attic to hold older versions of remapping code that are no longer used by MOM6. The subroutines is PosSumErrSignificant, remapByProjection, remapByDeltaZ and integrateReconOnInterval were moved to remapping_attic, where they can be tested by calling remapping_attic_unit_tests. The hard-coded old_algorithm logical module variable and the code it wraps were also eliminated. Also added a schematic description of the units of the real variables in the various routines in MOM_remapping and corrected some spelling errors. All answers are bitwise identical.
Moved interpolate_column and reintegrate_column (without changing anything) from MOM_diag_vkernels.F90 to MOM_remapping.F90 and incorporated the tests that had been in diag_vkernels_unit_tests into remapping_unit_tests. The entire MOM_diag_vkernels.F90 file was then removed. All answers are bitwise identical, although the module for two public routines was changed and a third was eliminated.
Remove missing_value arguments to interpolate_column and reintegrate_column, instead using 0 for the values in vanished cells. This change helps to address github.com/mom-ocean/issues/769. Also added comments schematically describing some of the argument units. Because 0 was already being used for the missing value (except in unit tests), all solutions are bitwise identical.
Added the new subroutine check_remapped_values with the duplicative error checking code in remapping_core_h and remapping_core_w, both to reduce code volume and promote code coverage, and to make the substance of these two routines easier to follow. All answers are bitwise identical.
- Adds `.gitlab/pipeline-ci-tool.sh` to enact most of the stages of the gitlab CI pipeline
- enables interactive/command-line reproduction of the pipeline
- `.gitlab/pipeline-ci-tool.sh` is documented in .gitlab/README.md with instructions
on how to use at the command line and what each function is doing
- All commands formerly in .gitlab-ci.yml are now one-line invocations of `.gitlab/pipeline-ci-tool.sh`
so .gitlab-ci.yml is now considerably smaller and easier to read with statements like
`.gitlab/pipeline-ci-tool.sh mrs-compile debug_gnu` or `.gitlab/pipeline-ci-tool.sh check-params gnu`
- Previously, all results were compared again the stored regression answers. This meant that
any error (e.g. layout) would show up as a fail of all types. We use the regression answers to
check the repro-symmetric mode and then compare everything else to repro-symmetric or other results
as appropriate. This allows us to distinguish between types of errors. The GH actions are doing it
this way, and we originally did this in the first forms of the pipeline, but in the last re-factor
I lazily switched to using the regression answers for everything.
- The subroutine categorize_axes cannot find the axes in ice restart files and gives warnings WARNING from PE 0: categorize_axes: Failed to identify x- and y- axes in the axis list (xaxis_1, yaxis_1, Time) of a variable being read from INPUT/ice_model.res.nc - This leads to an incorrect initializations and a subsequent sat.vap.press.overflow crash when using infra/FMS2 - Same experiment runs fine with infra/FMS1 - After the fix the infra/FMS1 and infra/FMS2 answers are bitwise identical
- The subroutine categorize_axes cannot find the axes in ice restart files and gives warnings
WARNING from PE 0: categorize_axes: Failed to identify x- and y- axes in the axis list (xaxis_1, yaxis_1, Time) of a variable being read from I
NPUT/ice_model.res.nc
- This leads to an incorrect initializations and a subsequent sat.vap.press.overflow crash when using
infra/FMS2
- Same experiment runs fine with infra/FMS1
- After the fix the infra/FMS1 and infra/FMS2 answers are bitwise
identical
Added a line initializing the string Cartesian to a blank string in categorize_axes, so that it not be uninitialized when it is used a few lines later.
Set the interpolation weights inside of interpolate_column to explicitly be the complement of one another, thereby saving an extra division at each point and reducing the number of variables that need to be stored, in preparation for the creation of a separate subroutine to find interface positions. This commit is mathematically equivalent to what was there before, and the extensive unit testing of interpolate_column is still passing, but it changes the value of some interpolated interface diagnostics at the level of roundoff (but not the MOM6 solutions themselves, as they do not depend on interpolate_column yet).
This patch introduces a new autoconf macro, AX_FC_CHECK_C_LIB, which confirms if the Fortran compiler can be linked to the netCDF C library. As with other netCDF tests, the nc-config tool is used if necessary (and available). This resolves some recent issues on platforms where netCDF and netCDF-Fortran are installed in separate locations, with different library directories (-L). It also resolves some false assumptions in configure.ac which presumed equivalent access by the configured C and Fortran compilers. Previously, we would test if the C compiler could be linked to netCDF, and then assume that the Fortran compiler shared the same relationship. We now use the Fortran compiler for both C and Fortran tests. This patch fixes many issues observed on MacOS systems, including some persistent problems on the GitHub Actions MacOS tests. For example, we can now use the default GCC 12 compilers, rather than forcing a rollback to GCC 11.
This patch fixes some issues with testing of C bindings in Fortran. Specifically, some tests are using a C compiler which may be unconfigured, causing unexpected errors. The autoconf script now uses the Fortran compiler to test these bindings, rather than using the C compiler to test for their existence. A new macro (AX_FC_CHECK_BIND_C) was added to run these tests. This achieves the actual goal (test of Fortran binding) on top of the original goal (availability of C function), while ensuring that the actual compiler of interest (FC) is used in the test. Two C-based tests are still present in the script for testing the size of jmp_buf and sigjmp_buf. The C compiler is now configured with the AX_MPI macro, and is only used to determine the size of these structs.
* Setup OBC segments for COBALT/OBGC tracers
- These are updates required to setup OBC segments for OBGC tracers.
- Since COBALT package has more than 50 tracers using the MOM6 table
mechanism for setting up OBC segments is not feasible. Rather, this
update delegates such setup to mechanims used in ocean_BGS tracers
leaving MOM6 mechanism for native tracers intact.
- Fixed issues caught by MOM6 githubCI
* Add capability to change obc segment update period
- COBALT tracers do not need as frequent segment bc updates and can
use a larger update period to speed up the model.
This commit introduces a new parameter DT_OBC_SEG_UPDATE_OBGC
that can be adjusted for obc segment update period.
- This commit applies the change only to BGC tracers but can easily
be changed to apply for all.
* Insert missing US%T_to_sec
- The unit conversion factor was missing causing a crash in a newer test.
* Updates from Andrew Ross
- Avoid low initial values in the tracer reservoirs
* Per Andrew Ross review
* corrected indentation per review
* Avoid using module vars per review request
- Reviewer asked to avoid using module variables with "save" attributes.
- This commit hides the module variables inside the existing OBC type.
* Coding style corrections per review
* Modification per review: do_not_log if .not.associated(CS%OBC)
Co-authored-by: Robert Hallberg <Robert.Hallberg@noaa.gov>
In this PR an option is added to use ice viscosity computed from the observed surface velocity, computed by the model and use a constant value (for debugging purposes). A new (char) parameter "ICE_VISCOSITY_COMPUTE" is introduced; its values can be "MODEL" (the ice viscosity computed by the model); "OBS" the ice viscosity is computed at the preprocessing step and read from a file (its name is defined by the parameter "ICE_STIFFNESS_FILE") into a variable with a name defined by "A_GLEN_VARNAME" parameter; "CONSANT" is a constant value defined by a parameter "A_GLEN". These changes are in MOM_ice_shelf_dynamics.F90. Minor changes are done to MOM_ice_shelf_initialize.F90 to correct units, scales.
Update from main
Added calls to get_param to set 12 input variable names in files via runtime parameters, including TIDEAMP_VARNAME, TEMP_COORD_VAR, SALT_COORD_VAR, THICKNESS_IC_VAR, INTERFACE_IC_RESCALE, TEMP_IC_VAR, SALT_IC_VAR, BASIN_VAR, TIDAL_DISSIPATION_VAR, ROUGHNESS_VARNAME, TIDEAMP_VARNAME and KH_BG_2D_VARNAME. Also added two new runtime parameters, THICKNESS_IC_RESCALE and INTERFACE_IC_RESCALE, to allow input thickness and interface height fields to be rescaled. A number of spelling errors in comments or output messages in the files that were being modified as a part of this commit, including changes in the documentation that appears in MOM_parameter_doc files. All answers are bitwise identical, but there are new entries and minor changes in many MOM_parameter_doc files.
Added calls to get_param to set 4 more input variable names in files via runtime, including U_IC_VAR, V_IC_VAR, OPEN_DY_CU_VAR and OPEN_DX_CV_VAR. Also added or amended comments describing internal variables to describe their units more consistently in MOM_shared_initialization. All answers are bitwise identical, but there may be new entries in some MOM_parameter_doc files.
Corrected a bug in converting depths read from an input file from units of cm to m when the ER03 version of tidal mixing is used. This commit will change answers when INT_TIDE_DISSIPATION = True, USE_CVMix_TIDAL = True, and TIDAL_ENERGY_TYPE = "ER03". There are no such configurations in the MOM6-examples pipeline tests, and it is not clear whether or where such a configuration has ever been used. This bug was introduced into dev/gfdl on Nov. 19, 2018 as a part of PR mom-ocean#883 in commit NOAA-GFDL@967e470, which was supposed to be a refactoring of this portion of the code without changing answers, but introduced this bug. This commit should restore solutions with impacted configurations to what they would have been before that earlier commit.
This patch removes the `build_{grid,data}.py` scripts from .testing's
tc4, along with the setup of the Python infrastructure used in the
.testing Makefile and GitHub Actions CI.
The Python scripts have been replaced with equivalent Fortran programs
which generate identical netCDF output.
A new rule (`preproc`) has been added to the .testing top Makefile for
generating the model input files.
The netCDF compiler dependenices are configured with autoconf, currently
duplicating the macros in `ac/configure.ac`. (NOTE: It may be possible
to share these with a common macro in ac/m4.
The configure script and Makefile are currently generated from
`configure.ac` and `autoreconf`. In the future, we could simply
pre-generate `configure` and add it to the repository.
This patch was motivated by the inability to install recent
netCDF-Python packages on systems with older gcc compilers, including
our main production machine. We could have possibly resolved this by
adding compiler configuration to pip, or perhaps reported the issue to
the netCDF-python project for them to resolve. But the costs of relying
on all this Python infrastructure is starting to exceed the benefits,
and I would recommend we excise it from our test suite.
GitLab CI includes the internal testing suite (.testing) and included an explicit setup of the Python environment (`make work/local-env`). The rule has since been removed, and the command now fails. This patch removes those steps, since we no longer use Python in the tests. It also slightly reworks the reporting of test output. Instead of re-running `make test`, it uses the `make test.summary` rule to report the final result.
Added a new logical argument to interpolate_column to specify whether the interpolated interface values outside of massless layers should be masked to zero. Also refactored the code in interpolate_column to separate out the determination of the interface position from the interpolation and the masking to facilitate the extension of this code to use higher order interpolation in planned subsequent changes. All answers are bitwise identical, although there is a new mandatory argument for a public interface.
Added ALE_remap_interface_vals and ALE_remap_vertex_vals to handle the interpolation of variables that are at the interfaces atop tracer cells or above the corners of the tracers cells from one grid to another. Because these are not yet used (but have been tested in calls that will be added with the next commit) all answers are bitwise identical, but there are two new publicly visible routines.
Added REMAP_AUXILIARY_VARS to control whether to remap the accelerations that are used in the predictor step of the split RK2 time stepping scheme. Also added the new routines remap_dyn_split_rk2_aux_vars, remap_OBC_fields and remap_vertvisc_aux_vars to do the remapping, and code to call these routines when REMAP_AUXILIARY_VARS is true. By default, REMAP_AUXILIARY_VARS is false, and all answers are bitwise identical, but the entire MOM6-examples regression suite has been run with this set to true, and they do appear to give physically plausible answers in all cases, partially addressing the issue noted at github.com/NOAA-GFDL/issues/203. New entries are added to the MOM_parameter_doc files, and there are three new publicly visible routines, but by default answers do not change.
* Adds the option to set the diffusivity KHTH as horizontally varying * Can be enabled via READ_KHTH = True, filename is provided by user via KHTH_FILE * Will return error if user sets both READ_KHTH = True and KHTH > 0
* full file path is now set as INPUTDIR/KHTH_FILE, where both INPUTDIR and KHTH_FILE are runtime parameters
thickness diffusivity --> isopycnal height diffusivity
Corrected the units written to the output files for 4 diagnostics (CAu_Stokes, CAv_Stokes, area_shelf_h and sfc_mass_flux) and added missing units arguments to the get_param calls for some (mostly unlogged) parameters. The logged calls where units are added include those for EKE_MAX, NDIFF_DRHO_TOL, NDIFF_X_TOL, and IMPULSE_SOURCE_TIME, while some unnecessary carriage returns were removed in the descriptions of some of these and closely related parameters. Also added units to the comment describing the AGlen argument to initialize_ice_AGlen. All answers are bitwise identical, but there can be minor changes in the metadata of some files, and some MOM_parameter_doc and available_diags files might exhibit minor changes.
Added a missing scale factor in the DENSE_WATER_EAST_SPONGE_SALT get_param call in dense_water_initialize_sponges, and added comments describing the local variables (and their units) throughout the dense_water_initialization module. The variable set by DENSE_WATER_SILL_HEIGHT was unused and it probably was always intended to be DENSE_WATER_SILL_DEPTH, which it now is. Units arguments were also added to two of the unlogged get_param calls in this module. Without this change, this test case would not reproduce with dimensional rescaling due to a scale factor that was omitted when salinity was being rescaled on May 3, 2022, which became a part of PR mom-ocean#122 to dev/gfdl, but answers should not change when dimensional rescaling of salinities is not used. All answers and output in the MOM6-examples test suite are bitwise identical.
Removed meaningless units arguments from 31 get_param calls for integer, character, or logical parameters. All answers are bitwise identical, but some entries in the various parameter_doc files are changed.
This patch fixes two issues in the preprocessing of tc4. * ocean_hgrid.nc is marked as a dependency of gen_data * Multiple ouputs are handled more safely in the Makefile Expanding on the second point: We were directing one rule to produce two output files, which resulted in the rule being run twice when invoked in parallel (make -j). This has been replaced with the recommended solution for handling concurrent outputs: use one to generate both, and connect the second to the first with a separate rule. I have also generalized the `make` command in the .testing Makefile. This should address (and hopefully fix) some intermittent errors in the .testing build on Gaea.
.testing: Fix concurrency errors in tc4 rules
Add units arguments to 30 unlogged get_param calls in 4 user modules (DOME2d_initialization, ISOMIP_initialization, Kelvin_initialization and seamount_initialization) to help detect inconsistent units and scaling factors. Also added comments describing many internal real variables and their units in the DOME2d_initialization module and seamount_initialize_temperature_salinity. All answers and output are bitwise identical.
Updated the default values of remap_answer_date as declared during the specification of the wave_speed_CS, remapping_CS and regridding_CS to 99991231 (to use the very latest version of the algorithms) if these control structures are used without further initialization via optional arguments to their various initialization routines (wave_speed_init(), initialize_remapping() and initialize_regridding() or set_regrid_params()). In testing both via the TC tests and the MOM6-examples regression suite, all answers are bitwise identical, because every instance where these are used either does have an appropriate call to the initialization routine, or because (as is the case for remapping diagnostics) they are hard-coded to use the PLM remapping scheme, which is not impacted by this answer_date. It is conceivable, if unlikely, however, that there could be cases outside of the standard MOM6 code that are using this remapping capability without proper initialization of this flag, in which case answers could change. In addition, a comment was updated in the interp_CS_type indicating the cases that would be altered by changing the hard coded default answer_date there, but the default in that case was not changed. All answers and output are identical in the MOM6-examples and TC test suites.
Completed the dimensional rescaling of internal variables in the Rossby_front_initialization and BFB_initialization modules, including the addition of comments describing numerous internal variables and their units. The BFB module was more extensively modified to place logging of the module name and version before the get_param calls for this module, following the MOM6 pattern, to correct spelling errors in get_param descriptions, and to use the grid mask in the grid file rather than a reexamination of the minimum depth to determine the land-mask. The internal subroutine write_BFB_log is no longer needed and has been folded into BFB_set_coord. All answers are bitwise identical, but there are minor changes in the MOM_parameter_doc files for the buoy_forced_basin test case.
The Gaea GitLab test configuration includes two no-library compilation tests for ocean-only and ice-ocean configurations. They are currently configured to run in the same directory. Recently, several unusual I/O errors suggest problems due to concurrent testing of the builds. In order to rule out this possibility, this patch moves the ice-ocean test to a separate directory.
This adds a new vertical viscosity parameterization as in Greatbatch and Lamb (1990), Ferreira & Marshall (2006) and Zhao & Vallis (2008), hereafter referred to as the GL90 vertical viscosity parameterization. This vertical viscosity scheme redistributes momentum in the vertical, and is the equivalent of the Gent & McWilliams (1990) parameterization, but in a TWA (thickness-weighted averaged) set of equations. The vertical viscosity coefficient nu is computed from kappa_GM via thermal wind balance, and the following relation: nu = kappa_GM * f^2 / N^2. The vertical viscosity del_z ( nu del_z u) is applied to the momentum equation with stress-free boundary conditions at the top and bottom. In the current implementation, kappa_GM is assumed either (a) constant or as (b) having an EBT structure. A third possible formulation of nu is depth-independent: nu = f^2 * alpha The latter formulation would be equivalent to a kappa_GM that varies as N^2 with depth. Currently, the GL90 parameterization is only implemented in stacked shallow water (SSW) mode, in which case we have 1/N^2 = h/g'. More specifically, this commit adds a new subroutine that computes the couping coefficient associated with GL90 via a_cpl_gl90 = nu / h = kappa_GM * f^2 / g' or a_cpl_gl90 = nu / h = f^2 * alpha / h. Further, a_cpl_gl90 is multiplied by a function (botfn), which is 0 within the GL90 bottom boundary layer, whose depth is set by Hbbl_gl90, and 1 otherwise. This modification is necessary to avlid fluxing momentum into vanished layers that ride over steep topography. Finally, a_cpl_gl90 is added to a_cpl, where the latter is the coupling coefficient associated with the remaining vertical stresses, used in the vertical viscosity solver. More information can be found in Loose et al. (https://www.essoar.org/doi/abs/10.1002/essoar.10512867.1), Appendix B.
|
Closing this because accidentally targeted PR at |
marshallward
pushed a commit
that referenced
this pull request
Apr 9, 2024
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.
GL90 parameterization
This adds a new vertical viscosity parameterization as in Greatbatch and Lamb (1990), Ferreira & Marshall (2006) and Zhao & Vallis (2008), hereafter referred to as the GL90 vertical viscosity parameterization. This vertical viscosity scheme redistributes momentum in the vertical, and is the equivalent of the Gent & McWilliams (1990) parameterization, but in a TWA (thickness-weighted averaged) set of equations. The vertical viscosity coefficient
nuis computed fromKD_GL90(corresponding tokappa_GM) via thermal wind balance, and the following relation:Currently, the GL90 parameterization is only implemented in stacked shallow water (SSW) mode, in which case we have
1/N^2 = h/g'.The vertical viscosity
del_z ( nu del_z u)is applied to the momentum equation with stress-free boundary conditions at the top and bottom.GL90 viscosity coefficient: Current options
In the current implementation,
KD_GL90is assumed either (a) constant or (b) horizontally varying. In both cases, (a) and (b), one can additionally impose an EBT structure in the vertical forKD_GL90. Another possible formulation of nu is depth-independent:The latter formulation would be equivalent to a
kappa_GMthat varies asN^2with depth.Figure: Vertical viscosity coefficient
Kv_vfor(upper left) baseline with GL90 parameterization turned off;
(upper right) GL90 parameterization turned on and spatially constant
KD_GL90;(lower left) GL90 parameterization turned on and
KD_GL90which is horizontally constant but varies in the vertical with and EBT structure;(lower right) GL90 parameterization turned on with vertically constant GL90 viscosity, achieved by setting
alpha_GL90 = 3 * 10^7 m^2 s.Note that in the last panel, the diagnostic
Kv_vis not exactly constant in the vertical becausewhere
h_vis the layer thickness used at a v-velocity points in the vertical viscosity scheme. The third factor is the composition of a harmonic and arithmetic mean, and is not equal to the identity.Implementation done via coupling coefficients
More specifically, this commit adds a new subroutine that computes the coupling coefficient associated with GL90 via
or
Further,
a_cpl_gl90is multiplied by a function (botfn), which is 0 within the GL90 bottom boundary layer, whose depth is set byHbbl_gl90, and 1 otherwise. This modification is necessary to avoid fluxing momentum into vanished layers that ride over steep topography.Finally,
a_cpl_gl90is added toa_cpl, where the latter is the coupling coefficient associated with the remaining vertical stresses, used in the vertical viscosity solver.More information can be found in Loose et al. (https://www.essoar.org/doi/abs/10.1002/essoar.10512867.1), Appendix B.