Hamiltonian model for coupled surface and internal waves in the presence of currents

TL;DR

This paper develops a Hamiltonian framework to analyze coupled surface and internal water waves influenced by depth-dependent currents, revealing that wave-current interactions are primarily affected by currents near the surface and interface.

Contribution

It introduces a Hamiltonian model for coupled surface and internal waves with currents, providing new insights into wave-current interactions in stratified fluids.

Findings

Wave-current interaction depends mainly on currents near the surface and interface.

Hamiltonian equations of motion are derived for the system.

Small amplitude and long-wave approximations are discussed.

Abstract

We examine a two dimensional fluid system consisting of a lower medium bounded underneath by a flatbed and an upper medium with a free surface. The two media are separated by a free common interface. The gravity driven surface and internal water waves (at the common interface between the media) in the presence of a depth-dependent current are studied under certain physical assumptions. Both media are considered incompressible and with prescribed vorticities. Using the Hamiltonian approach the Hamiltonian of the system is constructed in terms of 'wave' variables and the equations of motion are calculated. The resultant equations of motion are then analysed to show that wave-current interaction is influenced only by the current profile in the 'strips' adjacent to the surface and the interface. Small amplitude and long-wave approximations are also presented.

Equations

$$\Omega (\eta )=\{(x,y)\in {\mathbb{R}}^2:-h<y<\eta (x,t)\}.$$

$${\Omega }_1(\eta ,{\eta }_1)=\{(x,y)\in {\mathbb{R}}^2:\eta (x,t)<y<h_1+{\eta }_1(x,t)\}.$$

$$U(y)=\left\{\begin{array}[]{lcl}\sigma_{4}=\gamma_{1}h_{1}+\kappa_{1}<0\qquad y=h_{1},\\ \gamma_{1}y+\kappa_{1},\qquad m_{1}\leq y\qquad\mbox{ (layer V),}\\ \sigma_{3}=\gamma_{1}m_{1}+\kappa_{1}\qquad y=m_{1},\\ \sigma_{2}=\gamma l_{1}+\kappa\qquad y=l_{1},\\ \gamma y+\kappa,\qquad-l\leq y\leq l_{1}\qquad\mbox{ (layers II and III),}\\ \sigma_{1}=-\gamma l+\kappa\qquad y=-l,\\ 0\qquad y=-h\mbox{ (flatbed),}\end{array}\right.$$

$${\text{V}}^{\prime }:m_1\le y\le h_1+{\eta }_1(x,t),$$

$$-l\le \eta (x,t)\le l_1.$$

$$\begin{array}{rl}{\text{II}}^{\prime }:& -l\le y\le \eta (x,t),\\ {\text{III}}^{\prime }:& \eta (x,t)\le y\le l_1.\end{array}$$

$$\left\{\begin{array}[]{lcl}u_{1}=\varphi_{1,x}+\gamma_{1}y+\kappa_{1}\qquad\mbox{ (layer V${}^{\prime}$),}\\ u_{1}={\varphi}_{1,x}+\gamma y+\kappa\qquad\mbox{ (layer III${}^{\prime}$),}\\ v_{1}={\varphi}_{1,y},\\ u={\varphi}_{x}+\gamma y+\kappa\qquad\mbox{ (layer II${}^{\prime}$),}\\ v={\varphi}_{y}.\end{array}\right.$$

$$\rho >{\rho }_1.$$

$$\gamma =\frac{{\sigma }_2-{\sigma }_1}{l+l_1}.$$

$$\underset{|x|\to \infty }{\lim }\eta (x,t)=0,\underset{|x|\to \infty }{\lim }{\eta }_1(x,t)=0,$$

$$\underset{|x|\to \infty }{\lim }\phi (x,y,t)=0,\underset{|x|\to \infty }{\lim }{\phi }_1(x,y,t)=0.$$

$$-l\le \eta (x,t)<l_1,m_1\le h_1+{\eta }_1(x,t)$$

$$\bm{u}(x,y,t):=\left\{\begin{array}[]{ccc}u(x,y,t),&{\rm in}&\Omega(\eta),\\ u_{1}(x,y,t),&{\rm in}&\Omega_{1}(\eta,\eta_{1}),\end{array}\right.$$

$$\bm{v}(x,y,t):=\left\{\begin{array}[]{ccc}v(x,y,t),&{\rm in}&\Omega(\eta),\\ v_{1}(x,y,t),&{\rm in}&\Omega_{1}(\eta,\eta_{1}).\end{array}\right.$$

$${\mathbf{u}}_x+{\mathbf{v}}_y=0\text{in}\Omega \cup {\Omega }_1.$$

$$\bm{\psi}(t,x,y)=\left\{\begin{array}[]{ccc}\psi(t,x,y)&{\rm in}&\Omega,\\ \psi_{1}(t,x,y)&{\rm in}&\Omega_{1},\end{array}\right.$$

$$\left\{\begin{array}[]{lll}u=\psi_{y},&v=-\psi_{x},&{\rm in}\quad\Omega,\\ u_{1}=\psi_{1,y},&v_{1}=-\psi_{1,x},&{\rm in}\quad\Omega_{1}.\end{array}\right.$$

$$\bm{\varphi}(x,y,t):=\left\{\begin{array}[]{ccc}\varphi(x,y,t),&{\rm in}&\Omega(\eta),\\ \varphi_{1}(x,y,t),&{\rm in}&\Omega_{1}(\eta,\eta_{1}),\end{array}\right.$$

$$\mathbf{u}={\mathbf{\phi }}_x+U(y).$$

$$\mathbf{F}=2\omega \nabla \mathbf{\psi }=2\omega (-\mathbf{v},\mathbf{u}{)}^T$$

$$\left\{\begin{array}[]{lcl}\bm{u}_{t}+\bm{u}\bm{u}_{x}+\bm{v}\bm{u}_{y}+2\omega\bm{v}&=&-\frac{1}{\bm{\rho}^{*}}P_{x},\\ \bm{v}_{t}+\bm{u}\bm{v}_{x}+\bm{v}\bm{v}_{y}-2\omega\bm{u}&=&-\frac{1}{\bm{\rho}^{*}}P_{y}-g,\end{array}\right.$$

$$P=P_{atm}\text{o}ny={\eta }_1(x)+h_1,$$

$$v_1={\eta }_{1,t}+u_1{\eta }_{1,x}\mathrm{on}y={\eta }_1(x,t)+h_1,$$

$$\begin{array}[]{lll}v_{1}=\eta_{t}+u_{1}\eta_{x}&{\rm on}&y=\eta(x,t),\\ v=\eta_{t}+u\eta_{x}&{\rm on}&y=\eta(x,t),\end{array}$$

$$v=0\mathrm{on}y=-h.$$

$$\psi (t,x,\eta (x,t))={\psi }_1(t,x,\eta (x,t))=-{\int }_{-\infty }^x{\eta }_t(x^{\prime },t)d{x}^{\prime },$$

$$\chi (x,t)\equiv \psi (t,x,\eta (x,t))={\psi }_1(t,x,\eta (x,t)).$$

$${\chi }_1(x,t)\equiv {\psi }_1(t,x,{\eta }_1(x,t)+h_1)=-{\int }_{-\infty }^x{\eta }_{1,t}(x^{\prime },t)d{x}^{\prime }.$$

$${\eta }_{1,t}=({\phi }_{1,y}{)}_{s_1}-{\eta }_{1,x}[({\phi }_{1,x}{)}_{s_1}+{\gamma }_1(h_1+{\eta }_1)+{\kappa }_1]$$

$$\begin{array}[]{c}\eta_{t}=(\varphi_{1,y})_{s}-\eta_{x}[(\varphi_{1,x})_{s}+\gamma\eta+\kappa],\\ \eta_{t}=(\varphi_{y})_{s}-\eta_{x}[(\varphi_{x})_{s}+\gamma\eta+\kappa],\end{array}$$

Full text

††journal: Nonlinear Analysis: Real World Applications

1 Introduction

It has been known for many centuries that the ocean contains currents that flow along generally consistent paths. The Spanish galleons transporting gold and silver from Mexico to Spain made use of the Gulf Stream to help them return home. Since then, scientists have gained much more information on both where currents flow and why. In the oceans currents very often exist with undercurrents. The first undercurrent was discovered in 1951 by Townsend Cromwell who was investigating fishing techniques in the central Pacific Ocean. Undercurrents have since been found under most major currents. The equatorial region in the Pacific is characterised by a thin shallow layer of warm and less dense water over a much deeper layer of cold denser water. The two layers are separated by a sharp thermocline (where the temperature gradient has a maximum, it is very close to the pycnocline, where the pressure gradient has a maximum) at a depth, depending on the location, but usually at 100 – 200 m beneath the surface. For modelling purposes both layers are assumed homogeneous with a sharp boundary at the thermocline/pycnocline (see [36]).

The Equatorial Undercurrent (EUC) flows in a region that is roughly within 200 - 300 km (below 3∘ latitude) of the Equator, it is symmetric about the Equator and extends nearly across the whole length (more than 12000 km) of the Pacific Ocean basin [34]. With speeds in excess of 1 m/s, the EUC is one of the fastest permanent currents in the world.

The flow has nearly two-dimensional character, with small meridional variations. While at depths in excess of about 240 m there is, essentially, an abyssal layer of still water, the ocean dynamics near the surface is quite complex. In this region the wave motion typically comprises surface gravity waves with amplitudes of 1-2 m and oscillations with an amplitude of 10-20 m at the thermocline (of mean depth between 50 m and 150 m). These waves interact with the underlying currents. In that case the velocity is (anti-) parallel to the Earth’s angular speed ${\bf\omega},$ so their vector product is zero. This feature distinguishes the dynamics of the equatorial zone from the ocean dynamics at higher latitudes.

The strong stratification confines the wind-driven currents to a shallow near-surface region, less than 200 m deep. In the Atlantic and Pacific, the westward trade winds induce a westward surface flow at speeds of 25-75 cm/s, while a jet-like current – the Equatorial Undercurrent (EUC) – flows below it toward the East (counter to the surface current), attaining speeds of more than 1 m/s at a depth of nearly 100 m. The wind-generated equatorial current in the layer above the thermocline is with a strictly monotonic depth-dependence and exhibits flow-reversal, while beneath the thermocline the current simply decays with increasing depth, being irrelevant in the abyssal region.

While viscous theory is essential in explaining the generation of the equatorial current induced by wind forcing, inviscid theory is adequate for the study of non-turbulent wave-current interactions since the relevant Reynolds numbers are very large (see [43]).

For some general facts concerning the description of waves interacting with currents we refer to the following reviews and monographs [10, 49, 37, 52] and the references therein. The present study draws from previous single medium irrotational [54], [3], [46], [47], [48] and rotational [9], [11], [10], [12], [50], [17], [53], [42] studies as well as from studies of two-media systems such as [1], [2], [22], [21], [18], [19], [16], [15], [4], [5],[6], [7], [20],[27],[28], [29], [41], [44], [45].

The Hamiltonian approach to water waves dynamics has been put forward for the first time by Zakharov [54]. The Hamiltonian formulation describing the two-dimensional nonlinear interaction between coupled surface waves, internal waves, and an underlying current with piecewise constant vorticity, in a two-media fluid overlying a flat bed has been developed in [16], [15]. In the present study we will be following a similar approach, taking into account the shear current structure suggested in [18]. Related results for a flat surface (effectively rigid lid) has been studied in [4] [5],[6], [7], [8].

The model equations will be presented in a canonical Hamiltonian form and then small amplitude and long wave approximations will be derived.

2 Preliminaries

The system under study involves two-dimensional surface and internal gravity water waves and a depth dependent current as per Figure 1.

The medium underneath the internal wave is defined by the domain

[TABLE]

This medium is bounded at the bottom by an impermeable flatbed at a depth $-h$. The medium above the internal wave $y=\eta(x,t)$ is the domain

[TABLE]

This medium is regarded as being bounded on top by a surface wave at a hight $y=h_{1}+\eta_{1}(x,t)$ moving around the average level $y=h_{1}.$ Throughout the article the subscript $1$ will be used to mean evaluation for the upper medium $\Omega_{1}$, and no subscript means evaluation for the lower medium $\Omega$. Subscript $s$ will be used to denote evaluation at the common interface (thermocline/pycnocline), subscript $s_{1}$ - evaluation at the free surface.

The velocity field is denoted by ${\bf{V}}(x,y,z)=(u,v,0)$ in $\Omega$ and ${\bf{V}}_{1}(x,y,z)=(u_{1},v_{1},0)$ in $\Omega_{1}$. The function $\eta(x,t)$ describes the deviation of the internal wave from its average level $y=0$, i.e. $\int\limits_{\mathbb{R}}\eta(x,t)dx=0.$ Similarly, we define the mean of $\eta_{1}+h_{1}$ to be the unperturbed surface $y=h_{1},$ i.e. $\int\limits_{\mathbb{R}}\eta_{1}(x,t)dx=0.$

A depth dependent current $U(y)$ has the following structure:

[TABLE]

for constants $\sigma_{1}$, $\sigma_{2}$, $\sigma_{3}$, $\sigma_{4}<0,$ $\kappa$, $\kappa_{1},$ $l,$ $l_{1}$, $\gamma$ and $\gamma_{1}<0$, where $\kappa$ is the time-independent current velocity at $y=0;$ $\gamma$ and $\gamma_{1}$ are the non-zero constant vorticity for layers II, III and V, noting that the current is not specified explicitly in layers I and IV, and satisfies only the condition that in the whole fluid body $U(y)$ is a continuous function. This setup of the wave-current system is motivated in [18], see also [6, 7]. The wave motion on the surface is usually confined in the top layer V′, defined for depths $y$ such that

[TABLE]

and the internal wave is confined in the strip formed of layers II and III:

[TABLE]

We define also layers II′ and III′:

[TABLE]

We consider a velocity field, which is defined by the wave-related velocity potentials $\varphi(x,y,t)$ for the domain $\Omega$ and $\varphi_{1}(x,y,t)$ for the domain $\Omega_{1}$ as follows:

[TABLE]

We note that this representation separates the wave and current contributions to the velocity in layers II′, III′ and V′, so that the horizontal velocity field in $\Omega$ is nominally separated to a wave and current part, i.e. $u={{\varphi}}_{x}+U(y)$, and similarly $u_{1}={{\varphi}}_{1,x}+U(y)$ for $\Omega_{1},$ see also [15, 16].

The respective constant densities $\rho$ and $\rho_{1}$ of the lower and upper media satisfy the stability condition for immiscibility

[TABLE]

The rotationality of the layers II and III is given by constant vorticity

[TABLE]

The sign of $\gamma$ is not specified and the case $\gamma=0$ is of course possible.

We assume that for large $|x|$ the amplitude of $\eta$ attenuates and hence make the following assumptions

[TABLE]

[TABLE]

Moreover, we assume

[TABLE]

for all $x$ and $t$. This is an assumption that the surface wave motion takes place only in the strip of layer V′, and the internal wave motion is only in the strip of layers II and III.

3 Governing equations

The governing equations for an inviscid flow are the Euler’s equations modified with terms, taking into account the Coriolis force.

Denoting with $(\bm{u}(x,y,t),\bm{v}(x,y,t))$ the velocity field in $\Omega\cup\Omega_{1}$, the horizontal component of the velocity field is

[TABLE]

and the vertical component of the velocity field is

[TABLE]

The mass conservation is given by

[TABLE]

The equation of mass conservation (11) ensures the existence of a stream function

[TABLE]

determined up to an additive term that depends only on time, by

[TABLE]

For convenience we introduce also

[TABLE]

then

[TABLE]

Considering equatorial motion the following Coriolis forces per unit mass have to be taken into account:

[TABLE]

with $\omega$ being the rotational speed of the Earth. Then the Euler’s equations are

[TABLE]

where $P=P(x,y,t)$ denotes the pressure, $g$ is the gravitational acceleration and $\bm{\rho}^{*}$ denotes the density of the fluid which, in our case is assumed to be piecewise constant, equal to $\rho$ in $\Omega(\eta)$ and to $\rho_{1}$ in $\Omega_{1}(\eta,\eta_{1}).$

Complementing the equations of motion are the boundary conditions, of which

[TABLE]

(with $P_{atm}$ being the constant atmospheric pressure) decouples the motion of the water from that of the air. In addition to (16) we have the kinematic boundary conditions which refer to the flat bed, the interface $y=\eta(x,t)$, the free surface $y=h_{1}+\eta_{1}(x,t)$ and reflect the impermeability of these three surfaces. Thus, they read as

[TABLE]

[TABLE]

and

[TABLE]

From the kinematic boundary conditions (18) one can obtain

[TABLE]

i.e. $\bm{\psi}$ is a continuous function across $y=\eta.$ This also implies that the normal velocity field components are continuous across the interface $y=\eta(x,t)$.

We introduce for convenience

[TABLE]

Similarly

[TABLE]

With the velocity decompositions (3) the kinematic boundary conditions (17) and (18) can now be written as

[TABLE]

and respectively, as

[TABLE]

where the subscript $s_{1}$ means that we look at traces of the involved functions on the free surface $y=\eta_{1}(x,t)+h_{1}$, while the subscript $s$ denotes traces on the interface $y=\eta(x,t)$. The following notation will be used later on in the paper. Namely, we set

[TABLE]

The Dirichlet-Neumann operator $G(\eta)$ associated to the layer $\Omega(\eta)$ is defined as

[TABLE]

where $n$ denotes the unit outward normal vector to the layer $\Omega(\eta)$ along the interface $y=\eta(x)$. Recall that $\varphi_{1}$ is the solution of the boundary value problem

[TABLE]

The Dirichlet-Neumann operator $G_{1}(\eta,\eta_{1})$ associated to the upper media $\Omega_{1}(\eta,\eta_{1})$ is defined through

[TABLE]

where, we denote by $n_{1}$ the unit outward normal vector to $\Omega_{1}(\eta,\eta_{1})$ along the free surface $y=\eta_{1}(x)+h_{1}$. Of course, $G_{1}(\eta,\eta_{1})$ is a matrix-operator, for which we choose the notation

[TABLE]

Euler’s equations can be recast by means of the stream function and of the generalized velocity potential in the given layers in the form of Bernoulli conservation laws as follows:

[TABLE]

where $\tilde{f}_{1}(t)$ is an arbitrary function of $t$. This is related to the freedom to change $\varphi_{1}$, if necessary, by a time-dependent factor. Making use of (16) we can absorb the constant $P_{atm}$ in the arbitrary function and obtain on $y=h_{1}+\eta_{1}(x,t)$

[TABLE]

for another function

[TABLE]

Similarly, on the interface $y=\eta(x,t)$ the pressure is continuous function and thus it follows that

[TABLE]

where $f_{2}(t)$ is another function of $t$. The equality (32) can be written as

[TABLE]

Introducing the constants

[TABLE]

the variables (suggested as canonical variables in the Hamiltonian formulation in [1, 2] )

[TABLE]

and noting that

[TABLE]

we recast (30) and (33) in terms of the velocity potentials and the current parameters, using the relation between the stream function and the velocity potentials that follows from (12) and (3):

[TABLE]

[TABLE]

The balance of all quantities when $x\to\pm\infty$ gives

[TABLE]

Recall that in the absence of current and Coriolis force the system is Hamiltonian and can be represented in the form

[TABLE]

where $H_{0}(\eta,\eta_{1},\xi,\xi_{1})$ is the corresponding Hamiltonian. $H_{0}$ is evaluated in terms of the canonical variables, see [22, 21] for details:

[TABLE]

where the operator $B\equiv\rho_{1}G(\eta)+\rho G_{11}.$

Hence we can represent all terms not related to vorticity as variational derivatives of $H_{0}:$

[TABLE]

In addition, (23), (24) can be written as

[TABLE]

Now we are in a position to write the equations (44), (45) in the form

[TABLE]

where

[TABLE]

[TABLE]

is an integral with a constant Hamiltonian density and zero variational derivatives which keeps the overall Hamiltonian density decaying to zero at $x\to\pm\infty.$ The ghost terms with $f_{1,2}$ are un-physical (i.e. their values do not affect the measurable physical quantities such as velocities and elevations), since by definition

[TABLE]

Nevertheless their variational derivatives produce the $f_{1,2}$ terms in (44) which are also un-physical and can in principle be absorbed in the definition of the potentials $\xi,$ $\xi_{1},$ in our case, in such a way that $\xi,$ $\xi_{1}$ tend to zero when $x\to\pm\infty.$ From (42) it follows that $H$ does not contain terms, linear in the field variables. The lowest order terms are quadratic in the field variables and they produce the linearised equations.

The expression for the Hamiltonian (47) can be obtained alternatively by evaluation of the total energy of the fluid (up to a constant Hamiltonian density):

[TABLE]

Taking into account the stratification of the fluid the above expression can be rewritten as

[TABLE]

The computations follow the routine from [16, 6, 7]. We mention only that due to the two-dimensional character of the dynamics, the final expression (47) depends only on the variables on the surface and on the interface. It is an important feature, that there is no contribution from layers I and IV where the wave motion on $s$ and $s_{1}$ does not take place [6, 7]. The terms related to the internal wave are the same as in [7] where flat surface approximation is considered.

4 Hamiltonian dynamics

The ideas for the Hamiltonian formulation of water waves coupled to a flow with a constant vorticity originate from [17] followed by [53]. For internal waves with vorticity the problem is studied in [16, 15, 4, 5].

The evolution equations (46) can be written in a canonical Hamiltonian form in terms of the variables

[TABLE]

The system described by the phase space variables $\eta,\eta_{1},z,z_{1}$ is Hamiltonian. More precisely,

[TABLE]

[TABLE]

The proof follows the lines of the one in [16]. Since

[TABLE]

the system of equations (46) can be written as

[TABLE]

where for convenience $k=0,1$, $\eta_{0}\equiv\eta$, $\xi_{0}\equiv\xi,$ $\chi_{0}\equiv\chi$ and $\Gamma_{0}\equiv\Gamma.$ This is an equivalent Hamiltonian form, in terms of the original variables and with respect to the Poisson bracket (PB)

[TABLE]

i.e.

[TABLE]

Here we have to specify the spaces where the functionals $A,B$ may belong. The antisymmetry of the PB requires for example

[TABLE]

which is possible (e.g. due to integration by parts) if and only if

[TABLE]

Thus, at least one of the functionals in the PB should satisfy

[TABLE]

Since in (53) one of the functionals is always $H$ and

[TABLE]

5 The pressure in the body of the fluid

The pressure in the body of the fluid can be evaluated from the functions $\bm{\varphi}$ and $\bm{\psi}.$ They can be recovered from $\xi$ and $\xi_{1}(x,t).$ In addition, there is an interdependency between $\bm{\varphi}$ and $\bm{\psi}.$ since

[TABLE]

is an analytic function of the variable $z=x+iy$ in the domain $\Omega_{1}\cup\Omega.$ Thus, the corresponding analytic functions in $\Omega$ and $\Omega_{1}$ can be recovered from their values at the boundaries of $\Omega$ and $\Omega_{1},$ i.e. from $\Phi(x,t)$, $\Phi_{1}(x,t)$ and $\Phi_{2}(x,t).$

From the definition of the Dirichlet-Neumann operators and (24) we have that

[TABLE]

and

[TABLE]

Adding up the previous two relations we obtain

[TABLE]

Using (57), and $B:=\rho_{1}G(\eta)+\rho G_{11}$ and recalling that

[TABLE]

we can write $\Phi,\Phi_{1},\Phi_{2}$ in terms of the Hamiltonian variables $\xi$ and $\xi_{1}$ as follows

[TABLE]

In the layers with a fixed constant vorticity $\bm{\gamma}$ we have Bernoulli conservation laws

[TABLE]

From (31) and (42) for the layer V′ we have

[TABLE]

Another possible derivation of (62) is the following one. We can evaluate $\tilde{f}_{1}$ using the asymptotic values in (61) at $x\to\pm\infty$ where no wave motion takes place ($y=h_{1}$) and thus $\psi_{1,x}\to 0$ and $\psi_{1,y}\to\gamma_{1}h_{1}+\kappa_{1}$ are given by the velocity of the current. Moreover, $\psi_{1}$ on the surface is zero when $x\to\pm\infty$ due to (22). Hence

[TABLE]

We can apply the same approach in deriving the pressure in layers III*′* and II*′*. At $y=0$ and $x\to\pm\infty$ due to the lack of a wave motion the pressure is $P_{0}=P_{atm}+\rho_{1}gh_{1}.$ For the layer III′ we have

[TABLE]

For the layer II${}^{\prime},$

[TABLE]

6 Scales

Let us introduce non-dimensional variables (without bars) related to the dimensional (barred) as follows:

[TABLE]

The constant $a$ represents the average amplitude of the waves under consideration, $\varepsilon$ is a small parameter which will be used to separate the order of the terms in the model.

From $\bar{v}=\bar{\eta}_{\bar{t}}+\bar{u}\bar{\eta}_{\bar{x}}$ it follows

[TABLE]

therefore, if $\mathcal{O}(\eta_{t})=1$ then $\mathcal{O}(v)=\varepsilon$ and thus the dimensional expression with $\mathcal{O}(v)=1$ (and similar for $v_{1}$) should be

[TABLE]

Since $v$ is a $y$-derivative of the velocity potential, and with the adopted definitions $\bar{\varphi}=\varepsilon h_{1}\sqrt{gh_{1}}\varphi$ etc. and thus

[TABLE]

The scales for $u$, $u_{1}$ etc. do not change - their dominant parts are the vorticity and current components of order 1; only the ’wave’ component (which is $x$ derivative of $\varphi$) is of order $\varepsilon$. The Dirichlet-Neumann operators have the following structure (e.g. any of the introduced operators $G$, $G_{ij}$):

[TABLE]

where $\bar{G}^{(n)}\sim\bar{\eta}^{n}\partial_{\bar{x}}^{n+1}$ i.e. $\bar{G}^{(n)}=\frac{\varepsilon^{n}}{h_{1}}G^{(n)}$:

[TABLE]

With this scaling and ignoring the linear terms, whose average is 0, the Hamiltonian can be expanded as

[TABLE]

where

[TABLE]

where the leading order, $\mathcal{O}(1)$ expression for the operators is (does not depend on $\eta$, $\eta_{1}$)

[TABLE]

where $D=-i\partial_{x}$. The quadratic part produces the linearised equations. Similarly

[TABLE]

The order $\varepsilon$ terms $G^{(1)}$ of the operators are given in [22, 21].

7 Linearised equations

The Hamiltonian equations with a Hamiltonian $H^{(2)}$ are the linearised equations:

[TABLE]

where $\partial^{-1}$ is the inverse of $\partial_{x}$. We can change the coordinates via a linear transformation according to $\partial_{T}=\partial_{t}+\kappa\partial_{x}$:

[TABLE]

Note that the equations contain Coriolis terms dependent on $\omega.$ Usually $\kappa$ is of magnitude several m/s, $2\omega=1.46\times 10^{-4}$ s*-1* and hence $2\omega\kappa\ll g:$

[TABLE]

where we introduced the notation

[TABLE]

Next we search for solutions, proportional to

[TABLE]

where $\Omega_{0}(k)$ is the dispersion law for the wave, the wave speed is

[TABLE]

Such a solution would be an eigenfunction for $D$ with a corresponding eigenvalue $k$. From the system (78) one can express $\xi$, $\xi_{1}$ as

[TABLE]

The remaining two equations lead to a linear homogeneous system for $\eta$ and $\eta_{1}$. We introduce the notations

[TABLE]

The second equation of (78) gives

[TABLE]

Then the last equation of (78) leads to

[TABLE]

The compatibility of the two equations gives a 4-th order equation for $c_{0}(k)$:

[TABLE]

The last formula generalises the irrotational one from [21]. In addition, from (83), (82) one can determine if $\eta$ and $\eta_{1}$ have the same or an opposite polarity for each possible propagation speed $c_{0}$. The velocity $c_{0}(k)$ is relative to an observer moving together with the flow at $y=0$, i.e. with velocity $\kappa.$ The wave speed for a stationary observer therefore is

[TABLE]

8 Long waves approximation

In the long waves approximation the physical scales are measured by the dimensionless parameter like $\delta=\frac{h_{1}}{\bar{L}}$. We will study the equations under the additional approximation that the wavelengths $L$ are much bigger than $h$ and $h_{1}$. Since

[TABLE]

Thus the wave number $k=2\pi/L=2\pi\delta$ and $\mathcal{O}(k)=\delta$. We further assume that $\delta^{2}=\mathcal{O}(\varepsilon).$ Recall that the operator $D$ has an eigenvalue $k$, thus we shall keep in mind that $\mathcal{O}(D)=\delta.$ Moreover the $x$-derivative of the velocity potentials do not get an extra factor of $\delta$ since $\mathcal{O}(\bar{v})=\varepsilon$ remains unchanged and similarly $v_{1}$. In other words the ’wave’ component of $u$ let’s call it $\tilde{u}=\varphi_{x}$ is of order $\varepsilon$ and similarly $\tilde{u}_{1}=\varphi_{1,x}$ is of order $\varepsilon$. Despite the assumption $\delta^{2}=\mathcal{O}(\varepsilon)$ we will keep both scales $\delta$ and $\varepsilon$ in order to keep track of the origin of the various terms. We will keep track only of the scale variables $\varepsilon,\delta$ and not of the other dimensional factors. For example, $H^{(2)}$ itself contains the following type of terms:

[TABLE]

Here $H^{(2)}$ is given in terms of $\tilde{u},\eta,\tilde{u}_{1},\eta_{1}$ which are not canonical variables. The canonical are the variables

[TABLE]

It is more convenient to work however with the variables which are $x$-derivatives of $z$ and $z_{1}$:

[TABLE]

In terms of these new variables, the Hamiltonian structure changes accordingly,

[TABLE]

In terms of the old variables $\tilde{u},\eta,\tilde{u}_{1},\eta_{1}$

[TABLE]

Since $\mathcal{O}(\varepsilon^{2}\delta^{2})=\varepsilon^{3}$ then already $H^{(2)}$ in the long wave approximation produces terms of order $\varepsilon^{3}$. For this reason we will not need terms of order $\varepsilon^{3}\delta^{2}\sim\varepsilon^{4}$ from $H^{(3)}$. Therefore, the relevant terms from $H^{(3)}$ are

[TABLE]

The next assumption in our approximation is that the pair of canonical Hamiltonian variables $p_{1}$ and $\eta_{1}$ (i.e. with respect to the Hamiltonian structure in (89)) associated to the free surface are of smaller order,

[TABLE]

Clearly $\tilde{u}_{1}\rightarrow\varepsilon\tilde{u}_{1}$. The approximate Hamiltonian of the system is with terms of orders $\varepsilon^{2}$ and $\varepsilon^{3}$:

[TABLE]

For the sake of simplification, let us introduce the following notations for some constants:

[TABLE]

The equations are

[TABLE]

Due to the last equation we have

[TABLE]

and for the other variables

[TABLE]

The leading order linear equations for $\tilde{u}$ and $\eta$ are

[TABLE]

The wavespeed $c$ of the solutions, proportional to $e^{ik(x-ct)}$ satisfies the quadratic equation

[TABLE]

Introducing $c_{0}(k)=c(k)-\kappa$ and noting that

[TABLE]

we write the equation for $c_{0}(k)$ in the form

[TABLE]

The solution is

[TABLE]

There are right ($c_{0}>0$) and left ($c_{0}<0$) running waves. We notice that in this approximation $c(k)$ is $k$-independent, i.e. constant.

In the zero vorticity case

[TABLE]

Let us introduce also $c^{\prime}_{0}=c^{\prime}-\kappa.$ From (99) it also follows that in the leading order $\tilde{u}=(\rho c_{0}\eta)/h$, i.e.

[TABLE]

for some yet unknown quantity $q$. Nevertheless from (98) and (104) we determine

[TABLE]

Note that $\mathcal{O}(\tilde{u}_{1})=\varepsilon$, thus $\mathcal{O}(\varepsilon\tilde{u}_{1})=\varepsilon^{2}$ and such terms will be neglected in (96) and (97) . The most general form for $q$ is

[TABLE]

for some constants $b_{1},b_{2}.$ Now we are in a position to express everything in (96) and (97) only via the variable $\eta$. This way we obtain two evolutionary equations for $\eta$ which should coincide up to $\mathcal{O}(\varepsilon)$. The equality of their coefficients allows us to find

[TABLE]

$\eta$ satisfies the KdV equation [39], (see also [10, 35])

[TABLE]

which, in the limit when all vorticities are zero becomes

[TABLE]

From $\eta$ one can recover $\tilde{u}$ and as well as $\eta_{1}$ and $\tilde{u}_{1}$:

[TABLE]

We notice that in leading order

[TABLE]

and therefore both positive and negative relative polarities for $\eta$ and $\eta_{1}$ are possible. The KdV approximation for an internal wave coupled to a free surface for a different configuration of the currents is derived in [15].

The KdV equation represents a balance between a nonlinearity term $\varepsilon\eta\eta_{x}$, and dispersion term $\delta^{2}\eta_{xxx}.$ In the above considerations $\varepsilon$ and $\delta^{2},$ are of the same order and as a result we can have the stable soliton solutions of the KdV equation. However, there are various geophysical scales and other situations are possible, including $\delta\sim\varepsilon^{2}.$ In such case $\delta^{2}\sim\varepsilon^{4}\ll\varepsilon$ and instead of the KdV equation the relevant model is the dispersionless Burgers equation ($\partial_{\tau}=\partial_{t}+c\partial_{x}$)

[TABLE]

It is well known that the solutions of this equation always form a vertical slope and break. Such wave-breaking phenomenon is well known for internal waves in the ocean. This is a mechanism that causes mixing in the deep ocean, [40].

There are other integrable systems which provide an approximation of the equations in the Boussinesq regime, such as the Kaup-Boussinesq system investigated firstly by D.J. Kaup [38], see also [33]. Two-component integrable systems, that can match the model equations up to order $\delta^{2}$, are the 2-component Camassa-Holm system and the Zakharov-Ito system [14, 32, 31, 25, 23, 24].

9 Discussion and conclusions

We consider a two-media system of liquids with different densities, free surface and a free internal surface separating the liquids. The bottom of the system is considered horizontal and flat. We studied the surface and internal waves driven by gravity and Coriolis forces and interacting with a current. The underlying current is in the form of a shear flow with a specific velocity profile. The current has constant vorticity at the horizontal strips where the surface and internal waves are localised. The model is aimed at geophysical applications, where a typical configuration is the one of a thin shallow layer of warm and less dense water over a much deeper layer of cold denser water. The governing equations are written in a canonical Hamiltonian form, which gives rise to a systematic approach for possible approximations. In particular, small amplitude and long-wave regimes are studied. There are various geophysical scales, allowing for smooth solitons at the KdV regime as well as breaking waves in the very large wavelength regime, when the equations are asymptotically equivalent to the dispersionless Burgers equation. In the case of a free surface, even in the case of very small amplitudes, the internal wave is usually coupled to the surface wave. This has an impact on the propagation speeds and is observed in other related models [16, 15]. Other asymptotic regimes, e.g. related to the Nonlinear Schrödinger equation [51] remain to be studied. The stability of the waves interacting with currents is another important issue that needs to be addressed in the future. Related recent works in this connection are [13, 26, 30].

Acknowledgements

The author acknowledges Seed funding grant support from Dublin Institute of Technology for a project in association with the Environmental Sustainability and Health Institute, Dublin. The author is grateful to Prof. Adrian Constantin for many valuable discussions and to an anonymous referee for many important suggestions.