## Abstract

Theoretical understanding of the growth of wind-driven surface water waves has been based on two distinct mechanisms: growth due to random atmospheric pressure fluctuations unrelated to wave amplitude and growth due to wave coherent atmospheric pressure fluctuations proportional to wave amplitude. Wave-independent random pressure forcing produces wave growth linear in time, while coherent forcing proportional to wave amplitude produces exponential growth. While observed wave development can be parameterized to fit these functional forms and despite broad agreement on the underlying physical process of momentum transfer from the atmospheric boundary layer shear flow to the water waves by atmospheric pressure fluctuations, quantitative agreement between theory and field observations of wave growth has proved elusive. Notably, wave growth rates are observed to exceed laminar instability predictions under gusty conditions. In this work, a mechanism is described that produces the observed enhancement of growth rates in gusty conditions while reducing to laminar instability growth rates as gustiness vanishes. This stochastic parametric instability mechanism is an example of the universal process of destabilization of nearly all time-dependent flows.

## 1. Introduction

It has been found in a variety of field studies that the growth rate of wind-driven water waves depends on the stratification of the surface boundary layer and the related gustiness of the wind. For example, Kahma and Calkoen (1992) found that water wave growth rates differ depending on the boundary layer stratification and that the growth rate substantially increases under unstable conditions. They attributed this increase in growth rate to increased levels of gustiness in unstably stratified boundary layers. However, when Janssen and Komen (1985) investigated modal growth rates that result under conditions of unstable stratification, they obtained the surprising result that unstable stratification decreases the predicted growth rate of long water waves.

While this enhancement of growth rate in gusty conditions is now generally acknowledged, the mechanism responsible for the increase remains a matter of debate (Abdalla and Cavaleri 2002). Modal mechanisms for producing increased wave growth rates in fluctuating winds were investigated in the studies of Nikolayeva and Tsimring (1986), Komen et al. (1994; see Jansenn therein), Miles (1997), and Miles and Ierley (1998). In his investigation, Komen et al. (1994; see Jansenn therein) considered the effect of a normally distributed friction velocity on growth rate using the modal growth rate parameterization of Snyder et al. (1981). Because growth rate varies with friction velocity in this parameterization, being zero for negative friction velocities and positive for positive friction velocities, an increase in the average growth rate proportional to rms friction velocity variance results when friction velocity fluctuation symmetrically distributed about its mean value is assumed. Moreover, according to Cavaleri and Burgers (1992) the velocity fluctuations at the 10-m level rather than friction velocity are normally distributed. This leads to skewness in the distribution of friction velocity and to further enhancement of the growth rate. In general, such mechanisms produce increased growth rates under gusty conditions owing to the form of the dependence of modal growth rate on wind speed. Cavaleri and Burgers introduced these effects into a wave-height prediction model and found faster wave growth and higher waves in the fully developed sea, the increase being up to 30% in wave height. It was further found that model results improved when compared to observations with the introduction of this gustiness parameterization (Komen et al. 1994). Ponce and Ocampo-Torres (1998) introduced high-frequency wind observations made in the Gulf of California into a third-generation wave height prediction model and found an increase in the resulting wave height. Similar conclusions were obtained by Bauer and Weisse (2000). An extended investigation of the effect of gustiness on wave height was undertaken by Abdalla and Cavaleri (2002) who demonstrated that the increase in the resulting wave fields was larger in the Atlantic Ocean (where gustiness is more pronounced) than in the Mediterranean (where gustiness is comparatively reduced). Increase in growth rates in these models arises from the increase in the average of the instantaneous growth rate and also from the convexity of the exponential function relating wave growth rate to wave growth.

In this work we propose a new mechanism for enhancing wave growth in gusty conditions that involves the nonnormality of the full spectrum of the underlying dynamical operator of the coupled air–sea system. This mechanism is distinct from the aforementioned mechanisms for producing enhanced growth rate under gusty conditions based on the functional dependence of modal growth rate on wind speed (Cavaleri and Burgers 1992; Komen et al. 1994; Miles 1997; Janssen 2004). The enhanced growth mechanism that we propose can be viewed as an incoherent parametric instability, and it exemplifies the universal instability arising from the intrinsic nonnormality of time-dependent dynamical operators (Zel’Dovich et al. 1984; Farrell and Ioannou 1996b, 1999; Poulin et al. 2003; Hadjiioannou et al. 2006).

To put the stochastic mechanism in context, we now briefly review theories proposed to explain the phenomenon of wind-driven growth of surface water waves. These theories in general involve either wave incoherent stochastic forcing by random atmospheric pressure fluctuations (Eckart 1953; Phillips 1957) or wave coherent forcing by wave-induced atmospheric pressure fluctuations (Helmholtz 1868; Kelvin 1871; Jeffreys 1925, 1926; Miles 1957, 1959a, b). These mechanisms together with a parameterization of nonlinear interactions were incorporated by Hasselmann (1960) into a general prediction equation for the spatial and temporal distribution of the spectral energy of the waves. Various empirical versions of this prediction equation are presently used for operational wave forecasting purposes (WAMDI Group 1988; Komen et al. 1994; Tolman and Chalikov 1996; Booij et al. 1999; Janssen 2007). These prediction equations incorporate a parameterized variance growth linear in time to account for wave incoherent forcing and a parameterized growth exponential in time to account for wave-induced coherent forcing. Coherent forcing by wave-induced atmospheric pressure fluctuations proportional to wave amplitude is generally accepted as necessary to account for observed rates, but the growth rate is parameterized because growth rates obtained from laminar modal instability theory are too small by up to an order of magnitude for long waves (Abdalla and Cavaleri 2002; Janssen 2007). Despite recent progress and while these parameterized spectral prediction models provide useful predictions of the sea state, they remain loosely related to stability theory and are widely acknowledged to need improvement (Komen et al. 1994; Cavaleri 2006).

The surface water wave generation problem can be thought of as a shear stability problem in the presence of a flexible lower boundary. Perhaps the most familiar example of this class of stability problems is Kelvin–Helmholtz instability resulting from Bernoulli suction coherently 180° out of phase with surface elevation. However, Kelvin–Helmholtz instability requires that this suction exceed the gravitational restoring force, in turn requiring wind speeds in excess of those observed to be associated with wave generation. It follows that Kelvin–Helmholtz instability is not generally regarded as an important mechanism for water wave generation.

The instability mechanisms of Jeffreys (1925, 1926), Miles (1957, 1959a, b), and Belcher and Hunt (1993) result from the positive momentum flux from the atmosphere to the water wave that occurs when a component of atmospheric perturbation pressure lags surface water elevation by 90°. The Jeffreys theory postulates that turbulent flow separation produces the required atmospheric pressure to wave height phase lag, while Belcher and Hunt attribute the phase shift to a nonseparated sheltering effect. Similar to nonseparated sheltering, the Miles theory uses a linear laminar theory that has no direct role for turbulence in producing the phase lag required for instability.

Laminar instability theory takes account of the turbulent nature of atmospheric boundary layer pressure fluctuations only in that the fixed phase lag between atmospheric pressure fluctuations and wave height arising in the theory can be interpreted as the statistical average of the actual fluctuating pressure. The attraction of the laminar instability mechanism is that it predicts an average component of atmospheric pressure fluctuation systematically proportional to wave amplitude and, as a result, growth exponential in time. In contrast, wave incoherent turbulence produces pressure fluctuations unrelated to wave amplitude and therefore to variance growth linear in time. Even coherent resonant pressure fluctuations produce variance growth, at most, quadratic in time (Phillips 1957).

There is no doubt as to the existence in laminar flow over a flexible boundary of the instability described by Miles (1957, 1959a, b). However, theory and observation could be brought into better accord if there existed in a turbulent shear flow over a flexible boundary another instability mechanism exponential in time and substantially dominant in growth rate when compared with the laminar mechanism. A comparison of observations, laboratory experiments, numerical simulations, and analytic theories is shown in Fig. 1.^{1} Inspection of Fig. 1 suggests that there are two regimes of growth occurring in the wave generation problem: one in which observed growth rates roughly follow the predictions of the Miles theory, while lying in the mean above these predictions, and another branch in which the observed growth rate is higher and follows a functional form that appears unrelated to extant theory. When different regimes such as these are found in a physical problem, it is often indicative of a missing variable in the analysis. Indeed, it has been frequently remarked that the missing variable in this case is the stability of the airflow (Kahma and Calkoen 1992). The strongly growing wave regime is characteristically found in cold air flowing over warm water, which produces highly turbulent and gusty conditions (Barnett and Wilkerson 1967; Stewart and Teague 1980), while the more slowly growing wave regime is associated with less turbulent laboratory experiments, numerical simulations, and theoretical calculations and also observationally with more stable or, at most, neutrally stable tropical marine boundary layer flows (Sullivan and McWilliams 2002; Belcher and Hunt 1993; Snyder et al. 1981).

In this work we describe a mechanism that uses wave incoherent atmospheric turbulence to produce an increase over that predicted by laminar instability theory in the statistically coherent component of atmospheric pressure fluctuations proportional to and properly phased with wave height to produce exponential growth.

We begin by reviewing the formulation of the stability problem for the growth of water waves. Then we demonstrate that the stability problem is highly nonnormal and identify initial perturbations that produce water waves with greatest amplitude. The universal mechanism of destabilization time-dependent flows depends on this nonnormality of the fluctuating operator (Farrell and Ioannou 1996b, 1999). Having demonstrated the existence of high nonnormal growth in the wave generation problem, we proceed to the description of the time-dependent stability problem. Then we introduce the stochastic time-dependent instability mechanism and apply it to the wave generation problem for a simplified Kelvin–Helmholtz approximation to boundary layer flow and, also, for a more realistic logarithmic boundary layer profile.

## 2. The linear stability problem for surface water waves

Consider an inviscid incompressible atmosphere of constant density in which the mean wind, *U* (*z*, *t*), in the *x* direction varies in the vertical direction *z* and with time *t*. Assume that the mean wind at *z* = 0 is *U*_{0}(*t*). Harmonic perturbations in the atmosphere with *x* wavenumber *k* and streamfunction *ψ _{a}* (

*z*,

*t*)

*e*are governed by the equation

^{ikx}where *D* ≡ ∂/∂*z*, ∇^{2} ≡ *D*^{2} − *k*^{2} is the Laplacian operator, and *U* ″(*z*, *t*) ≡ *D*^{2}*U*(*z*, *t*) denotes the curvature of the mean wind profile. Equation (1) is assumed to be valid in the region *z* > 0 occupied by the atmosphere, which has density *ρ _{a}*. The streamfunction

*ψ*is assumed to decay to 0 as

_{a}*z*→ ∞.

A semi-infinite incompressible fluid of density *ρ _{w}* occupies the region

*z*< 0. This fluid models the water: it is assumed to have no motion other than that associated with irrotational small amplitude surface waves

^{2}with streamfunction form:

in which *ψ*^{0}_{w}(*t*) is the streamfunction of the water at the mean air–water interface, *z* = 0. The wave height, which is also assumed to take the harmonic form *η* = η̂(*t*)*e*^{ikx}, is a material boundary and satisfies in the small elevation limit the relations

where *ψ*^{0}_{a}(*t*) ≡ *ψ*_{a}(0, *t*). Subtracting the two conditions (3) and (4) we obtain that the streamfunctions in the air and water are given by

Streamfunction (2) satisfies the momentum equations for *z* < 0 with time dependence obtained by imposing (5) and continuity of normal stress, which on linearization at *z* = 0 takes the form

where *p _{w}*,

*are the values of the Fourier coefficients of the pressure on the water and atmospheric side of the interface at*

_{a}*z*= 0. In boundary condition (6) surface tension has been neglected.

Using the linearized momentum equation in the *x* direction in both the atmosphere and water layers, the normal stress boundary condition at *z* = 0 can be expressed as

where *α*(*t*) is the wind shear at *z* = 0 and *ε* = *ρ _{a}*/

*ρ*; for the water–air interface,

_{w}*ε*≈ 0.001.

## 3. Nonnormal perturbation growth for time-independent wind

We first consider the case of time-independent wind *U*(*z*) and *U*(0) = 0 at the water surface *z* = 0.

### a. Formulation of the normal-mode solution

With the assumption of time independence of the background wind, eigenstates of (1) with boundary conditions (3) and (7) are found by assuming solution form *ψ _{a}*(

*z*,

*t*) =

*ψ*(

*z*)

*e*

^{−}

*and solving the resulting eigenvalue problem:*

^{ikct}with boundary conditions *ψ* → 0 as *z* → ∞ and

at *z* = 0. This eigenproblem admits an analytic eigenfunction *ψ* with associated eigenvalue *c*, which can be identified with an unstable surface water wave. From the eigenvalue both growth rate, *kc _{i}* =

*k*ℑ(

*c*), and phase speed,

*c*= ℜ(

_{r}*c*), of the surface water wave can be determined. This eigenproblem can be solved using numerical integration and the shooting method.

Alternatively, following Miles (1957), the growth rate can be obtained by exploiting the smallness of the density ratio *ε* to expand the streamfunction and eigenvalue in a regular perturbation series:

Using these expressions the zero-order balance in (8) and (9) immediately yields the phase speed *c*^{2}_{0} = *g*/*k*, implying that the zero-order phase speed of the interface wave is the phase speed of surface gravity waves on a water surface bounded by a motionless atmosphere of vanishing density. The *O*(*ε*) correction to the phase speed is^{3}

in which *c*_{0} = *g*/*k* is the phase speed of surface gravity waves in the absence of an atmosphere, and *Dψ*_{0}/*ψ*_{0} at *z* = 0 is obtained by solving the Rayleigh equation

with specification that the regular singularity of the differential equation at the critical level *U*(*z*_{c}) = *g*/*k* be resolved by selecting the solution that results in the limit in which the phase speed *c*_{0} has a vanishing positive imaginary part (Drazin and Reid 1981). In that case the Wronskian

is constant except at the critical level *z _{c}* where it undergoes a jump. Because

*W*(∞) = 0, the value of the Wronskian at the interface can be obtained:

which then determines the *O*(*ε*) growth rate:

The last equality, obtained by Miles (1957), implies that all profiles with negative curvature (such as occur in the atmospheric boundary layer) are unstable with growth rate *O*(*ε*), even if the velocity profile is not inflectional. A physical explanation of the instability mechanism was given by Lighthill (1962). The unstable waves are prograde (moving in the direction of the flow) with phase speed *c* = *g*/*k* + *O*(*ε*). It is noteworthy that, according to the above solution, no exponentially growing water waves exist for mean winds with zero or positive curvature, such as Couette flow, and presumably there would be no surface waves generated by such a flow regardless of the wind speed. Moreover, the predicted growth rates can be shown to decrease under unstable boundary layer stratification in contrast to the observed increase (Janssen and Komen 1985).

### b. Normal-mode solution for the logarithmic boundary layer flow

We now present numerical calculations of the growth rate, phase speed, and structure of the eigenmode for the logarithmic boundary layer flow

shown in Fig. 2. In (17) *u*_{*} is the friction velocity, *K* = 0.42 is the von Kármán constant, and *z*_{0} = *κ**u*^{2}_{*}/*g* is the roughness length expressed in terms of the Charnock constant (*κ* = 0.0144), the friction velocity *u*_{*}, and the acceleration of gravity *g*. This profile is not inflectional but supports unstable prograde surface waves with the numerically obtained growth rates shown in Fig. 3 and phase speeds to *O*(*ε*) equal to *g*/*k*, the phase speed of the surface gravity wave in the absence of an atmosphere. The structure in the atmosphere associated with a typical unstable surface wave is shown in Fig. 4. The eigenfunction leans against the shear (see also Fig. 9), extracting energy from the mean flow by interaction mediated by Reynolds stress divergence. This energy is transferred to the water waves by pressure forces.

Departure of the numerical solutions from *O*(*ε*) perturbation theory is *O*(*ε*^{2}). For example, normalized growth rates *kc _{i}*/

*ε*and phase speeds

*c*are shown in Figs. 5 and 6, respectively, for various density ratios

_{r}*ε*. Perturbation theory predicts that the growth rates of the unstable waves should scale with

*ε*and that their phase speeds should be

*g*/

*k*to

*O*(

*ε*

^{2}): this is verified in Figs. 5 and 6. The shape of the growth rate curve predicted by the Miles equation (16) reflects the logarithmic boundary layer profile: notably, the rapid decrease in growth rates for longer waves shown in Fig. 5 disagrees with observations of wave growth (cf. Fig. 1).

### c. Optimal excitation of surface gravity waves

Although the solutions above verify the existence of an exponentially growing eigenstate that asymptotically dominates perturbation growth for stationary laminar winds, it does not follow that the growth of surface gravity waves proceeds exclusively from growth of this unstable mode even in the case of steady winds. It may be that a dominant role in wave growth results from transient growth occurring before the exponential asymptote is attained.

The Rayleigh equation has a continuous spectrum of nonanalytic modes that are neutral, and while each does not produce growth in isolation, because they are not orthogonal, they can cooperate to produce transient growth. Moreover, there exist initial conditions comprised of these neutral modes together with the unstable mode that excite the unstable mode optimally. In fact, because of the nonorthogonality of the continuous spectrum and the surface mode, the unit energy initial condition that optimally excites the unstable surface mode is not the mode itself, but rather the associated adjoint mode (Farrell 1988; Hill 1995; Farrell and Ioannou 1996a).

Consider introducing at *t* = 0 a perturbation *ψ*. This perturbation can be expanded in the basis of the eigenmodes *ψ _{i}* as

The adjoint *ϕ _{i}* of eigenfunction

*ψ*can be used to construct a biorthogonal basis such that (

_{i}*ϕ*,

_{i}*ψ*) = 0 when

_{j}*i*≠

*j*. Then perturbation

*ψ*excites eigenmode

*ψ*with amplitude

_{k}which is maximized when the initial perturbation *ψ* is chosen to be parallel to the adjoint *ϕ _{k}* of the target mode,

*ψ*. Consequently, the unit norm initial perturbation that optimally excites the unstable mode

_{k}*ψ*is its adjoint

_{u}*ϕ*:

_{u}which will excite the unstable mode magnified by the factor

The adjoint eigenfunctions of the Rayleigh equation satisfy the equation

with the boundary conditions *ϕ* → 0 at *z* → ∞ and

at the air–water interface *z* = 0.

The adjoint eigenfunction *ϕ _{i}* with eigenvalue

*λ*associated with the eigenfunction

_{i}*ψ*of the Rayleigh equation with eigenvalue

_{i}*c*has eigenvalue

_{i}*λ*

_{i}=

*c**

_{i}(the asterisk denotes complex conjugation). This adjoint eigenfunction is the same as the adjoint eigenfunction of the Rayleigh equation with solid boundary conditions (Drazin and Reid 1981):

The inner product for this adjoint is that associated with the square norm of total energy (kinetic plus potential). The inner product of eigenfunction *ϕ _{k}* with eigenvalue

*c**

_{k}and eigenfunction

*ψ*with eigenvalue

_{l}*c*is

_{l}With this inner product the analytic eigenfunctions and their adjoints form a biorthogonal basis, and the adjoint is the unit energy perturbation that maximally excites the mode.

### d. Amplification of wave height resulting from the adjoint perturbations

If the adjoint of the most unstable mode is introduced with unit energy, the corresponding normalized unstable mode is excited with amplitude |*α _{u}*|. Because the unstable mode dominates growth as

*t*→ ∞, this perturbation adjoint to the unstable mode is the perturbation that leads to optimal excitation of the surface wave as

*t*→ ∞ rather than the unstable mode itself.

This optimal amplification factor is shown as a function of wavenumber *k* for the logarithmic profile and for density ratios *ε* = 0.1, 0.001 in Fig. 7. The time needed to attain this large amplitude is *O*[1/(*kc _{i}*)] as seen in the example integrations shown in Fig. 8. The time development of the adjoint can be examined in numerical calculations in which a lid is placed at

*z*= 1 m so as to improve numerical conditioning of the initial value problem. For that case, the most unstable eigenmode and its adjoint are shown in Fig. 9, and the time development of the surface elevation is shown in Fig. 10. The adjoint perturbation in the neighborhood of the critical layer leans against the shear, consistent with extracting energy from the mean flow. The half-width of the adjoint perturbation variation in the critical-level neighborhood is seen from (24) to be

*c*/

_{i}*U*′(

*z*), where

_{c}*U*′(

*z*) is the background shear at the critical level. The growth rate

_{c}*kc*is proportional to

_{i}*ε*, so the width of the adjoint is also proportional to

*ε*. The adjoint perturbation can be idealized as a localized wave packet in the neighborhood of the critical level with central vertical wavenumber

*m*=

*πU*′(

*z*)/

_{c}*c*. This perturbation will shear over (Orr 1907; Kelvin 1887; Farrell 1982; Boyd 1983; Tung 1983; Farrell and Ioannou 1993) and will assume phase vertical orientation at a time

_{i}*T*=

*m*/[

*kU*′(

*z*)]. The maximum energy that can be attained by this perturbation is proportional to (

_{c}*m*/

*k*)

^{2}, and part of this energy will be transferred to the unstable surface wave exciting it with amplitude proportional to

*m*/

*k*, which is proportional to

*U*′(

*z*)/(

_{c}*kc*). This argument accounts for the scaling and time of development of the amplification factor obtained by the adjoint initial condition in Fig. 7.

_{i}A second factor important in the development of the adjoint is that the perturbation is concentrated near the critical layer where it is advected with the phase speed of the surface wave so that a form of resonance is induced with the pressure work building up the wave, as shown in Fig. 11. If a concentrated initial condition leaning against the shear is centered above or below the critical layer, then simulations confirm that, although the pressure work on the wave reaches the same maximum at the time the perturbation assumes a vertical position, it does not remain collocated with the surface wave and consequently will be unable to excite the wave to high amplitude.

### e. Discussion

In this section we have examined the problem of the growth of water waves from the point of view of generalized stability theory. It was demonstrated that in addition to the slowly growing exponential solutions, there exist initial perturbations not of modal form that produce transient growth rates far in excess of the growth rates produced by the unstable modes. In particular, unstable modes are shown to be highly suboptimal structures for the purpose of initiating wave growth. Nonmodal transient development results from atmospheric structures not of modal form but with the property of producing large pressure perturbations properly configured to produce wave growth.

These transiently growing structures no doubt play a role in instigating wave growth and certainly dominate the initial stages of growth. Nevertheless, it may be thought that to produce a fully developed wave field, initial perturbations must experience an extended period of exponential growth and that this implies that the role of transient growth is confined to exciting the unstable mode. However, it can be shown that transient growth can be sustained in turbulent flows so that the resulting growth is exponential in time (Farrell and Ioannou 1996b, 1999): this is the mechanistic basis for the stochastic parametric instability theory for the growth of surface water waves described below.

## 4. Exponential perturbation growth for time-dependent winds

### a. Gustiness-induced wave growth modeled as a generalized Kelvin–Helmholtz stability problem

Consider in (1) the simple time-dependent shear flow:

This flow is a generalization of the time-dependent version of the Kelvin–Helmholtz instability previously treated by Kelly (1965) with the crucial inclusion of a time-dependent shear in the atmosphere. The streamfunction in the air,

is consistent with the wave height, which is in turn determined by a single equation with the form of the harmonic oscillator equation,

in which appears the imaginary time-dependent coefficient in the place of damping,

and a complex restoring force:

where

is the square frequency of the surface gravity wave in the absence of flow. With the change of variable

the phase-shifted elevation *ζ* satisfies the harmonic oscillator equation:

with complex time-dependent square frequency:

In the absence of wind fluctuations instability occurs for Ω^{2} < 0, but this instability requires surface wind speeds *U*_{0} greater than observed to be associated with exponential growth of the longer water waves (Phillips 1977). The inclusion of sufficiently high shear extends the domain of the traditional Kelvin–Helmholtz instability without shear to longer waves, but the required surface winds remain unrealistically high.

While it is well known that harmonic modulations of Ω^{2}(*t*) lead to unstable growth and that this growth can primarily be traced to the second subharmonic of Ω^{2}(*t*) (Kelly 1965), it is a remarkable fact that even randomly modulated Ω^{2}(*t*), such that every realization of the oscillator is stable, can induce a positive Lyapunov exponent for wave height. This universal destabilization process arises when the instantaneous operators governing the dynamics produce transient growth and do not commute (Van Kampen 1976; Zel’Dovich et al. 1984; Farrell and Ioannou 1996b, 1999). It is instructive to demonstrate this process by assuming a constant shear *α*_{0} and time modulations of the velocity of the air–water interface of the form:

where *ξ*(*t*) is a zero mean random variable with unit variance. Assume further an air–water density ratio *ε* ≪ 1 and that the shear is sufficiently large that *εα*_{0} ≫ *kU*_{0}; with these assumptions, in the small noise limit, and retaining only the linear terms in *σ*, (33) is well approximated by

where

is the square frequency of the surface wave with mean surface velocity *U*_{0} and shear *α*_{0}.

Under these assumptions an estimate of the Lyapunov exponent can be obtained in the small *σ* limit (Arnold et al. 1986):

where

is the power of the random process at twice the unperturbed frequency, *ω _{d}*.

^{4}

At first glance, this *O*(*ε*^{2}) mean growth rate is disappointingly slow given that *ε* ≈ 0.001 for the air–water interface. However, the logarithmic layer that characterizes the mean wind profile near the surface has shears at the surface on the order of *u*^{2}_{*}/*ν* ≈ 10^{4} s^{−1}, where *u*_{*} is the friction velocity and *ν* is the molecular viscosity. Even as the surface becomes rough, shear remains large near the surface, being *u*_{*}/(*K z*_{0}) ≈ 10^{3} s^{−1} here, assuming a roughness length *z*_{0} of 1 mm and von Kármán constant *K* = 0.42 (Phillips 1977). It follows that the logarithmic boundary layer has shears that are at least *O*(1/*ε*) s^{−1} in the neighborhood of the air–water interface: because *U*_{0} is *O*(1) m s^{−1}, the growth rate is substantial.

The Lyapunov exponent *λ* is associated with the conventionally tabulated nondimensional growth rate per radian, *β*, through the relation:

For sufficiently long waves, *εα*_{0} ≫ *kU*_{0} and because the square of the surface wave frequency *ω*^{2}_{d} is nearly proportional to wavenumber *k*, the dependence of nondimensional growth rate on wavenumber is, to good approximation,

Correspondence of this estimate with observed growth rates for gusty boundary layer wind conditions is shown in Fig. 1 (dashed line) for friction velocity *u*_{*} = 0.2 m s^{−1}, *U*_{0} = 5*u*_{*}, normalized rms surface wind fluctuation *σ* = 0.3, mean normalized surface shears *εα* = 2 s^{−1}, and spectral density *ξ̂*_{1}(2*ω*_{d}) = 0.025 m^{2} s^{−1}.

This is a simple model of wave growth produced by a fluctuating boundary layer wind. This example demonstrates that variations in velocity and shear can destabilize a surface wave and that the expected growth rate is proportional to the product of the density ratio *ε* squared, the mean shear *α*_{0} squared, and the spectral energy content of the wind fluctuations at twice the wave frequency. In the next section we consider a more realistic time-varying wind and show that enhanced growth rates compared to growth rates obtained for the mean laminar wind also result in this more realistic boundary layer wind.

### b. Destabilization of surface water waves by atmospheric turbulence

The atmospheric boundary layer is typically turbulent and the associated characteristic wind field fluctuation is a phenomenon referred to as gustiness. From the point of view of time-dependent operator stability, this turbulent gustiness represents a structured time dependence of the operator: the structure being referred to confines fluctuations to terms in the operator representing the wind and its derivatives.

For time-dependent background winds we discretize (1) with the associated boundary conditions and write it symbolically as

where **Ψ** = [|η̂, *ψ̂*(*z*_{1}), *ψ̂*(*z*_{2}), · · · , *ψ̂*(*z*_{n})]^{T} is the state vector consisting of the surface elevation and the values of the streamfunction at the discretization points. Let *δ* be a small enough time interval that the operator can be approximated as constant on this interval. Over this interval the state **Ψ**(*t*) advances to exp[𝗔(*t*)*δ*]Ψ(*t*) and, consequently, the state at time *Nδ* is

The Lyapunov exponent that specifies the mean exponential growth rate of the wave height is

As discussed in the previous section, this exponent can be positive even if at each instant the operators 𝗔(*kδ*) have no unstable eigenvalues. Key to this behavior is the noncommutation of the operators at different time instants and the nonnormality of the operators 𝗔(*kδ*) (cf. Farrell and Ioannou 1999). This instability resulting from time dependence of the dynamical operator involves the full spectrum of the operator, and the resulting growth rate is therefore mechanistically distinct from that based solely on the dependence of the modal instability growth rate on the wind fluctuations (Komen et al. 1994; Miles 1997; Janssen 2004), although it reduces to these results for gusts with sufficiently long correlation times.

We now examine wave growth arising from gustiness in a boundary layer model with a logarithmic velocity profile:

with *K* = 0.42 the von Kármán constant. The gustiness can be parameterized using (45) by assuming a distribution for the wind at 10 m, *U*_{10}(*t*), and obtaining the consistent friction velocity *u*_{*}(*t*) and time-varying roughness length

expressed in terms of the Charnock constant (*κ* = 0.0144), the friction velocity *u**, and the acceleration of gravity *g*. We believe this simple model of gustiness to be adequate for the purposes of illustrating the mechanism of destabilization of waves by time-dependent winds. A more accurate viscous boundary layer model will be the subject of future study.

For *U*_{10}(*t*) we follow the wind model of Cavaleri and Burgers (1992) and Janssen (2004). There is good evidence that *U*_{10} is an approximately Gaussian random variable (Munn 1966; Smith et al. 1990; Cavaleri and Cardone 1994; Janssen 2004) with *U*_{10} (bar denotes time mean) normalized root-mean-square velocity variations, *σ*, that can exceed 0.3 for highly turbulent conditions (Cavaleri et al. 1981). Abdalla and Cavaleri (2002) plot the correlation between *σ* and the difference between water and air temperature and find that for a temperature difference between atmosphere and ocean of 10°C, the fluctuation level is *σ* ≈ 0.3: analysis of North Atlantic wind records indicate *σ* values in excess of 0.3. They report similar recorded values from the oceanographic tower of the Instituto Studio Dinamica Grandi Massi. The calculations that we report are made with *U*_{10} normalized rms velocity fluctuations of *σ* and assume

in which *ξ*(*t*) is a red noise process with zero mean, unit variance, and correlation time *τ*. Most calculations, presented here, use *τ* = 5 s. Calculations with different correlation times show that the results do not depend sensitively on *τ*, provided *τ* is not too long compared to the wave period. This is consistent with the analysis in the previous section in which the growth rate was shown to depend on the spectral energy content of the stochastically varying background flow at double the frequency of the surface. The stochastic forcing that we use is generated by the red noise process:

where *dW* is a Weiner process. This red noise process produces power at frequency *ω* with the following dependence on the correlation time *τ* :

In general, wind fluctuations that are sufficiently broadband will have spectral power at the required frequencies. The dependence of growth rate on correlation time for 6-m waves with *σ* = 0.25 and *ε* = 0.1 is verified in Fig. 12 to follow the theoretically predicted dependence given by (49).^{5}

As an example, consider wind fluctuations *U*_{10}(*t*) from (47) with *U*_{10} = 10 m s^{−1} corresponding to mean friction velocity *u*_{*} ≈ 0.4 and *U*_{10}-normalized rms variance *σ* = 0.15, which is the mean level of fluctuation observed when the temperature difference is just 6°C. As shown in Fig. 13, with this gustiness model the height of a 2.5-m wave grows exponentially with nondimensional growth rate *β* = 2.5, which is approximately three times greater than the corresponding rate for the constant logarithmic profile with *U*_{10} = *U*_{10}. Growth rates obtained for other wavenumbers and for various *u*_{*} are shown in Fig. 1 where for comparison growth rates predicted from the Miles laminar instability theory and growth rates from observations are also shown. The growth rate increases rapidly with the level of gustiness starting from the laminar instability value at *σ* = 0, as shown in Fig. 14. The calculation was done for a 6-m wave and *c*/*u*_{*} = 12.25 (*u*_{*} is the friction velocity of the time mean wind) and for *ε* = 0.1. After scaling, the growth rate for *ε* = 0.001 is approximately two orders of magnitude smaller. Note the rapid increase of the growth rate with the gustiness parameter *σ*. For this case the growth rate increases from the laminar growth rate with increasing gustiness approximated as Δ*β* = 0.1*σ*^{2} + 15*σ*^{4} + · · · , implying that for gustiness levels *σ* ≥ 0.35 the wave height growth rate is *O*(*σ*^{4}) larger than the laminar growth rate. This higher power dependence, obtained when the correlation time of the gusts is on the order of the wave period, contrasts with the linear dependence on *σ* that results from the modal quasi-static approximation of the effect of gustiness on wave height growth, valid when the correlation time of the gusts is considered to be long compared to the wave period (Janssen 2004).

### c. The lognormal distribution of wave heights and the occurrence of overshoots and rogue waves

Note that the wave height shown as a function of time in Fig. 13 fluctuates around its mean exponential growth rate, sometimes growing at a greater rate while at other times at a lesser rate. Realizations of wave height over 1500-s integrations for 3.14-m waves with density ratio *ε* = 0.001 are shown in Fig. 15. Realizations of wave height over 100-s integrations of a 6-m wave for density ratio *ε* = 0.1 sampled at 1-s intervals are shown in Fig. 16. These figures demonstrate the variation among realizations as well as the occurrence of large excursions in the amplitude of individual realizations of wave height. This provides an explanation for an important observed wave growth phenomenon referred to as overshooting by Liu and Ross (1980). It also provides an explanation for the occurrence of rogue waves with greater probability than would be predicted by a Gaussian fit to the ensemble of wave heights.

The emergence of overshoots, and subsequently of rogue waves, is a universal characteristic of multiplicative growth processes because such processes produce deviations from the mean that grow at a rate faster than the mean (Molčanov 1978; Arnold 1984; Kadanoff 2000; Farrell and Ioannou 2002). The rms fluctuation of the wave height, 〈*η*^{2}(*t*)〉, increases exponentially at rate *λ*_{2}, which is greater than the asymptotic growth rate *λ*_{1} of the magnitude of the wave height, 〈|*η*(*t*)|〉, which is in turn greater than the Lyapunov exponent (cf. Farrell and Ioannou 2002). In Fig. 16 the exponential growth of the wave height according to the Lyapunov growth rate, *λ* = 0.0982 s^{−1} for this example, is shown in a log plot of wave height. The rms wave-height growth, shown in Fig. 16 by the curve with squares, is dominated by extreme wave heights and grows exponentially at the much higher estimated growth rate.^{6} In Figure 16 the theoretical estimated exponential growth of the rms at rate *λ*_{2} = 0.2128 s^{−1} is also shown. This estimate is based on the accurate approximation that ensemble member growth rates over a sufficient interval have a normal distribution. In a plot such as Fig. 16 the points representing the wave heights of any finite number of simulations converge, almost surely, to growth at the Lyapunov rate as time progresses [this is the multiplicative ergodic theorem of Oseledec (Arnold 1998)]. For example, the wave height evolution of the realization (shown with crosses in Fig. 16) has outstanding growth initially but eventually regresses to growth at the mean Lyapunov rate. As time increases, the width of the distribution of the estimated Lyapunov growth rates narrows so that, in the limit *t* → ∞, its probability density function (pdf) becomes an infinitely narrow distribution about the Lyapunov exponent. This allows accurate calculation of the Lyapunov exponent from a finite number of ensemble simulations but makes it difficult to estimate with accuracy from a limited number of ensemble simulations the growth rate of higher moments of this distribution that increasingly depends, as the moment order increases, on extreme events. This is the reason why the rms growth of the ensemble (curve with squares) cannot capture, as time progresses, the theoretical rms growth from the accurate approximation of Gaussian distribution for the growth rate (dashed curve) in Fig. 16.

This model predicts that under gusty conditions local groups of waves that have experienced different realizations of gusts over time along their group path will have substantially different amplitudes tending toward a lognormal distribution with individual packets overshooting and then regressing to the mean, only to be overtaken by still larger overshoots in the limit producing the lognormal distribution. Although neighboring wave packets mix in the physical system to produce a mean observed wave field with less variation than that found in these idealized calculations, the implication of this calculation is that the probability of finding a rogue wave increases rapidly as the gustiness level increases.

## 5. Discussion

Excitation of surface water waves by wind is a familiar phenomenon of great theoretical and practical importance that continues to present challenges to our understanding (Liu et al. 2002; Cavaleri 2006; Janssen 2007). The atmospheric boundary layer is clearly turbulent, but incoherent pressure forcing from atmospheric pressure fluctuations produces variance growth linear in time (Phillips 1957) while growth exponential in time appears to be required by observations. The required exponential growth rate implies coherent and properly phased pressure fluctuations proportional to wave amplitude. The laminar critical layer instability mechanism of Miles (1957) produces exponential growth rates by maintaining this phase relationship but fails to correspond with observations of growth rates in general (Abdalla and Cavaleri 2002).

It has been found that the scatter in observed growth rates can be reduced by taking atmospheric stability into account (Kahma and Calkoen 1992). Data collected during strong offshore flow of cold air over warm water in winter (Barnett and Wilkerson 1967; Dobson 1971; Stewart and Teague 1980; Kahma 1981; Kahma and Calkoen 1992), which is expected to produce highly turbulent winds, exhibits enhanced wave growth rates compared with the growth rates found in data obtained in relatively stable tropical marine boundary layer flows (Snyder et al. 1981), which sustain much less turbulence. We argue here that this increased growth rate in unstable air arises from the enhancement of turbulent fluctuations in statically unstable boundary layers. In contrast to predictions of laminar stability theory in which growth rates rapidly decrease with wavelength, in this theory long waves can have high growth rates that increase as the gustiness level increases.

In this work we have demonstrated an essentially turbulent exponential instability, the fundamental explanation of which can be traced to the universal instability of time-dependent flows arising from the necessary nonnormality of time-dependent operators that are not at all times commuting. Returning to the classical mechanisms of wave growth: stochastic forcing by the turbulent surface pressure field associated with the boundary layer turbulent wind, and laminar modal instability; the stochastic parametric mechanism is essentially a result of the turbulence of the boundary layer, like the first mechanism. However, while the stochastic forcing enters additively in the Phillips theory, it enters multiplicatively in this theory and, thus, crucially produces exponential growth reducing to and extending the Miles theory as the turbulence level increases. This stochastic parametric instability dominates the laminar critical layer instability for sufficiently turbulent boundary layers and provides a mechanism for explaining the observed enhancement of surface water wave growth rate in gusty conditions as well as the occurrence of overshoots and rogue waves.

## Acknowledgments

This work was supported by NSF ATM-0123389. Petros Ioannou has been partially supported by an ELKE grant by the University of Athens and by the European Social Fund and National Resources (EPEAK II) PYTHAGORAS.

## REFERENCES

**.**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**.**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

## Footnotes

*Corresponding author address:* Brian Farrell, Department of Earth and Planetary Sciences, Geological Museum 452, Harvard University, 24 Oxford St., Cambridge, MA 02138. Email: farrell@seas.harvard.edu

^{1}

The data of Shemdin and Hsu (1967), Larson and Wright (1975), and Wu et al. (1979), and the Bight of Abaco data of Snyder et al. (1981), were taken from Plant (1982). The data of Kahma (1981), Barnett and Wilkerson (1967), Snyder and Cox (1966), Schule et al. (1971), and DeLeonibus and Simpson (1972) were taken from Kahma (1981). Here *u*_{*}/*U*_{10} = 0.0361 was used for the conversion of *c*/*U*_{10} to *c*/*u*_{*}. For the excellent Bragg scatterometer data of Stewart and Teague (1980), this conversion was made under the assumption that the wind profile was logarithmic.

^{2}

The assumption that the water motion arises only from the inviscid wave constitutes a limitation of our study. The velocities that arise from the coupling of the two layers by molecular viscosity is known to give rise to a new set of unstable modes (Hooper and Boyd 1983, 1987); in addition, inclusion of the extension of the turbulent boundary layer into the water (Bye 1995) may affect the stability properties of the mean flow.

^{3}

Miles (1957) argues that the linearized boundary condition should be formulated on the curvilinear surface of the wave and, consequently, that the shear term in the boundary condition (9) should be neglected, which has no effect on the growth rate of the waves at first order in *ε*. It does, however, result in predicted phase speeds in atmospheric boundary layer flows, for which *U* ′(0) is *O*(1/*ε*), that depart *O*(1) from the phase speed of the surface gravity waves, in disagreement with observations. The reason is that in (12) the real part of *Dψ*_{0}/*ψ*_{0} nearly cancels *U* ′(0)/*c*_{0}, giving a correction to the phase speed that is *O*(*ε*) even for large shears. If the term *U* ′(0)/*c*_{0} were not present in the boundary condition, then the phase speed correction would be *O*(1) for boundary layer flows for which *U* ′(0) is *O*(1/*ε*). That incorrect phase speeds are predicted if the mean shear term is not included in (12) can be verified by numerical solution of the eigenvalue problem.

^{4}

This estimate is a good asymptotic approximation as long as *ε**α*_{0}*kU*_{0}*σ*/*ω*^{2}_{d} ≪ 1.

^{5}

The stochastic fluctuations of the operator 𝗔(*t*) are dominated by the fluctuations of the curvature of the wind near the surface, which from (45) and (46) varies as *O*(1/*u*^{3}_{*}), resulting for the chosen level of *U*_{10} fluctuation in correlation time for the curvature-dominated wave forcing, approximately ¼ of the correlation time of *U*_{10}.

^{6}

*t*, approach a normal distribution with mean 〈

*λ*(

*t*)〉 that converges to the Lyapunov exponent and the rms variance,

*σ*(

*t*), which decreases as 1/

*t*. Knowing the pdf of

*λ*(

*t*) allows calculation of the

*p*th moment exponent: from which the growth rate of the rms wave height can be calculated as

*λ*

_{2}(

*t*) = 〈

*λ*(

*t*)〉 +

*tσ*

^{2}(

*t*), which asymptotically approaches a constant (cf. Farrell and Ioannou 2002). That the Lyapunov growth rates become normally distributed is verified in Fig. 16 in which the

*λ*(

*t*)

*t*for the ensemble is shown along with the curve of one standard deviation in growth that is seen itself to grow as

*t*.