A Conservative Model for Nonlinear Dynamics in a Stratified, Rotating Fluid

We present a set of equations describing the nonlinear dynamics of flows constrained by environmental rotation and stratification (Rossby numbers [ ] 0.1,0.5 Ro ∈ and Burger numbers of order unity). The fluid is assumed incompressible, adiabatic, inviscid and in hydrostatic balance. This set of equations is derived from the Navier Stokes equations (with the above properties), using a Rossby number expansion with second order truncation. The resulting model has the following properties: 1) it can represent motions with moderate Rossby numbers and a Burger number of order unity; 2) it filters inertia-gravity waves by assuming that the divergence of horizontal velocity remains small; 3) it is written in terms of a single function of space and time (pressure, generalized streamfunction or Bernoulli function); 4) it conserves total (Ertel) vorticity in a Lagrangian form, and its quadratic norm (potential enstrophy) at the model order in Rossby number; 5) it also conserves total energy at the same order if the work of pressure forces vanishes when integrated over the fluid domain. The layerwise version of the model is finally presented, written in terms of pressure. Integral properties (energy, enstrophy) are conserved by these layerwise equations. The model equations agree with the generalized geostrophy equations in the appropriate parameter regime. Application to vortex dynamics are mentioned.

). Modeling these flows has therefore been the subject of many studies, using the "primitive equations", i.e. the three-dimensional Navier-Stokes equations with hydrostatic balance, Boussinesq approximation and incompressibility. These primitive equations, describe both fast (wave-like) and slow (vortex-like) motions. These two subspaces of solutions of the primitive equations are independent only for linear dynamics. In nonlinear dynamics, fast and slow motions can interact : initially slow motions, free of fast waves, remain in the slow subspace, but nonlinear interaction of fast waves can contribute to slow motions. Therefore the slow motions subspace is called a quasi manifold.
The importance of quasi-2D jets and vortices in geophysical flows and laboratory experiments has motivated the search for asymptotic limits of the primitive equations on the slow quasimanifold (see for instance [4]). This regime is governed by the time evolution of potential vorticity (also called Ertel vorticity, a combination of vertical vorticity with vertical density gradient). For incompressible and inviscid fluids, potential vorticity is conserved by each fluid element (see [5]). To filter out fast (inertia-gravity) waves, zero or weak divergence of horizontal velocity is usually imposed. The most prominent example of slow regime is the horizontal balance between the Coriolis acceleration and the pressure gradient (for parallel flows). This equilibrium is known as geostrophic balance, for which horizontal velocity has zero divergence. The following model in complexity, with small horizontal velocity divergence, is the quasi-geostrophic model [6]. This model is derived from the primitive equations, by assuming dominant Coriolis acceleration compared with relative directly from the primitive equations for continuously stratified flows, contrary to the generalized geostrophic model, derived from the shallow-water equations (layerwise primitive equations). Then, the projection of these equations onto a set of homogeneous layers is shown to be in agreement with equations directly derived from the shallow-water equations. Such a model is important for the study of mesoscale vortex dynamics in the ocean; for these vortices, Ro is not very small but often remains smaller than 1. When modeling the stability or the interaction of such vortices, it is essential that the model preserves their potential vorticity, and their energy when both are physically conserved ( [15] [16] [17]).
Here, we present the primitive equations and the physical scaling (Section 2). In Section 3, the intermediate model is derived and its conservation properties at ( ) 3 O Ro are presented: Lagrangian invariance of Ertel vorticity, energy conservation in any volume for which the integral flux of the Bernoulli function across the boundary vanishes (as for primitive equations). The model equations are written in terms of a generalized streamfunction or of Bernoulli function. In Section 4, the layerwise model is presented and compared with previously derived models (frontal geostrophic equations, generalized geostrophic equations). Finally, conclusions are drawn.

Dimensional Primitive Equations
The fundamental equations for fluid motions strongly constrained by rotation and stratification are the primitive equations, derived from the three-dimensional Navier-Stokes equations with background rotation, under the following assumptions: • incompressibility (moderate particle velocities with respect to the speed of sound); • Boussinesq approximation (invariance of density in the horizontal momentum equations); • negligible vertical accelerations leading to hydrostatic balance; • moderate vertical velocities, implying that the Coriolis force mostly acts on the horizontal velocity; • background rotation a linear function of the y coordinate 0 f f y β = + (beta-plane approximation). On Earth, this corresponds to a projection of the motion on the tangential plane to Earth at a given latitude: 1 : where * * * , , u v w are the velocity components (functions of , , , x y z t ), g is the local gravity intensity; * T p and * T ρ are the total pressure and density; pressure is split into static and dynamic parts and density into mean and anomaly:

Scale Analysis
The horizontal velocity * * , u v is now scaled by U, the vertical velocity * Bu From these numbers and dynamical relations, we scale: • pressure via geostrophy, To illustrate this scaling, consider vortices in the ocean (at mid latitudes), for

~1
Ro , the primitive equations should be used. Moreover, as β is much smaller than Ro , the beta-plane approximation is justified.

Dimensionless Primitive Equations
With this scaling and dropping stars for dimensional variables, equation (QM x ) is: and Equation H is ( ), D wS z Dt ρ = (11) where the total (three-dimensional) derivative is written in scaled form: for which no truncation of , , u v w in Ro has yet been made. Open Journal of Marine Science Following [5], it is straightforward to show that these equations conserve total (potential) vorticity where z n is the unit vertical vector.

Generalized Vorticity Equation
We first rewrite equation H to extract w as an exact derivative: This defines X, a quantity physically related to the static and dynamic densities. This expression of the vertical velocity w is central to the subsequent algebra. We note that this expression is similar in form to the quasi-geostrophic expression of w, for which the Lagrangian derivative is simply two-dimensional: We also introduce the absolute vorticity Y, the sum of the relative vorticity and of the Coriolis parameter:
x y RoY y Ro v u The generalized vorticity equation (GV), that we derive hereafter, can be obtained by successive steps: • first by taking the curl of (QM x , QM y ), • then by substituting , u v via (QM x , QM y ), z p ∂ via (QM z ) and w via (H), • and by canceling the three-dimensional velocity divergence via equation (C) (without any approximation). This rather tedious algebra can be more elegantly replaced by the following symbolic combination leading to (GV): is a symbolic writing of Equations (7)-(9).
Physically, in this expression: • the first term on the first line is the absolute vorticity; • the rest of the first line represents the coupling of the baroclinicity vector with the vortex stretching; N. Filatoff, X. Carton Open Journal of Marine Science • the first term on the second line is the coupling of the velocity divergence with the absolute vorticity; • the second part of the second line is the coupling of that velocity divergence with baroclinicity; • finally, the third line contains the vortex stretching, its coupling with the absolute vorticity and with the baroclinicity vector.
At first order in Ro, the total vorticity is that of the quasi-geostrophic (QG) model, and at the following order in Ro all terms are nonlinear and represent couplings between the QG terms.
Following the formal calculation, the (GV) equation is obtained: We remark that Q has the form A B C * + , where A is the absolute vorticity, B is the "baroclinic" part and C is related to the baroclinicity vector. Neglecting the horizontal variations of vertical velocity, Q can be expressed as the product of the horizontal vorticity by a scalar gradient.
For very large Bu, the barotropic vorticity is retrieved and we note From here on, equation (C) is replaced by equation (GV). For consistency, remark that the 5 unknowns are still described by 5 equations. The following step relates the generalized vorticity Q and the Ertel vorticity.

Conservation of the Ertel Vorticity
In dimensionless variables, the Ertel vorticity is written Note that: We remark that A is a function of density only. The conservation of the generalized vorticity, Equation (15), implies that ( ) some straightforward algebra, with the use of relation (13), leads to: which is mass conservation, and thus immediately to: with the form of Π given by Equation (17). The Ertel vorticity is thus conserved as for the primitive equations, but here only up to ( ) 3 O Ro . We remark that, contrary to quasi-geostrophic motions where total vorticity is advected by a two-dimensional geostrophic velocity, here the Ertel vorticity is advected by the three-dimensional velocity.
With less complex algebra, the equivalent relation between Q and Π for quasi-geostrophic dynamics is: Both in the QG and GV models, Π contains both the static and dynamic parts of potential vorticity, while Q expresses the dynamic anomaly, relevant to the description of motions (Y. Morel, priv. comm.). A represents the new distribution of fluid masses induced by these motions.

Total Energy
Kinetic energy density is expressed in full form as: At the model order, the kinetic energy density is multiplied by Ro and the following form is kept: where the horizontal velocity expansion in Ro is kept up to second order for consistency.
Potential energy density is expressed as (see Appendix A): The energy conservation equation (E) is obtained via: The main steps of the calculation are: which expresses the conventional result that the variation of energy with time is due to the work of pressure forces, or in Eulerian form  The total energy is then exactly conserved.

Model Equations in Generalized Streamfunction
Here the three-dimensional Lagrangian derivative is expanded in terms of  (8)- (13). The product L of these operators will permute with the time derivative and the Jacobian. This specific property will help transform Equation

Vertical Derivative
The expression of the vertical velocity is given by equation (H), that is, the quantity is substracted from the total derivative: where we define the vertical operator z L as . Using Equation (9) and (13) as a second intermediate step.

Explicit Form of the Generalized Vorticity in Generalized Streamfunction
With the 3D Lagrangian derivative written in terms of ψ (or p) hereabove, the generalized vorticity Q can also be expressed in terms of ψ . This expression is obtained in two steps: first, the generalized vorticity ( ) , Q ζ ρ (Equation (16)) is given as a function of pressure p, using the expressions ( ) ( ) The quasi-geostrophic form of the total vorticity Q is: As all terms containing ψ in Q will be factor of Ro , the expression of ( ) p ψ given in Equation (27) We recognize in this expression the quasi-geostrophic terms at In generalized streamfunction, the conservation of generalized vorticity is thus: which, with expressions (27) and (30), is the final result for Section 3.4.
Obviously energy conservation could be derived directly from Equation (33)

An Alternative Form of the Dynamics in Bernoulli Function
We recall here that all tilded quantities are massless (i.e. are divided by s ρ ), as it was already defined in Equation (27) The two forms in h B  and Φ are equivalent (see Appendix B4).
As for the horizontal Lagrangian operator in Section 3.4.2, the vertical operator z L must now be applied to obtain the total Lagrangian derivative: where ( ) h Ro Bu L X has finally been replaced by ( ) Thus the Bernoulli formulation of the evolution equation for the generalized vorticity is:

Layer Model
The form of the relative vorticity and of the Lagrangian derivative in terms of pressure given in part 1 is projected onto homogeneous layers. The basic property of the layer model is the independence of the horizontal velocities (and of the horizontal vorticity ζ ) on the vertical coordinate in each homogeneous layer.
Obviously, energy and enstrophy are conserved to the same order in Ro as for the continuously stratified case.

Layerwise Conservation of Ertel Vorticity
To obtain the conservation of the Ertel vorticity in each layer, successive steps are followed: Rather than using potential vorticity, its inverse, potential thickness Indeed, in case of strong interface deformation, i h could locally vanish or become considerably small, whereas the relative vorticity remains positive since 0.5 Ro ≤ (and β is usually much smaller). Therefore, i Z cannot become infinite (inertial instability cannot arise).
In each layer, potential thickness conservation is written: In terms of pressure, potential thickness conservation in each layer is written: and the operator i L is, in terms of pressure: The term i i L Z has the following form in layers for the coupling with the bottom topography.
Numerical implementation of these equations in the N-layer form, with iterative solution of pressure from potential thickness is described and applied in (Filatoff et al., 1997).

Comparison with the Generalized Geostrophy Model
We show here that a simple restriction of the two-layer model yields the well-known generalized geostrophy equations ( [11]). We recall that these equations were derived in the context of a rigid lid and flat bottom ocean. For an easier comparison, we use similar notations to these authors' (beware that some might differ from ours in the previous sections). We define the interface elevation as: is not used.
The equations derived for continuous stratification have then been projected onto a discrete representation of vertical stratification. In the layer representation, the dynamical model equations conserve potential vorticity for each fluid particle and integral invariants for the whole fluid. For practical purposes, potential thickness was used rather than potential vorticity; this avoids singularity when layers become vanishingly thin. The model equations were compared successfully to the frontal and to the generalized geostrophic equations derived by [11]; this comparison also identified the terms which appear in the equations at higher order than that considered by these authors.
This model has been implemented in generalized streamfunction, both with pseudospectral and finite-difference horizontal schemes. It has been applied to various studies of vortex dynamics: • barotropic vortex merger, showing a cyclone/anticyclone asymmetry when Ro increases [18]; this asymmetry is similar in the SWE and GV models but does not exist in the QG model. • intrathermocline vortex merger (the merger of two vortices in the middle layer of a three-layer ocean); the merger of two distant anticyclones is possible both in the SWE and the GV models, but not in the QG model.

B4. Equivalence between the Expressions of Generalized Vorticity in Bernoulli Function and in Streamfunction
To prove the identity between the two expressions of the generalized vorticity