Merge main to dev/gfdl (2013-10-25)
#1612
Closed
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
The autoconf Python interpreter search was slightly modified to search for Python even if $PYTHON is set to an empty string. This is done by unsetting PYTHON if it is set but empty, then following the usual macro. This was required since `export PYTHON` in a Makefile will create the `PYTHON` variable but will assign it no value (i.e. empty string). This causes issues in some build environments. The backup `configure~` script was also added to the developer `ac-clean` cleanup rule.
Refactored 6 files in the ALE directory to calculate the nominal depth in thickness units in a single place (it is done in regridding_main now) and pass it to the various places where it is used, in a preparatory step to modify how this calculation is done in non-Boussinesq mode. There are new arguments to several publicly visible routines, including: - Add non_depth_H arguments to hybgen_regrid, build_zstar_grid, build_sigma_grid, build_rho_grid, build_grid_HyCOM1, build_grid_adaptive, build_adapt_column and build_grid_arbitrary - Add optional zScale arguments to build_zstar_grid and build_grid_HyCOM1 - Add unit_scale_type arguments to regridding_main, ALE_regrid_accelerated and ALE_offline_inputs Also eliminated an incorrect rescaling GV%Z_to_H facto when calculating the total column thickness from the layer thicknesses when an ice shelf is used with a Hycom grid. This would have caused dimensional consistency testing to fail. Added the new runtime parameters HYBGEN_H_THIN, HYBGEN_FAR_FROM_SURFACE HYBGEN_FAR_FROM_BOTTOM, and HYBGEN_DENSITY_EPSILON to set previously hard-coded dimensional parameters used in the Hybgen regridding code and store these values in new variables in hybgen_regrid_CS. Two of these are no longer passed to hybgen_column_regrid as separate parameters. By default these new runtime parameters recover the previous hard-coded values. Also eliminated an unused block of code in build_rho_column. Several comments documenting variables or their units were also added. All answers are bitwise identical, but there are 4 new runtime parameters that would appear in some MOM_parameter_doc files and there are changes to the arguments to 11 routines.
As described in issue #372, I would like to be able to create restart files that contain information about the particle location. These files will be written at the same time as other restart files. I cannot add these calls directly to the driver, because the driver does not have information about the particle location. We have added save_MOM6_internal_state as a subroutine in MOM.F90, and we added calls to this subroutine from each of the drivers. We hope this will allow for more new packages to write restart files in the future. Co-authored by Spencer Jones <[email protected]>
Added the integer valid_SpV_halo to the thermo_var_ptrs type to indicate whether the SpV_array has been updated and its valid halo size, to facilitate error detection and debugging in non-Boussinesq mode. Tv%valid_SpV_halo is set to the halo size in calc_derived_thermo or after a halo update is done to tv%SpV_avg, and it is set to a negative value right after calls that change temperatures and salinities (such as by ALE remapping) unless there is a call to calc_derived_thermo. Tests for the validity of tv%SpV_avg are added to the routines behind thickness_to_dz, with fatal errors issued if invalid arrays would be used, but more tests could perhaps be used in any parameterization routines where tv%SpV_avg is used directly. Handling the updates to tv%SpV_avg this way helps to avoid unnecessary calls to calc_derived_thermo, which in turn has equation of state calls that can be expensive, while also providing essential verification of new code related to the non-Boussinesq code. These tests can probably be commented out or removed for efficiency once there is a full suite of regression tests for the fully non-Boussinesq mode of MOM6. In addition, a new optional debug argument was added to calc_derived_thermo which can be used to triggers checksums for the variables used to calculate tv%SpV_avg. One call to calc_derived_thermo was also added just before the initialization call to ALE_regrid that will be needed with the next commit, but does not change answers yet. All answers are bitwise identical, but there is a new element in a transparent and widely used type and a new optional argument to a public interface.
Use RHO_KV_CONVERT instead of RHO_0 to set the non-Boussinesq version of GV%m_to_H, so that there is a mechanism for testing the independence of the fully non-Boussinesq mode from the Boussinesq reference density. With this change, GV%Z_to_H is not guaranteed to be equal to (GV%Z_to_m*GV%m_to_H), with the latter expression preferred when setting parameters. By default the two parameters are the same, and they will probably only ever differ in testing the code. All Boussinesq solutions are bitwise identical, but there are differences in the description of RHO_KV_CONVERT that will appear in MOM_parameter_doc files.
Add new arguments to 7 routines that will be needed for the non-Boussinesq capability, but do not use them yet, so that there will be fewer cross file dependencies as the various changes are being reviewed simultaneously. The impacted interfaces are MEKE_int, vertvisc_coef, sumSWoverBands, KPP_calculate, differential_diffuse_T_S, set_BBL_TKE, and apply_sponge In the three step_MOM_dyn_... routines and in calculateBuoyancyFlux1d, this change includes calls to thickness_to_dz to calculate the new vertical distance arrays that will be passed into vertvisc_coef or sumSWoverBands. The only place where the new arguments are actually used is in sumSWoverBands and set_opacity where the changes are particularly simple. All answers are bitwise identical, but there are new non-optional arguments to seven publicly visible routines.
* Restore functionality for reading slices from 3d volumes in MOM_io - The recent MOM_io modifications in support of FMS2_io accidentally removed support for reading on-grid data (same horizontal grid as model) k-slices. This is needed in some configurations in the model state initialization. * Add FMS1 interfaces * Additional patches to enable reading ongrid state initialization data - read local 3d volume rather than attempting to slice ongrid data vertically. - Related bugfixes in MOM_io
- We were reading KV_ML_INVZ2 without logging, then checking for KVML and finally logging based on a combination of the two. This had the side affect that we get warnings about not using KVML even if KVML was not present. - The fix checks for KVML first, and then changes the default so that when KVML=1e-4 is replaced by KV_ML_INVZ2=1e-4 we end up with no warnings and KVML can be obsoleted safely. Note: this commit alone does not remove all warnings from the MOM6-examples suite because we still need to fix the MOM_input that still use KVML - KVML needs to be unscaled since it is the default for KV_ML_INVZ2 - tc3 used KVML and has been corrected.
Use the new runtime parameter RHO_PGF_REF instead of RHO_0 to set the reference density that is subtracted off from the other densities when calculating the finite volume pressure gradient forces. Although the answers are mathematically equivalent for any value of this parameter, a judicious choice can reduce the impacts of roundoff errors by about 2 orders of magnitude. By default, RHO_PGF_REF is set to RHO_0, and all answers are bitwise identical. However, there is a new runtime parameter that appears in many of the MOM_parameter_doc.all files.
The message that a file is being created was issued as a WARNING when we all agree it should really be a NOTE. Depth_list.nc is read if it is present to avoid recomputing a sorted list, but the absence of the file is not an error and does not warrant a warning. Changes: - Changed WARNING to NOTE. - Removed MOM_mesg from imports since it wasn't being used.
The interpolation scheme for state-dependent diagnostic coordinates was incorrectly registering as the same parameter as the main model. This meant it was never possible to change the interpolation scheme from the default (which was not the same as the main model). Fix registers the generated parameter name which was always computed but not used. A typical example of the generated parameter is "DIAG_COORD_INTERP_SCHEME_RHO2".
Fixed a bug in which wave_speed_init was effectively discarding any values of mono_N2_depth passed to it via the optional argument mono_N2_depth, but also changed the default value of RESOLN_N2_FILTER_DEPTH, which was previously being discarded, to disable the monotonization and replicate the previous results. There were also clarifying additions made to the description how to disable RESOLN_N2_FILTER_DEPTH. This will change some entries in MOM_parameter_doc files, and it will change solutions in cases that set RESOLN_N2_FILTER_DEPTH to a non-default value and have parameter settings that use the resolution function to scale their horizontal mixing. There are, however, no known active simulations where the answers are expected to change.
Revised the calculation of gprime and the coordinate densities (GV%Rlay) in fully non-Boussinesq mode to use the arithmetic mean of adjacent coordinate densities in the denominator of the expression for g_prime in place of RHO_0. Also use LIGHTEST_DENSITY in place of RHO_0 to specify the top-level coordinate density in certain coordinate modes. Also made corresponding changes to the fully non-Boussinesq APE calculation when CALCULATE_APE is true, and eliminated an incorrect calculation of the layer volumes in non-Boussinesq mode using the Boussinesq reference density that was never actually being used when CALCULATE_APE is false. This commit will change answers in some fully non-Boussinesq calculations, and an existing runtime parameter is used and logged in some new cases, changing the MOM_parameter_doc file in those cases.
Refactored thickness_diffuse when in non-Boussinesq mode to avoid any dependencies on the Boussinesq reference density, and to translate the volume streamfunction into the mass streamfunction using an appropriately defined in-situ density averaged to the interfaces at velocity points. This form follows the suggestions of Appendix A.3.2 of Griffies and Greatbatch (Ocean Modelling, 2012) when in non-Boussinesq mode. Thickness_diffuse_full was also revised to work properly in non-Boussinesq mode (and not depend on the Boussinesq reference density) when no equation of state is used. As a part of these changes, the code now uses thickness-based streamfunctions and other thickness-based internal calculations in MOM_thickness_diffuse. For example, the overturning streamfunctions with this change are now in m3/s in Boussinesq mode, but kg/s in non-Boussinesq mode. These changes use a call to thickness_to_dz to set up a separate variable with the vertical distance across layers, and in non-Boussinesq mode they use tv%SpV_avg to estimate in situ densities. Additional debugging checksums were added to thickness_diffuse. The code changes are extensive with 15 new or renamed internal variables, and changes to the units of 9 other internal variables and 3 arguments to the private routine streamfn_solver. After this change, GV%Rho, GV%Z_to_H and GV%H_to_Z are no longer used in any non-Boussinesq calculations (12 such instances having been elimated). Because some calculations have to be redone with the separate thickness and dz variables, this will be more expensive than the original version. No public interfaces are changed, and all answers are bitwise identical in Boussinesq or semiBoussinesq mode, but they will change in non-Boussinesq mode when the isopycnal height diffusion parameterization is used.
* Salt data structures * First steps at brine plume: pass info from SIS2 * The brine plume parameterization, - including now passing the dimensional scaling tests. * Fix problem when running Tidal_bay case with gnu. * Avoiding visc_rem issues inside land mask. Tweaking the brine plume code. * Using the proper MLD in the brine plumes - it now works better on restart * Always including MLD in call to applyBoundary... - I could move it up and make it not optional. * Adding some OpenMP directives to brine plumes
This commit brings the drifters interface up-to-date with the current version of the drifters package, which requires h (layer thickness) to calculate the vertical movement of particles. The interfaces in the code and in config_src/external are updated to pass this information to the drifters package.
Pass dt_kappa_smooth to calc_isoneutral_slopes and vert_fill_TS in units of [H Z ~> m2 or kg m-1] instead of [Z2 ~> m2] for consistency with the units of other diffusivities in the code and to reduce the depenency on the Boussinesq reference density in non-Boussinesq configurations. In addition to the changes to the units of these two arguments, there is a new unit_scale_type argument to vert_fill_TS and MOM_calc_varT and a new verticalGrid_type argument to MOM_stoch_eos_init. The units of 4 vertical diffusivities in the control structures in 4 different modules are also changed accordingly. All answers are bitwise identical in Boussinesq mode, but they can change for some non-Boussinesq configurations. There are new mandatory arguments to three publicly visible routines.
Added a comment justifying the use of a fixed rescaling factor for the diffusivity used in vert_fill_TS. All answers and output are identical.
Added the new public interface find_ustar to extract the friction velocity from either a forcing type argument, or a mech_forcing_type argument, either directly or from tau_mag, and in non-Boussinesq mode by using the time-evolving surface specific volume. Find_ustar is an overloaded interface to find_ustar_fluxes or find_ustar_mech_forcing, which are the same but for the type of one of their arguments. For now, the subroutines bulkmixedlayer, mixedlayer_restrajt_OM4, mixedlayer_restrat_Bodner and mixedlayer_restrat_BML are calling find_ustar to avoid code duplication during the transition to work in fully non-Boussinesq mode, but it will eventually be used in about another half dozen other places. All Boussinesq answers are bitwise identical, but non-Boussinesq answers will change and become less dependent on the Boussinesq reference density, and there is a new publicly visible interface wrapping two subroutines.
Changed the units of the optional mono_N2_depth argument to wave_speed, wave_speed_init and wave_speed_set_param in thickness units instead of height units. Accordingly, the units of one element each in the diagnostics_CS and wave_speed_CS and a local variable in VarMix_init are also changed to thickness units. The unit descriptions of some comments describing diagnostics were also amended to also describe the non-Boussinesq versions. Because this is essentially just changing when the unit conversion occurs, all answers are bitwise identical, but there are changes to the units of an optional argument in 3 publicly visible routines.
Added the new runtime parameter BT_RHO_LINEARIZED to specify the density that is used to convert total water column thicknesses into mass in non-Boussinesq mode with linearized options in the barotropic solver or when estimating the stable barotropic timestep without access to the full baroclinic model state. The default is set to RHO_0 and answers do not change by default. This new parameter is used in non-Boussinesq mode with some options in btcalc and find_face_areas, when LINEARIZED_BT_CORIOLIS = True or BT_NONLIN_STRESS = False, and in the unit conversion of the ice strength with dynamic pressure. Also cancelled out factors of GV%Z_to_H in MOM_barotropic.F90 to simplify the code and reduce the dependence on the value of GV%Rho_0 in non-Boussinesq mode. This involved changing the units of 4 variables in the barotropic_CS type, 3 internal variables in btstep and an internal variable in barotropic_init to use thickness units. The rescaled internal variable mass_to_Z was also replaced with the equivalent GV%RZ_to_H. There are also 4 new debugging messages. Also modified the units of the gtot_est argument to match those of pbce. There is a new element in barotropic_CS. Because GV%Z_to_H is an exact power of 2 in Boussinesq mode, all answers are bitwise identical in that mode, but in non-Boussinesq mode this conversion involves multiplication and division by GV%Rho_0, so while all answers are mathematically equivalent, this change does change answers at roundoff in non-Boussinesq mode. Additionally there is a new runtime parameter that will appear in some MOM_parameter_doc files.
Use thickness_to_dz to convert layer thicknesses to depths in 4 tracer modules (DOME_tracer, dye_example, ideal_age_example and nw2_tracers) so that this conversion is done correctly in non-Boussineq mode, and there is no longer any dependency on the Boussinesq reference density in that mode. This change includes the addition of a thermo_var_ptrs argument to 5 routines (initialize_DOME_tracer, initialize_dye_tracer, dye_tracer_column_physics ideal_age_tracer_column_physics and count_BL_layers) and changes to the units of some internal variables, and the addition of 6 new 2-d or 3-d arrays with the vertical distance across layers. An unused param_file_type argument to initialize_DOME_tracer was also eliminated. Comments were also added to describe the units of 5 of the variables in the ideal age tracer control structure and 7 internal variables in that same module, and there was some minor cleanup of the formatting cf calls in tracer_flow_control_init. There was some minor refactoring in the ns2_tracers module to use SZK_(GV) instead of SZK_(G) to declare the vertical extent of some arrays, and the vertical indexing convention for interfaces in nw2_tracer_dist was revised from starting at 0 to start at 1 for consistency with all the other code in MOM6. Also moved the code to do halo updates for the physical model state variables and call calc_derived_thermo before calling tracer_flow_control_init, because some routines there are now using the layer average specific volume to convert between thicknesses and heights when in non-Boussinesq mode. All answers in Boussinesq mode are bitwise identical, but these passive tracer modules have slightly different answers in non-Boussinesq mode. There are changes to the non-optional arguments to 4 public interfaces.
Changed a recently added OMP directive for plume_flux from private to firstprivate to reflect how this variable is actually used. This bug was introduced with PR #401, but was causing sporadic failures in some of our pipeline tests with the intel compiler (essentially due to initialized memory when openMP is used) for subsequent commits.
This patch merges the internal `save_restart` function with the new `save_MOM6_internal_state` function into a new general MOM restart function. It also makes an effort to eliminate `MOM_restart` as a driver dependency, narrowing the required MOM API for existing and future drivers. Also removes the `restart_CSp` argument from `MOM_wave_interface_init`, since it appeared to be used for nothing.
MOM simulations typically abort of the restart directory (usually RESTART) are absent. This patch adds POSIX support for mkdir() and creates the directory if it is missing.
Using inquire() to check for directory existence is not possible, since at least one compiler (Intel) does not consider directories to be files. The inquire call is replaced with a C interface to the POSIX stat() function. We do not fully emulate the behavior of stat, but we use its return value to determine existence of directories. This provides a more reliable method for identifying the existence of the directory. This should resolve many of the observed problems with RESTART creation in coupled runs.
Cancelled out factors of GV%Z_to_H in MOM_hor_visc.F90 to simplify the code and reduce the dependence on the value of GV%Rho_0 in non-Boussinesq mode. This involved changing the units of 3 internal variables in horizontal_viscosity and one element in the hor_visc_CS type to use thickness units or their inverse. Because GV%Z_to_H is an exact power of 2 in Boussinesq mode, all answers are bitwise identical in that mode, but in non-Boussinesq mode this conversion involves multiplication and division by GV%Rho_0, so while all answers are mathematically equivalent, this change does change answers at roundoff in non-Boussinesq mode.
Revised the units of 12 vertvisc_type elements to be based on thicknesses, so that vertical viscosities (in [H Z T-1 ~> m2 s-1 or Pa s]) are stored as dynamic viscosites when in non-Boussinesq mode, with analogous changes to the diapycanl diffusivity (now in [H Z T-1 ~> m2 s-1 or kg m-1 s-1]). Similarly changed the units of the 2 Rayleigh drag velocity elements (Ray_u and Ray_v) of the vertvisc_type from vertical velocity units to thickness flux units and to more accurately reflect the nature of these fields. The bottom boundary layer TKE source element (TKE_BBL) was also revised to [H Z2 T-3 ~> m3 s-3 or W m-2]. This commit also adds required changes to the units of the viscosities or shear-driven diffusivities returned from KPP_calculate, calculate_CVMix_shear, calculate_CVMix_conv, Calculate_kappa_shear, Calc_kappa_shear_vertex, calculate_tidal_mixing and calculate_CVMix_tidal. Because GV%Z_to_H is an exact power of 2 in Boussinesq mode, all answers are bitwise identical in that mode, but in non-Boussinesq mode this conversion involves multiplication and division by GV%Rho_0, so while all answers are mathematically equivalent, this change does change answers at roundoff in non-Boussinesq mode unless GV%Rho_0 is chosen to be an integer power of 2.
Rescaled diapycnal diffusivities passed as arguments in non-Boussinesq mode, to be equivalent to the thermal conductivity divided by the heat capacity, analogously to the difference between a kinematic viscosity and a dynamic viscosity, so that the new units are [H Z T-1 ~> m2 s-1 or kg m-1 s-1]. This includes changing the units of 4 arguments to set_diffusivity; 3 arguments each to calculate_bkgnd_mixing, add_drag_diffusivity, add_LOTW_BBL_diffusivity, user_change_diff, calculate_tidal_mixing and add_int_tide_diffusivity; 2 arguments to KPP_calculate, calculate_CVMix_conv, compute_ddiff_coeffs, differential_diffuse_T_S, entrainment_diffusive, double_diffusion, add_MLrad_diffusivity, and calculate_CVMix_tidal; and one argument to energetic_PBL. The units of 36 internal variables were also changed, as were a total of 29 elements in various opaque types, including 8 elements in bkgnd_mixing_cs, 2 in diabatic_CC, 3 in tidal_mixing_diags type, 1 in user_change_diff_CS, 9 in set_diffusivity_CS type, and 6 elements in diffusivity_diags. Two new internal variables were added, and several redundant GV%H_to_Z conversion factors were also cancelled out, some using that GV%H_to_Z*GV%Rho0 = GV%H_to_RZ. Because GV%Z_to_H is an exact power of 2 in Boussinesq mode, all answers are bitwise identical in that mode, but in non-Boussinesq mode this conversion involves multiplication and division by GV%Rho_0, so while all answers are mathematically equivalent, this change does change answers at roundoff in non-Boussinesq mode unless GV%Rho_0 is chosen to be an integer power of 2.
Changed the units for TKE arguments to [H Z2 T-3 ~> m3 s-3 or W m-2] for find_TKE_to_Kd, add_drag_diffusivity, calculate_tidal_mixing and add_int_tide_diffusivity, with similar changes to the units of 21 diagnostics or internal variables in the same routines and in add_LOTW_BBL_diffusivity and add_MLrad_diffusivity. Dozens of unit conversion factors were also cancelled out with these changes, including using that GV%Z_to_H/GV%Rho_0 = GV%RZ_to_H and that GV%Rho_0*GV%H_to_Z = GV%H_to_RZ for both Boussinesq or non-Boussinesq configurations. Because GV%Z_to_H is an exact power of 2 in Boussinesq mode, all answers are bitwise identical in that mode, but in non-Boussinesq mode this conversion involves multiplication and division by GV%Rho_0, so while all answers are mathematically equivalent, this change does change answers at roundoff in non-Boussinesq mode unless GV%Rho_0 is chosen to be an integer power of 2.
Nested includes are tracked for the purpose of include flags (-I), but not with respect to the content within those files. This patch tracks the module usage statements (`use ...`) inside of any include files and adds them to the Makefile rules of the top-level file. This was implemented within the `nested_inc` function by adding a new argument.
Changed the default value of Z_INIT_REMAP_GENERAL to true for fully non-Boussinesq configurations. All existing fully non-Boussinesq cases that use INIT_LAYERS_FROM_Z_FILE = True use this setting, and it is likely that such cases will not work at all if this is false. This change reduces the parameters that need to be changed to go between equivalent Boussinesq and non-Boussinesq configurations to just BOUSSINESQ and SEMI_BOUSSINESQ. The previous default (false) is being retained for any Boussinesq or semi-Bousssinesq cases so that all answers and output are bitwise identical in those cases, but by default some non-Boussinesq solutions change answers, and there are changes to the MOM_parameter_doc files for those non-Boussinesq cases.
This commit revises lateral_mixing_coeffs to work in an appropriate mixture of thickness and vertical extent variables to avoid any dependence on the Boussinesq reference density in non-Boussinesq mode, while retaining the previous answers in Boussinesq mode. This commit adds the new runtime parameter FULL_DEPTH_EADY_GROWTH_RATE to indicate that the denominator of an Eady growth rate calculation should be based on the full depth of the water column, rather than the nominal depth of the bathymetry. The new option is only the default for fully non-Boussinesq cases. A primordial horizontal indexing bug was corrected in the v-direction slope calculation. Because it only applies for very shallow bathymetry, does not appear to impact any existing test cases and went undetected for at least 12 years, it was corrected directly rather than wrapping in another new runtime flag. However, this bug is being retained for now in a comment to help with review and debugging if the answers should change unexpectedly in some yet-to-be identified configuration. Two debugging checksums were added for the output variables calculated in calc_resoln_function. The case of some indices was corrected to follow the MOM6 soft convention using case to indicate the staggering position of variables. The previously incorrect units of one comment were also fixed. There is a new logical element in the VarMix_CS type. To accommodate these changes there are three new internal variables in calc_slope_functions_using_just_e. A total of 9 GV%H_to_Z conversion factors were eliminated with this commit. N2 is no longer calculated separately in calc_slope_functions_using_just_e, but this code is left in a comment as it may be instructive. This commit involved changing the units of one internal variable in calc_QG_Leith_viscosity to use inverse thickness units (as its descriptive comment already indicated). There are already known problems with calc_QG_Leith_viscosity as documented with a fatal error; this will be addressed in a subsequent commit. All answers are bitwise identical in the existing MOM6-examples test suite, but they will change when fully non-Boussinesq, and there is a new entry in some MOM_parameter_doc files.
Anonymous uploads to Codecov are throttled (though it is not clear if this is happening on the Codecov or the GitHub Actions side). Regardless, the advice to get around this throttling seems to be to use the Codecov token associated with the project. The .testing Makefile was modified to use this token if provided, and the GitHub Actions environment now attempts to fetch this from the secrets of the repository. Hopefully individual groups can set this for their own projects, and it will fall back to anonymous upload if unset, but I guess we'll have to see how this plays out.
Renamed OBC related variables to emphasize that they are sea surface heights, and not sea surface heights or total column mass depending on the Boussinesq approximation. Also changed the units of these renamed variables to [Z ~> m] instead of [H ~> m or kg m-2]. The renamed element of the OBC_segment_type is eta (which is now SSH). The renamed and rescaled elements of the BT_OBC_type are H_[uv] (which are now dZ_[uv]) and eta_outer_[uv] (which are now SSH_outer_[uv]). The internal variables H_[uv] and h_in in apply_velocity_OBCs are now dZ_[uv] and ssh_in. Tidal_elev in update_OBC_segment_data and cff_eta in tidal_bay_set_OBC_data were rescaled but not renamed. There is also a new vertical grid type argument to apply_velocity_OBCs for use in changing the scaling of the SSH variables. A total of 11 GV%Z_to_H or GV%H_to_Z rescaling factors were cancelled out as a result of these changes, while 16 new ones were added, but most of these will in turn be dealt with in a follow-on commit that enables the use of OBCs in non-Boussinesq mode. Because any cases that use Flather open boundary conditions with the non-Boussinesq mode issue a fatal error, all answers are bitwise identical, but there are changes to the names and units of elements in a transparent type.
Read in open boundary condition SSH and SSHamp segment data directly in rescaled units of [Z ~> m] instead of [H ~> m or kg m-2], which then has to undergo a subsequent conversion. All answers in Boussinesq mode are bitwise identical and non-Boussinesq mode does not work with OBCs yet.
Added the new subroutine find_col_avg_SpV to return the column-averaged specific volume based on the coordinate mode and the layer averaged specific volumes that are in tv%SpV_avg. All answers are bitwise identical, but there is a new publicly visible routine.
Get Flather open boundary conditions working properly in non-Boussinesq mode. This includes calculating the column-average specific volume in step_MOM_dyn_split_RK2 and passing it as a new argument to btstep. Inside of btstep, this is copied over into a wide halo array and then passed on to set_up_BT_OBC and apply_velocity_OBCs, where it is used to determine the free surface height or the vertical column extent (in [Z ~> m]) from eta for use in the Flather radiation open boundary conditions. In addition, there are several places in MOM_barotropic related to the open boundary conditions where the usual G%bathyT needed to be replaced with its wide-halo counterpart, CS%bathyT. Also, a test was added for massless OBC columns in apply_velocity_OBCs, which are then assumed to be dry rather than dividing by zero. A fatal error message that is triggered in the case of Flather open boundary conditions in non-Boussiesq mode was removed. With this change, all Boussinesq answers are bitwise identical, but non-Boussinesq cases with Flather open boundary conditions are now working and giving answers that are qualitatively similar to the Boussinesq cases.
Add the new element dZtot to the OBC_segment_type to hold the total vertical extent of the water column, and use thickness_to_dz in update_OBC_segment_data to convert the layer thicknesses to the vertical layer extents used to set dZtot. This change leads to the cancellation of 6 unit conversion factors. All answers are bitwise identical in Boussinesq mode, but they will change and become less dependent on the Boussinesq reference density in some Boussinesq cases with certain types of open boundary condition.
This adds logical tests to determine the global OBC data update patterns to allow for halo updates related to these updates without having the model hang up with inconsistently blocked message passing calls. This commit corrects a bug with some test cases that was introduced with the halo update after the call to thicknesses_to_dz in update_OBC_segment_data two commits previously.
Previously, ice-shelf Dirichlet boundary conditions only allowed u-velocity and v-velocity to be set to the same value (which is enforced by setting the value of parameters u_face_mask or v_face_mask to 3). This functionality is retained here, but now setting u_face_mask or v_face_mask to 5 enforces a Dirichlet boundary for u-velocity only along the respective cell face. Similarly, setting u_face_mask or v_face_mask to 6 enforces the Dirichlet boundary for v-velocity only along the respective cell face. This functionality is required for most ice-sheet modeling configurations, e.g. the idealized MISMIP+ configuration requires setting v-velocity only to 0 along its lateral boundaries, but with free-slip conditions enforced for u-velocity. Adding this capability required changes throughout the ice-shelf code. Further changes were needed for how driving stress and Neummann conditions at computational boundaries are calculated in subroutine calc_shelf_driving_stress, and calls to subroutine apply_boundary_values were eliminated because they were not justified and caused errors. The new boundary treatment was tested by comparing simulated and analytical solutions for 1-D ice shelf flow with free-slip lateral boundary conditions and positive u_velocity enforced at the western boundary. This required the addition of a new parameter ADVECT_SHELF, which if false (as in the 1-D test case), turns off ice-shelf thickness evolution. By default, ADVECT_SHELF=True. The new boundary treatment was also justified in 2-D by correctly simulating the expected MISMIP+ steady-state.
The codecov token was added to the tc* test upload, but not the unit test upload. This patch adds the token to the unit testing.
Corrected the non-Boussinesq calculation of the total depth used by the v-component of the Flather open boundary condition, making the v-component consistent with the u-component and correcting an oversight with a recent commit. This commit could change answers in some non-Boussinesq cases with Flather open boundary conditions.
Previously, when ice_viscosity_compute == MODEL, ice shelf viscosity was calculated using 1 quadrature point per cell; however, this quadrature point was not centered in the cell. -Added a routine bilinear_shape_fn_grid_1qp, which is used to instead calculate viscosity for a single cell-centered quadrature point -Added an option to define parameter ice_viscosity_compute = MODEL_QUADRATURE, where viscosity is calculated at the same four quadrature points used during the SSA solution (the typical approach in finite element codes). Note that when using one quadrature point (ice_viscosity_compute == MODEL), array ice_visc(:,:) is the cell-centered ice viscosity. When using four quadrature points (ice_viscosity_compute == MODEL_QUADRATURE), ice_visc(:,:) is the cell-centered ice viscosity divided by the effective stress, Ee, which varies between each quadrature point and is saved in a separate array CS%Ee(:,:,:). In the SSA, ice viscosity is calculated for a cell with indices i,j as ice_visc(i,j) * Ee, where Ee=1 if using one quadrature point for viscosity and Ee=CS%Ee(i,j,k) if using four quad points (where k the current quad point). Also note the post_data call for ice_visc when using four quad points, where Ice_visc is outputted for visualization from a cell as ice_visc(i,j)*Ee_av(i,j), where Ee_av(i,j) is the average Ee in the cell.
Added an option NONLIN_SOLVE_ERR_MODE=3 to check for convergence of the ice shelf SSA solution by testing the change of norm, i.e. 2*abs(|u_{t}|-|u_{t-1}|) / (|u_{t}|+|u_{t-1}|)
Use thickness_to_dz to convert thicknesses from thickness units to height units in dumbbell_initialize_sponges with the traditional (non-ALE) sponges. Boussinesq answers are identical, but non-Boussinesq answers with an equation of state will change to be less dependent on the value of RHO_0.
* Fix nudged OBCs for tracers * Fix index in tracer reservoirs * fix index in the function update_segment_tracer_reservoirs * Fix index bugs in OBC tracer reservoirs * Fix tracer index bugs at open boundaries
* +Add particle code option to advect with uhtr The particle code has so far used the same velocity has was used in the dynamics step. I would like to add the option for the particle code to use uhtr/h and vhtr/h, so that the velocities used to advect particles may include the effects of parameterized eddies. To make this work, I have added a flag that controls which velocity to use and moved the particles_run step to take place after uhtr and vhtr are defined. The interfaces in the code and in config_src/external are updated to pass this information to the drifters package.
* Added ice shelf Coulomb friction law (Schoof 2005, Gagliardini et al 2007) needed for MISMIP+ experiments (Asay-Davis et al 2016).
* +REMAP_AUX needs at least one more halo update. - This one is for CS%u_av, CS%v_av, which need to be updated coming into step_MOM_dyn_split_RK2. * +Next stab at fixing REMAP_AUX fallout. - This fixes the Bering ORLANSKI OBCs for differing processor counts. - This is either the wrong way to do group_pass for OBLIQUE OBC's or there is more wrong with them. * Adding a group pass, still not solving the problem - Problem is in tangential_vel at tile boundaries. It matches right at the boundary, but needs some halo points to match too. * +Fixing oblique OBCs - Without this, u_av and v_av don't update a wide enough halo to get answers to reproduce across different processor counts with oblique OBCs. * Fixed an oopsie with OBC * Getting rid of extra exchange (that didn't help)
Refactored diapyc_energy_req_test and diapyc_energy_req_calc to remove the dependence on the Boussinesq reference density when in non-Boussinesq mode. This includes changes to the scaled units of the Kd_int argument to diapyc_energy_req_calc and the Kd argument to diapyc_energy_req_calc and the addition of a new argument to diapyc_energy_req_calc. A call to thickness_to_dz is used for the thickness unit conversions. There are 5 new internal variables, and changes to the units of several others. These routines are not actively used in MOM6 solutions, but instead they are used for testing and debugging new code, so there are no changes to solutions, but the results of these routines can differ in fully non-Boussinesq mode.
Modified the MOM_diapyc_energy_req.F90 version of find_PE_chg to align more closely with the version in MOM_energetic_PBL.F90, including making PE_chg into a mandatory argument, changing the name of the ColHt_cor argument to PE_ColHt_cor, and modifying some variable descriptions in units. Also removed find_PE_chg_orig from MOM_diapyc_energy_req.F90 and the old_PE_calc code that calls it. Extra values were also added to Te, Te_a and Te_b and the equivalent salinity variables so that the logical branches at (K==2) and (K=nz) could be simplied out of diapyc_energy_req_calc. Because old_PE_calc had been hard-coded to .false., all answers are bitwise identical.
…equired for NW2 (#484) * Update of Zanna-Bolton-2020 closure: code optimization and features required in NW2 configuration * Resolving compilation errors and doxygen by Alistair * Remove force sync for clock * Change naming of functions according to MOM_Zanna_bolton module
- Several edits to the ice shelf solo driver so that it works with the rest of the current MOM6 - Added capability to initialize a surface mass balance (SMB) that is held contstant over time when running from the ice-shelf solo driver (see new subroutine initialize_ice_SMB). This is required for MISMIP+. A constant SMB can also be used from the MOM driver for coupled ice-shelf/ocean experiments (e.g. MISOMIP). - The new, constant SMB is passed into solo_step_ice_shelf, where change_thickness_using_precip is called - Added capability to save both non-time-stamped and time-stamped restart files when using the ice shelf solo driver. This is useful for debugging. - slight reorganization to when ice shelf post_data calls are made - Added safety checks to diag_mediator_end() so that it works with the ice shelf solo-driver, which now calls it instead of (now removed) solo_ice_shelf_diag_mediator_end() routine. Removed the runtime parameter SAVE_BOTH_RESTARTS from the ice shelf solo-driver, which is no longer needed.
-fixed a bug in change_thickness_using_precip (was missing a division by ice density) -optimized ice shelf pass_var calls with optional complete arguments -corrected the grid area to multiply with ice shelf driving stress before its post_data call -changed some order of operations by adding parentheses, with the hope that it would improve symmetry of the ice shelf solution during MISMIP+. There was no effect, but this version of the code was used for MISMIP+ and MISOMIP.
|
wrong repo... 🤦 |
Codecov Report
@@ Coverage Diff @@
## main #1612 +/- ##
==========================================
- Coverage 37.91% 37.54% -0.38%
==========================================
Files 269 270 +1
Lines 77155 79030 +1875
Branches 14164 14633 +469
==========================================
+ Hits 29254 29672 +418
- Misses 42623 43892 +1269
- Partials 5278 5466 +188
... and 1 file with indirect coverage changes 📣 We’re building smart automated test selection to slash your CI/CD build times. Learn more |
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.
Merged branch from
main, includes updates from NCAR.Gaea regression: https://gitlab.gfdl.noaa.gov/ogrp/MOM6/-/pipelines/21103 ✔️ 🟡
Do not merge this PR, it will be merged manually. This is for documentation purposes.