1.1. Primitive Equations
At resolutions larger than 1 km (more marginally above 100 m), The ocean is a fluid that can be described to a good approximation by the primitive equations. The PE equations are simplifications from the Navier-Stokes equations made from scale considerations, along with a nonlinear equation of state, which couples the two active tracers (temperature and salinity):
Hydrostatic hypothesis: the vertical momentum equation is reduced to a balance between the vertical pressure gradient and the buoyancy force (nonhydrostatic processes such as convection must be parametrized)
Boussinesq hypothesis: density variations are neglected except in their contribution to the buoyancy force
Incompressibility hypothesis (steming from the former): the three dimensional divergence of the velocity vector is assumed to be zero.
Spherical earth approximation: the geopotential surfaces are assumed to be spheres so that gravity (local vertical) is parallel to the earth’s radius
Thin-shell approximation: the ocean depth is neglected compared to the earth’s radius
Turbulent closure hypothesis: the turbulent fluxes (which represent the effect of small scale processes on the large-scale) are expressed in terms of large-scale features
By default (#undef NBQ), CROCO solves the primitive equations. But it has also the ability to relax the first 3 hypothesis (#define NBQ). When SOLV3D is not activated, CROCO can be used as a classical shallow water model.
1.1.1. Equations in Cartesian coordinates
\[ \begin{align}\begin{aligned}\frac{\partial u}{\partial t}
+ \vec{\bf \nabla} . \left ( \vec{\textbf v} u \right ) - f v =
- \frac{\partial \phi}{\partial x} + \mathcal{F}_u + \mathcal{D}_u\\\frac{\partial v}{\partial t}
+ \vec{\bf \nabla} . \left ( \vec{\textbf v} v \right ) + f u =
- \frac{\partial \phi}{\partial y} + \mathcal{F}_v + \mathcal{D}_v\end{aligned}\end{align} \]
Turbulent closure schemes are applied to parametrized subgrid-scale vertical fluxes.
The time evolution of a scalar concentration field, \(C(x,y,z,t)\) (e.g. salinity, temperature, or nutrients), is governed by the advective-diffusive equation :
\[\frac{\partial C}{\partial t}
+ \vec{\bf \nabla} . \left ( \vec{\textbf v} C \right ) =
\mathcal{F}_C + \mathcal{D}_C\]
\[\rho = \rho(T,S,P)\]
In the Boussinesq approximation, density variations are neglected in the momentum equations except in their contribution to the buoyancy force in the vertical momentum equation. Under the hydrostatic approximation, it is further assumed that the vertical pressure gradient balances the buoyancy force :
\[\frac{\partial \phi}{\partial z} = - \frac{\rho g}{\rho_0}\]
\[\vec{\bf \nabla} . \vec{\textbf v} =
\frac{\partial u}{\partial x} + \frac{\partial v}{\partial y}
+ \frac{\partial w}{\partial z} = 0\]
The variables used are :
\(\mathcal{D}_u, \mathcal{D}_v, \mathcal{D}_C\) : diffusive terms
\(\mathcal{F}_u, \mathcal{F}_v, \mathcal{F}_C\) : forcing terms
\(f(x,y)\) : Traditional Coriolis parameter \(2 \Omega sin \phi\)
\(g\) : acceleration of gravity
\(\phi(x,y,z,t)\) : dynamic pressure \(\phi=P/\rho_0\), with P the total pressure
\(\rho_0+\rho(x,y,z,t)\) : total in situ density
\(u,v,w\) : the (x,y,z) components of vector velocity \(\vec{\textbf v}\)
1.1.2. Equations in terrain following coordinates
We first introduce a generalized stretched vertical coordinate system (\(s\)), which sets the variable bottom flat at \(z=-h(x,y)\). \(s\) spans the range from -1 (bottom) to 0 (surface) and the transformation rules are:
\[\begin{split} \begin{aligned}
{\left(\frac{\partial}{\partial x}\right)_{z}=\left(\frac{\partial}{\partial x}\right)_{s}-\left(\frac{1}{H_{z}}\right)\left(\frac{\partial z}{\partial x}\right)_{s} \frac{\partial}{\partial s}} \\
{\left(\frac{\partial}{\partial y}\right)_{z}=\left(\frac{\partial}{\partial y}\right)_{s}-\left(\frac{1}{H_{z}}\right)\left(\frac{\partial z}{\partial y}\right)_{s} \frac{\partial}{\partial s}} \\
{\frac{\partial}{\partial z}=\left(\frac{\partial s}{\partial z}\right) \frac{\partial}{\partial s}=\frac{1}{H_{z}} \frac{\partial}{\partial s}} \\
\text{where } {H_{z} \equiv \frac{\partial z}{\partial s}}
\end{aligned}\end{split}\]
The vertical velocity in \(s\) coordinate is:
\[\begin{split}\begin{aligned} \Omega(x, y, s, t) &=\frac{1}{H_{z}}\left[w-(1+s) \frac{\partial \zeta}{\partial t}-u \frac{\partial z}{\partial x}-v \frac{\partial z}{\partial y}\right] \\
w &=\frac{\partial z}{\partial t}+u \frac{\partial z}{\partial x}+v \frac{\partial z}{\partial y}+\Omega H_{z} \end{aligned}\end{split}\]
\(\Omega=0\) at both surface and bottom.
Next, the requirement for a laterally variable grid resolution can also be met, for suitably smooth domains, by introducing an appropriate orthogonal coordinate transformation in the horizontal. Let the new coordinates be \(\xi(x,y)\) and \(\eta(x,y)\) where the relationship of horizontal arc length to the differential distance is given by:
\[ \begin{align}\begin{aligned}\qquad(d s)_{\xi}=\left(\frac{1}{m}\right) d \xi\\\qquad(d s)_{\eta}=\left(\frac{1}{n}\right) d \eta\end{aligned}\end{align} \]
Here \(m(\xi,\eta)\) and \(n(\xi,\eta)\) are the scale factors which relate the differential distances \((\Delta \xi, \Delta \eta)\) to the actual (physical) arc lengths.
\[\begin{split} \begin{aligned} \frac{\partial}{\partial t}\left(\frac{H_{z} u}{m n}\right)+\frac{\partial}{\partial \xi}\left(\frac{H_{z} u^{2}}{n}\right)+\frac{\partial}{\partial \eta}\left(\frac{H_{z} u v}{m}\right)+\frac{\partial}{\partial s}\left(\frac{H_{z} u \Omega}{m n}\right) \\-\left\{\left(\frac{f}{m n}\right)+v \frac{\partial}{\partial \xi}\left(\frac{1}{n}\right)-u \frac{\partial}{\partial \eta}\left(\frac{1}{m}\right)\right\} H_{z} v=\\-\left(\frac{H_{z}}{n}\right)\left(\frac{\partial \phi}{\partial \xi}+\frac{g \rho}{\rho_{o}} \frac{\partial z}{\partial \xi}+g \frac{\partial \zeta}{\partial \xi}\right)+\frac{H_{z}}{m n}\left(\mathcal{F}_{u}+\mathcal{D}_{u}\right) \\ \frac{\partial}{\partial t}\left(\frac{H_{z} v}{m n}\right)+\frac{\partial}{\partial \xi}\left(\frac{H_{z} u v}{n}\right)+\frac{\partial}{\partial \eta}\left(\frac{H_{z} v^{2}}{m}\right)+\frac{\partial}{\partial s}\left(\frac{H_{z} v \Omega}{m n}\right) \\+\left\{\left(\frac{f}{m n}\right)+v \frac{\partial}{\partial \xi}\left(\frac{1}{n}\right)-u \frac{\partial}{\partial \eta}\left(\frac{1}{m}\right)\right\} H_{z} u=\\-\left(\frac{H_{z}}{m}\right)\left(\frac{\partial \phi}{\partial \eta}+\frac{g \rho}{\rho_{o}} \frac{\partial z}{\partial \eta}+g \frac{\partial \zeta}{\partial \eta}\right)+\frac{H_{z}}{m n}\left(\mathcal{F}_{v}+\mathcal{D}_{v}\right) \\
\frac{\partial}{\partial t}\left(\frac{H_{z} T}{m n}\right)+\frac{\partial}{\partial \xi}\left(\frac{H_{z} u T}{n}\right)+\frac{\partial}{\partial \eta}\left(\frac{H_{z} v T}{m}\right)+\frac{\partial}{\partial s}\left(\frac{H_{z} \Omega T}{m n}\right) &=\frac{H_{z}}{m n}\left(\mathcal{F}_{T}+\mathcal{D}_{T}\right) \\ \frac{\partial}{\partial t}\left(\frac{H_{z} S}{m n}\right)+\frac{\partial}{\partial \xi}\left(\frac{H_{z} u S}{n}\right)+\frac{\partial}{\partial \eta}\left(\frac{H_{z} v S}{m}\right)+\frac{\partial}{\partial s}\left(\frac{H_{z} \Omega S}{m n}\right) &=\frac{H_{z}}{m n}\left(\mathcal{F}_{S}+\mathcal{D}_{S}\right) \\ \rho &=\rho(T, S, P) \\ \frac{\partial \phi}{\partial s} &=-\left(\frac{g H_{z} \rho}{\rho_{o}}\right) \\ \frac{\partial}{\partial t}\left(\frac{H_{z}}{m n}\right)+\frac{\partial}{\partial \xi}\left(\frac{H_{z} u}{n}\right)+\frac{\partial}{\partial \eta}\left(\frac{H_{z} v}{m}\right)+\frac{\partial}{\partial s}\left(\frac{H_{z} \Omega}{m n}\right) &=0 \end{aligned}\end{split}\]