Interaction of Tollmien-Schlichting Waves in the Air with the Sea Surface

TL;DR

This study analyzes the linear and non-linear stability of air-water flows to understand soliton generation by wind, deriving a non-linear Schrödinger equation and comparing different airflow profiles with sea observations.

Contribution

It introduces a comprehensive stability analysis of air-water interface flows, deriving a non-linear Schrödinger equation for soliton amplitude considering various airflow profiles.

Findings

Snake solitons observed for certain airflow profiles

Violent surface motion in turbulent boundary layer cases

Airflow nonlinearity dominates water interactions in wind-wave dynamics

Abstract

Linear stability of fully developed flows of air over water is carried out in order to study non-linear effects in the generation of solitons by wind. A linear stability analysis of the basic flow is made and the conditions at which solitons first begin to grow is determined. Then, following [10], the non-linear stability of the flow is examined and the quintic non-linear Schr\"{o}dinger equation is derived for the amplitude of disturbances. The coefficients of the non-linear Schr\"odinger equation are calculated from the eigenvalue problem which determines the stability of air-water interface. An asymptotic and a numerical stability analysis is carried out to determine the neutrally stable flow conditions for air-sea interface. Four different profiles are considered for the airflow blowing over the surface of the sea, namely, plane Couette flow (pCf), plane Poiseuille flow (pPf), laminar and turbulent boundary layer (L,TBL) profiles. For each of the above cases the shear flow counterpart in the water is assumed to be a pPf. It is shown that the above amplitude equation produces `snake' solitons [9] for pCf, pPf and LBL profiles, with striking similarities. On the other hand, for TBL we observe a very violent surface motion. For cases of pCf and LBL remarkable similarity is observed with observations made at sea. We conclude that the effect of nonlinearity in the airflow over the sea surface is much larger than nonlinear interactions in the water, and hence it is not possible to decouple the motion in the air and the water for finite amplitude wind-wave interactions, particularly in the case of wind-generated solitons in shallow waters.

Figures

Tables

Table 1. Table 1: Parameters for various models

Profile	$d_1$	$a_2$	$\kappa$	$\varpi$
pCf	$(0.895\times {10}^{-5},0.217\times {10}^{-2})$	$(0.895\times {10}^{-5},0.217\times {10}^{-2})$	$(-29.7,237)$	$(-5.94,47.4)$
pPf	$(0.404\times {10}^{-5},0.742\times {10}^{-3})$	$(0.557\times {10}^{-3},-0.17\times {10}^{-1})$	$(-9.95,912.3)$	$(-2.49,228)$
LBL	$(0.554\times {10}^{-5},0.275\times {10}^{-2})$	$(0.165\times {10}^{-2},-0.159)$	$(-0.221,20.61)$	$(-0.0737,6.87)$
TBL	$(0.222\times {10}^{-2},0.376\times {10}^{-3})$	$(0.042\times {10}^{-3},0.11\times {10}^{-1})$	$(-19.2,1142)$	$(-9.6,571)$

Equations

$$\frac{\partial A}{\partial \tau }-a_2\frac{{\partial }^{2}A}{\partial {\xi }^2}=\frac{d_1}{d_{1r}}A+\kappa A|A{|}^2+\varpi A|A{|}^4,$$

$$\eta =a\cos kx+\frac{1}{2}{a}^{2}k\cos 2kx,$$

$$\frac{\partial A}{\partial \tau }-a_2\frac{{\partial }^{2}A}{\partial {\xi }^2}=\frac{d_1}{d_{1r}}A+\kappa A|A{|}^2+\varpi A|A{|}^4,$$

$${\eta }_w=a\cos \left\{k_{c}x-k_c{c}_{r}t\left[\frac{1}{2}\alpha (a{k}_c{)}^2+\beta (a{k}_c{)}^2\right]\right\}$$

$$\displaystyle\left.\begin{array}[]{lll}U_{b}\rightarrow U_{\infty},&&\mbox{as}\hskip 7.11317pty\rightarrow D_{u}\\ l^{2}\left(\displaystyle\frac{dU_{b}}{dy}\right)^{2}=U_{*}^{2};\,\,l=\kappa(y+y_{0}),&&\mbox{in upper fluid}\\ u_{b}=U_{b},\hskip 7.11317pt\mu_{w}\displaystyle\frac{du_{b}}{dy}=\mu_{a}\displaystyle\frac{dU_{b}}{dy},&&\mbox{on}\hskip 7.11317pty=0\\ \displaystyle\frac{G}{\rho_{w}}+\nu_{w}\displaystyle\frac{{\rm d}^{2}u_{b}}{{\rm d}y^{2}}=0,&&\mbox{in lower fluid}\\ u_{b}=0,&&\mbox{on}\hskip 7.11317pty=-D_{\ell}\\ \end{array}\right\}$$

$$U_b(y)=U_1\ln \left(\frac{y+y_0}{y_0}\right)+U_0$$

$$u_b(y)=\frac{G}{2{\mu }_w}(D_{\ell }^2-y^2)+\frac{{\mu }_a}{{\mu }_w}\frac{U_1}{y_0}(D_{\ell }+y)$$

$$U_0=\frac{D_{\ell }}{{\mu }_w}\left(\frac{U_1{\mu }_a}{y_0}+\frac{G{D}_{\ell }}{2}\right),U_1=\frac{U_{\ast }}{\kappa }$$

$$\displaystyle\left.\begin{array}[]{l}\chi=x-\displaystyle\frac{a\cosh k(y\mp d)}{\sinh kd}\sin k(x-ct)\\ \\ \eta=y-\displaystyle\frac{a\sinh k(y\mp d)}{\sinh kd}\cos k(x-ct)\\ \\ \vartheta=z\end{array}\right\}$$

$$J=1-\frac{2ak}{\sinh kd}\cosh k(\eta \mp d)\cosh k(\chi -ct)+O(ak{)}^2$$

$$J=1+O(ak)$$

$$\displaystyle\left.\begin{array}[]{l}\displaystyle\frac{\partial\zeta}{\partial t}+U\displaystyle\frac{\partial\zeta}{\partial\chi}+V\displaystyle\frac{\partial\zeta}{\partial\eta}=W\\ \\ \displaystyle\frac{\partial\zeta}{\partial t}+u\displaystyle\frac{\partial\zeta}{\partial\chi}+v\displaystyle\frac{\partial\zeta}{\partial\eta}=w\\ \end{array}\right\}$$

$$\displaystyle\left.\begin{array}[]{l}\displaystyle\frac{\partial U}{\partial\eta}+\displaystyle\frac{\partial V}{\partial\chi}=\mu^{-1}\left(\displaystyle\frac{\partial u}{\partial\eta}+\displaystyle\frac{\partial v}{\partial\chi}\right)\\ \\ \displaystyle\frac{\partial W}{\partial\eta}+\displaystyle\frac{\partial V}{\partial\vartheta}=\mu^{-1}\left(\displaystyle\frac{\partial w}{\partial\eta}+\displaystyle\frac{\partial v}{\partial\vartheta}\right)\\ \\ \displaystyle\frac{\partial U}{\partial\vartheta}+\displaystyle\frac{\partial W}{\partial\chi}=\mu^{-1}\left(\displaystyle\frac{\partial u}{\partial\vartheta}+\displaystyle\frac{\partial w}{\partial\chi}\right)\\ \end{array}\right\}$$

$$-\rho P+\frac{2\rho }{R}\frac{\partial V}{\partial \eta }+\rho \mathcal{F}\zeta =-\rho p+\frac{2}{R\nu }\frac{\partial v}{\partial \eta }-\mathcal{T}{\nabla }_{\parallel }^2\zeta +\mathcal{F}\zeta$$

$$\displaystyle\mbox{\boldmath$U$}=[U_{b}(\eta),0,0]+\delta[U(\eta),V(\eta),W(\eta)]e^{i(k\chi+\ell\vartheta-kct)}+{\rm c.c.}$$

$$\displaystyle\mbox{\boldmath$u$}=[u_{b}(\eta),0,0]+\delta[u(\eta),v(\eta),w(\eta)]e^{i(k\chi+\ell\vartheta-kct)}+{\rm c.c.}$$

$$P=\frac{\chi U_b^{\prime \prime }}{\rho \mu R}+\delta P(\eta )e^{i(k\chi +\ell \vartheta -kct)}+c.c.$$

$$p=\frac{\chi u_b^{\prime \prime }}{\rho \mu R}+\delta p(\eta )e^{i(k\chi +\ell \vartheta -kct)}+c.c.$$

$$\zeta =\delta a{e}^{i(k\chi +\ell \vartheta -kct)}+c.c.,|\delta |\ll 1$$

$$\displaystyle\left.\begin{array}[]{lcl}U=V=W=0&&\mbox{on $\eta=1$}\\ \\ ik(U_{b}-c)U+VU_{b}^{\prime}=-ikP+R^{-1}[U^{\prime\prime}-(k^{2}+\ell^{2})U]&&\\ \\ ik(U_{b}-c)V=-P^{\prime}+R^{-1}[V^{\prime\prime}-(k^{2}+\ell^{2})V]&&\\ \\ ik(U_{b}-c)W=-i\ell P+R^{-1}[W^{\prime\prime}-(k^{2}+\ell^{2})W]&&\\ \\ ikU+V^{\prime}+i\ell W=0&&\\ \\ u+au^{\prime}_{b}=U+aU_{b}^{\prime},\,\,\,v=V,\,\,w=W&&\mbox{on $\eta=0$}\\ \\ ik(u_{b}-c)a=v,\,\,\,ik(U_{b}-c)a=V&&\mbox{on $\eta=0$}\\ \\ u^{\prime}+ikv=\mu(U^{\prime}+ikV),\,\,\,w^{\prime}+i\ell v=\mu(W^{\prime}+i\ell V)&&\mbox{on $\eta=0$}\\ \\ \rho p-{\mathscr{F}}a-2v^{\prime}/\nu R-{\mathscr{T}}a(k^{2}+\ell^{2})=\rho\left(P-{\mathscr{F}}a-2V^{\prime}/R\right)&&\mbox{on $\eta=0$}\\ \\ u=v=w=0&&\mbox{on $\eta=d^{-1}$}\\ \\ ik(u_{b}-c)u+vu_{b}^{\prime}=-ik\rho p+(\nu R)^{-1}[u^{\prime\prime}-(k^{2}+\ell^{2})u]&&\\ \\ ik(u_{b}-c)v=-\rho P^{\prime}+(\nu R)^{-1}[v^{\prime\prime}-(k^{2}+\ell^{2})v]&&\\ \\ ik(u_{b}-c)w=-i\ell\rho P+(\nu R)^{-1}[w^{\prime\prime}-(k^{2}+\ell^{2})w]&&\\ \\ iku+v^{\prime}+i\ell w=0&&\\ \end{array}\right\}$$

$$\displaystyle\left.\begin{array}[]{lll}U={\Phi}^{\prime},&&u={\phi}^{\prime}\\ &&\\ V=-ik\Phi,&&v=-ik\phi\end{array}\right\}$$

$$\displaystyle\left.\begin{array}[]{lll}\Phi={\Phi}^{\prime}=0,&&\mbox{on}\hskip 7.11317pt\eta=1\\ \\ ik(U_{b}-c)({\Phi}^{\prime\prime}-k^{2}\Phi)-ik{U}^{\prime\prime}_{b}\Phi=&&\\ &&\\ \hskip 85.35826ptR^{-1}(\Phi^{iv}-2k^{2}{\Phi}^{\prime\prime}+k^{4}\Phi)&&\mbox{in upper fluid}\\ \\ \left.\begin{array}[]{l}\phi=\Phi\\ \\ {\phi}^{\prime}-\displaystyle\frac{{u}^{\prime}_{b}\phi}{\sigma}={\Phi}^{\prime}-\displaystyle\frac{{U}^{\prime}_{b}\Phi}{\sigma}\\ \\ {\phi}^{\prime\prime}+k^{2}\phi=\mu({\Phi}^{\prime\prime}+k^{2}\Phi)\\ \\ \displaystyle\frac{1}{\nu R}({\phi}^{\prime\prime\prime}-3k^{2}{\phi}^{\prime})+ik(\phi u_{b}-{\phi}^{\prime}\sigma)+i\displaystyle\frac{k\phi}{\sigma}\displaystyle\frac{({\sf F}+k^{2}{\sf T})}{R^{2}}=\\ \\ \displaystyle\frac{\rho}{R}({\Phi}^{\prime\prime\prime}-3k^{2}{\Phi}^{\prime})+\rho ik(\Phi U_{b}-{\Phi}^{\prime}\sigma)+\rho i\displaystyle\frac{k\Phi}{\sigma}\displaystyle\frac{{\sf F}}{R^{2}}\\ \end{array}\right]&&\mbox{on}\hskip 7.11317pt\eta=0\\ \\ ik(u_{b}-c)({\phi}^{\prime\prime}-k^{2}\phi)-ik{u}^{\prime\prime}_{b}\phi=&&\\ &&\\ \hskip 85.35826pt(\nu R)^{-1}(\phi^{iv}-2k^{2}{\phi}^{\prime\prime}+k^{4}\phi)&&\mbox{in lower fluid}\\ \\ \phi={\phi}^{\prime}=0,&&\mbox{on}\hskip 7.11317pt\eta=-d^{-1}\\ \end{array}\right\}$$

$$\displaystyle\left.\begin{array}[]{l}c=c_{0}+kc_{1}+...\\ \\ \phi=\phi_{0}+k\phi_{1}+...\\ \\ \Phi=\Phi_{0}+k\Phi_{1}+...\\ \end{array}\right\}$$

$$\displaystyle\left.\begin{array}[]{lll}\Phi_{0}={\Phi}^{\prime}_{0}=0,&&\mbox{on}\hskip 7.11317pt\eta=1\\ \\ 0=R^{-1}\Phi_{0}^{iv}&&\mbox{in upper fluid}\\ \\ \\ \left.\begin{array}[]{l}\phi_{0}=\Phi_{0}\\ \\ {\phi}^{\prime}_{0}-\displaystyle\frac{{u}^{\prime}_{b}\phi_{0}}{\sigma}={\Phi}^{\prime}_{0}-\displaystyle\frac{{U}^{\prime}_{b}\Phi_{0}}{\sigma},\hskip 7.11317pt\sigma=u_{b}c_{0}\\ \\ {\phi}^{\prime\prime}_{0}=\mu{\Phi}^{\prime\prime}_{0}\\ \\ (\nu R)^{-1}{\phi}^{\prime\prime\prime}_{0}=\rho R^{-1}{\Phi}^{\prime\prime\prime}_{0}\\ \end{array}\right]&&\mbox{on}\hskip 7.11317pt\eta=0\\ \\ \\ 0=(\nu R)^{-1}\Phi_{0}^{iv}&&\mbox{in lower fluid}\\ \\ \phi_{0}={\phi}^{\prime}_{0}=0,&&\mbox{on}\hskip 7.11317pt\eta=-d^{-1}\\ \end{array}\right\}$$

$${\Phi }_0=(\eta -1{)}^2\left(\eta \frac{d}{2\mu }\frac{\mu -d^2}{1+d}+1\right)$$

$${\varphi }_0=(\eta +d^{-1}{)}^2\left(\eta \frac{d}{2}\frac{\mu -d^2}{1+d}+d^2\right).$$

$$\displaystyle c_{0}=u_{b}(0)+\displaystyle\frac{2\mu d(1+d)({U}^{\prime}_{b}-{u}^{\prime}_{b})}{4\mu d(1+d)^{2}+(\mu-d^{2})^{2}}\hskip 14.22636pt\mbox{provided $\mu\neq 1$}$$

$${\Phi }_0=A[\cosh k(\eta -1)-\cosh \lambda (\eta -1)]+B[\sinh \lambda (\eta -1)-\lambda k^{-1}\sinh k(\eta -1)]$$

$${\varphi }_0=a[\cosh k(\eta +d^{-1})-\cosh \gamma (\eta +d^{-1})]+b[\sinh \gamma (\eta +d^{-1})-\gamma k^{-1}\sinh k(\eta +d^{-1})]$$

$$\displaystyle\left.\begin{array}[]{lll}\Phi^{(1)}={\Phi^{(1)}}^{\prime}=0,&&\mbox{on}\hskip 7.11317pt\eta=1\\ \\ i(U_{b}-c_{0}){\Phi}^{\prime\prime}_{0}-i{U}^{\prime\prime}_{b}\Phi_{0}=R^{-1}\Phi^{iv}_{1}&&\mbox{in upper fluid}\\ &&\\ \left.\begin{array}[]{l}\phi_{0}=\Phi_{0}\\ \\ {\phi}^{\prime}_{1}-\displaystyle\frac{{u}^{\prime}_{b}}{\sigma}\left(\phi_{1}+\displaystyle\frac{c_{1}}{\sigma}\phi_{0}\right)=\\ \\ \hskip 78.24507pt{\Phi}^{\prime}_{1}-\displaystyle\frac{{U}^{\prime}_{b}}{\sigma}\left(\Phi_{1}+\displaystyle\frac{c_{1}}{\sigma}\Phi_{0}\right)\\ \\ {\phi}^{\prime\prime}_{1}=\mu{\Phi}^{\prime\prime}_{1}\\ \\ (\nu R)^{-1}{\phi}^{\prime\prime\prime}_{1}-i\left({\phi}^{\prime}_{1}\sigma-\phi_{0}{u}^{\prime}_{b}\right)+\displaystyle\frac{i\phi_{0}}{\sigma}\displaystyle\frac{{\sf F}}{R^{2}}=\\ \\ \hskip 28.45274pt\rho R^{-1}{\Phi}^{\prime\prime\prime}_{1}-\rho i\left({\Phi}^{\prime}_{1}\sigma-\Phi_{0}{U}^{\prime}_{b}\right)+\\ \\ \hskip 28.45274pt\displaystyle\frac{\rho i\Phi_{0}{\sf F}}{\sigma R^{2}}\\ \\ \end{array}\right]&&\mbox{on}\hskip 7.11317pt\eta=0\\ \\ i(u_{b}-c_{0}){\phi}^{\prime\prime}_{0}-i{u}^{\prime\prime}_{b}\phi_{0}=(\nu R)^{-1}\phi^{iv}_{1}&&\mbox{in lower fluid}\\ \\ \phi_{1}={\phi}^{\prime}_{1}=0,&&\mbox{on}\hskip 7.11317pt\eta=-d^{-1}\\ \end{array}\right\}$$

$$\displaystyle\left.\begin{array}[]{lll}\Psi={\Psi}^{\prime}=0,&&\mbox{on}\hskip 7.11317pt\eta=1\\ \\ 0=R^{-1}\Psi^{iv}&&\mbox{in upper fluid}\\ \\ \left.\begin{array}[]{l}\psi=\Psi,\hskip 14.22636pt{\psi}^{\prime}={\Psi}^{\prime},\hskip 14.22636pt{\psi}^{\prime\prime}=\mu{\Psi}^{\prime\prime}\\ \\ {\psi}^{\prime\prime\prime}-\displaystyle\frac{{u}^{\prime}_{b}}{\sigma}{\psi}^{\prime}=\mu\left({\Psi}^{\prime\prime\prime}-\displaystyle\frac{{U}^{\prime}_{b}}{\sigma}{\Psi}^{\prime}\right)\\ \\ \end{array}\right]&&\mbox{on}\hskip 7.11317pt\eta=0\\ \\ 0=(\nu R)^{-1}\psi^{iv}&&\mbox{in lower fluid}\\ \\ \psi={\psi}^{\prime}=0,&&\mbox{on}\hskip 7.11317pt\eta=-d^{-1}\\ \end{array}\right\}$$

$$(\nu R{)}^{-1}{\int }_{-d^{-1}}^0\psi {\varphi }_0^{iv}d\eta +\rho R^{-1}{\int }_0^1\Psi {\Phi }_0^{iv}d\eta =(\nu R{)}^{-1}{\int }_{-d^{-1}}^0{\varphi }_0{\psi }^{iv}d\eta +\rho R^{-1}{\int }_0^1{\Phi }_0{\Psi }^{iv}d\eta =0$$

$$\Psi =(\eta -1{)}^2\left[\eta \frac{\mu -d^2+4d(1+d)}{2(\mu +d)}+1\right]$$

Full text

An asymptotic and a numerical stability analysis is carried out to determine the neutrally stable flow conditions for air-sea interface. Four different profiles are considered for the airflow blowing over the surface of the sea, namely, plane Couette flow (pCf), plane Poiseuille flow (pPf), laminar and turbulent boundary layer (L,TBL) profiles. For each of the above cases the shear flow counterpart in the water is assumed to be a pPf.

A nonlinear stability analysis results in the nonlinear Schrödinger equation

[TABLE]

where $\xi$ and $\tau$ are local variables, and the amplitude of the surface wave is proportional to $A$. The complex constants $a_{2}$, $d_{1}$, $\kappa$ and $\varpi$ are evaluated from the linear stability of the two-fluid interface. The profile for the initial condition considered here is that of the Stokes wave

[TABLE]

where $a$ is the amplitude and $k$ is the wavenumber of the surface Stokes wave.

It is shown that the above amplitude equation produces ‘snake’ solitons [9] for pCf, pPf and LBL profiles, with striking similarities. On the other hand, for TBL we observe a very violent surface motion. For cases of pCf and LBL remarkable similarity is observed with observations made at sea.

We conclude that the effect of nonlinearity in the airflow over the sea surface is much larger than nonlinear interactions in the water, and hence it is not possible to decouple the motion in the air and the water for finite amplitude wind-wave interactions, particularly in the case of wind-generated solitons in shallow waters.

1 Introduction

The wind-generated solitons which moves nondispersively on the thin wind-driven drift layer of wind-excited surface wave waves are important because they are products of the processes under which surface waves are generated by wind. When the local atmospheric conditions are such that significant bubble densities are created, principally by breaking waves, the resulting bubble layer will dominate and mask these solitons, in underwater acoustic scattering. A qualitative mixed linear-nonlinear model of ocean-wave generation is briefly described, a key element of which is the generation of solitons by wind. These solitons are created, and destroyed on moving gravity waves, caused by waves being generated by wind along with accompanying capillary components. The gravity-capillary waves are dispersive, in contrast the solitons are nondispersive, that is to say, they move at constant speed, on the thin moving wind drift layer generated by the upper part of the water surface.

In this paper, we will examine the stability of laminar and turbulent flows of air over water surface where viscous effects are taken into account in both fluids. The flow that will be studied here is the fully developed flow of air over water. In contrast to [2], we shall restrict our analysis to only the case when the motion is driven solely by a constant pressure gradient.

A linear stability analysis will be performed to determine the neutrally stable flow conditions. This approach is different to [4], or indeed other similar theories such as [1], and most experimental work, where the emphasis is on determining wave growth rates rather than the conditions at which waves will first start to grow. Following [10], a non-linear analysis will be performed which will yield the non-linear quintic Schrödinger equation

[TABLE]

governing the wave amplitude. In equation (1) $A$ is proportional to amplitude of the non-linear surface wave. The scaled length and time, $\xi$ and $\tau$, and the complex constants $a_{2},d_{1}$ and $\kappa$ and $\varpi$ will be defined in the proceeding analysis. The constants $a_{2}$ and $d_{1}$ are properties of the flow which is obtained from the linear stability theory, while the effect of non-linear interactions is determined by $\kappa$ and $\varpi$.

It is to be noted that although equation (1) is strictly applicable to velocity profiles of laminar flow, there are experimental and theoretical suggestions that equations similar to that of (1) are also applicable to wind-wave generation situations. For example, the experimental work [6] suggests that wind-waves at fixed fetch and under steady wind conditions have properties that are similar to non-linear Stokes wavetrains. Indeed, the wind-waves given by the equilibrium amplitude solutions of (1) (for which $\kappa_{r}<0$ and $\varpi=0$) are similar to Stokes wavetrains

[TABLE]

where $k_{c}$ is the critical wavenumber, $\alpha$ and $\beta$ are the dimensionless energy transfer parameters [13], and $c$ is the wave phase speed. It has been shown [11] that such wavetrains with non-linear effects lead to changes in the phase speed of the wave and provides the sideband instability mechanism, similar to that of Benjamin-Fier, which will eventually result to the break-up of the wave-train.

The wind-generated solitons which moves nondispersively on the thin wind-driven drift layer of wind-excited surface wave waves are important because they are products of the processes under which surface waves are generated by wind. When the local atmospheric conditions are such that significant bubble densities are created, principally by breaking waves, the resulting bubble layer will dominate and mask these solitons, in underwater acoustic scattering. A qualitative mixed linear-nonlinear model of ocean-wave generation is briefly described, a key element of which is the generation of solitons by wind. These solitons are created, and destroyed on moving gravity waves, caused by waves being generated by wind along with accompanying capillary components. The gravity-capillary waves are dispersive, in contrast the solitons are nondispersive, that is to say, they move at constant speed, on the thin moving wind drift layer generated by the upper part of the water surface.

In this paper, we will examine the stability of laminar and turbulent flows of air over water surface where viscous effects are taken into account in both fluids. The flow that will be studied here is the fully developed flow of air over water. In contrast to [2], we shall restrict our analysis to only the case when the motion is driven solely by a constant pressure gradient.

2 Formulation of the problem

We consider the steady, flow of two super-imposed fluid layers confined between a stationary horizontal parallel plate, placed at $y=-D_{\ell}$, and open to atmosphere above at a distance $y=D_{u}\gg 1$. We assume the flow of the upper fluid is turbulent, satisfying Prandtl’s mixing-length model, and that of the lower one is a fully developed laminar flow. The velocity vector \mbox{\boldmathU}=[U_{b}(y),0,0], for the upper fluid, and \mbox{\boldmathu}=[u_{b}(y),0,0], for the lower fluid,111Hereafter, the upper and the lower symbols (Greek or Roman) denote the upper (air) and the lower (water) fluids, respectively. Also the subscripts $a$ and $w$ denote air and water, respectively. both satisfy the Navier-Stokes equations. In the present case we shall confine ourselves to the case when the motion of the lower fluid is driven solely by a constant pressure gradient. Thus, the equations of motion and their corresponding boundary conditions reduce to

[TABLE]

where $-G=\displaystyle\frac{dp}{dx}.$ The solutions of (8) are given by

[TABLE]

where

[TABLE]

$U_{*}$ is the friction velocity of the upper fluid, $\kappa$ is the von Kármán’s constant, and $y_{0}$ is an effective roughness length.

Having obtained the velocity fields governing the basic flow, we shall now formulate the non-dimensional equations that govern linear instability. The motion is governed by the three-dimensional Navier-Stokes equations in each layer of fluid. The boundary conditions are the no-slip conditions on the plates and continuity of stress and velocity across the unknown interface $y=\zeta(x,z,t)$, which will be determined from the proceeding non-linear analysis.

We will introduce the following curvilinear coordinates [5]

[TABLE]

where in (17), the minus sign is for the upper media and the plus sign for the lower media. The Jacobian of the transformation is

[TABLE]

Note that, in the linear stability theory we shall be restricted to infinitesimal disturbances for which the Jacobian (18) reduce to

[TABLE]

Thus, within this order of approximation, the governing equations are symbolically the same as that of Cartesian coordinates. We remark that under the transformation (17), the lines $y=d$ and $y=-d$ will map to lines $\eta=1$ and $\eta=-d^{-1}$, respectively, and the unknown interface $y=\zeta$ will then be at $\eta=0$.

The boundary conditions are:

on $\eta=1$ no-slip conditions $\hskip 34.14322ptU=V=W=0$,

on $\eta=0$ (the unknown interface);

(i) the continuty of velocities \hskip 25.6073pt\mbox{\boldmathu}=\mbox{\boldmathU}

(ii) the kinematic conditions

[TABLE]

(iii) the continuity of tangential stresses

[TABLE]

(iv) the dynamic condition

[TABLE]

where $\nabla_{\parallel}^{2}\equiv\partial^{2}/\partial\chi^{2}+\partial^{2}/\partial\eta^{2}$, and $d=D_{u}/D_{\ell},\mu=\mu_{a}/\mu_{w},\rho=\rho_{a}/\rho_{w},\nu=\nu_{a}/\nu_{w}.$ Also, on $\eta=-d^{-1}$ there is no-slip conditions $u=v=w=0$.

In (29) $R=U_{\infty}d/\nu_{a}$ is the Reynolds number, ${\mathscr{F}}=gd/U^{2}_{\infty}$ is the Froude number and ${\mathscr{T}}=S/\rho_{w}U_{\infty}^{2}d$ is the Weber number, where $U_{\infty}$ is a reference velocity, $S$ is the dimensional surface tension and $g$ is the acceleration due to gravity.

We next perturb the basic flow by infinitesimal disturbances, i.e.

[TABLE]

where c.c. denotes a complex conjugate, then the linearized equations of motion may be cast as

[TABLE]

Substituting in equations of motion and their corresponding boundary conditions, making use of the Squire transformation [14] $m^{2}=k^{2}+\ell^{2},\hskip 7.11317pt\hat{U}=m^{-1}(kU+\ell W),\hskip 7.11317pt\hat{V}=V,\hskip 7.11317pt\hat{c}=c,\hat{P}=(m/k)P,\hskip 7.11317pt\hat{a}=(k/m)a$, etc. (dropping the hat notation) and introducing the stream functions $\Phi(\eta)$ and $\phi(\eta)$ such that

[TABLE]

we obtain, after linearization and some manipulations,

[TABLE]

where $\sigma=u_{b}(0)-c=U_{b}(0)-c$. Note that, we have redefined the Froude and Weber number to be ${\sf F}=gD_{u}^{3}/\nu_{a}^{2}$ and ${\sf T}=SD_{u}/\rho_{w}\nu_{a}^{2}$, respectively.

Systems similar to that of (84) have been derived [7,15–17]. However, as was pointed out by Blennerhassett [2], in the derivations [7,17] the relevance of ${\mathscr{F}}$ and ${\mathscr{T}}$ to the stability of two-dimensional waves was not noted. He further remarks that in earlier derivations [15,16], for the stability of superposed fluid layers, the terms $u_{b}^{\prime}\phi/\sigma$ and $U_{b}^{\prime}\Phi/\sigma$ was incorrectly omitted from the second interface condition in the system (84).

We remark that in general the solution of (84), and determining the location of the neutral stability curve in $(R,k)$–plane is a computational task. The full details of the numerical methods is given in section 6 below.

3 The long-wave limit

Following [2], we consider the stability of the basic flow in the limit of long wavelength disturbances and expand the eigenfunctions and the eigenvalue in regular perturbation series in $k$ [7]. Thus we may write

[TABLE]

as $k\rightarrow 0$ for a fix value of $R,{\sf F},{\sf T}$. Substituting (90) in (84) we obtain up to $O(1)$ the following homogeneous system

[TABLE]

In the case of pPf, the system (109) has a eigensolution

[TABLE]

and

[TABLE]

The eigenvalue $c_{0}$ is given by

[TABLE]

Similarly, in the limit of slow basic flow, when the Reynolds number for the basic flow is small, we may expand $c\sim R^{-1}c_{0}+c_{1}+...$, $\phi\sim\phi_{0}+R\phi_{1}+...$ and $\Phi\sim\Phi_{0}+R\Phi_{1}+...$ as $R\rightarrow 0$ for fixed ${\mathscr{F}}$ and ${\mathscr{T}}$. Thus we obtain a system similar to (109) whose solution is

[TABLE]

where $\lambda^{2}=k(k-ic_{0})$ and $\gamma^{2}=k(k-i\nu c_{0})$, and the eigenvalues $c_{0}$ and the constants $A,B,a$ and $b$ are determined by applying the interface conditions at $\eta=0$. This yields a homogeneous $4\times 4$ set of equations whose non-trivial solutions can only be obtained for certain values of $c_{0}$.

We next consider $O(k)$. To this order the system (84) reduce to

[TABLE]

Since $\phi^{(0)},\Phi^{(0)}$ and $c^{(0)}$ are real, we can easily see that $c^{(1)}$ is purely imaginary. Note that the homogeneous part (134) is the same as (109). Thus, in order that the inhomogeneous system to have a solution we must satisfy certain solvability criterion. We introduce the adjoint functions [6], $\Psi$ and $\psi$, and write down the system which is adjoint to (134)

[TABLE]

with the adjoint condition 222Stuart’s solvability criterion states that if $L$ is a differential operator and $N$ is its adjoint operator, then if $L\Phi=f$ and its adjoint is given by $N\Psi=0$ the solvability condition states

$\int_{a}^{b}\Psi L\Phi\,d\eta=\int_{a}^{b}\Phi\underbrace{N\Psi}_{=0}\,d\eta$ which implies $\int_{a}^{b}\Psi f\,d\eta=0.$

[TABLE]

The eigenvalue $c_{0}$ for the adjoint system (148) is the same as that of (109), while the adjoint functions for pPf are

[TABLE]

When the solvability condition is applied to (109) we obtain the following expression for $c_{1}$:

[TABLE]

To $O(k)$, the condition for neutral stability is that $c_{1}=0$. Hence, from (149) we obtain the following relationship between $k$ and $R$:

[TABLE]

Equation (150) is a quadratic in $\sigma$ which may be expressed as

[TABLE]

with its solutions given by

[TABLE]

Since $\sigma=u_{b}(0)-c$, we obtain the following compact relation between the wavenumber $k$ and the Reynolds number $R$, i.e.

[TABLE]

4 The short-wave limit

We shall now focus our attention to conditions under which short-wavelengths disturbances may grow along the interface. We shall pose the following restrictions on the problem

[TABLE]

To study this, we will adopt an asymptotic analysis for pPf. However, before we commence with our analysis, it is convenient to express the interfacial boundary conditions (84) in more compact mathematical form. We do this by defining

[TABLE]

Thus, the second to the fourth interfacial conditions (84) may be written as

[TABLE]

We shall also scale the coordinates $(\chi,\eta)$ with the wavenumber $k$ such that

[TABLE]

and assume the following asymptotic expansions for $c,\phi$ and $\Phi$:

[TABLE]

Neglecting terms of $O(d^{-1})$ by invoking the restriction (a) in (24), replacing the independent variable by $\hat{\varphi}=-\varphi$ and defining ${\mathscr{Y}}=\hat{\varphi}/d$  $(0\leq{\mathscr{Y}}\leq 1)$, the eigenfunction for the lower fluid will satisfy the following Orr-Sommerfeld equation

[TABLE]

in which we anticipate $c\sim 1-O(k^{1/3})$.

Since the effective Reynolds number, $R\rho d^{3}/\mu^{2}$, is large, we may neglect the viscous terms in (177) and thus reduce it to the Rayleigh equation

[TABLE]

Equation (178) has a regular solution

[TABLE]

and a solution with logarithmic branch point at ${\mathscr{Y}}=0$ and ${\mathscr{Y}}=1$, namely

[TABLE]

where the superscript $o$ in (179) and (180) denotes the outer solution.

Now to recover the boundary layer near the interface, and following the usual practice in the theory of hydrodynamic stability, we introduce $\varepsilon=(iR\rho kd^{3}/\mu^{2})^{-1/3}$ and let ${\sf Y}={\mathscr{Y}}/\varepsilon$. Thus, (177) becomes

[TABLE]

Expanding $\phi$ in powers of $\varepsilon$ we obtain, at zero order in $\varepsilon$

[TABLE]

and at $O(\varepsilon)$ we get

[TABLE]

The solution of (182) is

[TABLE]

In (184), ${\mathscr{A}}_{k}({\sf Y},2),\,\,k=1,2$, denote generalized Airy function

[TABLE]

where ${\mathscr{L}}_{k}$ is the appropriate path of integrals shown in figure 1.

For ${\sf Y}\gg 0$ the generalized Airy function ${\mathscr{A}}_{2}({\sf Y},2)$ is exponentially large. Thus we set the constant $m_{4}=0$ in (184). The constants $m_{1}$ and $m_{2}$ are found by matching (184) with the inviscid solutions (179) and (180). However, to do this we require the next-order solutions to the constant and linear parts in (184).

From (184) we see that

[TABLE]

where the superscript $I$ represent the inner solution and ${\mathscr{B}}_{3}({\sf Y},2,1)$ is a generalized Airy function given by

[TABLE]

Now matching $\phi_{2}^{I}$ with $\phi_{2}^{o}$, and $\phi_{1}^{I}$ with $\phi_{1}^{o}$, we obtain, respectively

[TABLE]

where $\psi(2)$ is the psi function.

To evaluate the above integral we note that

[TABLE]

where Ai is the usual Airy function. Using the asymptotic representation of the generalized Airy function

[TABLE]

and the following relationship

[TABLE]

it can be shown that

[TABLE]

Note that, solutions of this kind are used as an initial guess for the full numerical integration.

5 Nonlinear stability

The governing equations and boundary conditions for the two-dimensional motion of two superposed fluids in curvilinear coordinates are

[TABLE]

where the upper and lower case symbols refer to fluids above and below the interface, respectively, and $\rho=\rho_{a}/\rho_{w}$; $\mu=\mu_{a}/\mu_{w}$; $\nu=\nu_{a}/\nu_{w}$ and $d=D_{a}/D_{w}$, being the ratio of the height above and below the interface. Here, $\rho$, $\mu$, and $\nu$ are dimensionless density, dynamic viscosity and kinematic viscosity, respectively, with understanding that the suffices $a$ and $w$ refer to air and water. Also, ${\mathscr{F}}=gD_{a}/U_{\infty}$, ${\mathscr{T}}=S/\rho_{w}U_{\infty}^{2}D_{a}$ and $R=U_{\infty}D_{a}/\nu_{a}$ are the Froude, Weber and Reynolds numbers respectively, with $g$ denoting the acceleration due to gravity, $S$ the coefficient of surface tension, $U_{\infty}$ the free stream speed of the air, and $\eta$ is the interface profile which is the Stokes wave.

For the disturbance near theocratical conditions we expand the frequency and eigenfunction as a Taylor series. Thus

[TABLE]

Now, since $d(kc)/dk|_{k_{c}}=-a_{1}$, then $-a_{1}$ is the group velocity of the disturbance; $a_{2}$ is related to the curvature of the nose of the neutral curve and $d_{1}$ is the exponential growth of the wave for $R>R_{c}$. Likewise we expand $\Phi$ in the form

[TABLE]

with a similar expansion for $\phi$ in the lower fluid.

It is well known that the basic flow becomes unstable to travelling wave disturbances when $R=R_{c}$ [3,8]. To examine the growth of these Tollmien-Schlichting waves in the vicinity of $R_{c}$ we introduce the small parameter

[TABLE]

and scaled variable

[TABLE]

The continuity equations in (222) are satisfied by the streamfunctions $\Lambda$ and $\lambda$ for the upper and lower fluids, respectively, thus

[TABLE]

Introducing the expansions

[TABLE]

where c.c. denotes complex conjugate, and ${\sf E}=\exp[ik_{c}(\chi-c_{r}t)].$ We further expand the functions in (227) as

[TABLE]

The pressure and interface displacement are expanded as

[TABLE]

where

[TABLE]

Substituting in (222) and equating powers of $O(\varepsilon^{\tiny\frac{1}{2}}{\sf E})$ we obtain the coupled Orr-Sommerfeld problem

[TABLE]

Note that the system (254) is the same as (84) with $R=R_{c}$, $k=k_{c}$ and $c=c_{r}$. The solution of (254) may be written as

[TABLE]

where $A(\tau,\xi)$ is the amplitude function to be determined.

The distortion of the mean flow occurs at $O(\varepsilon{\sf E}^{0})$, which is

[TABLE]

where

[TABLE]

and a tilde denotes complex conjugate. The system (254) has unique solution of the form

[TABLE]

where

[TABLE]

The constants $B,C,b$ and $c$ are found by applying the interface conditions.

At $O(\varepsilon{\sf{E}}^{2})$ we obtain the system governing the second harmonic, namely

[TABLE]

The solution of (325) may be expressed as

[TABLE]

with the assumption that the solution is unique.

We emphasize that the adjoint equation in the upper and he lower fluid layer is the adjoint of the Orr-Sommerfeld equation for a homogeneous fluid. The adjoint system possesses the following relationship

[TABLE]

where the operators ${\mathscr{L}}$ and ${\mathfrak{l}}$ are the same as that given in (254), and

[TABLE]

The system of equations for $O(\varepsilon{\sf E})$ is given by

[TABLE]

where UF and LF stands for upper and lower fluid, respectively. The solution of the $O(\varepsilon{\sf E})$ system is

[TABLE]

where $A_{2}$ is and arbitrary function of $\tau$ and $\xi$, and ${\it\Psi}_{10}$ and ${\it\psi}_{10}$ are the solutions of the following system

[TABLE]

since the adjoint to the system (84), at critical conditions, satisfies

[TABLE]

Next, at $O\left(\varepsilon^{\tiny\frac{3}{2}}{\sf E}\right)$, we obtain the system of equations

[TABLE]

where

[TABLE]

Proceeding as above we obtain the following set of equations at $O\left(\varepsilon^{\tiny\frac{5}{2}}{\sf E}^{2}\right)$.

[TABLE]

where

[TABLE]

and

[TABLE]

In (425), the remaining coefficients are given by Blennerhassett [2].

This system of equations has a unique solution provided the inhomogeneous terms are orthogonal to the adjoint function. Applying the solvability condition yields the amplitude equation (1).

The complex coefficients $\kappa,a_{2}$, $d_{1}$ and $\varpi$ in equation (1) are given by

[TABLE]

where $\Phi,\phi$ satisfies the Orr-Sommerfeld equation for the basic flow with $\Psi,\psi$ satisfying their adjoints.

It is important to note that in the formulation presented here quantities such as $\Phi$ and $\phi$, and their adjoints, are functions of $\eta$ alone, but, by virtue of (17) they are functions of both spacial coordinates $x$, $y$ and time $t$. Also, since the Jacobian of the transformation is $1+O(ak)$, the resulting equations expressed in curvilinear coordinates are, to first order in $ak$, symbollically identical to their Caretesian counterpart. Thus, the solutions for the basic flow, in curvilinear coordinates, are the same as those given by (10) and (489), in Cartesian coordinates, with $y$ simply replaced with $\eta$. The error in this approximation is $O(ak)^{2}$ and may be neglected since $ak\ll 1$.

6 Numerical Schemes

The complex cubic–quintic nonlinear Schrödinger equation (1) is solved numerically, using a second-order central difference for the spacial part and the third-order explicit Runge-Kutta method for the transient part, whose complex coefficients are determined from the system of ordinary differential equations (254), (329) and (329) which coupled via interface conditions which are solved numerically using the finite difference scheme as described below. First, we rewrite the system in a compact form as

[TABLE]

where $\Phi=({\it\Phi}_{11},{\it\Phi}_{22},\Psi)$, $\phi=(\phi_{11},\phi_{22},\psi)$ and the definitions of the coefficients $c$, $C$, $d$, $D$, $e$, $E$, $FL_{np}$, $FU_{np}$, $FLNH_{n}$, $FUNH_{n}$, $f$, $F$ are those given by system of equation (254), (329) and (395), respectively.

We approximate

[TABLE]

where $\delta$ is the usual central difference operator, i.e., $\delta g_{i}=g_{i+1/2}-g_{i-1/2}$. Then

[TABLE]

Here $h=\eta_{i}-\eta_{i-1}$ and $\mu$ is defined as $\mu g_{i}=\frac{1}{2}\left(g_{i+1/2}+g_{i-1/2}\right)$.

For $\phi^{\prime\prime\prime}$ we use the following approximation

[TABLE]

Using (474) and the appropriate equation of motion from (466), we have

[TABLE]

Hence we maintain $O\left(h^{4}\right)$ accuracy and a five point formula for the third derivative in the interface. The similar expressions result for the derivatives of $\Phi$. Substituting (471)–(473) in (466), we obtain the following finite difference approximation.

[TABLE]

where

[TABLE]

where the subscript $i$ implies evaluation at $\eta=\eta_{i}$. Similar expressions result for the approximation to $\Phi$. Thus, the finite difference approximation to the coupled ODEs system with the boundary conditions and the interface condition (466) leads to the matrix equation ${\bf{\mathscr{A}}x}={\bf b}$ with the following form:

[TABLE]

where the entries of the matrix ${\bf\mathscr{A}}$ and unknown vector ${\bf x}$ are complex numbers.

6.1 Complex cubic–quintic nonlinear Schrödinger equation

Let $h=(b-a)/n$ be the grid spacing and $x_{j}=a+jh,\,j=0,1,\cdots,n$ be the grid points. Define $A_{j}(t)$ as an approximation to $A(x_{j},t),\,j=0,1,\cdots,n$, ${\bf A}(t)=(A_{1}(t),\cdots,A_{n-1}(t))^{\sf T}$ and ${\bf A}_{0}(x)=(A_{0}(x_{1}),\cdots,A_{0}(x_{n-1}))^{\sf T}$. Let $F(A_{j}(t))=\displaystyle\frac{d_{1}}{d_{1r}}A_{j}(t)+\kappa A_{j}(t)|A_{j}(t)|^{2}+\varpi A_{j}(t)|A_{j}(t)|^{4}$. Using the central difference approximation to $\partial_{x}^{2}A(x_{j},t)$, we write the semi–discretization of the IBVP for equation (1) as the following system of ODEs

[TABLE]

where ${\bf F}({\bf A}(t))=\left(F(A_{1}(t)),\,\cdots,\,F(A_{n}(t))\right)^{\sf T}$. The upper dot in (481) indicates derivative with respect to $t$ and ${\mathscr{B}}$ is the finite difference matrix.

To construct an integration scheme to solve the ODE system (481)–(482), let $t_{n+1}=t_{n}+\Delta t$, and let ${\bf A}^{n}$ denote the value of the variable ${\bf A}$ at time $t_{n}$. Employing a low storage variant third–order Runge–Kutta scheme [12], we write the fully discrete system as

[TABLE]

where ${\bf G}({\bf A}^{n})={\mathscr{B}}{\bf A}^{n}+{\bf F}({\bf A}^{n})$.

The numerical scheme is parallelized for distributed memory clusters of processors or heterogeneous networked computers using the MPI (Message Passing Interface) library and implemented in FORTRAN.

7 Results

For the numerical simulation reported in this paper, we choose ${\mathscr{F}}=10^{3}k_{s}^{3}$ and ${\mathscr{T}}=10^{3}(1-\rho)k_{s}$ for the present case of wind blowing over waves. Here the quantity $k_{s}$ is defined as

[TABLE]

with $\rho_{a}=1.225\times 10^{-3}\,{\rm g/cm^{3}}$, $\rho_{w}=0.9991\,{\rm g/cm^{3}}$, $g=981\,{\rm cm/s^{2}}$, $S=73.5\,{\rm dyn/cm}$, and we choose $D_{\ell}=D_{u}$ such that $d=1$. The kinematic viscosities are taken to be $\nu_{a}=0.145\,{\rm cm^{2}/s}$ and $\nu_{w}=1.138\times 10^{-2}\,{\rm cm^{2}/s}$.

The velocity profiles for the basic flow, for the case of turbulent boundary layer, are those given by (8). The velocity profiles for pPf and pCf, for the upper fluid, satisfy

[TABLE]

The solution of this system is given by [2]

[TABLE]

where

[TABLE]

and

[TABLE]

Note that, for pPf $U_{p}=1$ which implies that he upper boundary is stationary, and the flow is driven by the pressure gradient. In contrast, for pCf $U_{p}=0$ and the motion is generated soley by the relative motion of the boundaries.

Computations were performed on a Linux cluster (zeus.db.erau.edu: 256 Intel Xeon 3.2GHz 1024 KB cache 4GB with Myrinet MX, GNU Linux) at Embry-Riddle University. In all simulations presented here, we use the spacial domain of $0\leq\xi\leq 3$ and temporal domain $0\leq\tau\leq 15$, the typical values of the parameters $a=0.3$, $l=0.5$ for the amplitude and the wavelength of the Stokes wave profile, respectively, and the other parameters of the model (1) are given in Table 1 below.

Figure 2 shows the surface structure for pCf (left) and pPf (right). For the case of pPf we observe the classical ‘snake’ type solitons [9], whereas for pCf we see a different structure similar to colliding solitons which, in time, bifurcates to weaker soliton structure.

In figure 3 the surface structure for LBL (left) and TBL (right) are depicted. As can be seen, the soliton structures for TBL resembles the chaotic solitons [9]. This is, of course, not surprising as turbulence generates more complex surface agitation compared to its laminar counterpart. We may term this as a violent surface motion.

The most intriguing results obtained is that shown in figure 3 for LBL. This results clearly confirms the generation of solitons by wind, which has commonly been observed is shallow waters, c.f. figure 4.

8 Conclusions

In this paper we presented an asymptotic theory for linear stability of non-linear air-sea interfaces by taking into account the motion in the air and the water. The two motions are then matched at the free surface via the usual kinematic, dynamic and free stress boundary conditions. The flow above water is assumed to be turbulent, with a mean logarithmic velocity profile. On the other hand, a plane Poiseuille flow is assumed for the motion in the water. Also presented, is the linear stability calculations of two fluid interface which suggest the flow become unstable for wavenumbers above that given by equation (23). Finally, the growth of surface waves by the shear wind blowing over is calculated. It is shown and extra term arises in the expression for the dimensionless energy transfer parameter which is due to the coupling between wind and surface waves.

In the second part of this paper [3], we use a non-linear hydrodynamic stability [7] to derive the quintic non-linear Schrödinger equation (1) with complex coefficients. We then use the above analysis to determine the complex amplitude of the equation (1).

Using a nonlinear stability analysis of a coupled air-sea system is performed. An amplitude equation, namely the quintic nonlinear Schrödinger equation with complex coefficients is derived. The complex coefficients are obtained from the linear stability of the same system in the neighborhood of the critical Reynolds number. Four different profiles for the shear flow for the airflow over the interface is considered, namely, plane Couette flow (pCf), plane Poiseuille flow (pPf), laminar and turbulent boundary layer (L,TBL) profiles. For each of the above cases the shear flow counterpart in the water is assumed to be a pPf.

It is shown that the above amplitude equation produces ‘snake’ solitons for pCf, pPf and LBL profiles, with striking similarities. On the other hand, for TBL we observe a very violent surface motion. For cases of pCf and LBL remarkable similarity is observed with observations made at sea.

We conclude that the effect of nonlinearity in the airflow over the sea surface is much larger than nonlinear interactions in the water, and hence it is not possible to decouple the motion in the air and the water for finite amplitude wind-wave interactions. We remark, under the conditions shown in this paper, the generation of solitons by wind is possible particularly in shallow waters.