2.2.1. Governing equations

Assuming that both the liquid and gas are initially irrotational, the non-dimensional velocity potentials of the liquid and gas phases are defined as ${\it\phi}_{l}$ and ${\it\phi}_{gj}$ , where $j=1$ represents the gas above the sheet and $j=2$ represents the gas below.

The mass conservation equations of the liquid and gas are

(2.2) $$\begin{eqnarray}\displaystyle & {\it\phi}_{l,xx}+{\it\phi}_{l,yy}=0\quad \text{for }-1+{\it\eta}_{2}\leqslant y\leqslant 1+{\it\eta}_{1}, & \displaystyle\end{eqnarray}$$

(2.3) $$\begin{eqnarray}\displaystyle & {\it\phi}_{gj,xx}+{\it\phi}_{gj,yy}=0\quad j=1,2\text{ for }y\geqslant 1+{\it\eta}_{1},y\leqslant -1+{\it\eta}_{2}. & \displaystyle\end{eqnarray}$$

(2.4) $$\begin{eqnarray}\displaystyle & {\it\phi}_{l,y}-{\it\eta}_{j,t}-{\it\phi}_{l,x}{\it\eta}_{j,x}=0\quad \text{at }y=(-1)^{j+1}+{\it\eta}_{j}, & \displaystyle\end{eqnarray}$$

(2.5) $$\begin{eqnarray}\displaystyle & {\it\phi}_{gj,y}-{\it\eta}_{j,t}-{\it\phi}_{gj,x}{\it\eta}_{j,x}=0\quad \text{at }y=(-1)^{j+1}+{\it\eta}_{j}. & \displaystyle\end{eqnarray}$$

(2.6) $$\begin{eqnarray}\displaystyle & & \displaystyle \frac{1}{2}+{\it\rho}{\it\phi}_{gj,t}-{\it\phi}_{l,t}+\frac{1}{2}{\it\rho}({\it\phi}_{gj,x}^{2}+{\it\phi}_{gj,y}^{2})-\frac{1}{2}({\it\phi}_{l,x}^{2}+{\it\phi}_{l,y}^{2})-\frac{(-1)^{j}{\it\eta}_{j,xx}}{\mathit{We}(1+{\it\eta}_{j,x}^{2})^{1.5}}=0\nonumber\\ \displaystyle & & \displaystyle \quad \text{at }y=(-1)^{j+1}+{\it\eta}_{j}.\end{eqnarray}$$

A second-order weakly nonlinear analysis is used to solve the nonlinear equations. Assuming that the nonlinearity is weak, all the variables are expressed by a second-order expansion in power series of ${\it\eta}_{0}$ :

(2.7) $$\begin{eqnarray}({\it\phi}_{l},{\it\phi}_{gj},{\it\eta}_{j})=(x,0,0)+(\text{}^{1}{\it\phi}_{l},\text{}^{1}{\it\phi}_{gj},\text{}^{1}{\it\eta}_{j}){\it\eta}_{0}+(\text{}^{2}{\it\phi}_{l},\text{}^{2}{\it\phi}_{gj},\text{}^{2}{\it\eta}_{j}){\it\eta}_{0}^{2}+O({\it\eta}_{0}^{3}).\end{eqnarray}$$

Note that the terms in the first bracket on the right-hand side of (2.7) represent the basic flow field. For the variables on the disturbed boundaries, their values are obtained by a Taylor series using their values on the undisturbed boundaries, for example,

(2.8) $$\begin{eqnarray}{{\it\phi}_{l,t}|}_{y=(-1)^{j+1}+{\it\eta}_{j}}={{\it\phi}_{l,t}|}_{y=(-1)^{j+1}}+{\it\eta}_{j}{{\it\phi}_{l,yt}|}_{y=(-1)^{j+1}}+\cdots \,.\end{eqnarray}$$

Substituting (2.7) and (2.8) into (2.1)–(2.6), collecting the terms with ${\it\eta}_{0}$ and ${\it\eta}_{0}^{2}$ , the first-order and the second-order governing equations are achieved.

The first-order equations are

(2.9) $$\begin{eqnarray}\displaystyle & {\text{}^{1}{\it\eta}_{j}|}_{t=0}=\hat{{\it\eta}}_{j}\exp (\text{i}kx)+\text{c.c.}, & \displaystyle\end{eqnarray}$$

(2.10) $$\begin{eqnarray}\displaystyle & \text{}^{1}{\it\phi}_{l,xx}+\text{}^{1}{\it\phi}_{l,yy}=0\quad \text{for }-1\leqslant y\leqslant 1, & \displaystyle\end{eqnarray}$$

(2.11) $$\begin{eqnarray}\displaystyle & \text{}^{1}{\it\phi}_{gj,xx}+\text{}^{1}{\it\phi}_{gj,yy}=0\quad j=1,2\text{ for }y\geqslant 1,y\leqslant -1, & \displaystyle\end{eqnarray}$$

(2.12) $$\begin{eqnarray}\displaystyle & \text{}^{1}{\it\phi}_{l,y}-\text{}^{1}{\it\eta}_{j,t}-\text{}^{1}{\it\eta}_{j,x}=0\quad \text{at }y=(-1)^{j+1}, & \displaystyle\end{eqnarray}$$

(2.13) $$\begin{eqnarray}\displaystyle & \text{}^{1}{\it\phi}_{gj,y}-\text{}^{1}{\it\eta}_{j,t}=0\quad \text{at }y=(-1)^{j+1}, & \displaystyle\end{eqnarray}$$

(2.14) $$\begin{eqnarray}\displaystyle & \displaystyle {\it\rho}\text{}^{1}{\it\phi}_{gj,t}-\text{}^{1}{\it\phi}_{l,t}-\text{}^{1}{\it\phi}_{l,x}-(-1)^{j}\frac{\text{}^{1}{\it\eta}_{j,xx}}{\mathit{We}}=0\quad \text{at }y=(-1)^{j+1}. & \displaystyle\end{eqnarray}$$

The second-order equations are

(2.15) $$\begin{eqnarray}\displaystyle & {\text{}^{2}{\it\eta}_{j}|}_{t=0}=0, & \displaystyle\end{eqnarray}$$

(2.16) $$\begin{eqnarray}\displaystyle & \text{}^{2}{\it\phi}_{l,xx}+\text{}^{2}{\it\phi}_{l,yy}=0\quad \text{for }-1\leqslant y\leqslant 1, & \displaystyle\end{eqnarray}$$

(2.17) $$\begin{eqnarray}\displaystyle & \text{}^{2}{\it\phi}_{gj,xx}+\text{}^{2}{\it\phi}_{gj,yy}=0\quad j=1,2\text{ for }y\geqslant 1,y\leqslant -1, & \displaystyle\end{eqnarray}$$

(2.18) $$\begin{eqnarray}\displaystyle & \text{}^{2}{\it\phi}_{l,y}-\text{}^{2}{\it\eta}_{j,t}-\text{}^{2}{\it\eta}_{j,x}=\text{}^{1}{\it\phi}_{l,x}\text{}^{1}{\it\eta}_{j,x}-\text{}^{1}{\it\eta}_{j}\text{}^{1}{\it\phi}_{l,yy}\quad \text{at }y=(-1)^{j+1}, & \displaystyle\end{eqnarray}$$

(2.19) $$\begin{eqnarray}\displaystyle & \text{}^{2}{\it\phi}_{gj,y}-\text{}^{2}{\it\eta}_{j,t}=\text{}^{1}{\it\phi}_{gj,x}\text{}^{1}{\it\eta}_{j,x}-\text{}^{1}{\it\eta}_{j}\text{}^{1}{\it\phi}_{gj,yy}\quad \text{at }y=(-1)^{j+1}, & \displaystyle\end{eqnarray}$$

(2.20) $$\begin{eqnarray}\displaystyle & & \displaystyle {\it\rho}\text{}^{2}{\it\phi}_{gj,t}-\text{}^{2}{\it\phi}_{l,t}-\text{}^{2}{\it\phi}_{l,x}-(-1)^{j}\frac{\text{}^{2}{\it\eta}_{j,xx}}{\mathit{We}}\nonumber\\ \displaystyle & & \displaystyle \quad =\text{}^{1}{\it\eta}_{j}(-{\it\rho}\text{}^{1}{\it\phi}_{gj,yt}+\text{}^{1}{\it\phi}_{l,yt}+\text{}^{1}{\it\phi}_{l,yx})-\frac{1}{2}{\it\rho}(\text{}^{1}{\it\phi}_{gj,x}^{2}+\text{}^{1}{\it\phi}_{gj,y}^{2})+\frac{1}{2}(\text{}^{1}{\it\phi}_{l,x}^{2}+\text{}^{1}{\it\phi}_{l,y}^{2})\nonumber\\ \displaystyle & & \displaystyle \qquad \text{at }y=(-1)^{j+1}.\end{eqnarray}$$

2.2.2. First-order solutions

The first-order equations are linear in $x$ and $t$ ; therefore, the sinuous mode and the varicose mode can be added together directly to give

(2.21) $$\begin{eqnarray}\displaystyle (\text{}^{1}{\it\phi}_{l},\text{}^{1}{\it\phi}_{gj},\text{}^{1}{\it\eta}_{j}) & = & \displaystyle (\text{}^{1,s}\hat{{\it\phi}}_{l},\text{}^{1,s}\hat{{\it\phi}}_{gj},\text{}^{1,s}\hat{{\it\eta}}_{j})\exp (\text{i}kx-\text{i}{\it\omega}_{1s}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{1,v}\hat{{\it\phi}}_{l},\text{}^{1,v}\hat{{\it\phi}}_{gj},\text{}^{1,v}\hat{{\it\eta}}_{j})\exp (\text{i}kx-\text{i}{\it\omega}_{1v}t)+\text{c.c.},\end{eqnarray}$$

where ${\it\omega}_{1s}$ and ${\it\omega}_{1v}$ are the first-order sinuous and varicose complex frequencies; their real parts represent the oscillation frequencies and their imaginary parts represent the temporal growth rates. The superscript ‘ $1,s$ ’ represents the linear sinuous mode and ‘ $1,v$ ’ represents the linear varicose mode. The symbol ‘ $\hat{~}$ ’ denotes the disturbance amplitudes, which are functions of $y$ only and, in particular, $\text{}^{1,s}\hat{{\it\eta}}_{j}$ and $\text{}^{1,v}\hat{{\it\eta}}_{j}$ are constants.

The two unstable modes are solved using standard linear analysis combining (2.21) and (2.10)–(2.14). Here only the results are presented. The upper and lower interface displacements are identical for the sinuous mode and opposite for the varicose mode:

(2.22a,b ) $$\begin{eqnarray}\text{}^{1,s}\hat{{\it\eta}}_{1}=\text{}^{1,s}\hat{{\it\eta}}_{2},\quad \text{}^{1,v}\hat{{\it\eta}}_{1}=-\text{}^{1,v}\hat{{\it\eta}}_{2}.\end{eqnarray}$$

The dispersion relations must be satisfied in linear instability to make the solution non-trivial:

(2.23a,b ) $$\begin{eqnarray}D_{s}({\it\omega}_{1s},k)=0,\quad D_{v}({\it\omega}_{1v},k)=0.\end{eqnarray}$$

Here $D_{s}$ and $D_{v}$ are defined according to the character equations for the linear sinuous and varicose modes:

(2.24a,b ) $$\begin{eqnarray}D_{s}({\it\omega},k)=-{\it\rho}{\it\omega}^{2}+\frac{k^{3}}{\mathit{We}}-({\it\omega}-k)^{2}\tanh k,\quad D_{v}({\it\omega},k)=-{\it\rho}{\it\omega}^{2}+\frac{k^{3}}{\mathit{We}}-({\it\omega}-k)^{2}\coth k.\end{eqnarray}$$

The velocity potentials are

(2.25) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle \text{}^{1,s}\hat{{\it\phi}}_{l}=\text{i}\left(1-\frac{{\it\omega}_{1s}}{k}\right)\frac{\sinh ky}{\cosh k}\text{}^{1,s}\hat{{\it\eta}}_{1},\quad \text{}^{1,s}\hat{{\it\phi}}_{gj}=(-1)^{j+1}\frac{\text{i}{\it\omega}_{1s}}{k}\exp [k+(-1)^{j}ky]\text{}^{1,s}\hat{{\it\eta}}_{1},\\ \displaystyle \text{}^{1,v}\hat{{\it\phi}}_{l}=\text{i}\left(1-\frac{{\it\omega}_{1v}}{k}\right)\frac{\cosh ky}{\sinh k}\text{}^{1,v}\hat{{\it\eta}}_{1},\quad \text{}^{1,v}\hat{{\it\phi}}_{gj}=\frac{\text{i}{\it\omega}_{1v}}{k}\exp [k+(-1)^{j}ky]\text{}^{1,v}\hat{{\it\eta}}_{1}.\end{array}\!\!\right\}\qquad & & \displaystyle\end{eqnarray}$$

Equations (2.22)–(2.25) are the well-known results of standard linear analysis and they are identical to those presented by Squire (Reference Squire1953) and Hagerty & Shea (Reference Hagerty and Shea1955). Note that (2.24) and (2.25) are valid only at $k>0$ ; $k<0$ renders different expressions for the solutions but the same instability behaviour governed by the real parts of the complex solutions; at $k=0$ , all the temporal frequencies and velocity potentials are zero. The unique solution of the linear homogeneous equations is determined by the initial condition. Substituting (2.21) and (2.22) into the initial condition (2.9), the disturbance amplitudes of the first-order sinuous and varicose modes are derived as

(2.26a,b ) $$\begin{eqnarray}\text{}^{1,s}\hat{{\it\eta}}_{j}={\textstyle \frac{1}{2}}(\hat{{\it\eta}}_{1}+\hat{{\it\eta}}_{2}),\quad \text{}^{1,v}\hat{{\it\eta}}_{j}=(-1)^{j+1}{\textstyle \frac{1}{2}}(\hat{{\it\eta}}_{1}-\hat{{\it\eta}}_{2}).\end{eqnarray}$$

To measure the proportion of initial sinuous amplitude in the linear part, we define

(2.27) $$\begin{eqnarray}r_{s}=\frac{|\text{}^{1,s}\hat{{\it\eta}}_{j}|}{|\text{}^{1,s}\hat{{\it\eta}}_{j}|+|\text{}^{1,v}\hat{{\it\eta}}_{j}|}.\end{eqnarray}$$

2.2.3. Second-order solutions

The second-order solutions are assumed as

(2.28) $$\begin{eqnarray}\displaystyle (\text{}^{2}{\it\phi}_{l},\text{}^{2}{\it\phi}_{gj},\text{}^{2}{\it\eta}_{j}) & = & \displaystyle (\text{}^{21,ss}\hat{{\it\phi}}_{l},\text{}^{21,ss}\hat{{\it\phi}}_{gj},\text{}^{21,ss}\hat{{\it\eta}}_{j})\exp (2\text{i}kx-2\text{i}{\it\omega}_{1s}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{21,s\overline{s}}\hat{{\it\phi}}_{l},\text{}^{21,s\overline{s}}\hat{{\it\phi}}_{gj},\text{}^{21,s\overline{s}}\hat{{\it\eta}}_{j})\exp (-\text{i}{\it\omega}_{1s}t+\text{i}\overline{{\it\omega}}_{1s}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{21,vv}\hat{{\it\phi}}_{l},\text{}^{21,vv}\hat{{\it\phi}}_{gj},\text{}^{21,vv}\hat{{\it\eta}}_{j})\exp (2\text{i}kx-2\text{i}{\it\omega}_{1v}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{21,v\overline{v}}\hat{{\it\phi}}_{l},\text{}^{21,v\overline{v}}\hat{{\it\phi}}_{gj},\text{}^{21,v\overline{v}}\hat{{\it\eta}}_{j})\exp (-\text{i}{\it\omega}_{1v}t+\text{i}\overline{{\it\omega}}_{1v}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{21,sv}\hat{{\it\phi}}_{l},\text{}^{21,sv}\hat{{\it\phi}}_{gj},\text{}^{21,sv}\hat{{\it\eta}}_{j})\exp (2\text{i}kx-\text{i}{\it\omega}_{1s}t-\text{i}{\it\omega}_{1v}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{21,s\overline{v}}\hat{{\it\phi}}_{l},\text{}^{21,s\overline{v}}\hat{{\it\phi}}_{gj},\text{}^{21,s\overline{v}}\hat{{\it\eta}}_{j})\exp (-\text{i}{\it\omega}_{1s}t+\text{i}\overline{{\it\omega}}_{1v}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{22,s}\hat{{\it\phi}}_{l},\text{}^{22,s}\hat{{\it\phi}}_{gj},\text{}^{22,s}\hat{{\it\eta}}_{j})\exp (2\text{i}kx-\text{i}{\it\omega}_{2s}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{22,\overline{s}}\hat{{\it\phi}}_{l},\text{}^{22,\overline{s}}\hat{{\it\phi}}_{gj},\text{}^{22,\overline{s}}\hat{{\it\eta}}_{j})\exp (-\text{i}{\it\omega}_{2\overline{s}}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{22,v}\hat{{\it\phi}}_{l},\text{}^{22,v}\hat{{\it\phi}}_{gj},\text{}^{22,v}\hat{{\it\eta}}_{j})\exp (2\text{i}kx-\text{i}{\it\omega}_{2v}t)\nonumber\\ \displaystyle & & \displaystyle +\,(\text{}^{22,\overline{v}}\hat{{\it\phi}}_{l},\text{}^{22,\overline{v}}\hat{{\it\phi}}_{gj},\text{}^{22,\overline{v}}\hat{{\it\eta}}_{j})\exp (-\text{i}{\it\omega}_{2\overline{v}}t)+\text{c.c.}\end{eqnarray}$$

The modes with superscript ‘21’ have the same exponent coefficients as the inhomogeneous terms formed by the first-order products in (2.18)–(2.20), and they represent the transfer of perturbation from the first order to the second order. Note that $\overline{{\it\omega}}_{1s}$ and $\overline{{\it\omega}}_{1v}$ represent the complex conjugates of ${\it\omega}_{1s}$ and ${\it\omega}_{1v}$ , but $\overline{s}$ and $\overline{v}$ are merely symbols indicating that the disturbances are partly generated by the complex conjugates of first-order disturbances. The modes with superscript 22 have second-order frequencies: ${\it\omega}_{2s}$ and ${\it\omega}_{2\overline{s}}$ are the second-order sinuous frequencies, and ${\it\omega}_{2v}$ and ${\it\omega}_{2\overline{v}}$ are the second-order varicose frequencies. Thus they represent the second-order inherent disturbances. Their corresponding wavenumbers are chosen to be $2k$ and 0 (identical to the ‘21’ modes) in order to satisfy the initial condition (2.15).

In (2.28), modes ‘ $21,ss$ ’ and ‘ $21,s\overline{s}$ ’ are generated from first-order sinuous waves, so they represent the first harmonics of the linear sinuous mode. Similarly modes ‘ $21,vv$ ’ and ‘ $21,v\overline{v}$ ’ represent the first harmonics of the linear varicose mode. Modes ‘ $21,sv$ ’ and ‘ $21,s\overline{v}$ ’ are the first harmonics generated by the first-order sinuous and varicose waves together; therefore, they represent the coupling between the linear sinuous and varicose modes in the second order. They are produced only when the linear sinuous and varicose waves are both present, so their existence signifies that the sinuous and varicose waves do not propagate independently when the nonlinearity is taken into account.

The solutions for ‘ $21$ ’ modes are derived by solving inhomogeneous equations. Take the mode ‘ $21,ss$ ’ as an example. Substituting (2.28) into the mass conservation equations (2.16) and (2.17), and collecting terms with coefficient $\exp (2\text{i}kx-2\text{i}{\it\omega}_{1s}t)$ , one arrives at

(2.29) $$\begin{eqnarray}\displaystyle & \text{}^{21,ss}\hat{{\it\phi}_{l}^{\prime \prime }}-4k^{2}(\text{}^{21,ss}\hat{{\it\phi}}_{l})=0\quad \text{for }-1\leqslant y\leqslant 1, & \displaystyle\end{eqnarray}$$

(2.30) $$\begin{eqnarray}\displaystyle & \text{}^{21,ss}\hat{{\it\phi}_{gj}^{\prime \prime }}-4k_{}^{2}(\text{}^{21,ss}\hat{{\it\phi}}_{gj})=0\quad j=1,2\text{ for }y\geqslant 1,y\leqslant -1. & \displaystyle\end{eqnarray}$$

$y\rightarrow \pm \infty$

(2.31) $$\begin{eqnarray}\displaystyle & \text{}^{21,ss}\hat{{\it\phi}}_{l}=A\exp (2ky)+B\exp (-2ky)\quad \text{for }-1\leqslant y\leqslant 1, & \displaystyle\end{eqnarray}$$

(2.32) $$\begin{eqnarray}\displaystyle & \text{}^{21,ss}\hat{{\it\phi}}_{gj}=C_{j}\exp [(-1)^{j}2ky]\quad j=1,2\text{ for }y\geqslant 1,y\leqslant -1, & \displaystyle\end{eqnarray}$$

$A$

$B$

$C_{j}~(j=1,2)$

$\exp (2\text{i}kx-2\text{i}{\it\omega}_{1s}t)$

(2.33) $$\begin{eqnarray}\displaystyle & \text{}^{21,ss}\hat{{\it\phi}}_{l}^{\prime }+2\text{i}({\it\omega}_{1s}-k)\text{}^{21,ss}\hat{{\it\eta}}_{j}=-k^{2}(\text{}^{1,s}\hat{{\it\phi}}_{l})\text{}^{1,s}\hat{{\it\eta}}_{j}-\text{}^{1,s}\hat{{\it\eta}}_{j}\text{}^{1,s}\hat{{\it\phi}}_{l}^{\prime \prime }\quad \text{at }y=(-1)^{j+1}, & \displaystyle\end{eqnarray}$$

(2.34) $$\begin{eqnarray}\displaystyle & \text{}^{21,ss}\hat{{\it\phi}}_{gj}^{\prime }+2\text{i}{\it\omega}_{1s}\text{}^{21,ss}\hat{{\it\eta}}_{j}=-k^{2}(\text{}^{1,s}\hat{{\it\phi}}_{gj})\text{}^{1,s}\hat{{\it\eta}}_{j}-\text{}^{1,s}\hat{{\it\eta}}_{j}\text{}^{1,s}{\hat{{\it\phi}}_{gj}}^{\prime \prime }\quad \text{at }y=(-1)^{j+1}, & \displaystyle\end{eqnarray}$$

(2.35) $$\begin{eqnarray}\displaystyle & & \displaystyle 2\text{i}[-{\it\rho}{\it\omega}_{1s}\text{}^{21,ss}\hat{{\it\phi}}_{gj}+{\it\omega}_{1s}\text{}^{21,ss}\hat{{\it\phi}}_{l}-k(\text{}^{21,ss}\hat{{\it\phi}}_{l})]+(-1)^{j}\frac{4k^{2}(\text{}^{21,ss}\hat{{\it\eta}}_{j})}{\mathit{We}}\nonumber\\ \displaystyle & & \displaystyle \quad =\text{i}\text{}^{1,s}\hat{{\it\eta}}_{j}[{\it\rho}{\it\omega}_{1s}\text{}^{1,s}\hat{{\it\phi}}_{gj}^{\prime }-{\it\omega}_{1s}\text{}^{1,s}{\hat{{\it\phi}}^{\prime }}_{l}+k(\text{}^{1,s}{\hat{{\it\phi}}^{\prime }}_{l})]-\frac{1}{2}{\it\rho}[\text{}^{1,s}\hat{{\it\phi}}_{gj}^{\prime 2}-k^{2}(\text{}^{1,s}\hat{{\it\phi}}_{gj}^{2})]\nonumber\\ \displaystyle & & \displaystyle \qquad +\,\frac{1}{2}[\text{}^{1,s}\hat{{\it\phi}}_{l}^{\prime 2}-k^{2}(\text{}^{1,s}\hat{{\it\phi}}_{l}^{2})]\quad \text{at }y=(-1)^{j+1}.\end{eqnarray}$$

Substituting (2.31) and (2.32) into (2.33)–(2.35), six equations in six unknown constants $A$ , $B$ , $C_{j}~(j=1,2)$ and $\text{}^{21,ss}\hat{{\it\eta}}_{j}~(j=1,2)$ are obtained.

The solvability of the inhomogeneous equations is discussed next according to the Fredholm alternative. In the second-order equations (2.16)–(2.20), the forms of the second-order disturbances, i.e. the linear parts on the left-hand side of the equations, are exactly the same as the forms of the first-order disturbances in the first-order equations (2.10)–(2.14), so the coefficient matrix of $A$ , $B$ , $C_{j}~(j=1,2)$ and $\text{}^{21,ss}\hat{{\it\eta}}_{j}~(j=1,2)$ , referred to as $\unicode[STIX]{x1D648}(2{\it\omega}_{1s},2k)$ , is the same as the coefficient matrix of the corresponding constants in the first-order linear analysis, $\unicode[STIX]{x1D648}({\it\omega}_{1s},k)$ , as long as $2{\it\omega}_{1s}$ and $2k$ are replaced by ${\it\omega}_{1s}$ and $k$ , respectively. The linear analysis indicates that the singularity of $\unicode[STIX]{x1D648}({\it\omega},k)$ , which guarantees the non-trivial solutions of the linear homogeneous equations, is satisfied when and only when the dispersion relation $D_{s}({\it\omega},k)D_{v}({\it\omega},k)=0$ is satisfied. The property of $D_{s}$ and $D_{v}$ governed by (2.24) determines that $D_{s}({\it\omega}_{1s},k)=0$ (the linear dispersion relation required by (2.23)) automatically results in $D_{s}(2{\it\omega}_{1s},2k)\neq 0$ and $D_{v}(2{\it\omega}_{1s},2k)\neq 0$ , and therefore $\unicode[STIX]{x1D648}(2{\it\omega}_{1s},2k)$ is non-singular and thus invertible. So according to the first Fredholm alternative, $A$ , $B$ , $C_{j}~(j=1,2)$ and $\text{}^{21,ss}\hat{{\it\eta}}_{j}~(j=1,2)$ have one and only one solution.

The solvability of the other inhomogeneous modes can be examined by the same procedure. Similar to mode ‘ $21,ss$ ’, the modes ‘ $21,vv$ ’, ‘ $21,sv$ ’ and ‘ $21,s\overline{v}$ ’ are solvable in the whole wavenumber range $k>0$ , but the modes ‘ $21,s\overline{s}$ ’ and ‘ $21,v\overline{v}$ ’ are solvable only at wavenumbers smaller than the cutoff wavenumber $k_{c}$ (the wavenumber corresponding to zero temporal growth rate). This is because, when $k>k_{c}$ , ${\it\omega}_{1s}$ and ${\it\omega}_{1v}$ become real and this makes their corresponding coefficient matrices become $\unicode[STIX]{x1D648}^{\prime }(0,0)$ . (The form of the coefficient matrix for zero wavenumber modes $\unicode[STIX]{x1D648}^{\prime }({\it\omega},0)$ is different from the non-zero wavenumber counterpart $\unicode[STIX]{x1D648}({\it\omega},k)$ because the velocity potential is in the form of $Ay+B$ rather than the exponential terms in (2.31) and (2.32).) The matrix $\unicode[STIX]{x1D648}^{\prime }(0,0)$ is singular and not invertible, and the vector of the inhomogeneous terms generated by the first order, referred to as $\boldsymbol{b}$ , does not guarantee the second Fredholm alternative, that is, $\boldsymbol{b}$ is orthogonal to all eigenvectors of $\unicode[STIX]{x1D648}^{\prime }(0,0)^{\text{T}}$ with zero eigenvalue, so the inhomogeneous equations are not solvable. It can be confirmed from (A 4) and (A 8) in appendix A that the solution of $\text{}^{21,s\overline{s}}\hat{{\it\phi}}_{gj}$ and $\text{}^{21,v\overline{v}}\hat{{\it\phi}}_{gj}$ does not exist when ${\it\omega}_{1s}$ and ${\it\omega}_{1v}$ are real, which happens at $k>k_{c}$ . Therefore, for the present study, the solvability condition requires that the wavenumber $k$ should be smaller than the cutoff wavenumber $k_{c}$ . This is also desirable in the current analysis because, in this wavenumber range, the temporal growth rates $\text{Im}({\it\omega}_{1s})$ and $\text{Im}({\it\omega}_{1v})$ are positive, representing the instability of the flow field. In order to extend the wavenumber range to $k>k_{c}$ , other measures, such as expanding the cutoff wavenumber $k_{c}$ in power series of ${\it\eta}_{0}$ , may be required (Yuen Reference Yuen1968), but that is not the focus of the present work.

Back to the solution of mode ‘ $21,ss$ ’. Knowing that the solution is unique, the following parity analysis is helpful. If the relations

(2.36a−c ) $$\begin{eqnarray}\text{}^{21,ss}\hat{{\it\eta}}_{1}=-\text{}^{21,ss}\hat{{\it\eta}}_{2},\quad A=B,\quad C_{1}=C_{2},\end{eqnarray}$$

exist, then the boundary conditions (2.33)–(2.35) will be identical at $y=1$ and $y=-1$ , and solving only those at $y=1$ is enough. Thus $\text{}^{21,ss}\hat{{\it\eta}}_{1}$ is solved as

(2.37) $$\begin{eqnarray}\text{}^{21,ss}\hat{{\it\eta}}_{1}=\frac{k[({\it\omega}_{1s}-k)^{2}(3-\tanh ^{2}k-4\tanh k\coth 2k)+2{\it\rho}{\it\omega}_{1s}^{2}]}{D_{v}(2{\it\omega}_{1s},2k)}\text{}^{1,s}\hat{{\it\eta}}_{1}^{2}.\end{eqnarray}$$

The solutions for $\text{}^{21,ss}\hat{{\it\phi}}_{l}$ and $\text{}^{21,ss}\hat{{\it\phi}}_{gj}$ are given in appendix A.

Using the same solution procedure for other ‘21’ modes, it is easy to deduce the relations between the upper and lower interface displacements via parity analysis:

(2.38) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\text{}^{21,s\overline{s}}\hat{{\it\eta}}_{1}=-\text{}^{21,s\overline{s}}\hat{{\it\eta}}_{2},\quad \text{}^{21,vv}\hat{{\it\eta}}_{1}=-\text{}^{21,vv}\hat{{\it\eta}}_{2},\quad \text{}^{21,v\overline{v}}\hat{{\it\eta}}_{1}=-\text{}^{21,v\overline{v}}\hat{{\it\eta}}_{2},\\ \text{}^{21,sv}\hat{{\it\eta}}_{1}=\text{}^{21,sv}\hat{{\it\eta}}_{2},\quad \text{}^{21,s\overline{v}}\hat{{\it\eta}}_{1}=\text{}^{21,s\overline{v}}\hat{{\it\eta}}_{2}.\end{array}\!\!\right\} & & \displaystyle\end{eqnarray}$$

The solutions for the upper interface are

(2.39) $$\begin{eqnarray}\displaystyle & \text{}^{21,s\overline{s}}\hat{{\it\eta}}_{1}=0, & \displaystyle\end{eqnarray}$$

(2.40) $$\begin{eqnarray}\displaystyle & \displaystyle \text{}^{21,vv}\hat{{\it\eta}}_{1}=\frac{k[({\it\omega}_{1v}-k)^{2}(3-\coth ^{2}k-4\coth k\coth 2k)+2{\it\rho}{\it\omega}_{1v}^{2}]}{D_{v}(2{\it\omega}_{1v},2k)}\text{}^{1,v}\hat{{\it\eta}}_{1}^{2}, & \displaystyle\end{eqnarray}$$

(2.41) $$\begin{eqnarray}\displaystyle & \text{}^{21,v\overline{v}}\hat{{\it\eta}}_{1}=0, & \displaystyle\end{eqnarray}$$

(2.42) $$\begin{eqnarray}\displaystyle \text{}^{21,sv}\hat{{\it\eta}}_{1} & = & \displaystyle \frac{2k}{D_{s}({\it\omega}_{1s}+{\it\omega}_{1v},2k)}\{\!({\it\omega}_{1s}-k)^{2}+({\it\omega}_{1v}-k)^{2}-({\it\omega}_{1s}+{\it\omega}_{1v}-2k)\nonumber\\ \displaystyle & & \displaystyle \times \,[({\it\omega}_{1s}-k)\tanh k+({\it\omega}_{1v}-k)\coth k]\tanh 2k+2{\it\rho}{\it\omega}_{1s}{\it\omega}_{1v}\!\}\text{}^{1,s}\hat{{\it\eta}}_{1}\text{}^{1,v}\hat{{\it\eta}}_{1},\qquad\end{eqnarray}$$

(2.43) $$\begin{eqnarray}\text{}^{21,s\overline{v}}\hat{{\it\eta}}_{1}=0.\end{eqnarray}$$

Their corresponding velocity potentials are given in appendix A. It can be seen that $\text{}^{21,ss}\hat{{\it\eta}}_{j}$ and $\text{}^{21,s\overline{s}}\hat{{\it\eta}}_{j}$ are identical to $\hat{{\it\eta}}_{j21}$ and $\hat{{\it\eta}}_{j22}$ respectively (the corresponding harmonics of linear sinuous mode for temporal instability) in the study of Yang et al. (Reference Yang, Wang, Fu, Du and Tong2013) when viscosity vanishes.

One important characteristic of the solutions is that the zero wavenumber interface displacements, i.e. $\text{}^{21,s\overline{s}}\hat{{\it\eta}}_{j}$ , $\text{}^{21,v\overline{v}}\hat{{\it\eta}}_{j}$ and $\text{}^{21,s\overline{v}}\hat{{\it\eta}}_{j}$ , are zero, indicating that up to the second order there is no unstable mode that makes the sheet dilate or attenuate uniformly everywhere. This is a special property of the temporal instability of planar liquid sheets, which can also be proved by the study of Clark & Dombrowski (Reference Clark and Dombrowski1972), and it is very different from the nonlinear instability of cylindrical jets (Yuen Reference Yuen1968). Another important characteristic is that the first harmonic of the linear sinuous mode $\text{}^{21,ss}\hat{{\it\eta}}_{j}$ and the first harmonic of the linear varicose mode $\text{}^{21,vv}\hat{{\it\eta}}_{j}$ are varicose, but the first harmonic generated by the coupling between the linear sinuous and varicose modes, $\text{}^{21,sv}\hat{{\it\eta}}_{j}$ , is sinuous.

Since $\text{}^{21,ss}\hat{{\it\eta}}_{1}$ , $\text{}^{21,vv}\hat{{\it\eta}}_{1}$ and $\text{}^{21,sv}\hat{{\it\eta}}_{1}$ are respectively in proportion to $\text{}^{1,s}\hat{{\it\eta}}_{1}^{2}$ , $\text{}^{1,v}\hat{{\it\eta}}_{1}^{2}$ and $\text{}^{1,s}\hat{{\it\eta}}_{1}\text{}^{1,v}\hat{{\it\eta}}_{1}$ , from which they are generated, it is useful define the ‘relative second-order disturbance amplitudes’ as

(2.44a−c ) $$\begin{eqnarray}\text{}^{2r,ss}\hat{{\it\eta}}=\frac{\text{}^{21,ss}\hat{{\it\eta}}_{1}}{\text{}^{1,s}\hat{{\it\eta}}_{1}^{2}},\quad \text{}^{2r,vv}\hat{{\it\eta}}=\frac{\text{}^{21,vv}\hat{{\it\eta}}_{1}}{\text{}^{1,v}\hat{{\it\eta}}_{1}^{2}},\quad \text{}^{2r,sv}\hat{{\it\eta}}=\frac{\text{}^{21,sv}\hat{{\it\eta}}_{1}}{\text{}^{1,s}\hat{{\it\eta}}_{1}\text{}^{1,v}\hat{{\it\eta}}_{1}}.\end{eqnarray}$$

Here $\text{}^{2r,ss}\hat{{\it\eta}}$ , $\text{}^{2r,vv}\hat{{\it\eta}}$ and $\text{}^{2r,sv}\hat{{\it\eta}}$ represent the ability of each linear mode to generate the corresponding nonlinear harmonic; their norms are a good indicator of the initial (at $t=0,$ so the effect of temporal growth rate is excluded) extent of nonlinearity, especially when $|\text{}^{1,s}\hat{{\it\eta}}_{j}|=|\text{}^{1,v}\hat{{\it\eta}}_{j}|$ .

The equations of modes ‘ $22$ ’ do not contain the first-order inhomogeneous products, so they are solved using the same method as standard linear analysis. The second-order temporal frequencies are determined by the dispersion relations

(2.45a−d ) $$\begin{eqnarray}D_{s}({\it\omega}_{2s},2k)=0,\quad D_{s}({\it\omega}_{2\overline{s}},0)=0,\quad D_{v}({\it\omega}_{2v},2k)=0,\quad D_{v}({\it\omega}_{2\overline{v}},0)=0.\end{eqnarray}$$

Similar to linear analysis, the sinuous and varicose disturbance amplitudes yield

(2.46a−d ) $$\begin{eqnarray}\text{}^{22,s}\hat{{\it\eta}}_{1}=\text{}^{22,s}\hat{{\it\eta}}_{2},\quad \text{}^{22,\overline{s}}\hat{{\it\eta}}_{1}=\text{}^{22,\overline{s}}\hat{{\it\eta}}_{2},\quad \text{}^{22,v}\hat{{\it\eta}}_{1}=-\text{}^{22,v}\hat{{\it\eta}}_{2},\quad \text{}^{22,\overline{v}}\hat{{\it\eta}}_{1}=-\text{}^{22,\overline{v}}\hat{{\it\eta}}_{2}.\end{eqnarray}$$

Substituting (2.28), (2.36a ), (2.38) and (2.46) into the second-order initial condition (2.15), $\text{}^{22,s}\hat{{\it\eta}}_{j}$ , $\text{}^{22,\overline{s}}\hat{{\it\eta}}_{j}$ , $\text{}^{22,v}\hat{{\it\eta}}_{j}$ and $\text{}^{22,\overline{v}}\hat{{\it\eta}}_{j}$ are solved as

(2.47) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\text{}^{22,s}\hat{{\it\eta}}_{j}=-\text{}^{21,sv}\hat{{\it\eta}}_{j},\quad \text{}^{22,\overline{s}}\hat{{\it\eta}}_{j}=-\text{}^{21,s\overline{v}}\hat{{\it\eta}}_{j}=0,\\ \text{}^{22,v}\hat{{\it\eta}}_{j}=-\text{}^{21,ss}\hat{{\it\eta}}_{j}-\text{}^{21,vv}\hat{{\it\eta}}_{j},\quad \text{}^{22,\overline{v}}\hat{{\it\eta}}_{j}=-\text{}^{21,s\overline{s}}\hat{{\it\eta}}_{j}-\text{}^{21,v\overline{v}}\hat{{\it\eta}}_{j}=0.\end{array}\!\!\right\} & & \displaystyle\end{eqnarray}$$

The solutions of the velocity potentials of the ‘22’ modes are given in appendix A.

2.2.4. Summary

Removing the modes whose interface displacements are zero and neglecting the remainder term $O({\it\eta}_{0}^{3})$ , the second-order approximate interface displacement is summarized as

(2.48) $$\begin{eqnarray}\displaystyle {\it\eta}_{j}=\text{}^{1,s}{\it\eta}_{j}+\text{}^{1,v}{\it\eta}_{j}+\text{}^{21,ss}{\it\eta}_{j}+\text{}^{21,vv}{\it\eta}_{j}+\text{}^{21,sv}{\it\eta}_{j}+\text{}^{22,s}{\it\eta}_{j}+\text{}^{22,v}{\it\eta}_{j}, & & \displaystyle\end{eqnarray}$$

with

(2.49a ) $$\begin{eqnarray}\displaystyle \text{}^{1,s}{\it\eta}_{j} & = & \displaystyle \text{}^{1,s}\hat{{\it\eta}}_{j}{\it\eta}_{0}\exp (\text{i}kx-\text{i}{\it\omega}_{1s}t)+\text{c.c.},\end{eqnarray}$$

(2.49b ) $$\begin{eqnarray}\displaystyle \text{}^{1,v}{\it\eta}_{j} & = & \displaystyle \text{}^{1,v}\hat{{\it\eta}}_{j}{\it\eta}_{0}\exp (\text{i}kx-\text{i}{\it\omega}_{1v}t)+\text{c.c.},\end{eqnarray}$$

(2.49c ) $$\begin{eqnarray}\displaystyle \text{}^{21,ss}{\it\eta}_{j} & = & \displaystyle \text{}^{21,ss}\hat{{\it\eta}}_{j}{\it\eta}_{0}^{2}\exp (2\text{i}kx-2\text{i}{\it\omega}_{1s}t)+\text{c.c.},\end{eqnarray}$$

(2.49d ) $$\begin{eqnarray}\displaystyle \text{}^{21,vv}{\it\eta}_{j} & = & \displaystyle \text{}^{21,vv}\hat{{\it\eta}}_{j}{\it\eta}_{0}^{2}\exp (2\text{i}kx-2\text{i}{\it\omega}_{1v}t)+\text{c.c.},\end{eqnarray}$$

(2.49e ) $$\begin{eqnarray}\displaystyle \text{}^{21,sv}{\it\eta}_{j} & = & \displaystyle \text{}^{21,sv}\hat{{\it\eta}}_{j}{\it\eta}_{0}^{2}\exp (2\text{i}kx-\text{i}{\it\omega}_{1s}t-\text{i}{\it\omega}_{1v}t)+\text{c.c.,}\end{eqnarray}$$

(2.49f ) $$\begin{eqnarray}\displaystyle \text{}^{22,s}{\it\eta}_{j} & = & \displaystyle \text{}^{22,s}\hat{{\it\eta}}_{j}{\it\eta}_{0}^{2}\exp (2\text{i}kx-\text{i}{\it\omega}_{2s}t)+\text{c.c.},\end{eqnarray}$$

(2.49g ) $$\begin{eqnarray}\displaystyle \text{}^{22,v}{\it\eta}_{j} & = & \displaystyle \text{}^{22,v}\hat{{\it\eta}}_{j}{\it\eta}_{0}^{2}\exp (2\text{i}kx-\text{i}{\it\omega}_{2v}t)+\text{c.c.}\end{eqnarray}$$

3.2.1. Effect of wavenumber

The complex temporal frequency as a function of wavenumber $k$ , i.e. the dispersion relation, is given in figure 9. For most of the wavenumber range, the linear sinuous mode has a much larger temporal growth rate $\text{Im}({\it\omega})$ , and the phase velocity represented by $\text{Re}({\it\omega})/k$ is different for the linear sinuous and varicose modes. Figure 9 also reveals that the linear varicose mode has greater temporal growth rates at very large wavenumbers, and that it has a larger cutoff wavenumber (the wavenumber corresponding to zero temporal growth rate) for this mode.

Figure 9.

Dispersion relation: $\mathit{We}=5$ , ${\it\rho}=0.1$ .

Figure 10.

Relative second-order disturbance amplitudes at different wavenumbers: $\mathit{We}=5$ , ${\it\rho}=0.1$ .

The relative second-order disturbance amplitudes, which indicate the second-order initial amplitudes, are given in figure 10. For most wavenumbers, $|\text{}^{2r,vv}\hat{{\it\eta}}|$ is the largest and $|\text{}^{2r,ss}\hat{{\it\eta}}|$ is the smallest, so the first harmonic of the linear varicose mode ‘ $21,vv$ ’ has the largest initial amplitude and the first harmonic of the linear sinuous mode ‘ $21,ss$ ’ has the smallest initial amplitude. Figure 10 also shows that, while $|\text{}^{2r,ss}\hat{{\it\eta}}|$ increases monotonically with increase of wavenumber $k$ , indicating that short waves have stronger nonlinear effects, both $|\text{}^{2r,vv}\hat{{\it\eta}}|$ and $|\text{}^{2r,sv}\hat{{\it\eta}}|$ have a peak at a certain wavenumber.

Our examination shows that, unfortunately, in the wavenumber range $0<k<0.39$ where both modes grow (required by the solvability condition discussed in § 2.2.3), only at approximately $k>0.24$ does the weakly nonlinear analysis agree well with numerical simulation when the sheet breaks up. For $0<k<0.1$ , the small values of $\text{Im}({\it\omega}_{1v})$ and $|\text{}^{2r,ss}\hat{{\it\eta}}|$ make the breakup time very long, resulting in extremely large wave amplitudes (dominated by the linear sinuous mode ‘ $1,s$ ’). This makes the nonlinearity significantly strong and the weakly nonlinear analysis is no longer appropriate to study the breakup. For $0.1<k<0.24$ , the large values of $|\text{}^{2r,sv}\hat{{\it\eta}}|$ and $|\text{}^{2r,vv}\hat{{\it\eta}}|$ bring about large initial second-order velocities of the weakly nonlinear analysis, making it less precise at $t=0$ . Consequently, the initial vortex sheet strengths of the numerical simulation given by (2.58) and (2.59) do not induce an initial velocity field that agrees well with that of the weakly nonlinear analysis, and this certainly brings about significant disparities between the analytical and numerical results in the subsequent disturbance evolution. Only at $k>0.24$ are all the first-order and second-order amplitudes not very large when the sheet breaks up, and the weakly nonlinear analysis agrees qualitatively well with numerical simulation. The present study only focuses on this range.

Figure 11.

Shapes of sheets at different wavenumbers: (a,b) the results of numerical simulation and (c,d) the results of weakly nonlinear analysis; $\mathit{We}=5$ , ${\it\rho}=0.1$ , ${\it\eta}_{0}=0.2$ ; (a,c)  $k=0.25$ , $t=57.5$ and (b,d)  $k=0.35$ , $t=60$ .

Table 1.

Disturbance amplitudes of various modes at different wavenumbers $k$ when sheets break up (predicted by weakly nonlinear analysis); $\mathit{We}=5$ , ${\it\rho}=0.1$ , ${\it\eta}_{0}=0.2$ .

The shapes of sheets at two wavenumbers $k=0.25$ and $k=0.35$ are given in figure 11. The most prominent characteristic is that an increase of wavenumber results in a reduction of wave amplitude. The reason is explained in the weakly nonlinear analysis: the amplitude of the linear sinuous mode ‘ $1,s$ ’ becomes significantly smaller at $k=0.35$ (see table 1). This can be easily understood from the dispersion relation shown in figure 9: the increase of $k$ from 0.25 to 0.35 leads to a significant reduction of sinuous growth rate $\text{Im}({\it\omega}_{1s})$ but a slight increase of varicose growth rate $\text{Im}({\it\omega}_{1v})$ . For the amplitudes of various modes given in table 1, while the linear varicose mode ‘ $1,v$ ’ is enhanced significantly due to the increase of $\text{Im}({\it\omega}_{1v})$ , modes ‘ $21,ss$ ’ and ‘ $21,vv$ ’ do not change very much, as the rise of $|\text{}^{2r,ss}\hat{{\it\eta}}|$ and the reduction of $|\text{}^{2r,vv}\hat{{\it\eta}}|$ with the increase of $k$ (see figure 10) offset the corresponding change in temporal growth rates. Therefore, the breakup mechanism does not vary a lot with the increase of wavenumber.

The cutoff wavenumbers are generally in proportion to ${\it\rho}\mathit{We}$ . In the following analysis, we only consider the case of $k=0.5{\it\rho}\mathit{We}$ . At this wavenumber, the instability behaviour is conspicuous and the weakly nonlinear analysis agrees well with numerical simulation when the sheet breaks up.

3.2.2. Effect of initial disturbance amplitude

The shapes of the sheets at two different initial disturbance amplitudes, ${\it\eta}_{0}=0.2$ and ${\it\eta}_{0}=0.01$ , are given in figure 12. It is easy to understand that the sheet with a smaller initial amplitude requires a longer time for the disturbance to grow to one that is large enough to cause breakup, so the breakup time for ${\it\eta}_{0}=0.01$ is significantly longer than that for ${\it\eta}_{0}=0.2$ . Apart from this obvious difference, the instability behaviour of ${\it\eta}_{0}=0.01$ exhibits many other different properties from that of ${\it\eta}_{0}=0.2$ which was studied in detail in § 3.1: it has a larger wave amplitude, its wave’s upgrade and downgrade have nearly the same slopes, and it has two locations of attenuation in one wavelength.

Figure 12.

Shapes of sheets at different initial disturbance amplitudes: (a,b) the results of numerical simulation and (c,d) the results of weakly nonlinear analysis; $\mathit{We}=5$ , ${\it\rho}=0.1$ , $k=0.25$ ; (a,c)  ${\it\eta}_{0}=0.2$ , $t=57.5$ , (b)  ${\it\eta}_{0}=0.01$ , $t=119$ , and (d)  ${\it\eta}_{0}=0.01$ , $t=114$ .

Figure 13.

Effect of initial amplitude on the temporal evolution of the amplitudes of various modes: $\mathit{We}=5$ , ${\it\rho}=0.1$ , $k=0.25$ ; (a)  ${\it\eta}_{0}=0.2$ and (b)  ${\it\eta}_{0}=0.01$ .

The reason for these differences can be explained by the weakly nonlinear analysis. As figure 13 shows, with ${\it\eta}_{0}$ reducing from 0.2 to 0.01, at the later stage of disturbance evolution, while the amplitudes of modes ‘ $1,v$ ’, ‘ $21,sv$ ’ and ‘ $21,vv$ ’ decrease significantly, the amplitude of ‘ $21,ss$ ’ increases extensively. The large amplitude of mode ‘ $1,s$ ’ is not covered by figure 13, but at the figure’s upper limit, ${\it\eta}_{a}=2$ , the slope of ${\it\eta}_{0}=0.01$ is much steeper, indicating that it grows much faster, and so, when the sheet breaks up, the amplitude of mode ‘ $1,s$ ’ is 8.019 for ${\it\eta}_{0}=0.01$ , much larger than 4.192 for ${\it\eta}_{0}=0.2$ . The increase of the ‘ $1,s$ ’ mode’s amplitude results in a larger wave amplitude; the decrease of the ‘ $21,sv$ ’ mode’s amplitude reduces its function in modulating the waves, making the slopes of upgrade and downgrade nearly identical. As for the attenuation, the domination of the ‘ $21,ss$ ’ mode makes the thinning interspaced by one half-wavelength (see figure 12 d), the one-wavelength mode ‘ $1,v$ ’ significantly reducing its amplitude and weakening its effect of eliminating one of the attenuations.

Figure 14.

Amplification factor of modes ‘ $21,ss$ ’ and ‘ $21,sv$ ’:  $\mathit{We}=5$ , ${\it\rho}=0.1$ , $k=0.25$ .

We give a further explanation of the dominance of modes ‘ $21,ss$ ’ and ‘ $1,s$ ’ at small initial amplitudes. The main reason for this phenomenon is that different modes grow at different speeds. Figure 14 gives the amplification factors of modes ‘ $21,ss$ ’ and ‘ $21,sv$ ’, denoted as $\text{}^{21,ss}f$ and $\text{}^{21,sv}f$ :

(3.8a,b ) $$\begin{eqnarray}\text{}^{21,ss}f=\exp [2\,\text{Im}({\it\omega}_{1s})t],\quad \text{}^{21,sv}f=\exp [\text{Im}({\it\omega}_{1s})t+\text{Im}({\it\omega}_{1v})t].\end{eqnarray}$$

At $t=57.5$ , the breakup time for ${\it\eta}_{0}=0.2$ , $\text{}^{21,ss}f$ and $\text{}^{21,sv}f$ have parallel values, but at a much longer breakup time, $t=114$ , rendered by a much smaller initial amplitude ${\it\eta}_{0}=0.01$ , $\text{}^{21,ss}f$ becomes far larger than $\text{}^{21,sv}f$ , the reason being that its dominance in exponential growth rate becomes extensive at longer times. This explains why the reduction of ${\it\eta}_{0}$ leads to the domination of mode ‘ $21,ss$ ’ over modes ‘ $21,sv$ ’ and ‘ $21,vv$ ’. For the first-order modes ‘ $1,s$ ’ and ‘ $1,v$ ’, this analysis is not applicable, as their dependences on ${\it\eta}_{0}$ are different: the first-order modes are proportional to ${\it\eta}_{0}$ , while the second-order modes are proportional to ${\it\eta}_{0}^{2}$ . The amplitudes of modes ‘ $1,s$ ’ and ‘ $21,ss$ ’, denoted as $\text{}^{1,s}A$ and $\text{}^{21,ss}A$ , have a simple relation: according to their expressions

(3.9a,b ) $$\begin{eqnarray}\text{}^{1,s}A=2{\it\eta}_{0}|\text{}^{1,s}\hat{{\it\eta}}_{j}|\exp [\text{Im}({\it\omega}_{1s})t],\quad \text{}^{21,ss}A=2{\it\eta}_{0}^{2}|\text{}^{21,ss}\hat{{\it\eta}}_{j}|\exp [2\,\text{Im}({\it\omega}_{1s})t],\end{eqnarray}$$

using the definition of $\text{}^{2r,ss}\hat{{\it\eta}}$ given in (2.44), we have

(3.10) $$\begin{eqnarray}\text{}^{1,s}A=\sqrt{\frac{2(\text{}^{21,ss}A)}{|\text{}^{2r,ss}\hat{{\it\eta}}|}}.\end{eqnarray}$$

Therefore, the increase of the amplitude of mode ‘ $21,ss$ ’ results in a corresponding amplification of mode ‘ $1,s$ ’, explaining the large amplitude of the linear sinuous mode at a small initial disturbance amplitude. For the linear varicose mode ‘ $1,v$ ’, such a simple relation does not exist. We can only explain a little vaguely that, because its growth rate $\text{Im}({\it\omega}_{1v})$ is small, the increment of its amplification factor $\text{}^{1,v}f=\exp [\text{Im}({\it\omega}_{1v})t]$ as a result of a longer breakup time is not enough to compensate the effect of the reduction in ${\it\eta}_{0}$ , leading to the decrease of its amplitude.

3.2.3. Effect of Weber number and gas-to-liquid density ratio

The effects of Weber number and gas-to-liquid density ratio have been extensively studied by previous researchers. It is now well known that an increase in Weber number or gas-to-liquid density ratio significantly increases temporal growth rates, thereby prominently enhancing instability, and that the instability is mainly caused by the interaction between liquid and gas. In the present study, while the typical discussion, i.e. the effect of $\mathit{We}$ and ${\it\rho}$ on the temporal growth rates, second-order amplitudes and breakup time, is briefly noted, we also reveal some interesting properties of unstable waves, and give explicit explanations for them through the weakly nonlinear analysis.

Figure 15.

Shapes of sheets at different Weber numbers: (a,b) the results of numerical simulation and (c,d) the results of weakly nonlinear analysis; ${\it\rho}=0.1$ , $k=0.5{\it\rho}\mathit{We}$ , ${\it\eta}_{0}=0.2$ ; (a,c)  $\mathit{We}=3$ , $t=124.8$ , (b)  $\mathit{We}=6$ , $t=46$ , and (d)  $\mathit{We}=6$ , $t=45$ .

Table 2.

Effect of Weber number on temporal growth rates and relative second-order disturbance amplitudes; ${\it\rho}=0.1$ , $k=0.5{\it\rho}\mathit{We}$ .

Table 3.

Effect of gas-to-liquid density ratio on temporal growth rates and relative second-order disturbance amplitudes; $\mathit{We}=3$ , $k=0.5{\it\rho}\mathit{We}$ .

Figure 16.

Shapes of sheets at different gas-to-liquid density ratios: (a,b) the results of numerical simulation and (c,d) the results of weakly nonlinear analysis; $\mathit{We}=3$ , $k=0.5{\it\rho}\mathit{We}$ , ${\it\eta}_{0}=0.2$ ; (a,c)  ${\it\rho}=0.1$ , $t=124.8$ , (b)  ${\it\rho}=0.13$ , $t=87$ , and (d)  ${\it\rho}=0.13$ , $t=83.5$ .

The effects of Weber number and gas-to-liquid density ratio on the major parameters, i.e. temporal growth rates and relative second-order disturbance amplitudes, are given in tables 2 and 3. An increase in Weber number and gas-to-liquid density ratio results in an extensive increase in the temporal growth rates, and this leads to a much shorter breakup time (see figures 15 and 16). For the second-order amplitudes, the amplitude of mode ‘ $21,ss$ ’ increases significantly, while the other two modes only change slightly. But this does not bring about the predominance of ‘ $21,ss$ ’ when the sheet breaks up, because the effect of this increase in second-order initial amplitude is somewhat marginal compared to the prominent effect of growth rates. A further examination shows that such a change in Weber number and gas-to-liquid density ratio does not result in a qualitative change in the breakup mechanism.

The numerical simulation shows that, as the Weber number and gas-to-liquid density ratio increase, the protuberance near the point of breakup becomes larger. This indicates that a more extensive interaction between gas and liquid can significantly enhance the local nonlinearity and distort the sheet more severely. The protuberances on the sheets have already been observed in the experiments of Asare et al. (Reference Asare, Takahashi and Hoffman1981) and Tammisola et al. (Reference Tammisola, Sasaki, Lundell, Matsubara and Söderberg2011); their experiments also showed that these protuberances are very likely to break up into smaller drops. But the mechanism of their formation and evolution is not very clear so far, and we look forward to future studies on this issue.

Next, we use the weakly nonlinear analysis to explain two interesting phenomena in figures 15 and 16. The first phenomenon is that, in figure 15, at $\mathit{We}=3$ , the sheet breaks up on its wave crest, but at $\mathit{We}=6$ , the sheet breaks up near its wave trough. The main reason is the shift of the linear varicose mode ‘ $1,v$ ’, as the discussion of figure 7 has indicated that the sheet breaks up near the contraction of mode ‘ $1,v$ ’ (see figure 15 c,d). At $\mathit{We}=3$ , the contraction of mode ‘ $1,v$ ’ is near the wave crest, but at $\mathit{We}=6$ , its contraction moves rightwards, resulting in the corresponding change in breakup location.

The second phenomenon is that, in figure 16, using $y=0$ as the origin, at ${\it\rho}=0.1$ , the depth of the wave trough is larger than the height of the wave crest, but at ${\it\rho}=0.13$ , the situation is the opposite, the height of the wave crest being larger than the depth of the wave trough. This is mainly caused by the first harmonic generated by the coupling between the linear sinuous mode and the varicose mode, i.e. mode ‘ $21,sv$ ’. At ${\it\rho}=0.1$ , the troughs of mode ‘ $21,sv$ ’ overlap both the peak and trough of the basic sinuous wave, shifting both of them downwards; while at ${\it\rho}=0.13$ , it is the crests of mode ‘ $21,sv$ ’ that are in the positions of the peak and trough of the basic sinuous wave, and thus shift both of them upwards. This phenomenon has again highlighted the importance of the nonlinear coupling between the linear sinuous and varicose modes.