Equations of Motion

The primitive equations in Cartesian coordinates are shown here. The momentum balance in the $\,\!x$- and $\,\!y$-directions are:

${\frac {\partial u}{\partial t}}+{\vec {v}}\cdot \nabla u-fv=-{\frac {\partial \phi }{\partial x}}+{F}_{u}+{D}_{u}$

${\frac {\partial v}{\partial t}}+{\vec {v}}\cdot \nabla v+fu=-{\frac {\partial \phi }{\partial y}}+{F}_{v}+{D}_{v}$

The time evolution of the potential temperature and salinity fields, $\,\!T(x,y,z,t)$ and $\,\!S(x,y,z,t)$, are governed by the advective-diffusive equations:

${\frac {\partial T}{\partial t}}+{\vec {v}}\cdot \nabla T={F}_{T}+{D}_{T}$

${\frac {\partial S}{\partial t}}+{\vec {v}}\cdot \nabla S={F}_{S}+{D}_{S}$

The equation of state is given by:

$\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 _{o}}}$

The final equation expresses the continuity equation for an incompressible fluid:

${\frac {\partial u}{\partial x}}+{\frac {\partial v}{\partial y}}+{\frac {\partial w}{\partial z}}=0$

For the moment, the effects of forcing and dissipation will be represented by the schematic terms $F\,\!$ and $D\,\!$, respectively. The horizontal and vertical mixing will be described more fully in Horizontal Mixing and Vertical Mixing Parameterizations. The variables used are shown here:

${D}_{u},{D}_{v},{D}_{T},{D}_{S}\,\!$	diffusive terms
${F}_{u},{F}_{v},{F}_{T},{F}_{S}\,\!$	forcing terms
$f(x,y)\,\!$	Coriolis parameter
$g\,\!$	acceleration of gravity
$h(x,y)\,\!$	bottom depth
$\nu ,\kappa \,\!$	horizontal viscosity and diffusivity
$K_{m},K_{T},K_{S}\,\!$	vertical viscosity and diffusivity
$P\,\!$	total pressure $P\approx -\rho _{o}gz$
$\phi (x,y,z,t)\,\!$	dynamic pressure $\phi =\left(P/\rho _{o}\right)$
$\rho _{o}+\rho (x,y,z,t)\,\!$	total in situ density
$S(x,y,z,t)\,\!$	salinity
$t\,\!$	time
$T(x,y,z,t)\,\!$	potential temperature
$u,v,w\,\!$	the ($x,y,z\,\!$) components of vector velocity ${\vec {v}}\,\!$
$x,y\,\!$	horizontal coordinates
$z\,\!$	vertical coordinate
$\zeta (x,y,t)\,\!$	the surface elevation

Vertical Boundary Conditions

The vertical boundary conditions can be prescribed as follows:

top ($z=\zeta (x,y,t))\,\!$:

$K_{m}\,{\frac {\partial u}{\partial z}}=\tau _{s}^{x}(x,y,t)$

$K_{m}\,{\frac {\partial v}{\partial z}}=\tau _{s}^{y}(x,y,t)$

$K_{T}\,{\frac {\partial T}{\partial z}}={Q_{T} \over \rho _{o}c_{P}}+{1 \over \rho _{o}c_{P}}{dQ_{T} \over dT}(T-T_{\rm {ref}})$

$K_{S}\,{\frac {\partial S}{\partial z}}={(E-P)S \over \rho _{o}}$

$w={\partial \zeta \over \partial t}$

and bottom ($z=-h(x,y)\,\!$):

$K_{m}\,{\frac {\partial u}{\partial z}}=\tau _{b}^{x}(x,y,t)$

$K_{m}\,{\frac {\partial v}{\partial z}}=\tau _{b}^{y}(x,y,t)$

$K_{T}\,{\frac {\partial T}{\partial z}}=0$

$K_{S}\,{\frac {\partial S}{\partial z}}=0$

$-w+{\vec {v}}\cdot \nabla h=0$

The new variables above are:

Variable	Description
$E-P\,\!$	evaporation minus precipitation
$\gamma _{1},\gamma _{2}\,\!$	linear and quadratic bottom stress coefficients
$Q_{T}\,\!$	surface heat flux
$\tau _{s}^{x},\tau _{s}^{y}\,\!$	surface wind stress
$\tau _{b}^{x},\tau _{b}^{y}\,\!$	bottom stress
$T_{\rm {ref}}\,\!$	surface reference temperature

The surface boundary condition variables are defined in the table above. Since $Q_{T}\,\!$ is a strong function of the surface temperature, it is also prudent to include a correction term for the change in $Q\,\!$ as the surface temperature drifts away from the reference temperature that was used in computing $Q_{T}\,\!$. On the variable bottom, $z=-h(x,y)\,\!$, the horizontal velocity components are constrained to accommodate a prescribed bottom stress which is a sum of linear and quadratic terms (actually one or the other, or a log layer, depending on the cpp options):

$\tau _{b}^{x}=(\gamma _{1}+\gamma _{2}{\sqrt {u^{2}+v^{2}}})u$

$\tau _{b}^{y}=(\gamma _{1}+\gamma _{2}{\sqrt {u^{2}+v^{2}}})v$

The vertical heat and salt flux may also be prescribed at the bottom, although they are usually set to zero.

Horizontal Boundary Conditions

As distributed, the model can easily be configured for a periodic channel, a doubly periodic domain, or a closed basin. Code is also included for open boundaries which may or may not work for your particular application. Appropriate boundary conditions are provided for $u,v,T,S,\,\!$ and $\zeta \,\!$. At every timestep the subroutines xxxx and yyyy are called to fill in the necessary boundary values.

The model domain is logically rectangular, but it is possible to mask out land areas on the boundary and in the interior. Boundary conditions on these masked regions are straightforward, with a choice of no-slip or free-slip walls.

If biharmonic friction is used, a higher order boundary condition must also be provided. The model currently has this built into the code where the biharmonic terms are calculated. The high order boundary conditions used for $u\,\!$ are ${\frac {\partial }{\partial x}}\left({\frac {h\nu }{mn}}{\frac {\partial ^{2}u}{\partial x^{2}}}\right)=0$ on the eastern and western boundaries and ${\frac {\partial }{\partial y}}\left({\frac {h\nu }{mn}}{\frac {\partial ^{2}u}{\partial y^{2}}}\right)=0\,$ on the northern and southern boundaries. The boundary conditions for $v,T,\,\!$ and $S\,\!$ are similar. These boundary conditions were chosen because they preserve the property of no gain or loss of volume-integrated momentum, temperature, or salt.