Slip-enhanced Rayleigh-Plateau instability of a liquid film on a fibre

Surface-tension-driven instability of liquid jets and their subsequent disintegration into droplets has been much investigated since the pioneering work of Plateau (1873) and Rayleigh (1878, 1892). This instability also plays an important role in the dynamics of a liquid film coated on a fibre with additional complexities introduced by liquid-solid interfaces, which has received considerable scientific interest (Quéré, 1999) because it is crucial to a range of technologies such as optical fibre manufacturing (Deng et al., 2011), droplet transport (Lee et al., 2022) and water collection through fog harvesting (Chen et al., 2018; Zhang et al., 2022).

The Rayleigh-Plateau (RP) instability analysis for annular films was first performed by Goren (1962), who showed that a film of outer radius $h_{0}$ is unstable to sufficiently long wavelength disturbances $\lambda>\lambda_{\mathrm{crit}}=2\pi h_{0}$ ($\lambda_{\mathrm{crit}}$ is the critical wavelength of the instability beyond which no growth of the instability occurs) and dominant (the most unstable/fast growing) modes depend on the ratio of $h_{0}$ to the fibre radius $a$, confirmed by experiments (Goren, 1964). Hammond (1983) found the same instability of a thin film inside a capillary tube with a nonlinear analysis for the appearance and growth of periodically spaced lobes. Under the lubrication approximation that the film thickness is much smaller than the fibre radius (i.e. $h_{0}-a\ll a$), Frenkel (1992) proposed a weakly nonlinear thin-film equation. Because of its simplicity and capability of predicting nonlinear behaviours, this lubrication equation and its higher-order versions (Craster & Matar, 2006; Ruyer-Quil et al., 2008) have been widely used to study thin-film dynamics flowing down a fibre theoretically (Kalliadasis & Chang, 1994; Yu & Hinch, 2013) and explain experimental findings (Quéré, 1990; Duprat et al., 2007; Craster & Matar, 2009; Ji et al., 2019). Recently, these lubrication models have been extended for more complicated situations with other physics taken into account, such as electric fields (Ding et al., 2014), heat transfer (Zeng et al., 2017), thermal fluctuations (Zhang et al., 2021) and Van der Waals forces (Tomo et al., 2022).

Most of the previous works mentioned employed the classical no-slip boundary condition at the liquid-solid interface, which is inaccurate in some practical cases. For example, silicone oil on a solid wall made of nylon, used in Quéré’s experiments (Quéré, 1990), has been found to exhibit non-negligible fluid slippage (Brochard-Wyart et al., 1994; Lauga et al., 2007). This slip effect has attracted much research attention recently (Secchi et al., 2016; Zhang et al., 2020; Kavokine et al., 2021, 2022) and been known to be crucial to dynamics of a variety of interfacial flows (Liao et al., 2014; Halpern et al., 2015; Martínez-Calvo et al., 2020; Zhao et al., 2022). For a cylindrical film, Ding & Liu (2011) proposed a lubrication equation incorporating slip to explore the instability of a film flowing down a porous vertical fibre. They found that the instability is enhanced by the fluid-porous interface, modelled as a slip boundary condition in the lubrication model. A two-equation lubrication model was then developed for the case at non-infinitesimal Reynold numbers with inertia taken into account (Ding et al., 2013). For an annular film inside a slippery tube, Liao et al. (2013) solved a lubrication equation with leading order terms numerically and found that the instability is exaggerated by a fractional amount of wall slip, resulting in much more rapid draining compared to the usual no-slip case (Hammond, 1983). Haefner et al. (2015) experimentally investigated the influence of slip on the RP instability for a film on a horizontal fibre. Similar to the findings of both Ding & Liu (2011) and Liao et al. (2013), the instability was found to be enhanced by the wall slip with larger growth rates of undulations obtained. Additionally, these experimental results were shown to match predictions of a slip-modified lubrication equation. Halpern & Wei (2017) later demonstrated that wall slip can enhance the drop formation in a film flowing down a vertical fibre, providing a plausible interpretation for the discrepancy between the experimentally-obtained and theoretically-predicted critical Bond numbers for drop formation. More recently, the slippage hydrodynamics of cylindrical films were further explored in a non-isothermal environment (Chao et al., 2018) and under the influence of intermolecular forces (Ji et al., 2019) via more complicated lubrication models.

Noticeably, Kliakhandler et al. (2001) experimentally demonstrated that the lubrication models for no-slip cases are only valid when the film thickness is smaller than the fibre radius, confirmed by following theoretical analysis performed by Craster & Matar (2006) based on the Stokes equations. However, when slip is considered, the accuracy of lubrication models for films on a fibre has never been examined carefully, despite its wide usages in theoretical studies (Ding & Liu, 2011; Liao et al., 2013; Halpern & Wei, 2017; Chao et al., 2018) and experimental investigations (Haefner et al., 2015; Ji et al., 2019). In fact, slip-modified lubrication models predict a dominant wavelength independent of slip length, which somehow is not consistent with experimental data (Haefner et al., 2015). Thus, it is necessary to provide a more general theoretical analysis to go beyond the lubrication paradigm and clarify the validities of slip-modified lubrication models. For example, it is unclear whether slip-modified lubrication models work for large-slip cases.

In this work, linear stability analysis of the axisymmetric Stokes equations is performed to investigate slip-enhanced RP instability of liquid films on a fibre. Direct numerical simulations of the Navier-Stokes (NS) equations are employed to validate the theoretical findings and provide more physical insights. We also compare our theory with available experimental results (Haefner et al., 2015). The article is laid out as follows. Non-dimensionalised governing equations for films on fibres are introduced in § 2. In § 3, a theoretical model incorporating slip is developed for the RP instability. Numerical simulations, performed in § 4, are compared with the predictions of the theoretical model for the influence of slip on the growth rate (§ 4.1) and the dominant wavelength (§ 4.2) of perturbations. In §5, the theoretically-predicted dominant wavelengths are compared with the experimental data of Haefner et al. (2015).

We consider a Newtonian liquid film coating on a horizontal fibre of the radius $a$ with the $x$-axis along the centre line (figure 1). The initial radius of the film measured from the $x$-axis is $r=h_{0}$.

The incompressible NS equations are employed to predict the dynamics of the flow inside the liquid film. To identify the governing dimensionless parameters, we non-dimensionalise the NS equations with the rescaling variables shown below:

where $\tilde{h},\,\,\tilde{t},\,\,\tilde{\mathbf{u}},\,\,\textrm{and}\,\,\tilde{p}$ represent the dimensional interface height, time, velocity, and pressure respectively (note that the dimensional material parameters are not given tildes). $\tilde{\bm{\tau}}$ is the shear stress, which is proportional to the strain rate in Newtonian fluids. $\mu$ is the liquid dynamic viscosity and $\gamma$ is surface tension of the liquid-gas interface. The dimensionless NS equations can be written as:

The non-dimensional quantity $\mathrm{Oh}=\mu/\sqrt{\rho\gamma h_{0}}$ is the Ohnesorge number, which relates the viscous forces to inertial and surface-tension forces. Note that, in the previous experiment for the slip-enhanced RP instability (Haefner et al., 2015), a high viscosity polymer liquid, entangled polystyrene, is coated on a fibre with radius $10\textrm{ - }25\,\mu\mathrm{m}$. The density of the entangled polystyrene $\rho=1.05\,\mathrm{g\,cm^{-3}}$ and its liquid-gas surface tension $\gamma=30.8\,\mathrm{mN\,m^{-1}}$. The liquid viscosity depends on temperature with the various values in the range $0.2\textrm{ - }10\,\mathrm{kg\,m^{-1}s^{-1}}$ (Haefner, 2015). So the Ohnesorge number of liquid film in the experiment is larger than 10 and $\mathrm{Oh}^{2}\gg 1$ in equation (3). It is reasonable to assume that the inertial terms are negligible compared to the viscous terms. With axisymmetric initial perturbations, equation (2) and (3) can be reduced to axisymmetric Stokes equations, written as:

where $u$ and $v$ represent the axial and radial velocity respectively.

Since the density of gas around the film is much smaller than that of liquid, the gas flow outside can be assumed to be dynamically passive to simplify the problem. The liquid-gas interface height $h(x,t)$ satisfies the kinematic boundary condition

The normal stress balances at the interface $r=h$ gives

where $\bm{n}$ is the outward normal and $\nabla\cdot\bm{n}$ represents the dimensionless Laplace pressure. The tangential force balance is

where $\bm{t}$ is the tangential vector. With $\bm{n}$ and $\bm{t}$ expressed in terms of the unit vectors in $x$-direction ($\hat{\bm{e}}_{x}$) and $r$-direction ($\hat{\bm{e}}_{r}$),

(8) and (9) explicitly give

for the normal forces, and

for the tangential forces. Here $\partial_{x}$ and $\partial_{x}^{2}$ refers to first and second spatial derivatives.

In terms of the boundary conditions at the fibre surface $r=\alpha$ ($\alpha$ is the dimensionless fibre radius, i.e. $\alpha=a/h_{0}$), we have the slip boundary condition in tangential direction and the no-penetration boundary condition in the normal direction such that

Here, $l_{\mathrm{s}}$ represents the dimensionless slip length rescaled by $h_{0}$.

The governing equations above can be further simplified under the long-wave approximation to give a one-dimensional lubrication equation, compared with experimental results in the previous work (Haefner et al., 2015). This lubrication equation is displayed here for the instability analysis appearing in § 3 (see Appendix A for the derivation). Its dimensionless format is

where $M(h)$ is the the mobility term

and $p$ is the Laplace pressure

In this section, linear instability analysis based on equations (4 - 12) is performed using the normal mode method, which has been widely used for the RP instability in different fluid configurations (Rayleigh, 1878; Tomotika, 1935; Craster & Matar, 2006; Liang et al., 2011). The dimensionless perturbed quantities are expressed as

where $\omega$ is the growth rate of perturbations and $k$ is the wavenumber. Substituting these quantities into the mass and momentum equations (4 - 6) leads to

With equations (16 - 18), we can eliminate $\hat{u}$ and $\hat{p}$ to get a fourth order ordinary partial differential equation for $\hat{v}$

The general solution for (19) can be obtained in terms of Bessel functions, written as

where $I$ and $K$ are modified Bessel function of first and second kind, The subscripts ‘$0$’ and ‘$1$’ represent order of the Bessel functions. $C_{1}\textrm{ - }C_{4}$ are four arbitrary constants to be determined by the boundary conditions and $r\in[\alpha,1]$. Substituting (20) into (16) and (17) gives us the solutions for $\hat{u}$ and $\hat{p}$

In a similar approach, the dimensionless perturbed quantities for $\hat{u},\hat{v},\hat{p}$, combined with $h(x,t)=1+\hat{h}e^{\omega t+ikx}$, are substituted into perturbed boundary equations (7 - 12). For the boundary conditions at the interface ($r=1$), their linearisation gives

And for the boundary conditions on the fibre surface ($r=\alpha$), their linearised forms are

According to equation (25), $\hat{h}$ in (24) can be eliminated to give the final four equations of the boundary conditions, i.e. (23), (24), (26) and (27). Substituting the Bessel functions (20 - 22) into these perturbed equations leads to a homogeneous system of linear equations for $C_{1}\textrm{ - }C_{4}$, which has a non-trivial solution only if the determinant of the coefficients vanishes. In this way, we have the final equation

where

Because $\omega$ only appears linearly in the first line of the determinant, the dispersion relation between $\omega$ and $k$ can be expressed explicitly

where

To compare the Stokes model (29) with the lubrication theory, equation (13) is linearised using $h=1+\hat{h}e^{\omega t+ikx}$, providing the dispersion relation

where $M=\left[3+\alpha^{4}-4\alpha^{2}+4\,\mathrm{ln}\,\alpha-4l_{\mathrm{s}}(1-\alpha^{2})^{2}/\alpha\right]/16.$

In figure 2, predictions of (29) and (46) in the limiting case ($l_{\mathrm{s}}=0$) are found to agree with Craster & Matar (2006)’s results. Note that the lubrication dispersion relation (dashed lines) matches the Stokes dispersion relation (solid lines) when $\alpha\geq 0.6$. However, the agreement deteriorates for the smaller $\alpha$, mainly due to the limitations of lubrication approximations for ‘thicker’ films. Despite the failure in predicting the dispersion relation for small $\alpha$, the Plateau instability criterion with critical wavenumber, $k_{\mathrm{crit}}=2\,\pi h_{0}/\lambda_{\mathrm{crit}}=1$ (Plateau, 1873) can be predicted by the lubrication model for all $\alpha$ values. The reason is that both the circumferential and tangential curvature for the Laplace pressure is kept in (15) to give the term $k^{2}-1$ in the dispersion relation.

Figure 3 illustrates the dispersion relations of different $l_{\mathrm{s}}$ values. Both the Stokes model and the lubrication model show that the RP instability is enhanced by slip with faster growing perturbations for larger $l_{\mathrm{s}}$, consistent with previous findings based on lubrication models (Liao et al., 2013; Haefner et al., 2015; Halpern & Wei, 2017). However, obvious discrepancies between the Stokes dispersion relation and the lubrication dispersion relation are found in most cases in figure 3. Even for thin films ($\alpha\geq 0.7$), where the lubrication model has been shown to perform well on the no-slip boundary condition (see the agreement of black lines in figure 3 (a,b)), remarkable discrepancies appear when $l_{\mathrm{s}}\geq 0.5$, indicating that the lubrication model (13) is not suitable for large-slip cases. The main reason is that the lubrication approximation imposes a restriction on $l_{\mathrm{s}}$, namely $l_{\mathrm{s}}\ll\lambda^{2}-1+\alpha$, where $\lambda$ is the perturbation wavelength (see the derivation of (73) and its corresponding explanation in Appendix A). Moreover, the discrepancy increases for larger $\alpha$, similar to the trend found in the no-slip cases in figure 2.

The $k_{\mathrm{crit}}$ values for the RP instability in figure 3 are found not to be affected by the slip. These values are determined by the competition between the two curvature terms for the Laplace pressure, shown on the right-hand side of (10), where the term of the circumferential curvature, $1/\left(h\sqrt{1+\left(\partial_{x}h\right)^{2}}\right)$ is the driving force and the term of tangential curvature, $\partial_{x}^{2}h/\left(1+\left(\partial_{x}h\right)^{2}\right)^{3/2}$ is the resisting force. The balance of the two forces gives the $k_{\mathrm{crit}}$, which is independent of $l_{\mathrm{s}}$. The dominant wavenumber $k_{\mathrm{max}}$ predicted by the Stokes model in figure 3 is shown to be modified by slip, while $k_{\mathrm{max}}$ from the lubrication model is unchanged with an analytical expression $k_{\mathrm{max}}=\sqrt{1/2}$, derived from (46). The variations of $k_{\mathrm{max}}$ with slip are further illustrated in figure 4 (a), where $k_{\mathrm{max}}$ declines with increasing $l_{\mathrm{s}}$, possibly leading to larger drop sizes after the film breakup. This trend holds for different fibre radii with more rapid decrease of $k_{\mathrm{max}}$ of thicker films (smaller $\alpha$). These findings indicate the limitations of the classical slip lubrication models (Liao et al., 2013; Haefner et al., 2015; Halpern & Wei, 2017; Chao et al., 2018; Ji et al., 2019), which can be severe with large slip, while the Stokes model proposed here has the capability of predicting the instability over a much larger range of $l_{\mathrm{s}}$. Additionally, the lubrication model shows that corresponding growth rates $\omega_{\mathrm{max}}$ of each $k_{\mathrm{max}}$ grow linearly with increasing $l_{\mathrm{s}}$, which can be described by $\omega_{\mathrm{max}}=-M/4$, while the Stokes model presents nonlinear growth with smaller $\omega_{\mathrm{max}}$.

In this section, direct numerical simulations are performed to support the theoretical findings in § 3 and provide more physical insights into the slip-enhanced RP instability of the cylindrical films on fibres. The open-source framework Basilisk (Popinet, 2014, 2018) is employed to solve axisymmetric incompressible Navier-Stokes equations with surface tension. The interface between the high-density liquid and the low-density ambient air is reconstructed by a volume-of-fluid method, which has been validated by different kinds of interfacial flows such as jet breakup (Deblais et al., 2018), bubble bursting (Berny et al., 2020) and breaking waves (Mostert & Deike, 2020).

The simulation domain is a rectangle (the section of a hollow cylinder in cylindrical coordinates) with a size $[\alpha,2]\times[0,L]$. Here $L$ is the length of the film/fibre and $\alpha$ is the radius of the fibre. The liquid film of initial radius $h_{0}=1$ is placed at the bottom of the domain with small perturbations at the liquid-gas interface. The left and right boundary conditions of the simulation domain are considered periodic. The top is a symmetry boundary and bottom is a slip-wall boundary, as described by equation (12). All simulations are performed in the dimensionless unit with the rescaling variables displayed in equation (1). The non-dimensional parameter $\mathrm{Oh}=10$, which comes from liquid properties displayed in the previous experiment (Haefner, 2015) with $h_{0}=10\,\mu\mathrm{m},\,\rho=1.05\,\mathrm{g\,cm^{-3}},\,\gamma=30.8\,\mathrm{mN\,m^{-1}}\,\,\mathrm{and}\,\,\mu=0.22\,\mathrm{kg\,m^{-1}s^{-1}}$. Two configurations with different film lengths and initial perturbations are simulated here to investigate the influence of slip on the growth rates of perturbations (§ 4.1) and the dominant wavelength (§ 4.2) respectively.

Besides the numerical validations, the Stokes model is further compared with the experimental results in the previous work done by Haefner et al. (2015). The experiments were done by measuring the evolutions of entangled polystyrene films with homogeneous initial thicknesses on no-slip and slip fibres of radius $a=9.6\,\mu$m.

Figure 12 displays the dimensional dominant wavelengths $\tilde{\lambda}_{\mathrm{max}}$ ($\lambda_{\mathrm{max}}h_{0}$) of liquid films with different initial thickness $h_{0}$, where the blue triangles and red dots represent experimental data of the films on no-slip and slip fibres respectively. Since $\lambda_{\mathrm{max}}$ predicted by the lubrication model (46) is the constant ($\lambda_{\mathrm{max}}=2\sqrt{2}\pi$), independent of both $l_{\mathrm{s}}$ and $\alpha$, $\lambda_{\mathrm{max}}h_{0}$ grows linearly with increasing $h_{0}$ (see the black dashed lines in figure 12). However, the lubrication model seemingly underestimates experimentally-obtained $\lambda_{\mathrm{max}}$. The possible reason is that most $\alpha$ values investigated in the experiment are smaller than $0.5$, outside the range of validity of the lubrication approximation. A series of dominant wavelengths with different $\alpha$ values are extracted from the Stokes model (29) to generate a blue solid line for no-slip cases and a red dash-dotted line for slip cases with a fitting $\tilde{l}_{\mathrm{s}}$ in figure 12. Note that the dimensionless slip length varies along the red dash-dotted curve with different $h_{0}$ values. Though both the blue and red curves still look linear due to small changes of $\lambda_{\mathrm{max}}$, obvious differences between the Stokes predictions and the lubrication ones are found for the thick films ($h_{0}>40\,\mu$m). Promisingly, despite dispersion of experimental data, the Stokes model can provide better predictions for dominant wavelengths of both no-slip cases and slip cases compared to the lubrication model.

In addition to the dominant wavelengths, we then make a further comparison between the experimental growth rates and theoretical ones, illustrated in figure 13, where the growth rate is non-dimensionalised by the ratio of the capillary velocity to the fibre radius, i.e., $\omega=\frac{a}{\gamma/\mu}\tilde{\omega}$. However, using the parameters from Haefner et al. (2015) ($\gamma/\mu=294\,\mu\mathrm{m\,min}^{-1}$ and $\tilde{l}_{\mathrm{s}}=0.3\,a$), the Stokes model is found to provide worse agreements with the experimental data than the lubrication model. Noticeably, predictions of the lubrication model match the experimental data even in the cases of thick films ($h_{0}/a>3$), which is unreasonable due to the thin-film assumption in the lubrication theory. The reason for the disagreement is currently unclear and should be the subject of future investigation.

In this work, to study the influence of liquid-sold slip on the RP instability of liquid films on a fibre, a theoretical model based on the linear stability analysis of axisymmetric Stokes equations is developed, which goes beyond the classical slip lubrication model. The lubrication model is found to overestimate the slip enhanced growth rates of perturbations and fail to predict the dominant wavelength, while the new model can provide more accurate predictions for the instability. Direct numerical simulations of the NS equations are then performed to validate the two theoretical findings via two fluid configurations: (i) a short film with a fixed wavelength for the slip-enhanced growth rate of the perturbation, and (ii) a long film with random initial perturbations for the slip-dependent dominant modes of perturbations. The dominant wavelengths are also confirmed by experimental data from the previous work (Haefner et al., 2015). Overall, the applicability of the slip lubrication model should be restricted to thin-film and small-slip flows. For more general problems, it is better to use the Stokes theory proposed here.

Promisingly, the wall slip has been realised experimentally for films on a fibre (Haefner et al., 2015; Ji et al., 2019) and the slip length can also be measured directly (Huang et al., 2006; Maali & Bhushan, 2012; Maali et al., 2016). It is our hope that the wall slip can be controlled experimentally to verify our new predictions, most notably the relation between drop sizes ($\lambda_{\mathrm{max}}$) and slip lengths. Potential extensions of the framework are numerous. For example, to include the influence of other physics such as electric fields (Ding et al., 2014) and intermolecular forces (Ji et al., 2019; Tomo et al., 2022), both of which are know to affect $k_{\mathrm{crit}}$ and $k_{\mathrm{max}}$ in the RP instability, and to explore related flow configurations, like a film flowing down a fibre (Liu & Ding, 2021), where an interesting open problem is to understand the influence of slip on dynamics of travelling waves in different flow regimes.

[Acknowledgements] Useful discussions with Dr. Mykyta Chubynsky, Dr. Mu Kai, Mr. Qiao Ran and Dr. Qin Jian are gratefully acknowledged.

[Funding]This work was supported by the National Natural Science Foundation of China (grant number $11621202$) and the Youth Innovation Promotion Association CAS (grant number $2018491$).

[Declaration of interests] The authors report no conflict of interest.

[Author ORCIDs]
Chengxi Zhao https://orcid.org/0000-0002-3041-0882
Yixin Zhang https://orcid.org/0000-0003-4632-3780
Ting Si https://orcid.org/0000-0001-9071-8646