Benjamin-Ono model of an internal wave under a flat surface

A two-layer fluid system separated by a pycnocline in the form of an internal wave is considered. The lower layer is infinitely deep, with a higher density than the upper layer which is bounded above by a flat surface. The fluids are incompressible and inviscid. A Hamiltonian formulation for the fluid dynamics is presented and it is shown that an appropriate scaling leads to the integrable Benjamin-Ono equation.


Introduction
Due to the ubiquitous presence of waves on the surface of oceans much research has been completed which attempts to explain the mechanisms responsible for such waves. However, despite many centuries of reporting, by mariners in particular, of unexplained phenomena beneath the surface, such as the "dead water" phenomenon coined by Fridtjof Nansen [30], it is relatively recently that some progress of a mathematical nature has been achieved in the description of internal waves. In many situations both the surface and the internal waves are moving in the presence of currents. Due to the nonlinearity the wave-current interactions are quite complex and necessitate the study of rotational fluids [7,9,25,31,33] .
The significant findings of Vladimir Zakharov in 1968 in [35] have established the use of the Hamiltonian approach for dynamic descriptions of wave motion. For single layer systems both irrotational [4,18,27,28] and rotational [10,11,14,16,17,32,34] set-ups have been examined within the Hamiltonian framework. Stratified systems which contain a pycnocline separating two layers in the form of an internal wave have been considered in an irrotational context in [2,3,19,20] and in a rotational context, which is of most interest, in [5,6,12,15].
We aim to establish a model for a system with a flat-bed, flat surface and internal wave. Our setting captures the nonlinear dynamics as influenced by both the internal wave behaviour and the presence of a depth-varying current. In physical terms the system can be thought of as a non-mixing oceanic environment consisting of two discrete fluid bodies at different densities in the presence of a current such as the case in the Pacific ocean where the EUC (Equatorial Under Current) influences a stratified region of the ocean [22].
Perturbative techniques can be used to develop such models. In [8] the chosen scaling regime leads to a KdV type model. We choose to adopt a similar approach, but with a scaling which leads to a Benjamin-Ono type approximation.

Set-up and governing equations
The flow we consider consists of two layers, Ω and Ω 1 , which have different densities, due to different salinity levels or temperatures. They are separated by a very thin layer termed as a pycnocline (whose thickness is neglected), as shown in Figure 1 where the horizontal axis is x and the vertical axis is y. The internal wave is formed at the pycnocline y = η(x, t). The mean of η is assumed to be zero, that is for all t so that η(x, t) measures the elevation of the internal wave with respect to the level y = 0. The domains Ω and Ω 1 are defined as and Throughout the paper we will use the subscript 1 to denote quantities pertaining to the upper layer, while quantities referring to the lower layer will appear without subscript. The domain Ω is infinitely deep and Ω 1 has a flat surface. The physical justification for this assumption of absence of surface motion is due to the small amplitude of the surface waves in comparison to the amplitude of the internal waves. The fluids are considered to be inviscid and incompressible with ρ and ρ 1 being the respective constant densities of the lower and upper media. Stable stratification requires The velocity fields are given by and the incompressibility implies divV = 0 and divV 1 = 0 or u x + v y = 0 and u 1,x + v 1,y = 0.
The Euler equations govern the fluid dynamics in each layer, where p, p 1 are the dynamic pressure terms in the corresponding layers and g = (0, −g) is the Earth's gravitational acceleration.
In addition, the kinematic boundary condition for the interface between the layers (the pycnocline) must be satisfied, that is Velocity potentials ϕ and ϕ 1 are introduced to capture the irrotational components of the velocity fields. These components in general include timeindependent constant currents so that the velocity potentials can be further decomposed as ϕ ≡ ϕ + κx and ϕ 1 ≡ ϕ 1 + κ 1 x where κ and κ 1 are the respective time-independent currents at y = 0. In this decomposition ϕ and ϕ 1 represent the potential components related to the wave motion. In the presence of a current with constant vorticity we can therefore represent the velocity field components in Ω via and similarly for Ω 1 via where γ = u y − v x and γ 1 = u 1,y − v 1,x are the constant non-zero vorticities (see [13], where the situation with a free surface is studied). The incompressibility of the fluid (4) allows us to introduce a stream function ψ in Ω via u = ψ y and v = −ψ x and a stream function ψ 1 in Ω 1 via The following assumptions will be made: η(x, t) < h 1 for all values of x and t; η(x, t) and ϕ 1 (x, y, t), are in the Schwartz class S(R) with respect to the x variable (for any y and t); ϕ(x, y, t) is in the Schwartz class S(R 2 ) with respect to both the x and y variables (for any t). Due to these assumptions for large absolute values of x the internal wave attenuates, meaning that Moreover, lim for all values of x and t. The physical reasoning for this is the absence of wave motion at infinite depth y → −∞.
Euler's equations (5), (6) in terms of the introduced variables are [8] and where p and p 1 are the dynamic pressure terms, ρ and ρ 1 are the constant densities, f and f 1 are some so far arbitrary functions of time. Their presence is related to the fact that the velocity potentials are determined up to an additive term whose gradient is zero. For further convenience we choose At the interface y = η(x, t) (denoted by a subscript "i") the dynamic pressure terms are equal giving the Bernoulli equation [24] ρ where χ and χ 1 are the interface stream functions. Furthermore, it could be shown [5] that χ = χ 1 , and so the Bernoulli equation becomes This form suggests the introduction of a single variable ρ ϕ − ρ 1 ϕ 1 which is going to provide one of the Hamiltonian coordinates (the momentum) in the next section. The other one (the coordinate) is the variable η(x, t), which satisfies the kinematic boundary condition (7)

The Hamiltonian formulation
The functional H, which describes the total energy of the system, can be written as the sum of the kinetic, K, and potential energy, V contributions. The potential part, must be However, the potential energy is always measured from some reference value, e.g. V (η = 0) which is the potential energy of the current (without wave motion). Therefore, the relevant part of the potential energy, contributing to the wave motion is In order to determine the kinetic energy of the wave motion, from the total kinetic energy of the fluid one should subtract again the constant, but infinite kinetic energy of the current which is In terms of the dependent variables η(x, t), ϕ(x, t) and ϕ 1 (x, t) this kinetic energy is The Hamiltonian is therefore The assumption that ϕ(x, y, t) is in the Schwartz class S(R 2 ) with respect to both the x and y variables gives lim x→±∞ 0 −∞ (γy + κ) ϕ(x, y, t)dy = 0 for any t, and furthermore We would like to write the Hamiltonian only in terms of the one dimensional variables pertaining to the interface. To this end the Dirichlet-Neumann operators are introduced, where n and n 1 are the unit exterior normals, 1 + (η x ) 2 is a normalisation factor and are the interface velocity potentials.
Since usually there is no jump in the current velocity, in what follows we take κ = κ 1 .
The Hamiltonian can therefore be written in terms of conjugate variables η(x, t) and ξ(x, t), following the procedure in [8], as and the operator B, as per [19], is introduced as We point out that due to the initial assumptions on η, ϕ and ϕ 1 and (25), (27), the Hamiltonian variables η and ξ are in S(R) with respect to the x variable for all t.

Scaling of the Hamiltonian
By the introduction of a small arbitrary constant parameter, δ, defined by where λ is the wavelength of the internal wave, meaning the scaled system will be considered as having long waves, the variables will be scaled according to η → δη, ξ → ξ and ∂ → δ∂ noting that the differential operators ∂ and D are related by and the wave number k := 2π/λ is an eigenvalue of D for the monochromatic linear waves in the form e ikx and therefore ∂ ∼ O(δ). The (unscaled) Dirichlet-Neumann operators can be expanded in terms of orders of η as [21] and so the expanded Dirichlet-Neumann operators are scaled as, noting from [19] that the constant term for the infinite lower layer is |D|, Using the expansion for the hyperbolic tangent it can be written that and therefore the Dirichlet-Neumann operators can be expanded further as The details of the calculations are provided in the Appendix. The final expression for the Hamiltonian (26) is where A := (ργ − ρ 1 γ 1 )κ + g(ρ − ρ 1 ).

The Benjamin-Ono approximation
The variable u = ξ x which assumes the role of momentum in the Hamiltonian approach (cf. [2,3]), in a similar fashion to η assuming the role of the generalised coordinate, is introduced.
Noting that we can rescale the Hamiltonian by a factor δ 2 by the choice of a proper time scale, the Hamiltonian (37) in terms of η and u is given by The equations of motion (29) can be rewritten as and so therefore and In the leading order (40) -(41) give where c = c(k) is the wave speed. The linearised equations of motion produce and − (c − κ)u = (−A + Γc)η.
Multiplying both equations by each other gives their compatibility condition which, in turn, gives the dispersion relation As usual there are two solutions for the two signs ±, which correspond to left and right-running waves with respect to the average flow with velocity κ. Moreover, the leading order approximation from (42) is for a known c. Our aim now is to find the next order approximation of the form for some, as yet, unknown values α and β, with a view to establishing a Benjamin-Ono type approximation. Equation (40) is hence written, to first order of δ, as Likewise, equation (41) is written as Noting that η t = −cη x + O(δ) It is noted that the coefficient of the η x term Comparing the |D|η x terms in (46) and (47) and so Comparing the ηη x terms in (46) and (47) 2 (50) The equation for η, from (46), is therefore given by which can be written as The second component u can be expressed with η by (45) where now α and β are given by (50) and (49).
We should keep in mind that there are always two sets of equations for the left and right running waves corresponding to the different choices of the sign in (44).
The obtained equation (52) is the well known Benjamin-Ono (BO) equation [1,29] which is an integrable equation whose solutions can be obtained by the Inverse Scattering method, see [23,26] and the references therein.

Solitary wave solution
The standard form of the BO equation is for which the one-soliton solution is known, where C 0 and X 0 are constants. Hence, the equation has a solution , which is the same solution, if one replaces the arbitrary constant σC 0 with C 0 . Transforming X using X → X − cT gives the equation with a solution , and hence after a rescaling of η, X and T , equation (52) has a solution where the amplitude η 0 and the initial displacement x 0 are arbitrary constants, and Expression (58) shows that the wavespeed of the soliton depends on its amplitude η 0 . Also, note that the sign of C 0 depends on η 0 and the parameters of the system.

Discussion
The illustrative one-soliton solution of the BO equation (55) suffers, however, from the following disadvantages. First, it is not in the Schwartz class in the x-variable, which is not a very serious disadvantage from the physical point of view. Second, it violates the assumption (1) for η since for the expression (55) R η(X, T )dX = π = 0. This is due to the fact that the initial condition η(x, 0) for this particular solution does not satisfy the mentioned assumptions. Therefore, extra care is needed when the inverse scattering, or other methods are applied to the BO equation when modelling internal waves.

Appendix
Writing the Dirichlet-Neumann operators as where the superscript identifies the order of δ the relevant terms that we will be using are: G [1] = δ|D|, G [2] = 0, G [3] = δ 3 DηD − δ 3 |D|η|D| Therefore noting that ρ 1 G(η; δ) = δρ 1 |D| + δ 3 ρ 1 DηD − δ 3 ρ 1 |D|η|D| + O(δ 5 ) and ρG 1 (η; δ) = δ 2 ρh 1 D 2 − δ 3 ρDηD + O(δ 5 ) the operator B is transformed to This can be written as and so Noting that D 2 |D| = |D| we can write Using the expansion the inverse of the operator B is given by Writing the operator as where, again, the superscript identifies the order of δ the relevant terms that we will be using are: The Hamiltonian can therefore be written, using components of the expanded operators as H(η, ξ) = 1 2 Replacing the expressions for G [j] , G