#### Abstract

We consider three inviscid, incompressible, irrotational fluids that are contained between the rigid walls and and that are separated by two free interfaces and . A generalized nonlocal spectral (NSP) formulation is developed, from which asymptotic reductions of stratified fluids are obtained, including coupled nonlinear generalized Boussinesq equations and -dimensional shallow water equations. A numerical investigation of the -dimensional case shows the existence of solitary wave solutions which have been investigated for different values of the characteristic parameters.

#### 1. Introduction

Since some early studies started in the 1950s (see [1, 2] and refs. therein), much interest has been devoted over the years to the flow of stratified fluids, both from a theoretical and from applied point of view. In particular, multifluid flows characterized by gravitationally stable density interfaces are a useful mathematical model in geophysics and engineering [3–5]. On the other hand, many analytical studies have been dedicated to the equations that describe the system of two ideal fluids, separated by a free interface [6]. Indeed, a great effort has been devoted to modeling the evolution of the internal waves in a two-fluid system and to deriving asymptotic reductions (see [7] for a recent review). Among them, of particular relevance are the Benjamin-Ono (BO) equation [8, 9] and the intermediate long wave (ILW) equation [10, 11]. Moreover, in [12] model equations governing the evolution of fully nonlinear long waves at the interface of two immiscible ideal fluids were derived, while the weakly nonlinear limit was treated in [13]; finally in [14] a Hamiltonian formulation of the two-fluid system was derived and through a perturbative theory a systematic analysis of the long waves scaling regimes was carried out. On the other hand, in recent years, a nonlocal formulation of water waves for both and dimensions was presented in [15], where the original equations with unknown boundary conditions are replaced by an integrodifferential equation and a nonlinear partial differential equation, both of which are formulated in a known domain. The nonlocal formulation obtained in [15] is derived from the general approach to studying boundary value problems for linear and nonlinear PDEs introduced in [16]. A crucial role in such approach is played by a nonlocal equation, called the global relation [16]. A generalization of the results obtained in [15] was presented in [17], where a nonlocal formulation was derived, governing two ideal fluids separated by a free interface and bounded above either by a rigid lid or by a free surface [17]. Due to the dependence on a free spectral parameter, the corresponding equations are usually called the nonlocal spectral (NSP) equations of the two-fluid system. The NSP equations were particularly useful for deriving asymptotic approximations; we wish to point out an asymptotically -dimensional generalization of the intermediate long wave (ILW) equation reported in [17] which includes the KP equation and the Benjamin-Ono equation as limiting cases. Numerical investigations indicated the existence of lump type solutions, with a speed versus amplitude relationship shown to be linear in the shallow, intermediate, and deep water regime.

However, to the best of our knowledge, the phenomenological models for more than two fluids mentioned at the beginning of this section have not been paralleled by any analytical study. This prompted us to develop a generalization of the NSP formulation to the case of three ideal fluids, separated by two free interfaces and limited above by a rigid lid. Namely, we consider three inviscid, incompressible, irrotational fluids that are confined between the rigid lids and and are separated by two free interfaces and . We derive in the following the NSP equations governing the evolution of the three-fluid system. Specifically, define the functions , , and , where , , and are, respectively, the velocity potentials in the lower, the intermediate, and upper layer (see Figure 1).

The NSP formulation is given by the following equations. *First layer (bottom)* *intermediate layer* *third layer (top)* *Bernoulli’s equations*In (1)–(5), , , , where is a free spectral parameter and the constants , , and denote gravity, surface tension of , and surface tension of , respectively.

In the next section, we derive the above equations, starting from the classic equations governing three ideal fluids separated by two free interfaces bounded above and below by rigid lids. Moreover, we derive conservation laws and integral identities for the three-fluid system for the NSP formulation.

Section 3 is devoted to the derivation of the weakly nonlinear equations: after a suitable nondimensionalization of the variables, we obtain the reduction to a system of shallow water equations in the weakly nonlinear limit.

In Section 4, under the assumption of maximal balance, we introduce travelling wave variables moving only to the right and we study only the -dimensional case. In terms of the new variables, we finally obtain a system of coupled nonlinear shallow water equations which we study numerically in terms of the parameters entering the theory and that we show to admit solitary wave solutions.

#### 2. A Nonlocal Spectral Formulation of the Classic Three-Fluid Equations

##### 2.1. A Weak Formulation of the Classic Three-Fluid Equations

We recall the classic equations governing three ideal fluids separated by two free interfaces and and bounded above and below by rigid lids. It is assumed that the lower fluid is of density , the intermediate fluid is of density , and the upper fluid is of density with . The equations are given in terms of the interface variables and the velocity potentials , , and associated with the lower, intermediate, and upper fluid domains, respectively, We also require that , and as . In the previous equations, the constants , , and denote gravity, the surface tension associated with the free interface , and the surface tension associated with the free surface . Equations (6), (9), and (12) express the fact that fluids are divergence and curl-free. Equations (15)-(16) express that the jump in pressure across an interface is balanced by surface tension.

We now obtain a weak formulation of (6)–(16), expressed in terms of , , , , , and . Let us define the domains For any bounded, harmonic functions , defined in , , and , respectively, and satisfying the following identities hold at any given time : Equations (19) and (21), together with the Bernoulli equations (15) and (16), constitute a nonlocal system of equations describing the three-fluid system.

##### 2.2. Derivation of the Nonlocal Spectral Equations

We derive from (19) the NSP equation (1) given in the Introduction. To motivate the derivation, suppose that satisfiesDefine by Then, conditions (22) turn into the following: The solution to the above ODE is given by reverting back to physical space, The previous equation shows that any harmonic function whose -derivative at vanishes can formally be written as a sum of functions . Therefore, by linearity it suffices to require that (19) holds for the parametrized family of functions Substituting into (19), we obtain We simplify (28) by noting thatUsing this identity in (28) and integrating by parts, we finally get the nonlocal spectral equation (1) as follows:

Along the same lines used for the first layer, let us derive (2)-(3) from (20).

We define the basic functions and putting into (20) givesWe simplify (32) by noting that Using these identities in (32) and integrating by parts, we finally get the nonlocal spectral equation (2).

Similarly, let us consider the basic functionsand put into (20), obtainingWe simplify (35) by noting that Using these identities in (35) and integrating by parts, we finally get the nonlocal spectral equation (3).

For the top layer, we require that (21) holds for the parametrized family of functions Substituting into (21), we obtain We simplify (38) by noting that Using this identity in (28) and integrating by parts, we finally get the nonlocal spectral equation (4).

##### 2.3. Conservation Laws and Integral Identities

In [17], conservation laws and integral identities for a two-fluid system were derived from the nonlocal spectral formulation. We now derive the analogous conservation laws for the three-fluid system for the NSP formulation. We start by expanding (2), (3) for small . At first order, by setting to zero the coefficient of and putting , we get Subtracting (41) from (40), we obtain the resultwhich corresponds to mass conservation in the intermediate fluid domain.

Similarly, by setting to zero the coefficient of , with , we getWhen we subtract (43) from (44), we obtainsimilarly, if we set , at the same order we get the corresponding equation:Equations (45) and (46) describe the evolution of the center of mass; the right-hand side is the momentum of the fluid.

Next, we obtain a last identity at third order, setting to zero the coefficient of :with ; the above is a virial type formula, analogous to the one obtained in [17].

#### 3. Weakly Nonlinear Equations

##### 3.1. Nondimensionalization of the NSP Equations

In order to derive weakly nonlinear equations, we nondimensionalize all physical variables in (1)–(5) according towhere is a characteristic wavelength, is a nondimensional parameter, and and are a characteristic amplitude and velocity, respectively.

Then, (1)–(5) become, after dropping primes and letting , *first layer (bottom)* *intermediate layer* *third layer (top)* *Bernoulli’s equations*In the above relations, is the* nonlinearity ratio* and is the so-called* aspect ratio* [12]; moreover, it is In the following, we assume , and , which corresponds to the case of* long waves*; moreover, we take the nonlinearity ratio to be .

The above assumptions imply that we are interested to derive asymptotic reductions of the NSP equations in the case of weakly nonlinear long waves.

##### 3.2. Derivation of the Equations

We first expand (49) and (53) in and . From (49), taking the Inverse Fourier Transform (IFT), we getwhile (53) givesOn the other hand, from (49), it is also useful to derive the following reduction: We now use (58) at leading order (L.O.) in (56) and, after some manipulations, we obtain from (57) the following equation:By (58) and (59), the* first generalized Boussinesq* equation obtains

We now turn our attention to the intermediate layer.

By equating the expansions of (50) and (51), we get the following expression for : which in turn gives where, from (58) and (69) at L.O., we have When (63) are used in (62), after integrating with respect to , we finally get We now take the -derivative of (57) and use (64), obtaining by equating (65) and (58), we get the* second generalized Boussinesq* equation:

We now expand (52) and take the IFT, gettingwhile from the expansion of (54) we get On the other hand, from (52), it is also useful to derive the following expansion: We use (68) at L.O. in (69) and get In order to get an expression for in (68), we now expand (50), getting By taking the -derivative of (71) and integrating with respect to , we get where the L.O. of (58) has also been used.

When we substitute (72) back into (68) and take the -derivative, we obtain By equating (70) and (73), we finally obtain the* third generalized Boussinesq* equation:

The system of 3 generalized Boussinesq equations (60), (66), and (74) will now be reduced to a system of 2 independent equations. We first take the shallow water limit, given by in (60) and in (66), obtaining the following reductions: Next, the combination of (75) and (76) gives which we substitute back into (74) and we obtain Rearranging (66), we get then the system of weakly nonlinear, shallow water equations is given by (78) and (79).

#### 4. Multiple Scale Derivation: Special Case

In order to obtain some interesting limiting equations, we now make the assumption of* maximal balance* (the small terms are of the same order) . This reflects a balance of weak nonlinearity and weak dispersion.

In this section, we study a* Special Case*, that is, when the waves velocities in the shallow water equations (78) and (79) are the same:

Here below, we assume the following asymptotic expansions for , , and : We also introduce new variables: which are* travelling waves variables* and describe the direction of waves propagation along the positive -axis: represents right moving waves and represents left moving waves.

In terms of the new variables, one obtains

Considering the assumption of maximal balance together with (80), (79) becomes Substituting the expansions (81), (83) into (86) and equating L.O. terms yield the wave equations: and , whose solutions, respectively, are and .

We also assume unidirectional waves and* only* work with the right moving waves, so that Then, keeping the terms, (86) becomes

In order to remove secular terms, the right-hand side of (88) needs to be set to zero (see [18]):

Following the same steps as above for (78), together with (89), and studying only the -dimensional case, we obtain a set of two equations:

In the following, we will study only the -dimensional case of the derivation.

Finally, define the functions , ; then, differentiating (90) with respect to , we can write the resulting equations in terms of and as Substitute (91) into (92) and obtain Our multiple scale derivation leads therefore to having a system of two equations: (91) and (93) constitute a system of coupled nonlinear* shallow water* equations.

The two coupled shallow water equations will be studied numerically as a function of the parameters entering the theory, in order to prove the existence of solitary waves and analyze their behaviour.

##### 4.1. Numerical Investigation of Traveling Wave Solutions of the Coupled Shallow Water Equations

We now investigate whether (91) and (93) possess traveling wave solutions. For convenience, take . Passing to a traveling coordinate system moving with velocity in the direction, (91) and (93) become, respectively,

To solve (94) numerically, we use the spectral renormalization (SPRZ) method, developed in [19]. We take the Fourier transform of the previous equations and get, upon rearranging,

In general, we cannot find a solution to (95)-(96) by naive iteration. Instead, we assume that and , where , are unknown parameters and , are unknown functions (this step is the renormalization part). Then, (95) can be written in terms of and as Note that by multiplying (97) by (where denotes the conjugate of ), rearranging and integrating the result we get Along the same line, (96) can be written in terms of and as By multiplying (99) by (where denotes the conjugate of ), rearranging and integrating the result we get

We use the above SPRZ scheme to solve for the modes and when , , , and .

Figure 2 displays the resulting speed versus amplitude relationship, which is nearly linear for (shallow water regime). Note that the horizontal axis in Figure 2 is .

Figure 3 shows the cross sections of and ; on the left there are the solutions starting from a Gaussian function and on the right the ones starting from a step function, both for a typical speed in the shallow water regime.

After obtaining these results, it would now be of particular interest from the applied point of view to address the issue of a local breakdown of the nonlinear internal waves propagating in the density stratified fluid [20]. Phenomena of this kind can induce localized turbulence in the stratified fluid and are relevant for their occurrence in the atmosphere and the oceans. We plan to address this issue in the future together with a generalization of our model to the -dimensional case.

#### Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

#### Acknowledgment

The authors wish to thank M. J. Ablowitz for useful suggestions.