mom_lateral_mixing_coeffs module reference

Variable mixing coefficients.

More…

Data Types

varmix_cs

Variable mixing coefficients.

Functions/Subroutines

calc_depth_function()

Calculates the non-dimensional depth functions.

calc_resoln_function()

Calculates and stores the non-dimensional resolution functions.

calc_slope_functions()

Calculates and stores functions of isopycnal slopes, e.g.

calc_visbeck_coeffs_old()

Calculates factors used when setting diffusivity coefficients similar to Visbeck et al., 1997.

calc_eady_growth_rate_2d()

Calculates the Eady growth rate (2D fields) for use in MEKE and the Visbeck schemes.

calc_slope_functions_using_just_e()

The original calc_slope_function() that calculated slopes using interface positions only, not accounting for density variations.

calc_qg_leith_viscosity()

Calculates the Leith Laplacian and bi-harmonic viscosity coefficients.

varmix_init()

Initializes the variables mixing coefficients container.

varmix_end()

Destructor for VarMix control structure.

Detailed Description

This module provides a container for various factors used in prescribing diffusivities, that are a function of the state (in particular the stratification and isoneutral slopes).

The resolution function

The resolution function is expressed in terms of the ratio of grid-spacing to deformation radius. The square of the resolution parameter is

\[R^2 = \frac{L_d^2}{\Delta^2} = \frac{ c_g^2 }{ f^2 \Delta^2 + c_g \beta \Delta^2 }\]

where the grid spacing is calculated as

\[\Delta^2 = \Delta x^2 + \Delta y^2 .\]

Todo

Check this reference to Bob on/off paper. The resolution function used in scaling diffusivities (Hallberg, 2010) is

\[r(\Delta,L_d) = \frac{1}{1+(\alpha R)^p}\]

The resolution function can be applied independently to thickness diffusion (module mom_thickness_diffuse()), tracer diffusion (mom_tracer_hordiff) lateral viscosity (), tracer diffusion (mom_tracer_hordiff) lateral viscosity (mom_hor_visc()).).

Robert Hallberg, 2013: Using a resolution function to regulate parameterizations of oceanic mesoscale eddy effects. Ocean Modelling, 71, pp 92-103. http://dx.doi.org/10.1016/j.ocemod.2013.08.007

Symbol

Module parameter

USE_VARIABLE_MIXING

RESOLN_SCALED_KH

RESOLN_SCALED_KHTH

RESOLN_SCALED_KHTR

\(\alpha\)

KH_RES_SCALE_COEF (for thickness and tracer diffusivity)

\(p\)

KH_RES_FN_POWER (for thickness and tracer diffusivity)

\(\alpha\)

VISC_RES_SCALE_COEF (for lateral viscosity)

\(p\)

VISC_RES_FN_POWER (for lateral viscosity)

GILL_EQUATORIAL_LD

Visbeck diffusivity

This module also calculates factors used in setting the thickness diffusivity similar to a Visbeck et al., 1997, scheme. The factors are combined in mom_thickness_diffuse::thickness_diffuse() but calculated in this module.but calculated in this module.

\[\kappa_h = \alpha_s L_s^2 S N\]

where \(S\) is the magnitude of the isoneutral slope and \(N\) is the Brunt-Vaisala frequency.

Visbeck, Marshall, Haine and Spall, 1997: Specification of Eddy Transfer Coefficients in Coarse-Resolution Ocean Circulation Models. J. Phys. Oceanogr. http://dx.doi.org/10.1175/1520-0485(1997)027%3C0381:SOETCI%3E2.0.CO;2

Symbol

Module parameter

USE_VARIABLE_MIXING

\(\alpha_s\)

KHTH_SLOPE_CFF (for mom_thickness_diffuse() module) module)

\(\alpha_s\)

KHTR_SLOPE_CFF (for mom_tracer_hordiff module)

\(L_{s}\)

VISBECK_L_SCALE

\(S_{max}\)

VISBECK_MAX_SLOPE

Vertical structure function for KhTh

The thickness diffusivity can be prescribed a vertical distribution with the shape of the equivalent barotropic velocity mode. The structure function is stored in the control structure for thie module (varmix_cs()) but is calculated using subroutines in ) but is calculated using subroutines in mom_wave_speed()..

Symbol

Module parameter

KHTH_USE_EBT_STRUCT

Type Documentation

type mom_lateral_mixing_coeffs/varmix_cs

Variable mixing coefficients.

Type fields:
  • % id_sn_u [integer] :: Diagnostic identifier.

  • % id_sn_v [integer] :: Diagnostic identifier.

  • % id_l2u [integer] :: Diagnostic identifier.

  • % id_l2v [integer] :: Diagnostic identifier.

  • % id_res_fn [integer] :: Diagnostic identifier.

  • % id_n2_u [integer] :: Diagnostic identifier.

  • % id_n2_v [integer] :: Diagnostic identifier.

  • % id_s2_u [integer] :: Diagnostic identifier.

  • % id_s2_v [integer] :: Diagnostic identifier.

  • % id_dzu [integer] :: Diagnostic identifier.

  • % id_dzv [integer] :: Diagnostic identifier.

  • % id_dzsxn [integer] :: Diagnostic identifier.

  • % id_dzsyn [integer] :: Diagnostic identifier.

  • % id_rd_dx [integer] :: Diagnostic identifier.

  • % id_kh_u_qg [integer] :: Diagnostic identifier.

  • % id_kh_v_qg [integer] :: Diagnostic identifier.

  • % diag [type( diag_ctrl ),pointer] :: A structure that is used to regulate the timing of diagnostic output.

  • % initialized [logical] :: True if this control structure has been initialized.

  • % use_variable_mixing [logical] :: If true, use the variable mixing.

  • % resoln_scaling_used [logical] :: If true, a resolution function is used somewhere to scale away one of the viscosities or diffusivities when the deformation radius is well resolved.

  • % resoln_scaled_kh [logical] :: If true, scale away the Laplacian viscosity when the deformation radius is well resolved.

  • % resoln_scaled_khth [logical] :: If true, scale away the thickness diffusivity when the deformation radius is well resolved.

  • % depth_scaled_khth [logical] :: If true, KHTH is scaled away when the depth is shallower than a reference depth.

  • % resoln_scaled_khtr [logical] :: If true, scale away the tracer diffusivity when the deformation radius is well resolved.

  • % interpolate_res_fn [logical] :: If true, interpolate the resolution function to the velocity points from the thickness points; otherwise interpolate the wave speed and calculate the resolution function independently at each point.

  • % use_stored_slopes [logical] :: If true, stores isopycnal slopes in this structure.

  • % resoln_use_ebt [logical] :: If true, uses the equivalent barotropic wave speed instead of first baroclinic wave for calculating the resolution fn.

  • % khth_use_ebt_struct [logical] :: If true, uses the equivalent barotropic structure as the vertical structure of thickness diffusivity.

  • % kdgl90_use_ebt_struct [logical] :: If true, uses the equivalent barotropic structure as the vertical structure of diffusivity in the GL90 scheme.

  • % calculate_cg1 [logical] :: If true, calls

  • % calculate_rd_dx [logical] :: If true, calculates Rd/dx and populate CSRd_dx_h. This parameter is set depending on other parameters.

  • % calculate_res_fns [logical] :: If true, calculate all the resolution factors. This parameter is set depending on other parameters.

  • % calculate_depth_fns [logical] :: If true, calculate all the depth factors. This parameter is set depending on other parameters.

  • % calculate_eady_growth_rate [logical] :: If true, calculate all the Eady growth rates. This parameter is set depending on other parameters.

  • % use_stanley_iso [logical] :: If true, use Stanley parameterization in MOM_isopycnal_slopes.

  • % use_simpler_eady_growth_rate [logical] :: If true, use a simpler method to calculate the Eady growth rate that avoids division by layer thickness. This parameter is set depending on other parameters.

  • % full_depth_eady_growth_rate [logical] :: If true, calculate the Eady growth rate based on an average that includes contributions from sea-level changes in its denominator, rather than just the nominal depth of the bathymetry. This only applies when using the model interface heights as a proxy for isopycnal slopes.

  • % cropping_distance [real] :: Distance from surface or bottom to filter out outcropped or incropped interfaces for the Eady growth rate calc [Z ~> m].

  • % h_min_n2 [real] :: The minimum vertical distance to use in the denominator of the bouyancy frequency used in the slope calculation [H ~> m or kg m-2].

  • % sn_u [real(:,:),allocatable] :: S*N at u-points [T-1 ~> s-1].

  • % sn_v [real(:,:),allocatable] :: S*N at v-points [T-1 ~> s-1].

  • % l2u [real(:,:),allocatable] :: Length scale^2 at u-points [L2 ~> m2].

  • % l2v [real(:,:),allocatable] :: Length scale^2 at v-points [L2 ~> m2].

  • % cg1 [real(:,:),allocatable] :: The first baroclinic gravity wave speed [L T-1 ~> m s-1].

  • % res_fn_h [real(:,:),allocatable] :: Non-dimensional function of the ratio the first baroclinic deformation radius to the grid spacing at h points [nondim].

  • % res_fn_q [real(:,:),allocatable] :: Non-dimensional function of the ratio the first baroclinic deformation radius to the grid spacing at q points [nondim].

  • % res_fn_u [real(:,:),allocatable] :: Non-dimensional function of the ratio the first baroclinic deformation radius to the grid spacing at u points [nondim].

  • % res_fn_v [real(:,:),allocatable] :: Non-dimensional function of the ratio the first baroclinic deformation radius to the grid spacing at v points [nondim].

  • % depth_fn_u [real(:,:),allocatable] :: Non-dimensional function of the ratio of the depth to a reference depth (maximum 1) at u points [nondim].

  • % depth_fn_v [real(:,:),allocatable] :: Non-dimensional function of the ratio of the depth to a reference depth (maximum 1) at v points [nondim].

  • % beta_dx2_h [real(:,:),allocatable] :: The magnitude of the gradient of the Coriolis parameter times the grid spacing squared at h points [L T-1 ~> m s-1].

  • % beta_dx2_q [real(:,:),allocatable] :: The magnitude of the gradient of the Coriolis parameter times the grid spacing squared at q points [L T-1 ~> m s-1].

  • % beta_dx2_u [real(:,:),allocatable] :: The magnitude of the gradient of the Coriolis parameter times the grid spacing squared at u points [L T-1 ~> m s-1].

  • % beta_dx2_v [real(:,:),allocatable] :: The magnitude of the gradient of the Coriolis parameter times the grid spacing squared at v points [L T-1 ~> m s-1].

  • % f2_dx2_h [real(:,:),allocatable] :: The Coriolis parameter squared times the grid spacing squared at h [L2 T-2 ~> m2 s-2].

  • % f2_dx2_q [real(:,:),allocatable] :: The Coriolis parameter squared times the grid spacing squared at q [L2 T-2 ~> m2 s-2].

  • % f2_dx2_u [real(:,:),allocatable] :: The Coriolis parameter squared times the grid spacing squared at u [L2 T-2 ~> m2 s-2].

  • % f2_dx2_v [real(:,:),allocatable] :: The Coriolis parameter squared times the grid spacing squared at v [L2 T-2 ~> m2 s-2].

  • % rd_dx_h [real(:,:),allocatable] :: Deformation radius over grid spacing [nondim].

  • % slope_x [real(:,:,:),allocatable] :: Zonal isopycnal slope [Z L-1 ~> nondim].

  • % slope_y [real(:,:,:),allocatable] :: Meridional isopycnal slope [Z L-1 ~> nondim].

  • % ebt_struct [real(:,:,:),allocatable] :: Vertical structure function to scale diffusivities with [nondim].

  • % real (* kh_v_qg [*, *) :: Laplacian metric-dependent constants [L3 ~> m3].

  • % real :: Laplacian metric-dependent constants [L3 ~> m3].

  • % real :: QG Leith GM coefficient at u-points [L2 T-1 ~> m2 s-1].

  • % real :: QG Leith GM coefficient at v-points [L2 T-1 ~> m2 s-1].

  • % use_visbeck [logical] :: Use Visbeck formulation for thickness diffusivity.

  • % varmix_ktop [integer] :: Top layer to start downward integrals.

  • % visbeck_l_scale [real] :: Fixed length scale in Visbeck formula [L ~> m], or if negative a scaling factor [nondim] relating this length scale squared to the cell area.

  • % eady_gr_d_scale [real] :: Depth over which to average SN [Z ~> m].

  • % res_coef_khth [real] :: A coefficient [nondim] that determines the function of resolution, used for thickness and tracer mixing, as: F = 1 / (1 + (Res_coef_khth*Ld/dx)^Res_fn_power)

  • % res_coef_visc [real] :: A coefficient [nondim] that determines the function of resolution, used for lateral viscosity, as: F = 1 / (1 + (Res_coef_visc*Ld/dx)^Res_fn_power)

  • % depth_scaled_khth_h0 [real] :: The depth above which KHTH is linearly scaled away [Z ~> m].

  • % depth_scaled_khth_exp [real] :: The exponent used in the depth dependent scaling function for KHTH [nondim].

  • % kappa_smooth [real] :: A diffusivity for smoothing T/S in vanished layers [H Z T-1 ~> m2 s-1 or kg m-1 s-1].

  • % res_fn_power_khth [integer] :: The power of dx/Ld in the KhTh resolution function. Any positive integer power may be used, but even powers and especially 2 are coded to be more efficient.

  • % res_fn_power_visc [integer] :: The power of dx/Ld in the Kh resolution function. Any positive integer power may be used, but even powers and especially 2 are coded to be more efficient.

  • % visbeck_s_max [real] :: Upper bound on slope used in Eady growth rate [Z L-1 ~> nondim].

  • % use_qg_leith_gm [logical] :: If true, uses the QG Leith viscosity as the GM coefficient.

  • % use_beta_in_qg_leith [logical] :: If true, includes the beta term in the QG Leith GM coefficient.

  • % wave_speed [type( wave_speed_cs )] :: Wave speed control structure.

  • % pass_cg1 [type(group_pass_type)] :: For group halo pass.

  • % debug [logical] :: If true, write out checksums of data for debugging.

Function/Subroutine Documentation

subroutine mom_lateral_mixing_coeffs/calc_depth_function(G, CS)

Calculates the non-dimensional depth functions.

Parameters:
  • g :: [in] Ocean grid structure

  • cs :: [inout] Variable mixing control structure

Call to:

mom_error_handler::mom_error

Called from:

mom::step_mom mom::step_offline

subroutine mom_lateral_mixing_coeffs/calc_resoln_function(h, tv, G, GV, US, CS)

Calculates and stores the non-dimensional resolution functions.

Parameters:
  • g :: [inout] Ocean grid structure

  • gv :: [in] Vertical grid structure

  • h :: [in] Layer thickness [H ~> m or kg m-2]

  • tv :: [in] Thermodynamic variables

  • us :: [in] A dimensional unit scaling type

  • cs :: [inout] Variable mixing control structure

Call to:

mom_error_handler::mom_error mom_diag_mediator::query_averaging_enabled mom_wave_speed::wave_speed

Called from:

mom::step_mom mom::step_offline

subroutine mom_lateral_mixing_coeffs/calc_slope_functions(h, tv, dt, G, GV, US, CS, OBC)

Calculates and stores functions of isopycnal slopes, e.g. Sx, Sy, S*N, mostly used in the Visbeck et al. style scaling of diffusivity.

Parameters:
  • g :: [inout] Ocean grid structure

  • gv :: [in] Vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • h :: [inout] Layer thickness [H ~> m or kg m-2]

  • tv :: [in] Thermodynamic variables

  • dt :: [in] Time increment [T ~> s]

  • cs :: [inout] Variable mixing control structure

  • obc :: Open boundaries control structure

Call to:

calc_eady_growth_rate_2d mom_isopycnal_slopes::calc_isoneutral_slopes calc_slope_functions_using_just_e calc_visbeck_coeffs_old mom_error_handler::mom_error mom_diag_mediator::query_averaging_enabled

subroutine mom_lateral_mixing_coeffs/calc_visbeck_coeffs_old(h, slope_x, slope_y, N2_u, N2_v, G, GV, US, CS)

Calculates factors used when setting diffusivity coefficients similar to Visbeck et al., 1997. This is on older implementation that is susceptible to large values of Eady growth rate for incropping layers.

Parameters:
  • g :: [inout] Ocean grid structure

  • gv :: [in] Vertical grid structure

  • h :: [in] Layer thickness [H ~> m or kg m-2]

  • slope_x :: [in] Zonal isoneutral slope [Z L-1 ~> nondim]

  • n2_u :: [in] Buoyancy (Brunt-Vaisala) frequency at u-points [L2 Z-2 T-2 ~> s-2]

  • slope_y :: [in] Meridional isoneutral slope [Z L-1 ~> nondim]

  • n2_v :: [in] Buoyancy (Brunt-Vaisala) frequency at v-points [L2 Z-2 T-2 ~> s-2]

  • us :: [in] A dimensional unit scaling type

  • cs :: [inout] Variable mixing control structure

Call to:

mom_error_handler::mom_error mom_diag_mediator::query_averaging_enabled

Called from:

calc_slope_functions

subroutine mom_lateral_mixing_coeffs/calc_eady_growth_rate_2d(CS, G, GV, US, h, e, dzu, dzv, dzSxN, dzSyN, SN_u, SN_v)

Calculates the Eady growth rate (2D fields) for use in MEKE and the Visbeck schemes.

Parameters:
  • cs :: [inout] Variable mixing coefficients

  • g :: [in] Ocean grid structure

  • gv :: [in] Vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • h :: [in] Interface height [Z ~> m]

  • e :: [in] Interface height [Z ~> m]

  • dzu :: [in] dz at u-points [Z ~> m]

  • dzv :: [in] dz at v-points [Z ~> m]

  • dzsxn :: [in] dz Sx N at u-points [Z T-1 ~> m s-1]

  • dzsyn :: [in] dz Sy N at v-points [Z T-1 ~> m s-1]

  • sn_u :: [inout] SN at u-points [T-1 ~> s-1]

  • sn_v :: [inout] SN at v-points [T-1 ~> s-1]

Called from:

calc_slope_functions

subroutine mom_lateral_mixing_coeffs/calc_slope_functions_using_just_e(h, G, GV, US, CS, e, calculate_slopes)

The original calc_slope_function() that calculated slopes using interface positions only, not accounting for density variations.

Parameters:
  • g :: [inout] Ocean grid structure

  • gv :: [in] Vertical grid structure

  • h :: [inout] Layer thickness [H ~> m or kg m-2]

  • us :: [in] A dimensional unit scaling type

  • cs :: [inout] Variable mixing control structure

  • e :: [in] Interface position [Z ~> m]

  • calculate_slopes :: [in] If true, calculate slopes internally otherwise use slopes stored in CS

Call to:

mom_error_handler::mom_error

Called from:

calc_slope_functions

subroutine mom_lateral_mixing_coeffs/calc_qg_leith_viscosity(CS, G, GV, US, h, k, div_xx_dx, div_xx_dy, vort_xy_dx, vort_xy_dy)

Calculates the Leith Laplacian and bi-harmonic viscosity coefficients.

Parameters:
  • cs :: [inout] Variable mixing coefficients

  • g :: [in] Ocean grid structure

  • gv :: [in] The ocean’s vertical grid structure.

  • us :: [in] A dimensional unit scaling type

  • h :: [inout] Layer thickness [H ~> m or kg m-2]

  • k :: [in] Layer for which to calculate vorticity magnitude

  • div_xx_dx :: [in] x-derivative of horizontal divergence (d/dx(du/dx + dv/dy)) [L-1 T-1 ~> m-1 s-1]

  • div_xx_dy :: [in] y-derivative of horizontal divergence (d/dy(du/dx + dv/dy)) [L-1 T-1 ~> m-1 s-1]

  • vort_xy_dx :: [inout] x-derivative of vertical vorticity (d/dx(dv/dx - du/dy)) [L-1 T-1 ~> m-1 s-1]

  • vort_xy_dy :: [inout] y-derivative of vertical vorticity (d/dy(dv/dx - du/dy)) [L-1 T-1 ~> m-1 s-1]

Called from:

mom_hor_visc::horizontal_viscosity

subroutine mom_lateral_mixing_coeffs/varmix_init(Time, G, GV, US, param_file, diag, CS)

Initializes the variables mixing coefficients container.

Parameters:
  • time :: [in] Current model time

  • g :: [in] Ocean grid structure

  • gv :: [in] The ocean’s vertical grid structure

  • us :: [in] A dimensional unit scaling type

  • param_file :: [in] Parameter file handles

  • diag :: [inout] Diagnostics control structure

  • cs :: [inout] Variable mixing coefficients

Call to:

mom_error_handler::mom_error mom_wave_speed::wave_speed_init

subroutine mom_lateral_mixing_coeffs/varmix_end(CS)

Destructor for VarMix control structure.