# Energy transfers in internal tide generation, propagation and dissipation in the deep ocean

Jochem Willem Floor, Francis Auclair, Patrick Marsaleix

► To cite this version:

Jochem Willem Floor, Francis Auclair, Patrick Marsaleix. Energy transfers in internal tide generation, propagation and dissipation in the deep ocean. Ocean Modelling, 2011, 38 (1-2), pp.22-40. 10.1016/j.ocemod.2011.01.009. hal-01000133

HAL Id: hal-01000133

https://hal.science/hal-01000133

Submitted on 5 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Energy transfers in internal tide generation, propagation and dissipation in the deep ocean

J.W. Floor*, F. Auclair, P. Marsaleix

Laboratoire d'Aérologie, University of Toulouse/CNRS UMR 5560, 14 avenue Edouard Belin, F-31400 Toulouse, France

The energy transfers associated with internal tide (IT) generation by a semi-diurnal surface tidal wave impinging on a supercritical meridionally uniform deep ocean ridge on the f-plane, and subsequent IT-propagation are analysed using the Boussinesq, free-surface, terrain-following ocean model Symphonie. The energy diagnostics are explicitly based on the numerical formulation of the governing equations, permitting a globally conservative, high-precision analysis of all physical and numerical/artificial energy transfers in a sub-domain with open lateral boundaries. The net primary energy balances are quantified using a moving average of length two tidal periods in a simplified control simulation using a single time-step, minimal diffusion, and a no-slip sea floor. This provides the basis for analysis of enhanced vertical and horizontal diffusion and a free-slip bottom boundary condition. After a four tidal period spin-up, the tidally averaged (net) primary energy balance in the generation region, extending ±20 km from the ridge crest, shows that the surface tidal wave loses approximately C = 720 W/m or 0.3% of the mean surface tidal energy flux (2.506 × 10⁵ W/m) in traversing the ridge. This corresponds mainly to the barotropic-to-baroclinic energy conversion due to stratified flow interaction with sloping topography. Combined with a normalised net advective flux of baroclinic potential energy of 0.9 × C this causes a net local baroclinic potential energy gain of 0.72 × C and a conversion into baroclinic kinetic energy through the baroclinic buoyancy term of 1.18 × C. Tidally averaged, about 1.14 × C is radiated into the abyssal ocean through the total baroclinic flux of internal pressure associated with the IT- and background density field. This total baroclinic pressure flux is therefore not only determined by the classic linear surface-to-inter nal tide conversion, but also by the net advection of baroclinic (background) potential energy, indicating the importance of local processes other than linear IT-motion. In the propagation region (PR), integrated over the areas between 20 and 40 km from the ridge crest, the barotropic and baroclinic tide are decou pled. The net incoming total baroclinic pressure flux is balanced by local potential energy gain and out-ward baroclinic flux of potential energy associated with the total baroclinic density. The primary net energy balances are robust to changes in the vertical diffusion coefficient, whereas relatively weak hor izontal diffusion significantly reduces the outward IT energy flux. Diapycnal mixing due to vertical diffu sion causes an available potential energy loss of about 1% of the total domain-averaged potential energy gain, which matches $\frac{k_m-1}{k_m} \rho_0 K^V N^2$ to within 0.5%, for $k_m$ linearly distributed grid-levels and constant back ground density $\rho_0$, vertical diffusivity ($K^V$) and buoyancy frequency ($N$).

## 1. Introduction

The world's oceans dissipate about 3.6 TW of the tidal energy in the earth-moon-sun system, of which 2.54 TW is associated with the semi-diurnal (M2) surface tide (Cartwright and Ray, 1991; Egbert and Ray, 2001). By interaction with sloping sea floors, roughly 1 TW of M2 surface tidal energy is converted into M2 internal tides (Egbert and Ray, 2001; Wunsch and Ferrari, 2004), which are observed to propagate for thousands of kilometres from

their generation regions, such as the Hawaiian ridge (Ray and Mitchum, 1997). Thus, internal tides (IT) appear to play a fundamental role in the meridional overturning circulation, providing the abyssal ocean with half of the small-scale mixing energy required to maintain the observed global stratification (Munk and Wunsch, 1998; Wunsch and Ferrari, 2004). The required energy cascade from large to small-scale internal waves may occur through the formation of higher harmonics due to wave-wave interactions, parametric sub-harmonic instability (Gerkema et al., 2006a), reflection of IT beams (Gerkema et al., 2006b) and fission of low-mode internal waves into solitons (e.g. Gerkema and Zimmerman, 1995; Shaw et al., 2009). Small-scale internal wave breaking, resulting from a multitude of processes reviewed by

<footnote>* Corresponding author. Tel.: +33 6 52 48 11 77.
E-mail address: jochem.w.floor@gmail.com (J.W. Floor).</footnote>
Thorpe (2005), causes overturning of isopycnal surfaces and irreversible diapycnal transfer of heat and salinity.

Even with the advent of dedicated large-scale campaigns (Lien and Gregg, 2001; Rudnick et al., 2003), existing in situ measurements cannot cover the IT energy transfers through the wide range from basin scale tidal dynamics to millimetric mixing. Altimetry and SAR observations provide ample data on propagation of low-mode internal tides and solitons, but cannot tell us much about weaker, high-mode internal tide motion (Egbert and Ray, 2001; Apel, 1987; Hyder et al., 2005; Niwa and Hibiya, 2001).

Time-averaged tidal barotropic-to-baroclinic energy conversion is often assumed to match the local internal tide energy flux divergence, which is modelled as

$$ \langle \vec{\nabla} \cdot (\vec{v}'p') \rangle_T = -\rho_0 \langle b' \bar{v}_z \rangle_T, \quad (1.1) $$

where $\vec{v}', p'$, and $b'$ are the baroclinic velocity, pressure and buoyancy anomaly and $\bar{v}_z$ the vertical velocity induced by barotropic tidal flow over sloping topography, while $\rho_0$ is a constant background density and $\langle \cdot \rangle_T$ indicates the temporal average (Niwa and Hibiya, 2001; Holloway and Merrifield, 1999; Gerkema et al., 2004; Di Lorenzo et al., 2006). This formulation is inappropriate for large-amplitude internal tides, which would require inclusion of (non-linear) advective energy fluxes (Scotti et al., 2006; Lamb, 2007) and does not explicitly represent the energy lost by the surface tidal wave.

Recent theoretical estimates are constrained to two-dimensional motion, the rigid-lid or infinite-depth assumptions, inviscid, non-diffusive and (weakly non-) linear cases (Baines, 1982; Bühler and Muller, 2007; Pétrélis et al., 2006; Khatiwala, 2003; Llewellyn Smith and Young, 2002; Lamb, 2007). Viscosity is essential to internal tide theory, which otherwise contains small-scale singularities in the internal tide beams (Pétrélis et al., 2006). Peacock et al. (2008) showed good correspondence between laboratory experiments and weakly viscous, linear theory, but were limited to either sub-critical or knife-edge topography of small height compared to the fluid depth. Non-linearity, viscosity and inclusion of a free surface to explicitly calculate tidal energy conversion render the theory intractable.

A global-scale model that accurately represents fine-scale mixing is not currently within reach, although low-resolution global simulations point to primary IT-production and conversion sites (Simmons, 2008). Recent regional-scale efforts focused on M2-tidal conversion and internal tide energy flux (Holloway and Merrifield, 1999; Munroe and Lamb, 2005; Di Lorenzo et al., 2006; Katsumata, 2006; Carter et al., 2008). Carter et al. (2008) used POM to estimate the M2 internal tide energetics around the Hawaiian ridge, employing the sum of kinetic and linearised available potential energy, which is appropriate for linear stratification and small-amplitude internal waves, discretised on the horizontal and vertical mid-point of the model grid-cells. Although they concluded that local dissipation is a non-negligible factor in the IT generation zone, they observed an error in both the barotropic and baroclinic global energy balance on the order of 10% of the primary energy conversions.

We propose a complete energetics analysis of the internal tide from generation to mixing, explicitly based on the numerical scheme of the energy-conserving, Boussinesq, hydrostatic, free-surface, terrain-following (σ-) coordinate ocean model Symphonie, which is similar to POM and described in detail by Marsaleix et al. (2008). We evaluate all physical and numerical energy transfers, permitting quantification, control and motivation of the artificial energy transfers due to model choices, for instance temporal diffusion and discretisation onto the C-grid (Arakawa and Lamb, 1977). Expanding upon the global formulation for closed basins by Marsaleix et al. (2008), we also consider lateral boundary fluxes

in sub-domains, e.g. to estimate internal tide energy radiation from the generation region.

We analyse M2 internal tide generation by a barotropic surface wave that impinges on a supercritical Gaussian ridge in the rotating stratified deep ocean, a representative case related to the Hawaiian Island chain, previously studied in detail by Holloway and Merrifield (1999), Munroe and Lamb (2005) and Lamb (2007). We aim to quantify numerically (a) the primary energy balances in internal tide generation and in IT propagation away from topography and (b) the associated energy expended in diapycnal diffusion. The energy lost to enhanced diapycnal mixing associated with the internal tide is evaluated following Winters et al. (1995), using a novel adiabatic redistribution algorithm adapted to terrain-following coordinates and a free surface. Finally, we underline the importance of respecting the discretised model formulation for precise energy diagnostics.

The paper is organised as follows: In Section 2, the model setup is introduced and following parameter space analysis the resulting internal tide field is characterised. In Section 3, the energy evolution equations are presented, which are adapted to internal tide generation in Section 4. In Section 5, we show the energy balance is closed, analyse the numerical energy transfers and the primary energy transfers in internal tide generation and propagation and evaluate the effects of enhanced vertical and horizontal diffusion and a free-slip bottom boundary condition. We discuss the results and conclude and present perspectives in Section 6.

## 2. Numerical simulation of the internal tide

The energy transfers inherent in internal tide generation by free-surface tidal flow over a Gaussian ridge are studied using the hydrostatic approximation, on the traditional f-plane. The model equations are presented in Appendix A.1 and discretised on the C-grid in the horizontal and on the Lorenz grid in the vertical, using second order accurate schemes for advection of momentum and tracers and an explicit leapfrog scheme in time, using a single time-step. Vertical diffusion is calculated using an implicit scheme. A Robert-Asselin filter is applied to the momentum and tracer equations to limit high-frequency numerical noise associated with the leapfrog scheme (Robert, 1966; Asselin, 1972; Marsaleix et al., 2008, Section 3.2.2). A control simulation is defined in Section 2.1, aimed primarily at the analysis of mechanical energy conversion from the barotropic to the baroclinic tide. Therefore, numerical and physical viscosity and tracer diffusivity are minimised, while maintaining model stability.

### 2.1. Model setup

The numerical model is set up in the zonal/vertical Oxz-plane and is uniform and cyclic in the Oy-direction, permitting the Coriolis effect (cf. Appendix A.2), with $f = 10^{-4}$ s$^{-1}$, corresponding to a latitude of 43.4°N. The numerical domain extends 1200 km zonally, with ambient depth $H = 5$ km and a Gaussian ridge defined by

$$ h = h_0 \exp \left\{ - \frac{(x - x_0)^2}{a^2} \right\}, \quad (2.1) $$

centred at $x_0 = 600$ km from the western boundary, and of characteristic height $h_0 = 1500$ m and width $a = 6.45$ km. The horizontal resolution is $\Delta x = 1$ km and in the vertical 40 linearly distributed $\sigma$-layers are used.

Initially, the fluid is salt-stratified with a constant buoyancy frequency $N = 10^{-3}$ s$^{-1}$. With the linear equation of state (A5), the heat and salinity Eqs. (A4a) and (A4b) reduce to a single evolution equation for the density anomaly
$$ \dot{\rho} = \rho - \rho_0, \qquad (2.2) $$

where $\rho_0 = 1.029056$ kg m$^{-3}$ is an (arbitrary) constant background value. In the control simulation only constant vertical turbulent viscosity and diffusivity are used, $K^V = 10^{-6}$ m$^2$ s$^{-1}$, for momentum and scalars.

Focusing on tidal dynamics, all free surface fluxes are (arbitrarily) set to zero. The bottom boundary condition is no-slip, so that at $z = -H + h(x)$ the velocity vanishes, while temperature and salinity fluxes through the sea floor are zero.

An M2-tidal surface wave of frequency $\omega = 1.4075 \times 10^{-4}$ s$^{-1}$ and period $T = 12.4$ h is forced on the west boundary, at $x = 0$, by specifying the depth-averaged zonal current as $\bar{v}_x^F = -U_0 \sin(\omega t)$, with $U_0 = 0.025$ ms$^{-1}$ and assuming the free surface elevation behaves as the linear, inviscid shallow water solution $\eta = \bar{v}_x \sqrt{(1 - \mu)H/g}$, where $\mu = f^2/\omega^2$ represents the strength of Coriolis dispersion for (arbitrary) frequency $\omega$.

The method of characteristics (e.g. Johns et al., 1983; Marsaleix et al., 2009a) is used to radiate the surface tidal wave out of the domain at the eastern boundary and (small) surface wave anomalies with respect to the forcing at the western boundary (cf. Appendix A.2). The simulation is terminated before any reflected internal wave might reach the domain of interest. The model characteristics are summarised in Table 2.1.

## 2.2. Parameter space

The internal tide is characterised by the surface tidal amplitude $U_0$ and frequency $\omega$, the Coriolis frequency $f$, the strength of the stratification through $N$, the characteristic height $h_0$ and width $a$ of the topography and the ambient fluid depth, $H$ (e.g. Garrett and Kunze, 2007; Legg and Huijts, 2006). These can be organised into the following set of non-dimensional parameters: The tidal excursion parameter, $U_0/(\omega a)$, the topographic slope $h_0/a$, the angle of inclination $\theta_T = \tan^{-1}(h_0/a)$, the relative height of the topography with respect to the surrounding fluid depth $h_0/H$, and the vertical Froude number $Fr = U_0/(Nh_0)$. The final non-dimensional parameter is the slope of the internal wave beams at the tidal frequency, $\sqrt{(\omega^2 - f^2)/(N^2 - \omega^2)} = \tan \theta$ in the non-hydrostatic case, where $\theta$ is the angle with respect to the horizontal of the internal wave beam, along which the internal tide energy propagates. If the internal waves are not scattered or dissipated in the ocean interior, they

<table><caption>**Table 2.1**
Summary of physical and numerical model characteristics.</caption><thead><tr><th>Parameter</th><th>Symbol</th><th>Control simulation</th></tr></thead><tbody><tr><td colspan="3"><strong>(a) Physical parameters</strong></td></tr><tr><td>Domain length</td><td>L (km)</td><td>1.2 × 10<sup>3</sup></td></tr><tr><td>Ambient depth</td><td>H (m)</td><td>5 × 10<sup>3</sup></td></tr><tr><td>Ridge height</td><td>h<sub>0</sub> (m)</td><td>1.5 × 10<sup>3</sup></td></tr><tr><td>Ridge e-folding width</td><td>a (m)</td><td>6.45 × 10<sup>3</sup></td></tr><tr><td>Coriolis parameter</td><td>f (s<sup>-1</sup>)</td><td>10<sup>-4</sup></td></tr><tr><td>Reference density</td><td>ρ<sub>0</sub> (kg m<sup>-3</sup>)</td><td>1.029056 × 10<sup>3</sup></td></tr><tr><td>Buoyancy frequency</td><td>N (s<sup>-1</sup>)</td><td>10<sup>-3</sup></td></tr><tr><td>Vertical diffusion coefficient</td><td>K<sup>V</sup> (m<sup>2</sup> s<sup>-1</sup>)</td><td>10<sup>-6</sup></td></tr><tr><td>Horizontal diffusion coefficient</td><td>K<sup>H</sup> (m<sup>2</sup> s<sup>-1</sup>)</td><td>0 (for v<sub>z</sub>, T, and S)</td></tr><tr><td>Tidal velocity amplitude</td><td>U<sub>0</sub> (ms<sup>-1</sup>)</td><td>2.5 × 10<sup>-2</sup></td></tr><tr><td>Tidal period</td><td>T (h)</td><td>12.4</td></tr><tr><td colspan="3"><strong>(b) Numerical parameters</strong></td></tr><tr><td>Zonal grid-cells</td><td>M</td><td>1201</td></tr><tr><td>Meridional grid-cells</td><td>N</td><td>3</td></tr><tr><td>Mass vertical levels</td><td>k<sub>m</sub></td><td>40</td></tr><tr><td>Velocity vertical levels</td><td>k<sub>v</sub></td><td>41</td></tr><tr><td>Horizontal grid-scale</td><td>Δx (m)</td><td>1000</td></tr><tr><td>Vertical grid-scale</td><td>Δσ</td><td>1/40</td></tr><tr><td>Time-step</td><td>Δt (s)</td><td>1.437</td></tr><tr><td>Asselin coefficient</td><td>A</td><td>0.1</td></tr></tbody></table>

are reflected at the free surface or ocean floor so that the energy propagates horizontally along the oceanic wave-guide.

In principle, internal waves are generated at the tidal and its harmonic frequencies (Bell, 1975), but for small tidal excursion compared to the width of the topography the tidal frequency dominates, which is the case here because $U_0/(\omega_0 a) = 0.03 \ll 1$. The relative height of the topography is $h_0/H = 0.3$. The non-linearity of the tidal dynamics can be characterised by the product of the tidal excursion parameter and the relative height of the topography. Here, $U_0 h_0/(a\omega_0 H) = 0.008 \ll 1$ and it is safe to say that non-linearity will not play an important role. The internal tide will be predominantly generated at the M2 tidal forcing frequency.

For large Froude numbers, the flow is relatively unaffected by the sloping topography, whereas if $Fr \ll 1$ the topography acts as a barrier to the flow. At present, blocking effects may play a role because $Fr = 0.02$, and the internal tide generation by oscillating tidal flow perpendicular to the ridge should be relatively strong.

The ratio of the maximum topographic angle with respect to the horizontal to that of the internal tide beam indicates whether the topography is sub-critical ($\theta_T/\theta < 1$) or super-critical ($\theta_T/\theta > 1$). Internal tide generation is significantly stronger for super-critical than for sub-critical topography because there is always a (critical) point on the slope where free internal gravity waves can resonate with the tidal frequency. For super-critical topography, internal tide beams are observed to emanate from the critical points, where the topographic and beam slope are equal (e.g. Gerkema, 2006c). Here, the greatest topographic angle of inclination, $\theta_T = 11.3^\circ$, exceeds the angle with respect to the horizontal of the internal tide beam, $\theta = 5.7^\circ$, so that the topography is supercritical and the generated internal tide should be relatively strong (e.g. Munroe and Lamb, 2005). Regimes ranging from strongly sub-critical to strongly super-critical were studied by Floor (2009).

For internal waves, non-hydrostatic effects are likely to be small if $\omega^2/N^2 \ll 1$. In this limit the slope of internal tide beams is $\sqrt{(\omega^2 - f^2)/N^2} = \tan \theta$. The low tidal frequency leads to $\omega^2/N^2 = 0.014$, which implies that the hydrostatic approximation used in Symphonie is appropriate for the M2 internal tide. Hydrostatic tidal beams have a slope that is 1% smaller than in the non-hydrostatic case. Non-hydrostatic effects are important in small-scale processes such as wave breaking, which we do not consider explicitly (see e.g. Lamb, 2007; Auclair et al., 2010).

## 2.3. Velocity decomposition

To distinguish (barotropic) motion associated with the free surface from that due to internal (baroclinic) motion, the horizontal velocity is decomposed as

$$ v_x = \bar{v}_x + v'_x, \qquad (2.3) $$

where horizontal direction is indicated by the subscript $\alpha = x, y$ and $v'_x$ is the anomaly with respect to the vertically averaged horizontal velocity

$$ \bar{v}_x = \frac{1}{D} \int_{-H}^{H} v_\alpha dz, \qquad (2.4) $$

where $z = \eta(x, y, t)$ is the equation of the free surface elevation with respect to its static equilibrium position ($z = 0$), and $D \equiv \eta + H - h$ the total depth.

Symphonie is formulated in terrain-following Oxyσ-coordinates (Marsaleix et al., 2008; Johns et al., 1983), where the vertical coordinate is the relative height of a fluid parcel with respect to total depth,

$$ \sigma = \frac{z + H}{D}, \qquad (2.5) $$

so that in the present model $\sigma = 0$ at the sea floor and $\sigma = 1$ at the free surface. With the same number of levels, the model grid has en-
hanced resolution over shallow bathymetry and, due to surface waves, the σ-levels move with respect to the absolute vertical coordinate, z. Some useful relations for partial derivatives are shown in Appendix A.1.

Using (2.3)-(2.5), the vertical velocity can also be split up into barotropic and baroclinic contributions, so that

$$v_z = \bar{v}_z + v'_z, \qquad (2.6)$$

where respectively

$$\bar{v}_z = \sigma \left( \frac{\partial \eta}{\partial t} + \bar{v}_x \frac{\partial \eta}{\partial x_x} \right) + (\sigma - 1) \bar{v}_x \frac{\partial H}{\partial x_x}, \qquad (2.7)$$

$$v'_z = \sigma v'_x \frac{\partial \eta}{\partial x_x} + (\sigma - 1) v'_x \frac{\partial H}{\partial x_x} + v_\sigma. \qquad (2.8)$$

Repeated dummy indices α (or β) indicate summation, while $v_σ$ is the velocity component defined by $v_σ = Ddσ/dt$, which equals zero at sea-floor and surface. The barotropic component of $v_z$, (2.7), depends linearly only on local vertical sea-surface motion, the flow of the depth-averaged current along sloping topography, which is important for the IT-generation, and along sloping sea-surface. In the absence of baroclinic motion, the relative height σ of the fluid particles is conserved (Pedlosky, 1987, Section 3.3). The baroclinic component (2.8) comprises the anomalies with respect to barotropic motion.

## 2.4. Internal tide field

Fig. 1a shows the vertical velocity field four tidal periods after the tidal surface wave first impinges on the ridge. Wave fronts corresponding to baroclinic normal modes in uniform stratification (Pedlosky, 2003) propagate horizontally away from the ridge in both directions at a phase speed corresponding to the non-rotating limit, analogous to the propagation of transients away from a sea-level discontinuity undergoing geostrophic adjustment (Gill, 1982, Section 7.3). The fronts of modes 1, 2, and 3 are visible near 270, 130 and 85 km from the ridge summit. Far behind each front, after

a period of order T near the topography, the internal wave modes have generation frequency ω and feel the effect of rotation. The rotating phase speed corresponds to $c_{p,n} = c_n/\sqrt{1-\mu}$ and e.g. the first mode wavelength is $λ_1 = 101$ km (not shown). After four tidal periods, when about five baroclinic modes have travelled 50 km away from the ridge, within this distance the familiar internal tide beam takes shape closely corresponding to the theoretical angle with respect to the horizontal of $θ = 5.7°$ and a stationary state ensues.

The difference between the external ($\bar{v}_z$) and internal ($v'_z$) components of the vertical velocity $v_z$ is evident in Panels 1b and c. Being primarily associated with the impingement of the surface tidal wave on the slope, $\bar{v}_z$ is concentrated over the topographic slope, whereas $v'_z$ is significant throughout the domain and associated with the propagating baroclinic modes. Over the topography $v'_z$ is similar in strength and in the opposite direction to $\bar{v}_z$, which represents the topographic coupling of the barotropic and baroclinic wave modes, and thus the conversion from barotropic to baroclinic wave energy.

## 3. Energy evolution equations of a water column

Under the Boussinesq approximation, the mechanical and thermodynamical energy are essentially decoupled. The mechanical energy density is the sum of kinetic ($E_K$) and potential energy density ($E_p$) and is conserved for inviscid, non-diffusive flow in a closed domain.

Our approach is explicitly based on the discretised governing equations of Symphonie, using separate discrete diagnostic equations in terrain-following Oxyσ-coordinates for the barotropic and baroclinic $E_K$- and $E_p$-evolution. This permits the evaluation of consistency of the reversible exchange between $E_K$ and $E_p$ due to vertical mass-flux, an essential feature of gravity wave motion. This consistency requirement constrains the discretised formulation of the surface and internal pressure gradient forces in the momentum Eqs. (A1) and (A2) and the advective density (tracer) flux divergence

<figure><img src="image_2.png" /><figcaption>Fig. 1. Main characteristics of the linear internal tide after four tidal periods, in a vertical section extending to ±300 km from the ridge crest. (a) $v_z$, fronts of internal wave modes 1, 2, 3 are visible near 270, 130 and 85 km from the ridge crest, (b) $\bar{v}_z$, and (c) $v'_z$.</figcaption></figure>
in the tracer evolution Eqs. (A3) and (A4), as discussed by Marsaleix et al. (2008) for a closed domain. Their Fig. 9 shows the discrete energy transfers within and between the different compartments. The discretised formulation of energy transfers in a domain with open lateral boundaries is discussed in Appendix B.

Previous studies of internal tide energetics considered an evolution equation for the sum of the kinetic energy and available potential energy (i.e. the pseudo-energy, Lamb, 2007) or linearised available potential energy (Carter et al., 2008). To obtain a rigorously closed numerical energy balance, we choose to analyse the complete baroclinic potential energy balance and use the complete baroclinic density $\bar{\rho}$, which includes the background stratification, rather than an internal wave 'density anomaly' as is commonly used in rigid lid models of IT-generation (e.g. Lamb, 2007). As a consequence, $\bar{\rho}$ may include effects due to free-surface motion, non-linearity, (sub-) harmonics and might in a realistic context further contain geostrophic and wind-induced circulations, etc. Since, wave motion is most fundamentally described in terms of pseudo-energy (Shepherd, 1993), we also consider the available potential energy loss due to diapycnal mixing in Section 4.3, following Winters et al. (1995).

In the following, we present evolution equations for kinetic and potential energy per unit horizontal area, i.e. integrated vertically over the depth of the fluid. We distinguish between (a) external (barotropic) energy associated with surface motion and (b) internal (baroclinic) energy associated with internal wave motion, using the density decomposition (2.2) and velocity decomposition (2.3). In the following, $F_q$ indicates horizontal flux divergences or net vertical boundary flux (from the vertical integral of vertical flux divergence), while local energy conversions are indicated by $\phi_q$. The subscript q indicates the type of conversion under consideration. For reference, Table A.1 summarises the energy transfers and their abbreviations.

The Boussinesq, hydrostatic kinetic energy density is defined as $E_K = \frac{1}{2} \rho_0 v_\alpha v_\alpha$ and contains contributions from horizontal velocity only. Denoting vertical integration by a hat and employing the velocity decomposition (2.3), the kinetic energy tendency per unit area in the Oxy-plane of a fixed fluid column is

$$ \frac{\partial \hat{E}_K}{\partial t} = \int_0^1 \frac{\partial E_K D}{\partial t} d\sigma = \underbrace{\rho_0 \bar{v}_\alpha \frac{\partial D \bar{v}_\alpha}{\partial t} - \bar{E}_K \frac{\partial D}{\partial t}}_{\hat{E}_K/\partial t} + \underbrace{\rho_0 \int_0^1 v'_\alpha \frac{\partial D v'_\alpha}{\partial t} d\sigma - \frac{\partial D}{\partial t} \int_0^1 E'_K d\sigma}_{\hat{E}'_K/\partial t} \quad (3.1) $$

where barotropic kinetic energy density is $\bar{E}_K = \frac{1}{2} \rho_0 \bar{v}_\alpha \bar{v}_\alpha$ and the baroclinic part $E'_K = \frac{1}{2} \rho_0 v'_\alpha v'_\alpha$. Here, the first two right-hand side (rhs)-terms represent the local barotropic kinetic energy tendency per unit area, $\partial \bar{E}_K / \partial t$, while the last two terms constitute the baroclinic kinetic energy tendency per unit area, $\partial \hat{E}'_K / \partial t$.

Firstly, vertically integrating the momentum equations for a Boussinesq fluid (A1, 2) and the continuity equation (A3) and using (3.1), the evolution equation for barotropic kinetic energy per unit horizontal area can be written as

$$ \begin{aligned} \frac{\partial \hat{E}_K}{\partial t} + \underbrace{\frac{\partial D \bar{v}_\alpha \bar{E}_K}{\partial x_\alpha}}_{\hat{F}_K} + g \rho_0 \underbrace{\frac{\partial D \bar{v}_\alpha \eta}{\partial x_\alpha}}_{\hat{F}_T} + \underbrace{\int_0^1 \frac{\partial D \bar{v}_\alpha P}{\partial x_\alpha} d\sigma}_{\hat{F}_X} \\ = -D \underbrace{\int_0^1 \bar{\rho} g \bar{v}_\alpha d\sigma}_{\phi_\alpha} + g \rho_0 \eta \underbrace{\frac{\partial D \bar{v}_\alpha}{\partial x_\alpha}}_{-\phi_T} - \underbrace{\rho_0 \bar{v}_\alpha \int_0^1 \frac{\partial D v'_\beta v'_\alpha}{\partial x_\beta} d\sigma}_{\phi_{3D2D}} \\ + \rho_0 \bar{v}_\alpha \int_0^1 \underbrace{\frac{\partial}{\partial x_\beta} \left( D K^H \frac{\partial v_\alpha}{\partial x_\beta} \right)}_{\zeta_H} d\sigma + \rho_0 \bar{v}_\alpha \underbrace{\left[ \frac{K^V}{D} \frac{\partial v_\alpha}{\partial \sigma} \right]_0^1}_{\zeta_V} \end{aligned} \quad (3.2) $$

Eq. (3.2) involves the following energy transfers: $F_{E_K}$ is the divergence of the barotropic flux of barotropic kinetic energy, while

the sum ($F_T + F_X$) represents the divergence of the pressure flux due to the barotropic tidal current. $F_T$ is the divergence of the linear energy flux corresponding to small-amplitude, shallow water surface waves (e.g. Gill, 1982, Section 5.7; Floor, 2009, Section B.4). $F_X$ is the divergence of the barotropic flux of internal pressure, associated with the ambient stratification and the internal wave pressure anomaly. The energy conversion between barotropic kinetic energy and baroclinic potential energy is $\phi_\alpha$, the vertical integral of the barotropic buoyancy flux $\bar{\rho} g \bar{v}_\alpha$, and is hereinafter referred to as the barotropic buoyancy term (similarly for the total and baroclinic buoyancy terms). Furthermore, $\phi_T$ represents a sink of $\bar{E}_K$ equal to the local gain of barotropic potential energy. The term $\phi_{3D2D}$ represents non-linear interaction between baroclinic and barotropic motion. Finally, barotropic kinetic energy is dissipated through the mixing terms $\bar{\zeta}_H$ and $\bar{\zeta}_V$, including mixing parameterisation, numerical diffusion and bottom friction, which we will turn to in Section 5.5.

Secondly, the tendency of baroclinic kinetic energy per unit area is similarly expressed as

$$ \begin{aligned} \frac{\partial \hat{E}'_K}{\partial t} + \underbrace{\int_0^1 \frac{\partial D v'_\beta E'_K}{\partial x_\beta} d\sigma}_{F_{E'_K}} + \underbrace{\int_0^1 \frac{\partial D \bar{v}_\beta E'_K}{\partial x_\beta} d\sigma}_{F_{X1}} + \underbrace{\int_0^1 \frac{\partial D v'_\alpha P}{\partial x_\alpha} d\sigma}_{F_{IW}} \\ + \rho_0 \underbrace{\int_0^1 \frac{\partial D v'_\beta v'_\alpha \bar{v}_\alpha}{\partial x_\beta} d\sigma}_{F_{X2}} - \underbrace{D \int_0^1 \bar{\rho} g v'_z d\sigma}_{\phi'_z} + \rho_0 \bar{v}_\alpha \underbrace{\int_0^1 \frac{\partial D v'_\beta v'_\alpha}{\partial x_\beta} d\sigma}_{\phi_{3D2D}} \\ + \rho_0 \underbrace{\int_0^1 v'_\alpha \frac{\partial}{\partial x_\beta} \left( D K^H \frac{\partial v_\alpha}{\partial x_\beta} \right)}_{\zeta'_H} d\sigma + \rho_0 \underbrace{\int_0^1 v'_\alpha \frac{\partial}{\partial \sigma} \left( \frac{K^V}{D} \frac{\partial v_\alpha}{\partial \sigma} \right)}_{\zeta'_V} d\sigma. \end{aligned} \quad (3.3) $$

Here, $F_{E'_K}$ is the divergence of the baroclinic advective flux of baroclinic kinetic energy ($E'_K$), whereas $F_{X1}$ is the barotropic $E'_K$-flux divergence. $F_{IW}$ is the divergence of the flux of internal pressure due to baroclinic motion, closely related to the classical definition of the internal tide flux divergence (see Eq. (1.1)) but here it is important to stress that $F_{IW}$ also includes the background stratification. The term $F_{X2}$ represents non-linear interaction between barotropic and baroclinic motion. The energy conversion between baroclinic kinetic and potential energy is $\phi'_z$, which is the vertical integral of the buoyancy flux due to baroclinic vertical motion ($\bar{\rho} g \bar{v}_z$). The term $\phi_{3D2D}$ represents an advective exchange with barotropic kinetic energy, where it appears with opposite sign. Finally, $\bar{\zeta}_H$ and $\bar{\zeta}_V$ are the energy conversions associated with horizontal and vertical diffusion of baroclinic momentum, where the latter includes sea-floor friction effects.

Using the density decomposition (2.2), the potential energy tendency per unit area of a fixed fluid column is

$$ \frac{\partial \hat{E}_P}{\partial t} = \rho_0 g \eta \frac{\partial \eta}{\partial t} + \underbrace{\int_0^1 D \bar{\rho} g \frac{\partial z}{\partial t} d\sigma}_{\hat{\bar{E}}_P/\partial t} + \underbrace{\int_0^1 gz \frac{\partial D \bar{\rho}}{\partial t} d\sigma}_{\hat{\bar{E}}_P/\partial t} \quad (3.4) $$

where the barotropic potential energy density purely associated with surface motion is $\bar{E}_P = \rho_0 g z$ and the baroclinic part is given by $\bar{E}_P = gz \bar{\rho}$, which includes the background stratification and density anomalies associated with internal wave motion. This decomposition corresponds exactly to the formulation of the numerical model and is required to close the energy balance including diffusive effects.

Using the vertically integrated continuity equation, it can be shown that the local barotropic potential energy gain corresponds to
$$ \frac{\partial \hat{E}_p}{\partial t} = \phi_T, \qquad (3.5) $$

the rhs of which appears in the barotropic kinetic energy equation (3.2). The two equivalent terms in Eqs. (3.2) and (3.5) are maintained to emphasise that they originate from different governing equations, respectively the barotropic momentum and continuity equation, which is essential in the derivation of the equivalent discrete energy balances.

Finally, the baroclinic potential energy tendency equation can be written in a form similar to (3.2), namely

$$ \begin{aligned} & \frac{\partial \hat{E}_p}{\partial t} + \underbrace{\int_0^1 \frac{\partial D v_\beta \tilde{E}_p}{\partial x_\beta} d\sigma}_{F_{E_p} = \bar{F}_{E_p} + F'_{E_p}} - \underbrace{\int_0^1 \frac{\partial}{\partial x_\beta} \left( gzDK_\beta^H \frac{\partial \tilde{\rho}}{\partial x_\beta} \right) d\sigma}_{F_{m,H}} - \underbrace{\left[ \frac{gz K^V}{D} \frac{\partial \tilde{\rho}}{\partial \sigma} \right]_0^1}_{F_{m,V}} \\ & = D \underbrace{\int_0^1 \tilde{\rho} g v_z d\sigma}_{\phi_z = \phi_z + \phi'_z} - \underbrace{\int_0^1 g DK_\beta^H \frac{\partial z}{\partial x_\beta} \frac{\partial \tilde{\rho}}{\partial x_\beta} d\sigma}_{\phi_{m,H}} - \underbrace{\int_0^1 g K^V \frac{\partial \tilde{\rho}}{\partial \sigma} d\sigma}_{\phi_{m,V}}. \end{aligned} \qquad (3.6) $$

Here, $F_{E_p}$ is the horizontal divergence of the sum of the horizontal baroclinic and barotropic flux of potential energy, while $F_{m,H}$ ($F_{m,V}$) represents the horizontal (vertical) divergence of the density diffusion flux. Local exchanges with kinetic energy are the buoyancy term ($\phi_z$), and 'horizontal' diffusion along $\sigma$-surfaces and vertical diffusion ($\phi_{m,H}$ and $\phi_{m,V}$). The latter are discussed in Section 6.

$K^H$ corresponds to diffusion along $\sigma$-surfaces and primarily acts over sloping topography. Although strictly an abuse of language, we use 'horizontal diffusion' in this paper for convenience. Mellor and Blumberg (1985) suggest that preferably $K^H = 0$, to preserve temperature and salinity fronts, and should only be non-zero to limit noise. Non-zero $K^H$ is nonetheless employed regularly so it is interesting to study its effect on internal tide generation. In particular, the 'horizontal' diffusion of density along sloping $\sigma$-surfaces may change the potential energy $\bar{E}_p$ even if the density stratification only depends on depth, since in case of a sloping sea surface or -floor the 'horizontal' derivative of both depth, $\partial z/\partial x_\alpha$, and density, $\partial \bar{\rho}/\partial x_\alpha$, is non-zero.

Let us summarise this section by listing the shorthand evolution equations for vertically integrated energy per unit area,

$$ \frac{\partial \hat{E}_K}{\partial t} + F_{\bar{E}_K} + F_T + F_X = -\bar{\phi}_z - \phi_T - \phi_{3D2D} + \bar{\xi}_H + \bar{\xi}_V, \qquad (3.7a) $$

$$ \frac{\partial \hat{E}_p}{\partial t} = \phi_T, \qquad (3.7b) $$

$$ \frac{\partial \hat{E}_p}{\partial t} + F_{\bar{E}_p} - F_{m,H} - F_{m,V} = \bar{\phi}_z + \phi'_z - \phi_{m,H} - \phi_{m,V}, \qquad (3.7c) $$

$$ \frac{\partial \hat{E}'_K}{\partial t} + F'_{\bar{E}_K} + F_{X1} + F_{IW} + F_{X2} = -\phi'_z + \phi_{3D2D} + \xi'_H + \xi'_V. \qquad (3.7d) $$

Clearly, the baroclinic potential energy equation may act as an intermediary between the barotropic and baroclinic kinetic energy equations, through the barotropic ($\bar{\phi}_z$) and baroclinic buoyancy term ($\phi'_z$). Carter et al. (2008) considered energy balances for the sum of kinetic and (linearised available) potential energy, of the form

$$ \frac{\partial \hat{E}}{\partial t} + F_{\bar{E}} + F_T + F_X = -\bar{\phi}_z - \phi_{3D2D} + \bar{\xi}_H + \bar{\xi}_V, \qquad (3.8a) $$

$$ \begin{aligned} & \frac{\partial \hat{E}'}{\partial t} + F'_{\bar{E}} + F_{X1} + F_{IW} + F_{X2} - F_{m,H} - F_{m,V} \\ & = \bar{\phi}_z + \phi_{3D2D} + \xi'_H + \xi'_V - \phi_{m,H} - \phi_{m,V} \end{aligned} \qquad (3.8b) $$

which are simpler in form than (3.7), but do not allow the analysis of the redistribution of energy in the potential energy compartment.

Here, $\hat{\bar{E}} = \hat{\bar{E}}_K + \hat{\bar{E}}_p$ and $\hat{\bar{E}}' = \hat{\bar{E}}'_K + \hat{\bar{E}}'_p$, while the barotropic buoyancy conversion term ($\bar{\phi}_z$) and the advective term $\phi_{3D2D}$ directly carry energy between the barotropic and baroclinic energy compartments.

The evolution equations of total kinetic (3.9a) and potential energy per unit area (3.9b) are obtained by summing (3.7a) and (3.7d), respectively (3.7b) and (3.7c), giving

$$ \frac{\partial \hat{E}_K}{\partial t} + F_{\bar{E}_K} + F_T + (F_X + F_{IW}) = -\phi_z - \phi_T + \xi_H + \xi_V, \qquad (3.9a) $$

$$ \frac{\partial \hat{E}_p}{\partial t} + F_{\bar{E}_p} - F_{m,H} - F_{m,V} = \phi_z + \phi_T - \phi_{m,H} - \phi_{m,V}, \qquad (3.9b) $$

where $E_K = \bar{E}_K + E'_K$ and $E_p = \bar{E}_p + \tilde{E}_p$ and the advective flux divergence of kinetic energy is $F_{\bar{E}_K} = F'_{\bar{E}_K} + F'_{\bar{E}_p} + F_{X1} + F_{X2}$. Let us compare Eq. (3.9) to the mechanical energy equations in Cartesian Oxyz-coordinates (e.g. Gill, 1982, Section 6.7). Without detailing viscosity and diffusion for simplicity, assuming that the pressure at the free surface equals zero and that the sea-floor is fixed and no-slip, in Cartesian coordinates, it can be shown (cf. Appendix A.4) that the evolution equations for total kinetic and potential energy integrated over the depth of the fluid, from $z = -H$ to $z = \eta$, are

$$ \begin{aligned} & \frac{\partial}{\partial t^+} \int_{-H}^{\eta} E_K dz + \frac{\partial}{\partial x_\alpha^+} \int_{-H}^{\eta} (v_\alpha (E_K + \rho_0 g \eta + P)) dz \\ & = - \int_{-H}^{\eta} \tilde{\rho} g v_z dz - g \rho_0 \eta \frac{\partial D \bar{v}_\alpha}{\partial x_\alpha^+} + \xi^*, \end{aligned} \qquad (3.10a) $$

$$ \begin{aligned} & \frac{\partial}{\partial t^+} \int_{-H}^{\eta} E_p dz + \frac{\partial}{\partial x_\alpha^+} \int_{-H}^{\eta} v_\alpha \tilde{E}_p dz - F_m^+ \\ & = \int_{-H}^{\eta} \tilde{\rho} g v_z dz + g \rho_0 \eta \frac{\partial D \bar{v}_\alpha}{\partial x_\alpha^+} - \phi_m^+. \end{aligned} \qquad (3.10b) $$

Asterisks denote partial derivatives holding $z$ constant (cf. (A6)) and total viscous and diffusion terms are denoted $\xi^*$, respectively $\phi_m^*$ and $F_m^*$. Since horizontal derivatives of z-independent variables are identical in $\sigma$- and $z$-coordinates, Eq. (3.10) correspond exactly to Eq. (3.9). The kinetic energy density may increase by a local convergence of the kinetic energy density flux ($F_{\bar{E}_K} = \frac{\partial}{\partial x_\alpha^+} \int_{-H}^{\eta} v_\alpha E_K dz$) and/or a local convergence of the flux associated with internal pressure. Because the total horizontal velocity $v_\alpha$ is considered, the divergence of the vertically integrated flux of total internal pressure ($\frac{\partial}{\partial x_\alpha^+} \int_{-H}^{\eta} v_\alpha P dz$) in (3.10a) corresponds to the sum ($F_X + F_{IW}$) in (3.9a). The buoyancy conversion term $\phi_z$ in (3.9a) and (3.9b) is exactly equivalent to the term $\int_{-H}^{\eta} \tilde{\rho} g v_z dz$ in (3.10a) and (3.10b): upward (downward) motion of fluid parcels causes a sink (source) of kinetic energy density to (from) potential energy. The term $\phi_T$ is the standard form of the conversion between the barotropic part of kinetic and potential energy (e.g. Pedlosky, 1987, Section 3.6) appearing as $g \rho_0 \eta \partial D \tilde{v}_\alpha / \partial x_\alpha^+$ in (3.10a) and (3.10b). The divergence of potential energy ($F_{\bar{E}_p}$) includes only the baroclinic potential energy associated with $\tilde{\rho}$, because the divergence of barotropic potential energy is incorporated in $\phi_T$ (cf. Floor, 2009, Section B.4).

## 4. Energetics of internal tide generation and propagation

For a flow dominated by periodic wave motion, the sign of the instantaneous reversible energy transfers depends on the wave phase. When averaged over several tidal periods, relatively small but irreversible energy transfers such as the energy conversions from one wave-type to another that radiates away, or associated with diffusion, may be significant compared to large amplitude reversible transfers. We integrate the energy transfers horizontally
over a sub-domain and apply a moving average over two tidal periods. Thus, we consider in the meridionally uniform case

$$ \langle \phi_q \rangle = \frac{1}{2T} \int_{t-T}^{t+T} \int_{X_L}^{X_R} \phi_q(x, t') dx dt', $$

where $\langle \cdot \rangle$ indicates spatial integration from the western lateral boundary at $X_L$ to the eastern one at $X_R$ in addition to a time-average over $t-T \le t < t+T$. In what follows, $\langle F_q \rangle$ denotes a net outward flux through the boundaries of the domain. We also use $\langle \cdot \rangle_X$ to indicate spatial integration only.

### 4.1. Exchange between barotropic and baroclinic wave energy

We are interested in the tidal average of the conversion mechanisms described in Section 3, integrated over a finite-width domain centred on sloping topography defined by (2.1). These govern the redistribution of incoming barotropic tidal energy in the form $-\langle F_T + F_X \rangle$, into local $\bar{E}_K$, $\bar{E}_P$, $\bar{E}_P$ and $E'_K$ as well as net outward baroclinic and cross-fluxes of energy and pressure, such as the potential energy flux $\langle F_{\bar{E}_P} \rangle$, and the baroclinic pressure flux, $\langle F_{\bar{I}W} \rangle$. Applying the $\langle \cdot \rangle$-operator to the energy Eq. (3.7) the tidally averaged energy equations are written

$$ \left\langle \frac{\partial \hat{E}_K}{\partial t} \right\rangle + \left\langle F_{\bar{E}_K} \right\rangle + \langle F_T + F_X \rangle = -\langle \bar{\phi}_z \rangle - \langle \phi_T \rangle - \langle \phi_{3D2D} \rangle + \langle \bar{\xi}_H + \bar{\xi}_V \rangle, \quad (4.1a) $$

$$ \left\langle \frac{\partial \hat{E}_P}{\partial t} \right\rangle + \left\langle F_{\bar{E}_P} \right\rangle + \langle F_{m,H} \rangle + \langle F_{m,V} \rangle = \langle \bar{\phi}_z \rangle + \langle \phi'_z \rangle - \langle \phi_{m,H} \rangle - \langle \phi_{m,V} \rangle, \quad (4.1b) $$

$$ \left\langle \frac{\partial \hat{E}'_K}{\partial t} \right\rangle + \left\langle F_{\bar{E}_K} \right\rangle + \langle F_{X1} \rangle + \langle F_{\bar{I}W} \rangle + \langle F_{X2} \rangle = -\langle \phi'_z \rangle + \langle \phi_{3D2D} \rangle + \langle \bar{\xi}'_H \rangle + \langle \bar{\xi}'_V \rangle. \quad (4.1c) $$

The term $-\langle F_T + F_X \rangle$ on the left-hand side (lhs) of (4.1a) represents the net local gain of barotropic energy in the domain, lost by the surface wave as it passes through the domain, which may cause a local change of $\bar{E}_K$, while $\langle \phi_T \rangle$ corresponds to the net local exchange with $\bar{E}_P$ associated with the propagation of a linear shallow water wave (e.g. Gill, 1982, Section 5.7). In the absence of stratification, topography and viscosity, the tidally averaged energy flux divergence $\langle F_T \rangle_T$ is theoretically zero, and the surface tidal wave propagates through the domain unhindered causing no net local energy gain. With non-zero stratification, the surface tidal wave also does work against the baroclinic (background plus internal tide) pressure field, so that the barotropic pressure flux divergence is $\langle F_T + F_X \rangle$. Besides that, $\langle \phi_{3D2D} \rangle$ represents non-linear energy exchange between $\bar{E}_K$ and $E'_K$.

The buoyancy term $\langle \phi_z \rangle = \langle \bar{\phi}_z \rangle + \langle \phi'_z \rangle$ is the fundamental reversible energy conversion in gravity wave motion, and redistributes the incoming barotropic energy between the compartments of $\bar{E}_K$, $\bar{E}_P$ and $E'_K$. When the hydrostatic approximation is made, the Cartesian vertical velocity $v_z$ is a diagnostic variable, which must be derived from the vertical integral of the continuity equation (A3), which is a fundamental constraint in Boussinesq models. Because $\phi_z$ requires $v_z$, in hydrostatic models the continuity equation plays a pivotal role in determining the conversion between kinetic and potential energy. A similar argument applies to the barotropic and baroclinic constituents of vertical velocity, $v_z$ and $v'_z$, and the buoyancy term, $\bar{\phi}_z$ and $\phi'_z$.

Using Eq. (2.7), the barotropic buoyancy term can be decomposed as

$$ \bar{\phi}_z = \underbrace{\int_0^1 D \hat{\rho} g \sigma \frac{\partial \eta}{\partial t} d\sigma}_{\phi_S} + \underbrace{\int_0^1 D \hat{\rho} g \bar{v}_\alpha \frac{\partial z}{\partial x_\alpha} d\sigma}_{\phi_B}, \quad (4.2a) $$

consisting of the baroclinic potential energy gain (loss) due to a raising (lowering) of fluid parcels by free-surface motion, $\phi_S$, and vertical motion due to the barotropic flow along sloping $\sigma$-surfaces, respectively $\phi_B$, depending on both the free-surface and sea-floor slope. Similarly, using Eq. (2.8) the baroclinic buoyancy term consists of

$$ \phi'_z = \underbrace{\int_0^1 D \hat{\rho} g v_\sigma d\sigma}_{\phi_\sigma} + \underbrace{\int_0^1 D \hat{\rho} g v'_\alpha \frac{\partial z}{\partial x_\alpha} d\sigma}_{\phi_B}, \quad (4.2b) $$

the baroclinic potential energy gain (loss) due to a raising (lowering) of fluid parcels by flow across $\sigma$-surfaces, $\phi_\sigma$, and baroclinic flow along sloping $\sigma$-surfaces, $\phi'_B$. These decompositions are particular to the Oxy$\sigma$-coordinates, but are also convenient to the present study in facilitating the analysis of energy transfers between the barotropic and baroclinic waves, due to the clear-cut distinction between the contribution to $\phi_z$ by the free surface wave (4.2a) and internal wave motion (4.2b).

### 4.2. Generation and propagation region

The domain can be split up into two physically distinct regions, the 'generation region' (GR) where the bottom topography varies significantly and the 'propagation region' (PR) where the sea floor is flat. The GR reflects internal tide generation and radiation, because there is significant coupling between barotropic and baroclinic motion primarily due to the sloping sea-floor. Where the topographic slope greatly exceeds the free surface slope, over the ridge slopes in the present model set-up $\langle \phi_B \rangle$ corresponds to the classical barotropic conversion term appearing on the rhs of Eq. (1.1). We define the GR to extend 20 km on either side of the ridge crest.

The PR involves energy transfers between the potential and kinetic compartments associated primarily with propagation of the internal and surface tide. The flat sea-floor of the PR permits propagation of the barotropic and baroclinic tidal waves without significant mutual interaction, provided the ratio of the internal and surface gravity wave celerity, respectively $c_n$ and $c_0$, satisfies $c_n^2/c_0^2 \ll 1$ (Pedlosky, 2003, p. 85), which applies in the present case study. The energy transfers in the PR are defined as the difference between those in a domain extending 40 km on either side of the ridge crest and the energy transfers in the GR.

### 4.3. Diapycnal mixing

In order to obtain a closed energy balance, the energy analysis is explicitly based on the formulation of the numerical model, in particular considering the total baroclinic potential energy.

The potential energy of a finite-width fluid domain can be defined as the sum of dynamically available ($\bar{E}_A$) and background potential energy ($\bar{E}_B$). Regarded as global quantities of a fluid volume, the former may be defined as the maximum amount of potential energy available for conversion into kinetic energy under any adiabatic redistribution of mass (Lorenz, 1955), so that simultaneously $\bar{E}_B$ be the minimum potential energy attainable by adiabatically rearranging the fluid in the volume under consideration (Winters et al., 1995). The energy of a wave field is best characterised by the pseudo-energy, the sum of kinetic energy and $\bar{E}_A$ (Shepherd, 1993). In a dynamic fluid, the energy conversion due to diffusive mixing (e.g. $F_{m,H}$, $\phi_{m,H}$) may contain both an adiabatic (reversible) part, due to diffusion along isopycnals, and an
irreversible diabatic part due to diapycnal diffusion, by which $\hat{E}_A$ is lost from the wave-field by increasing the background potential energy, $\hat{E}_B$.

To isolate the energy expended in irreversible mixing, following Winters et al. (1995), we examine the gain in background potential energy in a domain of finite width, through local diapycnal mixing, $\langle\phi_d\rangle$, and through the domain boundaries, $\langle F_d \rangle$. The instantaneous definition per unit horizontal area is

$$F_d + \phi_d = \frac{\int_0^1 \frac{\partial}{\partial x_\alpha} \left( gZ^* DK_\alpha^H \frac{\partial \tilde{\rho}}{\partial x_\alpha} \right) d\sigma + \underbrace{\left[ gZ^* \frac{K^V}{D} \frac{\partial \tilde{\rho}}{\partial \sigma} \right]}_{F_{d,H}} \underbrace{\int_0^1 gDK_\alpha^H \frac{\partial Z^*}{\partial x_\alpha} \frac{\partial \tilde{\rho}}{\partial x_\alpha} d\sigma - \int_0^1 \underbrace{gK^V \frac{\partial Z^*}{\partial \sigma} \frac{\partial \tilde{\rho}}{\partial \sigma}}_{\phi_{d,V}} d\sigma}{\phi_{d,H}} \quad (4.3)$$

where $\phi_{d,H}$ and $\phi_{d,V}$ are local conversions, and $F_{d,H}$ and $F_{d,V}$ represent the divergence of the diapycnal mixing flux in an infinitesimally wide fluid column. Note that in the present case, $F_{d,V} \equiv 0$. Here, $Z^* = Z^*(\tilde{\rho})$ can be identified as the depth in the background state of a fluid parcel of density $\tilde{\rho}$. The local $\hat{E}_A$-change due to diapycnal mixing is defined as the difference between the total potential energy gain and the $\hat{E}_B$-gain due to diffusion, denoted e.g. $-\phi_{lV} = -(\phi_{m,V} - \phi_{d,V})$.

In non-uniformly discretised coordinate systems, numerically adjacent grid-cells may have different physical volume and altitude. The calculation of background potential energy by adiabatic rearrangement of the density field should be carefully considered in order to ensure conservation of volume. This issue was not addressed in the algorithms for Boussinesq fluids proposed by Winters et al. (1995), Huang (1998), or Tseng and Ferziger (2001). To fill this gap, following similar principles we propose an algorithm adapted to $\sigma$-coordinates, guaranteeing volume conservation and applicable to free-surface flow, which is described in Appendix C.

# 5. Numerical energy diagnostics of internal tide generation over a deep-ocean ridge

## 5.1. Global energy balance

The results in this sub-section correspond to the internal tide GR, since it experiences the strongest energy transfers, showing characteristic limits of precision with respect to the smallest energy transfers.

<caption>**Table 5.1**</caption>

Magnitude in $10^6 \log(\text{power/m})$ of the LHS ($M_{lhs}$) of the energy equations (3.2), (3.3), (3.5), (3.6), and the pressure-gradient decompositions (A8), (A9). Also shown are magnitudes of the largest and smallest term in each balance ($M_{max}$, $M_{min}$), and that of the maximum absolute error ($\epsilon_{max}$), calculated over 20 tidal periods for the control simulation. KE refers to total kinetic energy, PO and PI, respectively to barotropic and baroclinic potential energy, PGO and PGI to barotropic and baroclinic pressure gradient decomposition, and $\phi_z$ refers to the buoyancy conversion term between KE and PI.

<table><thead><tr><th>Balance</th><th>Equation</th><th>M<sub>lhs</sub></th><th>M<sub>max</sub></th><th>M<sub>min</sub></th><th>ε<sub>max</sub></th></tr></thead><tbody><tr><td>KE</td><td>(3.2)<br/>+(3.3)</td><td>4</td><td>4</td><td>-3</td><td>-9</td></tr><tr><td>PO</td><td>(3.5)</td><td>4</td><td>4</td><td>4</td><td>-12</td></tr><tr><td>PGO</td><td>(3.2), (A8)</td><td>3</td><td>3</td><td>-3</td><td>-9</td></tr><tr><td>PI</td><td>(3.6)</td><td>4</td><td>4</td><td>0</td><td>-5</td></tr><tr><td>PGI</td><td>(3.3), (A8)</td><td>3</td><td>4</td><td>2</td><td>-9</td></tr><tr><td>φ<sub>z</sub></td><td>(4.2ab)</td><td>4</td><td>4</td><td>3</td><td>-5</td></tr></tbody></table>

### 5.1.1. Numerical precision

Table 5.1 shows the magnitude of the instantaneous energy transfers during 20 tidal periods, integrated over the GR, for the left-hand side of the energy evolution equations (3.2), (3.3), (3.5) and (3.6), and the pressure gradient decompositions (A7) and (A8). The magnitude of the largest and smallest term, and the maximum absolute difference betweenlhs and rhs,

$$\epsilon_{\max} = \max(|lhs - rhs|), \quad (5.1)$$

which should be zero theoretically, are shown for each balance. The instantaneous barotropic kinetic and potential energy balance, as well as the pressure gradient decompositions, are in equilibrium up to 13 significant digits. The error occurs due to numerical rounding. The baroclinic potential energy, however, has up to nine significant digits, which is also reflected in the buoyancy term. This reduced numerical precision is due to the summation of large terms which nearly cancel, associated with potential energy flux divergence in particular. The smallest term in all balances remains of order $O[10^5]$ W/m larger than $\epsilon_{max}$, supporting the utility of model-oriented energetics analysis. An overview of the instantaneous and tidally-averaged complete kinetic and potential energy balances and the decomposition into the primary terms of the tidal energy balance is shown in Table 4.2 of Floor (2009).

### 5.1.2. Coriolis energy flux

In theory, the Coriolis effect does not do work (e.g. Gill, 1982, Section 4.6) and integrated over a closed domain the sum of the discrete version of the Coriolis terms in the kinetic energy balance $-v_x f v_y$ and $v_y f v_x$ vanishes identically. However, due to the spatial discretisation of the momentum equations on the C-grid, this is not the case for a domain with open lateral boundaries, as shown in Appendix B.4.

When the internal tide field is well-established, in the GR the net energy flux due to spatial discretisation of the Coriolis force primarily acts as a small source of baroclinic kinetic energy, on the order of 1% of the fluxes in the predominant net energy balance between the barotropic pressure flux, $-\langle F_T + F_X \rangle$, the buoyancy term, $\langle \phi_z \rangle$, and the baroclinic IT-pressure flux, $-\langle F_{lW} \rangle$. The net barotropic Coriolis flux is several orders of magnitude smaller and negligible. In the control simulation, the net baroclinic Coriolis flux is an order of magnitude larger than the net potential energy change due to vertical mixing, $-\phi_{m,V}$. This emphasises the importance of accurate energy diagnostics.

### 5.1.3. Temporal diffusion

The energy transfer associated with the Robert–Asselin (R–A) filter resembles a Laplacian diffusion of energy in time (cf. Marsaleix et al., 2008; Floor, 2009, for details). Its strength is expressed by the coefficient $0 \le A < 0.5$. If $A=0$, a pure leapfrog scheme is used so that the odd and even trajectories evolve separately and may diverge. If $A=0.5$, the two trajectories are averaged each iteration. The appropriate value for A is determined case-by-case. In the control simulation $A=0.1$ is used. The filter chiefly acts during periods of model spin-up and ideally only negligibly affects the model physics. Since the R–A filter is used only to enhance model stability, only two weaker values of A are tested in a simulation otherwise corresponding to the control run. Table 5.2 shows the

<table><caption>**Table 5.2**
Diffusion tests.</caption><thead><tr><th>Parameter</th><th>Control</th><th>Tests</th></tr></thead><tbody><tr><td>A</td><td>0.1</td><td>0.05, 0.01</td></tr><tr><td>K<sup>V</sup> (m<sup>2</sup> s<sup>-1</sup>)</td><td>10<sup>-6</sup></td><td>Unforced/tidal motion, 10<sup>-5</sup>, 10<sup>-4</sup>, 10<sup>-3</sup></td></tr><tr><td>K<sup>H</sup> (m<sup>2</sup> s<sup>-1</sup>)</td><td>0</td><td>Unforced/tidal motion, 0.5, 1, 2</td></tr></tbody></table>
values of the diffusion coefficients and A used in the control simulation as well as the values of A tested while keeping the remaining parameters constant.

For all three cases, during the first tidal period, the R-A filter acts relatively strongly yet intermittently, at O[10³] W/m times A for durations of 0.001 T, due to the strong velocity shear associated with the surface wave front that initially rapidly propagates through the domain. In the stationary state after about 4 T, the R-A filter causes a kinetic energy loss of about -1 W/m times A. Its direct effect on potential energy is negligible. This is of the same order as $\phi_{m,v}$ in the control simulation. Reducing the coefficient A negligibly enhances the other (net) energy transfers. Since its purpose is to enhance model stability and both test simulations with smaller R-A coefficients are stable over 20 tidal periods, the R-A coefficient should be reduced accordingly.

## 5.2. Instantaneous and local energy transfers in GR and PR

The barotropic buoyancy term $\bar{\phi}_z$ depends on $\bar{\phi}_z$ and therefore predominantly acts over the topographic slopes in the GR, with mean 720 W/m and amplitude O[10³] W/m. In the PR $\bar{\phi}_z$ has zero mean and small amplitude. This confirms decoupling of barotropic and baroclinic wave motion, so we distinguish a fortiori the internal tide GR and the PR, where both tidal wave types interact little. Since $v'_z$ is significant throughout the domain, so is the baroclinic buoyancy term, $\phi'_z$, with an amplitude of O[10³] W/m.

For zero stratification and/or topography, the modelled time-mean of the eastward barotropic energy density flux $Fl_T = g\rho_0 D \bar{v}_x \eta$ per unit width is

$$ \langle Fl_T \rangle_T \approx 2.506 \times 10^5 \text{ W/m}, \quad (5.2) $$

which is matched in the control simulation far west of the ridge. This is slightly smaller than the theoretical time-mean energy flux $2.529 \times 10^5 \text{ W/m}$ that corresponds to the linear forcing (cf. Gill, 1982, Section 8.2). Far east of the ridge, the tidal energy flux is $\langle Fl_T \rangle_T = 2.499 \times 10^5 \text{ W/m}$, so that the tidal wave loses approximately 0.3% of its energy in traversing the ridge. Nearer the ridge, it is the sum $\langle Fl_T + Fl_X \rangle_T \approx 2.50 \times 10^5 \text{ W/m}$ where $Fl_X = \int_0^1 D\bar{v}_x P d\sigma$, that corresponds to the barotropic tidal energy flux, since the tidal wave feels the internal tide pressure perturbations.

Consistency in the energy diagnostics is indeed important, since a small error on the scale of the instantaneous tidal energy flux may lead to large uncertainties in the net energy transfers of O[10²] W/m after 10 T, which are discussed below.

## 5.3. Primary energy transfers in the generation region

Fig. 2 illustrates the primary net energy transfers involved in internal tide generation by a barotropic tide impinging on a ridge. The initial spin-up period lasts about 3.5 T. In Fig. 2c, an initial build-up of baroclinic kinetic energy is evident, associated with the establishment of the internal tide field. Fig. 2b shows that the local baroclinic potential energy tendency remains positive, reflecting the fact that baroclinic potential energy is being advected into the water column. After 4T, some 12 normal modes propagating at the group speed have reached the edge of the GR, which corresponds to the maximum wave number that can be reasonably represented at the present resolution, requiring a minimum of three grid-points/levels per wavelength. For waves at the tidal frequency, the vertical resolution is the limiting factor here.

Fig. 2a shows that the net barotropic pressure flux into the domain, $-(\langle F_T + F_X \rangle)$, of 720 W/m is balanced by the barotropic buoyancy term $\langle \bar{\phi}_z \rangle$. The former becomes stationary after approximately 4 T, whereas $\langle \bar{\phi}_z \rangle$ shows only a slight increase from the outset. The net work done by the barotropic current matches the barotropic

buoyancy term to ±1% after 4 T and less than ±0.05% after 6 T. After a spin-up period of 4 T, the running mean barotropic kinetic energy equation (4.1a) therefore reduces to

$$ \langle F_T + F_X \rangle \approx -\langle \bar{\phi}_z \rangle. \quad (5.3) $$

The barotropic buoyancy term matches the classical barotropic to baroclinic tide conversion term, i.e. $\langle \bar{\phi}_z \rangle \approx \langle \bar{\phi}_B \rangle$ with a margin of 0.01%, because the surface-contribution $\langle \phi_S \rangle$ is negligible.

Fig. 2b shows that the primary balance in the potential energy equation (4.1b) involves all terms except diffusion, i.e.

$$ \left\langle \frac{\partial \tilde{E}_p}{\partial t} \right\rangle + \left\langle F_{E_p} \right\rangle \approx \langle \bar{\phi}_z \rangle + \langle \phi'_z \rangle, \quad (5.4) $$

which is satisfied with a margin less than 0.05%. The barotropic buoyancy term $\langle \bar{\phi}_z \rangle$, and the potential energy flux into the fluid column $-(\langle F_{E_p} \rangle)$, act as sources of baroclinic potential energy and a transfer into baroclinic kinetic energy through the baroclinic buoyancy term $\langle \phi'_z \rangle$.

Until about 4 T the baroclinic current contributes about two-thirds of the potential energy advection, compared to one third by the barotropic part $\langle F_{E_p} \rangle$. After 4T, the latter is reduced to a few percent of the total baroclinic potential energy flux, while the baroclinic part $\langle F'_{E_p} \rangle$, dominates. The local baroclinic potential energy gain $\langle \partial \tilde{E}_p / \partial t \rangle$ is about 80% of the net advective flux $\langle F_{E_p} \rangle$. The remainder plus $\langle \bar{\phi}_z \rangle$ is converted into baroclinic kinetic energy by $\langle \phi'_z \rangle$, which is thus about 13% larger than the barotropic buoyancy term.

Finally, during the quasi-stationary period, in the baroclinic kinetic energy balance the internal wave pressure flux, $\langle F_{IW} \rangle \approx 823 \text{ W/m}$, closely matches the energy provided by the baroclinic buoyancy term $-\langle \phi'_z \rangle$. Thus, Fig. 2c illustrates that to within a margin of 3%

$$ \langle F_{IW} \rangle \approx -\langle \phi'_z \rangle. \quad (5.5) $$

The remaining 3% mainly consists of the local gain in baroclinic kinetic energy and the baroclinic numerical Coriolis-flux described in Section 5.1.2, which represents approximately 0.3% of the total running mean kinetic energy tendency. Effects of advection and viscosity in the kinetic energy equations are several orders of magnitude smaller. These results are neither sensitive to the initial phase of the tidal forcing, nor to the length of the integration region, although the outward internal tide energy flux reaches steady state later in a longer domain (cf. Lamb, 2007, Fig. 4).

It is remarkable that the baroclinic pressure flux, $\langle F_{IW} \rangle$, exceeds the energy lost by the barotropic tide, $-\langle F_T + F_X \rangle$, by a factor 1.14. The excess is attributed to part of the baroclinic flux of potential energy $\langle F'_{E_p} \rangle$ into the GR. This significant contribution may be due to the inclusion of the background density field in the baroclinic density $\bar{\rho}$ and pressure $P$ in our model-oriented approach.

The strength of the time-averaged baroclinic flux associated with the background stratification can be estimated by a simple scale analysis from the expression $Dv'_x \bar{\rho} gz$, averaged over time and integrated over the fluid depth.

Assuming a sinusoidal dependence of the free-surface motion and considering the first baroclinic normal mode horizontal velocity, $v'_{x,1} = -u'_1 \cos(m_1 z) \sin(k_1 x \pm \omega t)$, with amplitude $u' = O[10^{-2}] \text{ms}^{-1}$, using the physical parameters in Table 2.1a and using $\eta = \sqrt{\frac{\bar{H}}{g}} U_0 \sin(Kx - \omega t)$ for the time-varying part of the total depth D, we find

$$ \frac{1}{2T} \int_{-T}^{T} \int_{0}^{1} D v'_{x} \bar{\rho} gz d\sigma dt = - \frac{\rho_0 N^2}{2T} \sqrt{\frac{\bar{H}}{g}} U_0 \int_{0}^{1} \left\{ z^2 \sum_{n} u'_{n} \cos(m_n z) \times \int_{-T}^{T} \sin(Kx - \omega t) \sin(k_n x - \omega t) dt \right\} d\sigma = O[10] \text{ W/m,} $$
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Primary volume-integrated energy transfers in internal tide generation, integrated over a 2D column extending 20 km on either side of the ridge. The moving average over two tidal periods is shown. (a) The barotropic pressure-flux $-\langle F_T + F_X \rangle$ is balanced by the barotropic buoyancy term, $\langle \phi_z \rangle \approx \langle \phi_B \rangle$. (b) The barotropic buoyancy term $\langle \phi_z \rangle$ and the baroclinic flux of baroclinic potential energy $-\langle F_{E_p} \rangle$ cause a local gain in potential energy $\langle \partial \hat{E}_p / \partial t \rangle$, and a conversion into baroclinic kinetic energy by the baroclinic buoyancy term $-\langle \phi'_z \rangle$. The term $\langle \text{sum} \rangle$ represents the sum of the terms in Eq. (5.5) that together balance $\langle \partial \hat{E}_p / \partial t \rangle$. (c) The baroclinic buoyancy term $-\langle \phi'_z \rangle$ is radiated by the internal waves, represented by the baroclinic pressure flux, $\langle F_{IW} \rangle$.</figcaption></figure>

into the GR for eastward travelling internal waves. Both east- and westward internal waves transfer energy into the GR at order $O[10]$ W/m. This term is non-zero for any normal mode when $k_n \neq K$, which is the case here since the surface wavelength is much greater than the internal wavelength, so that $k_n \gg K$. This interaction term between free-surface and internal tide is non-existent in a rigid-lid model such as used by Lamb (2007) and particular to the free-surface formulation used here. In a future study, the internal tide signature may be separated from the background density (or potential energy) field using the pseudo-energy formulation of Winters et al. (1995), or the combined wavelet-eof analysis of Pairaud and Auclair (2005).

Munroe and Lamb (2005) numerically estimated the internal tide energy flux in nearly identical physical circumstances, but at latitude 20°N and with much stronger vertical diffusivity/viscosity, with the pressure flux formulation $F_{IW}^{ML} = \int_H^\eta p_I v_x dz$, where $p_I$ is the approximate hydrostatic pressure anomaly associated with the internal tide field, rather than the complete internal pressure $P$. They found a net eastward baroclinic flux of internal pressure of $F_{IW}^{ML} \approx 400$ W/m. Multiplying this by 2, this does correspond closely to our estimate of $\langle F_{IW} \rangle \approx 823$ W/m of the sum of east- and westward baroclinic pressure flux, despite the fact that $\langle F_{IW} \rangle$ includes both the IT-pressure and the pressure associated with the background density field.

### 5.4. Primary energy transfers in the propagation region

Fig. 3 illustrates the primary energy transfers integrated over the PR, between 20 and 40 km on either side of the ridge. Initially, as the internal tide front arrives, the internal wave pressure flux into the PR causes a local increase in kinetic energy and is partially converted into baroclinic potential energy, so that

$$ \left\langle \frac{\partial \hat{E}'_K}{\partial t} \right\rangle + \langle F_{IW} \rangle \approx -\langle \phi'_z \rangle, \qquad (5.6) $$

as shown in Fig. 3a.

After about 6 T, internal wave modes up to mode 12 have travelled over 40 km away from the ridge, and the local kinetic energy gain is reduced to zero. Therefore, the net baroclinic kinetic energy balance reduces to

$$ \langle F_{IW} \rangle \approx -\langle \phi'_z \rangle, \qquad (5.7) $$

which amounts to about 360 W/m from 6 T onwards, showing a gradual decline towards about 308 W/m after nearly 10 T.

Fig. 3b shows that until about 2 T, the baroclinic buoyancy term causes a local increase in baroclinic potential energy. Simultaneously, the advective flux of baroclinic potential energy $\langle F_{E_p} \rangle$, out of the PR, strongly increases. The approximate baroclinic potential energy balance is

$$ \left\langle \frac{\partial \hat{E}_p}{\partial t} \right\rangle + \langle F'_{E_p} \rangle \approx \langle \phi'_z \rangle, \qquad (5.8) $$

because the net barotropic flux of potential energy, $\langle \bar{F}_{E_p} \rangle$, converges to zero (not shown). After 6T, the baroclinic potential energy gain is of order $O[10]$ W/m and the sum of the balances (5.7) and (5.9) then approximately reduces to

$$ \langle F_{IW} \rangle \approx -\langle F'_{E_p} \rangle, $$

i.e. the net total energy gain in the PR is small. The strength of the constituents of balance (5.9) depends on the length of the domain of propagation. When the PR length matches the first normal mode wavelength, the net transfers in (5.9) reduce to about 10% of the baroclinic flux of baroclinic pressure from the GR, while the local baroclinic potential energy gain remains of $O[10]$ W/m. As in the GR, our results (5.7) and (5.9) for the PR may be attributed to the combined use of a free surface and the total baroclinic density, $\bar{\rho}$, that includes the background stratification in accordance with the present model-oriented approach.
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Primary volume-integrated energy transfers in internal tide propagation, between 20 and 40 km from ridge crest. (a) The internal wave pressure flux $\langle F_{IW} \rangle$ is balanced by buoyancy term $-\langle \phi_z' \rangle$; (b) The buoyancy term $-\langle \phi_z' \rangle$ is balanced by a potential energy gain $\langle \partial \tilde{E}_P / \partial t \rangle$ and radiated out of the propagation region by the baroclinic potential energy flux characterised by $\langle F_{Ep} \rangle$, while the barotropic potential energy flux is negligible. The term $\langle \text{sum} \rangle$ represents the sum of the terms in Eq. (5.9) that together balance $\langle \partial \tilde{E}_P / \partial t \rangle$.</figcaption></figure>

## 5.5. Vertical and horizontal diffusion and bottom friction

In this section, we address the following questions: how strong are the diffusive potential energy transfers? How sensitive are the primary energy balances to the parameterisation of sub-grid processes? How much energy do the dynamics lose through diapycnal mixing? What is the effect of a free-slip bottom boundary condition? The different tests for $K^H$ and $K^V$ keeping the remaining parameters constant are summarised in Table 5.2.

### 5.5.1. Diffusive energy transfers

In the absence of tidal motion both vertical and 'horizontal' (along $\sigma$-surfaces) diffusion of density associated with sub-grid scale processes may already change potential energy. Therefore we also analyse an unforced case labelled STILL. The relative effect of the tidal dynamics on the diffusive energy transfers is labelled DYN, so that e.g. $\phi_{m,V}^{DYN} = \phi_{m,V} - \phi_{m,V}^{STILL}$.

The STILL domain-averaged potential energy gain due to vertical density diffusion is for arbitrary constant $K^V$ and $N$

$$-\overline{\langle \phi_{m,V}^{STILL} \rangle} = \frac{k_m - 1}{k_m} \rho_0 K^V N^2, \quad (5.9)$$

where $k_m = 40$ is the number of vertical mass levels in the fluid interior, $\overline{\langle \cdot \rangle}$ indicates the volume average and the factor $\frac{k_m-1}{k_m} = 0.975$, stems from summation by parts in the discretised formulation to obtain the integrand of $\phi_{m,V}$, using vertically uniform spacing $\Delta\sigma$. The initial GR and PR volumes are respectively, 187.9 km³ and 200.0 km³. In the control simulation in the GR the initial diffusive potential energy gain is $-\langle \phi_{m,V}^{STILL} \rangle = 0.18848$ W/m, or about 0.02% of $-\langle F_T + F_X \rangle$.

In the STILL case, non-zero 'horizontal' diffusion along $\sigma$-surfaces causes a potential energy gain in the GR of

$$-\overline{\langle \phi_{m,H}^{STILL} \rangle} = 3.23 \times 10^{-3} \rho_0 K^H N^2 \text{ Wm}^{-3},$$

which is non-zero because of the sloping sea floor. In the GR, e.g. with $K^H = 1 \text{ m}^2 \text{ s}^{-1}$ this gives $-\langle \phi_{m,H}^{STILL} \rangle = 6.24 \times 10^2 \text{ W/m}$ or a significant fraction (86%) of $-\langle F_T + F_X \rangle$. In the PR, the initial horizontal diffusion is zero, because isopycnal and $\sigma$-surfaces coincide. After 10 T, the local impact of $K_H$ on the potential energy tendency in the PR is two orders of magnitude smaller than in the GR.

For several values of $K^V$ and $K^H$, Table 5.3 shows the relative potential energy change associated with tidal dynamics due to vertical and horizontal density diffusion in the GR after 10 T, normalised

<caption>**Table 5.3**</caption>

Potential energy gain relative to the STILL case, for various strengths of vertical and horizontal diffusion.

<table><thead><tr><th colspan="3">(a) Instantaneous potential energy gain in the generation region due to vertical diffusion at t = 10T, relative to STILL case. For K<sup>H</sup> ≠ 0, the vertical diffusion equals K<sup>V</sup> = 10<sup>-6</sup> m<sup>2</sup> s<sup>-1</sup></th></tr><tr><th>K<sup>V</sup> (m<sup>2</sup> s<sup>-1</sup>)</th><th>&lang;&phi;<sub>m,V</sub><sup>DYN</sup>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;</th><th>-&lang;&phi;<sub>I,V</sub>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;</th></tr></thead><tbody><tr><td>10<sup>-6</sup></td><td>3.32 &times; 10<sup>-3</sup></td><td>-1.40 &times; 10<sup>-2</sup></td></tr><tr><td>10<sup>-5</sup></td><td>3.33 &times; 10<sup>-3</sup></td><td>-1.38 &times; 10<sup>-2</sup></td></tr><tr><td>10<sup>-4</sup></td><td>3.18 &times; 10<sup>-3</sup></td><td>-1.36 &times; 10<sup>-2</sup></td></tr><tr><td>10<sup>-3</sup></td><td>1.71 &times; 10<sup>-3</sup></td><td>-1.16 &times; 10<sup>-2</sup></td></tr><tr><th>K<sup>H</sup> (m<sup>2</sup> s<sup>-1</sup>)</th><th>&lang;&phi;<sub>m,V</sub><sup>DYN</sup>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;</th><th>-&lang;&phi;<sub>I,V</sub>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;/&lang;&phi;<sub>m,V</sub><sup>STILL</sup>&rang;</th></tr><tr><td>0.5</td><td>3.08 &times; 10<sup>-3</sup></td><td>-1.22 &times; 10<sup>-2</sup></td></tr><tr><td>1</td><td>2.79 &times; 10<sup>-3</sup></td><td>-1.11 &times; 10<sup>-2</sup></td></tr><tr><td>2</td><td>2.20 &times; 10<sup>-3</sup></td><td>-9.45 &times; 10<sup>-3</sup></td></tr></tbody><thead><tr><th colspan="3">(b) Instantaneous potential energy gain in the generation region due to diffusion along &sigma;-surfaces at t = 10T, relative to STILL case. The vertical diffusion equals K<sup>V</sup> = 10<sup>-6</sup> m<sup>2</sup> s<sup>-1</sup></th></tr><tr><th>K<sup>H</sup> (m<sup>2</sup> s<sup>-1</sup>)</th><th>&lang;&phi;<sub>m,H</sub><sup>DYN</sup>&rang;/&lang;&phi;<sub>m,H</sub><sup>STILL</sup>&rang;/&lang;&phi;<sub>m,H</sub><sup>STILL</sup>&rang;</th><th>-&lang;&phi;<sub>I,H</sub>&rang;/&lang;&phi;<sub>m,H</sub><sup>STILL</sup>&rang;/&lang;&phi;<sub>m,H</sub><sup>STILL</sup>&rang;</th></tr></thead><tbody><tr><td>0.5</td><td>-1.63 &times; 10<sup>-2</sup></td><td>3.22 &times; 10<sup>-4</sup></td></tr><tr><td>1</td><td>-1.94 &times; 10<sup>-2</sup></td><td>5.75 &times; 10<sup>-3</sup></td></tr><tr><td>2</td><td>-2.53 &times; 10<sup>-2</sup></td><td>1.38 &times; 10<sup>-2</sup></td></tr></tbody></table>

respectively, by $-\langle \phi_{m,V}^{STILL} \rangle_X$ and $-\langle \phi_{m,H}^{STILL} \rangle_X$, which increases nearly linearly from zero.

Fig. 4a shows the time-integral over $5T \le t \le 10T$ of the integrand of $-\phi_{m,V}$ in the control simulation. Along the internal tide beams, the local $\tilde{E}_P$-change due to vertical diffusion is significantly modulated, by up to ±5%. This applies to all tested $K^V$. However, as shown in Table 5.3, after 10 T integrated over the GR the dynamics cause an increase of $-\langle \phi_{m,V} \rangle_X$ of less than 0.5% compared to $-\langle \phi_{m,V}^{STILL} \rangle_X$, while $-\langle \phi_{m,H} \rangle_X$ decreases by about 1.5–2.5% compared to $-\langle \phi_{m,H}^{STILL} \rangle_X$.

### 5.5.2. Sensitivity of primary energy transfers to diffusivity

Enhanced $K^V$ increases the local potential energy gain by 5% for $K^V = 10^{-3} \text{ m}^2 \text{ s}^{-1}$, but otherwise negligibly affects the primary energy transfers. Fig. 5 shows the effect of enhanced horizontal diffusion on the primary energy transfers in the GR. The tested 'horizontal' diffusion coefficients are weak compared to commonly used values for a horizontal resolution of 1 km, such as $K^H = 10 \text{ m}^2 \text{ s}^{-1}$. However, the primary baroclinic energy transfers in the GR are strongly reduced: the baroclinic buoyancy term
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Vertical section of local tendency of the potential energy density due to vertical diffusion with $K^V = 10^{-6}$ m$^2$ s$^{-1}$, integrated from $t = 5$ T to $t = 10$ T and zoomed into ±300 km from the ridge crest. The colour bar is cut off at 1.98 and 2.48 J/kg. (a) Change in total potential energy density due to diffusion, integrand of $-\phi_{m,V}$. (b) Change in background potential energy density due to diffusion, integrand of $-\phi_{d,V}$.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Primary volume-integrated energy transfers in internal tide generation, with horizontal diffusion, super-imposed for $K^V = 10^{-6}$ with $K^H = 0, 0.5, 1$, and $2$ m$^2$ s$^{-1}$. (a) Barotropic energy transfers are hardly affected; (b) $\langle \partial \tilde{E}_p / \partial t \rangle$ increases, while $\langle F_{Ep} \rangle$ and $\langle \phi'_z \rangle$ decrease with $K^H$. The term $\langle \text{sum} \rangle$ represents the sum of the terms in Eq. (5.5) that together balance $\langle \partial \tilde{E}_p / \partial t \rangle$ and (c) both $\langle F_{IW} \rangle$ and $\langle \phi'_z \rangle$ strongly decrease with $K^H$.</figcaption></figure>

$\langle \bar{\phi}_z \rangle$ and the baroclinic flux of internal pressure $\langle F_{IW} \rangle$ decrease by a factor 2 for $K^H = 2$ m$^2$ s$^{-1}$.

### 5.5.3. Diapycnal mixing

In the STILL case, all vertical diffusion is diapycnal because the isopycnals are horizontal, giving $\phi_{m,V}^{\text{STILL}} = \phi_{d,V}^{\text{STILL}}$. Mixing along $\sigma$-surfaces initially acts as a sink of available potential energy in the GR, given by

$$-\langle \phi_{1,H}^{\text{STILL}} \rangle = -0.69 \times 10^{-5} \rho_0 K^H N^2 \text{ Wm}^{-3}.$$

Table 5.3 shows the relative available potential energy change associated with tidal dynamics due to vertical and horizontal density diffusion in the GR after $10$ T, normalised respectively by

$-\langle \phi_{m,V}^{\text{STILL}} \rangle$ and $-\langle \phi_{m,H}^{\text{STILL}} \rangle$. For vertical diffusion, the dynamic background potential energy gain exceeds the total potential energy gain, so that diapycnal diffusion acts as a sink of available potential energy. Fig. 4b shows the time-integral over $5T \le t \le 10T$ of the integrand of $-\phi_{d,V}$ in the control simulation, but a similar picture applies to all tested $K^V$. The integrand of $-\phi_{d,V}$ is more strongly modulated than that of $-\phi_{m,V}$, but remains close to $-\langle \phi_{m,V}^{\text{STILL}} \rangle$ when averaged over the GR. After about $t = 3T$, the available potential energy gain associated with 'horizontal diffusion', $-\langle \phi_{I,H} \rangle$, starts increasing linearly from zero.

### 5.5.4. Free slip sea-floor friction

In the control simulation, the bottom boundary condition is no-slip, with zero tracer flux. This unrealistically neglects bottom
boundary layer turbulence and uses weak vertical viscosity $K^V = 10^{-6} m^2 s^{-1}$. Therefore, we consider a free-slip condition, where at the mid-level of the bottom grid-cell the bottom stress is parameterised as $\tau_x = \rho_0 C_D \sqrt{\bar{v}_\beta \bar{v}_\beta} \nu_x$, with drag coefficient $C_D = 2.5 \times 10^{-3}$.

In the free-slip case, the barotropic tidal flux is enhanced by about 1% compared to the control simulation, to $\langle F_{IT} + F_{IX} \rangle_T \approx 2.535 \times 10^5 W/m$. In the stationary state after about 4T the free-slip (FS) bottom boundary condition causes less than 1% change in the amplitudes of the instantaneous energy transfers. In both the GR and PR, the energy expended by bottom friction in the no-slip case is $O[10^{-3}] W/m$, whereas in the FS case it is about 21 W/m. A simple scale analysis gives the order of magnitude of no-slip kinematic bottom stress, $K^V U_0 / \Delta z^2 = O[10^{-12}] m s^{-2}$. In the FS case, the kinematic bottom stress is of order $C_D U_0^2 = O[10^{-7}] m s^{-2}$. This explains the large difference in frictional energy dissipation.

Let us focus on the GR. The wave energy deposited there, lost from the barotropic tidal pressure flux, $\langle F_{IT} + F_{IX} \rangle_T$, is altered both by the changed bottom boundary condition 'upstream' and by local energy transfers. Compared to the no-slip case, the net barotropic tidal pressure flux into the GR $-\langle F_T + F_X \rangle$ gained 3% due to weakening of the net surface pressure flux, $-\langle F_T \rangle$, while the conversion into baroclinic potential energy $\langle \phi_z \rangle$, gained 2%. The difference between $-\langle F_T + F_X \rangle$ and $\langle \phi_z \rangle$ is attributed to enhanced $\tilde{E}_K$-dissipation by local bottom friction. Both the barotropic, $\langle \tilde{E}_{E_p} \rangle$, and baroclinic flux of potential energy $\langle F'_{E_p} \rangle$, increased in strength by about 70 W m$^3/kg$, which is a respective amplification of 100% and 12%, reflecting that $\langle \tilde{E}_{E_p} \rangle$ remains small in comparison. The internal potential energy tendency $\langle \partial \tilde{E}_p / \partial t \rangle$ increased by about 20%. Being of the same order of magnitude, the baroclinic buoyancy term $\langle \phi'_z \rangle$, increased by only 7%. The baroclinic flux of internal pressure, $\langle F_{IW} \rangle$, of the same order as the baroclinic buoyancy term $\langle \phi'_z \rangle$, increased by about 5%, while the baroclinic kinetic energy tendency $\langle \partial \tilde{E}_K / \partial t \rangle$, diminished even further compared to the control simulation. The difference $\langle \phi'_z \rangle - \langle F_{IW} \rangle$ is matched by the energy lost in bottom friction.

Increasing the vertical diffusion coefficient to a more realistic $K^V = 10^{-4} m^2 s^{-1}$ in the FS case changes the previous primary GR energy transfers by less than 1%.

# 6. Discussion and conclusion

We analysed the energetics of internal tide generation by stratified barotropic tidal flow impinging on a supercritical Gaussian ridge and its subsequent propagation, using the energy-conserving, Boussinesq, hydrostatic, free-surface ocean model Symphonie. The originality of our approach lies in the closed, i.e. numerically consistent, global energy balance, evaluated explicitly in terms of all physical and numerical energy transfers in the model. The complete, high-precision energy diagnostics allowed the quantification of the energy transfers associated with physical diffusion and numerical choices, such as the discretisation of the Coriolis acceleration.

## 6.1. Basic energy balance of internal tide generation

After a spin-up of four tidal periods, within 50 km from the ridge crest the simulated internal tide beams correspond well to linear hydrostatic plane internal waves on the f-plane, at the tidal frequency. The generation region (GR) extending 20 km on either side of the ridge crest is characterised by strong topographic slopes over which the barotropic tidal flow generates strong $\bar{v}_z$ and the barotropic buoyancy term associated with the topographic slope, $\bar{\phi}_b$, identified as the primary barotropic-to-baroclinic energy conversion term.

The stationary primary net energy transfers associated with internal tide generation are summarised in Fig. 6a and to a good approximation given by

$$ \langle F_T + F_X \rangle \approx -\langle \bar{\phi}_z \rangle, \qquad (6.1a) $$

$$ \langle \frac{\partial \tilde{E}_p}{\partial t} \rangle + \langle F'_{E_p} \rangle \approx \langle \phi_z \rangle - \langle F_{IW} \rangle, \qquad (6.1b) $$

where the outward baroclinic pressure flux $-\langle F_{IW} \rangle$ was substituted for the baroclinic buoyancy term $\langle \phi'_z \rangle$ in (6.1b). About 720 W/m of barotropic energy is lost by the surface tide, $-(F_T + F_X)$, in the GR, which is converted into internal tide energy through the barotropic buoyancy term $\langle \bar{\phi}_z \rangle$. This corresponds to about 0.3% of the total surface tide energy flux $\langle F_{IT} + F_{IX} \rangle_T \approx 2.506 \times 10^5 W/m$. The internal tide energy flux out of the GR, $\langle F_{IW} \rangle$, amounts to 823 W/m and exceeds $\langle \phi_z \rangle$ by 13%. Thus, $\langle F_{IW} \rangle$ is fed by both $\langle \phi_z \rangle$ and the difference between the local baroclinic potential energy gain and the net baroclinic potential energy flux into the GR $-\langle F'_{E_p} \rangle$. The significance of $\langle F'_{E_p} \rangle$ is most likely associated with the density decomposition (2.2) that was used for consistency with the numerical formulation: the baroclinic density $\bar{\rho}$ includes the background stratification (cf. Section 5.3). Furthermore, $\langle F'_{E_p} \rangle$ and $\langle F_{IW} \rangle$ may include effects of trapped waves and tidal rectification (e.g. Pérenne et al., 2000), besides the presence of (sub-) harmonics and non-linearity.

These results are not very sensitive to changes in vertical diffusion, which primarily affects the background density field. The second objective of this work was to quantify the energy transfer due to diapycnal mixing. In all tested cases, diapycnal diffusion acts as an APE-sink around 1% of the STILL potential energy gain. However, to determine the deep-ocean mixing energy budget an appropriate horizontal diffusion formulation is crucial, because it strongly affects the baroclinic energy conversions in the GR: the internal wave pressure flux to the abyssal ocean, $\langle F_{IW} \rangle$, is reduced by 40% for weak $K_H = 2 m^2 s^{-1}$. Adopting a free-slip bottom boundary formulation causes significant local energy dissipation due to seafloor friction since it parameterises bottom boundary layer turbulence. However, the associated change in net flux of energy and pressure into and out of the GR is only a few percent of the energy lost in the GR by the barotropic tide, $-(F_T + F_X)$.

Kinetic energy advective terms play a negligible role in the primary energy balance (6.1) in the GR of the present simulation, whereas $-\langle F'_{E_p} \rangle$ is important due to interaction between the surface and internal tide and the background stratification. Non-linearity can transfer energy between the surface and internal tide, to (sub-) harmonics or rectified tides (Pérenne et al., 2000). As noted by Carter et al. (2008) it is possible that non-linear interactions are underestimated here, because the simulation was started from rest.

## 6.2. Basic energy balance of internal tide propagation

In the propagation region (PR), extending between 20 and 40 km from the ridge crest, the surface and internal waves propagate rather independently because the barotropic vertical velocity, $\bar{v}_z$, is negligible. Vertical motion is primarily due to the internal tide and the buoyancy term is dominated by $\phi'_z$.

In the PR, when after about 6 T a stationary state is reached, the primary net energy balance is

$$ -\langle F_{IW} \rangle \approx \langle \phi'_z \rangle \approx \langle \partial \tilde{E}_p / \partial t \rangle + \langle F'_{E_p} \rangle, \qquad (6.2) $$

which is summarised in Fig. 6b, amounting to about 308 W/m in the control simulation. The net baroclinic pressure flux, $-\langle F_{IW} \rangle$, into the PR is approximately balanced by the baroclinic flux of potential energy out of the PR, $\langle F'_{E_p} \rangle$, resulting in a small local baroclinic
<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Conceptual scheme of net internal tide energy transfers and balance at arbitrary time $t = 8.925 T$, accurate to <1%. (a) Generation region: the mean eastward barotropic tidal energy flux $\langle F_{IX} + F_{IT} \rangle_T$ is shown west and east of the ridge. (b) Propagation region.</figcaption></figure>

potential energy gain. The internal and surface tide are decoupled, while the barotropic energy flux divergence is negligible in comparison with the terms in (6.2). The values of $-\langle F_{IW} \rangle$, $\langle \phi'_z \rangle$ and $\langle F'_{E_p} \rangle$ depend on the length of the PR and reduce to about 10% of the baroclinic flux of internal pressure from the GR if it equals integer wavelengths of the first (and all other) internal wave normal modes at the tidal frequency.

As mentioned in the introduction, small-amplitude internal tides are often described using a rigid lid and density anomaly with respect to the background stratification (e.g. Lamb (2007), illustrative example 2). In the linear description of propagating small-amplitude sinusoidal internal wave modes, the net tidally averaged IT-pressure flux and baroclinic potential energy flux into the PR are zero. As in the GR, our results (6.2) for the PR may be attributed to the combined use of a free surface and the total baroclinic density, $\bar{\rho}$, that includes the background stratification in accordance with the present model-oriented approach. Additional effects that may play a role include non-linearity, the presence of (sub-) harmonics and in a realistic setting would incorporate the geostrophic

circulation and wind induced currents etc., which may be separated using wavelet-eof analysis (Pairaud and Auclair, 2005).

Except for the extreme case with $K^V = 10^{-3} \text{ m}^2 \text{ s}^{-1}$, in the PR vertical diffusion reduces the energy transfers in (6.2) by less than 0.5%. The potential energy tendency in the PR due to 'horizontal' diffusion remains negligible.

### 6.3. Numerical evaluation of energy transfers

The energy diagnostics are explicitly based on the numerical scheme of the *Symphonie* model equations for two main reasons. Firstly, this permits closing the global energy balance to high precision, in the sense that all physical and numerical sources and sinks in the individual energy compartments of $\bar{E}_K$, $\bar{E}_P$, $\bar{E}_P$ and $E'_K$ and their mutual exchanges are quantified. This allowed a detailed analysis of energy transfers associated with the 'physical' diffusion coefficient $K^V$ and the bottom boundary condition (Section 5.5). We also found that the staggered grid causes a net baroclinic energy boundary flux due to the Coriolis effect, on the order of 1% of the
primary energy fluxes. Secondly, consistency in the energy diagnostics is needed because a small error on the scale of the instantaneous local $E_p$-flux would lead to significant uncertainties in the net energy transfers and obscure the diffusive energy transfers, $\phi_{m,v}$. Such uncertainties may have arisen in previous studies, which singled out specific (approximate) energy transfers (Munroe and Lamb, 2005), or discretised an approximate balance directly onto the model grid (Carter et al., 2008).

Katsumata (2006) and Carter et al. (2008) partitioned the energy equations into barotropic and baroclinic components, without distinguishing between kinetic and potential energy in the results. Carter et al. (2008) found a difference of 5% between the barotropic and baroclinic estimate of the buoyancy conversion term. In fact, making the distinction between kinetic and potential energy allowed us to show that the buoyancy term is conservative in Symphonie, i.e. no energy is lost in the instantaneous exchange between potential and kinetic energy that is fundamental to gravity wave motion. Additionally, we confirmed explicitly that the barotropic tidal energy lost in the generation region, $-\langle F_T + F_X \rangle$, is indeed primarily converted into baroclinic potential energy, which was previously assumed implicitly.

Comparing simulations at resolutions of $\Delta x = 4$ km and $\Delta x = 1$ km, Carter et al. (2008) concluded tidal conversion may be underestimated if the horizontal resolution is too low. In a higher resolution simulation the baroclinic flux of internal pressure $\langle F_{IW} \rangle$ might be larger due to an improved representation of previously unresolved steep sea-floor features, but also because a larger number of internal wave modes can be represented. In the present case, the vertical resolution is the limiting factor for the number of modes that can be represented. The ridge is a smooth Gaussian and its features are well-resolved with a resolution $\Delta x = 1$ km. We performed additional tests using 80 and 160 vertical levels as well as $\Delta x = 0.5$ km (not shown), but did not find an appreciable change of the primary energy transfers (6.1) and (6.2).

Pressure gradient truncation errors (PGE) associated with the terrain-following coordinates and the initialisation of the density field can cause spurious circulation and energy transfers (Marsalleix et al., 2009b, Section 2.2). However, analysis of the STILL case indicates that the weak PGE-induced circulation is unlikely to have affected the tidal dynamics and energy transfers significantly. Separating the physical energy transfers from the PGE-induced transfers is beyond the scope of the present study, but may prove useful for taller and/or steeper topography, realistic non-linear stratification and/or vertical discretisation.

In conclusion, we have analysed the energetics of internal tide generation and propagation in detail, based explicitly on the discretised governing equations of Symphonie. The model formulation was deliberately straightforward, both numerically and physically to enforce consistency in the buoyancy term, $\phi_z$, and limiting vertical diffusion to a minimum in the control simulation to focus on mechanical energy conservation. We presented a schematic overview of the primary energy balances in the generation region (over sloping topography) and the propagation region, Fig. 6. The complete energy balance permitted the quantification of physical and numerical energy transfers and the analysis of the impact of diffusion and viscosity. Finally, using a novel adiabatic rearrangement algorithm adapted to free-surface flow and terrain-following coordinates, we showed that in the generation region the tidal dynamics cause an additional background potential energy gain equal to about 1% of the total diffusive potential energy gain.

Progress in distinguishing the surface and internal tide energy fluxes can, especially in the baroclinic potential energy balance, be made by decomposing density into constant, vertically varying 'background' and internal wave anomaly components. Numerous free-surface studies have used such an approach (e.g. Munroe and Lamb, 2005; Carter et al., 2008). Although the background den-

sity is straightforwardly defined in rigid-lid models (e.g. Lamb, 2007), in sigma-coordinates with a free surface that is non-trivial and will be the subject of future work.

The present approach is neither limited to the present model set-up, nor to the particular discretised formulation. A generalisation to other models, using different coordinates (Cartesian, spherical) and/or grids (A, B, unstructured), of the present energy diagnostics approach is envisageable but would require careful inspection of the particular discretised formulation. However, the present approach becomes complex in case of e.g. multi-step integration methods, or non-linear diffusion schemes.

## Acknowledgements

This project benefited from LEFE-IDAO (2006) project grant "Ondes Internes" and the ANR "PIWO" grant. The first author was supported by the French MENESR scholarship 21859-2006. The authors thank Cyril Nguyen and the Laboratoire d'Aérologie IT-team, Serge Prieur, Laurent Cabanas, Jérémy Leclercq, Didier Gagen, and Juan Escobar for their support.

## Appendix A. Model formulation

### A.1. Model equations

The governing equations under the Boussinesq and hydrostatic approximation in terrain-following coordinates are summarised as follows. The horizontal momentum equations are given by

$$ 
\begin{aligned}
\frac{\partial D v_x}{\partial t} + \frac{\partial D v_\beta v_x}{\partial x_\beta} + \frac{\partial v_\sigma v_x}{\partial \sigma} + \varepsilon_{\alpha\beta\gamma} D \Omega_\beta v_\gamma &= -Dg \frac{\partial \eta}{\partial x_\alpha} - \frac{D}{\rho_0} \frac{\partial P}{\partial x_\alpha} + \frac{1}{\rho_0} \frac{\partial P}{\partial \sigma} \frac{\partial z}{\partial x_\alpha} + \frac{\partial}{\partial x_\beta} \left( D K_\beta^H \frac{\partial v_\alpha}{\partial x_\beta} \right) + \frac{\partial}{\partial \sigma} \frac{K^V}{D} \frac{\partial v_\alpha}{\partial \sigma}.
\end{aligned}
\quad (A1) $$

The horizontal velocity is $v_\alpha$ where $\alpha = x, y$ indicates direction and $v_\sigma$ is the velocity perpendicular to $\sigma$-surfaces, which is classically called $\omega$ (Blumberg and Mellor, 1987). Repeated dummy indices $\alpha, \beta$ indicate summation, while $\varepsilon_{\alpha\beta\gamma}$ selects the vertical component of the Coriolis effect. $\Omega_\beta$ is the $\beta$-component of the Coriolis acceleration, where $\beta = x, y, z$. Furthermore, the total depth is $D = H - h(x) + \eta$, where $h$ is the sea floor elevation with respect to ambient depth $H$ and $\eta$ is surface elevation, $g$ is the gravitational acceleration, $\rho_0$ is a constant reference density, $P$ is pressure due to the density anomaly $\tilde{\rho}$, $z$ is depth, $K_\beta^H$ and $K^V$, respectively are the 'horizontal' (along $\sigma$-surfaces) and vertical coefficients for both viscosity and diffusion, which implies that the turbulent Prandtl number equals 1.

In the vertical, the hydrostatic equation holds, so that

$$ 0 = -\frac{\partial P}{\partial \sigma} - D \tilde{\rho} g. \quad (A2) $$

Furthermore, the continuity equation is

$$ \frac{\partial \eta}{\partial t} + \frac{\partial D v_\alpha}{\partial x_\alpha} + \frac{\partial v_\sigma}{\partial \sigma} = 0 \quad (A3) $$

and the heat and salinity equation are respectively given by

$$ \frac{\partial D \tilde{T}}{\partial t} + \frac{\partial D v_\alpha \tilde{T}}{\partial x_\alpha} + \frac{\partial v_\sigma \tilde{T}}{\partial \sigma} = \frac{\partial}{\partial x_\alpha} D K_\alpha^H \frac{\partial \tilde{T}}{\partial x_\alpha} + \frac{\partial}{\partial \sigma} \left( \frac{K^V}{D} \frac{\partial \tilde{T}}{\partial \sigma} \right), \quad (A4a) $$

$$ \frac{\partial D \tilde{S}}{\partial t} + \frac{\partial D v_\alpha \tilde{S}}{\partial x_\alpha} + \frac{\partial v_\sigma \tilde{S}}{\partial \sigma} = \frac{\partial}{\partial x_\alpha} D K_\alpha^H \frac{\partial \tilde{S}}{\partial x_\alpha} + \frac{\partial}{\partial \sigma} \left( \frac{K^V}{D} \frac{\partial \tilde{S}}{\partial \sigma} \right). \quad (A4b) $$

A linear equation of state is used, so that
$$ \tilde{\rho} = -\alpha_T \tilde{T} + \beta_S \tilde{S}, \quad (A5) $$

where $\alpha_T$ is the thermal expansion coefficient, while $\beta_S$ is the saline contraction coefficient.

Partial derivatives in Cartesian coordinates (indicated by an asterisk, except for the vertical) are related to those in terrain-following coordinates by

$$ \frac{\partial}{\partial t^*} = \frac{\partial}{\partial t} - \frac{1}{D} \frac{\partial z}{\partial t} \frac{\partial}{\partial \sigma}; \quad \frac{\partial}{\partial x_*} = \frac{\partial}{\partial x_*} - \frac{1}{D} \frac{\partial z}{\partial x_*} \frac{\partial}{\partial \sigma}; \quad \frac{\partial}{\partial z} = \frac{1}{D} \frac{\partial}{\partial \sigma}. \quad (A6) $$

Using these relations, it is straightforward to recover the governing equations in Cartesian coordinates, except for the viscous and diffusive terms in (A1) and (A4a), (A4b), which correspond to the formulation of Mellor and Blumberg (1985).

## A.2. Boundary conditions

Cyclic boundary conditions are employed in the meridional direction, that is all northward (scalar, volume, energy-) fluxes at the northern model boundary are recycled at the southern boundary, while respecting the C-grid discretisation. The wave propagation is in the zonal direction. Anomalies with respect to the forcing are defined as $\bar{v}_x^A = \bar{v}_x - \bar{v}_x^F$. In case of linear dynamics and monochromatic waves, the total barotropic flow at the western boundary $\bar{v}_x = \bar{v}_x^F + \bar{v}_x^A$ satisfies

$$ \frac{\partial \bar{v}_x}{\partial t} = |c_{p,0}| \frac{\partial \bar{v}_x - \bar{v}_x^F}{\partial x} + \frac{\partial \bar{v}_x^F}{\partial t}, $$

while at the eastern boundary

$$ \frac{\partial \bar{v}_x}{\partial t} = -|c_{p,0}| \frac{\partial \bar{v}_x}{\partial x}, $$

where $|c_{p,0}| = \sqrt{gH/(1 - \mu)}$ is the surface wave phase speed in the rotating regime. The remaining model variables satisfy unforced boundary conditions, e.g. at the western boundary

$$ \frac{\partial \phi}{\partial t} = |c_{p,0}| \frac{\partial \phi}{\partial x}. $$

In this setup, partial surface wave reflections occur at the eastern and western open boundary due to non-linearity. This leads to a weak initial 'sloshing' relative to the tidal wave that would occur in the absence of boundaries, associated with the propagation of the initial surface wave front through the domain. Nonetheless, for a flat sea floor, after 3 T throughout the domain the tidal wave phase is indistinguishable from that of the linear, inviscid prediction $\eta(x, t) = \eta^F(t - x/c_{p,0})$ with phase velocity $c_{p,0} = \sqrt{gH/(1 - \mu)}$ and matches its amplitude to within 1%.

## A.3. Decomposition of the pressure work rate

For reference, the energy transfers defined in Section 3 are listed in Table A.1. The barotropic (baroclinic) buoyancy term is derived from the vertically integrated work done by the barotropic (baroclinic) current against the internal pressure gradient force, using the chain-rule of partial differentiation, the continuity equation (A3) and hydrostatic balance (A2). The barotropic and baroclinic pressure work rate decomposition can thus be expressed as

$$ -D \bar{v}_x \int_0^1 \left( \frac{\partial P}{\partial x_*} + \tilde{\rho} g \frac{\partial z}{\partial x_*} \right) d\sigma = -F_x - \bar{\phi}_z, \quad (A7) $$

$$ -D \int_0^1 v_x^* \left( \frac{\partial P}{\partial x_*} + \tilde{\rho} g \frac{\partial z}{\partial x_*} \right) d\sigma = -F_{IW} - \phi_z^*. \quad (A8) $$

The magnitude and accuracy of the equivalent discretised decompositions PGO and PGI integrated over the GR are illustrated in Table 5.1 and by Floor (2009, Table 4.2).

## A.4. Energy equations in Oxyz-coordinates

Neglecting viscosity and diffusivity, the kinetic and potential energy density evolution equations are given by Gill (1982, Section 6.7)

$$ \frac{\partial E_K}{\partial t^*} + \frac{\partial}{\partial x_*} (v_x(E_K + g\rho_0\eta + P)) + \frac{\partial}{\partial z} (v_z(E_K + P)) = -\tilde{\rho}g v_z - g\rho_0\eta \frac{\partial v_x}{\partial x_*^*}, \quad (A9) $$

$$ \frac{\partial E_p}{\partial t^*} + \frac{\partial}{\partial x_*^*} (v_z E_p) + \frac{\partial}{\partial z} (v_z \tilde{E}_p) = \tilde{\rho} g v_z + g \rho_0 \eta \frac{\partial v_x}{\partial x_*^*}, \quad (A10) $$

where $E_K = \frac{1}{2}\rho_0 v_x v_x$ and $E_p = \tilde{E}_p + \bar{E}_p$. Integrating the kinetic energy density over the fluid depth, from $z = -H$ to $z = \eta$ and taking the partial time-derivative, we find the tendency of kinetic energy per unit area in the Oxy-plane of a fluid column, namely:

$$ \frac{\partial}{\partial t^*} \int_{-H}^{H} E_K dz = \int_{-H}^{H} \frac{\partial E_K}{\partial t^*} dz + E_K(z = \eta) \frac{\partial \eta}{\partial t^*} + E_K(z = -H) \frac{\partial H}{\partial t^*}, $$

where Leibniz' rule was invoked to obtain the rhs. Thus, the energy tendency per unit area of a fluid column depends on the vertical integral of the local energy equations (A9) and (A10), and two terms accounting for the variations of free surface and sea floor. Similar terms arise from the advective terms in (A9) and a similar argument applies to potential energy. These boundary terms cancel partially with terms in the vertical integral of $\partial(\cdots)\partial z$ in (A9) and partially when the free-surface pressure is zero and the sea-floor height is time-independent, impenetrable and no-slip. The result is shown in the energy equations per unit horizontal area of a fluid column, (3.10a) and (3.10b) in Section 3.

# Appendix B. Model-oriented energy analysis

## B.1. Local transfers and boundary fluxes of energy

The discretisation of the global energy evolution equations imposes the discretised formulation of the boundary fluxes of both kinetic ($F_X$, $F_T$, $F_{IW}$ and $F_{E_K}$) and potential energy ($F_{E_p}$, $F_{m,H}$) required to obtain a closed energy balance. Let us neglect variations along $Oy$ for simplicity. In general, kinetic energy conversion associated with a spatial gradient along the $Ox$-axis of arbitrary variable $A^i$, representing pressure, momentum or diffusive fluxes, integrated over a box running over $(i1 + 1)\Delta x \le x \le i2\Delta x$ can be expressed on the domain $(i1, i2) \in N$ in the form

$$ -\sum_{i=i1+1}^{i2} v_x^{i-\frac{1}{2}} (A^i - A^{i-1}) = \sum_{i=i1}^{i2-1} A^i (t_x^{i-\frac{1}{2}} - v_x^{i-\frac{1}{2}}) - v_x^{i2-\frac{1}{2}} A^{i2} \\ + v_x^{i+\frac{1}{2}} A^{i1}, \quad (B1) $$

where the sum on the rhs represents local energy conversion. The last two terms are respectively the $A^i$-flux at the eastern and western boundary, where each involves the velocity on the interior grid-point of each box boundary. For expressions involving pressure, surface elevation, longitudinal momentum flux or diffusion, one substitutes for $A^i$, respectively $P^i$, $\eta^i$, $F_{m,L}^i$, or finally the diffusive momentum flux $F_{d,L}^i$ (cf. B.2). In case of $P^i$ the discretised continuity equation allows the rhs-sum to be written in terms of the buoyancy term.

Similarly, potential energy conversions associated with divergence of scalar-flux or -diffusive flux $B^i$ (see B.3) integrated over the mass-points $(i1)\Delta x \le x \le (i2)\Delta x$, can be expressed in the form

$$ -\sum_{i=i1}^{i2} z^i \left( B_{i\frac{1}{2}}^{i+1} - B_{i\frac{1}{2}}^{i-1} \right) = \sum_{i=i1}^{i2-1} B_{i\frac{1}{2}}^{i+1} (z^{i+1} - z^i) - B_{i\frac{1}{2}}^{i2+\frac{1}{2}} z^{i2} + B_{i\frac{1}{2}}^{i1-\frac{1}{2}} z^i, \quad (B2) $$

where the boundary fluxes are defined using $B$ exterior to the boundary points $x = i1\Delta x$ and $x = i2\Delta x$. This is required for
<table><caption>**Table A.1**
List of energy transfers.</caption><thead><tr><th>Symbol</th><th>Expression</th><th>Meaning</th></tr></thead><tbody><tr><td>E<sub>K</sub></td><td><sup>1</sup>&frasl;<sub>2</sub> ρ<sub>0</sub> v<sub>x</sub> v<sub>x</sub></td><td>Total kinetic energy density</td></tr><tr><td>E&#773;<sub>K</sub></td><td><sup>1</sup>&frasl;<sub>2</sub> ρ<sub>0</sub> v&#773;<sub>x</sub> v&#773;<sub>x</sub></td><td>Barotropic kinetic energy density</td></tr><tr><td>E'<sub>K</sub></td><td><sup>1</sup>&frasl;<sub>2</sub> ρ<sub>0</sub> v'<sub>x</sub> v'<sub>x</sub></td><td>Baroclinic kinetic energy density</td></tr><tr><td>E<sub>P</sub></td><td>ρgz</td><td>Total potential energy density</td></tr><tr><td>E&#773;<sub>P</sub></td><td>ρ<sub>0</sub>gz</td><td>Barotropic potential energy density</td></tr><tr><td>E&#773;<sub>P</sub></td><td>gzρ</td><td>Baroclinic potential energy density</td></tr><tr><td>F<sub>T</sub></td><td>gρ<sub>0</sub>Dv<sub>x</sub>η</td><td>Surface tide energy flux</td></tr><tr><td>F<sub>X</sub></td><td>∫<sub>0</sub><sup>1</sup> Dv<sub>x</sub>Pdσ</td><td>Barotropic internal pressure flux</td></tr><tr><td>(F<sub>T</sub> + F<sub>X</sub>)</td><td>gρ<sub>0</sub> ∂v<sub>x</sub>η / ∂v<sub>x</sub> + ∫<sub>0</sub><sup>1</sup> ∂v<sub>x</sub>P / ∂v<sub>x</sub> dσ</td><td>Divergence of barotropic flux of surface plus internal pressure</td></tr><tr><td>φ<sub>T</sub></td><td>gρ<sub>0</sub>η ∂v<sub>x</sub> / ∂η</td><td>Local E&#773;<sub>K</sub> to E&#773;<sub>P</sub> conversion</td></tr><tr><td>φ<sub>z</sub> = φ&#773;<sub>z</sub> + φ'<sub>z</sub></td><td>D ∫<sub>0</sub><sup>1</sup> ρg v<sub>x</sub> dσ = D ∫<sub>0</sub><sup>1</sup> ρ&#773;g (v&#773;<sub>x</sub> + v'<sub>x</sub>) dσ</td><td>Total equals barotropic plus</td></tr><tr><td>φ<sub>B</sub> = φ&#773;<sub>B</sub> + φ'<sub>B</sub></td><td>D ∫<sub>0</sub><sup>1</sup> ρg v<sub>x</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub> dσ = D ∫<sub>0</sub><sup>1</sup> ρ&#773;g (v&#773;<sub>x</sub> + v'<sub>x</sub>) ∂v<sub>x</sub> / ∂v<sub>x</sub> dσ</td><td>Baroclinic buoyancy conversion term between E&#773;<sub>K</sub>, E'<sub>K</sub> and E&#773;<sub>P</sub></td></tr><tr><td>φ<sub>S</sub>, φ<sub>σ</sub></td><td>D ∫<sub>0</sub><sup>1</sup> ρg (σ ∂η / ∂η) dσ, D ∫<sub>0</sub><sup>1</sup> ρ&#773;g v<sub>σ</sub> dσ</td><td>Total equals barotropic plus baroclinic buoyancy term due to surface and topographic slope</td></tr><tr><td>E&#773;<sub>H</sub></td><td>ρ<sub>0</sub> v&#773;<sub>x</sub> ∫<sub>0</sub><sup>1</sup> (DK<sup>H</sup><sub>β</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>) dσ</td><td>Buoyancy term due to local vertical surface motion and v<sub>σ</sub></td></tr><tr><td>E&#773;<sub>V</sub></td><td>ρ<sub>0</sub> v<sub>x</sub> [K<sup>g</sup><sub>D</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>]<sup>1</sup> ∫<sub>0</sub><sup>1</sup></td><td>E&#773;<sub>K</sub>-dissipation of hor. mixing</td></tr><tr><td>F<sub>E&#773;<sub>K</sub></sub></td><td>∂v<sub>x</sub>E&#773;<sub>K</sub> / ∂v<sub>x</sub></td><td>E&#773;<sub>K</sub>-dissipation by vertical mixing</td></tr><tr><td>φ<sub>3D2D</sub></td><td>ρ<sub>0</sub> v<sub>x</sub> ∫<sub>0</sub><sup>1</sup> ∂v'<sub>x</sub>v'<sub>x</sub> / ∂v<sub>x</sub> dσ</td><td>Divergence of barotropic advective E&#773;<sub>K</sub>-flux</td></tr><tr><td>F<sub>IW</sub></td><td>∫<sub>0</sub><sup>1</sup> ∂v<sub>x</sub>P / ∂v<sub>x</sub> dσ</td><td>Non-lin. interaction of E&#773;<sub>K</sub> and E'<sub>K</sub></td></tr><tr><td>E&#773;<sub>H</sub></td><td>ρ<sub>0</sub> ∫<sub>0</sub><sup>1</sup> v'<sub>x</sub> ∂(DK<sup>H</sup><sub>β</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>) dσ</td><td>Divergence of baroclinic flux of internal pressure</td></tr><tr><td>E&#773;<sub>V</sub></td><td>ρ<sub>0</sub> ∫<sub>0</sub><sup>1</sup> v'<sub>x</sub> ∂(K<sup>g</sup><sub>D</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>) dσ</td><td>E&#773;<sub>K</sub>-dissipation of horizontal mixing</td></tr><tr><td>F<sub>E&#773;<sub>K</sub></sub></td><td>∫<sub>0</sub><sup>1</sup> ∂v'<sub>x</sub>E&#773;<sub>K</sub> / ∂v<sub>x</sub> dσ</td><td>E&#773;<sub>K</sub>-dissipation of vertical mixing</td></tr><tr><td>F<sub>X1</sub></td><td>∫<sub>0</sub><sup>1</sup> ∂v<sub>x</sub>E&#773;<sub>K</sub> / ∂v<sub>x</sub> dσ</td><td>Divergence of barotropic E&#773;<sub>K</sub>-flux</td></tr><tr><td>F<sub>X2</sub></td><td>ρ<sub>0</sub> ∫<sub>0</sub><sup>1</sup> ∂v'<sub>x</sub>v'<sub>x</sub>E&#773;<sub>K</sub> / ∂v<sub>x</sub> dσ</td><td>Divergence of baroclinic (v'<sub>x</sub>, v&#773;<sub>x</sub>)-flux</td></tr><tr><td>F<sub>E<sub>P</sub></sub> = F&#773;<sub>E<sub>P</sub></sub> + F'<sub>E<sub>P</sub></sub></td><td>∫<sub>0</sub><sup>1</sup> ∂v<sub>x</sub>E<sub>P</sub> / ∂v<sub>x</sub> dσ = ∫<sub>0</sub><sup>1</sup> ∂(E&#773;(v<sub>x</sub> + v'<sub>x</sub>))E<sub>P</sub> / ∂v<sub>x</sub> dσ</td><td>Divergence of total equals barotropic plus baroclinic E&#773;<sub>P</sub>-flux.</td></tr><tr><td>F<sub>m,H</sub></td><td>∫<sub>0</sub><sup>1</sup> ∂(gzDK<sup>H</sup><sub>β</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>) dσ</td><td>Horizontal divergence of diffusive E&#773;<sub>P</sub>-flux</td></tr><tr><td>F<sub>m,V</sub></td><td>[gzK<sup>g</sup><sub>D</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>]<sup>1</sup> ∫<sub>0</sub><sup>1</sup></td><td>Top and bottom vertical diffusive E&#773;<sub>P</sub>-flux</td></tr><tr><td>φ<sub>m,H</sub></td><td>∫<sub>0</sub><sup>1</sup> gDK<sup>H</sup><sub>β</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub> dσ</td><td>Local E&#773;<sub>P</sub>-transfer due to horizontal density diffusion</td></tr><tr><td>φ<sub>m,V</sub></td><td>∫<sub>0</sub><sup>1</sup> gK<sup>V</sup> ∂v<sub>x</sub> / ∂v<sub>x</sub> dσ</td><td>Local E&#773;<sub>P</sub>-transfer due to vertical density diffusion</td></tr><tr><td>F<sub>d,H</sub></td><td>∫<sub>0</sub><sup>1</sup> (gzDK<sup>H</sup><sub>β</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>) dσ</td><td>Horizontal divergence of diffusive E&#773;<sub>P</sub>-flux</td></tr><tr><td>F<sub>d,V</sub></td><td>[gz<sup>g</sup><sub>D</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub>]<sup>1</sup> ∫<sub>0</sub><sup>1</sup></td><td>Surface and bottom vert. diffusive E&#773;<sub>P</sub>-flux</td></tr><tr><td>φ<sub>d,H</sub></td><td>∫<sub>0</sub><sup>1</sup> gDK<sup>H</sup><sub>β</sub> ∂v<sub>x</sub> / ∂v<sub>x</sub> dσ</td><td>Local E&#773;<sub>P</sub>-transfer due to horizontal diffusion</td></tr><tr><td>φ<sub>d,V</sub></td><td>∫<sub>0</sub><sup>1</sup> gK<sup>V</sup> ∂v<sub>x</sub> / ∂v<sub>x</sub> dσ</td><td>Local E&#773;<sub>P</sub>-transfer due to vertical diffusion</td></tr><tr><td>φ<sub>m,H</sub> - φ<sub>d,H</sub> = ∫<sub>0</sub><sup>1</sup> gDK<sup>H</sup><sub>β</sub> ∂(v<sub>x</sub> - v<sub>x</sub>) / ∂v<sub>x</sub> dσ</td><td>Local E&#773;<sub>A</sub>-transfer due to vertical diffusion</td></tr><tr><td>φ<sub>m,V</sub> - φ<sub>d,V</sub> = ∫<sub>0</sub><sup>1</sup> gK<sup>V</sup> ∂(v<sub>x</sub> - v<sub>x</sub>) / ∂v<sub>x</sub> dσ</td><td>Local E&#773;<sub>A</sub>-transfer due to vertical diffusion</td></tr></tbody></table>

consistency in the local buoyancy conversion term, such as due to flow along $\sigma$-planes, which is found by substitution of $B^{i+\frac{1}{2}} = \frac{1}{2}g(\hat{\rho}^{i+1} + \hat{\rho}^i)v_x^{i+\frac{1}{2}}$ in (B2). Its counterpart is retrieved from the kinetic energy equations by vertical integration of the rmsum of (B1) with $A^i = P^i$ and requires the hydrostatic and continuity equations (cf. Marsaleix et al., 2008, Section 3.3.1).

## B.2. (Diffusive) momentum fluxes

The momentum fluxes are defined as in Blumberg and Mellor (1987), so that on a mass-point at vertical level $k$ the longitudinal flux of zonal momentum is

$$F_{m,x,L}^{ij} = \frac{1}{4} (\hat{v}_x^{i+\frac{1}{2}j} + \hat{v}_x^{i-\frac{1}{2}j}) (v_x^{i+\frac{1}{2}j} + v_x^{i-\frac{1}{2}j}),$$

where $\hat{v}_x^{ij} = D^{ij}v_x^{jj} \Delta\sigma^{ij}$ and the subscript $m, x, L$ indicates a zonal (i.e. longitudinal) flux of zonal momentum. In contrast, a transversal

momentum flux is for example the meridional flux of zonal momentum,

$$F_{m,x,T}^{ij} = \frac{1}{4} (\hat{v}_y^{i+\frac{1}{2}} + \hat{v}_y^{i-\frac{1}{2}})^{j-\frac{1}{2}} (v_x^j + v_x^{j-1})^i.$$

Integrated over a horizontal area at level $k$, the latter leads to the energy flux decomposition

$$-\sum_{i=1}^{i2} \sum_{j=1}^{j2} v_x^{i-\frac{1}{2}} (F_{mT}^{i-\frac{1}{2}j+\frac{1}{2}} - F_{mT}^{i-\frac{1}{2}j-\frac{1}{2}}) = \sum_{i=1}^{i2} \sum_{j=1}^{j2-1} F_{mT}^{i-\frac{1}{2}j+\frac{1}{2}} (v_x^{i-\frac{1}{2}j+1} - v_x^{i-\frac{1}{2}j}) \\ - \sum_{i=1}^{i2} v_x^{i-\frac{1}{2}j2} F_{mT}^{i-\frac{1}{2}j2+\frac{1}{2}} \\ + v_x^{i-\frac{1}{2}j1} F_{mT}^{i-\frac{1}{2}j1-\frac{1}{2}},$$

where the boundary flux similarly has velocity on the interior grid-point, but the summation limits of the local conversion term differ
from the form presented in paragraph B.1. due to the discretisation onto the staggered C-grid.

In case of a constant diffusion coefficient, the diffusive flux of momentum along a σ-level in the Ox-direction is of the form $Fd_{m,x,L}^{i,j,k} = \frac{K^H}{\Delta x} (D\Delta\sigma)^{i,j,k,t} (v_x^{i+\frac{1}{2}} - v_x^{i-\frac{1}{2}})$.

Note that $(D\Delta\sigma)^{i,j,k,t}$ is averaged onto mass-points.

### B.3. Horizontal (diffusive) flux of scalars

The potential energy flux terms along Ox labelled $B^{i+\frac{1}{2},j,k,t}$ are firstly the horizontal scalar flux,

$$-v_x^{i,j,k} \overline{D_x \Delta\sigma^x} (\overline{\psi^{i,j,k}})^x,$$

where $\psi$ represents $\tilde{T}$ or $\tilde{S}$. Secondly, the horizontal scalar diffusive flux is defined as

$$\frac{K^H}{\Delta x} (\psi^{i,j,k;t-\Delta t} - \psi^{i-1,j,k;t-\Delta t}) D^{i,j;t-\Delta t} \Delta\sigma^{i,j,k},$$

here $\tilde{T}$ and $\tilde{S}$ are evaluated at the previous time-step for numerical stability reasons (Marsaleix et al., 2008).

### B.4. Coriolis effect

Although no net work is done in the inner domain (Arakawa and Lamb, 1977), the Coriolis terms at the lateral boundaries pose a non-vanishing energy source or sink, which is illustrated by the integral of the Coriolis-contribution to the kinetic energy equations over a horizontal domain at an arbitrary σ-level, of the form

$$\sum_{i=i_1 j=j_1+1}^{i_2 j_2} v_y^{i,j-\frac{1}{2}} \overline{(F(v_x)^x)^{y:i_2 j-\frac{1}{2}}} - \sum_{i=i_1+1 j=j_1}^{i_2 j_2} v_x^{i-\frac{1}{2}j} \overline{(F(v_y)^y)^{x:i_1-\frac{1}{2}j}} \\ = \sum_{j=j_1+1}^{j_2-1} \left( F_{i_1 j} v_x^{i_1-\frac{1}{2}j} \overline{(v_y)^{y:i_1 j}} + F_{i_2 j} v_x^{i_2+\frac{1}{2}j} \overline{(v_y)^{y:i_2 j}} \right) \\ - \sum_{i=i_1+1}^{i_2-1} \left( F_{i_1 j_1} v_y^{i_1-\frac{1}{2}j_1} \overline{(v_x)^{x:i_1 j_1}} + F_{i_2 j_1} v_y^{i_2+\frac{1}{2}j_1} \overline{(v_x)^{x:i_2 j_1}} \right) \\ + F_{i_1 j_1} \left( v_x^{i_1-\frac{1}{2}j_1} v_y^{i_1 j_1+\frac{1}{2}} - v_x^{i_1+\frac{1}{2}j_1} v_y^{i_1 j_1-\frac{1}{2}} \right) \\ + F_{i_2 j_1} \left( v_x^{i_2+\frac{1}{2}j_1} v_y^{i_2 j_1+\frac{1}{2}} - v_x^{i_2-\frac{1}{2}j_1} v_y^{i_2 j_1-\frac{1}{2}} \right) \\ + F_{i_1 j_2} \left( v_x^{i_1-\frac{1}{2}j_2} v_y^{i_1 j_2+\frac{1}{2}} - v_x^{i_1+\frac{1}{2}j_2} v_y^{i_1 j_2-\frac{1}{2}} \right) \\ + F_{i_2 j_2} \left( v_x^{i_2+\frac{1}{2}j_2} v_y^{i_2 j_2-\frac{1}{2}} - v_x^{i_2-\frac{1}{2}j_2} v_y^{i_2 j_2+\frac{1}{2}} \right). \quad (B3)$$

Here, $F^{ij} = D^{ijf}$ and the following averaging operators were used, for arbitrary variable $\psi$ and grid-point $i, j$ and fixed $\sigma$, namely $\overline{(\psi)^{x,i,j}} = \frac{1}{2} (\psi^{i+\frac{1}{2}j} + \psi^{i-\frac{1}{2}j})$ along Ox and $\overline{(\psi)^{y,i,j}} = \frac{1}{2} (\psi^{i,j+\frac{1}{2}} + \psi^{i,j-\frac{1}{2}})$ along Oy. The sums on the right-hand side in (B3) represent energy fluxes through the domain boundaries, while the remaining terms correspond to each of the four corners of the domain. In a closed domain the velocities perpendicular to the boundaries at $i_1-1, i_2+1, j_1-1$ and $j_2+1$ are set to zero, so that in effect $v_x^{i_1-\frac{1}{2}j} = v_x^{i_2+\frac{1}{2}j} = 0$ and $v_y^{i_1 j_1-\frac{1}{2}} = v_y^{i_2 j_2+\frac{1}{2}} = 0$, and therefore (B3) vanishes and kinetic energy is indeed conserved. It is straightforward to define the energy fluxes due to the Coriolis effect associated with the external and internal mode, by replacing $v_x^{ij}$ by respectively $\overline{v_x^{ij}}$ and $v_x^{ij}$.

## Appendix C. Free-surface, sigma-coordinate sorting algorithm for background state computation

The computation of the diapycnal mixing at location $(i,j,k)$ requires knowledge of $Z^*(\bar{\rho}^{i,j,k})$, the depth of a fluid parcel of density $\bar{\rho}^{i,j,k}$ in the adiabatically reorganised “background” state. This

computation is done based on a sorting algorithm (Winters et al., 1995). The originality of the present implementation is that it is applied to a σ-coordinate, free surface model. The minimisation of the potential energy implies in particular that the background state free surface be levelled. To this end, we propose a five-step algorithm.

Step 1: in the original domain $(R_\eta)$, mass parcels $(\bar{\rho}\Delta V)^{i,j,k}$ are sorted based on their density, into a volume distribution $\Delta V_\eta(\bar{\rho})$.

Step 2: neighbouring mass parcels in $\Delta V_\eta(\bar{\rho})$ are gathered together if their densities are 'too close', giving the (possibly lower resolution) distribution $\Delta V_2(\bar{\rho})$. By 'too close', it is meant that the model vertical resolution imposes a minimum density difference that can be evaluated, such as $\Delta \bar{\rho} \approx \frac{1}{K}(\bar{\rho}_{max} - \bar{\rho}_{min})$, where $\bar{\rho}_{min}$ and $\bar{\rho}_{max}$ are respectively, the minimum and maximum density in the domain, and $K$ is the number of vertical levels.

Step 3: a 'levelled' domain $R_L$ is defined using the domain average of the original free surface elevation. The volume parcels in $R_L$ are sorted based on their depth, into a volume distribution $\Delta V_L(Z^*)$. It is important to note that the total volume in $R_L$ equals that of the original domain $R_\eta$.

Step 4: the volume distribution $\Delta V_L(Z^*)$ is now filled from the bottom up with mass parcels from $\Delta V_2(\bar{\rho})$, starting with the densest parcel, giving the vertical density distribution $\bar{\rho}(Z^*)$. Thus, the volume in $R_L$ below the vertical level $Z^*(\bar{\rho})$ contains parcels with density $\rho > \bar{\rho}(Z^*)$.

Step 5: The resulting depth-distribution $Z^*(\bar{\rho})$ of density $\bar{\rho}$ in domain $R_L$ is finally interpolated onto the parcels $\Delta V^{i,j,k}$ in $R_\eta$ with corresponding $\bar{\rho}$, giving the spatial distribution of $Z^*(\bar{\rho}^{i,j,k})$.

## References

Apel, J.R., 1987. Principles of Ocean Physics. Academic Press, London.

Arakawa, A., Lamb, V.R., 1977. Computational design of the basic dynamical processes of the UCLA general circulation model. Methods Comput. Phys. 17, 173-265.

Asselin, R., 1972. Frequency filter for time integrations. MWR 100 (6), 487-490.

Auclair, F. et al., 2010. A non-hydrostatic algorithm for free-surface ocean modelling. Ocean Modell. 36, 49-70.

Baines, P.G., 1982. On internal tide generation models. DSR 29 (3A), 307-338.

Bell Jr., T., 1975. Topographically generated internal waves in the open ocean. JGR 80, 320-327.

Blumberg, A.F., Mellor, G.L., 1987. A description of a three-dimensional coastal ocean circulation model. In: Heaps, N.S. (Ed.), Three-Dimensional Coastal Ocean Models. American Geophysical Union, Washington DC, pp. 1-16.

Bühler, O., Muller, C.J., 2007. Instability and focusing of internal tides in the deep ocean. J. Fluid Mech. 558, 1-28.

Carter, G.S. et al., 2008. Energetics of M$_2$ barotropic-to-baroclinic tidal conversion at the Hawaiian islands. JPO 28, 2205-2223.

Cartwright, D.E., Ray, R.D., 1991. Energetics of global ocean tides from Geosat altimetry. JGR 96, 16897-16912.

Di Lorenzo, E. et al., 2006. Numerical and analytical estimates of M$_2$ tidal conversion at steep oceanic ridges. JPO 36, 1072-1084.

Egbert, G.D., Ray, R.D., 2001. Estimates of M$_2$ tidal energy dissipation from TOPEXPoseidon altimeter data. JPO 106, 22475-22502.

Floor, J.W., 2009. Energetics of Internal Tide Generation, Propagation and Dissipation. PhD Thesis, l'Université Toulouse III - Paul Sabatier, Toulouse.

Garrett, C., Kunze, E., 2007. Internal tide generation in the deep ocean. Annu. Rev. Fluid Mech. 39, 57-87.

Gerkema et al., 2004. Internal tides in the Bay of Biscay: conversion rates and seasonal effects. DSR-II 51, 2995-3008.

Gerkema et al., 2006a. Decay of semi-diurnal internal-tide beams due to subharmonic resonance. GRL 33, L08604.

Gerkema et al., 2006b. Non-linear effects in internal-tide beams, and mixing. Ocean Modell. 12, 302-318.

Gerkema, T., 2006c. Internal-wave reflection from uniform slopes: higher harmonics and Coriolis effects. Nonlinear Process. Geophys. 13, 265-273.

Gerkema, T., Zimmerman, J.T.F., 1995. Generation of nonlinear internal tides and solitary waves. JPO 25 (6), 1081-1094.

Gill, A.E., 1982. Atmosphere-Ocean Dynamics. Academic Press, London.

Holloway, P.E., Merrifield, M.A., 1999. Internal tide generation by seamounts, ridges and islands. JGR 104 (C11), 25937-25951.
Huang, R.X., 1998. Mixing and available potential energy in a Boussinesq ocean. JPO 28, 669-678.

Hyder et al., 2005. Observations and predictability of internal solitons in the northern Andaman Sea. AOR 27, 1-11.

Johns, B. et al., 1983. Simulation of storm surges using a three-dimensional numerical model: an application to the 1977 Andhra Cyclone. Quart. J. Roy. Meteorol. Soc. 109, 211-224.

Katsumata, K., 2006. Tidal stirring and mixing on the Australian North West Shelf. Mar. Freshwater Res. 57, 243-254.

Khatiwala, S., 2003. Generation of internal tides in an ocean of infinite depth: analytical and numerical calculations. DSR-I 50, 3-21.

Lamb, K.G., 2007. Energy and pseudoenergy flux in the internal wave field generated by tidal flow over topography. CSR 27, 1208-1232.

Legg, S., Huijts, K.M.H., 2006. Preliminary simulations of internal waves and mixing generated by finite amplitude tidal flow over isolated topography. DSR-II 53, 140-156.

Lien, R.-C., Gregg, M.C., 2001. Observations of turbulence in a tidal beam and across a coastal ridge. JPO 106, 4575-4591.

Llewellyn Smith, S.G., Young, W.R., 2002. Conversion of the barotropic tide. JPO 32, 1554-1566.

Lorenz, E.N., 1955. Available potential energy and the maintenance of the general circulation. Tellus VII (2), 157-167.

Marsaleix, P. et al., 2008. Energy conservation issues in sigma-coordinate free-surface ocean models. Ocean Modell. 20, 61-89.

Marsaleix, P. et al., 2009a. Open boundary conditions for internal gravity waves modelling using polarization relations. Ocean Modell. 29, 27-42.

Marsaleix, P. et al., 2009b. Low-order pressure gradient schemes in sigma coordinate models: the seamount test revisited. Ocean Modell. 30, 169-177.

Mellor, G.L., Blumberg, A.F., 1985. Modeling vertical and horizontal diffusivities with the sigma coordinate system. MWR 113 (8), 1379-1383.

Munk, W., Wunsch, C., 1998. Abyssal recipes II: energetics of tidal and wind mixing. DSR-I 45, 1977-2010.

Munroe, J.R., Lamb, K.G., 2005. Topographic amplitude dependence of internal wave generation by tidal forcing over idealized three-dimensional topography. JGR 110, C02001. doi:10.1029/2004JC002537.

Niwa, Y., Hibiya, T., 2001. Numerical study of the spatial distribution of the M₂ internal tide in the Pacific Ocean. JPO 106, 22441-22449.

Pairaud, I., Auclair, F., 2005. Combined wavelet and principal component analysis (WEOF) of a scale oriented model of coastal ocean gravity waves. Dynam. Atmos. Oceans 40 (4), 254-282.

Peacock, T., Echeverri, P., Balmforth, N.J., 2008. An experimental investigation of internal tide generation by two-dimensional topography. JPO 38, 235-242.

Pedlosky, J., 1987. Geophysical Fluid Dynamics, second ed. Springer-Verlag, New York.

Pedlosky, J., 2003. Waves in the Ocean and Atmosphere – Introduction to Wave Dynamics. Springer-Verlag, Berlin.

Pérenne, N. et al., 2000. A numerical study of stratified tidal rectification. CSR 20, 37-68.

Pétrélis, F. et al., 2006. Tidal conversion at a submarine ridge. JPO 36, 1053-1071.

Ray, R.D., Mitchum, G.T., 1997. Surface manifestation of internal tides in the deep ocean: observations from altimetry and island gauges. Prog. Oceanogr. 40, 135-162.

Robert, A.J., 1966. The integration of a low order spectral form of the primitive meteorological equations. J. Meteorol. Soc. Jpn. Ser. 44 (5), 237-245.

Rudnick et al., 2003. From tides to mixing along the Hawaiian ridge. Science 301, 355-357.

Scotti, A. et al., 2006. On the interpretation of energy and energy fluxes of nonlinear internal waves: an example from Massachusetts Bay. J. Fluid. Mech. 561, 103-112.

Shaw, P.-T., Ko, D.S., Chao, S.-Y., 2009. Internal solitary waves induced by flow over a ridge: with applications to the northern South China Sea. JGR 114, C02019. doi:10.1029/2008JC005007.

Shepherd, T.G., 1993. A unified theory of available potential energy. Atmos. Ocean 31, 1-26.

Simmons, H.L., 2008. Spectral modification and geographic redistribution of the semi-diurnal internal tide. Ocean Modell. 21, 126-138.

Thorpe, S.A., 2005. The Turbulent Ocean. Cambridge University Press, Cambridge.

Tseng, Y., Ferziger, J., 2001. Mixing and available potential energy in stratified flows. Phys. Fluids 13 (5), 1281-1293.

Winters, K.B. et al., 1995. Available potential energy and mixing in density-stratified fluids. J. Fluid Mech. 289, 115-128.

Wunsch, Ferrari, 2004. Vertical mixing, energy, and the general circulation of the oceans. Annu. Rev. Fluid Mech. 36, 281-314.