2 Experiments

2.1 Methods

Figure 1 shows a schematic of the set-up we used to experimentally study the characteristics of a liquid film flowing down a vertical string. We use a programmable syringe pump to introduce a working liquid into the nozzle and generate flows. The high-speed camera is mounted on a translation stage for focusing and positioning and is operated at a frame rate of 1000 frames per second.

Figure 1.

Experimental set-up.

The working fluid is a well-wetting liquid of low surface energy, Rhodorsil silicone oil v50, with density $\unicode[STIX]{x1D70C}=963~\text{kg}~\text{m}^{-3}$ , kinematic viscosity $\unicode[STIX]{x1D708}=50~\text{mm}^{2}~\text{s}^{-1}$ , surface tension $\unicode[STIX]{x1D70E}=20.8~\text{mN}~\text{m}^{-1}$ at $20\,^{\circ }\text{C}$ , and capillary length $l_{c}=1.5~\text{mm}$ . We use stainless-steel nozzles with the nozzle inner diameter ( $ID$ ) ranging from 0.6 mm to 2.5 mm and the wall thickness ranging from 0.1 mm to 0.2 mm. The experiments are performed using 0.6 m long nylon strings of diameter $0.2$ and 0.43 mm, both smaller than the capillary length. A weight is attached to vertically align the string. Two $X$ – $Y$ stages are used to centre the string with respect to the nozzle. The liquid mass flow rate, monitored using a weight scale connected to a computer, is varied from $0.02~\text{g}~\text{s}^{-1}$ to $0.09~\text{g}~\text{s}^{-1}$ .

2.2 Data analysis

We use image processing tools to extract parameters such as the average bead speed, period $L^{\ast }$ , maximum film thickness $h_{m}^{\ast }$ , and mass contained in each bead. In order to characterize the liquid bead profile, we use the longitudinal distance between the fibre and the maximum curvature point along the liquid. We utilize the colour contrast on the image to extract the contour of the fluid and the fibre to thereby determine the liquid film thickness $h^{\ast }$ . A local least-squares smoothing is performed to account for pixelation noise. From this profile, we determine the average period $L^{\ast }$ as the distance between two adjacent maxima. Integrating the profile over the chosen domain, we find the mass constraint value $M_{0}$ as defined later in § 3. Figure 2 shows a representative experimental frame and the film profile produced after the data are processed.

Figure 2.

Experimental frame (a) and processed data analysis (b) for fibre radius $R^{\ast }=0.1~\text{mm}$ , nozzle inner diameter $ID=0.8~\text{mm}$ , and flow rate $Q_{m}^{\ast }=0.04~\text{g}~\text{s}^{-1}$ , shown at a distance 8–10 mm from the nozzle. The red dots superimposed on the experimental frame correspond to the extracted film profile and the black dots on the bottom plot correspond to the locations of the maxima. The average distance between two maxima is $L^{\ast }=5.92~\text{mm}$ and the film thickness between the drops ranges from 0.095 to 0.120 mm.

The uncertainties in the parameters obtained from our image processing are estimated in terms of a single pixel scaling and the selection of the colour value for the contour. The estimated uncertainty is $\pm 0.08~\text{mm}$ for the streamwise length and $\pm 0.10~\text{mm}$ for the bead height. Using a total of 1000 images for each run, we measure the averages for maximum fluid height, mass, bead speed and distance between two subsequent beads.

3 Model formulation

We consider a flow of two-dimensional axisymmetric Newtonian fluid down a vertical cylinder of radius $R^{\ast }$ (see figure 3). The liquid properties, including the surface tension $\unicode[STIX]{x1D70E}$ , density $\unicode[STIX]{x1D70C}$ and kinematic viscosity $\unicode[STIX]{x1D708}$ , are all assumed to be constant. This model formulation does not include the nozzle size as a system parameter because we are interested in the flow downstream where the nozzle does not affect the dynamics. Instead, we consider the scales for frequency and mass of the droplets, which depend on the nozzle size (Sadeghpour et al. Reference Sadeghpour, Zeng and Ju2017). We will discuss this dependence later in appendix B. We review below the derivation of the governing equations and boundary conditions by Ruyer-Quil et al. (Reference Ruyer-Quil, Treveleyan, Giorgiutti-Dauphiné, Duprat and Kalliadasis2008) and Craster & Matar (Reference Craster and Matar2006) and discuss our inclusion of additional physics related to slip, curvature and wetting properties.

Figure 3.

Schematic of a thin liquid film flowing down a vertical cylindrical fibre.

The dimensional Navier–Stokes equations for axisymmetric flows are

(3.1a ) $$\begin{eqnarray}\displaystyle & \displaystyle u_{t^{\ast }}^{\ast }+u^{\ast }u_{x^{\ast }}^{\ast }+v^{\ast }u_{y^{\ast }}^{\ast }=-\frac{1}{\unicode[STIX]{x1D70C}}(p_{x^{\ast }}^{\ast }+\unicode[STIX]{x1D6F1}_{x^{\ast }}^{\ast })+g+\unicode[STIX]{x1D708}\left(u_{x^{\ast }x^{\ast }}^{\ast }+\frac{u_{y^{\ast }}^{\ast }}{y^{\ast }}+u_{y^{\ast }y^{\ast }}^{\ast }\right), & \displaystyle\end{eqnarray}$$

(3.1b ) $$\begin{eqnarray}\displaystyle & \displaystyle v_{t^{\ast }}^{\ast }+v^{\ast }v_{y^{\ast }}^{\ast }+u^{\ast }v_{x^{\ast }}^{\ast }=-\frac{p_{y^{\ast }}^{\ast }}{\unicode[STIX]{x1D70C}}+\unicode[STIX]{x1D708}\left(\frac{v_{y^{\ast }}^{\ast }}{y^{\ast }}+v_{y^{\ast }y^{\ast }}^{\ast }-\frac{v^{\ast }}{y^{\ast 2}}+v_{x^{\ast }x^{\ast }}^{\ast }\right), & \displaystyle\end{eqnarray}$$

where $t^{\ast }$ represents the time, $u^{\ast }$ and $v^{\ast }$ represent the axial and radial components of the velocity, $p^{\ast }$ is the pressure and $g$ is the gravitational acceleration. We adopt the functional form of the disjoining pressure used in Reisfeld & Bankoff (Reference Reisfeld and Bankoff1992) and introduce the film stabilization term $\unicode[STIX]{x1D6F1}^{\ast }(h^{\ast })$ ,

(3.1c ) $$\begin{eqnarray}\unicode[STIX]{x1D6F1}^{\ast }(h^{\ast })=-\frac{A^{\ast }}{{h^{\ast }}^{3}},\end{eqnarray}$$

where $A^{\ast }$ is a stabilization parameter. The equation of continuity is given by

(3.1d ) $$\begin{eqnarray}v_{y^{\ast }}^{\ast }+\frac{v^{\ast }}{y^{\ast }}+u_{x^{\ast }}^{\ast }=0.\end{eqnarray}$$

Along the fibre, at the interface between the solid substrate and the fluid $y^{\ast }=R^{\ast }$ , we impose the Navier slip and no-penetration boundary conditions:

(3.1e ) $$\begin{eqnarray}v^{\ast }=0,\quad u^{\ast }=\unicode[STIX]{x1D706}^{\ast }u_{y^{\ast }}^{\ast }\quad \text{at }y^{\ast }=R^{\ast },\end{eqnarray}$$

where $\unicode[STIX]{x1D706}^{\ast }>0$ is the slip length in standard slip models (Münch, Wagner & Witelski Reference Münch, Wagner and Witelski2005; Haefner et al. Reference Haefner, Benzaquen, Bäumchen, Salez, Peters, McGraw, Jacobs, Raphaël and Dalnoki-Veress2015). The no-slip boundary condition corresponds to $\unicode[STIX]{x1D706}^{\ast }=0$ . Typical slip lengths for polymeric liquids such as silicone oil range from $1$ to $10~\unicode[STIX]{x03BC}\text{m}$ (Quéré Reference Quéré1990; Halpern & Wei Reference Halpern and Wei2017). The normal and shear stress balances on the free surface $y^{\ast }=R^{\ast }+h^{\ast }$ are given by

(3.1f ) $$\begin{eqnarray}\displaystyle & \displaystyle p^{\ast }=\frac{2\unicode[STIX]{x1D707}}{1+h_{x^{\ast }}^{\ast 2}}(h_{x^{\ast }}^{\ast 2}u_{x^{\ast }}^{\ast }-h_{x^{\ast }}^{\ast }(v_{x^{\ast }}^{\ast }+u_{y^{\ast }}^{\ast })+v_{y^{\ast }}^{\ast })+\frac{\unicode[STIX]{x1D70E}}{(1+h_{x^{\ast }}^{\ast 2})^{3/2}}\left(\frac{1+h_{x^{\ast }}^{\ast 2}}{R^{\ast }+h^{\ast }}-h_{x^{\ast }x^{\ast }}^{\ast }\right),\quad & \displaystyle\end{eqnarray}$$

(3.1g ) $$\begin{eqnarray}\displaystyle & \displaystyle (1-h_{x^{\ast }}^{\ast 2})(v_{x^{\ast }}^{\ast }+u_{y^{\ast }}^{\ast })+2h_{x^{\ast }}^{\ast }(v_{y^{\ast }}^{\ast }-u_{x^{\ast }}^{\ast })=0, & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D707}$ is the dynamic viscosity, and $\unicode[STIX]{x1D70E}$ scales the total curvature which consists of a destabilizing azimuthal curvature term and a stabilizing axial curvature term. We then complete the system by including the kinematic boundary condition on the free surface $y^{\ast }=R^{\ast }+h^{\ast }$ ,

(3.1h ) $$\begin{eqnarray}h_{t^{\ast }}^{\ast }+u^{\ast }h_{x^{\ast }}^{\ast }=v^{\ast }\quad \text{at }y^{\ast }=R^{\ast }+h^{\ast }.\end{eqnarray}$$

The next step is to choose the appropriate dimensionless parameters in order to non-dimensionalize the model. Following Duprat et al. (Reference Duprat, Ruyer-Quil and Giorgiutti-Dauphiné2009) we choose the scales for the system as follows: the length scale in the radial direction $y$ is ${\mathcal{H}}$ , and the length scale in the streamwise direction $x$ is ${\mathcal{L}}={\mathcal{H}}/\unicode[STIX]{x1D716}$ . The scale ratio $\unicode[STIX]{x1D716}$ is set by the balance between the surface tension term $\unicode[STIX]{x1D70E}h_{x^{\ast }x^{\ast }x^{\ast }}^{\ast }$ , arising from $p_{x^{\ast }}^{\ast }$ , and the gravity $g$ , and is given by $\unicode[STIX]{x1D716}=(\unicode[STIX]{x1D70C}g{\mathcal{H}}^{2}/\unicode[STIX]{x1D70E})^{1/3}$ . This scale ratio is small (around $0.4$ ) in typical experiments and can also be rewritten as $\unicode[STIX]{x1D716}=We^{-1/3}$ , where the Weber number $We=(l_{c}/{\mathcal{H}})^{2}$ compares the capillary length $l_{c}=\sqrt{\unicode[STIX]{x1D70E}/(\unicode[STIX]{x1D70C}g)}$ to the radial length scale ${\mathcal{H}}$ . Then the characteristic streamwise velocity is ${\mathcal{U}}=(g{\mathcal{H}}^{2})/\unicode[STIX]{x1D708}$ , and the pressure and time scales are given by $\unicode[STIX]{x1D70C}g{\mathcal{L}}$ and $(\unicode[STIX]{x1D708}{\mathcal{L}})/(g{\mathcal{H}}^{2})$ , respectively. With these scales, we drop the star superscript in (3.1) and write the non-dimensional Navier–Stokes equations using dimensionless variables in the form

(3.2a ) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D716}^{2}Re(u_{t}+uu_{x}+vu_{y})=-p_{x}-\unicode[STIX]{x1D6F1}_{x}+1+\frac{u_{y}}{y}+u_{yy}+\unicode[STIX]{x1D716}^{2}u_{xx}, & \displaystyle\end{eqnarray}$$

(3.2b ) $$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D716}^{4}Re(v_{t}+vv_{y}+uv_{x})=-p_{y}+\unicode[STIX]{x1D716}^{2}\left(\frac{v_{y}}{y}+v_{yy}-\frac{v}{y^{2}}\right)+\unicode[STIX]{x1D716}^{4}v_{xx}, & \displaystyle\end{eqnarray}$$

where the Reynolds number $Re={\mathcal{U}}{\mathcal{L}}/\unicode[STIX]{x1D708}$ . The dimensionless continuity equation is identical to (3.1d ) in form. The balances of normal and tangential stresses at $y=h+R$ are expressed as

(3.2c ) $$\begin{eqnarray}(1-\unicode[STIX]{x1D716}^{2}h_{x}^{2})(\unicode[STIX]{x1D716}^{2}v_{x}+u_{y})+2\unicode[STIX]{x1D716}^{2}h_{x}(v_{y}-u_{x})=0,\end{eqnarray}$$

(3.2d ) $$\begin{eqnarray}\displaystyle p & = & \displaystyle \frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D716}^{2}(1+\unicode[STIX]{x1D6FC}h)(1+\unicode[STIX]{x1D716}^{2}h_{x}^{2})^{1/2}}-\frac{h_{xx}}{(1+\unicode[STIX]{x1D716}^{2}h_{x}^{2})^{3/2}}\nonumber\\ \displaystyle & & \displaystyle +\,\frac{\unicode[STIX]{x1D716}^{4}}{1+\unicode[STIX]{x1D716}^{2}h_{x}^{2}}[h_{x}^{2}u_{x}-h_{x}v_{x}+\unicode[STIX]{x1D716}^{2}(-h_{x}u_{y}+v_{y})],\end{eqnarray}$$

where the dimensionless parameter $\unicode[STIX]{x1D6FC}={\mathcal{H}}/R^{\ast }$ is the aspect ratio of the characteristic radial length scale and the fibre radius, and the dimensionless fibre radius $R=R^{\ast }/{\mathcal{H}}$ . The slip and no-penetration boundary conditions on $y=R$ are

(3.2e ) $$\begin{eqnarray}v=0,\quad u=\unicode[STIX]{x1D706}u_{y},\end{eqnarray}$$

where the non-dimensional slip length $\unicode[STIX]{x1D706}=\unicode[STIX]{x1D706}^{\ast }/{\mathcal{H}}$ . The kinematic boundary condition remains unaltered in its form.

We next simplify the above set of governing equations by following Ruyer-Quil et al. (Reference Ruyer-Quil, Treveleyan, Giorgiutti-Dauphiné, Duprat and Kalliadasis2008). Under the lubrication approximation, the inertial contributions can be neglected since $Re=O(1)$ and $\unicode[STIX]{x1D716}\ll 1$ . Note that we do not assume the small ratio $\unicode[STIX]{x1D6FC}$ between the liquid film thickness and the fibre radius. Omitting the terms of order $O(\unicode[STIX]{x1D716}^{2})$ , we rewrite the leading-order non-dimensional reduced momentum and continuity equations for the velocity field $(u,v)$ and the dynamic pressure $p$ as

(3.3a ) $$\begin{eqnarray}\displaystyle & \displaystyle 1-\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}x}-\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F1}}{\unicode[STIX]{x2202}x}+\frac{\unicode[STIX]{x2202}^{2}u}{\unicode[STIX]{x2202}y^{2}}+\frac{u_{y}}{y}=0, & \displaystyle\end{eqnarray}$$

(3.3b ) $$\begin{eqnarray}\displaystyle & \displaystyle -\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}y}=0, & \displaystyle\end{eqnarray}$$

(3.3c ) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}u}{\unicode[STIX]{x2202}x}+\frac{\unicode[STIX]{x2202}v}{\unicode[STIX]{x2202}y}+\frac{v}{y}=0. & \displaystyle\end{eqnarray}$$

The balance of tangential stresses at the free surface $y=h+R$ is reduced to

(3.3d ) $$\begin{eqnarray}u_{y}=0.\end{eqnarray}$$

The destabilizing azimuthal curvature and the stabilizing streamwise curvature terms in (3.2d ) are both important throughout the liquid film. For $\unicode[STIX]{x1D716}\ll 1$ a formal expansion of the azimuthal curvature term in (3.2d ) yields

(3.3e ) $$\begin{eqnarray}\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D716}^{2}(1+\unicode[STIX]{x1D6FC}h)\sqrt{1+\unicode[STIX]{x1D716}^{2}h_{x}^{2}}}=\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D716}^{2}(1+\unicode[STIX]{x1D6FC}h)}\left(1-\frac{1}{2}\unicode[STIX]{x1D716}^{2}h_{x}^{2}+O(\unicode[STIX]{x1D716}^{4})\right),\end{eqnarray}$$

which shows that this term is of order $O(\unicode[STIX]{x1D6FC}/\unicode[STIX]{x1D716}^{2})$ , for both cases $\unicode[STIX]{x1D6FC}h\ll 1$ and $\unicode[STIX]{x1D6FC}h=O(1)$ . At the steep front of a moving droplet, we have $h_{x}=O(1)$ , and the second term in the expansion contributes an additional $O(\unicode[STIX]{x1D6FC})$ term to the dynamic pressure. Similarly, for the streamwise curvature term in (3.2d ), we have

(3.3f ) $$\begin{eqnarray}\frac{h_{xx}}{(1+\unicode[STIX]{x1D716}^{2}h_{x}^{2})^{3/2}}=h_{xx}\left(1-\frac{3}{2}\unicode[STIX]{x1D716}^{2}h_{x}^{2}+O(\unicode[STIX]{x1D716}^{4})\right).\end{eqnarray}$$

At the front of the moving droplet, the $O(\unicode[STIX]{x1D716}^{2})$ term introduced by the nonlinear curvature is negligible. Therefore, we only keep a fully nonlinear azimuthal curvature term and use the linearized curvature in the streamwise direction. The balance of normal stresses at the free surface $y=h+R$ is then reduced to

(3.3g ) $$\begin{eqnarray}p=\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D716}^{2}(1+\unicode[STIX]{x1D6FC}h)\sqrt{1+\unicode[STIX]{x1D716}^{2}h_{x}^{2}}}-\frac{\unicode[STIX]{x2202}^{2}h}{\unicode[STIX]{x2202}x^{2}}.\end{eqnarray}$$

It is important to note that the first term on the right-hand side of (3.3g ) accounts for the balance between the azimuthal and axial scales characterized by $\unicode[STIX]{x1D6FC}$ and $\unicode[STIX]{x1D716}$ . We will further discuss appropriate forms of this term for different cases later.

To derive the evolution equation for $h$ from above, we first consider a uniform Nusselt flow without any interfacial instabilities. The velocity field of this flow is obtained by balancing the viscosity and gravity acceleration in (3.2). That is, the streamwise velocity $u_{0}$ of the Nusselt flow satisfies the reduced system

(3.4a ) $$\begin{eqnarray}1+\frac{\unicode[STIX]{x2202}^{2}u_{0}}{\unicode[STIX]{x2202}y^{2}}+\frac{1}{y}\frac{\unicode[STIX]{x2202}u_{0}}{\unicode[STIX]{x2202}y}=0,\end{eqnarray}$$

with boundary conditions

(3.4b ) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}u_{0}}{\unicode[STIX]{x2202}y}=0\quad \text{at }y=h+R,\quad u_{0}=\unicode[STIX]{x1D706}\frac{\unicode[STIX]{x2202}u_{0}}{\unicode[STIX]{x2202}y}\quad \text{at }y=R.\end{eqnarray}$$

Solving (3.4) for $u_{0}$ leads to

(3.5) $$\begin{eqnarray}u_{0}(y)=-\frac{1}{4}(y^{2}-R^{2})+\frac{1}{2}(h+R)^{2}\ln \left(\frac{y}{R}\right)+h\unicode[STIX]{x1D706}\left(\frac{h}{2R}+1\right).\end{eqnarray}$$

The rescaled kinematic boundary conditions together with the slip and no-penetration boundary conditions (3.2e ) and the shear stress condition (3.3d ) lead to the mass conservation equation

(3.6) $$\begin{eqnarray}(1+\unicode[STIX]{x1D6FC}h)\frac{\unicode[STIX]{x2202}h}{\unicode[STIX]{x2202}t}+\frac{\unicode[STIX]{x2202}q}{\unicode[STIX]{x2202}x}=0,\quad \text{where }q=\frac{1}{R}\int _{R}^{h+R}uy\,\text{d}y.\end{eqnarray}$$

Following Ruyer-Quil et al. (Reference Ruyer-Quil, Treveleyan, Giorgiutti-Dauphiné, Duprat and Kalliadasis2008) where an approach based on a projection of the velocity field of a test function is used, we take the inner product of (3.3a ) with the Nusselt uniform solution $u_{0}$ , and obtain

(3.7) $$\begin{eqnarray}\int _{R}^{h+R}\left(1+u_{yy}+\frac{u_{y}}{y}\right)u_{0}y\,\text{d}y=\int _{R}^{h+R}(p_{x}+\unicode[STIX]{x1D6F1}_{x})u_{0}y\,\text{d}y.\end{eqnarray}$$

From equation (3.3b ), we see that the pressure $p$ is a function of $x$ only, and the right-hand side of the above equation becomes $p_{x}\bar{R}q_{0}$ . Note that the Nusselt solution $u_{0}$ satisfies (3.4). By integrating by parts and applying the slip and no-penetration boundary conditions (3.2e ) and the shear stress condition (3.3d ), we obtain

(3.8) $$\begin{eqnarray}q=\left(1-\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}x}-\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F1}}{\unicode[STIX]{x2202}x}\right)q_{0}.\end{eqnarray}$$

Here $q_{0}$ , the flow rate per unit circumference length for the uniform Nusselt layer, is

(3.9) $$\begin{eqnarray}q_{0}=\frac{1}{R}\int _{R}^{R+h}u_{0}y\,\text{d}y=\frac{h^{3}}{3}\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC}h)+\frac{h^{2}}{4}(\unicode[STIX]{x1D6FC}h+2)^{2}\unicode[STIX]{x1D706},\end{eqnarray}$$

where the shape factor $\unicode[STIX]{x1D719}$ is a function defined by

(3.10) $$\begin{eqnarray}\unicode[STIX]{x1D719}(X)=\frac{3}{16X^{3}}[(1+X)^{4}(4\log (1+X)-3)+4(1+X)^{2}-1].\end{eqnarray}$$

Given a dimensional volumetric flow rate $Q_{m}^{\ast }$ and fibre radius $R^{\ast }$ , we define the volumetric flow rate per circumference unit $q_{0}^{\ast }$ as $q_{0}^{\ast }=Q_{m}^{\ast }/(2\unicode[STIX]{x03C0}\unicode[STIX]{x1D70C}R^{\ast })$ . From (3.9) the characteristic axial length scale ${\mathcal{H}}$ for a uniform Nusselt flow can then be obtained. For the no-slip case $\unicode[STIX]{x1D706}=0$ , equation (3.9) gives the Nusselt solution $h\equiv {\mathcal{H}}$ used in Ruyer-Quil et al. (Reference Ruyer-Quil, Treveleyan, Giorgiutti-Dauphiné, Duprat and Kalliadasis2008) and Duprat et al. (Reference Duprat, Ruyer-Quil and Giorgiutti-Dauphiné2009).

Finally, by rescaling the time scale $t\rightarrow t/\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})$ and combining (3.6), (3.9) and (3.8), we obtain the dimensionless governing equation

(3.11) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}t}\left(h+\frac{\unicode[STIX]{x1D6FC}}{2}h^{2}\right)+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}\left[{\mathcal{M}}(h)\left(1-\frac{\unicode[STIX]{x2202}p}{\unicode[STIX]{x2202}x}-\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F1}}{\unicode[STIX]{x2202}x}\right)\right]=0,\end{eqnarray}$$

where the mobility function takes the form

(3.12) $$\begin{eqnarray}{\mathcal{M}}(h;\unicode[STIX]{x1D706},\unicode[STIX]{x1D6FC})=\frac{h^{3}}{3}\frac{\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC}h)}{\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})}+\frac{h^{2}(\unicode[STIX]{x1D6FC}h+2)^{2}\unicode[STIX]{x1D706}}{4\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})}.\end{eqnarray}$$

This choice of time scale leads to a normalized mobility function such that ${\mathcal{M}}=1/3$ for $h=1$ and $\unicode[STIX]{x1D706}=0$ . Moreover, in the limit of the aspect ratio $\unicode[STIX]{x1D6FC}\rightarrow 0$ , the mobility function becomes

(3.13) $$\begin{eqnarray}{\mathcal{M}}(h;\unicode[STIX]{x1D706},\unicode[STIX]{x1D6FC})=\frac{1}{3}h^{2}(h+3\unicode[STIX]{x1D706})+\frac{\unicode[STIX]{x1D6FC}}{3}h^{2}(h+3\unicode[STIX]{x1D706})(h-1)+O(\unicode[STIX]{x1D6FC}^{2}),\end{eqnarray}$$

and its leading-order term agrees with the mobility function used in the slip model proposed by Halpern & Wei (Reference Halpern and Wei2017). Substituting (3.3g ) into (3.11) yields an evolution equation for the film thickness $h(x,t)$ :

(3.14)

where ${\mathcal{Z}}(h)$ includes the destabilizing azimuthal curvature of the film and the film stabilization term. Comparing (3.14) and (3.11), and using (3.3g ) with a scaling parameter $\unicode[STIX]{x1D702}=\unicode[STIX]{x1D716}^{2}$ , we write the complete form of ${\mathcal{Z}}(h)$ as

(3.15) $$\begin{eqnarray}{\mathcal{Z}}(h)=p+h_{xx}+\unicode[STIX]{x1D6F1}=\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D702}(1+\unicode[STIX]{x1D6FC}h)\sqrt{1+\unicode[STIX]{x1D702}h_{x}^{2}}}+\unicode[STIX]{x1D6F1}(h).\end{eqnarray}$$

Equation (3.14) is a fourth-order nonlinear partial differential equation (PDE) for the thickness $h(x,t)$ . This model accounts for the surface tension, gravity and azimuthal instabilities, but neglects inertia and streamwise viscous dissipation (Ruyer-Quil et al. Reference Ruyer-Quil, Trevelyan, Giorgiutti-Dauphiné, Duprat and Kalliadasis2009). When the Rayleigh–Plateau instability dominates over the inertia and streamwise viscous dissipation, our model (3.14) is expected to provide a good agreement with experimental data.

Different forms of ${\mathcal{Z}}(h)$ appear in the literature to represent the azimuthal curvature of the film. For instance, Yu & Hinch (Reference Yu and Hinch2013) assumed that the film is much thinner than the fibre radius, ${\mathcal{H}}\ll R^{\ast }$ (which corresponds to $\unicode[STIX]{x1D6FC}\ll 1$ ), and used a simple form ${\mathcal{Z}}(h)=-(\unicode[STIX]{x1D6FC}/\unicode[STIX]{x1D716})^{2}h$ where the leading-order constant term in the expansion of the azimuthal curvature is neglected. For thicker films that are of interest to the present study, $\unicode[STIX]{x1D6FC}=O(1)$ and $\unicode[STIX]{x1D6FC}h$ is at most $O(1)$ so we use (3.16) as the first form of ${\mathcal{Z}}$ ,

(3.16) $$\begin{eqnarray}{\mathcal{Z}}_{CM}(h)=\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D702}(1+\unicode[STIX]{x1D6FC}h)},\end{eqnarray}$$

following Craster & Matar (Reference Craster and Matar2006), Sisoev et al. (Reference Sisoev, Craster, Matar and Gerasimov2006) and Ruyer-Quil et al. (Reference Ruyer-Quil, Treveleyan, Giorgiutti-Dauphiné, Duprat and Kalliadasis2008, Reference Ruyer-Quil, Trevelyan, Giorgiutti-Dauphiné, Duprat and Kalliadasis2009). With ${\mathcal{Z}}_{{\mathcal{C}}{\mathcal{M}}}(h)$ in (3.16) and $\unicode[STIX]{x1D706}=0$ , the evolution equation (3.14) reduces to

(3.17) $$\begin{eqnarray}(\text{CM})\qquad \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}t}\left(h+\frac{\unicode[STIX]{x1D6FC}}{2}h^{2}\right)+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}\left[\frac{h^{3}}{3}\frac{\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC}h)}{\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})}\left(1-\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}x}\left[\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D702}(1+\unicode[STIX]{x1D6FC}h)}-h_{xx}\right]\right)\right]=0,\end{eqnarray}$$

which is consistent with the evolution equation derived by Craster & Matar (Reference Craster and Matar2006) except for a scaling difference. For the rest of this paper, we will refer to the model (3.17) as the Craster & Matar (CM) model.

To incorporate the slip effects under the framework of the CM model, we use the evolution equation (3.14) with the non-dimensional slip length $\unicode[STIX]{x1D706}>0$ , the mobility function (3.12) and the azimuthal curvature term (3.16). This will be referred to as the Slip Craster & Matar (SCM) model. We will show that this slip model promotes droplet formation and leads to an increased speed of propagation.

With the aspect ratio $\unicode[STIX]{x1D6FC}$ being of $O(1)$ , in the limit as $\unicode[STIX]{x1D702}\rightarrow 0$ , (3.16) becomes singular. To address this problem, in this paper we will consider two additional forms of ${\mathcal{Z}}(h)$ to reflect the balance between the azimuthal and axial scales under different flow conditions.

As the second form of ${\mathcal{Z}}$ , we define

(3.18) $$\begin{eqnarray}{\mathcal{Z}}_{FC}(h)=\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D702}(1+\unicode[STIX]{x1D6FC}h)\sqrt{1+\unicode[STIX]{x1D702}h_{x}^{2}}}.\end{eqnarray}$$

This fully nonlinear azimuthal curvature term provides an $O(\unicode[STIX]{x1D6FC}/\unicode[STIX]{x1D702}^{3/2})$ contribution to the dynamic pressure near the advancing edge of the moving beads. It has been shown in many applications (Snoeijer Reference Snoeijer2006; Lopes, von Borries & Hazel Reference Lopes, von Borries and Hazel2018) that using the full expression for the curvature term can provide better accuracy for lubrication models, and this issue was recently reviewed in Thiele (Reference Thiele2018). We will show that the inclusion of the fully nonlinear curvature yields an increased speed of propagation of the travelling beads and partially improves agreement with our experimental results. Compared with the thick-film model proposed in Kliakhandler et al. (Reference Kliakhandler, Davis and Bankoff2001), where both curvature terms are fully nonlinear, our model has only a nonlinear azimuthal curvature term. The stabilizing axial curvature term is linearly approximated. For the remainder of this paper, we will refer to the evolution equation (3.14) with the mobility function (3.12) and the azimuthal curvature term (3.18) as the Full curvature model (FCM).

As the third form of ${\mathcal{Z}}$ , we include the film stabilization term to deal with the unbalanced azimuthal curvature term. Thus we replace (3.16) with

(3.19) $$\begin{eqnarray}{\mathcal{Z}}_{FS}(h)=\frac{\unicode[STIX]{x1D6FC}}{\unicode[STIX]{x1D702}(1+\unicode[STIX]{x1D6FC}h)}-\frac{A}{h^{3}}.\end{eqnarray}$$

The last term of (3.19) is motivated by the functional form of the long-range attractive part of the well-known apolar van der Waals model (de Gennes Reference de Gennes1985; Oron, Davis & Bankoff Reference Oron, Davis and Bankoff1997; Oron & Bankoff Reference Oron and Bankoff2001) for the well-wetting liquids used in our experiments. The stabilization parameter $A=(\unicode[STIX]{x1D6FC}^{2}\unicode[STIX]{x1D716}_{p}^{4})/(3\unicode[STIX]{x1D702}(\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D716}_{p}+1)^{2})>0$ is expressed in terms of a dimensionless thickness $\unicode[STIX]{x1D716}_{p}$ , below which a thin uniform fluid layer on a fibre is stable. We will discuss this stability criterion further in § 4. The film stabilization term is found to significantly improve agreement with our experimental data as discussed in § 5. The model (3.14) with ${\mathcal{M}}(h)$ from (3.12) and the azimuthal curvature with film stabilization term (3.19) will be referred to as the film stabilization model (FSM).

Table 1.

A summary of the models using (3.14).

In summary, this paper considers four versions of the model (3.14) listed in table 1 to incorporate different physical effects. When contrasted with the CM model, which corresponds to the asymptotic model studied in Craster & Matar (Reference Craster and Matar2006), the other three new models take into consideration the slip effects, the nonlinear curvature and the film stabilization term, respectively.

4 Film stabilization model and stability analysis

In this section, we derive the stabilization parameter $A$ in (3.19) and study its effects. Using a stability analysis study, we show that $A>0$ is important to reproduce the stable travelling wave solutions (TWS) observed in experiments.

We begin by examining the linear stability of the uniform Nusselt solution. Following the approach in Craster & Matar (Reference Craster and Matar2006), we perturb the uniform film ( $h_{0}=1$ ) by an infinitesimal Fourier mode,

(4.1) $$\begin{eqnarray}h=h_{0}+\unicode[STIX]{x1D6FF}\text{e}^{\text{i}kx+\unicode[STIX]{x1D6EC}t},\end{eqnarray}$$

where $k$ is the wavenumber, $\unicode[STIX]{x1D6EC}$ describes the growth rate of the perturbation and $\unicode[STIX]{x1D6FF}$ $(\ll 1)$ is the initial amplitude. Expanding PDE (3.14) then gives the dispersion relation:

(4.2) $$\begin{eqnarray}\unicode[STIX]{x1D6EC}=-\text{i}kc_{k}(\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D706})+\frac{k^{2}}{3(1+\unicode[STIX]{x1D6FC})}\left(\frac{\unicode[STIX]{x1D6FC}^{2}}{(1+\unicode[STIX]{x1D6FC})^{2}\unicode[STIX]{x1D702}}-k^{2}-3A\right)\left(1+\frac{3(\unicode[STIX]{x1D6FC}+2)^{2}\unicode[STIX]{x1D706}}{4\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})}\right),\end{eqnarray}$$

where $c_{k}$ is the speed of linear kinematic wave solutions of (3.14) for small wavenumbers,

(4.3) $$\begin{eqnarray}c_{k}(\unicode[STIX]{x1D6FC},\unicode[STIX]{x1D706})=\frac{1}{\unicode[STIX]{x1D6FC}+1}\left(1+\frac{\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D719}^{\prime }(\unicode[STIX]{x1D6FC})}{3\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})}\right)+\frac{(\unicode[STIX]{x1D6FC}+2)\unicode[STIX]{x1D706}}{\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})},\end{eqnarray}$$

and it increases linearly with the presence of slip. The real part of $\unicode[STIX]{x1D6EC}$ is the effective growth rate, and the thickness $h_{0}=1$ is long-wave unstable with respect to perturbations with $0<k<k_{c}$ , where the critical wavenumber $k_{c}=[\unicode[STIX]{x1D6FC}^{2}/((1+\unicode[STIX]{x1D6FC})^{2}\unicode[STIX]{x1D702})-3A]^{1/2}$ . For $A=0$ this cut-off wavenumber corresponds to the classical Rayleigh–Plateau (RP) mode for the capillary instability of a viscous jet, where $Re(\unicode[STIX]{x1D6EC})=0$ at $k=k_{c}$ . The most unstable mode occurs at $k_{m}=[\unicode[STIX]{x1D6FC}^{2}/(2(1+\unicode[STIX]{x1D6FC})^{2}\unicode[STIX]{x1D702})-(3A/2)]^{1/2}$ , where the largest growth rate $\unicode[STIX]{x1D6EC}_{m}$ is attained. Equation (4.2) also shows that the slip effects always enhance the RP instability, which agrees with the work by Halpern & Wei (Reference Halpern and Wei2017).

Figure 4.

Dispersion relation plots of $Re(\unicode[STIX]{x1D6EC}_{0})$ against the wavenumber $k$ with the constant film thickness (a) $h_{0}=1$ and (b) $h_{0}=0.2$ for the CM model, the SCM with slip $\unicode[STIX]{x1D706}=0.1$ , and the FSM with $\unicode[STIX]{x1D716}_{p}=0.2$ , which corresponds to $A=0.027$ from (4.5). Other system parameters are $\unicode[STIX]{x1D6FC}=4.96$ and $\unicode[STIX]{x1D702}=0.23$ .

We repeat the above analysis for a uniform thin layer $h=h_{0}~({<}1)$ . The effective growth rate of disturbances with wavenumber $k$ is given by

(4.4) $$\begin{eqnarray}Re(\unicode[STIX]{x1D6EC}_{0})=\frac{k^{2}}{3h_{0}(\unicode[STIX]{x1D6FC}h_{0}+1)}\left(\frac{h_{0}^{4}\unicode[STIX]{x1D6FC}^{2}}{(\unicode[STIX]{x1D6FC}h_{0}+1)^{2}\unicode[STIX]{x1D702}}-k^{2}h_{0}^{4}-3A\right)\left(\frac{\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC}h_{0})}{\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})}+\frac{3(\unicode[STIX]{x1D6FC}h_{0}+2)^{2}}{4h_{0}\unicode[STIX]{x1D719}(\unicode[STIX]{x1D6FC})}\unicode[STIX]{x1D706}\right).\end{eqnarray}$$

Figure 4 gives the dispersion relation (4.4) for the uniform Nusselt solution $h_{0}=1$ and a thinner undisturbed layer of $h_{0}=0.2$ under the CM, SCM and FSM models. The Nusselt solution $h_{0}=1$ is unstable for small wavenumbers for all three models. The finite slippage increases the growth rates of the unstable modes. For the CM and SCM models, the thinner undisturbed liquid layer ( $h_{0}=0.2$ ) is unstable over larger ranges of the wavenumber. In stark contrast, the FSM model saturates the unstable modes, rendering the thinner undisturbed layer linearly stable for all wavenumbers (see the solid curve in figure 4 b) for certain values of $A$ .

In the literature, $A$ is typically expressed as a Hamaker constant in terms of microscopic quantities (de Gennes Reference de Gennes1985). We instead choose the values of $A$ based on stable undisturbed liquid layers. Empirical observations of the thin-film layers between droplets indicate that these undisturbed liquid layers are comparable to the corresponding fibre radius $R^{\ast }$ . Therefore, we pick a coating thickness $\unicode[STIX]{x1D716}_{p}^{\ast }\approx R^{\ast }$ , obtained from experimental measurements. Using the dimensionless undisturbed layer thickness $\unicode[STIX]{x1D716}_{p}=\unicode[STIX]{x1D716}_{p}^{\ast }/{\mathcal{H}}$ and the dispersion relation in (4.4), we derive a formula for $\tilde{A}$ ,

(4.5) $$\begin{eqnarray}\tilde{A}=\frac{\unicode[STIX]{x1D6FC}^{2}\unicode[STIX]{x1D716}_{p}^{4}}{3\unicode[STIX]{x1D702}(\unicode[STIX]{x1D6FC}\unicode[STIX]{x1D716}_{p}+1)^{2}}.\end{eqnarray}$$

For $A=\tilde{A}$ , any thin flat film of thickness less than the threshold value $\unicode[STIX]{x1D716}_{p}$ is linearly stable, i.e. $Re(\unicode[STIX]{x1D6EC}_{0})<0$ , for all wavenumbers.

We next examine the influences of the film stabilization terms (3.19) on the profiles and speeds of propagation of the travelling wave patterns governed by (3.14). We consider the model over a periodic domain $0\leqslant x\leqslant L$ and introduce a change of variables to the reference frame of the travelling wave,

(4.6a-c ) $$\begin{eqnarray}\unicode[STIX]{x1D709}=x-ct,\quad s=t,\quad h(x,t)=\tilde{h}(\unicode[STIX]{x1D709},s),\end{eqnarray}$$

where $c$ is the speed of the travelling wave. Then $\tilde{h}(\unicode[STIX]{x1D709},s)$ satisfies the PDE

(4.7) $$\begin{eqnarray}\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}s}\left(\tilde{h}+\frac{\unicode[STIX]{x1D6FC}}{2}\tilde{h}^{2}\right)-c\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\unicode[STIX]{x1D709}}\left(\tilde{h}+\frac{\unicode[STIX]{x1D6FC}}{2}\tilde{h}^{2}\right)+\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\unicode[STIX]{x1D709}}\left[{\mathcal{M}}(\tilde{h})\left(1-\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}\unicode[STIX]{x1D709}}[{\mathcal{Z}}(\tilde{h})-\tilde{h}_{\unicode[STIX]{x1D709}\unicode[STIX]{x1D709}}]\right)\right]=0.\end{eqnarray}$$

The travelling wave solution $H(\unicode[STIX]{x1D709})$ is a steady state of the PDE (4.7) and satisfies the fourth-order ordinary differential equation (ODE)

(4.8) $$\begin{eqnarray}c\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}\left(H+\frac{\unicode[STIX]{x1D6FC}}{2}H^{2}\right)=\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}\left[{\mathcal{M}}(H)\left(1-\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}\left[{\mathcal{Z}}(H)-H_{\unicode[STIX]{x1D709}\unicode[STIX]{x1D709}}\right]\right)\right].\end{eqnarray}$$

This is a nonlinear eigenvalue problem, where the speed of propagation $c$ corresponds to the eigenvalue. The effects of slippage $\unicode[STIX]{x1D706}$ and full nonlinear curvature (3.18) will be studied in § 5.3. Newton’s method is used to solve the nonlinear ODE (4.8) where the speed $c$ is treated as an unknown variable. In order to achieve local uniqueness, we impose a constraint of mass conservation,

(4.9) $$\begin{eqnarray}\int _{0}^{L}H+\frac{\unicode[STIX]{x1D6FC}}{2}H^{2}\,\text{d}\unicode[STIX]{x1D709}=M_{0},\end{eqnarray}$$

for given values of mass $M_{0}$ and wavelength $L$ extracted from our experimentally obtained liquid film profiles. We also set $H(\unicode[STIX]{x1D709}_{0})=H_{0}$ for $0\leqslant \unicode[STIX]{x1D709}_{0}\leqslant L$ .

Numerical investigations reveal that including the film stabilization term in (3.19) enhances the moving speed of travelling wave solutions. In figure 5(a), we plot the profiles of travelling wave solutions with a varying stabilization parameter $A$ and identical mass constraint $M_{0}$ and period $L$ . Since a larger value of $A$ corresponds to a stronger wetting potential in the lubrication model, the near-flat coating thickness increases with increasing $A$ , as expected. Correspondingly, the height of the beads decreases as the coating thickness increases under the fixed overall mass constraint. That is, the presence of the film stabilization term saturates capillary instabilities and generates smaller moving beads. Figure 5(b) shows that the predicted velocity $c$ increases with the parameter $A$ superlinearly for small $A$ , and increases linearly as $A$ becomes larger. In § 5.1, we will show that, compared to the CM model (3.17), the FSM model with parameter $A$ given by (4.5) significantly improves predictions of the travelling wave velocity against experimental observations.

Figure 5.

Travelling wave solutions of (4.8) with a varying stabilization parameter $A=0,0.01,0.02,0.03$ showing that (a) stronger regularization yields smaller droplet heights and thicker precursor layers and (b) the moving speed of drops is significantly increased with larger values of $A$ .

Earlier works have not reported a detailed stability analysis of the travelling wave. Therefore, we take a closer look at the linear stability of the travelling wave solutions in the film stabilization model. We consider a positive periodic travelling wave solution $H(\unicode[STIX]{x1D709})$ over the domain $0\leqslant \unicode[STIX]{x1D709}\leqslant L$ , and perturb it by setting $\tilde{h}(\unicode[STIX]{x1D709},s)=H(\unicode[STIX]{x1D709})+\unicode[STIX]{x1D6FF}\unicode[STIX]{x1D6F9}(\unicode[STIX]{x1D709})\text{e}^{\unicode[STIX]{x1D6EC}s}$ , where $\unicode[STIX]{x1D6FF}\ll 1$ and $\unicode[STIX]{x1D6F9}(\unicode[STIX]{x1D709})$ is also $L$ -periodic. Here we only focus on perturbations of the same period since the dynamics in both RP and isolated droplet regimes have fixed spatial periods. More complicated droplet dynamics (Kalliadasis & Chang Reference Kalliadasis and Chang1994) in the convective regime is not investigated here. We linearize equation (4.7) around the steady state $H(\unicode[STIX]{x1D709})$ and obtain the $O(\unicode[STIX]{x1D6FF})$ equation

(4.10) $$\begin{eqnarray}\unicode[STIX]{x1D6EC}\unicode[STIX]{x1D6F9}=\mathscr{L}\unicode[STIX]{x1D6F9},\end{eqnarray}$$

where the linear operator $\mathscr{L}$ is

(4.11) $$\begin{eqnarray}\displaystyle \mathscr{L}\unicode[STIX]{x1D6F9} & \equiv & \displaystyle c\left(\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}+\frac{\unicode[STIX]{x1D6FC}H^{\prime }}{1+\unicode[STIX]{x1D6FC}H}\right)\unicode[STIX]{x1D6F9}+\frac{1}{1+\unicode[STIX]{x1D6FC}H}\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}\left({\mathcal{M}}(H)\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}[{\mathcal{Z}}^{\prime }(H)\unicode[STIX]{x1D6F9}-\unicode[STIX]{x1D6F9}^{\prime \prime }]\right)\nonumber\\ \displaystyle & & \displaystyle -\,\frac{1}{1+\unicode[STIX]{x1D6FC}H}\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}\left({\mathcal{M}}^{\prime }(H)\left[1-\frac{\text{d}}{\text{d}\unicode[STIX]{x1D709}}({\mathcal{Z}}(H)-H^{\prime \prime })\right]\unicode[STIX]{x1D6F9}\right).\end{eqnarray}$$

Here $\unicode[STIX]{x1D6F9}(\unicode[STIX]{x1D709})$ is a normalized eigenmode associated with an eigenvalue $\unicode[STIX]{x1D6EC}$ and $\Vert \unicode[STIX]{x1D6F9}\Vert _{2}=1$ . If there are any eigenvalues $\unicode[STIX]{x1D6EC}$ with $Re(\unicode[STIX]{x1D6EC})>0$ to the problem (4.10), then the periodic travelling wave solution $H(\unicode[STIX]{x1D709})$ is unstable. For simplicity we only focus on the case where a one-period solution fits in the domain, and numerically calculate the spectrum and corresponding eigenfunctions for the eigenproblem.

Figure 6.

(a) The dependence of dominant eigenvalues of $H(\unicode[STIX]{x1D709})$ on the stabilization parameter $A$ for the eigenproblem (4.10) where the unstable eigenvalues are plotted in dashed lines and the stable eigenvalues are in solid lines. (b) Corresponding travelling wave profile $H$ and unstable eigenmodes $\unicode[STIX]{x1D6F9}_{1}$ , $\unicode[STIX]{x1D6F9}_{2}$ , $\unicode[STIX]{x1D6F9}_{3}$ for the CM model ( $A=0$ ). The parameters $(L,M_{0})=(8.86,22.16)$ correspond to the experiment with flow rate $Q_{m}^{\ast }=0.04~\text{g}~\text{s}^{-1}$ , fibre radius $R^{\ast }=0.1~\text{mm}$ and nozzle diameter $1.06~\text{mm}$ .

Figure 6(b) shows the unstable modes predicted by the CM model ( $A$ = 0). These instabilities contradict the stable TWS from the corresponding experiment and generate small wavy patterns in the flat film connecting the moving beads. However, these instabilities can be saturated by the FSM model with an appropriate stabilization parameter $A$ . For a typical travelling wave solution $H(\unicode[STIX]{x1D709})$ to the FSM model, its dominant eigenvalues with a varying $A$ are plotted in figure 6(a). It shows that, without the film stabilization term (i.e. $A=0$ ), the TWS is unstable. In addition to the translational eigenmode $\unicode[STIX]{x1D6EC}_{T}=0$ , the dominant unstable eigenvalues are given by complex conjugate pairs $\unicode[STIX]{x1D6EC}_{1}=0.0283\pm 0.756\text{i}$ , $\unicode[STIX]{x1D6EC}_{2}=0.0264\pm 1.006\text{i}$ , and $\unicode[STIX]{x1D6EC}_{3}=0.0097\pm 0.503\text{i}$ . Increasing the parameter $A$ yields three pairs of complex conjugate eigenvalues crossing the imaginary axis, suggesting three Hopf bifurcations occur at these crossing points. At each of these bifurcation points, a branch of time-periodic solution emerges corresponding to typical dynamics in the isolated droplet regime.

For $A>A_{c}$ all the eigenvalues of $H(\unicode[STIX]{x1D709})$ satisfy $Re(\unicode[STIX]{x1D6EC})<0$ and the travelling wave is stabilized. In § 5.2, we will show that the film stabilization term helps better arrest the stability transition from the RP regime to isolated droplet regime in experiments with different nozzle diameters.