mom_meke module reference

Implements the Mesoscale Eddy Kinetic Energy framework with topographic beta effect included in computing beta in Rhines scale.

More…

Data Types

meke_cs

Control structure that contains MEKE parameters and diagnostics handles.

Functions/Subroutines

step_forward_meke()

Integrates forward-in-time the MEKE eddy energy equation.

meke_equilibrium()

Calculates the equilibrium solution where the source depends only on MEKE diffusivity and there is no lateral diffusion of MEKE.

meke_equilibrium_restoring()

meke_lengthscales()

Calculates the eddy mixing length scale and \(\gamma_b\) and \(\gamma_t\) functions that are ratios of either bottom or barotropic eddy energy to the column eddy energy, respectively.

meke_lengthscales_0d()

Calculates the eddy mixing length scale and \(\gamma_b\) and \(\gamma_t\) functions that are ratios of either bottom or barotropic eddy energy to the column eddy energy, respectively.

meke_init()

Initializes the MOM_MEKE module and reads parameters.

meke_alloc_register_restart()

Allocates memory and register restart fields for the MOM_MEKE module.

meke_end()

Deallocates any variables allocated in MEKE_alloc_register_restart.

Detailed Description

The Mesoscale Eddy Kinetic Energy (MEKE) framework

The MEKE framework accounts for the mean potential energy removed by the first order closures used to parameterize mesoscale eddies. It requires closure at the second order, namely dissipation and transport of eddy energy.

Monitoring the sub-grid scale eddy energy budget provides a means to predict a sub-grid eddy-velocity scale which can be used in the lower order closures.

MEKE equations

The eddy kinetic energy equation is:

\[\partial_{\tilde{t}} E = \overbrace{ \dot{E}_b + \gamma_\eta \dot{E}_\eta + \gamma_v \dot{E}_v }^\text{sources} - \overbrace{ ( \lambda + C_d | U_d | \gamma_b^2 ) E }^\text{local dissipation} + \overbrace{ \nabla \cdot ( ( \kappa_E + \gamma_M \kappa_M ) \nabla E - \kappa_4 \nabla^3 E ) }^\text{smoothing}\]

where \(E\) is the eddy kinetic energy (variable MEKE) with units of m 2s -2, and \(\tilde{t} = a t\) is a scaled time. The non-dimensional factor \(a\geq 1\) is used to accelerate towards equilibrium.

The MEKE equation is two-dimensional and obtained by depth averaging the the three-dimensional eddy energy equation. In the following expressions \(\left< \phi \right> = \frac{1}{H} \int^\eta_{-D} \phi \, dz\) maps three dimensional terms into the two-dimensional quantities needed.

MEKE source terms

The source term \(\dot{E}_b\) is a constant background source of energy intended to avoid the limit \(E\rightarrow 0\).

The “GM” source term

\[\dot{E}_\eta = - \left< \overline{w^\prime b^\prime} \right> = \left< \kappa_h N^2S^2 \right> \approx \left< \kappa_h g\prime |\nabla_\sigma \eta|^2 \right>\]

equals the mean potential energy removed by the Gent-McWilliams closure, and is excluded/included in the MEKE budget by the efficiency parameter \(\gamma_\eta \in [0,1]\).

The “frictional” source term

\[\dot{E}_{v} = \left< \partial_i u_j \tau_{ij} \right>\]

equals the mean kinetic energy removed by lateral viscous fluxes, and is excluded/included in the MEKE budget by the efficiency parameter \(\gamma_v \in [0,1]\).

MEKE dissipation terms

The local dissipation of \(E\) is parameterized through a linear damping, \(\lambda\), and bottom drag, \(C_d | U_d | \gamma_b^2\). The \(\gamma_b\) accounts for the weak projection of the column-mean eddy velocty to the bottom. In other words, the bottom velocity is estimated as \(\gamma_b U_e\). The bottom drag coefficient, \(C_d\) is the same as that used in the bottom friction in the mean model equations.

The bottom drag velocity scale, \(U_d\), has contributions from the resolved state and \(E\):

\[U_d = \sqrt{ U_b^2 + |u|^2_{z=-D} + |\gamma_b U_e|^2 } .\]

where the eddy velocity scale, \(U_e\), is given by:

\[U_e = \sqrt{ 2 E } .\]

\(U_b\) is a constant background bottom velocity scale and is typically not used (i.e. set to zero).

Following Jansen et al., 2015, the projection of eddy energy on to the bottom is given by the ratio of bottom energy to column mean energy:

\[\gamma_b^2 = \frac{E_b}{E} = \gamma_{d0} + \left( 1 + c_{b} \frac{L_d}{L_f} \right)^{-\frac{4}{5}} ,\]
\[\gamma_b^2 \leftarrow \max{\left( \gamma_b^2, \gamma_{min}^2 \right)} .\]

MEKE smoothing terms

\(E\) is laterally diffused by a diffusivity \(\kappa_E + \gamma_M \kappa_M\) where \(\kappa_E\) is a constant diffusivity and the term \(\gamma_M \kappa_M\) is a “self diffusion” using the diffusivity calculated in the section Diffusivity derived from MEKE. \(\kappa_4\) is a constant bi-harmonic diffusivity.

Diffusivity derived from MEKE

The predicted eddy velocity scale, \(U_e\), can be combined with a mixing length scale to form a diffusivity. The primary use of a MEKE derived diffusivity is for use in thickness diffusion (module mom_thickness_diffuse()) and optionally in along isopycnal mixing of tracers (module ) and optionally in along isopycnal mixing of tracers (module mom_tracer_hor_diff()). The original form used (enabled with MEKE_OLD_LSCALE=True):). The original form used (enabled with MEKE_OLD_LSCALE=True):

\[\kappa_M = \gamma_\kappa \sqrt{ \gamma_t^2 U_e^2 A_\Delta }\]

where \(A_\Delta\) is the area of the grid cell. Following Jansen et al., 2015, we now use

\[\kappa_M = \gamma_\kappa l_M \sqrt{ \gamma_t^2 U_e^2 }\]

where \(\gamma_\kappa \in [0,1]\) is a non-dimensional factor and, following Jansen et al., 2015, \(\gamma_t^2\) is the ratio of barotropic eddy energy to column mean eddy energy given by

\[\gamma_t^2 = \frac{E_t}{E} = \left( 1 + c_{t} \frac{L_d}{L_f} \right)^{-\frac{1}{4}} ,\]
\[\gamma_t^2 \leftarrow \max{\left( \gamma_t^2, \gamma_{min}^2 \right)} .\]

The length-scale is a configurable combination of multiple length scales:

\[l_M = \left( \frac{\alpha_d}{L_d} + \frac{\alpha_f}{L_f} + \frac{\alpha_R}{L_R} + \frac{\alpha_e}{L_e} + \frac{\alpha_\Delta}{L_\Delta} + \frac{\delta[L_c]}{L_c} \right)^{-1}\]

where

\[\begin{split}\begin{eqnarray*} L_d & = & \sqrt{\frac{c_g^2}{f^2+2\beta c_g}} \sim \frac{ c_g }{f} \\\\ L_R & = & \sqrt{\frac{U_e}{\beta^*}} \\\\ L_e & = & \frac{U_e}{|S| N} \\\\ L_f & = & \frac{H}{c_d} \\\\ L_\Delta & = & \sqrt{A_\Delta} . \end{eqnarray*}\end{split}\]

\(L_c\) is a constant and \(\delta[L_c]\) is the impulse function so that the term \(\frac{\delta[L_c]}{L_c}\) evaluates to \(\frac{1}{L_c}\) when \(L_c\) is non-zero but is dropped if \(L_c=0\).

\(\beta^*\) is the effective \(\beta\) that combines both the planetary vorticity gradient (i.e. \(\beta=\nabla f\)) and the topographic \(\beta\) effect, with the latter weighed by a weighting constant, \(c_\beta\), that varies from 0 to 1, so that \(c_\beta=0\) means the topographic \(\beta\) effect is ignored, while \(c_\beta=1\) means it is fully considered. The new \(\beta^*\) therefore takes the form of

\[\beta^* = \sqrt{( \partial_xf - c_\beta\frac{f}{D}\partial_xD )^2 + ( \partial_yf - c_\beta\frac{f}{D}\partial_yD )^2}\]

where \(D\) is water column depth at T points.

Viscosity derived from MEKE

As for \(\kappa_M\), the predicted eddy velocity scale can be used to form a harmonic eddy viscosity,

\[\kappa_u = \gamma_u \sqrt{ U_e^2 A_\Delta }\]

as well as a biharmonic eddy viscosity,

\[\kappa_4 = \gamma_4 \sqrt{ U_e^2 A_\Delta^3 }\]

Limit cases for local source-dissipative balance

Note that in steady-state (or when \(a>>1\)) and there is no diffusion of \(E\) then

\[\overline{E} \approx \frac{ \dot{E}_b + \gamma_\eta \dot{E}_\eta + \gamma_v \dot{E}_v }{ \lambda + C_d|U_d|\gamma_b^2 } .\]

In the linear drag limit, where \(U_e << \min(U_b, |u|_{z=-D}, C_d^{-1}\lambda)\), the equilibrium becomes \(\overline{E} \approx \frac{ \dot{E}_b + \gamma_\eta \dot{E}_\eta + \gamma_v \dot{E}_v }{ \lambda + C_d \sqrt{ U_b^2 + |u|^2_{z=-D} } }\).

In the nonlinear drag limit, where \(U_e >> \max(U_b, |u|_{z=-D}, C_d^{-1}\lambda)\), the equilibrium becomes \(\overline{E} \approx \left( \frac{ \dot{E}_b + \gamma_\eta \dot{E}_\eta + \gamma_v \dot{E}_v }{ \sqrt{2} C_d \gamma_b^3 } \right)^\frac{2}{3}\).

MEKE module parameters

Symbol

Module parameter

USE_MEKE

\(a\)

MEKE_DTSCALE

\(\dot{E}_b\)

MEKE_BGSRC

\(\gamma_\eta\)

MEKE_GMCOEFF

\(\gamma_v\)

MEKE_FrCOEFF

\(\lambda\)

MEKE_DAMPING

\(U_b\)

MEKE_USCALE

\(\gamma_{d0}\)

MEKE_CD_SCALE

\(c_{b}\)

MEKE_CB

\(c_{t}\)

MEKE_CT

\(\kappa_E\)

MEKE_KH

\(\kappa_4\)

MEKE_K4

\(\gamma_\kappa\)

MEKE_KHCOEFF

\(\gamma_M\)

MEKE_KHMEKE_FAC

\(\gamma_u\)

MEKE_VISCOSITY_COEFF_KU

\(\gamma_4\)

MEKE_VISCOSITY_COEFF_AU

\(\gamma_{min}^2\)

MEKE_MIN_GAMMA2

\(\alpha_d\)

MEKE_ALPHA_DEFORM

\(\alpha_f\)

MEKE_ALPHA_FRICT

\(\alpha_R\)

MEKE_ALPHA_RHINES

\(\alpha_e\)

MEKE_ALPHA_EADY

\(\alpha_\Delta\)

MEKE_ALPHA_GRID

\(L_c\)

MEKE_FIXED_MIXING_LENGTH

\(c_\beta\)

MEKE_TOPOGRAPHIC_BETA

MEKE_KHTH_FAC

MEKE_KHTR_FAC

Symbol

Model parameter

\(C_d\)

CDRAG

References

Jansen, M. F., A. J. Adcroft, R. Hallberg, and I. M. Held, 2015: Parameterization of eddy fluxes based on a mesoscale energy budget. Ocean Modelling, 92, 2841, http://doi.org/10.1016/j.ocemod.2015.05.007 .

Marshall, D. P., and A. J. Adcroft, 2010: Parameterization of ocean eddies: Potential vorticity mixing, energetics and Arnold first stability theorem. Ocean Modelling, 32, 188204, http://doi.org/10.1016/j.ocemod.2010.02.001 .

Type Documentation

type mom_meke/meke_cs

Control structure that contains MEKE parameters and diagnostics handles.

Type fields
  • % id_meke [integer] :: Diagnostic handles.

  • % id_ue [integer] :: Diagnostic handles.

  • % id_kh [integer] :: Diagnostic handles.

  • % id_src [integer] :: Diagnostic handles.

  • % id_ub [integer] :: Diagnostic handles.

  • % id_ut [integer] :: Diagnostic handles.

  • % id_gm_src [integer] :: Diagnostic handles.

  • % id_mom_src [integer] :: Diagnostic handles.

  • % id_gme_snk [integer] :: Diagnostic handles.

  • % id_decay [integer] :: Diagnostic handles.

  • % id_khmeke_u [integer] :: Diagnostic handles.

  • % id_khmeke_v [integer] :: Diagnostic handles.

  • % id_ku [integer] :: Diagnostic handles.

  • % id_au [integer] :: Diagnostic handles.

  • % id_le [integer] :: Diagnostic handles.

  • % id_gamma_b [integer] :: Diagnostic handles.

  • % id_gamma_t [integer] :: Diagnostic handles.

  • % id_lrhines [integer] :: Diagnostic handles.

  • % id_leady [integer] :: Diagnostic handles.

  • % id_meke_equilibrium [integer] :: Diagnostic handles.

  • % equilibrium_value [real(:,:),pointer] :: The equilbrium value of MEKE to be calculated at each time step [L2 T-2 ~> m2 s-2].

  • % meke_frcoeff [real] :: Efficiency of conversion of ME into MEKE [nondim].

  • % meke_gmcoeff [real] :: Efficiency of conversion of PE into MEKE [nondim].

  • % meke_gmecoeff [real] :: Efficiency of conversion of MEKE into ME by GME [nondim].

  • % meke_damping [real] :: Local depth-independent MEKE dissipation rate [T-1 ~> s-1].

  • % meke_cd_scale [real] :: The ratio of the bottom eddy velocity to the column mean eddy velocity, i.e. sqrt(2*MEKE). This should be less than 1 to account for the surface intensification of MEKE.

  • % meke_cb [real] :: Coefficient in the

  • % meke_min_gamma [real] :: Minimum value of gamma_b^2 allowed [nondim].

  • % meke_ct [real] :: Coefficient in the

  • % visc_drag [logical] :: If true use the vertvisc_type to calculate bottom drag.

  • % meke_geometric [logical] :: If true, uses the GM coefficient formulation from the GEOMETRIC framework (Marshall et al., 2012)

  • % meke_geometric_alpha [real] :: The nondimensional coefficient governing the efficiency of the GEOMETRIC thickness diffusion.

  • % meke_equilibrium_alt [logical] :: If true, use an alternative calculation for the equilibrium value of MEKE.

  • % meke_equilibrium_restoring [logical] :: If true, restore MEKE back to its equilibrium value, which is calculated at each time step.

  • % gm_src_alt [logical] :: If true, use the GM energy conversion form S^2*N^2*kappa rather than the streamfunction for the MEKE GM source term.

  • % rd_as_max_scale [logical] :: If true the length scale can not exceed the first baroclinic deformation radius.

  • % use_old_lscale [logical] :: Use the old formula for mixing length scale.

  • % use_min_lscale [logical] :: Use simple minimum for mixing length scale.

  • % cdrag [real] :: The bottom drag coefficient for MEKE [nondim].

  • % meke_bgsrc [real] :: Background energy source for MEKE [L2 T-3 ~> W kg-1] (= m2 s-3).

  • % meke_dtscale [real] :: Scale factor to accelerate time-stepping [nondim].

  • % meke_khcoeff [real] :: Scaling factor to convert MEKE into Kh [nondim].

  • % meke_uscale [real] :: MEKE velocity scale for bottom drag [L T-1 ~> m s-1].

  • % meke_kh [real] :: Background lateral diffusion of MEKE [L2 T-1 ~> m2 s-1].

  • % meke_k4 [real] :: Background bi-harmonic diffusivity (of MEKE) [L4 T-1 ~> m4 s-1].

  • % khmeke_fac [real] :: A factor relating MEKEKh to the diffusivity used for MEKE itself [nondim].

  • % viscosity_coeff_ku [real] :: The scaling coefficient in the expression for viscosity used to parameterize lateral harmonic momentum mixing by unresolved eddies represented by MEKE.

  • % viscosity_coeff_au [real] :: The scaling coefficient in the expression for viscosity used to parameterize lateral biharmonic momentum mixing by unresolved eddies represented by MEKE.

  • % lfixed [real] :: Fixed mixing length scale [L ~> m].

  • % adeform [real] :: Weighting towards deformation scale of mixing length [nondim].

  • % arhines [real] :: Weighting towards Rhines scale of mixing length [nondim].

  • % africt [real] :: Weighting towards frictional arrest scale of mixing length [nondim].

  • % aeady [real] :: Weighting towards Eady scale of mixing length [nondim].

  • % agrid [real] :: Weighting towards grid scale of mixing length [nondim].

  • % meke_advection_factor [real] :: A scaling in front of the advection of MEKE [nondim].

  • % meke_topographic_beta [real] :: Weight for how much topographic beta is considered when computing beta in Rhines scale [nondim].

  • % meke_restoring_rate [real] :: Inverse of the timescale used to nudge MEKE toward its equilibrium value [s-1].

  • % meke_advection_bug [logical] :: If true, recover a bug in the calculation of the barotropic transport for the advection of MEKE, wherein only the transports in the deepest layer are used.

  • % fixed_total_depth [logical] :: If true, use the nominal bathymetric depth as the estimate of the time-varying ocean depth. Otherwise base the depth on the total ocean mass per unit area.

  • % kh_flux_enabled [logical] :: If true, lateral diffusive MEKE flux is enabled.

  • % initialize [logical] :: If True, invokes a steady state solver to calculate MEKE.

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

  • % diag [type(diag_ctrl),pointer] :: A type that regulates diagnostics output.

  • % id_clock_pass [integer] :: Clock for group pass calls.

  • % pass_meke [type(group_pass_type)] :: Group halo pass handle for MEKEMEKE and maybe MEKEKh_diff.

  • % pass_kh [type(group_pass_type)] :: Group halo pass handle for MEKEKh, MEKEKu, and/or MEKEAu.

Function/Subroutine Documentation

subroutine mom_meke/step_forward_meke(MEKE, h, SN_u, SN_v, visc, dt, G, GV, US, CS, hu, hv)

Integrates forward-in-time the MEKE eddy energy equation. See MEKE equations.

Parameters
  • meke :: MEKE data.

  • g :: [inout] Ocean grid.

  • gv :: [in] Ocean vertical grid structure.

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

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

  • sn_u :: [in] Eady growth rate at u-points [T-1 ~> s-1].

  • sn_v :: [in] Eady growth rate at v-points [T-1 ~> s-1].

  • visc :: [in] The vertical viscosity type.

  • dt :: [in] Model(baroclinic) time-step [T ~> s].

  • cs :: MEKE control structure.

  • hu :: [in] Accumlated zonal mass flux [H L2 ~> m3 or kg]

  • hv :: [in] Accumlated meridional mass flux [H L2 ~> m3 or kg]

Call to

meke_equilibrium meke_equilibrium_restoring meke_lengthscales mom_error_handler::mom_error

subroutine mom_meke/meke_equilibrium(CS, MEKE, G, GV, US, SN_u, SN_v, drag_rate_visc, I_mass, depth_tot)

Calculates the equilibrium solution where the source depends only on MEKE diffusivity and there is no lateral diffusion of MEKE. Results is in MEKEMEKE.

Parameters
  • g :: [inout] Ocean grid.

  • gv :: [in] Ocean vertical grid structure.

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

  • cs :: MEKE control structure.

  • meke :: A structure with MEKE data.

  • sn_u :: [in] Eady growth rate at u-points [T-1 ~> s-1].

  • sn_v :: [in] Eady growth rate at v-points [T-1 ~> s-1].

  • drag_rate_visc :: [in] Mean flow velocity contribution to the MEKE drag rate [L T-1 ~> m s-1]

  • i_mass :: [in] Inverse of column mass [R-1 Z-1 ~> m2 kg-1].

  • depth_tot :: [in] The depth of the water column [Z ~> m].

Call to

meke_lengthscales_0d

Called from

step_forward_meke

subroutine mom_meke/meke_equilibrium_restoring(CS, G, US, SN_u, SN_v, depth_tot)
Parameters
  • g :: [inout] Ocean grid.

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

  • cs :: MEKE control structure.

  • sn_u :: [in] Eady growth rate at u-points [T-1 ~> s-1].

  • sn_v :: [in] Eady growth rate at v-points [T-1 ~> s-1].

  • depth_tot :: [in] The depth of the water column [Z ~> m].

Called from

step_forward_meke

subroutine mom_meke/meke_lengthscales(CS, MEKE, G, GV, US, SN_u, SN_v, EKE, depth_tot, bottomFac2, barotrFac2, LmixScale)

Calculates the eddy mixing length scale and \(\gamma_b\) and \(\gamma_t\) functions that are ratios of either bottom or barotropic eddy energy to the column eddy energy, respectively. See MEKE equations.

Parameters
  • cs :: MEKE control structure.

  • meke :: MEKE data.

  • g :: [inout] Ocean grid.

  • gv :: [in] Ocean vertical grid structure.

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

  • sn_u :: [in] Eady growth rate at u-points [T-1 ~> s-1].

  • sn_v :: [in] Eady growth rate at v-points [T-1 ~> s-1].

  • eke :: [in] Eddy kinetic energy [L2 T-2 ~> m2 s-2].

  • depth_tot :: [in] The depth of the water column [Z ~> m].

  • bottomfac2 :: [out] gamma_b^2 [nondim]

  • barotrfac2 :: [out] gamma_t^2 [nondim]

  • lmixscale :: [out] Eddy mixing length [L ~> m].

Call to

meke_lengthscales_0d

Called from

step_forward_meke

subroutine mom_meke/meke_lengthscales_0d(CS, US, area, beta, depth, Rd_dx, SN, EKE, bottomFac2, barotrFac2, LmixScale, Lrhines, Leady)

Calculates the eddy mixing length scale and \(\gamma_b\) and \(\gamma_t\) functions that are ratios of either bottom or barotropic eddy energy to the column eddy energy, respectively. See MEKE equations.

Parameters
  • cs :: MEKE control structure.

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

  • area :: [in] Grid cell area [L2 ~> m2]

  • beta :: [in] Planetary beta = \(\nabla f\) [T-1 L-1 ~> s-1 m-1]

  • depth :: [in] Ocean depth [Z ~> m]

  • rd_dx :: [in] Resolution Ld/dx [nondim].

  • sn :: [in] Eady growth rate [T-1 ~> s-1].

  • eke :: [in] Eddy kinetic energy [L2 T-2 ~> m2 s-2].

  • bottomfac2 :: [out] gamma_b^2

  • barotrfac2 :: [out] gamma_t^2

  • lmixscale :: [out] Eddy mixing length [L ~> m].

  • lrhines :: [out] Rhines length scale [L ~> m].

  • leady :: [out] Eady length scale [L ~> m].

Called from

meke_equilibrium meke_lengthscales

function mom_meke/meke_init(Time, G, US, param_file, diag, CS, MEKE, restart_CS) [logical]

Initializes the MOM_MEKE module and reads parameters. Returns True if module is to be used, otherwise returns False.

Parameters
  • time :: [in] The current model time.

  • g :: [inout] The ocean’s grid structure.

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

  • param_file :: [in] Parameter file parser structure.

  • diag :: [inout] Diagnostics structure.

  • cs :: MEKE control structure.

  • meke :: MEKE-related fields.

  • restart_cs :: Restart control structure for MOM_MEKE.

Call to

mom_error_handler::mom_error mom_error_handler::mom_mesg

subroutine mom_meke/meke_alloc_register_restart(HI, param_file, MEKE, restart_CS)

Allocates memory and register restart fields for the MOM_MEKE module.

Parameters
  • hi :: [in] Horizontal index structure

  • param_file :: [in] Parameter file parser structure.

  • meke :: A structure with MEKE-related fields.

  • restart_cs :: Restart control structure for MOM_MEKE.

Call to

mom_error_handler::mom_error mom_error_handler::mom_mesg mom_io::var_desc

subroutine mom_meke/meke_end(MEKE)

Deallocates any variables allocated in MEKE_alloc_register_restart.

Parameters

meke :: [inout] A structure with MEKE-related fields.