Vertical Mixing Parametrizations

These equations are as given in the ROMS documentation

The type of vertical mixing parametrization can be chosen in the inputs file with the remora.vertical_mixing_type option.

Analytic Mixing

By default, the vertical mixing coefficients \(K_v\) and \(K_t\) (vec_Akv and vec_Akt in the code) are defined as a custom function in init_analytic_vmix in the problem’s prob.cpp file. This function can be defined to be a function of simulation variables, most commonly depth \(z_w\). It is re-evaluated at every time step.

Generic Length Scale (GLS)

The Generic Length Scale (GLS) model is a two-equation turbulence closure scheme which can be tuned to behave like many of the traditional schemes such as Mellor and Yamada 2.5. This class of schemes add additional prognostic equations for the turbulent kinetic energy \(k=\frac{q^2}{2}\) and the product \(kl\), where \(l\) is a length scale. The evolution equation for TKE is:

\[\frac{D}{Dt}\left(\frac{q^2}{2}\right)-\frac{\partial}{\partial z} \left[K_q \frac{\partial}{\partial z} \left(\frac{q^2}{2}\right)\right] = P_s + P_b + \epsilon\]

where \(P_s\) is the shear production, \(P_b\) is the buoyant production, and \(\epsilon\) is the dissipation. In model coordinates, this evolution equation becomes:

\[\begin{split}\frac{\partial}{\partial t} \left(\frac{H_z q^2}{mn}\right) &+ \frac{\partial}{\partial \xi} \left(\frac{H_z u q^2}{n}\right) + \frac{\partial}{\partial \eta}\left(\frac{H_z v q^2}{m}\right) + \frac{\partial}{\partial s}\left(\frac{H_z \Omega q^2}{mn}\right) - \frac{\partial}{\partial s}\left(\frac{K_q}{mnH_z}\frac{\partial q^2}{\partial s}\right)\\ &= \frac{2 H_z}{mn}\left(P_s + P_b + \epsilon\right).\end{split}\]

The terms on the right-hand side are:

\[\begin{split}P_s &= K_v \left[\left(\frac{\partial u}{\partial z}\right)^2 + \left(\frac{\partial v}{\partial z}\right)^2\right]\\ P_b &= -K_t N^2\\ \epsilon &= \left(c_mu^0\right)^{3+p/n} k^{3/2+m/n} \psi^{-1/n}.\end{split}\]

The parameter \(\psi\) is used to extablished the turbulence length scale. Its evolution equation is:

\[\frac{D\psi}{Dt} = \frac{\partial}{\partial z} \left(K_{\psi} \frac{\partial\psi}{\partial z}\right) + \frac{\psi}{k} \left(c_1 P_s + c_3 P_b - c_2 \epsilon F_{\mathrm{wall}}\right).\]

The coefficients \(c_1\) and \(c_2\) are chosen to be consistent with observations of decaying homogeneous, isotropic turbulence. The other parameter \(c_3\) has different values of stable (\(c_3^+\)) and unstable (\(c_3^-\)) stratification. Also,

\[\begin{split}\psi &= \left(c_{\mu}^0\right)^p k^m l^n\\ l &= \left(c_{\mu}^0\right)^3 k^{3/2}\epsilon - 1\end{split}\]

The indices \(p\), \(m\), and \(n\) as well as coefficients are set in the inputs file. The default values correspond to the \(k-\epsilon\) turbulence model.

The equations for \(q\) and \(\psi\) are evolved much like the model tracer equations, including an implicit solve for vertical operations. We use a predictor-corrector scheme in which the predictor step only computes advection. The mixing coefficients are calculated from \(q\) and \(l\):

\[\begin{split}K_v &= q l S_m + K_{v,\mathrm{background}}\\ K_t &= q l S_h + K_{t,\mathrm{background}}\\ K_k &= q l S_m / \sigma_k + K_{k,\mathrm{background}}\\ K_{\psi} &= q l S_m / \sigma_{\psi} + K_{\psi,\mathrm{background}}\end{split}\]

The constants \(\sigma_k = K_v / K_k\) and \(\sigma_{\psi} = K_v / K_{\psi}\) are also set in the inputs file. The stability coefficients, \(S_m\) and \(S_h\) are calculated using either the Galperin, Canuto-A or Canuto-B schemes.