# Implementation of a time-dependent bathymetry in a free-surface ocean model: Application to internal wave generation

Francis Auclair, Lucie Bordois, Yvan Dossmann, Thomas Duhaut, Claude Estournel, Jochem Willem Floor, Patrick Marsaleix, Cyril Nguyen, Alexandre Paci, Caroline Ulses

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

## To cite this version:

Francis Auclair, Lucie Bordois, Yvan Dossmann, Thomas Duhaut, Claude Estournel, et al. Implementation of a time-dependent bathymetry in a free-surface ocean model: Application to internal wave generation. Ocean Modelling, 2014, 80, pp.1-9. 10.1016/j.ocemod.2014.03.002. hal-02347345

HAL Id: hal-02347345

https://hal.science/hal-02347345

Submitted on 28 Sep 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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Implementation of a time-dependent bathymetry in a free-surface ocean model: Application to internal wave generation

F. Auclair<sup>a,*</sup>, L. Bordois<sup>a,b</sup>, Y. Dossmann<sup>b,a</sup>, T. Duhaut<sup>a</sup>, C. Estournel<sup>a</sup>, J.W. Floor<sup>a</sup>, P. Marsaleix<sup>a</sup>, C. Nguyen<sup>a</sup>, A. Paci<sup>b</sup>, C. Ulses<sup>a</sup>

We show that a time-varying bathymetry can be implemented in a σ-coordinate free-surface ocean model with only slight modifications to the original numerical algorithm. A consistent choice of variables for the height of the water column is proposed. The resulting algorithm creates new opportunities in ocean model configurations in which bathymetry changes at time scales comparable with those of ocean dynamics (e.g., near shore bed motion and seismic or volcanic catastrophes) or in direct modeling of laboratory experiments with a moving bottom. The new numerical implementation is carefully compared with a laboratory experiment involving internal wave generation by an oscillating, ridge-shaped topography.

An energy formulation of the time-varying bathymetry is proposed, and energy fluxes are computed to characterize and quantify the cascade of energy leading to the generation of internal waves by an oscillating ridge.

## 1. Introduction

Most ocean models postulate that bottom topography is time independent. Indeed, topography variations associated with slow geological events usually occur on timescales much longer than the timescales of ocean dynamics (hundreds to millions of years for geological events compared with a few decades at most for ocean dynamics). There are several circumstances in which this assumption is no longer valid. High frequency variations in ocean bathymetry can occur as a result of the motion of sand banks (Williams et al., 2000) or of sediment re-suspension or deposition (Le Hir et al., 2011). The bathymetry of the prodelta region of the Amazon shelf can change from one season to another (Le Bars et al., 2010). Tides or storm surges generate gravitational variations and can cause deformations of ocean bathymetry with amplitudes of a few tens of millimeters (Le Provost et al., 1998; Fratepietro et al., 2006). On even shorter timescales, seismic motion, volcanic eruptions or even asteroid impacts can induce dramatic changes in the water column depth, generating tsunami waves (Heinrich et al., 1999; Pavlov et al., 2009; Weiss et al., 2006).

At the laboratory scale, the bottom surface in water tanks can also be varied to set water in motion (Dossmann et al., 2011). For example, using “horizontally oscillating topography” can be an efficient way of studying tide-topography interactions at the laboratory scale. Gerkema and Zimmerman (2005) demonstrated that for a linear Boussinesq regime, such a (laboratory-based) modeling setup can be considered “equivalent to the original problem” of internal tide generation observed from an accelerated frame of reference co-moving with the surface tide.

We show hereafter that a time-varying bathymetry can be implemented in a free-surface, non-hydrostatic Boussinesq model with few modifications. The free-surface assumption can be associated with several types of vertical coordinates ranging from σ-vertical coordinates (Stacey et al., 1995) to pressure-coordinate or z-coordinate models (Griffies et al., 2010; Marshall et al., 2004; Losch et al., 2004). Whatever the vertical coordinate used, an anomaly of the free surface results in a stretching of the entire water column. To allow variations in the bathymetry, this stretching can be generalized to the uplift or to the subsidence of the sea floor. As for the evolution of the free surface, the water column is separated into two components: one is time-invariant (the mean depth), and the other is time-variant and allows for a varying sea floor. This latter component is the sum of the surface and bathymetry anomalies.

<footnote><sup>a</sup>SIROCCO Group, Laboratoire d'Aérologie, Université de Toulouse, CNRS, Toulouse, France</footnote>

<footnote><sup>b</sup>CNRM-GAME, GMEI, SPEA, CNRS, Météo-France, Toulouse, France</footnote>

<footnote>* Corresponding author. Tel.: +33 561332775; fax: +33 561332790.</footnote>

<footnote>E-mail address: francis.auclair@aeroobs-mip.fr (F. Auclair).</footnote>
As a consequence, the transformation of a free-surface bathymetry model to a free-surface, time-varying bathymetry model requires a few modifications of the original algorithm: (i) a new definition of the water column anomaly variable (it must now include both the free-surface and bathymetry anomalies), (ii) the vertical positioning of the calculation grid, (iii) the bottom boundary conditions, and (iv) the computation of the barotropic component of the pressure gradient.

The resulting algorithm is implemented in the non-hydrostatic algorithm proposed by Auclair et al. (2011) and is applied to a direct numerical simulation of a laboratory experiment of internal tide generation in a CNRM-GAME tank.

The proposed implementation of time-varying bathymetry is detailed in Section 2. In Section 3, a direct simulation of a laboratory configuration of internal tides generated by a moving ridge is proposed and compared with Synthetic Schlieren measurements.

## 2. Implementation of a time-varying bathymetry

### 2.1. Dynamical equations

The non-hydrostatic kernel proposed by Auclair et al., 2011 has been extended to allow time variations in the bathymetry, but the proposed algorithm can be implemented in any type of free surface model.

We designated $D(t)$ as the total depth of the water column ($D(t) > 0$). If $\zeta(t)$ and $\Delta H(t)$ are the surface and bottom anomalies of the water column, respectively, and $H$ is the time-independent mean depth, we can write $D(t) = H + \zeta(t) + \Delta H(t)$. The $\sigma$-coordinates can then be defined by the depth of the $\sigma$-level:

$$z(x, y, \sigma; t) = \sigma D(x, y; t) - H(x, y) - \Delta H(x, y; t) \quad (1)$$

At the sea surface ($z = \zeta$), $\sigma$ equals 1, whereas at the bottom ($z = -(H + \Delta H)$), $\sigma$ vanishes. The true vertical velocity $v_z$ can be written as follows:

$$\begin{aligned} v_z(x, y, \sigma; t) &\equiv \frac{dz}{dt} = v_{\sigma} + \left. \frac{dz}{dt} \right|_{\sigma} \\ &= v_{\sigma} + \sigma \left. \frac{d\zeta}{dt} \right|_{\sigma} + (\sigma - 1) \frac{d(H + \Delta H)}{dt} \end{aligned} \quad (2)$$

where the subscript “$\sigma$” refers to the derivation along a $\sigma$-surface. Because of the time-dependency of the bottom bathymetry, the last term of this expression now can be defined as follows:

$$\frac{d(H + \Delta H)}{dt} = \frac{\partial \Delta H}{\partial t} + v_{\sigma} \frac{\partial (H + \Delta H)}{\partial x_{\sigma}} \quad (3)$$

where $\alpha = 1$ (for $x$) or 2 (for $y$). The density field can be written as $\rho_0 + \rho$, with $\rho_0$ representing a reference density and with $\rho$ representing the anomaly satisfying the equation of state, $\rho = \rho(\theta, S, P)$, with $\theta, S$ and $P$ representing potential temperature, salinity and “total” pressure, respectively. Potential temperature and salinity satisfy the heat and salt evolution equations:

$$\left. \frac{1}{D} \frac{\partial D \theta}{\partial t} \right|_{\sigma} = A_{\theta}^{\sigma} \quad (4)$$

$$\left. \frac{1}{D} \frac{\partial D S}{\partial t} \right|_{\sigma} = A_{S}^{\sigma} \quad (5)$$

On the right hand side of equations 4 and 5 above, “A” includes diffusion and advection (Auclair et al., 2011). Based on the decomposition of the density and on the hydrostatic balance, the total pressure can be further decomposed as follows:

$$p_0 + P = p_0 + p + q \quad (6)$$

where $p_0$ and $p$ are the hydrostatic components of the pressure field, satisfying $\left. \frac{\partial(p_0 + P)}{\partial z} \right|_z = \frac{\partial(p_0 + P)}{\partial \sigma} = -g(\rho_0 + \rho)$ and $\left. \frac{\partial p_0}{\partial x_\sigma} \right|_z = \rho_0 g \frac{\partial \zeta}{\partial x_\sigma}$. The remaining non-hydrostatic component of the pressure is given as $q$.

The horizontal and vertical momentum equations can then be written:

$$\left. \frac{1}{D} \frac{\partial D v_{\sigma}}{\partial t} \right|_{\sigma} = -g \left. \frac{\partial \zeta}{\partial x_{\sigma}} - \frac{1}{\rho_0} \left. \frac{\partial P}{\partial x_{\sigma}} \right|_{z} \right|_{\sigma} + A_{\sigma}^{\sigma} \quad (7)$$

$$\left. \frac{1}{D} \frac{\partial D v_z}{\partial t} \right|_{\sigma} = - \frac{1}{\rho h \sigma_0} \left. \frac{\partial P}{\partial z} - \frac{\rho g}{\rho_0} \right|_{\sigma} + A_3^{\sigma} \quad (8)$$

The variable “A” in Eqs. (7) and (8) includes advection and diffusion together with the Coriolis pseudo-force (for details, see Auclair et al. (2011)). Under the Boussinesq approximation, the continuity equation simplifies to the following equation:

$$\left. \frac{\partial v_{\sigma}}{\partial x_{\sigma}} \right|_{z} + \frac{\partial v_z}{\partial z} = 0 \quad (9)$$

Using “$\sigma$-coordinates”, this equation becomes:

$$\left. \frac{\partial D v_{\sigma}}{\partial x_{\sigma}} \right|_{\sigma} - \frac{\partial}{\partial \sigma} \left( \left. \frac{\partial z}{\partial x_{\sigma}} \right|_{\sigma} v_{\sigma} \right) + \frac{\partial v_z}{\partial \sigma} = 0 \quad (10)$$

The bottom and surface boundary conditions are given by the kinematic conditions:

$$v_z(\sigma = 0) = -\frac{\partial \Delta H}{\partial t} - v_{\sigma}(\sigma = 0) \frac{\partial(H + \Delta H)}{\partial x_{\sigma}} \quad (11)$$

$$v_z(\sigma = 1) = \frac{\partial \zeta}{\partial t} + v_{\sigma}(\sigma = 1) \frac{\partial \zeta}{\partial x_{\sigma}} \quad (12)$$

In addition to the bathymetry anomaly ($\Delta H$), both the horizontal and vertical components of the velocity of the moving bathymetry are thus specified as boundary conditions (see Auclair et al. (2011) for details). This velocity and the bathymetry anomaly must satisfy Eq. (11).

Integrating the continuity equation (10) from bottom to top and substituting the kinematic conditions, Eqs. (11) and (12), lead to the free-surface evolution equation:

$$\frac{\partial D \bar{v}_{\sigma}}{\partial x_{\sigma}} = - \frac{\partial (\zeta + \Delta H)}{\partial t} \quad (13)$$

where the overbar signifies the vertical average of the horizontal velocity: $\bar{v}_{\sigma} = \frac{1}{D} \int_{-H+\Delta H}^{H+\Delta H} v_{\sigma} dz$. Eqs. (13) and (10) provide exact volume conservation (rather than mass conservation) in the context of the Boussinesq assumption. In regard to the free-surface anomaly ($\zeta$), the consequences of high frequency variations of the bottom perturbation ($\Delta H$) on the evolution of the water column volume are dealt with in the barotropic mode of the mode-splitting algorithm presented in Auclair et al. (2011). Atmospheric pressure anomalies are neglected in the present study and the surface pressure vanishes:

$$P(\sigma = 1) = p(\sigma = 1) = q(\sigma = 1) = 0 \quad (14)$$

### 2.2. Formulation of time-varying bathymetry

In a free-surface, terrain-following coordinate algorithm, the total depth of the water column can only change as a result of the free surface anomaly. Several choices of variables may permit the implementation of a time-varying bathymetry. The best choice of variables is the one leading to an algorithm as close as possible (at least formally) to the original constant-bathymetry algorithm at the lowest possible computing cost.

In the original algorithm, the depth of the water column can be accurately described using two variables: ($\zeta, H$). The second
column of Table 1 summarizes the set of non-hydrostatic equations described by Auclair et al. (2011) in which the surface anomaly and the mean depth ($\zeta$, $H$) explicitly appear.

Several choices of variables can be proposed to implement time-varying bathymetry. An obvious choice is obtained by adding the bathymetry anomaly ($\Delta H$) to the original set of variables ($\zeta$, $\Delta H$, $H$). Another choice of variables is ($\eta$, $\Delta H$, $H$) = ($\zeta + \Delta H$, $\Delta H$, $H$). The third and fourth columns of Table 1 show the corresponding formulation of the resulting algorithm. Using the ($\zeta$, $\Delta H$, $H$) formulation, the surface pressure force is the only term or equation that does not need to be modified. The computation of the depth of the $\sigma$-layers, the total depth and thickness of the levels, the continuity equation, the bottom kinematic conditions and the free-surface anomaly evolution equation necessitate new numerical treatment as a result of the introduction of a time-varying bathymetry anomaly ($\Delta H$).

Another possible formulation of the algorithm is based on the three variables ($\eta$, $\Delta H$, $H$) = ($\zeta + \Delta H$, $\Delta H$, $H$) and leads to a simpler numerical implementation. The computation of the depth of the $\sigma$-layers needs to be slightly modified (the bathymetry anomaly $\Delta H$ must be subtracted), and the surface pressure force needs slight modifications as well. If a no-flux boundary condition is applied, the bottom kinematic condition can be satisfied automatically without changing the original formulation of Auclair et al. (2011) by summing the continuity and free-surface evolution equations from bottom to surface. Consequently, the implementation of this equation does not require any further modifications.

In conclusion, the second set of variables ($\eta$, $\Delta H$, $H$) leads to a small number of insignificant modifications associated with the depth of the $\sigma$-layer, with the surface pressure force and with the bottom kinetic condition. To switch from the original algorithm to this set of variables, free-surface anomaly variables must be associated with total water-column anomalies (i.e., surface and

bottom anomalies). Based on the advantages listed above, this second set of variables has been chosen in the present implementation. The original model variable for the surface anomaly is thus modified to become the total water-depth anomaly (i.e., $\eta = \zeta + \Delta H$). Such an implementation of time-varying bathymetry does not require additional computations. Thus, the numerical cost of the resulting algorithm remains (approximately) constant.

## 2.3. Evolution of the kinetic energy content

Once the free-surface $\sigma$-coordinate algorithm has been extended to allow time-varying bathymetry, the modification of the simulated dynamics can be analyzed in detail. Indeed, the possibility of bathymetry changing with time can be viewed as a new forcing in the ocean model. Tsunami generation by a violent modification of the submarine topography is a dramatic example of such forcing.

The motion of the topography can consequently either reduce or enhance the kinetic energy of the water column. The numerical version of the evolution equation of both the potential and kinetic energies has been discussed in Floor et al. (2011) regarding the generation of internal waves by the impingement of a barotropic tide over a time-invariant ridge. In Section 3, the evolution equation of the kinetic energy is considered for the generation of internal waves by an oscillating ridge. The evolution of the kinetic energy can be obtained by multiplying the momentum equations (7) and (8) by the corresponding velocities. Two terms are of particular interest in regard to the consequences of the variations of the bottom topography: the rate of work produced by the pressure force and the buoyancy flux. In this study, the investigation is limited to the depth integral of the energy exchanges to simplify the discussion and to avoid the analysis of most of the $\sigma$-coordinate specific contributions.

<table><caption>**Table 1**</caption><thead><tr><th></th><th>Time-independent<br>bathymetry<br>(ζ, H)</th><th>Time-varying bathymetry<br>(ζ, ΔH, H)</th><th>Time-varying bathymetry<br>(η, ΔH, H) = (ζ + ΔH, ΔH, H)</th></tr></thead><tbody><tr><td>Depth of<br>σ-layers</td><td>z = σ(H + ζ) - H</td><td>z = σ(H + ζ + ΔH) - H - ΔH</td><td>z = σ(H + η) - H - ΔH</td></tr><tr><td>Total depth and<br>Thickness of<br>levels</td><td>D = H + ζ<br>Δz = Δσ(H + ζ)</td><td>D = H + ζ + ΔH<br>Δz = Δσ(H + ζ + ΔH)</td><td>D = H + η<br>Δz = Δσ(H + η)</td></tr><tr><td>Surface pressure<br>force</td><td>-g(∂ζ/∂x<sub>α</sub>)</td><td>-g(∂ζ/∂x<sub>α</sub>)</td><td>-g(∂η - ΔH/∂x<sub>α</sub>)</td></tr><tr><td>Cross-σ vertical<br>velocity<br>continuity<br>equation</td><td>(∂hv<sub>α</sub>/∂x<sub>α</sub>) + (∂v<sub>σ</sub>/∂σ) + (∂ζ/∂t) = 0</td><td>(∂hv<sub>α</sub>/∂x<sub>α</sub>) + (∂v<sub>σ</sub>/∂σ) + (∂ζ/∂t) + (∂ΔH/∂t) = 0</td><td>(∂hv<sub>α</sub>/∂x<sub>α</sub>) + (∂v<sub>σ</sub>/∂σ) + (∂η/∂t) = 0</td></tr><tr><td>Bottom kinematic<br>condition</td><td>w(σ = 0) = -v<sub>α</sub>(σ = 0)(∂H/∂x<sub>α</sub>)</td><td>w(σ = 0) = -∂ΔH/∂t - v<sub>α</sub>(σ = 0)(∂H + ΔH/∂x<sub>α</sub>)</td><td>w(σ = 0) = -∂ΔH/∂t - v<sub>α</sub>(σ = 0)(∂H + ΔH/∂x<sub>α</sub>)</td></tr><tr><td>Free-surface<br>anomaly<br>evolution</td><td>(∂Dv̄<sub>α</sub>/∂x<sub>α</sub>) + (∂ζ/∂t) = 0</td><td>(∂Dv̄<sub>α</sub>/∂x<sub>α</sub>) + (∂ζ + ΔH/∂t) = 0</td><td>(∂Dv̄<sub>α</sub>/∂x<sub>α</sub>) + (∂η/∂t) = 0</td></tr></tbody></table>
Using the chain rule several times, the work rate of the depth-integrated horizontal component of the pressure force associated with the total anomaly P can be written as follows:

$$ 
\begin{aligned}
\int_0^1 -D v_\alpha \frac{\partial P}{\partial x_\alpha} |_z d\sigma &= \int_0^1 -D v_\alpha \frac{\partial P}{\partial x_\alpha} |_\sigma d\sigma + \int_0^1 v_\alpha \frac{\partial z}{\partial x_\alpha} |_\sigma \frac{\partial P}{\partial \sigma} d\sigma \\
&= \int_0^1 -\frac{\partial D P v_\alpha}{\partial x_\alpha} |_\sigma d\sigma + \int_0^1 P \frac{\partial D v_\alpha}{\partial x_\alpha} |_\sigma d\sigma \\
&\quad + \int_0^1 v_\alpha \frac{\partial P}{\partial \sigma} \frac{\partial z}{\partial x_\alpha} |_\sigma d\sigma \\
&= F + \int_0^1 P \frac{\partial D v_\alpha}{\partial x_\alpha} |_\sigma d\sigma \\
&\quad - \int_0^1 P \frac{\partial}{\partial \sigma} (v_\alpha \frac{\partial z}{\partial x_\alpha} |_\sigma) d\sigma + [P v_\alpha \frac{\partial z}{\partial x_\alpha} |_\sigma]_0^1
\end{aligned} 
\quad (15) $$

where $F = \int_0^1 -\frac{\partial D P v_\alpha}{\partial x_\alpha} |_\sigma d\sigma$ is the flux rate resulting from the pressure force though the open boundary. The depth integral of the vertical component is

$$ \int_0^1 -v_z \frac{\partial P}{\partial \sigma} d\sigma = -[P v_z]_0^1 + \int_0^1 P \frac{\partial v_z}{\partial \sigma} d\sigma \quad (16) $$

Eqs. (15) and (16) can be summed, leading to the following equation:

$$ 
\begin{aligned}
\int_0^1 -D v_\alpha \frac{\partial P}{\partial x_\alpha} |_z d\sigma + \int_0^1 -v_z \frac{\partial P}{\partial \sigma} d\sigma \\
&= F + \left[ P \left( -v_z + v_\alpha \frac{\partial z}{\partial x_\alpha} \right) \right]_0^1 \\
&\quad + \int_0^1 P \left( \frac{\partial D v_\alpha}{\partial x_\alpha} \Big|_\sigma - \frac{\partial}{\partial \sigma} \left( v_\alpha \frac{\partial z}{\partial x_\alpha} \Big|_\sigma \right) + \frac{\partial v_z}{\partial \sigma} \right) d\sigma
\end{aligned} 
$$

The last term on the right-hand side vanishes due to the continuity equation (Eq. (10)), and the expression for the vertical velocity (Eq. (2)) can be substituted in the brackets, leading to the equation below:

$$ 
\begin{aligned}
\int_0^1 -D v_\alpha \frac{\partial P}{\partial x_\alpha} |_z d\sigma + \int_0^1 -v_z \frac{\partial P}{\partial \sigma} d\sigma \\
&= F + \left[ P \left( -v_\sigma - \sigma \frac{\partial \zeta}{\partial t} - (\sigma - 1) \frac{\partial \Delta H}{\partial t} \right) \right]_0^1
\end{aligned} 
\quad (17) $$

As there is no water flux through the top and bottom σ-surfaces ($\sigma = 0$ and $\sigma = 1$), the rate of work resulting from the pressure force can consequently be decomposed into 3 terms:

$$ 
\begin{aligned}
\int_0^1 -D v_\alpha \frac{\partial P}{\partial x_\alpha} |_z d\sigma + \int_0^1 -v_z \frac{\partial P}{\partial \sigma} d\sigma &= F - P(1) \frac{\partial \zeta}{\partial t} - P(0) \frac{\partial \Delta H}{\partial t} \\
&= F - P(0) \frac{\partial \Delta H}{\partial t}
\end{aligned} 
\quad (18) $$

where $(-P(1)\frac{\partial \zeta}{\partial t})$ and $(-P(0)\frac{\partial \Delta H}{\partial t})$ are the work rates of the pressure force, induced by the free-surface and bottom-bathymetry anomalies, respectively. The first term vanishes in the present study because $P(1) = 0$.

The implementation of a time-varying bathymetry also has an impact on the buoyancy flux, obtained as the work rate of the force of gravity:

$$ \phi_z = \int_{-(H+\Delta H)}^{\zeta} -\rho g v_z dz \quad (19) $$

The buoyancy flux is responsible for the transfer of energy between the kinetic and potential energy compartments (Winters et al., 1995). The analytical derivation of this term is rather straightforward, whereas a consistent evaluation of this term in a numerical model is not necessarily as clear. The hydrostatic component of

the vertical pressure force and the force of gravity traditionally cancel out. As a result, the buoyancy force must be recovered in the kinetic energy evolution equation via the horizontal component of the hydrostatic pressure force and the continuity equation (see, for instance, Floor et al. (2011)). In the evolution equation of potential energy, and for a linear incompressible equation of state, this buoyancy flux is recovered exactly (with a minus sign) from the advection of potential temperature and salinity and from the vertical motion of the σ-surfaces associated with the motion of the free surface using the (linear) equation of state. Therefore, these schemes must be consistent with each other to avoid any numerical leakage during transfers between the energy compartments (Marsaleix et al., 2008). This need for consistency is particularly true when these transfers are associated with internal gravity waves and are consequently periodic.

Substituting for the vertical velocity (Eq. (2)), the buoyancy flux (Eq. (19)) can be rewritten:

$$ \phi_z = \int_0^1 -D \rho g \left( v_\sigma + \sigma \frac{\partial \zeta}{\partial t} + (\sigma - 1) \frac{\partial \Delta H}{\partial t} + v_\alpha \frac{\partial z}{\partial x_\alpha} \Big|_\sigma \right) d\sigma \quad (20) $$

The modification of the water column depth enhances the transfer of energy between the kinetic and potential energy compartments through $\int_0^1 -D \rho g (\sigma - 1) \frac{\partial \Delta H}{\partial t} d\sigma$. Integrating by parts, Eq. (20) can be rewritten in terms of the hydrostatic component of the pressure:

$$ 
\begin{aligned}
\int_0^1 -D \rho g \left( \sigma \frac{\partial \zeta}{\partial t} + (\sigma - 1) \frac{\partial \Delta H}{\partial t} \right) d\sigma \\
&= \frac{\partial \zeta}{\partial t} (p(1) - \bar{p}) + \frac{\partial \Delta H}{\partial t} (p(0) - \bar{p})
\end{aligned} 
\quad (21) $$

where $\bar{p} = \int_0^1 p d\sigma$ and $p(1) = 0$. The buoyancy flux can thus be decomposed into several terms, two of which depend explicitly on the time-variation of the free-surface and bottom-bathymetry anomalies:

$$ 
\begin{aligned}
\phi_z &= \int_0^1 -D \rho g v_z d\sigma \\
&= \int_0^1 -D \rho g \left( v_\sigma + v_\alpha \frac{\partial z}{\partial x_\alpha} \Big|_\sigma \right) d\sigma - \frac{\partial \zeta}{\partial t} \bar{p} + \frac{\partial \Delta H}{\partial t} (p(0) - \bar{p})
\end{aligned} 
\quad (22) $$

Relations (18) and (20) show that a moving bottom enhances direct energy exchanges in two ways. A source term $(-P(0)\frac{\partial \Delta H}{\partial t})$ appears in the equation of evolution of the kinetic energy. This term is directly associated with the modification of the water column depth and is proportional to the total (i.e., hydrostatic and non hydrostatic) bottom pressure ($P(0)$). The time-varying bathymetry also enhances energy transfers between the kinetic and potential energy compartments through a contribution in the buoyancy flux: $\int_0^1 -D \rho g (\sigma - 1) \frac{\partial \Delta H}{\partial t} d\sigma$.

Relations (18) and (22) additionally show the symmetry between the transfers associated with the free-surface and bottom-bathymetry anomalies: in both relations, their variation with time induces energy transfers.

Several terms in the expressions of the pressure force (Eq. (18)) and buoyancy flux (Eq. (22)) are not specific to the time-varying bathymetry. The pressure force (Eq. (18)) includes a flux term (F) whose amplitude depends both on the chosen "control volume" over which this term is integrated and on the processes directly induced by the bathymetry motion. The main contribution to the flux (F) is the energy radiated by internal waves (Section 3). For a closed domain, this flux term vanishes.

The form of the first term in the expression of the buoyancy flux (Eq. (22)), $(\int_0^1 -D \rho g (v_\sigma + v_\alpha \frac{\partial z}{\partial x_\alpha} |_\sigma) d\sigma)$, is specific to the choice of coordinate system (namely, to the expression of the coordinate system). In σ-coordinates, the first term in the buoyancy flux Eq. (22) is the sum of the flux caused by the vertical velocity across
the σ-layers ($v_σ$) and by the slope of these σ-layers with respect to the horizontal. This term is associated with the vertical transfer of energy from the time-varying bathymetry to the water column and can also be associated with internal waves. The second term in the expression of the buoyancy flux (Eq. (22)) ($-\frac{\partial x}{\partial t}\bar{p}$) is a direct consequence of the displacement of the free surface. The amplitude of this term strongly depends on the studied configuration as well. In the real ocean, a pressure disturbance generated by a modification of the bathymetry is, for instance, radiated by acoustic waves to the surface, where the resulting anomaly can generate several types of surface waves (see Nosov, 1999). On longer time scales, the response of the free surface can indirectly depend on the amplitude of the forcing because of the time-varying bathymetry and because of the possible occurrence of a hydraulic control (Baines, 2004, Section 2). In the case of the oscillating ridge studied in the following section, the contribution of the free-surface anomalies can be neglected in front of the energy exchanged with the time-varying bottom ($\frac{\partial\Delta H}{\partial t}(p(0) - \bar{p})$).

### 3. Modeling the generation of internal waves by an oscillating ridge

#### 3.1. Oscillating ridge configuration

Internal gravity waves can be generated in a stratified ocean by any vertical displacement of the isopycnal surfaces from equilibrium. In the ocean, such displacements can be induced by the impingement of surface tides on bathymetry variations such as ridges, seamounts or shelf breaks, by wind-induced motions and, more generally, by encountering cross-isopycnal velocities of significant amplitudes. These displacements can be reproduced in stratified water tanks by oscillating objects (Mowbray and Rarity, 1967) or by moving bottoms (Dossmann et al., 2011). A direct numerical simulation of this type of flow can be obtained with the algorithm derived in this study and is now presented.

Internal gravity waves are generated by moving a Gaussian ridge back and forth in a salt-stratified tank in the CNRM-GAME (Meteo-France/CNRS) hydraulic laboratory. These internal waves propagate from the oscillating ridge and reflect at the closed-tank walls. With the exception of the lateral boundaries, the experimental configuration is 2D-vertical and is consequently simulated numerically as a vertical section in the (x-z) plane.

The Gaussian ridge oscillates at a period of $T_0 = 10.05$ s over a flat bottom with depth $H = 394$ mm. For $x \in [-2m, 2m]$,

$$ \Delta H(x; t) = -h_0 \exp \left( - \frac{(x - x_0(t))^2}{L_0^2} \right) \quad (23) $$

where the center of the ridge is located at $x_0(t) = A_0 \sin(\frac{2\pi t}{T_0})$ and where $L_0$ is the horizontal characteristic length scale of the bathymetry. The characteristics of the experimental configuration and the model parameters are given in Table 2.

The tank stratification is obtained experimentally and mimicked numerically by generating salt gradients at a constant temperature. The horizontally and vertically averaged resolutions are approximately 1 mm, allowing for the use of the molecular viscosity and salt diffusivity.

#### 3.2. Corresponding large-scale ocean configuration

In regard to the topography, the aspect ratio of the tank bathymetry is given by $\varepsilon_{b,tank} = 100/394 \sim 0.25$. This ratio corresponds to a ridge height of approximately 1000 m over an average abyssal plain height $H_{ocean} \sim 4000$ m, which is the same order of magnitude as the height of the Pacific ridge in the region of Hawaii. The Pacific ridge (Pinkel and Rudnick, 2006) and the Luzon Strait

<caption>**Table 2**</caption>

Main characteristics of the experimental and numerical internal wave configurations in a closed tank.

<table><tr><td>Size of the stratified tank (mm)</td><td></td><td>(4000 × 500 × 394)</td></tr><tr><td>Gaussian ridge height and e-folding width (mm)</td><td>(h<sub>0</sub>, L<sub>0</sub>)</td><td>(100, 36.86)</td></tr><tr><td>Number of points of model grid in x, and y direction.</td><td>(i<sub>max</sub>, j<sub>max</sub>)</td><td>(4096, 3)</td></tr><tr><td>Horizontal model grid (mm)</td><td>Δx = Δy</td><td>1</td></tr><tr><td>σ-levels</td><td>k<sub>max</sub></td><td>400 equally-spaced levels</td></tr><tr><td>Internal mode time step (ms)</td><td>Δt<sub>int</sub></td><td>0.647</td></tr><tr><td>External mode time step (ms)</td><td>Δt<sub>ext</sub></td><td>0.162</td></tr><tr><td>Brunt-Väisälä frequency (rad s<sup>-1</sup>)</td><td>N</td><td>0.92461</td></tr><tr><td>Equation of state</td><td></td><td>Incompressible and linear</td></tr><tr><td>Viscosity (m<sup>2</sup> s<sup>-1</sup>)</td><td>ν</td><td>10<sup>-6</sup></td></tr><tr><td>Diffusivity (m<sup>2</sup> s<sup>-1</sup>)</td><td>K<sub>D</sub></td><td>10<sup>-9</sup></td></tr><tr><td>Bottom boundary condition</td><td></td><td>No-slip</td></tr><tr><td>Amplitude and period of forcing</td><td>(A<sub>0</sub>, T<sub>0</sub>)</td><td>(1 mm, 10.05 s)</td></tr></table>

(Zhang et al., 2011) are indeed some of the most well-known generators of internal tides in the Pacific Ocean. A characteristic length scale for the tank's ridge is $L \sim 2L_0 = 73$ mm, whereas its motion can be associated with a pulsation $\omega \sim 2\pi/T_0 = 0.6$ rad s<sup>-1</sup> and with a velocity $U \sim L/T_0 = 7$ mm/s.

The dynamical regime is thus characterized by a rather large non-hydrostatic parameter $\delta_{NH} = \omega/N \sim 0.6$. Such a regime can be observed in the ocean as a result of M2 tidal influence and with an average stratification Brunt Väisälä pulsation $N_{ocean} \sim 2.510^{-4}$ rad/s<sup>-1</sup>. The initial density stratification is given by

$$ \rho(z) = \rho_0 \left( 1 - \frac{N^2 z}{g} \right) \quad \text{with } \rho_0 = 1029 \text{ kg/m}^3 $$

The ratio of the internal wave to the surface wave speeds, $\delta_{speed} = NH/\sqrt{gH}$, is approximately equal to 0.014 in the ocean configuration and 0.2 in the tank configuration, indicating in both cases weak non-linear interactions between surface and internal waves (Pedloskly, 2003, p. 85; Floor, 2009, p. 130).

This rapid evaluation of the main characteristics of the dynamical regime shows that the tank configuration is non-hydrostatic and that it can reasonably provide insights into the generation of internal waves by a barotropic tide over an oceanic ridge. This evaluation also illustrates some of the limits of this approach as the ratio of the internal to the surface wave speeds, $\delta_{speed}$, is small but not insignificant. The tank configuration will not, however, provide significant information about the non-linear interactions between the barotropic tides and the internal tides they induce. In a similar configuration, Gerkema and Zimmerman (1995) also showed that the tidal motion over a ridge could be considered as equivalent to a moving ridge if the topography ratio $\varepsilon_b$ and the parameter of nonlinearity $\varepsilon_{NL} = U\varepsilon_b/\omega L$ were both small ($\varepsilon_b, \varepsilon_{NL} \ll 1$). In the present tank configuration, these parameters are $\varepsilon_{b,tank} \sim 0.25$ and $\varepsilon_{NL,tank} \sim 210^{-2}$. Although rather small, $\varepsilon_{NL}$ is lower than 1 but not insignificant, indicating possible discrepancies between internal tides generated by tidal motion compared with the generation of internal waves by the motion of the ridge. Because the tank is static with respect to the surface of the earth, the only rotation to which the tank flow is submitted is the earth's rotation. As a result, its characteristic Rossby number, $R_0 \sim U/(10^{-4}L)$, is several orders of magnitude smaller than typical ocean Rossby numbers at the mid-latitudes. Consequently, none of the rotation-induced characteristics of real-ocean internal tides can be studied in the present tank configuration.
### 3.3. Characteristics of the internal wave

A vertical section of the square Brunt-Väisälä frequency anomaly after ten periods is presented in Fig. 1. Fig. 1a shows the anomaly obtained by Synthetic Schlieren visualization, and Fig. 1b shows the same anomaly simulated numerically. The computed Brunt-Väisälä frequency anomaly is a diagnostic variable obtained from the density field.

In both figures, internal wave rays form a Saint-Andrew's cross, with a characteristic angle $\phi$ equal to the angle predicted by the dispersion equation ($\sin \phi = \omega/N$). The rays reflect at the bottom and at the free surface and weaken during the reflection process. In both fields, the motion of the ridge creates a low-amplitude, permanent circulation, forming a horizontal straight line at the level of the ridge top.

Two primary types of discrepancies can be found between the experiment and the numerical simulation. The first discrepancy concerns the reflection of the ray. The reflection occurs a few centimeters below the free surface in the experiment, whereas the reflection is located right at the free surface in the numerical simulation. This discrepancy can be attributed to the formation of a well-mixed layer along the boundary in the experiment before the ridge motion is started. A first consequence of this discrepancy is that the reflected rays are slightly translated in the experiments. A second consequence is that energy losses are more significant in the experiment because part of the ray energy appears to be trapped in the mixed layers, enhancing the “mixing process” and reducing the amplitude of the reflected ray. Because of the presence of the solid wall at the bottom boundary, Schlieren measurements are available only a few centimeters above the bottom, and

direct comparisons are more difficult. The conditions at the bottom are not identical in both simulations because the entire bottom surface moves with the ridge in the numerical simulation, whereas only the ridge moves in the experiment. However, unlike rays in the surface layer, both the locations of the reflection points and the amplitudes of the reflected rays are quite similar in both simulations. Therefore, no further refinement of the numerical configuration was necessary.

Another difference between the experiment and the numerical simulation concerns the ray itself, which appears much more regular in the numerical simulation than in the experiment. The differences in ray reflection might partially explain such discrepancies. In Fig. 1a, the first upward ray shows two lighter regions corresponding to a higher amplitude of the Brunt-Väisälä frequency, separated by a darker region (approximately halfway along the ray's path). These differences are consistent with a reflection of the ray at different locations, with the lighter regions associated with two reflected rays in the tank. Further comparisons of the ray structure between the experiment and the simulation can then be made by plotting the evolution of the Brunt-Väisälä frequency anomaly along segments transverse to the rays. These segments are shown in Fig. 2 at four locations along the first upward ray. A global comparison is satisfactory for the amplitude of the ray (approximately 0.05 rad$^2$/s$^2$) and for its e-folding width (approximately 6 cm), even if in numerical simulations rays have slightly smaller amplitudes and are less focused than in the tank experiment. Globally, the amplitude decreases and the thickness increases along the ray, moving away from the ridge (from $S_1$ to $S_4$). At the third location ($S_3$ and $S'_3$), ray shapes are rather similar. The largest differences are found at the same locations: near and far from the ridge (i.e., in Segments $S_1-S'_1$ and $S_4-S'_4$, respectively). We note in particular that at the fourth location, the amplitude of the ray significantly increases in the experiment while it continues to decrease in the model field. A possible explanation, already proposed for the corresponding differences in Fig. 1, is the presence of a reflected ray in the experiment.

In regard to the occurrence of a mixed layer along the free surface, the objective of the present simulation was to study the generation of internal rays by the oscillations of a Gaussian-shaped ridge in a linearly stratified fluid, rather than the study of reflections or of the general properties of the internal tide field away from the ridge. Namely, no specific effort has been made to initialize the numerical simulation with a salinity field that corresponds to that of the experiment.

<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Vertical section of square Brunt-Väisälä frequency anomaly (rad$^2$s$^{-2}$) after 10 forcing periods (a) measured by Synthetic Schlieren, (b) simulated numerically. Black line segments across internal wave rays indicate positions of cross radius Brunt-Väisälä frequency evolution shown in Fig. 2. Segments are labeled ($S_1$) to ($S_4$) for Synthetic Schlieren and ($S'_1$) to ($S'_4$) for numerical simulation.</figcaption></figure>

### 3.4. The internal wave generation process

A comparison between experimentally and numerically simulated internal waves shows close correspondence between the two types of waves. Both types of simulations lead to internal waves, but as shown in Section (3.2), the generation mechanism of these simulated waves is different from the generation of internal waves (the so-called internal tides) by tidal waves over a ridge (Floor et al., 2011).

All of the energy fluxes considered in the present section are integrated in space over the computational domain, and their running mean over a forcing period (written $\langle\rangle_{T_0}$) is plotted in Fig. 3. Fig. 4 gives a schematic representation of the corresponding exchanges. The density profile is given in Section 3.2 and is chosen so that at the surface, $\rho(1) = \rho_0$.

The running mean of the work rate of the time-dependent bathymetry against the bottom hydrostatic pressure (relation 18, triangles in Fig. 3a) initially leads to an increase in kinetic energy: $\langle -p(0) \frac{\partial \bar{H}}{\partial t} \rangle_{T_0} > 0$. The averaged non-hydrostatic contribution to this work rate (squares) is initially negative, indicating a net loss
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Evolution of Brunt–Väisälä frequency along black segments shown in Fig. 1, Left: measured by Synthetic Schlieren in Segments (S₁) to (S₄), Right: simulated numerically with SNH model in Segments (S₁) to (S₄). Curves are plotted every second for $t \in [10T_0, 11T_0]$.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Running mean of the various transfers of energy ($W/\rho_0$) over a forcing period ($T_0$). All energy transfers are integrated over the computational domain and normalized by $\rho_0$. (a) Transfers associated with the motion of the bottom, black line: «$\phi_z$» = $\langle\phi_z\rangle_{T_0} = \langle\int_0^1 -D\rho g v_z d\sigma\rangle_{T_0}$, circles: «$\phi_z$: bottom» = $\langle\int_0^1 -\frac{\partial\Delta H}{\partial t}(p(0)-\bar{p})\rangle_{T_0}$, dashed line: «$\phi_z$: hydro» = $\langle\int_0^1 -\frac{\partial\Delta H}{\partial t}p(0)\rangle_{T_0}$, dashed-dotted line: «$\phi_z$: mean» = $\langle-\frac{\partial\Delta H}{\partial t}\bar{p}\rangle_{T_0}$, plus signs: «PG: bottom» = $\langle-P(0)\frac{\partial\Delta H}{\partial t}\rangle_{T_0}$, triangles: «PG: hydro» = $\langle-p(0)\frac{\partial\Delta H}{\partial t}\rangle_{T_0}$, squares: «PG: NH» = $\langle-q(0)\frac{\partial\Delta H}{\partial t}\rangle_{T_0}$, crosses: «KE: bottom» = $\langle-\frac{\partial\Delta H}{\partial t}(q(0)+\bar{p})\rangle_{T_0}$, (b): Several contributions to the buoyancy flux, dashed line: «$\phi_z$» = $\langle\phi_z\rangle_{T_0} = \langle\int_0^1 -D\rho g v_z d\sigma\rangle_{T_0}$, circles: «$\phi_z$: $\Delta z/\Delta t$» = $\langle-\frac{\partial\Delta H}{\partial t}\bar{p} + \frac{\partial\Delta H}{\partial t}(p(0)-\bar{p})\rangle_{T_0}$, dashed-dotted line: «$\phi_z$: $u.\Delta z/\Delta x + \omega$» = $\langle\int_0^1 -D\rho g(v_\sigma + v_\zeta \frac{\partial z}{\partial x}|_\sigma)d\sigma\rangle_{T_0}$, black line: «$\phi_z$: $\omega$» = $\langle\int_0^1 -D\rho g v_\sigma d\sigma\rangle_{T_0}$, plus signs: «$\phi_z$: $u.\Delta z/\Delta x$» = $\langle\int_0^1 -D\rho g v_\zeta \frac{\partial z}{\partial x}|_\sigma d\sigma\rangle_{T_0}$.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Fig. 4.** Schematic representation of the running-mean of the energy fluxes associated with the motion of the bottom bathymetry over one forcing period ($T_0$) (solid arrow). Both kinetic energy (KE) and potential energy (PE) compartments are represented. The running mean of the internal wave radiative flux ($\langle \bar{F} \rangle_{T_0} = \langle \int_0^1 -\frac{\partial \bar{P}_x}{\partial x_x} \partial_z d\sigma \rangle_{T_0}$) is represented with a dashed arrow (this flux vanishes when integrated over the entire tank).</figcaption></figure>

of kinetic energy: $\langle -q(0) \frac{\partial H}{\partial t} \rangle_{T_0} < 0$. However, the overall pressure anomaly contribution remains initially positive: $\langle -P(0) \frac{\partial H}{\partial t} \rangle_{T_0} > 0$. The averaged contribution of the bottom bathymetry to the buoyancy flux $\langle \frac{\partial H}{\partial t} (p(0) - \bar{p}) \rangle_{T_0}$ (relation 22, circles) remains negative, indicating that kinetic energy is lost while potential energy is enhanced. The total averaged buoyancy flux (black line) remains positive.

The sum of the running mean of the contributions of the moving bottom (buoyancy flux plus pressure gradient force) and the kinetic energy balance (crosses) is eventually negative: overall, kinetic energy is lost through direct motion of the bathymetry, while potential energy is enhanced. The overall loss of kinetic energy is, however, much smaller than the gain through the total buoyancy flux (black line), and kinetic energy is thus enhanced overall as a result of the motion of the bottom “ridge”.

To further study the evolution of the buoyancy flux, Fig. 3b details several contributions to the running mean of the buoyancy flux, based on Eq. (22). The running mean of the buoyancy flux over one forcing period $T_0$ is recalled (black curve in Fig. 3b). The total contribution of the time variations of the $\sigma$-layer (circles) is rather small compared with the total buoyancy flux and is mostly induced by the time variations of the bathymetry: $|\langle -\frac{\partial \bar{p}}{\partial t} \rangle_{T_0}| \ll |\langle \frac{\partial H}{\partial t} (p(0) - \bar{p}) \rangle_{T_0}|$. The barotropic potential energy associated with the free-surface anomaly and the associated energy transfers can be neglected in the present configuration providing that the maximal amplitude of the free-surface anomaly never exceeds $10^{-5}$ m in the numerical configuration. This maximum amplitude is confirmed experimentally as no surface wave or free-surface anomalies have been observed in the CNRM-GAME tank. These very small free surface anomalies lead to negligible barotropic energy transfers.

The remaining portion of the buoyancy flux is approximately equal to the baroclinic contribution ($\phi'_z$) defined by Floor et al. (2011) (plus-signs in Fig. 3b):

$$ \left\langle \int_0^1 -D\rho g \left( v_\sigma + v_\alpha \frac{\partial z}{\partial x_\alpha} \right) d\sigma \right\rangle_{T_0} \approx \langle \phi'_z \rangle_{T_0} \\ = \left\langle \int_0^1 -D\rho g \left( v_\sigma + v'_\alpha \frac{\partial z}{\partial x_\alpha} \right) d\sigma \right\rangle_{T_0} \quad (24) $$

where $v'_\alpha = v_\sigma - \bar{v}_\alpha$.

The individual contribution of the cross-$\sigma$ vertical velocity $\langle \int_0^1 -D\rho g v_\sigma d\sigma \rangle_{T_0}$ (black curve) and of the sloping of the $\sigma$-layers $\langle \int_0^1 -D\rho g v_\alpha \frac{\partial z}{\partial x_\alpha} d\sigma \rangle_{T_0}$ (dashed curve) are approximately anti-correlated.

The remaining “barotropic” contribution ($\langle \bar{\phi}_z \rangle$, see Eq. (3.2) of Floor et al., 2011) is approximately equal to the contribution of the time-varying bathymetry (red curves):

$$ \left\langle \frac{\partial \Delta H}{\partial t} (p(0) - \bar{p}) \right\rangle_{T_0} \approx \langle \bar{\phi}_z \rangle_{T_0} \equiv \left\langle \int_0^1 -D\rho g \left[ \sigma \left( \frac{\partial z_\zeta}{\partial t} + \bar{v}_\alpha \frac{\partial z_\zeta}{\partial x} \right) + (\sigma - 1) \left( \frac{\partial \Delta H}{\partial t} + \bar{v}_\alpha \frac{\partial (H + \Delta H)}{\partial x} \right) \right] d\sigma \right\rangle_{T_0} \quad (25) $$

This contribution is smaller than the baroclinic component $\langle \phi'_z \rangle_{T_0}$.

Oscillations at half the forcing period can be observed over both the buoyancy flux and its “baroclinic” component. Period halving is a result of the correlated variations of the velocity and density anomalies in ($\phi'_z$). The oscillations of the “baroclinic component” at half the forcing period confirm the generation of internal waves by the ridge. After 4 to 5 forcing periods, internal waves have propagated over a large portion of the tank, allowing for their contribution to the overall buoyancy flux to become significant. The energy flux associated with these waves reaches a quasi-stationary state after 10 to 11 forcing periods (not shown) and, on average, potential energy is transformed into kinetic energy because the baroclinic contribution to the buoyancy flux remains positive.

As a conclusion (Fig. 4), the running mean of the buoyancy flux over one forcing period can be written approximately as follows:

$$ \langle \phi_z \rangle_{T_0} = \langle \phi'_z \rangle_{T_0} + \langle \bar{\phi}_z \rangle_{T_0} \\ \approx \left\langle \int_0^1 -D\rho g \left( v_\sigma + v'_\alpha \frac{\partial z}{\partial x_\alpha} \right) d\sigma \right\rangle_{T_0} \\ + \left\langle \frac{\partial \Delta H}{\partial t} (p(0) - \bar{p}) \right\rangle_{T_0} \quad (26) $$

The motion of the bathymetry initiates the generation of internal waves by initially enhancing kinetic energy through the varying-bathymetry contribution to the work rate of the pressure force ($\langle -P(0) \frac{\partial H}{\partial t} \rangle_{T_0}$). Kinetic energy is then reduced, and potential energy is enhanced through the contribution of the time-varying
bathymetry to the “barotropic component of the buoyancy flux” ($\bar{\phi}_z$). These results are rather similar to those obtained by Floor et al. (2011) for the generation of internal waves by barotropic tides over a ridge. In the configuration studied by Floor et al. (2011), energy is directly supplied by tides, i.e., by long surface waves and isopycnal surfaces are displaced by tidally induced diapycnal currents. Floor et al. (2011) further showed that potential energy is then enhanced by the barotropic component of the buoyancy flux (Eq. (25)), whereas in the present configuration, surface waves play no significant role. In both cases (time-varying versus constant bathymetry), potential energy is eventually reduced on average by the propagation of internal waves through the “baroclinic component of the buoyancy flux” ($\phi_z$).

## 4. Discussion and conclusion

A simple implementation of a time-varying bathymetry is proposed for a $\sigma$-coordinate free-surface ocean model, offering new perspectives for the modeling of ocean configurations presenting topography variations. We have shown, in particular, that a careful choice of model variables can notably reduce the modifications to the numerical model algorithm. Considering the total change in the water column depth ($\eta = \zeta + \Delta H$), the numerical equations associated with volume conservation under the Boussinesq approximation in particular remain advantageously unchanged. The resulting implementation of time-varying bathymetry requires only a few additional operations and has thus (nearly) no effect on the global numerical cost.

The new time-varying bathymetry ocean model has been compared to a laboratory configuration of internal tide generation by a moving ridge in a CNRM-GAME tank. The resulting internal tide rays have been shown to be in close agreement with Synthetic Schlieren measurements. The “new” forcing caused by the variation of the bathymetry has been characterized by deriving the mechanical energy equations and numerically computing the resulting mechanical energy transfers associated with the pressure force and buoyancy flux. The proposed mechanism for the energy cascading leading to the generation of internal waves by a moving “ridge” in a stratified tank is found to be similar to the one proposed by Floor et al. (2011) for the generation of internal tides.

## Acknowledgments

We gratefully thank the computer team of the Laboratoire d'Aérologie for "NUWA", our laboratory cluster, together with CALMIP (HYPERION) and CINES (JADE) for computer resources. The present study has been supported by the ANR projects PIWO (contract number ANR-08-BLAN-0113) and COMODO (ANR-11-MONU-005) and the LEFE-IDAO "Ondes internes" project.

## References

Auclair, F., Estournel, C., Floor, J., Herrmann, M., Nguyen, C., Marsaleix, P., 2011. A non-hydrostatic algorithm for free-surface ocean modelling. Ocean Modell. 36, 49-70.

Baines P.G., 2004. Topographic effects in stratified flows. Cambridge Monogrphs on Mechanics, p. 482.

Dossmann, Y., Paci, A., Auclair, F., Floor, J.W., 2011. Simultaneous velocity and density measurements for an energy-based approach to internal waves generated over a ridge. Exp. Fluids 51 (4), 1013-1028.

Floor J.W., 2009. Transferts d'Energie associés à la Génération, Propagation et Dissipation de la Marée Interne (Ph.D. thesis). Université de Toulouse Paul Sabatier - Toulouse III.

Floor, J.W., Auclair, F., Marsaleix, P., 2011. Energy transfers in internal tide generation, propagation and dissipation in the deep ocean. Ocean Modell. 38 (1-2), 22-40.

Fratepietro, F., Baker, T.F., Williams, S.D.P., Van Camp, M., 2006. Ocean loading deformations caused by storm surges on the northwest European shelf. Geophys. Res. Lett. 33, L06317.

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

Griffies, S.M., Schmidt, M., Herzfeld, M., 2010. Element of MOM4. GFDL Ocean Group Tenical, Report #6, 444 pages.

Heinrich, P., Guibourg, S., Mangeney, A., Roche, R., 1999. Numerical modeling of a landslide generated tsunami following a potential explosion of the Montserrat volcano. Phys. Chem. Earth (A) 24 (2), 163-168.

Le Bars, Y., Lyard, F., Jeandel, C., Dardengo, L., 2010. The AMANDES tidal model for the Amazon estuary and shelf. Ocean Modell. 31, 132-149.

Le Hir, P., Cayocca, F., Waeles, B., 2011. Dynamics of sand and mud mixtures: a multiprocess-based modelling strategy. Cont. Shelf Res. 31, S135-S149.

Le Provost, C., Lyard, F., Molines, J.M., Genco, M.L., Rabilloud, F., 1998. A hydrodynamic ocean tide model improved by assimilating a satellite altimeter-derived data set. J. Geophys. Res. 103 (C3), 5513-5529.

Losch, M., Adcroft, A., Campin, J.M., 2004. How sensitive are coarse general circulation models to fundamental approximations in the equations of motion? J. Phys. Oceanogr. 34, 306-319.

Marsaleix, P., Auclair, F., Floor, J.W., Herrmann, M.J., Estournel, C., Pairaud, I., Ulses, C., 2008. Energy conservation issues in sigma-coordinate free-surface ocean models. Ocean Modell. 20 (1), 61-89.

Marshall, J., Adcroft, A., Campin, J.M., Hill, C., 2004. Atmospheric-ocean modelling exploiting fluid isomorphisms. Monthly Weather Rev. 132, 2882-2894.

Mowbray, D.E., Rarity, B.S.H., 1967. A theoretical and experimental investigation of the phase configuration of internal waves of small amplitude in a density stratified liquid. J. Fluid Mech. 28 (01), 1-16.

Nosov, M.A., 1999. Tsunami generation in compressible ocean. Phys. Chem. Earth (B) 24 (5), 437-441.

Pavlov, V.I., Tromp, J., Tito, E.P., 2009. Tsunami generation by ocean floor rupture front propagation: Hamiltonian description. Nat. Hazards Earth Syst. Sci. 9, 217-227.

Pedlosky, J., 2003. Waves in the Ocean and Atmosphere. Introduction to Waves Dynamics. Springer-Verlag, 260 pages.

Pinkel, R., Rudnick, D., 2006. Hawaii ocean mixing experiment (home), editorial. J. Phys. Oceanogr. 36, 965-966.

Weiss, R., Wünnemann, K., Bahlburg, H., 2006. Numerical modelling of generation, propagation and run-up of tsunamis caused by oceanic impacts: model strategy and technical solutions. Geophys. J. Int. 167, 77-88.

Williams, J.J., MacDonald, N.J., O'Connor, B.A., Pan, S., 2000. Offshore sand bank dynamics. J. Mar. Syst. 24, 153-173.

Winters, K.B., Lombard, P.N., Riley, J.J., D'Asaro, E.A., 1995. Available potential energy and mixing in density-stratified fluids. J. Fluid Mech. 289, 115-128.

Zhang, Z., Fringer, O.B., Ramp, S.R., 2011. Three-dimensional, nonhydrostatic numerical simulation of nonlinear internal wave generation and propagation in the south china sea. J. Geophys. Res. 116, C05022.