Freak Waves as Nonlinear Stage of Stokes Wave Modulation Instability

V.E. Zakharov^{a,b,1} A.I. Dyachenko^{b,∗} A.O. Prokofiev^{b}

aDepartment of Mathematics, University of Arizona, Tucson, AZ, 857201, USA bLandau Institute for Theoretical Physics, 2 Kosygin str., Moscow, 119334, Russia

Abstract

Numerical simulation of evolution of nonlinear gravity waves is presented. Simulation is done using two-dimensional code, based on conformal mapping of the fluid to the lower half-plane. We have considered two problems: i) modulation instability of wave train and ii) evolution of NLSE solitons with different steepness of carrier wave. In both cases we have observed formation of freak waves.

Key words: free surface, freak wave, modulation instability PACS: 02.60.Cb, 47.15.Hg, 92.10.Dh, 92.10Hm

1 Introduction

Waves of extremely large size, alternatively called freak, rogue or giant waves are a well-documented hazards for mariners (see, for instance Smith (1976), Dean (1990), Chase (2003)). These waves are responsible for loss of many ships and many human lives. Freak waves could appear in any place of the world ocean (see Earle (1975), Mori et al. (2002), Divinski et al. (2004)); however, in some regions they are more probable than in the others. One of the regions where freak waves are especially frequent is the Agulhas current of the South-East coast of South Africa (see Gerber (1996), Gutshabash et al. (1986), Irvine and Tilley (1988), Lavrenov (1998), Mallory (1974)). In the paper by Peregrine (1976) it was suggested that in areas of strong current such as the Agulhas, giant waves could be produced when wave action

∗ Corresponding author.

Email addresses: zakharov@math.arizona.edu (V.E. Zakharov),

alexd@landau.ac.ru (A.I. Dyachenko), alexpro@itp.ac.ru (A.O. Prokofiev).

1 Waves and Solitons LLC, W. Sereno Dr., Gilbert, AZ, 85233, USA

Preprint submitted to Elsevier Science 1 May 2006

is concentrated by reflection into a caustic region. According to this theory, a variable current acts analogously to an optic lens to focus wave action. The caustic theory of freak waves was supported since that time by works of many authors. Among them Smith (1976), Gutshabash et al. (1986), Irvine and Tilley (1988), Sand et al. (1990), Gerber (1987), Gerber (1993), Kharif and Pelinovsky (2003). The statistics of caustics with application to calculation of the freak wave formation probability was studied in the paper of White and Fornberg (1998).

On our opinion, a connection between freak wave generation and caustics for swell or wind-driven sea is the indisputable fact. However, this is not the end of the story. Focusing of ocean waves by an inhomogeneous current is a pure linear effect. Meanwhile, no doubts that freak waves are essentially nonlinear objects. They are very steep. In the last stage of their evolution, the steepness becomes infinite, forming a ”wall of water”. Before this moment, the steepness is higher than one for the limiting Stokes wave. Moreover, a typical freak wave is a single event (see, for instance Divinski et al. (2004). Before breaking it has a crest, three-four (or even more) times higher than the crests of neighbor waves. The freak wave is preceded by a deep trough or ”hole in the sea”. A characteristic life time of a freak wave is short - ten of wave periods or so. If the wave period is fifteen seconds, this is just few minutes. Freak wave appears almost instantly from a relatively calm sea. Sure, these peculiar features of freak waves cannot be explained by a linear theory. Focusing of ocean waves creates only preconditions for formation of freak waves, which is a strongly nonlinear effect.

It is natural to associate appearance of freak waves with the modulation instability of Stokes waves. This instability is usually called after Benjamin and Feir, however, it was first discovered by Lighthill (1965). The theory of instability was developed independently by Benjamin and Feir (1967) and by Zakharov (1966). Feir (1967) was the first who observed the instability experimentally in 1967.

Slowly modulated weakly nonlinear Stokes wave is described by nonlinear Shrödinger equation (NLSE), derived by Zakharov (1968). This equation is integrable (see Zakharov and Shabat (1972)) and is just the first term in the hierarchy of envelope equations describing packets of surface gravity waves. The second term in this hierarchy was calculated by Dysthe (1979), the next one was found a few years ago by Trulsen and Dysthe (1996). The Dysthe equation was solved numerically by Ablovitz and his collaborates (see Ablovitz et al. (2000 and 2001)).

Since the first work of Smith (1976), many authors tried to explain the freak wave formation in terms of NLSE and its generalizations, like Dysthe equation. A vast scientific literature is devoted to this subject. The list presented below is long but incomplete: Peregrine (1983), Peregrine et al. (1988), Tanaka (1990), Trulsen and Dysthe (1996), Trulsen and Dysthe (1997), Trulsen (2000), Trulsen et al. (2000), Ablovitz et al. (2000 and 2001), Onorato et al. (2000a), Onorato et al. (2000b), Onorato et al. (2001), Onorato et al. (2002).

One cannot deny some advantages achieved by the use of the envelope equations. Results of many authors agree in one important point: nonlinear development of modulation instability leads to concentration of wave energy in a small spatial region. This is a ”hint” regarding possible formation of freak wave. On the other hand, it is clear that the freak wave phenomenon cannot be explained in terms of envelope equations. Indeed, NLSE and its generalizations are derived by expansion in series on powers of parameter λ ≃ \frac{1}{Lk}, where k is a wave number, L is a length of modulation. For real freak wave λ ∼ 1 and any ”slow modulation expansion” fails. However, the analysis in the framework of the NLS-type equations gives some valuable information about formation of freak waves.

Modulation instability leads to decomposition of initially homogeneous Stokes wave into a system of envelope solitons (more accurately speaking - quasisolitons (Zakharov and Kuznetsov (1998), Zakharov et al. (2004)). This state can be called ”solitonic turbulence”, or, more exactly ”quasisolitonic turbulence”. In the framework of pure NLSE, solitonic turbulence is ”integrable”. Solitons are stable , they scatter on each other elastically. However, even in this simplest scenario, spatial distribution of wave energy displays essential intermittency. More exact Dysthe equation is not integrable. In this model solitons can merge, this effect increases spatial intermittency and leads to establishing of chaotic intense modulations of energy density. So far this model cannot explain formation of freak waves with λ ∼ 1.

This effect can be explained if the envelope solutions of a certain critical amplitude are unstable, and can collapse. In the framework of 1-D focusing NLSE solitons are stable, thus solitons instability and the collapse must have a certain threshold in amplitude. Instability of a soliton of large amplitude and further collapse could be a proper theoretical explanation of the freak wave origin.

This scenario was observed in numerical experiment on the heuristic onedimensional Maida-McLaughlin Tabak (MMT) model (see Majda et al. (1997)) of one-dimensional wave turbulence Zakharov et al. (2004). At a proper choice of parameters this model mimics gravity waves on the surface of deep water. In the experiments described in the cited paper instability of a moderate amplitude monochromatic wave leads first to creation of a chain of solitons, which merge due to inelastic interaction into one soliton of a large amplitude. This soliton sucks energy from neighbor waves, becomes unstable and collapse up to λ ∼ 1, producing the freak wave. We believe that this mechanism of freak wave formation is universal.

The most direct way to prove validity of the outlined above scenario for freak wave formation is a direct numerical solution of Euler equation, describing potential oscillations of ideal fluid with a free surface in a gravitational field. This solution can be made by the methods published in several well-known articles, Dommermuth et al. (1987); West et al. (1987); Clamond and Grue (2001). Here we use another method, based on conformal mapping. It shoul be mention that idea to exploit comformal mapping for unsteady flows was presented in Ovsyannikov (1973), and later in Meison et al. (1981); Chalikov and Sheinin (1998). Method used in this article has origin in (Dyachenko et al. (1996), has been using in Zakharov et al. (2002), and was finally formulated in Dyachenko (2001). This method is applicable in 1 + 1 geometry, it includes conformal mapping of fluid bounded by the surface to the lower half-plane together with ”optimal” choice of variables, which guarantees well-posedness of the equations (Dyachenko (2005)) and existence of smooth, unique solution of the equations for a finite time (Shamin (2006)). Here we would like to stress that one of the main goal of this paper is to demonstrate effectiveness of the conformal variables to simulate exact 2D potential flow with a free boundary. Earlier, fully nonlinear numerical experiments regarding wave breaking, freak wave formation, comparison with weakly nonlinear model (such as Nonlinear Shredinger equation) were done in the papers Dold and Peregrine (1986); Tanaka (1990); Banner and Tian (1998); Henderson et al. (1999); Clamond and Grue (2002). From the other hand, using conformal approach we have studied in the paper Zakharov et al. (2002), the nonlinear stage of modulation instability for Stokes waves of steepness μ = ka ≃ 0.3 and μ = 0.1.

In the present article we perform similar experiment for waves of steepness μ ≃ 0.15. This experiment could be considered as a simulation of a realistic situation. If a typical steepness of the swell is 0.06 ÷ 0.07, in caustic area it could easily be two-three times more. In the new experiment, we start with the Stokes waver train, perturbed by a long wave with twenty time less amplitude. We observe development of modulation instability and finally, the explosive formation of the freak wave that is pretty similar to waves observed in natural experiments.

2 Basic equations

Suppose that incompressible fluid covers the domain

−∞ < y < η(x, t). (2.1)

The flow is potential, hence

V = ∇ϕ,  ∇V = 0,  ∇^{2}ϕ = 0.  (2.2)

Let ψ = ϕ|_{y=η} be the potential at the surface and H = T + U be the total energy. The terms

T = −\frac{1}{2} ∫_{−∞}^{∞} ψϕ_{n}dx, (2.3)

U = \frac{g}{2} ∫_{−∞}^{∞} η^{2}(x, t)dx, (2.4)

are correspondingly kinetic and potential parts of the energy, g is a gravity acceleration and ϕ_{n} is a normal velocity at the surface. The variables ψ and η are canonically conjugated; in these variables Euler equation of hydrodynamics reads

\frac{∂η}{∂t} = \frac{δH}{δψ} ,  \frac{∂ψ}{∂t} = −\frac{δH}{δη} .  (2.5)

One can perform the conformal transformation to map the domain that is filled with fluid,

−∞ < x < ∞,  −∞ < y < η(x, t),  z = x + iy

in z-plane to the lower half-plane

−∞ < u < −∞,  −∞ < v < 0,  w = u + iv

in w-plane. Now, the shape of surface η(x, t) is presented by parametric equations

y = y(u, t), x = x(u, t),

where x(u, t) and y(u, t) are related through Hilbert transformation

y = Ĥ (x(u, t) − u) ,  x(u, t) = u − Ĥy(u, t).  (2.6)

Here

Ĥ(f(u)) = P.V.\frac{1}{π} ∫_{−∞}^{∞} \frac{f(u′)du′}{u′ − u} .

Equations (2.5) minimize the action,

S = ∫ Ldt, (2.7)

L = ∫ ψ \frac{∂η}{∂t} dx − H. (2.8)

Lagrangian L can be expressed as follows,

L = ∫_{−∞}^{∞} ψ(y_{t}x_{u} − x_{t}y_{u})du + \frac{1}{2} ∫_{−∞}^{∞} ψĤψ_{u}du − \frac{g}{2} ∫_{−∞}^{∞} y^{2}x_{u}du+

+ ∫_{−∞}^{∞} f (y − Ĥ(x − u)) du .  (2.9)

Here f is the Lagrange multiplier which imposes the relation (2.6). Minimization of action in conformal variables leads to implicit equations (see Dyachenko et al. (1996))

y_{t}x_{u} − x_{t}y_{u} = −Ĥψ_{u}
ψ_{t}y_{u} − ψ_{u}y_{t} + gyy_{u} + Ĥ(ψ_{t}x_{u} − ψ_{u}x_{t} + gyx_{u}) = 0.  (2.10)

System (2.10) can be resolved with respect to the time derivatives. Omitting the details, we present only the final result

Z_{t} = iUZ_{u},
Φ_{t} = iUΦ_{u} − B + ig(Z − u).  (2.11)

Here

Φ = 2P̂ψ (2.12)

is a complex velocity potential, U is a complex transport velocity:

U = 2P̂ (\frac{−Ĥψ_{u}}{|z_{u}|^{2}} ) (2.13)

and

B = P̂ (\frac{|Φ_{u}|^{2}}{|z_{u}|^{2}} ) = P̂ (|Φ_{z}|^{2}) . (2.14)

In (2.12), (2.13) and (2.14) P̂ is the projector operator generating a function that is analytical in a lower half-plane

P̂(f) = \frac{1}{2} (1 + iĤ) f.

In the equations (2.11)

z(w) → w,  Φ(w) → 0,  at v → −∞.

All functions z, Φ, U and B are analytic ones in the lower half-plane v < 0.

Recently we found that equations (2.11) were derived in Ovsyannikov (1973), and we call them here Ovsyannikov’s equations, OE. Implicit equations (2.10) were not known until 1994, so we call them DKSZ-equations.

Note, that equation (2.10) can be used to obtain the Lagrangian description of surface dynamics. Indeed, from (2.10) one can get

Ψ = ∂^{−1}Ĥ(y_{t}x_{u} − x_{t}y_{u}) (2.15)

Plugging (2.15) into (2.8) one can express Lagrangian L only in terms of surface elevation. This result was independently obtained by Balk (1996). In Dyachenko (2001) equations (2.11) were transformed to a simple form, which is convenient both for numerical simulation and analytical study. By introducing of new variables

R = \frac{1}{Z_{w}} ,  and V = iΦ_{z} = i\frac{Φ_{w}}{Z_{w}}  (2.16)

one can transform system (2.11) into the following one

R_{t} = i(UR_{w} − RU_{w}),
V_{t} = i(UV_{w} − RB_{w}) + g(R − 1).  (2.17)

Now complex transport velocity U and B

U = P̂(V R̄ + V̄ R)
B = P̂(V V̄).  (2.18)

Thereafter, we will call equations (2.17), (2.18) Dyachenko equations, DE.

Both DKSZ-equations (2.10) and OE (2.11) have the same constants of motion

H = − ∫_{−∞}^{∞} ΨĤΨ_{u}du + \frac{g}{2} ∫_{−∞}^{∞} y^{2}x_{u}dy, (2.19)

the same total mass of fluid

M = ∫_{−∞}^{∞} yx_{u}du, (2.20)

and the same horizontal momentum

P_{x} = ∫_{−∞}^{∞} Ψy_{u}du. (2.21)

The Dyachenko equations (2.17), (2.18) have the same integrals. To express them in terms of R and V, one has to perform the integration

Z = u + ∫_{−∞}^{u} \frac{du}{R},  Φ = −i ∫_{−∞}^{u} \frac{V}{R}du.  (2.22)

3 Freak waves as a result of modulation instability

The Stokes wave is unstable with respect to long-scale modulation. This remarkable fact was first established in Lighthill (1965), who calculated a growth-rate of instability in the limit of long-wave perturbation. As far as Lighthill’s growth-rate coefficient was proportional to the wave number of perturbation length, the result was in principle incomplete: somewhere at short scales the instability must be arrested. The complete form of the growth-rate coefficient was found independently in Zakharov (1966), Zakharov (1968)) and in Benjamin and Feir (1967).

The presented technique based on the conformal mapping makes possible to study modulation instability in a very compact way. It is convenient to use the Dyachenko equation (2.17), (2.18). Let g = 1, k = 1. To study instability of the Stokes wave, propagating with the velocity c > 1, one has to go to the moving reference frame by the following change of variables:

u → u − ct, τ = t, R = 1 − \frac{iV}{c} + r. (3.1)

Then the Dyachenko equations take the form:

\frac{∂}{∂t} (r − \frac{iV}{c} ) + cr′ = i(Ũr′ − rŨ′)

\frac{∂V}{∂t} = i(V V′ − B′) − \frac{V}{c} B′ − \frac{ig}{c} V + gr + iŨV′  (3.2)

Ũ = p̄(r̄V + rV̄)

For the stationary progressive wave (with subscript 0) the following relation is valid:

R_{0} = 1 − \frac{iV_{0}}{c_{0}} . (3.3)

For the perturbation δV and δr one can derive linear system against the stationary solution

\frac{∂}{∂t} (r − i\frac{δV}{c} ) + cr_{u} = 0,  (3.4)

(\frac{∂}{∂t} + c\frac{∂}{∂u}) δV = i\frac{∂}{∂u}(V_{0}δV − δB) − \frac{V_{0}}{c} δB_{u} − \frac{ig}{c} δV + gr.

System (3.4) contains all information about stability of the Stokes wave.

The modulation instability is described by a perturbation presented as a sum of following harmonics:

δV, r ≡ e^{−iκu}, e^{(1±κ)u−inκ}, n = 1, …, κ < 1.

In the leading order of nonlinearity one can put

r = p_{1}e^{−i(1+κ)u} + p_{2}e^{−i(1−κ)u},
V = q_{1}e^{−i(1+κ)u} + q_{2}e^{−i(1−κ)u}.  (3.5)

Plugging (3.5) to (3.4) one obtains closed equations to p_{1}, p_{2}, q_{1}, q_{2}:

ṗ_{1} − \frac{i}{c} q̇_{1} = ic(1 + κ)p_{1},

˙p̄_{2} + \frac{i}{c} ˙ q̄_{2} = −ic(1 − κ)p_{2},  (3.6)

q̇_{1} − i (\frac{1}{c} − c(1 + κ)) q_{1} − p_{1} = V_{2}(1 + κ)q̄_{2},

˙q̄_{2} − i (\frac{1}{c} − c(1 − κ)) q̄_{2} − p̄_{2} = V_{2}(1 − κ) q_{1}.

Here V_{2} is the amplitude of second harmonics. Assuming p_{1}, q_{1}, p̄_{2}, q̄_{2} ≃ e^{i(Ω+κc)t}, one gets the following equation for Ω:

[(Ω − c)^{2} − 1 − κ] [(Ω + c)^{2} − 1 + κ] = (c^{2} − Ω^{2})(\frac{1}{c} − c)^{2}(1 − κ^{2}) (3.7)

To obtain this equation we put

|V_{2}|^{2} = (\frac{1}{c} − c)^{2} . (3.8)

This condition appears from the natural physical requirement: if κ = 0, then Ω = 0 is a solution of (3.7).

After simple calculations one can obtain dispersion relation for Ω:

Ω^{2} = \frac{1}{8} (−A^{2}κ^{2} + \frac{1}{8}κ^{4}) . (3.9)

Here A is the amplitude of the first Fourier harmonic of the Stokes wave train:

r = Ae^{−iu} + …

The result that was obtained by Zakharov (1966, 1968), and by Benjamin and Feir (1967). Lighthill in 1965 found long-wave asymptotic of the instability growth-rate,

Ω^{2} = −\frac{1}{8} A^{2} κ^{2}, (3.10)

with the maximum value of the growth-rate,

Ω^{2} = −\frac{1}{4} A^{2} κ^{2}, (3.11)

achieved at

κ^{2} = 4 A^{2}. (3.12)

The technique developed above makes possible to study the modulation and other instabilities with any arbitrary accuracy.

Here we study modulation instability of uniform wave train of Stokes wave. Question of great interest is the nonlinear stage of the instability. Here and everywhere below we do simulation in periodic domain L = 2π and

g = 1.

Wavetrain of the amplitude a with wavenumber k_{0} is unstable with respect to large scale modulation δk. Growth rate of the instability γ is

γ = \frac{ω_{0}}{2} ((\frac{δk}{k_{0}} )^{2} (ak_{0})^{2} − \frac{1}{4} (\frac{δk}{k_{0}} )^{4})^{\frac{1}{2}} . (3.13)

Here ω_{0} is the linear dispersion relation for gravity wave

ω_{0} = \sqrt{gk_{0}}.

• The shape of Stokes progressive wave is given by:

y = \frac{c^{2}}{2g} (1 − \frac{1}{|Z_{u}|^{2}} ),

while Φ is related to the surface as

Φ = −c(Z − u), V = ic(R − 1).

The amplitude of the wave \frac{h}{L} is the parameter for initial condition. (For the sharp peaked limiting wave \frac{h}{L} ≃ 0.141)

• Put 100 such waves with small perturbation in the periodic domain of 2π.

In such a way we prepared initial wave train with the steepness μ ≃ 0.095 Main Fourier harmonic of this wave train is k = 100. Similar problem was studied in Song and Banner (2002). But instead of long wavetrain they studied evolution of small group of waves.

For perturbation small value for Fourier harmonic with k_{p} = 1 was set. So, that

R_{k} = Runperturbed_{k} + 0.05R_{100} exp^{−ik_{p}u} .

Surface profile of this initial condition is shown in Figure 1

Fourier spectrum of the initial condition is shown in Figure 2 and Figure 3.

After sufficient large time, which is more than 1300 wave periods one can observe freak wave formation, as it is shown in Figure 4. Freak wave grows from mean level of waves to its maximal value for several wave periods, than vanishes or breaks.

Fig. 1. Initial profile of the wave train

Fig. 2. Fourier coefficients |R_{k}| for initial condition (μ ≃ 0.095).

Detailed view at the freak wave at the moment of maximal amplitude is shown in Figure 5. This set of experiments is similar to that of Dold and Peregrine (1986); Tanaka (1990). The difference is that we were able to increase the accuracy of the simulation, and consider much longer wavetrains. Also (due to using conformal mapping) we can simulate breaking with multivalued surface profile. Acccuracy in the simulation is very important because of the freak wave appears in a very subtle manner on the phase relations between

Fig. 3. Fourier coefficients |V_{k}| for initial condition (μ ≃ 0.095).

Fig. 4. Freak wave on the surface profile. T = 802.07

Fourier harmonics of the surface. Moreover, for shorter wavetrains threshold of modulation instability increases, and breaking doesn’t happen even for large steepness. In our experiments we have observed threshold of steepness for wave breaking a little less than in Tanaka (1990), but above μ = 0.1. Still, surface profile from Tanaka (1990) (Fig.5) is very similar to the picture in Figure 5 with μ = 0.095.

Fig. 5. Zoom in surface profile at T = 802.07

During numerical simulation of the final stage of freak wave formation, resolution must be increased to resolve high curvature of the surface profile. To do this we have ben increasing number of Fourier harmonics, which reached 2^{20} at the end (T = 802.07). Fourier coefficients of R_{k} are shown in Figure 6.

Fig. 6. Fourier coefficients |R_{k}| at T = 802.07.

If amplitude of the wave train is large, than freak wave may eventually break. Such a picture is presented in the Figure 7, which corresponds to the other numerical simulation with the initial steepness μ ≃ 0.14.

Fig. 7. Profile of breaking wave.

4 Exact equations and nonlinear Shrödinger approximation

Evolution of weakly nonlinear Stokes wavetrain can be described by nonlinear Shrödinger equation (NLSE), derived by Zakharov (1968). This equation is integrable (see Zakharov and Shabat (1972)) and is just the first term in the hierarchy of envelope equations describing packets of surface gravity waves. The second term in this hierarchy was calculated by Dysthe (1979), the next one was found a few years ago by Trulsen and Dysthe (1996). The Dysthe equation was solved numerically by Ablovitz and his collaborates (see Ablovitz et al. (2000 and 2001)).

Since the first work of Smith (1976), many authors tried to explain the freak wave formation in terms of NLSE and its generalizations, like Dysthe equation. A vast scientific literature is devoted to this subject. The list presented below is long but incomplete: Ablovitz et al. (2000 and 2001), Onorato et al. (2000a), Onorato et al. (2000b), Onorato et al. (2001), Onorato et al. (2002), Peregrine (1983), Peregrine et al. (1988), Trulsen and Dysthe (1996), Trulsen and Dysthe (1997), Trulsen (2000), Trulsen et al. (2000), Clamond and Grue (2002).

One cannot deny some advantages achieved by the use of the envelope equations. Results of many authors agree in one important point: nonlinear development of modulation instability leads to concentration of wave energy in a small spatial region. This is a ”hint” regarding possible formation of freak wave. On the other hand, it is clear that the freak wave phenomenon cannot be explained in terms of envelope equations. Indeed, NLSE and its generalizations are derived by expansion in series on powers of parameter λ ≃ \frac{1}{Lk}, where k is a wave number, L is a length of modulation. For real freak wave λ ∼ 1 and any ”slow modulation expansion” fails. At this point interesting question rises: what happens to NLSE approximation when increasing the steepness of the carrier wave? In particular, we study ”exact” soliton solutions for NLSE placed in the exact equations (2.17).

Such type of problem was considered in the Henderson et al. (1999), but with low resolution, and small length of periodic carrier. Also in Clamond and Grue (2002) numerical solutions for envelope equation was compared with ”almost” exact equations.

For the equations (2.17) NLSE model can be derived for the envelope of R.

R = 1 + R_{1}e^{−ik_{0}u−ω_{0}t} + …

iR_{1}_{t} + \frac{1}{8} \frac{ω_{0}}{k_{0}^{2}} R_{1}_{uu} + \frac{1}{2}ω_{0}k_{0}^{2}|R_{1}|^{2}R_{1} = 0.

Initial conditions consist of ”linear wave carrier” e^{−ik_{0}u}, modulated in accordance with soliton solution for NLSE:

R(u) = 1 + s_{0} \begin{matrix}e & −ik_{0}u \\ cosh & (λk_{0}u) \end{matrix},
V (u) = −ic_{0}s_{0} \begin{matrix}e & −ik_{0}u \\ cosh & (λk_{0}u) \end{matrix}.  (4.14)

Here s_{0} is the steepness of the carrier wavetrain, c_{0} - phase velocity of the carrier.

First comparison of fully nonlinear model for water wave with NLSE was done in Clamond and Grue (2002) for the wave carrier with the steepness μ ≃ 0.091. For such steepness there was a good agreement between two models, but only for the short time. After finite time weakly nonlinear model (NLSE) ceases to be valid.

In our work we want to study the situation with larger and smaller stepness, to find out how NLSE approximation breaks.

4.1 Small steepness

First experiment was intended to observe how NLSE works. In the initial conditions (4.14) we used

s_{0} ≃ 0.07, λ = 0.1, k_{0} = 100.

Initial surface of fluid is shown in Figure 8.

Fig. 8. Initial surface profile like for NLSE soliton with μ ≃ 0.07.

After couple of thousands wave periods soliton changes a little, as it is seen in Figure 9:

Fig. 9. Surface profile like for NLSE soliton with μ ≃ 0.07 at T=1500.

Also in the Figure 10 and Figure 11 Fourier spectra of the soliton at both moments of time are presented.

Fig. 10. Fourier harmonics of the initial soliton with μ ≃ 0.07.

Fig. 11. Fourier harmonics of the soliton with μ ≃ 0.07 at T=1500.

So, one can can see that for the steepness μ ≤ 0.07 NLSE model is quite reasonable.

Another numerical experiment showing effective simulation with equations (2.17) along with applicability NLSE model for moderate steepness, μ ≃ 0.085, is the collision of two solitons.

In the Figure 12 initial condition is shown:

Fig. 12. Initial surface profile of two NLSE solitons with μ ≃ 0.085.

Moment of collision is shown in the Figure 13:

Fig. 13. Two NLSE solitons with μ ≃ 0.085. collide at T=30.8

and detailed view showing carrier wavetrain under the envelope is in the Figure 14.

After second collision (recall that boundary conditions are periodic) solitons are plotted in the Figure 15:

Fig. 14. Detailed view of two colliding NLSE solitons with μ ≃ 0.085 at T=30.8.

Fig. 15. Two NLSE solitons with μ ≃ 0.085. after two collisions at T=250.0.

Fourier spectra of these two solitons at the moments of time T = 0.05, 30.8, 250.0 are shown in Figure 16, 17, 18.

Fig. 16. Fourier spectrum of the initial surface profile of two NLSE solitons with μ ≃ 0.085.

Fig. 17. Fourier spectrum of two colliding NLSE solitons with μ ≃ 0.085. at T = 30.8

4.2 Large steepness

Now let’s turn to the higher steepness of the carrier,

μ = 0.1.

Fig. 18. Fourier spectrum of two NLSE solitons with μ ≃ 0.085. at T = 250.0

In the Figure 19 there is initial condition:

Fig. 19. Initial surface profile like for NLSE soliton with μ ≃ 0.10.

Again, after couple of thousands wave periods soliton changes a little, as it is seen in Figure 20:

In the Figure 21 and Figure 22 Fourier spectra of the soliton at both moments of time are presented.

Fig. 20. Surface profile like for NLSE soliton with μ ≃ 0.10 at T=2345.

Fig. 21. Fourier harmonics of the initial soliton with μ ≃ 0.10.

From this pictures one can see that for steepness μ ≃ 0.10 some corrections to the NLSE model are desirable. Dysthe equations are exactly intended for that situation.

But what happens when further increasing the steepness? Below we consider the case of the steepness of the carrier

μ = 0.14.

Fig. 22. Fourier harmonics of the soliton with μ ≃ 0.10 at T=2345.

In the Figure 23 there is initial condition:

Fig. 23. Initial surface profile like for NLSE soliton with μ ≃ 0.14.

Very fast, after couple of dozen wave periods soliton drastically changes, as it is seen in Figure 24: One can see freak wave at the surface (in Figure 25):

Fig. 24. Surface profile like for NLSE soliton with μ ≃ 0.14 at T=38.4.

Fig. 25. Zoomed surface profile near freak wave μ ≃ 0.14 at T=38.4.

In the Figure 26 and Figure 27 Fourier spectra of the soliton at both moments of time are presented. They demonstrate the quality of the numerical simulation.

From the last case, with the steepness μ = 0.14, one can see that envelope approximation completely fails. Such event as one single crest (freak wave) can not be described in terms of wave envelope.

Fig. 27. Fourier harmonics of the soliton with μ ≃ 0.14 at T=38.4.

5 Do freak waves appear from quasisolitonic turbulence?

Let us summarize the results of our numerical experiments. Certainly, they reproduce the most apparent features of freak waves: single wave crests of very high amplitude, exceeding the significant wave height more than three times, appear from ”nowhere” and reach full height in a very short time, less than ten periods of surrounding waves. The singular freak wave is proceeded by the area of diminished wave amplitudes. Nevertheless, the central question about the physical mechanism of freak waves origin is still open.

In our experiments, the freak wave appears as a result of development of modulation instability, and it takes a long time for the onset of instability to create a freak wave. Indeed, the level of perturbation in our last experiment is relatively high. The two-three inverse growth-rate is enough to reach the state of full-developed instability, when the initial Stokes wave is completely decomposed. Meanwhile, the freak wave appears only after fifteenth inverse growth-rates of instability. What happens after developing of instability but before formation of freak wave?

During this relatively long period of time, the state of fluid surface can be characterized as quasisolitonic turbulence, that consists of randomly located quasi-solitons of different amplitudes moving with different group velocities. Numerical study of interaction of envelope soliton was done in Clamond and Grue (2002). Such interaction leads to formation of wave with large amplitude. Here we can think in term of quasisolitonic turbulence. Such turbulence was studied in the recent work of Zakharov, Dias and Pushkarev (Zakharov et al. (2004)) in a framework of so-called defocusing MMT model.

i\frac{∂Ψ}{∂t} = |\frac{∂}{∂x}|^{1/2} Ψ + |\frac{∂}{∂x}|^{3/4} (||\frac{∂}{∂x}|^{3/4} Ψ|^{2} |\frac{∂}{∂x}|^{3/4} Ψ) (5.1)

This is a heuristic model description of gravity surface waves in deep water. In this model, quasi-solitons of small amplitude are stable, interact inelastically and can merge. Above some critical level quasi-solitons of large amplitude are unstable. They collapse in finite time forming very short wave pulses, which can be considered as models of freak waves. Equation 5.1 has the exact solution:

Ψ = A e^{ikx−iωt}
ω = k^{1/2} (1 + k^{5/2} A^{2}) .  (5.2)

This solution can be constructed as a model of the Stokes wave and is unstable with respect to modulation instability. Development of this instability was studied numerically. On the first stage, the unstable monochromatic wave decomposes to a system of almost equal quasi-solitons. Then, the quasisolitonic turbulence is formed: quasi-solitons move chaotically, interact with each other, and merge. Finally they create one large quasi-soliton, which exceeds threshold of instability and collapses, creating a freak wave.

One can think that a similar scenario of freak wave formation is realized in a real sea. We like to stress that the key point in this scenario is the quasisolitonic turbulence and not the Stokes wave. The Stokes wave is just a ”generator” of this turbulence. The quasisolitonic turbulence can appear as a result of instability of narrow spectral distributions of gravity waves.

The formulated above concept is so far a hypothesis, which has to be confirmed by future numerical experiments.

Acknowledgments

This work was supported by the US Army Corps of Engineers Grant W912HZ-05-P-0351, by ONR Grant N00014-03-1-0648, NSF Grant DMS 0404577, RFBR Grant 06-01-00665, the Program ”Fundamental Problems in Nonlinear Dynamics” from the RAS Presidium, and Grant ”Leading Scientific Schools of Russia”.

References

Ablovitz, M.I., Hammack, D., Henderson, J. and Scholder, C.M., Modulated periodic Stokes waves in deep water, Phys. Rev. Lett. 2000, 887-890; Ablovitz, M.I., Hammack, D., Henderson, J. and Scholder, C.M., Long-time dynamics of the modulational instability of deep water, Physica D 152-153 (2001) 416-433.

Balk, A., A Lagrangian for water waves, Phys. Fluids, 8 (1996) 416-420. Banner, M. and Tian, X., On the determination of the onset of wave breaking

for modulating surface gravity waves, J. Fluid Mech. 367 (1998) 107-137. Benjamin, T.B. and Feir, J.E., The disintegration of wave trains on deep water.

Part 1. Theory. J. Fluid Mech. 27 (1967) 417-430. http://bell.mma.edu/ achase/NS-221-Big-Wave.html Chalikov, D. and Sheinin, D., Direct modeling of one-dimensional nonlinear

potential waves, in Nonlinear Ocean Waves, ed. W. Perrie, Advances in Fluid Mechanics, Vol 17, Computation Mechanics Publications, (1998) 89-110.

D. Clamond, J. Grue, A fast method for fully nonlinear-wave computations,

J. Fluid Mech. 447 (2001) 337-355. D. Clamond, J. Grue, Interaction between envelope solitons as a model for

freak wave formation. Part I: Long time interaction, C.R.Mecanique 330 (2002) 575-580.

Dean, R.G., Freak Waves: A possible explanation. In Water Wave Kinetics

(ed. A. Torum and O. T. Gudmestad), Kluwer, (1990), pp. 609-612. Divinsky, V.D., Levin, B.V., Lopatikhin, L.I., Pelinovsky, E.N. and Slyunyaev, A.V., A freak wave in the Black sea: Observations and simulation. Doklady Earth Science, 395 (2004) 438-443.

Dold, J.W., An efficient surface integral algorithm applied to unsteady gravity

waves. J. Computat. Phys. 103 (1992) 90-115. Dold, J.W. and Peregrine, D.H. Water-wave modulation, In Proc. 20th Intl

Conf. on Coastal Engng. Vol. 1, Chap. 13, 1986, pp. 163-175. Dommermuth, D.G. and Yue, D.K.P. A high-order spectral method for the

study of nonlinear gravity waves, J. Fluid Mech. 184 (1987) 267-288. Dyachenko, A.I., E. Kuznetsov, E.A., Spector, M.D. and Zakharov, V.E., An-

alytical description of the free surface dynamics of an ideal fluid (canonical formalism and conformal mapping, Phys. Lett. A, 221 (1996) 73-79.

Dyachenko, A.I. On the dynamics of an ideal fluid with free surface, Doklady

Mathematics, 63 (2001) 115-118. Dyachenko, A.I. and Zakharov, V.E., Modulation Instability of Stokes Wave

→ Freak Wave, JETP Letters, 81(6) (2005) 255-259. Dyachenko, A.I., Instability in the Numeric Models of Free Surface Hy-

drodynamics, International conference ”Frontiers of Nonlinear Physics”, St.Petersburg-Nizhny Novgorod, Russia, August 2-9, 2005.

Dysthe, K.B., Note on a modification to the nonlinear Schrodinger equation

for application to deep water waves, Proc. Roy. Ser. A 369 (1979) 105-114. Earle, M.D., Extreme wave conditions during Hurricane Camille, J. Geophys.

Res. 80 (1975) 377-379. Feir, J.E. Discussion: Some results from wave pulse experiments, Proc. R. Soc.

Lond. A 299 (1967) 54-58. Gerber, M., The Benjamin-Feir instability of a deep water Stokes wave packet

in the presence of a non-uniform medium, J. Fluid Mech. 176 (1987) 311-332.

Gerber, M., The interaction of deep water gravity waves and an annular cur-

rent: linear theory, J. Fluid Mech. 248 (1993) 153-172. Gerber, M., Giant waves and the Agulhas current, Deep sea Res. 1996. Gutshabash, Y.S. and Lavrenov, I.V., Swell transformation in the Cape Ag-

ulhas current, Izv. Atmos. Ocean Phys. 22 (1986) 494-497. Henderson, K.L. D., Peregrine, D.H. and Dold, J.W. Unsteady water wave

modulations: fully nonlinear solutions and comparison with nonlinear Schrodinger equation, Wave Motion, 29 (1999) 341-361.

Irvine, D.E., Tilley, D.G., Ocean wave directional spectra and wave-current interaction in the Agulhas from the shuttle imaging radar-B synthetic aperture radar. J. Geophys. Res. 93 (1988) 15389-15401.

Kharif, C., and Pelinovsky, E., Physical mechanisms of the rogue wave phe-

nomena, Europ. J. Mech. B/Fluids, 22 (2003) 603-634. Lavrenov, I.V., The wave energy concentration at the Agulhas current of South

Africa, Natural Hazards 17 (1998) 117-127. Lighthill, M.J., Contribution to the theory of waves in nonlinear dispersive

systems, J. Inst. Math. Appl. 1 (1965) 269-306. Majda, A., McLaughlin, D. and Tabak, E., A One-Dimensional Model for

Dispersive Wave Turbulence, J. Nonlinear Science 7 (1997) 9-44. Mallory, J.K., Abnormal waves on the south-east of South Africa, Inst. Hydrog.

Rev. 51 (1974) 89-129. Meison, D., Orzag, S. and Izraely, M. Applications of numerical conformal

mapping, J. Computational Physics, 40 (1981) 345-360. Mori, N., Liu, P.C.,and Yasuda, T., Analysis of freak wave measurements in

the sea of Japan, Ocean Engrg. 29 (2002) 1399-1414. Onorato, M., Osborne, A.R., Serio, M. and Damiani, T., Occurence of freak

waves from envelope equation in random ocean wave simulations, in Rogue waves 2000: Brest, France, November 2000, eds. M. Olagnon and G.A. Athanassoulis, Ifremer, 2001, 181-192.

Onorato, M., Osborne, A.R. and Serio, M., The nonlinear dynamics of rogue waves and holes in deep-water gravity wave trains, Phys. Lett. A 275 (2000) 386-393.

Onorato, M., Osborne, A.R., Serio, M. and Bertone, S., Freak waves in random

oceanic sea states, Phys. Rev. Letters 86 (2001) 5831-5834. Onorato, M., Osborne, A.R. and Serio, M., Extreme wave events in directional,

random oceanic sea states, Phys. of Fluids, 14 (2002) L25-L28. Ovsyannikov, L.V. Dinamika sploshnoi sredy, M.A. Lavrent’ev Institute of

Hydrodynamics Sib. Branch USSR Ac. Sci., No.15, (1973) pp. 104-125, (in Russian).

Peregrine, D.H., Interaction of water waves and currents, Adv. Appl. Mech.

16 (1976) 9-117. Peregrine, D.H., Water waves, nonlinear Schrodinger equations and their so-

lutions, J. Austral. Math. Soc. B 25 (1983) 16-43. Peregrine, D.H., Skyner, D., Stiassnie, M. and Dold, N., Nonlinear effects on

focused water waves, Proc. 21th Intl Conf. on Coastal Engng. Vol. 1, Chap. 54, 1988, pp. 732-742.

Sand, S.E., Hansen, N.E., Klinting, P., Gudmestad, O.T., Sterndorff, M.J., Freak wave kinematics, in Water Wave Kinematics, eds. A. Torum, O. T. Gudmestad, Kluwer Academic, Doedrecht, 1990, pp. 535-549.

Shamin, R.V., On the Existence of Smooth Solutions to the Dyachenko Equations Governing Free-Surface Unsteady Ideal Fluid Flows, Doklady Mathematics 73(1) (2006) 112-113.

Smith, R., Giant Waves, Fluid Mech. 77 (1976), 417-431. Song, J. and Banner, M.L., On determining the onset and strength of breaking

for deep water waves. Part 1: Unforced irrotational wave groups. J. Phys. Oceanogr. 32 (2002) 2541-2558.

Tanaka, M., Maximum amplitude of modulated wavetrain, Wave motion 12

(1990) 559-568. Trulsen, K. and Dysthe, K.B., A modified nonlinear Schrodinger equation for

broader bandwidth gravity waves on deep water, Waves motion 24 (1996) 281-289.

Trulsen, K. and Dysthe, K.B., Freak waves - a three-dimensional wave simulation, in Proc. 21st Symposium on Naval Hydrodynamics, 1997, 550-558;

http://www.nap.edu/books/0309058791/html/550.html Trulsen, K., Simulating the spatial evolution of a measured time series of a

freak wave, in Rogue waves 2000: Brest, France, November 2000, eds. M. Olagnon and G.A. Athanassoulis, Ifremer, 2001, 265-274.

Trulsen, K., Kliakhandler, I., Dysthe, K.B. and Velarde, M.G. On weakly nonlinear modulation of waves on deep water, Phys. Fluids 12 (2000) 2432-2437.

West, B., Brueckner, K., Janda, R., Milder, D. and Milton, R. (1987). A new numerical method for surface hydrodynamics. Journal of Geophysical Research C 92(11) (1987) 11803-11824.

White, B.S., Fornberg, B., On the chance of freak waves at sea, J. Fluid Mech.

355 (1998) 113-138. Zakharov, V.E., Stability of nonlinear waves in dispersive media, J. Teor. Prikl.

Fiz. 51 (1966) 668-671 (in Russian); Sov.Phys. JETP 24 (1967) 455-459. Zakharov, V.E., Stability of periodic waves of finite amplitude on a surface of

deep fluid, J. Appl. Mech. Tech. Phys. 9 (1968) 190-194. Zakharov, V.E. and Shabat, A.B., Exact theory of two-dimensional self-

focusing and one-dimensional self-modulation of waves in nonlinear media, Soviet Physics JETP 34 (1972) 62-69.

Zakharov, V.E. and Kuznetsov, E.A., Optical solitons and quasisolitons, JETP

86 (1998) 1035-1046. Zakharov, V.E., Dyachenko A.I. and Vasilyev, O.A., New method for numer-

ical simulation of nonstationary potential flow of incompressible fluid with a free surface, Eur. J. of Mech. B-Fluids, 21 (2002) 283-291.

Zakharov, V.E., Dias, F. and Pushkarev, A.N., One-Dimensional Wave Turbulence, Phys. Reports 398 (2004) 1-65.