1 Introduction

The Rayleigh-Taylor instability (RTI) can be applied to any interface between fluids, where the acceleration is directed from the dense fluid to the light fluid. The linear stability analyses for a number of systems with two semi-infinite domains have been outlined^[1]. The corresponding study has been extended to finite-thickness viscous layers^[2-4], including surface tension^[5], compressible effects^[6], and nonlinear effects^[7-8]. The dripping of a pendant film provides an ideal model to study the combined contributions of these different effects.

RTI and RTI mixing can be suppressed by high-frequency vibrations^[9], rotation^[10-11], Marangoni effects^[12-14], rough surface^[15], and electric fields^[16]. The film under an inclined wall was investigated^[17] as well, and the experimental results showed that no drops formed for the inclination greater than a critical value. The dripping of a pendent film will occur when the RTI overcomes the combined effects of the viscosity and surface tension. Different from the previous studies, it is shown numerically and theoretically in this paper that the dripping process of a viscous thin film can be retarded by rough ceilings.

2 Model and numerical method

In order to study the boundary effect on the RTI, two different ceilings are considered in this paper (see Fig. 1). A layer of liquid with the viscosity μ₁, the density ρ₁, the surface tension σ, and the thickness h(x, t) coats the underside of a horizontal flat wall (see Fig. 1(a)) at y=0, where h₀, λ, and A denote the averaged film thickness, the wave length of the perturbed interface, and the amplitude of the interface variation, respectively. The surrounding fluid has the viscosity μ₂ and the density ρ₂. The whole system is subject to the gravitational acceleration g. In Fig. 1(b), a layer of liquid coats the underside of a horizontal ceiling h₁(x) corrugated in-phase with the initial interface perturbation, where δ is the initial uniform film thickness.

The two-dimensional Navier-Stokes equations for an incompressible two-phase flow are

(1)

(2)

where u, p, ρ, and μ denote the velocity, the pressure, the density, and the viscosity, respectively. f_st is the volumetric surface tension force. The continuous surface force (CSF) model is used to discretize the surface tension force

where κ, δ_I, and represent the interface curvature, the Dirac delta function, and the unit interface normal, respectively.

The governing equations are solved with an adaptively refined projection method^[18] based on a variable-density fractional-step scheme^[19], and the volume of fluid (VOF) method is used to track the fluid interface. The computational domain is taken as a L× L square, where L equals the wavelength λ of the interface perturbation. No-slip boundary conditions are used at the ceiling, and the domain bottom and periodic boundary conditions are applied on the left and right boundaries of the domain, respectively. The initial velocity is uniform zero, and the initial interface disturbance is

In the adaptive mesh-refinement process, 512 grid points per wavelength are used as the maximum grid resolution. It is reached mostly near the fluid interface and the wall, and is found to be fine enough to describe the initial stage of the interface evolution. Gravity and surface tension are two important factors for the surface instability, and their ratio is the Bond number defined by

where Δρ is the density difference between the fluids. As shown in Fig. 2(a), the critical Bo in the simulations is Bo_c≈39.23. It agrees very well with 4π², the classical theoretical threshold obtained through the balance between the static pressure and the surface tension. In addition, when the initial amplitude is large enough, e.g., A₀=9%, the interface becomes unstable at Bo=36 (see Fig. 2(b)), which is a value less than 4π². It indicates a nonlinear subcritical instability.

The amplitude suppression factor is defined by

(3)

where A_corr and A_flat denote the amplitudes of the interface perturbations with the corrugated ceiling and the flat wall, respectively. In the following dimensionless analysis, the characteristic length and time scales are

Consequently, the Bond number may be looked as the square of a dimensionless wavelength λ/l_c.

3 Numerical results and mechanism

In this section, the numerical simulations of the supercritical and subcritical instabilities of a pendent water film in air are conducted, and the following parameters are used:

3.1 Retardation of supercritical instability

In order to compare the dripping processes, the minimum initial film thickness δ underneath the flat wall is assumed to be the same as the initial thickness of the film coating a rough ceiling corrugated in-phase with the initial disturbed interface (see Fig. 1). At supercritical states, it is shown in Fig. 3(a) that the interface amplitude grows more slowly than the case with a flat ceiling, and the value of the suppression factor η_A at the beginning of the dripping process is higher than that at the later stage of the dripping process (see Fig. 3(b)). Note that according to the definition of η_A, smaller η_A corresponds to stronger suppression effects. η_A decreases first and then increases with the increase in δ (see Fig. 3(c)). Consequently, there is an optimized initial thickness δ corresponding to the maximum suppression effect (or a minimum of η_A). It is noted that, during the dripping process, different instants have different optimized initial thicknesses (see Fig. 3(c)). Since Bo represents the ratio of the gravity to the surface tension, a larger Bo represents a smaller surface tension effect. Therefore, η_A decreases with the increase in Bo (see Fig. 3(d)). The retardation mechanism will be discussed in detail in Subsection 3.3.

3.2 Retardation of subcritical instability

As shown in Fig. 2(b), the interface is unstable with the disturbance amplitude A₀=0.09λ at a subcritical state Bo=36 < 4π². The growth of such a nonlinear instability can be retarded as well with an in-phase corrugated ceiling. The temporal evolution of the amplitude A and suppression factor η_A are shown in Figs. 4(a) and 4(b) and the variations of suppression factor η_A with respect to the thickness δ and Bo are given in Figs. 4(c) and 4(d), respectively. Similar to that in the supercritical case, the suppression factor η_A first decreases and then increases with respect to the film thickness δ, and a minimum is obtained at an intermediate value in the explored δ range. In addition, η_A decreases as well with respect to Bo. The suppression factor η_A clearly shows that the amplitude growth of the subcritical instability is significantly retarded by the in-phase corrugated ceiling.

3.3 Retardation mechanism

Considering a pendent liquid film with the viscosity μ in an invscid ambient fluid with the constant pressure P₀. The interface lies at y=-h(x, t). Then, its vertical velocity is

Integrating the continuity equation and substituting v into it, we have

(4)

where h₁ is the height of the corrugated wall (see Fig. 1(b)).

The lubrication approximation^[20] is applied, where the inertia terms in the momentum equations are ignored, h≪λ. Then, the interface curvature can be approximated by κ=-∇²h, and Eq. (1) can be simplified as follows:

(5)

(6)

where p is the reduced pressure, and the total pressure is reduced by P₀. Integrating Eq. (6) in y and considering the surface tension effect, we have

Substituting p into Eq. (5) and integrating it with the boundary conditions at y=-h^[5] and u=0 at the top wall, we obtain the following equations for a flat wall:

(7)

(8)

and the following equations for a corrugated ceiling:

(9)

(10)

For a thin film with the initial uniform thickness beneath an in-phase corrugated ceiling, i.e., h(x)-h₁(x) =δ at any x, Eq. (10) can be approximated as follows:

(11)

It is shown in Fig. 5 that the solutions of Eqs. (10) and (11) coincide with each other at the initial stage of the instability. Drops will be formed around the interface of the local maximum of h, e.g., h(x=0) (see Fig. 1), where the first term on the right-hand side of Eq. (8) vanishes because ∇(ρgh+σ∇²h)=0 and thus Eq. (8) turns to be

(12)

Comparing Eqs. (11) and (12), we can conclude that the amplitude growth will be retarded by the rough ceiling if h³ -δ³ > 0. Note that h reaches its maximum value at

which is positive and increases with δ. It indicates that the growth of an unstable interface can be retarded by the in-phase corrugated ceiling or the suppression factor η_A decreases with the increase in δ (see Figs. 3(c) and 4(c)). When the film is thick or δ is large, the interface becomes far from the ceiling and then the retardation effect is weakened. Consequently, there should be a local minimum of η_A at a specific thickness δ (see Figs. 3(c) and 4(c)).

Furthermore, substituting the normal mode h'=A₀e^{wt+ iαx} into Eq. (11), we can get the growth rate as follows:

(13)

Therefore, the growth rate of a uniform-thickness thin film underneath a corrugated ceiling is proportional to the power of three of its initial thickness, i.e., w~δ³. This scaling law agrees with the numerical simulations (see Fig. 6).

4 Conclusions

It is shown numerically in this paper that a pendent film may experience supercritical or subcritical instabilities, depending on the perturbation amplitude and the Bond number, and both the two kinds of instabilities can be retarded by the in-phase corrugated ceilings. There is an optimized initial thickness of the thin film, where the suppression factor η_A reaches its minimum value, representing the maximum retardation effect. Based on the lubrication approximation, the retardation effect of the corrugated ceiling is explained in term of the reduced growth rate, whose power-law relation with the initial film thickness is revealed by the linear stability analysis and confirmed by numerical simulations.