diff --git a/docs/zotero.bib b/docs/zotero.bib
index 01fe2c6185..eb2423de02 100644
--- a/docs/zotero.bib
+++ b/docs/zotero.bib
@@ -1,5 +1,34 @@
% This file is generated by zotero. Manual edits will be lost!
+@article{agarwal2023,
+ title={Impact of Stochastic Ocean Density Corrections on Air-Sea Flux Variability},
+ author={Agarwal, Niraj and Small, R Justin and Bryan, Frank O and Grooms, Ian and Pegion, Philip J},
+ journal={Geophysical Research Letters},
+ volume={50},
+ number={13},
+ pages={e2023GL104248},
+ year={2023},
+}
+
+@article{kenigson2022,
+ title={Parameterizing the impact of unresolved temperature variability on the large-scale density field: 2. Modeling},
+ author={Kenigson, JS and Adcroft, A and Bachman, SD and Castruccio, F and Grooms, I and Pegion, P and Stanley, Z},
+ journal={Journal of Advances in Modeling Earth Systems},
+ volume={14},
+ number={3},
+ pages={e2021MS002844},
+ year={2022},
+}
+
+@article{stanley2020,
+ title={Parameterizing the Impact of Unresolved Temperature Variability on the Large-Scale Density Field: {Part 1.} Theory.},
+ author={Stanley, Z and Grooms, I and Kleiber, W and Bachman, SD and Castruccio, F and Adcroft, A},
+ journal={Journal of Advances in Modeling Earth Systems},
+ volume={12},
+ number={12},
+ pages={e2020MS002185},
+ year={2020},
+}
@article{redi1982,
title = {Oceanic {Isopycnal} {Mixing} by {Coordinate} {Rotation}},
diff --git a/src/core/MOM_PressureForce_FV.F90 b/src/core/MOM_PressureForce_FV.F90
index aaabab3500..e55376d568 100644
--- a/src/core/MOM_PressureForce_FV.F90
+++ b/src/core/MOM_PressureForce_FV.F90
@@ -2030,8 +2030,6 @@ subroutine PressureForce_FV_init(Time, G, GV, US, param_file, diag, CS, ADp, SAL
type(tidal_forcing_CS), intent(in), target, optional :: tides_CSp !< Tides control structure
! Local variables
- real :: Stanley_coeff ! Coefficient relating the temperature gradient and sub-gridscale
- ! temperature variance [nondim]
integer :: default_answer_date ! Global answer date
logical :: use_temperature ! If true, temperature and salinity are used as state variables.
logical :: use_EOS ! If true, density calculated from T & S using an equation of state.
@@ -2043,6 +2041,7 @@ subroutine PressureForce_FV_init(Time, G, GV, US, param_file, diag, CS, ADp, SAL
# include "version_variable.h"
character(len=40) :: mdl ! This module's name.
logical :: use_ALE ! If true, use the Vertical Lagrangian Remap algorithm
+ logical :: stoch_eos ! Can't use Stanley param here unless stoch_eos is true
integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, nz
isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed ; nz = GV%ke
@@ -2182,16 +2181,15 @@ subroutine PressureForce_FV_init(Time, G, GV, US, param_file, diag, CS, ADp, SAL
"boundary cells is extrapolated, rather than using PCM "//&
"in these cells. If true, the same order polynomial is "//&
"used as is used for the interior cells.", default=.true.)
+ call get_param(param_file, mdl, "STOCH_EOS", stoch_eos, &
+ default=.false., do_not_log=.true.)
call get_param(param_file, mdl, "USE_STANLEY_PGF", CS%use_stanley_pgf, &
"If true, turn on Stanley SGS T variance parameterization "// &
"in PGF code.", default=.false.)
if (CS%use_stanley_pgf) then
- call get_param(param_file, mdl, "STANLEY_COEFF", Stanley_coeff, &
- "Coefficient correlating the temperature gradient and SGS T variance.", &
- units="nondim", default=-1.0, do_not_log=.true.)
- if (Stanley_coeff < 0.0) call MOM_error(FATAL, &
- "STANLEY_COEFF must be set >= 0 if USE_STANLEY_PGF is true.")
-
+ if (.not.stoch_eos) then
+ call MOM_error(FATAL, "PressureForce_FV_init: USE_STANLEY_PGF requires STOCH_EOS")
+ endif
CS%id_rho_pgf = register_diag_field('ocean_model', 'rho_pgf', diag%axesTL, &
Time, 'rho in PGF', 'kg m-3', conversion=US%R_to_kg_m3)
CS%id_rho_stanley_pgf = register_diag_field('ocean_model', 'rho_stanley_pgf', diag%axesTL, &
diff --git a/src/core/MOM_stoch_eos.F90 b/src/core/MOM_stoch_eos.F90
index 909c2e9a6a..8db35c72e7 100644
--- a/src/core/MOM_stoch_eos.F90
+++ b/src/core/MOM_stoch_eos.F90
@@ -4,7 +4,7 @@ module MOM_stoch_eos
! This file is part of MOM6. See LICENSE.md for the license.
use MOM_diag_mediator, only : register_diag_field, post_data, diag_ctrl
use MOM_error_handler, only : MOM_error, FATAL
-use MOM_file_parser, only : get_param, param_file_type
+use MOM_file_parser, only : get_param, log_version, param_file_type
use MOM_grid, only : ocean_grid_type
use MOM_hor_index, only : hor_index_type
use MOM_isopycnal_slopes, only : vert_fill_TS
@@ -62,27 +62,32 @@ logical function MOM_stoch_eos_init(Time, G, GV, US, param_file, diag, CS, resta
type(MOM_restart_CS), pointer :: restart_CS !< A pointer to the restart control structure.
! local variables
+ ! This include declares and sets the variable "version".
+# include "version_variable.h"
integer :: i,j
MOM_stoch_eos_init = .false.
CS%seed = 0
+ call log_version(param_file, "MOM_stoch_eos", version, "")
call get_param(param_file, "MOM_stoch_eos", "STOCH_EOS", CS%use_stoch_eos, &
- "If true, stochastic perturbations are applied "//&
- "to the EOS in the PGF.", default=.false.)
+ "If true, computes stochastic perturbations that can be applied "//&
+ "to the EOS in various places.", default=.false.)
call get_param(param_file, "MOM_stoch_eos", "STANLEY_COEFF", CS%stanley_coeff, &
"Coefficient correlating the temperature gradient "//&
"and SGS T variance.", units="nondim", default=-1.0)
+ if ((CS%stanley_coeff < 0.0) .and. CS%use_stoch_eos) call MOM_error(FATAL, &
+ "STANLEY_COEFF must be set >= 0 if STOCH_EOS is true.")
call get_param(param_file, "MOM_stoch_eos", "STANLEY_A", CS%stanley_a, &
"Coefficient a which scales chi in stochastic perturbation of the "//&
"SGS T variance.", units="nondim", default=1.0, &
- do_not_log=((CS%stanley_coeff<0.0) .or. .not.CS%use_stoch_eos))
+ do_not_log=.not.CS%use_stoch_eos)
call get_param(param_file, "MOM_stoch_eos", "KD_SMOOTH", CS%kappa_smooth, &
"A diapycnal diffusivity that is used to interpolate "//&
"more sensible values of T & S into thin layers.", &
units="m2 s-1", default=1.0e-6, scale=GV%m2_s_to_HZ_T, &
- do_not_log=(CS%stanley_coeff<0.0))
+ do_not_log=.not.CS%use_stoch_eos)
! Don't run anything if STANLEY_COEFF < 0
if (CS%stanley_coeff >= 0.0) then
@@ -258,4 +263,30 @@ subroutine MOM_calc_varT(G, GV, US, h, tv, CS, dt)
endif
end subroutine MOM_calc_varT
+!> \namespace mom_stoch_eos
+!!
+!! This module provides the foundation of the Stanley parameterization (\cite stanley2020) for correcting the
+!! computation of density. Density is not a prognostic variable in MOM6; it is computed for various purposes
+!! in various places. The correction to this calculation provided by this module has been implemented
+!! in some places where density is used, but not all.
+!!
+!! To use the correction, first set STOCH_EOS=True. Then, choose the constant c from (25) of
+!! \cite stanley2020. This is controlled using STANLEY_COEFF. Setting a negative value will
+!! result in an error. \cite stanley2020 found a value of 0.2 offline, coarsening from 0.1 to 1 degree
+!! resolution. \cite kenigson2022 proposed a value of 0.5 in a 2/3 degree resolution model.
+!!
+!! Whether the correction is deterministic or stochastic can be controlled using the variable
+!! STANLEY_A. Setting this to 0.0 uses the deterministic version, while a value of 1.0 produces
+!! the stochastic version. Reducing from 1 to 0 smoothly transitions from stochastic to deterministic.
+!!
+!! To turn the correction on in various parts of the code, use
+!! - USE_STANLEY_PGF=True for the pressure gradient force (cf. \cite kenigson2022)
+!! - USE_STANLEY_ISO=True to correct the computation of isopycnal slopes (used in many places)
+!! - USE_STANLEY_GM=True to use the parameterization within GM (cf. \cite agarwal2023)
+!! - USE_STANLEY_ML=True to use the parameterization within the mixed-layer restratification
+!! parameterization. It applies to both the OM4 and Bodner schemes. (cf. \cite agarwal2023)
+!!
+!! For ensemble simulations, the random number generator seed can be controlled using the parameter
+!! SEED_STOCH_EOS
+
end module MOM_stoch_eos
diff --git a/src/parameterizations/lateral/MOM_lateral_mixing_coeffs.F90 b/src/parameterizations/lateral/MOM_lateral_mixing_coeffs.F90
index f2f476b0c8..374cd2ce6c 100644
--- a/src/parameterizations/lateral/MOM_lateral_mixing_coeffs.F90
+++ b/src/parameterizations/lateral/MOM_lateral_mixing_coeffs.F90
@@ -1335,9 +1335,8 @@ subroutine VarMix_init(Time, G, GV, US, param_file, diag, CS)
! scaled by the resolution function.
logical :: better_speed_est ! If true, use a more robust estimate of the first
! mode wave speed as the starting point for iterations.
- real :: Stanley_coeff ! Coefficient relating the temperature gradient and sub-gridscale
- ! temperature variance [nondim]
logical :: om4_remap_via_sub_cells ! Use the OM4-era ramap_via_sub_cells for calculating the EBT structure
+ logical :: stoch_eos ! Can't use Stanley param here unless stoch_eos is true
! This include declares and sets the variable "version".
# include "version_variable.h"
character(len=40) :: mdl = "MOM_lateral_mixing_coeffs" ! This module's name.
@@ -1461,15 +1460,13 @@ subroutine VarMix_init(Time, G, GV, US, param_file, diag, CS)
call get_param(param_file, mdl, "DEBUG", CS%debug, default=.false., do_not_log=.true.)
+ call get_param(param_file, mdl, "STOCH_EOS", stoch_eos, &
+ default=.false., do_not_log=.true.)
call get_param(param_file, mdl, "USE_STANLEY_ISO", CS%use_stanley_iso, &
"If true, turn on Stanley SGS T variance parameterization "// &
"in isopycnal slope code.", default=.false.)
- if (CS%use_stanley_iso) then
- call get_param(param_file, mdl, "STANLEY_COEFF", Stanley_coeff, &
- "Coefficient correlating the temperature gradient and SGS T variance.", &
- units="nondim", default=-1.0, do_not_log=.true.)
- if (Stanley_coeff < 0.0) call MOM_error(FATAL, &
- "STANLEY_COEFF must be set >= 0 if USE_STANLEY_ISO is true.")
+ if (CS%use_Stanley_ISO .and. .not.stoch_eos) then
+ call MOM_error(FATAL, "VarMix_init: USE_STANLEY_ISO requires STOCH_EOS")
endif
if (CS%Resoln_use_ebt .or. CS%khth_use_ebt_struct .or. CS%kdgl90_use_ebt_struct &
diff --git a/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 b/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90
index 00b3e0e616..0f8fe55e73 100644
--- a/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90
+++ b/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90
@@ -863,6 +863,9 @@ subroutine mixedlayer_restrat_Bodner(CS, G, GV, US, h, uhtr, vhtr, tv, forces, d
"To use the Bodner et al., 2023, MLE parameterization, either MLE_USE_PBL_MLD or "// &
"Bodner_detect_MLD must be True.")
endif
+ if (CS%use_Stanley_ML .and. .not.GV%Boussinesq) call MOM_error(FATAL, &
+ "MOM_mixedlayer_restrat: The Stanley parameterization is not"//&
+ "available without the Boussinesq approximation.")
if (associated(bflux)) &
call pass_var(bflux, G%domain, halo=1)
@@ -1644,8 +1647,6 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS,
real :: flux_to_kg_per_s ! A unit conversion factor for fluxes. [kg T s-1 H-1 L-2 ~> kg m-3 or 1]
real :: omega ! The Earth's rotation rate [T-1 ~> s-1].
real :: ustar_min_dflt ! The default value for RESTRAT_USTAR_MIN [Z T-1 ~> m s-1]
- real :: Stanley_coeff ! Coefficient relating the temperature gradient and sub-gridscale
- ! temperature variance [nondim]
integer :: default_answer_date ! The default setting for the various ANSWER_DATE flags
! This include declares and sets the variable "version".
character(len=200) :: inputdir ! The directory where NetCDF input files
@@ -1653,6 +1654,7 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS,
character(len=128) :: mle_fl_file ! Data containing MLE front-length scale. Used
! when reading from file.
character(len=32) :: fl_varname ! Name of front-length scale variable in mle_fl_file.
+ logical :: stoch_eos ! Can't use Stanley param here unless stoch_eos is true
# include "version_variable.h"
integer :: i, j
@@ -1689,6 +1691,14 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS,
"This sets the default value for the various _ANSWER_DATE parameters.", &
default=99991231, do_not_log=.true.)
call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".")
+ call get_param(param_file, mdl, "STOCH_EOS", stoch_eos, &
+ default=.false., do_not_log=.true.)
+ call get_param(param_file, mdl, "USE_STANLEY_ML", CS%use_Stanley_ML, &
+ "If true, turn on Stanley SGS T variance parameterization "// &
+ "in ML restrat code.", default=.false.)
+ if (CS%use_Stanley_ML .and. .not.stoch_eos) then
+ call MOM_error(FATAL, "mixedlayer_restrat_init: USE_STANLEY_ML requires STOCH_EOS")
+ endif
call openParameterBlock(param_file,'MLE') ! Prepend MLE% to all parameters
if (GV%nkml==0) then
call get_param(param_file, mdl, "USE_BODNER23", CS%use_Bodner, &
@@ -1751,9 +1761,6 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS,
"Fraction by which to extend the mixed-layer restratification "//&
"depth used for a smoother stream function at the base of "//&
"the mixed-layer.", units="nondim", default=0.0)
- call get_param(param_file, mdl, "USE_STANLEY_TVAR", CS%use_Stanley_ML, &
- "If true, turn on Stanley SGS T variance parameterization "// &
- "in ML restrat code.", default=.false.)
call get_param(param_file, mdl, "USE_CR_GRID", CS%Cr_grid, &
"If true, read in a spatially varying Cr field.", default=.false.)
call get_param(param_file, mdl, "USE_MLD_GRID", CS%MLD_grid, &
@@ -1811,17 +1818,6 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS,
"geostrophic kinetic energy or 1 plus the square of the "//&
"grid spacing over the deformation radius, as detailed "//&
"by Fox-Kemper et al. (2011)", units="nondim", default=0.0)
- ! These parameters are only used in the OM4-era version of Fox-Kemper
- call get_param(param_file, mdl, "USE_STANLEY_ML", CS%use_Stanley_ML, &
- "If true, turn on Stanley SGS T variance parameterization "// &
- "in ML restrat code.", default=.false.)
- if (CS%use_Stanley_ML) then
- call get_param(param_file, mdl, "STANLEY_COEFF", Stanley_coeff, &
- "Coefficient correlating the temperature gradient and SGS T variance.", &
- units="nondim", default=-1.0, do_not_log=.true.)
- if (Stanley_coeff < 0.0) call MOM_error(FATAL, &
- "STANLEY_COEFF must be set >= 0 if USE_STANLEY_ML is true.")
- endif
call get_param(param_file, mdl, 'VON_KARMAN_CONST', CS%vonKar, &
'The value the von Karman constant as used for mixed layer viscosity.', &
units='nondim', default=0.41)
diff --git a/src/parameterizations/lateral/MOM_thickness_diffuse.F90 b/src/parameterizations/lateral/MOM_thickness_diffuse.F90
index cf06be4ed2..e07b95cdca 100644
--- a/src/parameterizations/lateral/MOM_thickness_diffuse.F90
+++ b/src/parameterizations/lateral/MOM_thickness_diffuse.F90
@@ -97,9 +97,6 @@ module MOM_thickness_diffuse
!! When this is true, it breaks rotational symmetry.
logical :: use_GM_work_bug !< If true, use the incorrect sign for the
!! top-level work tendency on the top layer.
- real :: Stanley_det_coeff !< The coefficient correlating SGS temperature variance with the mean
- !! temperature gradient in the deterministic part of the Stanley parameterization.
- !! Negative values disable the scheme. [nondim]
logical :: read_khth !< If true, read a file containing the spatially varying horizontal
!! isopycnal height diffusivity
logical :: use_stanley_gm !< If true, also use the Stanley parameterization in MOM_thickness_diffuse
@@ -2195,6 +2192,7 @@ subroutine thickness_diffuse_init(Time, G, GV, US, param_file, diag, CDp, CS)
! available.
logical :: use_meke = .false. ! If true, use the MEKE formulation for the thickness diffusivity.
integer :: default_answer_date ! The default setting for the various ANSWER_DATE flags.
+ logical :: stoch_eos ! Can't use Stanley param here unless stoch_eos is true
integer :: i, j
CS%initialized = .true.
@@ -2322,15 +2320,13 @@ subroutine thickness_diffuse_init(Time, G, GV, US, param_file, diag, CDp, CS)
"streamfunction formulation, expressed as a fraction of planetary "//&
"rotation, OMEGA. This should be tiny but non-zero to avoid degeneracy.", &
default=1.e-15, units="nondim", scale=US%Z_to_L, do_not_log=.not.CS%use_FGNV_streamfn)
+ call get_param(param_file, mdl, "STOCH_EOS", stoch_eos, &
+ default=.false., do_not_log=.true.)
call get_param(param_file, mdl, "USE_STANLEY_GM", CS%use_stanley_gm, &
"If true, turn on Stanley SGS T variance parameterization "// &
"in GM code.", default=.false.)
- if (CS%use_stanley_gm) then
- call get_param(param_file, mdl, "STANLEY_COEFF", Stanley_coeff, &
- "Coefficient correlating the temperature gradient and SGS T variance.", &
- units="nondim", default=-1.0, do_not_log=.true.)
- if (Stanley_coeff < 0.0) call MOM_error(FATAL, &
- "STANLEY_COEFF must be set >= 0 if USE_STANLEY_GM is true.")
+ if (CS%use_Stanley_GM .and. .not.stoch_eos) then
+ call MOM_error(FATAL, "thickness_diffuse_init: USE_STANLEY_GM requires STOCH_EOS")
endif
call get_param(param_file, mdl, "OMEGA", omega, &
"The rotation rate of the earth.", &