Phys. Fluids 26, 022101 (2014); https://doi.org/10.1063/1.4863846 26, 022101
© 2014 AIP Publishing LLC.
On the stabilizing effect of a liquid film on a
cylindrical core by oscillatory motions
Cite as: Phys. Fluids 26, 022101 (2014); https://doi.org/10.1063/1.4863846
Submitted: 29 May 2013 . Accepted: 20 January 2014 . Published Online: 03 February 2014
Wilko Rohlfs, Matthias Binz, and Reinhold Kneer
ARTICLES YOU MAY BE INTERESTED IN
The interaction between viscous fingering and wrinkling in elastic-walled Hele-Shaw cells Physics of Fluids 26, 022102 (2014); https://doi.org/10.1063/1.4864188
Rayleigh-Taylor instability in thin liquid films subjected to harmonic vibration Physics of Fluids 29, 052105 (2017); https://doi.org/10.1063/1.4984082
Dynamic stabilization of the Rayleigh-Taylor instability of miscible liquids and the related “frozen waves”
On the stabilizing effect of a liquid film on a cylindrical core
by oscillatory motions
Wilko Rohlfs, Matthias Binz, and Reinhold Kneer
Institute of Heat and Mass Transfer, RWTH Aachen University, Augustinerbach 6, 52056 Aachen, Germany
(Received 29 May 2013; accepted 20 January 2014; published online 3 February 2014)
Liquid films on cylindrical bodies like wires or fibres disintegrate if their length exceeds a critical size (Plateau-Rayleigh instability). Stabilization can be achieved by an axial oscillation of the solid core provided that a suitable combination of forcing amplitude and frequency is given. To investigate the stabilizing effect, direct numerical simulations (DNS) of the axisymmetric problem are conducted with a height function based solver. It is found that the mechanism of film stabilization is caused by the interaction between an inertia dominated region (high film thickness) and a viscosity dominated region (low film thickness). Replenishing of the thin film region is thereby supported while depleting is suppressed, finally leading to a stable film flow on an oscillating cylinder. To the end, a systematic variation of the main system parameters, e.g., the Weber number, the ratio between the radius of the inner core and the average film coating thickness, and the oscillation frequency is presented and the influence of the parameters discussed.C 2014 AIP Publishing LLC.
[http://dx.doi.org/10.1063/1.4863846]
I. INTRODUCTION
From Rayleigh’s theory for jets1 and cylindrical surfaces,2 it is known that an axisymmetric liquid cylinder above a critical length is unstable against infinitesimal disturbances leading to a spontaneous disintegration into droplets (Plateau-Rayleigh instability). The break-up is caused by surface tension effects and can be explained by a balance of the surface tension forces in radial direction resulting from the curvature in axial and angular direction. The same force balance is also valid for the case of a liquid layer surrounding a solid cylinder which causes non-uniform film thickness distributions and a possible break-up of film coatings.3–5Even for cases where the fluid
layer is in motion, Joseph et al.6 found this rather simple force balance to be valid for interfacial
perturbations at leading order. In various industrial production processes, e.g., in case of an insulating coat of lacquer on a wire, a uniform coating thickness is requested and a deformation of the liquid layer before solidification has to be avoided.
A way of stabilizing liquid layers (with an unperturbed film thickness h0) is to subject the
solid core of radius R to harmonic axial oscillations (of amplitude A and angular frequencyω) as illustrated in Fig.1.
It is important to note that the oscillations do not stabilize the unperturbed cylindrical shape of the liquid coating leading to a flat film solution. However, the oscillations applied to the solid cylinder coated by a liquid layer result in a saturation of the instability, preventing film rupture and thus ensuring a continuous coating. The effect of saturation is well-known in the literature for different kinds of interfacial instabilities which will lead to a break-up of the fluid in absence of oscillatory motion. For the Rayleigh-Taylor instability, where a denser fluid rests initially on top of a less dense fluid between two plates, the effects of oscillations in vertical7 and horizontal8
direction have been investigated in many studies. Halpern and Frenkel9found a nonlinear saturation
due to horizontal oscillations of the upper plate by using a weakly nonlinear evolution equation for the disturbances. Consequently, the Rayleigh-Taylor instability fails to rupture the film in a certain window of frequencies and amplitudes. For small frequencies, the authors found stabilization, 1070-6631/2014/26(2)/022101/13/$30.00 26, 022101-1 C 2014 AIP Publishing LLC
z solid core uz(t) = Aω sin(ω t) r R h0 λ
FIG. 1. Case description.
however, with a chaotic and non-periodic surface topology. Contrary, higher frequencies result in a time-periodic solution with only one single but fast moving wave inside the given domain. For a two-phase flow which is bounded by vertically oscillating walls, Blyth10found that the mechanism
of the Rayleigh-Taylor instability is intensified by these oscillations. For a Couette flow between two horizontal plates, which is generally unstable if the layer of the more viscous liquid is thinner, Coward and Papageorgiou11found that the effect of oscillations, when superimposed to the main flow, can
either be stabilizing or destabilizing, depending on the parameters. In a further study, Coward and Renardy12considered the additional effect of a streamwise pressure gradient, consisting of a steady and an unsteady part. In Halpern et al.,13 such an oscillating pressure gradient was imposed on a core-annular flow, for which the capillary force due to the azimuthal curvature of the interface is a main destabilizing mechanism, similar to the present study. Including thermo-capillary forces, Thiele et al.14investigated the effect of vibrations on the long-wave Marangoni instability. Using lubrication theory, the authors found that vibrations have in general a stabilizing influence for different flow configurations such as a flat film flowing on a horizontal or inclined substrate. External vibrations applied to falling liquid films have also been found to have a significant influence.15The flow rate of
falling films on vertical cylindrical surfaces, for instance, can be controlled by ultrasound forcing16
by increasing the forcing amplitude. Surpassing a critical amplitude, the film flow can be completely stopped. For a partially wetted vertically vibrating inclined plate, drops have been found to move upwards (against the direction of gravity).17 This phenomenon has been related to either contact
angle hysteresis, introducing a nonlinearity into the force-speed relation or to an asymmetric area of contact between the fluid and the plate.
Haimovich and Oron18have recently derived a nonlinear, one-dimensional (in space) evolution
equation from the incompressible Navier-Stokes equations for the case of a falling film on a horizontal cylinder subjected to harmonic axial oscillations under isothermal and non-isothermal conditions19 as well as for double frequency forcing20 using the methods of long-wave theory. The authors performed numerical simulations of this equation and found a correlation for the critical amplitude,
¯ A=hA
0, of the oscillation
¯
Ac= 0.1 · −1We−0.038H0.33(1+ H)−1.3L¯1.3 (1)
which is needed to stabilize the film. This critical value depends on the Womersly number = ωh20 ν , a
geometric parameter H =h0
R, the dimensionless length of the domain ¯L= hλ0, and the Weber number
We= σ h0
ρν2. In these equationsν, ρ, σ denote the fluid properties: kinematic viscosity, density, and
surface tension. It was found that the critical amplitude of the oscillation decreases proportionally with an increase of the oscillation’s frequency. According to the authors, the correlation is valid for 2138≤ We ≤ 12833, 36π ≤ ¯L ≤ 120π, and 0.2 ≤ H ≤ 0.4.
Furthermore, the authors associate the mechanism of rupture prevention with the advection induced by the substrate oscillation. Thus, a further thinning of the film interface’s troughs is suppressed by fluid which is carried due to the vibrating surface. They further found a short-time
film thinning in the trough which becomes localized in time and replenishes during the periodical oscillation.
The present study picks up the problem using a direct numerical approach for the solution of the full set of incompressible Navier-Stokes equations. With the direct numerical approach, additional forces, such as an interaction between the fluid and the gaseous phase or higher order effects can be accounted for. The original contribution of this paper is two-fold. First, a precise analysis of the numerical simulation results revealed a detailed understanding of the saturation that prevents film rupture. Second, a variation of the main influencing parameters allowed us to verify the results of the long-wave theory and to test its accuracy outside the prescribed range of validity.
II. PROBLEM FORMULATION
Picking up the case from Haimovich and Oron,18 this study considers a horizontal circular
cylinder in zero-gravity (g = 0) with the radius R coated by an incompressible liquid layer. The average film thickness of this layer is h0, the relevant fluid properties are: density, ρ, dynamic
viscosity,μ, and surface tension, σ. The generally unstable behavior of the liquid layer is influenced by an axial harmonic oscillation (of amplitude A and frequencyω). The direct numerical simulations are performed for a cylindrical domain with radial, axial, and azimuthal coordinates (r, z,θ).
Within both phases (the liquid and the ambient gaseous phase) the incompressible Navier-Stokes equations, the mass conservation equation
1 r ∂(rur) ∂r + ∂uz ∂z = 0, (2)
the momentum conservation equation in radial ρ ∂u r ∂t + ur ∂ur ∂r + uz ∂ur ∂z = −∂p ∂r + μ 1 r ∂ ∂r r∂ur ∂r +∂2ur ∂z2 − ur r2 , (3)
and axial direction
ρ ∂uz ∂t + ur∂u z ∂r + uz∂u z ∂z = −∂p∂z + μ 1 r ∂ ∂r r∂uz ∂r +∂2uz ∂z2 , (4)
are assumed to hold, where p is the static pressure.
The boundary conditions at the surface of the axially oscillating cylinder (at r= R) are
ur = 0, uz(t)= Aω sin(ωt), (5)
assuming no slip and no penetration. From the mathematical point of view, the boundary conditions at the free surface, r= R + h(z) are the kinematic condition,
∂h
∂t = ur− uz
∂h
∂z, (6)
and the balance of the radial and axial components of the interfacial stress (dynamic boundary conditions) 2μ∂ur ∂r − μ ∂uz ∂r + ∂ur ∂z ∂h ∂z − p + σκ = 0, (7) μ ∂u z ∂r + ∂ur ∂z − 2μ∂uz ∂z ∂h ∂z + p ∂h ∂z − σκ ∂h ∂z = 0. (8)
These boundary conditions can be directly applied in the derivation of the nonlinear evolution equation.18 However, for the direct numerical approach that is used in this study, those boundary
conditions are indirectly applied by solving for the advection of the two fluids and by introducing a surface tension force in the interfacial region (see also Sec.III). In the equations above,κ denotes the mean curvature of the interface given by
κ = 1+∂h∂z2 R+ h − ∂2h ∂z2 1+∂h ∂z 2−3/2 . (9)
The stability of this system for the case of R= 0 and a symmetry boundary condition (ur= 0,
∂uz/∂z = 0, at r = 0, z) has been intensively studied and is well known asRayleigh’s theory (1878).
This case will be used in the present study to validate the numerical procedures.
III. NUMERICAL METHOD
The described two-phase problem is solved by a direct numerical simulation of the equation system(2)–(4)using the volume of fluid (VOF) approach21and the continuum surface force method22
implemented in the open-source code OpenFOAM.23 Due to axisymmetry, the problem is solved
in a two-dimensional domain (see Fig.2), which is specified as a wedge with an angle of 1◦. The length of the domain isλ and height is rmax− R. Periodic boundary conditions have been applied
in axial direction enclosing the computational domain. At r= R, a moving wall boundary condition is applied according to Eq.(5). The pressure gradient in wall normal direction is set to zero, which follows from the no-slip condition in the governing equation. At the upper boundary of the domain at r= rmax, the total pressure (p0= p + u2/2) is set to a constant value in order to allow for an inflow and
outflow of the gaseous phase. The spatial discretization is done by hexagonal cells with a uniform cell size in axial direction and a wall refinement in radial direction. For the temporal discretization, an adjustable time-step size with a constant Courant-Number of 0.01 is applied. The central differences scheme is used for all spatial discretizations and a first order, bounded implicit scheme is chosen for temporal discretization. For the surface tension, the surface curvature is calculated based on a height function approach and Eq.(9)which yields to a better approximation of the pressure distribution compared to the originally implemented method of OpenFOAM’s multiphase solver.23This solver
calculates the height of the fluid layer by integrating the fluid column in wall normal direction (which agrees with the surface definition given in Rohlfs et al.,24Eq.(10)), yielding a stencil for the fluid
height of each cell column in axial direction. Due to the predefined direction of integration, this approach is limited to reasonable surface deformations, such as investigated in the present study or, for instance, in falling films.
The classical Plateau-Rayleigh problem is simulated for a validation of the modified solver. Quantitatively, the perturbations growth behavior of different wave numbers, k= 2π/λ, have been examined and compared to analytical values. The dimensionless numerical growth rate at the time step n is computed by ¯ wn= τln εn+1 max − lnεn−1 max tn+1− tn−1 with:εmax= h − h0, (10) rmax R h0 λ periodic b.c.
periodic b.c. moving wall variable inlet/outlet r
z
analytical
numerical
1 0.8
0 0.2 0.4 0.6
dimensionless growth rate
ω
0 0.1 0.2 0.3
dimensionless wave number k
FIG. 3. Comparison between analytical and numerical solution of the growth rate for different wave numbers with ¯k= kh0.
whereτ is the characteristic time-scale of breakup, τ =
ρlh30
γ , andεmaxis the maximum amplitude
of the perturbation.25
The results are compared to the analytic growth rate according to Rayleigh (1878), ¯ w = −¯k1− ¯k2 I1 ¯k I0 ¯k , (11)
whereρldenotes the density of the liquid and I0and I1 the modified Bessel functions of first and
second kind, respectively.
Figure 3 compares the results of the numerical simulation with the analytic solution of the dimensionless growth rate ¯ω in dependence of the dimensionless wave number for an initial distur-bance ofε0= 0.01 · h0. An adequate agreement with Rayleigh’s theory is found proving that the
solver modifications are well suited for simulating the Plateau-Rayleigh problem.
IV. RESULTS
The stabilizing effect and mechanism of the harmonic oscillations is investigated by imposing sub- and supercritical values for the amplitude according to the correlation function of Haimovich and Oron.18A reference case is chosen according to the parameters studied therein, however, with
the limitation of a reduced domain/wavelength of ¯L = 3 ¯Lcrit= 6πHH+1 = 113.097. This limitation
results from the high computational effort related with the direct numerical simulation of the proposed problem. The values of the dimensionless parameters are: We= 2138, H = 0.2, and = 2π.
The effect of film stabilization is illustrated in Fig.4for which the critical amplitude is ¯Ac= 2.57
(Eq.(1)). This corresponds to an angular frequency ofω ≈ 628.321/s(equal to an excitation frequency
400
0 100 200 300
dimensionless time t 1
0.1
dimensionless film thickness h
min
[log] A = 0A = 1.59
A = 3.18
FIG. 4. Temporal evolution of the minimum non-dimensional film thickness for various forcing amplitudes; the critical amplitude is ¯Ac= 2.57 according to correlation(1).
1 0 0.4 0.6 0.8 dimensionless coordinate z 0.2 4 3 0 2 1 t = 0 t = 10 t = 500 t = 20 t = 25 t = 30
dimensionless film thickness h
FIG. 5. Temporal evolution of the film topology for a dimensionless amplitude of ¯A= 3.5.
of f= 100 Hz). The diagram shows the temporal evolution (¯t = τt withτ = ρlh20
μ ) of the minimum
film thickness ( ¯hmin= hhmin0 ) for a fixed core ( ¯A= 0), a moving core that oscillates with a sub-critical
amplitude ( ¯A= 1.59) and a supercritical amplitude ( ¯A = 3.18). The ratio between the densities and the viscosities in the liquid and the gaseous phase isρl/ρg= 692 and νl/νg= 44, respectively,
representing a low viscous oil in an air-like atmosphere.
For all three cases, an initially moderately decreasing film thickness is found for ¯t< 25. After the start-up phase, the film thickness decreases rapidly, however, with a decreasing gradient. Compared to the static case, a moving cylindrical core has an ambiguous impact on the rate of film thinning. For ¯t< 55, the rate of film thinning is higher if oscillations are present. Contrary, for higher values of ¯t film thinning is faster for the static case. If the amplitude is sufficiently high, the minimal film thickness converges to a constant level and a stabilization of the film is achieved, as it is found for the supercritical amplitude ¯A= 3.18 in Fig.4.
Figure 5 shows the topology of the film at different times for a dimensionless amplitude of ¯
A= 3.5. This value has been chosen because it provides a stable film for all variations investigated in the framework of the parameter study which is presented below. The solutions are shown at a point in time where the cylinder has reached its maximum displacement in negative direction (¯z= −A) and is motionless. From the beginning of the simulation up to a dimensionless time of ¯t= 10, the shape of the perturbed surface does not change significantly, but the trough and crest regions become more distinct. At ¯t= 20, a new crest emerges at the position of the former trough and forms a satellite droplet.5,18 From ¯t= 20 to ¯t = 25, both crests grow, fed from fluid of the intermediate
troughs which thus decrease and spread, leading to steep flanks of the droplets’ edges. At ¯t= 30, a third and a fourth satellite droplet seem to emerge at the centres of both former valleys. However, these small droplets vanish, while the film approaches its stable quasi-steady position, shown here for ¯t= 500.
In Fig.5, the evolution of the film topology was shown for a fixed position of the cylindrical oscillation (cylinder position at ¯z= −A). However, periodical changes of the film topology within one period have not been considered. Those periodical changes of the film surface during a forcing period ¯T = Tτ = 1 are shown in Fig.6for equidistantly spaced points in time within 500≤ ¯t ≤ 501 (quasi-steady) in the laboratory frame of reference. The top and bottom contour lines show the film topology for a maximal amplitude in the negative direction (¯z = −A) and a motionless cylinder. The four plots above the bottom one show the topology for a cylinder movement in positive coordinate direction, while the four plots below the top one show the topology for a cylinder movement in negative coordinate direction (returning to the initial position). Additionally, contours of constant height (over time) are given for ¯h = 0.1025 and ¯h = 1.5. From the plot, a significant movement of the lower film thickness contour can be seen, while the movement of the higher film thickness contour is comparatively low. Further, both droplets show an asymmetric shape, which is characterized by a steep wave front in upstream direction and a shallow front in downstream direction (referred to the
0 0.2 0.4 0.6 0.8 1 0 5 10 15 dimensionless coordinate z
dimensionless film thickness h
(shifted vertically by 1.2 per
Δ t = 0.1) h = 0.1025 h = 1.5 t = 501 t = 500
FIG. 6. Film topology during one period of excitation (500≤ ¯t ≤ 501) for a dimensionless amplitude of ¯A = 3.5.
direction of the cylinder movement). Such asymmetric waves shapes can also be seen in horizontally oscillating flows as shown experimentally by Jalikop.26
A. The stabilizing mechanism of oscillatory motions
A different representation of the film movement in time is given in Fig.7, showing a color-coded (grey-scaled) space-time plot. In this plot, the ordinate represents the dimensionless time and the color-coding the logarithmic value of the film thickness. A phase shift between the motion of the thin film region and the motion of the crest region is found. It is worth noting that the simulation results show an asymmetric behavior in which the droplet in the center of the periodic domain is smaller than the other droplet. Additionally, the left trough is larger in size compared to the right trough.
Information about the movement in the laboratory frame of reference makes a physical under-standing of the stabilization mechanisms rather difficult. Figure8therefore depicts the film thickness distribution in a moving frame of reference showing a constant film thickness in the trough region ( ¯h< 0.1). The vertical isolines illustrate this stationary behavior. Contrary, a significant displace-ment of the droplet is found in the moving frame of reference, while the droplet is nearly stationary in the laboratory frame (Fig.7). The phase-shifted oscillatory behavior of the smaller droplet becomes more distinct in Fig.8, and also a phase-shifted movement of the larger drop becomes obvious.
The displacement of the droplet in the moving frame of reference is associated with a net volume flow, ¯q = hqρ
0μ, on the cylinder. This flow, i.e., the volume flow integrated over the height ¯h,
is presented in Fig.9. In this plot, a flow in negative coordinate direction is illustrated by dark colors and a flow in positive coordinate direction by light colors. The two thick contour lines crossing
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 dimensionless coordinate z
film thickness distribution [log] 0
1/2 3/2 2
dimensionless phase
0 0.2 0.4 0.6 0.8 1 0 dimensionless phase Velocity Filmthicknessdistribution[log] n o i t c e l f e D dimensionless coordinate z 0 dimensionless phase z´ = 0.031
FIG. 8. Film topology during one period of excitation (500≤ ¯t ≤ 501) for a dimensionless amplitude of ¯A = 3.5; Color-coded space-time plot (0.0806 < ¯h < 2.946).
the plot in horizontal direction mark the position in space and time for which the velocity in axial direction is zero. The two thick black arrows show the general flow direction.
For the trough region, the movement of the flow corresponds directly to the acceleration of the cylinder which follows a cosine function (similar to the cylinder position, but with an opposite sign). For the case that the acceleration is zero (¯t= 1/2π and 3/2π), the velocity in this region becomes zero, while for the case of maximal acceleration (¯t= 0 and π), the velocity is maximal. Because of the direct connection between acceleration and flow behavior, we associate this region as viscosity dominated.
The maximal volume flow at the wave crests, i.e., in the droplet region, occurs near ¯t= 1/2π and ¯t= 3/2π, when the cylinder displays its maximal velocity and the deflection of the drop is zero. However, the direction of the volume flow is opposite to the direction of the cylinder movement. In the laboratory frame of reference, both velocities cancel each other out for the most part, leading to an almost stationary droplet. Because of this behavior, we associate this region as inertia dominated. For the causing effect of stabilization, i.e., the prevention of film rupture, the continuous flooding/replenishing of the troughs has been found to play an important role.18 However, the
external harmonic oscillation has to generate a somehow asymmetric flow behavior in order to
replenishing depleting 0 1 2 1 2 1 2 1 2 dimensionless coordinate z 0 0.2 0.4 0.6 0.8 1
dimensionless film thickness h
Net volume flow viscosity dominated inertia dominated q = 0 q = 0 dimensionless coordinate z 0 0.2 ~ h 0.4 0.6 1 min ~ hmin 0 1/2 3/2 2 dimensionless phase t = 1.88 t = 1.72 t = 0.84 t = 0.72
FIG. 9. Volume flow during one period of excitation (500≤ ¯t ≤ 501) for a dimensionless amplitude of ¯A = 3.5; Color-coded space-time plot (−5.43 · 103< ¯q < 5.21 · 103).
suppress the constantly acting surface tension force (in de-wetting direction). In the left plot in Fig.9, two vertical lines shown at ¯z= 0.28 and ¯z = 0.75 indicate the location of minimal film thickness ( ¯hmin). The flow towards and away from the position of minimal film thickness is shown
for four instances of time in the four right plots. The respective point in time is marked in the left plot by the white horizontal dotted lines. For the bottom plot (¯t= 0.72π), the cylinder moves from the left to the right. While the viscosity dominated region follows the cylinder movement directly, the flow in the inertia dominated region is still moving in negative coordinate direction. The small white and black arrows indicate the mechanistic effect of the volume flow in the trough region. If the arrow points in the direction of the minimal film thickness (black arrow), the flow acts replenishing. Otherwise, the flow acts depleting. In the lower two plots (¯t= 0.72π and ¯t = 0.84π), the replenishing zone is on the left side of the position of minimal film thickness. In the upper two plots (¯t= 1.72π and ¯t = 1.88π), the replenishing zone is on the right side. In all four plots it becomes obvious, that – compared to the depleting zone – the replenishing zone is larger in size. This asymmetry is caused by the droplet in the inertia dominated zone which suppresses the depleting due to the counter-current flow direction. Contrary, the flooding mechanism in the viscosity dominated region is supported.
B. Influence of system parameters
Having identified the underlying mechanism for the stabilization of the film flow, the influence of the main system parameters, e.g., the Weber number, the geometric parameter H, and the oscillation frequency is investigated in this section. Special attention is drawn to system conditions at the edge and outside of the validity range of the asymptotic model.18
Following the empirical correlation (1), the critical value of the product ¯Ac· , a critical
velocity, is constant for a given system. As a result, a variation in frequency and amplitude, keeping the product of the two parameters constant, should not influence the stability. Figure10shows the temporal evolution of the minimum film thickness for a lower frequency (reduced by a factor of two compared to the reference case presented in Sec.IV A) and for cases with a ten- and hundred-fold increase in frequency. The temporal evolution shows for all four cases a decay of the minimum film thickness, however it converges toward a constant value in the long-term perspective. Although the prevention of a film rupture for an infinite time frame cannot be guaranteed as the simulations were conducted within a finite time frame, a general saturation behavior of the Rayleigh-Plateau instability is obvious. Figure11shows the spatial distribution of the film thickness covering the oscillating cylinder. The low frequency case is characterized by a single-droplet, contrary to the other cases in which two droplets are formed. The size of the single droplet is larger due to the accumulation of almost all liquid in this one drop. The two cases with an increased frequency are very similar to each other. Compared to the reference case, the difference in size of the two droplets
0 100 200 300 400 500 0.1 1 = , A = 6.36 = 2 , A = 3.18 = 20 , A = 0.318 = 200 , A = 0.0318
dimensionless film thickness h
min
dimensionless time t
FIG. 10. Temporal evolution of the minimum non-dimensional film thickness for various frequencies and amplitudes. The product ¯Ac· remains constant in value.
0 0.2 0.4 0.6 0.8 1 0 0.5 1 1.5 2 2.5 3 3.5 4 = , A = 6.36 = 2 , A = 3.18 = 20 , A = 0.318 = 200 , A = 0.0318
dimensionless film thickness h
dimensionless coordinate z
FIG. 11. Quasi-steady film topology on an oscillating cylinder moving with various frequencies and amplitudes. The product ¯
Ac· remains constant in value.
formed is reduced, while the minimal film thickness decreases with an increase in frequency. Another interesting observation is the high degree of symmetry of the film topology found for the cases with higher excitation frequencies.
In order to confirm the weak dependency of the stability problem on the Weber number found by Haimovich and Oron,18numerical simulations have been performed for 534.5< We < 8552. The
temporal evolution of the minimum film thickness (Fig.12) reveals a saturation independent from the value of the Weber number. However, the minimal film thickness is found to strongly decrease with an increase in Weber number, and thus with an increasing surface tension force. Additionally, the Weber number significantly influences the temporal evolution, leading to a faster decay for high Weber numbers. The spatial distribution of the film thickness is shown in Fig. 13. Similar to a reduced excitation frequency, low Weber numbers result in the formation of one large and a second smaller drop, while for higher Weber numbers, the sizes of the two drops approach each other.
Finally, the geometric parameter H has been varied, which describes the ratio between the average film thickness h0 and the radius of the inner solid core R. The value of the average film
thickness h0is not varied, with the result that all dimensionless variables remain constant. However,
a modification of the cylinder radius R has an additional effect on the effective wavelength of the domain, e.g., the critical wavelength Lcrit. The dimensional length of the domain has therefore been
modified in order to remain L = 3 · Lcrit. Figure14shows the temporal evolution of the minimum
film thickness for the reference case and four cases with a smaller diameter of the inner core. Due to the reduced ratio between inner core and film coating, the boundary layer model is not valid, resulting
0 100 200 300 400 500
0.1 1
dimensionless time t
dimensionless film thickness h
min We = 1069 We = 2138 We = 8552 We = 4276 We = 3207 We = 534.5
0 0.2 0.4 0.6 0.8 1 0 0.5 1 1.5 2 2.5 3 3.5 4 dimensionless coordinates z
dimensionless film thickness h
We = 1069 We = 2138 We = 8552 We = 4276 We = 3207 We = 534.5
FIG. 13. Quasi-steady film topology on an oscillating cylinder for various Weber numbers.
0 100 200 300 400 500
0.01 0.1 1
dimensionless film thickness h
min dimensionless time t H = 0.2 H = 0.25 H = 0.5 H = 0.75 H = 1
FIG. 14. Temporal evolution of the minimum non-dimensional film thickness for various values of the geometric parameter H. 0 0.2 0.4 0.6 0.8 1 0 0.5 1 1.5 2 2.5 3 3.5 4 H = 0.2 H = 0.25 H = 0.5 H = 0.75 H = 1
dimensionless film thickness h
dimensionless coordinate z
in a limitation of the correlation of Haimovic and Oron18to H< 0.4. The direct numerical simulations are not subject to this limitation, so that values up to H= 1 (equal thickness of inner core and film coating) have been investigated. For those higher values of H, the numerical simulations reveal a continuous reduction in film thickness up to the point where the film ruptures and the simulation stops. This is found to occur before ¯t= 500 for H ≥ 0.5. For the case H = 0.25, a decrease in film thickness is found up to the maximal time simulated, not showing any sign on saturation. Contrary to the results of the direct numerical simulation which suggest an increase in critical amplitude for higher values of H, the correlation function (see Eq.(1)) reveals a dependency of Ac∝ H−0.977for
0.2≤ H ≤ 0.4 and a constant value of L/Lcrit. The reason for the contrary behavior is still unknown
and should be subject of further research. For this, a direct comparison of the film thickness profiles obtained by the two different modeling approaches with respect to the geometric parameter H in the validity range of the evolution equation is necessary. The film thickness distributions in spatial direction are shown in Figure15. Similar to what has been found for the high frequency cases, an increase in the parameter H results in a higher symmetry of the film coating. The size of both, the smaller and the larger droplet, decreases with larger values of H due to the axisymmetrical effects and the reduced amount of liquid.
V. SUMMARY
In conclusion, the stabilizing effect of a liquid film on a cylinder by harmonic oscillatory motions was investigated by way of direct numerical simulations. The main stabilizing mechanism is found to be the replenishing flow in the thin film region. Due to dominating viscous forces in this region, the flow is directly coupled to the acceleration of the cylinder. Contrary, the region of high film thickness (droplet region) is dominated by inertia forces for which the fluid motion is phase-shifted with respect to the motion of the cylinder. This phase shift of the flow in the two regions supports the replenishing mechanism while suppressing depleting, finally stabilizing the film flow.
A systematic variation of the main system parameters, i.e., the Weber number, the geometric parameter H, and the oscillation frequency, reveals the following: First, the oscillation frequency and a reciprocal modification of the forcing amplitude as well as a variation of the Weber number do not influence the stability of the system, which corresponds to the results of the nonlinear evolution equation.18Second, an increase of the geometric parameter H (the ratio between the inner core and
the film coating) prevents the minimal film thickness from saturating, finally leading to film rupture. These findings are contrary to the results of the nonlinear evolution equation, especially for large values of H, for which the validity of the model is not given.
ACKNOWLEDGMENTS
The authors thank Alexander Oron for his constructive comments on a preliminary version of this article. This work was financially supported by the Deutsche Forschungsgemeinschaft (Grant No. DFG KN 764/3-2).
1W. S. Rayleigh, “On the instability of jets,”Proc. London Math. Soc.s1-10, 4–13 (1878). 2W. S. Rayleigh, “On the instability of cylindrical fluid surfaces,”Philos. Mag.34, 177–180 (1892). 3S. L. Goren, “The instability of an annular thread of fluid,”J. Fluid Mech.12, 309 (1962). 4S. L. Goren, “The shape of a thread of liquid undergoing breakup,”J. Colloid Sci.19, 81–86 (1964).
5F. Mashayek and N. Ashgriz, “Instability of liquid coatings on cylindrical surfaces,” Phys. Fluids 7, 2143–2153
(1995).
6D. D. Joseph, K. Nguyen, and G. Beavers, “Non-uniqueness and stability of the configuration of flow of immiscible fluids
with different viscosities,”J. Fluid Mech.141, 319–345 (1984).
7G. H. Wolf, “Dynamic stabilization of interchange instability of a liquid-gas interface,”Phys. Rev. Lett.24, 444–446
(1970).
8E. Talib and A. Juel, “Instability of a viscous interface under horizontal oscillation,”Phys. Fluids19, 092102 (2007). 9D. Halpern and A. L. Frenkel, “Saturated Rayleigh-Taylor instability of an oscillating Couette film flow,” J. Fluid Mech.
446, 67–93 (2001).
10M. Blyth, “Effect of pulsation on two-layer channel flow,”J. Eng. Math.59, 123–137 (2007).
11A. V. Coward and D. T. Papageorgiou, “Nonlinear stability of oscillatory core-annular flow: A generalized
12A. V. Coward and Y. Y. Renardy, “Small-amplitude oscillatory forcing on two-layer plane channel flow,”J. Fluid Mech. 334, 87–109 (1997).
13D. Halpern, J. A. Moriarty, and J. B. Grotberg, “Capillary-elastic instabilities with an oscillatory forcing function,” in
IUTAM Symposium on Non-linear Singularities in Deformation and Flow(Springer Science + Business Media Dordrecht,
1999), pp. 243–255.
14U. Thiele, J. M. Vega, and E. Knobloch, “Long-wave Marangoni instability with vibration,”J. Fluid Mech.546, 61–87
(2005).
15S. P. Lin and J. N. Chen, “The mechanism of surface wave suppression in film flow down a vertical plane,”Phys. Fluids 10, 1787–1792 (1998).
16L. Moldavsky, M. Fichman, and A. Oron, “Dynamics of thin liquid films falling on vertical cylindrical surfaces subjected
to ultrasonic forcing,”Phys. Rev. E76, 045301 (2007).
17P. Brunet, J. Eggers, and R. D. Deegan, “Vibration-induced climbing of drops,”Phys. Rev. Lett.99, 144501 (2007). 18O. Haimovich and A. Oron, “Nonlinear dynamics of a thin liquid film on an axially oscillating cylindrical surface,”Phys.
Fluids22, 032101 (2010).
19O. Haimovich and A. Oron, “Nonlinear dynamics of a thin nonisothermal liquid film on an axially oscillating cylindrical
surface,”Phys. Rev. E84, 061605 (2011).
20O. Haimovich and A. Oron, “Nonlinear dynamics of a thin liquid film on an axially oscillating cylindrical surface subjected
to double-frequency forcing,”Phys. Rev. E87, 052403 (2013).
21C. W. Hirt and B. D. Nichols, “Volume of fluid (VOF) method for the dynamics of free boundaries,”J. Comput. Phys.39,
201–225 (1981).
22J. Brackbill, D. Kothe, and C. Zemach, “A continuum method for modeling surface tension,”J. Comput. Phys.100,
335–354 (1992).
23H. Rusche, “Computational fluid dynamics of dispersed two-phase flows at high phase fractions,” Ph.D. thesis (Imperial
College, University of London, 2002).
24W. Rohlfs, F. G. Dietze, D. Haustein, and R. Kneer, “Two-phase electrohydrodynamic simulations using a volume-of-fluid
approach: A comment,”J. Comput. Phys.231, 4454–4463 (2012).
25J. Delteil, S. Vincent, A. Erriguible, and P. Subra-Paternault, “Numerical investigations in Rayleigh breakup of round liquid
jets with VOF methods,”Comput. Fluids50, 10–23 (2011).
26S. V. Jalikop and A. Juel, “Oscillatory transverse instability of interfacial waves in horizontally oscillating flows,”Phys. Fluids24, 044104 (2012).
