Plateau–Rayleigh instability of a soft layer coated on a rigid cylinder

Abstract

We study the Plateau–Rayleigh instability of a soft viscoelastic solid layer coated on a rigid cylinder i.e., a soft fibre with a rigid core. The onset of instability is examined using a linear stability analysis. We find that increasing the rigid cylinder radius or the stiffness of the layer reduces the growth rate of the fastest growing mode. For each rigid cylinder radius, a critical elastocapillary number is found below which all wavelengths of disturbances are stable. The critical value for a soft fibre with a thick rigid cylindrical core can be several orders of magnitudes larger than that for a totally soft fibre (no rigid core). This highlights the strong stabilizing effect of the rigid core on the system. Increasing the relaxation timescale of the viscoelastic material also slows down the growth of disturbances, but has no effect on the critical elastocapillary number. Interestingly, the wavelength of the fastest growing mode is independent of the rigid cylinder radius for the purely elastic case.

---

1 Introduction

Significant stress is usually required for a solid to deform, and as such the effects due to surface tension have often been ignored. Soft solids such as elastomers and gels,¹ and biological materials such as protein condensates,² cells and tissues,³ on the other hand, have elastic moduli ranging between kPa and MPa, which means that they can deform much easier. In recent years, there has been significant attention to how capillary effects can lead to soft solid deformations, and fascinating elastocapillary phenomena have been discovered.^4–8 The Plateau–Rayleigh instability (PRI),^9,10 namely the instability driven by surface tension and often illustrated by the breakup of liquid jets into droplets,¹¹ has currently been examined for fibres made of soft solids.^12–21

Studies of the PRI of soft fibres have focused on soft materials that demonstrate elastic or viscoelastic responses.^12–21 An experimental study using soft agar gel fibres by Mora et al.¹⁴ has shown that the instability occurs when the elastocapillary length γ/μ ≥ 6R, where R is the radius of the soft fibre, γ is the solid surface tension and μ is the shear modulus. Other studies of the PRI on a soft fibre include, for example, the formation of bead-on-string structures.^12,16,18,22 Another interesting factor that might significantly modify the instability is having an inner rigid core, i.e. a fibre consists of a soft-layer-coated on a rigid cylinder. Such a kind of setup, but with a coated liquid film instead of a soft solid layer, has been studied extensively^23–32 since the early works in the 1960s by Goren.^33,34 Recent studies on this problem have addressed the effects of liquid slip on the solid surface^28,32 and the dynamics of the droplets formed on the fibre.²⁵ However, investigations on situations in which the coated layer is a soft solid are still missing. How viscoelastic properties and the rigid core influence the PRI of a soft-layer-coated fibre remain unclear, which is addressed in this article.

Soft solids have recently been used in many applications, for example in 3D bioprinting,³⁵ mimicking muscle tissues in biomedicine³⁶ and water harvesting.³⁷ In living organisms, soft fibrous-shaped structures are often found in cellular tubes and compartments of cells. Instability of these biological ingredients have been rationalised through the PRI.^38–42 One example is a recent study on the undulation and droplet formation of a layer of condensed protein TPX2 on a microtubule.⁴² Hence, investigating the PRI of fibrous soft solids is becoming more important for both understanding the fundamental physical problems and the development of new technologies. In this study, we examine the onset of PRI of a soft layer coated on a rigid cylinder through a linear stability analysis.

2 Formulation

We consider a rigid cylinder of radius r₀ coated with a layer of incompressible viscoelastic material of density ρ and uniform thickness H at an undeformed state. The whole fibre consisting of the rigid cylinder and the coated soft layer has a total radius R = r₀ + H as shown in Fig. 1. Consider that the fiber is thin enough such that the Bond number Bo ≡ ρgR²/γ ≪ 1, where g is the gravitational acceleration, which allows us to neglect the effect of gravity. By only consider viscoelastic stresses, the equation of motion of a soft material element is given bywhere σ_ij(x,t) and u_i(x,t) are respectively the stress tensor and the displacement vector in index notation at a position vector x and time t. Assuming a linear viscoelastic response of the incompressible soft material, the relation between σ_ij and the strain tensor ε_ij iswhere μ(t) is the shear relaxation function, p is the pressure, δ_ij is the Kronecker delta, and the strain tensor ε_ij is related to the displacement as

Fig. 1 Schematic representation of an infinitely long, rigid cylinder of radius r₀ coated with a layer of incompressible viscoelastic material of uniform thickness H at an undeformed state. A random disturbance of the interface of the viscoelastic layer is decomposed into sinusoidal deformations of wavenumber k. In cylindrical coordinate systems, r and z are respectively the radial coordinate and the longitudinal coordinate.

There are different models that have been used to describe viscoelastic materials. The simplest models are the Kelvin–Voigt model or the Maxwell model, which consists of a spring and a viscous damper connected in parallel or in series respectively. In this study, we consider the soft material to behave as a gel described by the Chasset–Thirion model assuming a power law response given by

where μ₀ is the static shear modulus, τ is the relaxation timescale of the viscoelastic response, Γ is the gamma function and n is a parameter typically smaller than or equal to unity.

To compute the growth rate of disturbance on the viscoelastic layer, we follow the approach delineated by previous studies.^17,32 We decompose the time-dependent variables into normal modes e^st, where s is the growth rate. The amplitude of the normal mode of a function f(t) is obtained by the Laplace transform defined as

We apply the Laplace transform to the constitutive relation (2) and obtain

where [small sigma, Greek, tilde]_ij, ũ_i and p̃ are respectively the Laplace transform of σ_ij, u_i and p, and [small mu, Greek, circumflex](s) is the shear modulus in Laplace space defined as

Thus the governing eqn (1) in Laplace space can be written as

Next, we scale the lengths with the radius of the whole fibre R, the time with the inertial-capillary time and the stresses by μ_o. We define the following dimensionless variables as

Note that variables in Laplace space have an additional dimension of time. The dimensionless form of the governing equation (eqn (8)) in Laplace space is

whereis the inverse of the dimensionless shear modulus in Laplace space which is a function of s̄ and consists of 2 parameters: Σ ≡ γ/(μ₀R) and [small tau, Greek, macron] ≡ τ/t_c. In general, the expression of β depends on the viscoelastic model. The parameter Σ is called the elastocapillary number which compares the capillary pressure to the elastic stress. When Σ → ∞, capillarity dominates over elasticity. The other parameter [small tau, Greek, macron] is the ratio of the viscoelastic relaxation timescale to the inertial-capillary time and can be interpreted as the Deborah number since the material responds purely elastically when [small tau, Greek, macron] = 0. When [small tau, Greek, macron] → ∞, viscoelastic effects are expected to dominate over inertial effects. To sum up, the dimensionless control parameters of the system are: r̄₀, Σ, [small tau, Greek, macron] and n.

We consider only the longitudinal disturbance and neglect the azimuthal disturbance as azimuthal normal modes always increase the surface energy.¹¹ The deformation of the soft layer is axisymmetric. We hence use the cylindrical coordinate system (r,ϕ,z) and the corresponding unit vectors are denoted as (r̂,[small phi, Greek, circumflex],ẑ). The displacement vector is denoted as u(r,z,t) = u_r(r,z,t)r̂ + u_z(r,z,t)ẑ. Note that u_ϕ = 0 due to axisymmetry.

The governing eqn (11) can be solved by applying the Helmholtz decomposition of the displacement as described in ref. 17 and 32. The general solutions are given as

p̄(s̄) = −s̄²[A₁I₀(k̄r̄) +A₃K₀(k̄r̄)]e^ik̄z̄, (13)

ū_r = [A₁k̄I₁(k̄r̄) − A₂(ik̄)I₁(αr̄) − A₃k̄K₁(k̄r̄) − A₄(ik̄)K₁(αr̄)]e^ik̄z̄, (14)

and

ū_z = [A₁(ik̄)I₀(k̄r̄) + A₂αI₀(αr̄) + A₃(ik̄)K₀(k̄r̄) − A₄αK₀(αr̄)]e^ik̄z̄ (15)

where k̄ ≡ kR is the dimensionless wavenumber, I and K are respectively the modified Bessel functions of the first and second kind, is a function of s̄ and k̄, and is introduced here to simplify the expressions of equations, and A_m(m = 1,2,3,4) are the undermined coefficients.

Next we impose the boundary conditions. At r̄ = r̄₀ ≡ r₀/R, there is no penetration of material, thus

ū_r|_r̄=r̄0 = 0. (16)

For the z-direction, we impose a no-slip condition

ū_z|_r̄=r̄0 = 0. (17)

At r̄ = 1, we assume the slope of the deformed interface to be small. Balancing the Laplace pressure due to solid surface tension and the viscoelastic stresses gives

in the r-direction, and gives the vanishing shear stress

[small sigma, Greek, macron]_rz|_r̄=1 = 0 (19)

in the z-direction.

3 Dispersion relation

Using the expressions of the general solutions (eqn (13)–(15)) for the boundary conditions (eqn (16)–(19)) yields a set of linear equations for the unknowns A_m which is given in Appendix 6.1. The solvability condition for non-trivial solutions of the linear equations with unknowns A_m is that the determinant is vanishing, which gives the following dispersion relationwhere
when r̄₀ → 0, the modified Bessel functions I₀(r̄₀) and I₁(r̄₀) are bounded, but K₀(r̄₀) and K₁(r̄₀) diverge. These unbounded terms are eliminated for the case of r̄₀ = 0.¹⁷ To examine the feature of the dispersion relation for small r̄₀, we expand the determinant of the dispersion relation (20) in series of r̄₀. The result for the leading order terms is given in Appendix 6.2 which consists of logarithmic terms of r̄₀.

4 Results

The validation of the dispersion relation (eqn (20)) is presented in Appendix 6.3 in which we compare our results in the Newtonian fluid limit with previous studies in the literature. In this section, we first consider the purely elastic solid limit and then we focus on the effects of viscoelasticity.

4.1 Purely elastic solid limit

When the viscoelastic relaxation timescale τ is small compared to the inertial-capillary time , i.e. taking the limit [small tau, Greek, macron] → 0, our viscoelastic model reduces to the purely elastic model. In this limit, and β = Σ in the dispersion relation (eqn (20)). We examine how the instability depends on the dimensionless parameters Σ and r̄₀.

4.1.1 Dependence on the rigid core radius and the elastocapilliary number. We show the dispersion relation (eqn (20)) by plotting the dimensionless growth rate s̄ as a function of the dimensionless wavenumber k̄ for different values of Σ and a fixed r̄₀ = 0.1 in Fig. 2(a). We see that for each curve, the unstable modes (i.e. s̄ > 0) lie within a certain range of dimensionless wavenumbers, i.e. k̄₁ < k̄ < k̄₂, where k̄₁ and k̄₂ are defined as s̄(k̄₁) = 0 and s̄(k̄₂) = 0. Within each range of unstable modes, there is a dimensionless wavenumber k̄ = k̄_m that corresponds to the fastest growing mode with a maximum dimensionless growth rate s̄_m ≡ s̄(k̄ = k̄_m). To study how the characteristic quantities depend on the control parameters, we plot s̄_m, k̄_m, and k̄₁ (and k̄₂) as a function of Σ respectively in Fig. 2(b)–(d) for three different dimensionless rigid core radii, i.e. r̄₀ = 10⁻⁹, r̄₀ = 0.1 and 0.9. We also add the result of Tamim et al.¹⁷ for the situation of a soft fibre without a rigid cylindrical core (i.e. r̄₀ = 0). We see in Fig. 2(b) that for a fixed value of Σ, s̄_m is smaller for a larger r̄₀. For all four different r̄₀, the dimensionless growth rate of the fastest growing mode s̄_m decreases with decreasing Σ. When Σ is reduced to a critical value Σ_c, s̄_m drops to zero. There is no positive solution of s̄ for Σ < Σ_c. This means that the coated elastic layer is stable under disturbance of any wavelength when Σ < Σ_c. We would also point out that, firstly, the critical value Σ_c for r̄₀ = 0.9 is orders of magnitude larger than that for r̄₀ = 0.1. This can be understood as follows: when the coated layer becomes thinner, the strain is larger and hence a softer layer can still generate strong enough elastic stress to resist the growth of disturbance. Secondly, even for the rigid cylindrical core radius as small as r̄₀ = 10⁻⁹, s̄_m starts to deviate significantly from that for the case of r̄₀ = 0 when Σ is close to Σ_c. Moreover, the critical value Σ_c for r̄₀ = 10⁻⁹ (Σ_c = 8.7) is not the same as that for r̄₀ = 0 (Σ_c = 6)^14,17 but larger. These features will be discussed in details in Section 4.1.2 where a stability phase diagram of Σ_cvs. r̄₀ will be examined.

Fig. 2 (a) The dimensionless growth rate s̄ as a function of the dimensionless wavenumber k̄ for different values of Σ and a fixed r̄₀ = 0.1 and [small tau, Greek, macron] = 0 (purely elastic). (b)–(d): The characteristic quantities s̄_m in (b), k̄_m in (c), and k̄₁ and k̄₂ in (d) as a function of Σ for four different dimensionless rigid core radius, i.e. r̄₀ = 0 (from Tamim et al.¹⁷), r̄₀ = 10⁻⁹, r̄₀ = 0.1 and 0.9. Other parameters: [small tau, Greek, macron] = 0.

Interestingly, as shown in Fig. 2(c), we find that the dimensionless wavenumber of the fastest growing mode k̄_m is independent of the dimensionless radius of the rigid core r̄₀, and decreases with reducing Σ. It is thus expected that elastic beads formed on a softer coated layer are closer together. When Σ → ∞, we find k̄_m = 0.7, which agrees with the dimensionless wavenumber of the fastest growing mode for the classical PRI of inviscid fluid threads.^10,11 Another remarkable point is that k̄_m starts to drop significantly with reducing Σ only when Σ ⪅ 10². For cases with Σ_c > 10², e.g. for r̄₀ = 0.9, the dimensionless wavenumber of the fastest growing mode is always close to the asymptotic value, i.e. k̄_m ≈ 0.7. Regarding the range of unstable modes, as we can see in Fig. 2(a) and (d), it shrinks when Σ is reduced. In contrast to the r̄₀ = 0 case for which k̄₁ remains zero, for non-zero r̄₀ cases, k̄₁ increases when Σ is reduced. This also means that the transition from unstable fibers to stable fibers (i.e. at Σ = Σ_c) occurs at k̄_m ≠ 0, which needs to be verified experimentally. It is expected that the bead separation formed on a soft fiber with a rigid core is finite (k̄_m ≠ 0) at the transition, instead of diverging for the r̄₀ = 0 case.

4.1.2. Stability phase diagram. To show how the onset of instability depends on Σ and r̄₀, we plot Σ_c (in log scales) as a function of r̄₀ in a phase diagram, see Fig. 3(a). We see that increasing the radius of the rigid cylindrical core or decreasing the elastocapillary number can make the coated film more stable. Hence both the stiffness and the rigid core have a stabilizing effect. In Fig. 3(b), we plot Σ_c as a function of the thickness of the coated layer, i.e. H̄ = 1 − r̄₀, in log–log scales. We find that when r̄₀ ⪆ 0.2, the curve follows a power law Σ_c = 21.5 H̄^λ where λ ≈ −2.75. Hence the critical value Σ_c for a thin coated elastic film is orders of magnitudes larger than a thick film. In the opposite limit of r̄₀ → 0, we find that our results can be approximately described by a logarithmic relation given aswhere Σ_co = 6 is the critical value for a soft elastic fibre without a rigid cylindrical core obtained by Mora et al.¹⁴ or Tamim et al.,¹⁷ and a = 37 is a fitting parameter. The comparison is shown in Fig. 3(c) for small r̄₀. This result is consistent with the series expansion of the dispersion relation (eqn (20)) for small r̄₀ given in Appendix 6.2.

Fig. 3 (a) Stability phase diagram. Line: the critical elastocapillary number Σ_c (log scale) as a function of r̄₀. (b) The critical elastocapillary number Σ_c as a function of 1 − r̄₀ in log–log scales. (c) A linear-log plot of the critical elastocapillary number Σ_c as a function of r̄₀.

4.2. Viscoelastic effects

We study the instability when [small tau, Greek, macron] ≠ 0. We show results for two different rigid cylindrical core radii, namely r̄₀ = 0.1 (a thick coated layer) and r̄₀ = 0.9 (a thin coated layer). Since the critical elastocapillary numbers for r̄₀ = 0.1 and r̄₀ = 0.9 are different by orders of magnitude, we take a fixed value of Σ = 50 for the case of r̄₀ = 0.1 and a fixed value of Σ = 15 000 for the case of r̄₀ = 0.9. In Fig. 4(a) and (b), we present the dispersion relation from eqn (20) for different [small tau, Greek, macron], and in Fig. 4(c) and (d) for different n. We see that both k̄₁ and k̄₂ are independent of [small tau, Greek, macron] and n. The dimensionless growth rate of the fastest growing mode s̄_m decreases when [small tau, Greek, macron] is enhanced. This means that the viscoelastic relaxation of the material slows down the growth of disturbance. We also see that when [small tau, Greek, macron] is increased, the dimensionless wavenumber of the fastest growing mode becomes smaller for r̄₀ = 0.1 but remains constant for r̄₀ = 0.9. When varying the other viscoelastic parameter n, the dimensionless growth rate of the fastest growing mode s̄_m decreases with decreasing n̄. Hence smaller values of n have a larger damping effect. The change of s̄_m with varying n̄ is more sensitive for r̄₀ = 0.9 than for r̄₀ = 0.1 as we can see in Fig. 4(c) and (d). The dimensionless wavenumber of the fastest growing mode k̄_m decreases slightly when n is increased from n = 0.2 to n = 1 for r̄₀ = 0.1 but remains constant for r̄₀ = 0.9.

Fig. 4 The dimensionless growth rate s̄ as a function of the dimensionless wavenumber k̄ for different values of [small tau, Greek, macron] and a fixed value of exponent n = 1 with r̄₀ = 0.1 and Σ = 50 in (a), and r̄₀ = 0.9 and Σ = 15 000 in (b). The dimensionless growth rate s̄ vs. the dimensionless wavenumber k̄ for different values of exponent n and a fixed value of [small tau, Greek, macron] = 1 with r̄₀ = 0.1 and Σ = 50 in (c), and r̄₀ = 0.9 and Σ = 15 000 in (d).

Next, we take [small tau, Greek, macron] = 100 and n = 0.5, and show the results for s̄_m and k̄_m as a function of Σ for three different r̄₀ = 0, 0.1 and 0.9 in Fig. 5. We see that the behaviors of s̄_m and k̄_m for [small tau, Greek, macron] = 100 are similar to those shown in Fig. 2(b) and (c) for the purely elastic case, i.e. [small tau, Greek, macron] = 0, except that when Σ is close to the critical value Σ_c. When Σ is approaching Σ_c, the curves for both s̄_m and k̄_m bend concavely, and the value of k̄_m becomes dependent on r̄₀. However, varying the values of [small tau, Greek, macron] or n does not change the critical values Σ_c, namely Σ_c remains the same as for the purely elastic case ([small tau, Greek, macron] = 0). When Σ → ∞, viscoelastic effects disappear, and the results are the same as for the PRI of an inviscid fluid.¹¹

Fig. 5 The characteristic quantities s̄_m in (a) and k̄_m in (b) as a function of Σ for three different dimensionless rigid core radii, i.e. r̄₀ = 0 (from Tamim et al.¹⁷), r̄₀ = 0.1 and 0.9. Other parameters: [small tau, Greek, macron] = 100 and n = 0.5.

It is also interesting to investigate the situation when keeping Σ finite, and taking the limit that the Deborah number [small tau, Greek, macron] → ∞. It is then expected that viscoelastic effects dominate over inertial effects. We then introduce a dimensionless growth rate using the viscoelastic timescale τ (instead of the inertial-capillary time t_c) for the rescaling, i.e. s̄_τ ≡ sτ = s̄[small tau, Greek, macron]. The expression of β can be written as β = Σ/(1 + s̄ⁿ_τ). In the limit of [small tau, Greek, macron] → ∞, s̄ = s̄_τ/[small tau, Greek, macron] → 0, and thus α = k̄. Hence the dispersion relation (eqn (20)) consists of only the quantities r̄₀, β and k̄. Inertia plays no role in this limit. A plot of β vs. k̄ from eqn (20) for r̄₀ = 0.1 is shown in Fig. 6(a). The dimensionless maximum growth rate s̄_τ,m(k̄_m) can be computed from the value of β_m, defined as the minimum of β. Namely, s̄_τ,m = [Σ/β_m − 1]^1/n. We can then find the critical elastocapillary number by putting s̄_τ,m = 0, and thus obtain Σ_c = β_m. For r̄₀ = 0.1, we find that Σ_c = 21.7 (see Fig. 6(a)), which is the same as for the purely elastic case.

Fig. 6 (a) The relation between β and k̄ for the limiting case of [small tau, Greek, macron] → ∞. (b) and (c) The characteristic quantities k̄_m in (b) and s̄_m in (c) as a function of Σ for different values of [small tau, Greek, macron]. In (c) for [small tau, Greek, macron] → ∞, the dashed-dotted curve represents s̄_τ,m/100 vs. Σ. Other parameters: r̄₀ = 0.1 and n = 0.5.

Lastly, we compare the results when varying [small tau, Greek, macron]. In Fig. 6(b), we plot k̄_mvs. Σ for [small tau, Greek, macron] = 0, 100 and when [small tau, Greek, macron] → ∞. We can see that the curve of k̄_m for [small tau, Greek, macron] = 100 becomes flat when Σ is approaching Σ_c, which collapses with the constant value k̄_m = 0.586 for the case of [small tau, Greek, macron] → ∞. In Fig. 6(c), we plot s̄_mvs. Σ for [small tau, Greek, macron] = 0 and 100. To compare the results of [small tau, Greek, macron] = 100 and the case of [small tau, Greek, macron] → ∞, we divide s̄_τ,m for the case of [small tau, Greek, macron] → ∞ by the value [small tau, Greek, macron] = 100 since s̄ = s̄_τ/[small tau, Greek, macron], and add a curve for s̄_τ,m/100 vs. Σ in Fig. 6(c). We can see that for Σ close to Σ_c, s̄_m for τ = 100 and s̄_τ,m/100 collapse. Hence for [small tau, Greek, macron] ≠ 0, when Σ is approaching Σ_c, s̄_m ≪ 1 and viscoelastic effects dominate over inertial effects.

5 Conclusions

We investigate the onset of the PRI of a soft layer coated on a rigid cylinder by analytically deriving the dispersion relation using the linear stability analysis. We implement the Chasset–Thirion model for the viscoelastic response of the soft layer. We find that the stiffness (characterized by 1/Σ) and the rigid cylindrical core (characterized by r̄₀) have a stabilizing effect. The dimensionless growth rate of the fastest growing mode s̄_m decreases with decreasing Σ or increasing r̄₀. Importantly, there exists a critical elastocapillary number Σ_c for each r̄₀ such that the coated layer is stable for any Σ < Σ_c. The critical value depends strongly on r̄₀. For example, for r̄₀ = 0.9, the soft layer becomes unstable only when the fibre is very soft, namely when Σ ≥ 12 100. While for a soft fibre without a rigid core (r̄₀ = 0), the onset of instability occurs at Σ ≈ 6, which is three to four orders of magnitudes smaller. This remarkable result would be interesting to be verified by performing experiments for different soft coatings on a fibre. Regarding the dimensionless wavenumber of the fastest growing mode k̄_m, we find that k̄_m decreases with decreasing Σ. Interestingly, k̄_m is independent of r̄₀ for the purely elastic case.

Regarding the roles of the viscoelastic parameters [small tau, Greek, macron] and n, we find that increasing the relaxation timescale of viscoelastic material [small tau, Greek, macron] or reducing the power n can slow down the growth of disturbance. Changing either [small tau, Greek, macron] or n has no effect on the critical elastocapillary number. However, for [small tau, Greek, macron] ≠ 0, the curves for both s̄_mvs. Σ and k̄_mvs. Σ bend concavely when Σ is approaching Σ_c, and viscoelastic effects become dominant over inertial effects. The value of k̄_m becomes dependent on r̄₀ when Σ is around Σ_c. On the other hand, viscoelastic effects disappear in the limit of Σ → ∞.

Our linear stability analysis predicts the threshold of instability¹⁴ and the fastest growing modes. When the undulation of the interface grows, elastic stresses become larger and could stop the growth.^14,18 To study the evolution of the undulation and the equilibrium shape of the deformed soft layer, nonlinear elastic or viscoelastic models are required as large deformation is expected.¹⁸

6 Appendices

6.1 Expressions of the boundary conditions

Using the expressions of the general solutions (eqn (13)–(15)) for the boundary conditions (eqn (16)–(19)) yields the following set of linear equations for the unknowns A_m.

A₁k̄I₁(r̄₀k̄) − A₂ik̄I₁(r̄₀α) − A₃k̄K₁(r̄₀k̄) − A₄ik̄K₁(r̄₀α) = 0, (22)

A₁(ik̄)I₀(r̄₀k̄) + A₂αI₀(r̄₀α) + A₃(ik̄)K₀(r̄₀k̄) − A₄αK₀(r̄₀α) = 0, (23)

and

A₁2ik̄²I₁(k̄) + A₂(k̄² + α²)I₁(α) − A₃2ik̄²K₁(k̄) + A₄(k̄² + α²)K₁(α) = 0 (25)

6.2 Dispersion relation for small core radii

We expand the determinant D(k̄,s̄,r̄₀,[small tau, Greek, macron],n) of the dispersion relation (eqn (20)) in series of r̄₀ and the leading order terms are given as follows:

D(k̄,s̄;r̄₀,k̄,[small tau, Greek, macron],n) = F₁ + F₂ln(r̄₀) + …, (26)

where

F₁ = 2(C₂D₁ − C₁D₂)γ_E(α − k̄)(α + k̄) − 2iα²(B₂D₁ − B₁D₂)k̄K₀(α) + i(B₂D₁−B₁D₂)k̄(α²+k̄²)K₀(k̄) + α[2α²(B₂C₁ − B₁C₂) − i(2 + β)(B₂D₁ − B₁D₂)k̄ + 2(B₂C₁ − B₁C₂)k̄² + iβ(B₂D₁ − B₁D₂)k̄³]K₁(α) + ik̄²[(2 + β)(B₂D₁ − B₁D₂) + 4i(B₂C₁ − B₁C₂)k̄ + β(−B₂D₁ + B₁D₂)k̄²]K₁(k̄), (27)

F₂ = 2α²ln(α/2)(C₂D₁ − C₁D₂) + 2(−C₂D₁ + C₁D₂)k̄²ln(k̄/2) (28)

and

B₁ = (ik̄)I₀(r̄₀k̄),

B₂ = αI₀(r̄₀α),

D₁ = 2ik̄²I₁(k̄),

D₂ = (k̄² + α²)I₁(α),

γ _E is Euler's constant ≈ 0.57721.

6.3 Validation of our model: the Newtonian fluid limit

Our viscoelastic model reduces to the Newtonian fluid when n = 1 and in the limits of μ_o → 0 and μ_oτ → η, where η is the dynamic viscosity of the fluid. In terms of the dimensionless parameters, it means Σ → ∞ and , where Oh is called the Ohnesorge number. In this limiting situation, the dimensionless governing eqn (11) is reduced towhere we have introduced a new rescaled growth rate ŝ ≡ ηRs/γ = Ohs̄. The dispersion relation is obtained by substituting and β = 1/ŝ in eqn (20). The dimensionless control parameters are Oh and r̄₀.

6.3.1 Stokes flow case. Taking further the limit that Oh → ∞, eqn (29) reduces to the Stokes equation in Laplace space. The remaining dimensionless control parameter is only r̄₀. We validate the expression of our dispersion relation (eqn (20)) in this limiting case by comparing our results with that of Zhao et al.³² in which the dispersion relation is obtained using the normal mode method to solve the Stokes equations. Fig. 7 shows the comparison for the cases of very small fibre radius (r̄₀ = 10⁻⁹), thick coated liquid film (r̄₀ = 0.1) and thin coated liquid film (r̄₀ = 0.9).

Fig. 7 The dispersion relation between the growth rate ŝ and the wavenumber k̄ in the Stokes flow limit for r̄₀ = 10⁻⁹ in (a), r̄₀ = 0.1 in (b) and r̄₀ = 0.9 in (c). Solid lines: our results when substituting α = k̄ and β = 1/ŝ for eqn (20). Dashed lines: results of Zhao et al.³²

Data availability

Data for this article, including the codes used to compute the data shown in the figures, are available at https://doi.org/10.11582/2024.00099.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

The authors gratefully acknowledge financial support from the Research Council of Norway (Project No. 315110).

References

1. Y. Guo, J. Bae, Z. Fang, P. Li, F. Zhao and G. Yu, Chem. Rev., 2020, 120, 7642–7707.
2. Z. G. Qian, S. C. Huang and X. X. Xia, Nat. Chem. Biol., 2022, 18, 1330–1340.
3. I. Levental, P. C. Georges and P. A. Janmey, Soft Matter, 2007, 3, 299–306.
4. B. Roman and J. Bico, J. Phys.: Condens. Matter, 2010, 22, 493101.
5. A. Jagota, D. Paretkar and A. Ghatak, Phys. Rev. E, 2012, 85, 051602.
6. M. L. Wei, Europhys. Lett. (EPL), 2014, 106, 14002.
7. R. W. Style, A. Jagota, C. Y. Hui and E. R. Dufresne, Annual Rev. Cond. Matter Phys., 2017, 8, 99–118.
8. J. Bico, É. Reyssat and B. Roman, Annu. Rev. Fluid Mech., 2018, 50, 629–659.
9. J. A. F. Plateau, Acad. Sci. Bruxelles Mem., 1843, 16, 3.
10. L. Rayleigh, Proc. Roy. Soc. London, 1879, 10, 4.
11. J. Eggers and E. Villermaux, Rep. Progr. Phys., 2008, 71, 036601.
12. M. Taffetani and P. Ciarletta, Phys. Rev. E, 2015, E91, 032413.
13. C. Xuan and J. Biggins, Phys. Rev. E, 2017, 95, 1–8.
14. S. Mora, T. Phou, J. M. Fromental, L. M. Pismen and Y. Pomeau, Phys. Rev. Lett., 2010, 105, 1–4.
15. C. Lestringant and B. Audoly, Proc. R. Soc. A, 2020, 476, 20200337.
16. Y. Fu, L. Jin and A. Goriely, J. Mech. Phys. Solids, 2021, 147, 104250.
17. S. I. Tamim and J. B. Bostwick, Soft Matter, 2021, 17, 4170–4179.
18. A. Pandey, M. Kansal, M. A. Herrada, J. Eggers and J. H. Snoeijer, Soft Matter, 2021, 17, 5148–5161.
19. B. Dortdivanlioglu and A. Javili, Extreme Mech. Lett., 2022, 55, 101797.
20. G. Yang, C. F. Gao and C. Q. Ru, Int. J. Solids Struct., 2022, 241, 111491.
21. C. Q. Ru, Mech. Res. Communicat., 2022, 124, 103959.
22. P. P. Bhat, S. Appathurai, M. T. Harris, M. Pasquali, G. H. McKinley and O. A. Basaran, Nat. Phys., 2010, 6, 625–631.
23. D. Quéré, Annu. Rev. Fluid Mech., 1999, 31, 347–384.
24. H.-C. Chang and E. A. Demekhin, J. Fluid Mech., 1999, 380, 233–255.
25. I. L. Kliakhandler, S. H. Davis and S. G. Bankoff, J. Fluid Mech., 2001, 429, 381–390.
26. C. Ruyer-quil, P. Treveleyan, F. Giorgiutti-dauphiné, C. Duprat and S. Kalliadasis, J. Fluid Mech., 2008, 603, 431–462.
27. Y. Zheng, H. Bai, Z. Huang, X. Tian, F.-Q. Nie, Y. Zhao, J. Zhai and L. Jiang, Nature, 2010, 463, 640–643.
28. S. Haefner, M. Benzaquen, O. Bäumchen, T. Salez, R. Peters, J. D. McGraw, K. Jacobs, E. Raphaël and K. Dalnoki-Veress, Nat. Commun., 2015, 6, 7409.
29. A. Sadeghpour, Z. Zeng and Y. S. Ju, Langmuir, 2017, 33, 6292–6299.
30. H. Chen, T. Ran, Y. Gan, J. Zhou, Y. Zhang, L. Zhang, D. Zhang and L. Jiang, Nat. Mater., 2018, 17, 935–942.
31. H. Ji, C. Falcon, A. Sadeghpour, Z. Zeng, Y. S. Ju and A. L. Bertozzi, J. Fluid Mech., 2019, 865, 303–327.
32. C. Zhao, Y. Zhang and T. Si, J. Fluid Mech., 2023, 954, 1–22.
33. S. L. Goren, J. Fluid Mech., 1962, 12, 309–319.
34. S. L. Goren, J. Colloid Sci., 1964, 19, 81–86.
35. Y. Zhang, S. T. Ellison, S. Duraivel, C. D. Morley, C. R. Taylor and T. E. Angelini, Bioprinting, 2021, 21, e00121.
36. M. M. Fitzgerald, K. Bootsma, J. A. Berberich and J. L. Sparks, Biomacromolecules, 2015, 16, 1497–1505.
37. D. K. Nandakumar, Y. Zhang, S. K. Ravi, N. Guo, C. Zhang and S. C. Tan, Adv. Mater., 2019, 31, 1806730.
38. R. Bar-Ziv and E. Moses, Phys. Rev. Lett., 1994, 73, 1392–1395.
39. E. Hannezo, J. Prost and J.-F. M. C. Joanny, Phys. Rev. Lett., 2012, 109, 018101.
40. D. Gonzalez-Rodriguez, S. Sart, A. Babataheri, D. Tareste, A. I. Barakat, C. Clanet and J. Husson, Phys. Rev. Lett., 2015, 115, 088102.
41. M. R. King and S. Petry, Nat. Commun., 2020, 11, 270.
42. S. U. Setru, B. Gouveia, R. Alfaro-Aco, J. W. Shaevitz, H. A. Stone and S. Petry, Nat. Phys., 2021, 17, 493–498.

---