PHYSICAL REVIEW E 89, 053020 (2014)

Transient dynamics of confined liquid drops in a uniform electric field Shubhadeep Mandal, Kaustav Chaudhury, and Suman Chakraborty* Department of Mechanical Engineering, Indian Institute of Technology, Kharagpur, Kharagpur, 721302, India (Received 7 January 2014; revised manuscript received 10 April 2014; published 22 May 2014) We analyze the effect of confinement on the transient dynamics of liquid drops, suspended in another immiscible liquid medium, under the influence of an externally applied uniform dc electric field. For our analysis, we adhere to an analytical framework conforming to a Newtonian-leaky-dielectric liquid model in the Stokes flow regime, under the small deformation approximation. We characterize the transient relaxation of the drop shape towards its asymptotic configuration, attributed by the combined confluence of the charge-relaxation time scale and the intrinsic shape-relaxation time scale. While the former appears due to the charge accumulation process on the drop surface over a finite interval of time, the genesis of the latter is found to be intrinsic to the hydrodynamic situation under consideration. In an unbounded condition, the intrinsic shape-relaxation time scale is strongly governed by the viscosity ratio, defined as the ratio of dynamic viscosities of the droplet and the background liquid. However, when the wall effects are brought into consideration, the combined influence of the relative extent of the confinement and the intrinsic viscosity effects, acting in tandem, alter this time scale in a rather complicated and nontrivial manner. We reveal that the presence of confinement may dramatically increase the effective viscosity ratio that could have otherwise been required in an unconfined domain to realize identical time-relaxation characteristics. We also bring out the alterations in the streamline patterns because of the combinations of transient and confinement effects. Thus, our results reveal that the extent of fluidic confinement may provide an elegant alternative towards manipulating the transient dynamics of liquid drops in the presence of an externally applied electric field, bearing far-ranging consequences towards the design and functionalities of several modern-day microfluidic applications. DOI: 10.1103/PhysRevE.89.053020

PACS number(s): 47.65.−d

I. INTRODUCTION

The electromechanics of liquid drops has been a problem of immense interest to the scientific and industrial communities, as is attributable to its far-ranging consequences in several practical applications. Such applications range from natural processes such as the disintegration of rain drops in a thunderstorm [1], electric breakdown of insulating dielectric liquids on account of the presence of tiny water droplets [2], industrial processes such as enhanced coalescence and demixing in emulsions [3], electrohydrodynamic atomization [4], ink-jet printing [5,6], and biomicrofluidic applications such as electrowetting [7] and protein transfection into cells by interdroplet collision [8]. In all these applications, external electric fields may act as an effective means of manipulating the drop by inducing an interfacial stress effective over the relevant spatiotemporal scales [9–14]. The electrohydrodynamics of liquid drops is classically addressed by the leaky dielectric model, pioneered by Taylor [15] and subsequently advanced by Melcher [16] and others [13,17,18]. An underlying consideration behind this model is to assume that the two fluid phases have finite electrical conductivities, which allows for free charge accumulation at the interface, leading to the possibility of a net interfacial electrical shear [15–17,19]. The second assumption is based on the fact that the time scale of charge relaxation on account of electrical conduction from the bulk to the drop surface turns out to be significantly shorter than the convective time scale [15–17,19], so that the electric field equations may be decoupled from the momentum equations [13,15–17,20]. On

*

Corresponding author: [email protected]

1539-3755/2014/89(5)/053020(17)

application of an electric field, free charges accumulate on the droplet surface, which is attributable to a mismatch between the dielectric properties of the fluids. This induces a tangential component of electrostatic stress, in addition to a normal component. Triggered by this force, a flow field is induced in both the droplet phase and the bulk phase, which often has toroidal vortical structures. Accordingly, a viscous stress is induced, in an effort to balance the tangential electrical stress [13,15–17]. Typically, in response, the drop deforms either into a prolate (the droplet elongates in the direction of the applied electric field) or an oblate (the droplet elongates in the direction perpendicular to the applied electric field direction) spheroid shape, based on the distinctive electrical and hydrodynamic properties of the droplet and the bulk phases. In his pioneering effort, Taylor solved the steady state axisymmetric electrohydrodynamics inside and outside spherical drops in an unbounded domain subjected to a uniform electric field over the creeping flow regime [15]. His studies fundamentally established two important electrical parameters, namely the conductivity ratio R = σ1 /σ2 (where σ is the electrical conductivity, subscript 1 denotes the droplet phase and subscript 2 denotes the suspending medium) and permittivity ratio S = ε1 /ε2 (where ε is the electrical permittivity). He demonstrated that the electric field creates a rotational flow inside the drop, consisting of four vortices of equal strengths that are matched by corresponding vortices in the ambient fluid. He showed that for R < S, the ambient flow directs from the pole to the equator, whereas the flow direction is the opposite for R > S. For R = S, fluid flow ceases because of the charge-free interface. He further introduced a characteristic function T (also called the Taylor discriminating function) to predict the nature of drop deformation. The evolution of the drop from spherical shape

053020-1

©2014 American Physical Society

MANDAL, CHAUDHURY, AND CHAKRABORTY

PHYSICAL REVIEW E 89, 053020 (2014)

to oblate spheroid shape turns out to be a hallmark of T < 0, whereas its evolution to prolate spheroid shape turns out to be a scenario characterizing T > 0. For T = 0, the zero deformation of the drop signifies an intricate interplay of the electrical and hydrodynamic stresses. Since the pioneering work of Taylor [15] several studies have been reported concerning the dynamics of leaky dielectric drops [15,16,19,21–27]. However, a majority of these studies considered the steady state electrohydrodynamics in an unbounded domain. On the other hand, it has been well recognized, thanks to several contemporary applications [6,28] in the field of microfluidics, that the deformation time history of a leaky dielectric drop may turn out to be of immense practical consequence. The genesis of this significance primarily lies in the fact that in several applications of these kinds, information on the relative extent of the pertinent time scales of the electromechanics involved as compared to the characteristic system time scale as well as the transience of the flow development may act as a key towards the optimal design and functionalities of the concerned devices and systems. Therefore, it is not surprising that several studies have been reported in the literature, spanning over the last five decades [11,12,14,29–32], to bring out the implications of the transient electrohydrodynamics of liquid drops. It may be noted here that under transient conditions, the relaxation of the drop to oblate or prolate shapes occurs monotonically in an exponential function of time, as evidenced from the time history of the Taylor deformation parameter [11,12,14,26]. Such an evolution may be characterized by an appropriate relaxation time scale [11,12,29,30]. These monotonic behaviors prevail only when charging of the droplet interface takes place instantaneously. In the presence of finite charge relaxation time, droplet deformation behavior may turn out to be nonmonotonic and a lag in deformation relaxation could be notable [32]. Even in the presence of inertia, overshooting in deformation and decaying oscillation before reaching a final steady state deformation have been reported [32–34]. While both steady state electrohydrodynamics as well as unsteady state electrohydrodynamics of a liquid drop in an unbounded domain are reasonably well understood, significant insights have not yet been obtained about the effect of the domain confinement on the transient electrohydrodynamics of a drop. This consideration, however, finds its relevance in several microfluidics-based applications of contemporary relevance in which characteristic length scales of the confinement may approach the order of the drop size [6,28]. Although Taylor’s classical solution [15] may still be used as an estimate of the minimum distance over which confinement effects become consequential, wall effects must be introduced into the solution to bring out the consequences of confinement from a comprehensive quantitative perspective. Motivated by this proposition, Esmaeeli and Behjatian [20], in a recent study, solved the electrohydrodynamic equations for a spherical liquid drop in a confined domain. However, the effects of the characteristic transiences, which are important for contemporary microfluidic applications, were not addressed in their work. Another important effect that has been considered recently [14,32,35] in related contexts concerns the influence of finite charge relaxation time, which can produce a nonmonotonic deformation behavior of the droplet, although

the pertinent implications under confinement effects are yet to be brought out. In fact, a comprehensive review of the literature in the field reveals that no study has been reported in the literature that simultaneously considers the effects of unsteadiness and confinement on the electrohydrodynamics of the deformation of liquid drops in the presence of finite charge relaxation time. Here we report an analytical study to bring out the effect of confinement on the transient electrohydrodynamics of a drop subjected to a uniform electric field. In the small deformation regime, we investigate the effects of domain confinement in the presence of finite charge-relaxation time on the transient droplet deformation characteristic. Our analysis, thus, aptly retains the contribution of the charge-relaxation time scale and the intrinsic shape-relaxation time scale in governing the temporal behavior of the droplet as well as the underlying hydrodynamics. The central findings from our studies reveal that confinement effects hold the ability to endorse an elevation in the intrinsic shape-relaxation time scale of the drop, and can well surpass the effect of the charge-relaxation time scale in a highly confined environment, even though the interfacial charging process consumes a much longer time in the unbound domain. Therefore, any nontrivial implications of the combined influence of these two time scales occurs at a moderate extent of confinement, as shown herein. Subsequently, in the limit of instantaneous interface charging, we provide an explicit dependency of the intrinsic shape-relaxation time scale on the relative extent of the confinement and viscosity ratio; the relationship is found to be far from trivial. We consider a combination of confinement and viscous effects to draw equivalence with the case of an unconfined domain, albeit with an enhanced effective viscosity ratio. This equivalence demonstrates that to match the elevated intrinsic shape-relaxation time scale, realized otherwise in a confined domain, the droplet viscosity needs to be enormously large in an otherwise unbounded condition. Such high viscosities, however, are often impractical. Thus, fluidic confinement may provide an easy alternative towards manipulating the time-relaxation properties of liquid drops under transient conditions. We also demonstrate the time history of the streamline patterns under confined conditions. II. MATHEMATICAL MODELING AND THEORETICAL ANALYSIS

We consider a spherical liquid drop of radius Rd , immersed into another immiscible liquid domain, which is confined within a rigid, spherical container of radius Rb , as shown in Fig. 1. The system is otherwise quiescent, but is affected  of strength with a uniform-unidirectional dc electric field (E) E0 , as highlighted in Fig. 1. The droplet is neutrally buoyant and both the liquid phases are considered as leaky dielectric with constant electrohydrodynamic properties. Due to the application of the electric field, shear and normal stresses are induced over the drop surface, which results in electrodeformation of the droplet and subsequent fluid flow in and around the droplet. Here we consider that both the droplet and the suspending medium remain concentric throughout. Moreover, the flow field is assumed to be axisymmetric with respect to the direction of the applied electric field (z axis) [11,12,14,15,20].

053020-2

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . . z ϕ

r

θ Rd

E = E0 iˆz

ρ1 , σ 1,ε 1 μ1 ρ2 , σ 2,ε 2

Rb

μ2

FIG. 1. Schematic of a suspended droplet of radius Rd in a spherical confined geometry of radius Rb . The situation is essentially depicted in a spherical coordinate system (r,θ,ϕ).

This assumption may not be true always; nonaxisymmetric flow can be observed due to the application of a strong electric field [36]. However, in the present analysis we restrict ourselves to the condition that the applied electric field is small enough so that the flow field remains axisymmetric. Due to this symmetry, the spherical coordinate system (r,θ,ϕ) is an obvious choice, with no ϕ dependence; the origin of the coordinate system coincides with the center of the droplet and the fluidic confinement.

A. Governing equations and boundary conditions

Throughout our analysis we represent the governing equations and boundary conditions in their respective nondimensional forms, normalized using the conventional scales [11,14]: length lc = Rd , velocity uc = ε2 E02 Rd /μ2 , time

(e1) (e2) (e3) (e4)

PHYSICAL REVIEW E 89, 053020 (2014)

tc = μ2 Rd /γ , pressure pc = μ2 uc /Rd , electric potential φc = E0 Rd , and surface charge density qc = ε2 E0 . Here γ is the coefficient of the surface tension. Additionally, the system under consideration can be characterized by the following property ratios [10–12,14–18,20,25,36–41]:R = σ1 /σ2 , S = ε1 /ε2 , and λ = μ1 /μ2 , where σ , ε, and μ represent electrical conductivity, electrical permittivity, and dynamic viscosity, respectively (subscript 1 is used for the droplet phase and subscript 2 is used for the outer liquid phase throughout this article). Before going into the details of our present theoretical analysis, we first discuss the choices of the different physical electrohydrodynamic properties (conductivity, permittivity, viscosity, etc.), considered for the present analysis. Following the reported literature [11,12,14,20,32,35,36,42,43], a wide range of physical parameters may be plausible. From this wide range, we choose some representative parametric values, as shown in Table I. While other values of physical properties might well be adopted, those are unlikely to result in any qualitative influence on the results reported here, as verified by exhaustive theoretical simulations. 1. Governing equation for electric field distribution

In the present study, we consider the applied electric field  This can be attributed to the  × E = 0). E to be irrotational (∇ fact that the magnetic field induced due to the application of the electric field is negligible for a small dynamic electric current [16,17]. The irrotationality condition implies E =  where φ represents the scalar electric potential. Now, in −∇φ, the absence of bulk free charges and for the constant electrical properties in their respective phases, the electric potential in both phases (inside and outside of the droplet) satisfy the Laplace equation ∇ 2 φi = 0

(i = 1,2).

Here φ1 and φ2 represent the electric potential inside and outside the droplet, respectively. Equation (1) is supplemented with the following boundary and interfacial conditions:

φ1 (0,θ ) should be bounded φ1 (1,θ ) = φ2 (1,θ ) cos θ , where α = Rd /Rb is the confinement ratio φ2 (1/α,θ ) = − α  · n ˆ = q(t,θ ), where q(t,θ ) is surface charge density, −Q∇φ

with nˆ being the unit outward normal to the drop surface. Note that the dimensionless permittivity Q = 1 and Q = S in the suspending medium and within the droplet, respectively. Here F  = F2 − F1 denotes a jump in quantity F across the drop interface. One important thing to note regarding the boundary and interfacial conditions is that we are applying the interfacial conditions at r = 1. However, due to deformation, the droplet interface is not perfectly spherical but a function of space and time, which is not known a priori. Nevertheless, here we are considering the small deformation scenario, owing to the small

(1)

(2)

strength of the applied electric field E0 . In general, the drop deformation is quantified by the Taylor deformation parameter D = (L − L⊥ )/(L + L⊥ ), where L and L⊥ are the length of the principle axes of the deformed droplet, oriented parallel and normal to the direction of the applied electric field, respectively. The small deformation approximation necessitates |D|  1 [11,12,14,15,20]. Therefore, we can use the domain perturbation method by taking D as the perturbation parameter. By performing a formal domain perturbation analysis it can be shown that the above mentioned boundary conditions are correct for determining the electric potential to the leading

053020-3

MANDAL, CHAUDHURY, AND CHAKRABORTY

PHYSICAL REVIEW E 89, 053020 (2014)

TABLE I. Physical properties of droplet and suspending medium phases. Here ε0 = 8.854 × 10−12 [F m−1 ] is the permittivity of free space. The properties are taken from Refs. [32,36]. σ [S m−1 ]

ε/ε0

μ [kg m−1 s−1 ]

ρ [kg m−3 ]

Rd [mm]

γ [N m−1 ]

System A Droplet Suspending medium

1.2 × 10−12 4.5 × 10−11

3.0 5.3

0.0484 2.9040

960 960

4.0

3.2 × 10−3

System B Droplet (Silicone oil) Suspending medium (Castor oil)

1.25 × 10−12 4.5 × 10−11

3.0 5.3

0.1 0.69

960 961

0.5

4.5 × 10−3

order approximation (for further details see the Appendix in Ref. [20]). Now, in connection to the jump in electric displacement due to the charging of the droplet interface [criterion (e4)], it is imperative to estimate the distribution of the surface charge density q (t,θ ) by following the charge conservation law at the drop interface [17,32,42] Sa ∂q  s q − q nˆ · (nˆ · ∇)  u] = −β ∇φ  · n ˆ , + ReE [ u·∇ 2       Oh  ∂t

3

2 (3)  − n(  represents the surface gradient ops = ∇ ˆ nˆ · ∇) where ∇ erator, β is the dimensionless conductivity which takes a value of 1 in the suspending √ medium and R within the droplet, Sa = ε2 ν2 /σ2 Rd2 , Oh = μ2 / ρ2 Rd γ (the Ohnesorge number), and ReE = ε2 uc /σ2 Rd (the electric Reynolds number). In Eq. (3), while the transient term 1 characterizes the displacement current, the convection term 2 contains the contribution of both the convection of surface charges due to the tangential  s q), and the convection fluid velocity at the interface ( u·∇  u). Finally, due to the dilation of the interface (q nˆ · (nˆ · ∇) the term 3 , is reminiscent of the Ohmic conduction [14,32]. For a critical assessment of the significances of the different contributing factors in Eq. (3), we consider the order of magnitude analysis as follows: charge relaxation time scale (τe ) Sa ε2 /σ2 = = 2 μ2 Rd /γ capillary time scale (τc ) Oh ε2 /σ2 charge relaxation time scale (τe ) . = O ( 2 ) = ReE = Rd /uc convection time scale (τf )

O( 1 ) =

(4) It is important to mention that Ha and Yang [43] presented the values of charge-relaxation time scales τe (=ε2 /σ2 ) ranging from 0.02 to 4.68 s. Very small values of τe implicate instantaneous charging of the droplet interface. However, for most of the leaky dielectric liquids the conductivity is less, which finally leads to a finite charge-relaxation time scale. Thus, it is imperative to consider the finite charge-relaxation time of the drop and suspending liquids, as is considered here. In terms of a dimensionless number, the consideration of finite charge-relaxation time scale degenerates to a nonzero value of Sa. Now, referring to the physical parameters as depicted in Table I, it follows that O( 2 )  O( 3 ).

Sa ∂q  · n, ˆ = −β ∇φ Oh2 ∂t which is used for the subsequent analysis.

(6)

2. Hydrodynamic modeling and electrohydrodynamic coupling

1

O( 2 )  O( 1 ) and

Under the purview of criterion (5), we neglect the contribution of the charge convection from Eq. (3), without sacrificing the essential physics of interest. With these arguments in the background, Eq. (3) transforms to

(5)

Considering previously specified normalizing scales, the nondimensional form of the Navier-Stokes equation with the condition of incompressibility reads as Re ∂ u  u = −∇p  + χ ∇ 2 u + fe + Re u·∇ Ca ∂t  · u = 0, ∇

and (7)

where χ is the normalized viscosity which takes a value of 1 in a suspending medium and λ inside the droplet. Here the Reynolds number (Re = ρ2 uc Rd /μ2 , with ρ2 being the density of the suspending medium) signifies the importance of nonlinear inertia, whereas the ratio of the Reynolds number and Capillary number (Ca = μ2 uc /γ ) signifies the importance of transient inertia. Note that Re/Ca can be expressed in terms of Oh using the relationship Re/Ca = Oh−2 [11], following √ the definition Oh ≡ μ2 / ρ2 Rd γ . Now, Eq. (7) can be written using Oh as ∂ u  · u = −∇p  + χ ∇ 2 u + fe . + Re u∇ (8) ∂t The body force fe , as appears in Eqs. (7) and (8), represents the electric force density which can be found by taking the divergence of the Maxwell stress tensor (T M ) [11,17,20,36] Oh−2

 ∇ε,  · T M = ρe E − 1 (E · E)  fe = ∇ 2

(9)

where ρe is the bulk free charge density. The first term on the right-hand side of Eq. (9) is the Coulomb force and the second term is due to polarization stress. The leaky dielectric model [11,12,14–17] assumes free charge density in the bulk to be zero. Further, taking dielectric properties as constant, we finally obtain fe = 0 in the bulk. However, electric charges are present on the interface. Thus, electrical stress across the interface does not vanish, which results in electromechanical coupling at the interface [11,12,14–17]. For the present study, we are considering the low strength of the applied electric field (E0 ), which yields Re  1.

053020-4

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . .

PHYSICAL REVIEW E 89, 053020 (2014)

Using the physical properties given in Table I, we obtain Re = 3 × 10−6 and 4 × 10−6 for systems A and B, respectively. Therefore, the contribution of the nonlinear inertia term  u) in Eq. (8) is negligible, as compared to the other ( u·∇ effects. The simplified form of Eq. (8), thus becomes Oh−2

∂ u  + χ ∇ 2 u. = −∇p ∂t

(11)

The same kind of asymptotic expansion can be taken for p as well. Next, if we substitute the above kind of asymptotic expansion for all dependent variables in Eq. (10), we get ζ

∂ u(1) ∂ u(2) ∂ u(0) + ζ2 + ζ3 ∂t ∂t ∂t  (1) − ζ 2 ∇p  (2)  (0) − ζ ∇p = −∇p + χ ∇ 2 u(0) + ζ χ ∇ 2 u(1) + ζ 2 χ ∇ 2 u(2) + O(ζ 3 ).

(f1) u1r (0,θ ) and u1θ (0,θ ) should be bounded (f2) u1θ (1,θ ) = u2θ (1,θ ) ,  1 dξ  (f3) u1r (1,θ ) = u2r (1,θ ) = , Ca dt r=1 (16)  H   E (f4) Trθ  + Trθ  = 0,  H   E κ (f5) −p + Trr  + Trr  = , and Ca (f6) u2r (1/α,θ ) = u2θ (1/α,θ ) = 0.

(10)

Since in most scenarios the momentum diffusion time scale is very small as compared to the capillary time scale, one may consider Oh to be large. In fact, considering different practical cases, Ha and Yang [43] have discussed seven different systems of the droplet and suspending fluid, from which one can obtain Oh ranging from 40 to 130, implying Oh 1. Accordingly, here we are seeking a solution for large Oh or small ζ , where ζ ≡ Oh−2 . To get an approximate solution of Eq. (10), u can be expressed for fe = 0, ζ  1 under the framework of the regular perturbation approach [44] as u = u(0) + ζ u(1) + ζ 2 u(2) + O(ζ 3 ).

problem can be given by

On arriving at the above boundary conditions we consider that for small drop deformation the normal and tangent directions on a drop surface can be given by the r and θ directions, respectively, at the leading order of approximation. Here uir and uiθ represents the radial and the angular components of the velocity field in the ith phase. TrrE and TrθE denote the normal and the shear stresses, respectively, due to the electric field. TrrH and TrθH denote the hydrodynamic normal and shear stresses, respectively. Note that, for imposing the continuity of radial components of velocities [boundary condition (f3)], it is convenient to introduce the drop shape function ξ for the deformed droplet interface. The same is realized for defining the corresponding interfacial curvature κ, while imposing the normal stress balances [boundary condition (f5)]. Details regarding the drop shape consideration and subsequent solution of electrohydrodynamic equations are delineated post priori.

(12) B. Solution procedure

From the above equation, we can obtain the leading order approximation of Eq. (10), pertinent for the present study, as  (0) = χ ∇ 2 u(0) . ∇p

(13)

For small ζ or large Oh, we are considering only the zeroth order solution. For further analysis, we are omitting the subscript zero from the velocity and pressure terms of Eq. (13). In general, a velocity field can be represented by [44] u =  where  is a differentiable scalar function and  +∇  × A, ∇  · A = 0) vector function. A is a differentiable solenoidal (∇  · u = 0) gives ∇ 2  = 0. The incompressibility condition (∇ Putting the above stated form of u in Eq. (13) and using the incompressibility condition, one can show [44] ∇ 4 A = 0.

(14)

In terms of the spherical polar coordinates, one may write [44] A = (r sin θ )−1 ψ(r,θ )iˆϕ , where iˆϕ represents the unit vector in the azimuthal direction (ϕ). Substituting A in Eq. (14) one can obtain [44] 4

ψ = 0,

where ψ is the stream function and the operator 2

1. Distribution of electric potential, surface free charge density, and electrical stresses in a confined domain

To obtain the distribution of the electric potential, surface free charge density, and electrical stresses, we begin with the solution of the electric potential distribution. The general solution of Eq. (1) in axisymmetric spherical coordinates can be given by

∞  bn n an r + n+1 Pn (cos θ ) φ1 (r,θ ) = r n=0

∞  dn n cn r + n+1 Pn (cos θ ) , φ2 (r,θ ) = r n=0

where Pn (cos θ ) is the nth order Legendre polynomial. Imposing the first three boundary conditions [(e1) to (e3)] we obtain

(15)

2

∞ 

dn (1 − α 2n+1 )r n Pn (cos θ )

n=2

φ2 (r,θ) = [−(1 + d1 α 3 )r + d1 r −2 ]P1 is

≡ ∂ /∂r + (sin θ /r )[∂/∂θ {(1/sin θ )∂/∂θ}]. The given by boundary and interfacial conditions for the hydrodynamic 2

φ1 (r,θ) = [−r + d1 (1 − α 3 )r]P1 +

2

(17)

2

053020-5

+

∞  n=2

(18)

(−dn α 2n+1 r n + dn r −n−1 )Pn (cos θ) .

MANDAL, CHAUDHURY, AND CHAKRABORTY

PHYSICAL REVIEW E 89, 053020 (2014)

Imposing these forms of φ1 and φ2 into condition (e4) of Eq. (2) and invoking the surface charge conservation Eq. (6), the equations for the rest of the unknown constants can be obtained as d Sa [S + 2 − α 3 (S − 1)] (d1 ) = (R − 1) − [R + 2 − α 3 (R − 1)]d1 , and 2 dt Oh (19) d Sa 2n+1 2n+1 2n+1 2n+1 [nS(1 − α ) + (nα + n − 1)] (dn ) = −[nR(1 − α ) + (nα + n + 1)]dn , ∀ n  2. dt Oh2 Solving Eq. (19) with the initial condition q(t = 0,θ ) = 0, (implying an initially charge-free droplet interface), we obtain

R−1 3R S (S − R) −t /τ1 R + , dn = 0, ∀ n  2, (20) d1 = e (R + 2) (S + 2) R+2 where R = (R + 2)/[(R + 2) − α 3 (R − 1)] and S = (S + 2)/[(S + 2) − α 3 (S − 1)] are the correction factors due to domain confinement. Here τ1 is the nondimensional charge relaxation time scale of the form

 (S + 2) − α 3 (S − 1) Sa . (21) τ1 = (R + 2) − α 3 (R − 1) Oh2 Now, the closed form expressions for the distribution of the electric potential can be obtained as  S (S − R) (1 − α 3 ) −t/τ1 3R r cos θ 1− e φ1 (r,θ) = − R+2 S+2 



R−1 1 1 3S (S − R) 3 −t /τ1 φ2 (r,θ) = −R cos θ r − rα − 2 e + . (R + 2) (S + 2) R + 2 r2 r

(22)

 assuming the applied electric field is For obtaining the distribution of the electric field, we invoke the criterion E = −∇φ, irrotational. With this argument, the electric field inside and outside the droplet can be obtained as 

 S (S − R) 1 − α 3 −t τ1 3R 1− E1 = e / (cos θ iˆr − sin θ iˆθ ) R+2 S+2 

2 (R − 1) 1 3S (S − R) 2 3 −t /τ1 ˆ  (23) ir E2 = R cos θ 1 + + α + 3 e (R + 2) (S + 2) R + 2 r3 r 

R−1 1 3S (S − R) 1 3 −t /τ1 ˆ iθ , − R sin θ 1 − + α − 3 e (R + 2) (S + 2) R + 2 r3 r where iˆr and iˆθ are the unit vectors in the r and θ directions, respectively. Using the electric field distribution, we obtain the surface charge distribution q(t,θ ) across the interface as

3R (R − S)  ˆ ⇒ q(t,θ ) = cos θ [1 − et/τ1 ]. (24) q (t,θ ) = −Q∇φ · n (R + 2) With the knowledge of the electric potential and surface charge distribution, it is now imperative to unveil the jump in electrical stresses across the interface. It is the stress jump that eventually contributes to the electrodeformation of the droplet, and thereby drives the subsequent events. The jumps in the normal and the tangential electrical stresses, across the interface, are [35,42]     jump in normal stress : TrrE  = 12 Q Er2 − Eθ2  and (25)   jump in shear stress : TrθE  = qE1θ . Using Eq. (23) we can obtain the closed form expressions for the jumps in normal and tangential stresses across the drop interface, as   E  9R2 cos2 θ S (S − R) T  = (R 2 + 1 − 2S) + {2R(2 + α 3 ) + 2(α 3 − 1) + 4S(1 − α 3 )}e−t /τ1 rr (S + 2) 2 (R + 2)2  2 (S − R)2 3 2 3 2 3 2 −2t /τ1 + S {(2 + α ) + (α − 1) − 2S(1 − α ) }e (S + 2)2   9R2 (S − 1) 9R2 S (S − R) (α 3 − 1) (S + 1) −t /τ1 9R2 S2 (S − R)2 (α 3 − 1)2 (S + 1) −2t /τ1 + + + e (26) e (R + 2)2 (S + 2) (R + 2)2 (S + 2)2 2 (R + 2)2 and

    E  9R2 (S − R) sin 2θ S (S − R) (1 − α 3 ) −t /τ1 S (S − R) (1 − α 3 ) −2t/τ1 T  = 1− 1+ . e e + rθ S+2 S+2 2 (R + 2)2 053020-6

(27)

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . .

PHYSICAL REVIEW E 89, 053020 (2014)

2. Drop shape function: Small deformation approximation

For small deformation (|D|  1), we are considering the drop shape function ξ that has been used previously by [11,13,15,20,32] of the following form  2D ξ = 1+ (3 cos2 θ − 1) . (28) 3 For approximating the curvature across the interface, =r −ξ ≡ we can define the modified shape function as 0, and thereby define the local curvature of the interface  · (∇  /|∇  |). Following the small deformation as κ = ∇ approximation (|D|  1) one can conduct a formal domain perturbation analysis using D as the perturbation parameter to approximate the curvature along with the use of Eq. (28). To the leading order approximation, the curvature approximates the form [11,13,15,20,32]  8D 2 (3 cos θ − 1) . κ = 2+ (29) 3

Ca (C1 + D1 ) =

The flow field can be solved using Eq. (15) in conjunction with the boundary and interfacial conditions given in Eq. (16). Boundary condition (f4) indicates the form of the stream function as ψ(r,θ) = r n sin2 θ cos θ . Using this form we can determine the possible values of n = 0, −2, 3, and 5. Now, we can write the stream functions inside and outside of the droplet as [11,20]

ψ2 (r,θ) = (A2 + B2 r

+ C2 r + D2 r ) sin θ cos θ, 5

2

(31)

respectively, where A1 − D2 are eight constants to be determined using the above mentioned boundary conditions. For boundary condition (f1) to be true we must have A1 = B1 = 0 and the stream function inside the droplet gets simplified to ψ1 (r,θ) = (C1 r 3 + D1 r 5 ) sin2 θ cos θ.

−2B2 α 7 + 3C2 α 2 + 5D2 = 0.

(37)

1 ∂ψ1 = [C1 r + D1 r 3 ](3 cos2 θ − 1) r 2 sin θ ∂θ (38)

and u1θ (r,θ ) = −

1 ∂ψ1 = −[3C1 r 2 + 5D1 r 4 ] sin θ cos θ. r sin θ ∂r (39)

Similarly, the velocity field outside the drop can be given 1 ∂ψ2 r 2 sin θ ∂θ = [A2 r −2 + B2 r −4 + C2 r + D2 r 3 ](3 cos2 θ − 1)

u2r (r,θ) =

(40) and 1 ∂ψ2 r sin θ ∂r = −[−2B2 r −4 + 3C2 r + 5D2 r 3 ] sin θ cos θ. (41)

u2θ (r,θ ) = −

Using the above velocity fields, TrθH  and TrrH  can be obtained as  H T  = [(−6A2 − 16B2 − 6C2 − 16D2 ) rθ

+ λ (6C1 + 16D1 )] sin θ cos θ

3C1 + 5D1 = −2B2 + 3C2 + 5D2 ,

(33)

2 dD , 3 dt

(34)

(42)

and  H T  = [(4A2 + 8B2 − 2C2 − 6D2 ) + λ (2C1 + 6D1 )] rr × (1 − 3 cos2 θ).

(32)

Using boundary conditions (f2), (f3), and, (f6), we get the following criteria for the remaining constants:

Ca(A2 + B2 + C2 + D2 ) =

(36)

u1r (r,θ) =

(30)

and 3

A2 α 5 + B2 α 7 + C2 α 2 + D2 = 0,

as

3. Distribution of the flow field

−2

(35)

To impose the rest of the boundary conditions (f4) and (f5), first we have to determine TrθH , TrrH , and p. Before determining the above mentioned expressions, we need to obtain the expression of the velocity field in and out of the droplet using stream functions. The velocity field inside the drop can be given as

Here ξ and κ are represented in their respective nondimensional forms.

ψ1 (r,θ) = (A1 + B1 r −2 + C1 r 3 + D1 r 5 ) sin2 θ cos θ

2 dD , 3 dt

(43)

To determine p, we first substitute the velocity field into Eq. (13) and then after integrating it we obtain the pressure field inside and outside of the droplet. Following this approach, the jump in pressure across the interface can be obtained as p = [(2A2 + 7D2 ) − 7λD1 ](3 cos2 θ − 1) + X,

(44)

where X is an integration constant. Now, we can use the boundary conditions (f4) and (f5) to obtain the following relations:

[6A2 + 16B2 + 6C2 + 16D2 ] − λ [6C1 + 16D1 ]    9R2 (S − R) S (S − R) (1 − α 3 ) −t /τ1 S (S − R) (1 − α 3 ) −2t /τ1 + − 1− 1+ = 0, e e 2(R + 2)2 S+2 S+2 053020-7

(45)

MANDAL, CHAUDHURY, AND CHAKRABORTY

and

PHYSICAL REVIEW E 89, 053020 (2014)

 9R2 S (S − R) [−18A2 − 24B2 + 6C2 − 3D2 ] + λ[−6C1 + 3D1 ] + (R 2 + 1 − 2S) + {2R(2 + α 3 ) + 2(α 3 − 1) 2 2(R + 2) (S + 2)  2 (S − R)2 8D 3 2 3 2 3 2 −2t/τ1 . (46) {(2 + α ) + (α − 1) − 2S(1 − α ) }e + 4S(1 − α 3 )}e−t/τ1 + S = (S + 2)2 Ca

In the normal stress balance equation [Eq. (46)], we have only retained the terms having coefficient cos2 θ because of the fact that the normal stress balance is valid at the interface for each values of θ . 4. Shape evolution of the droplet

Using Eqs. (33) to (37) and (45) we obtain C1 , D1 , A2 , B2 , C2 , and D2 in terms of D and dD/dt. Substituting those expressions into Eq. (46) we finally obtain the equation of the evolution of the drop deformation parameter (D) as dD = g1 e−t /τ1 + g2 e−2t /τ1 + g3 D + g4 , (47) dt where g1 − g4 , are functions of R,S,λ, and α (the complete expressions are given in the Appendix). Imposing a spherical droplet criterion (D = 0, at t = 0) as the initial condition, we obtain the solution of Eq. (47), which gives the expression for the temporal evolution of the droplet deformation parameter, hence the shape of the droplet, as D(t) = Dss [1 − 1 e−t/τ2 − 2 e−2t/τ1 − 3 e−t/τ1 ]. (48) Here Dss is the steady state droplet deformation parameter in a confined domain, given as Dss =

9Ca R2 , 16 (R + 2)2

(49)

where (R,S,λ,α) is the discriminating function which decides the sense of the drop deformation, and is given as  = (R 2 − 2S + 1) + 35 (R − S) ,

(50)

where  = [3λf1 (α) + 2f2 (α)]/[λf1 (α) + f3 (α)] is a combined function of the viscosity ratio and confinement ratio, aiding the correction to the Taylor discriminating function. In addition to the charge-relaxation time scale τ1 , as defined in Eq. (21), now we have another relaxation time scale τ2 =

19λ2 f1 (α) + λf4 (α) + f5 (α) , 20λf1 (α) + 20f3 (α)

(51)

which, in tandem, governs the shape relaxation of the droplet, as is evident from Eq. (48). The expressions of f1 − f5 , and 1 − 3 are given in the Appendix; nevertheless, they are functions of α, λ, R, and S, which are considered as the input parameters. It is worth mentioning that τ2 is intrinsic to the hydrodynamic situation under consideration. Even in the limit τ1 → 0 (the instantaneous charging of the droplet interface), τ2 alone governs the shape relaxation. Thus, without any loss of generality, we delineate τ2 as the intrinsic shape-relaxation time scale. Equations (47) and (48) are the pivotal equations for the analysis of the electrohydrodynamics of confined droplets. From this equation, the exponential nature of the temporal

evolution of the drop deformation is noteworthy. With this knowledge of D(t) and dD/dt, we can also evaluate the constants C1 , D1 , A2 , B2 , C2 , and D2 , leading towards obtaining the velocity field. Because of the very exponential nature of D(t), the velocity fields are also developed in a similar exponential fashion. With this theoretical framework in the background, we now endeavor in describing the transient dynamics of confined droplets in uniform electric fields. III. RESULTS AND DISCUSSION A. Comparison to previous works

From Eq. (48) it is evident that D/Dss follows an exponential behavior with time. In Fig. 2, we compare this notion to the previously reported theoretical works [11,32]. Esmaeeli and Sharifi [11] have obtained the droplet deformation parameter as a function of time in an unbounded domain of the form 0 D0 = Dss [1 − e−t/τ ], 0

0 Dss

(52)

where represents the steady state droplet deformation and τ 0 represents the time scale that governs the shape relaxation of the droplet for an unbounded domain (α = 0). Equation (52) clearly depicts that in an unbounded domain the droplet relaxes monotonically to its final steady state shape. In arriving at this notion, they have made a critical assumption: charging of the interface occurs instantaneously. This approximation is valid when the charge-relaxation time scale is very small as compared to the capillary time scale. In terms of the dimensionless number this situation is denoted by Sa/Oh2 → 0. However, for many of the physical situations [43] Sa/Oh2 is not small, even the charge-relaxation time scale can be larger than the capillary time scale, indicating that the charging process of the droplet interface takes place over a finite time span. This issue was specifically addressed by Lanauze et al. [32]. If we set Sa/Oh2 = 0 or more specifically the dimensionless charge-relaxation time scale τ1 = 0 and α = 0 in Eq. (48), we recover Eq. (52), as depicted in Fig. 2. Note that under this condition, the depiction by Esmaeeli and Sharifi [11] agrees well with the notion of Lanauze et al. [32]. For finite Sa/Oh2 (=0.2873), however, in addition to the lag in the deformation behavior of the droplet, a transition from prolate to oblate shape of the droplet is noteworthy, as is evident from the sign change in the D/Dss pattern in Fig. 2(a). It is worth mentioning that this is in good agreement with the prediction by Lanauze et al. [32]. To rationalize the genesis of the lag and nonmonotony in the temporal evolution of the drop deformation parameter, it needs to be emphasized that initially the interface of the droplet was charge free, which implies that there is only polarization stress acting on the interface that leads towards the prolate shape of the droplet. In due course, charges reach the interface and in later time the dynamics

053020-8

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . .

PHYSICAL REVIEW E 89, 053020 (2014)

(a)

ss

(b)

Esmaeeli and Sharifi [11] Lanauze et al. [32]

Esmaeeli and Sharifi [11]

Present model Sa/Oh2 = 0 Sa/Oh2 = 0.2873

Present model Sa/Oh2 = 0 Sa/Oh2 = 13.6

t FIG. 2. (Color online) Temporal evolution of normalized droplet deformation parameter (D/Dss ) in unbounded domain for (a) system A and (b) system B given in Table I.

is governed by the combined influence of the polarization stress and Coulomb force, leading towards the shape transition of the droplet, whereas the finite time taken to completely charge the interface is responsible for the lag in relaxing the shape of the droplet towards the final steady state. Both these behaviors would be more noticeable for higher values of Sa/Oh2 , as is shown in Fig. 2(b). Here, we have considered a situation with Sa/Oh2 = 13.6 which endorses a noticeable amount of nonmonotony in the shape relaxation behavior (reminiscent of the shape transition), along with a significant amount of lag in attaining the steady state shape. It is important to mention in context that such behavior is in line with the findings of recently reported studies [32]. Next, we compare our analytical results to the numerical simulations of Haywood et al. [12], as demonstrated in Fig. 3. In their simulation, Haywood et al. [12] have considered the behavior of liquid droplets suspended into other immiscible liquids and subjected to a uniform electric field. Though they have traversed a wide range of other parameters, the interfacial surface charging is considered to be instantaneous (Sa/Oh2 → 0). From Fig. 3(a) it is apparent that the present estimation matches very well for Oh = 200. Even for a smaller Oh(=20) we obtain a good agreement between the present theory and their simulation. Therefore, with the leading order

approximation of Eq. (10), the present theoretical approach seems to be satisfactory. Moreover, from Fig. 3(b) one can appreciate that even for large variation of λ, over three orders from 10−1 to 102 , the present theoretical framework is sufficient to capture the essential physics of interest, without sacrificing the analytical tractability. B. Steady state droplet deformation

Figure 4 demonstrates the variation of the steady state value of deformation parameter Dss with capillary number Ca, for a different extent of the confinement ratio (α), using the physical properties specified in the figure caption. For the sake of comprehensiveness, in the same figure, we also plot the theoretical prediction of Taylor [15] and the experimental results of Salipante and Vlahovska [36]. The capillary number Ca plays a very significant role in determining the drop dynamics in an electric field [10–12,14–18,20,25,36–41]. In tune with Eq. (49), it is obvious that the Dss − Ca relationship is linear in nature, as is also noteworthy from the linear trends in Fig. 4. Since the analytical expression given by Eq. (49) recovers the same steady state deformation value as obtained from the Taylor theory [15] under an unbounded condition, variations of Dss with Ca are the same for the Taylor theory [15] and Eq. (49) with α = 0. Additionally, the prediction

(a)

Simulation results of Haywood et al. [12] for different Ohnesorge number Oh = 200 Oh = 20 Oh = 2 ss

(b) Simulation results of Haywood et al. [12] for different Ohnesorge number λ = 100 λ =1 λ = 0.1

Present model

Present model

t τ2 FIG. 3. (Color online) Temporal variation of normalized droplet deformation parameter (D/Dss ) for different (a) Ohnesorge number Oh and (b) viscosity ratio λ; the properties are taken from Ref. [12]. 053020-9

MANDAL, CHAUDHURY, AND CHAKRABORTY

PHYSICAL REVIEW E 89, 053020 (2014)

that this kind of shape reversal phenomenon, due to domain confinement, has already been reported elsewhere [20].

Present model α =0 α = 0.5 α = 0.8

C. Transient droplet deformation ss

Experimental results of Salipante and Vlahovska [36]

Taylor Theory [15]

Ca

FIG. 4. (Color online) Variation of steady state deformation parameter (Dss ) with Capillary number (Ca) for different confinement ratio (α); the properties are taken from Ref. [36].

from the present theoretical model is in satisfactory agreement with the experimental results of Salipante and Vlahovska [36], as shown in Fig. 4. The confinement ratio α, on the other hand, modulates Dss in a rather complicated manner, as is evident from the contribution of R and  in Eq. (49). In general, the electrohydrodynamics of a suspended drop, in an otherwise quiescent medium, is governed by the ratio of the electrical conductivity, permittivity, and viscosity of the drop to that of the suspending liquid. Under the condition R  S, the drop assumes a prolate shape, defined by the orientation of its major axis parallel to the direction of the applied electric field. For R < S the drop may assume either an oblate or prolate shape, characterized by Dss < 0 and Dss > 0, respectively [10–12,14–18,20,25,36–41]. These features are typical when the drop is in an unbounded domain (confinement ratio α = 0) [11,12]. Under the strong influence of confinement (0 < α < 1), a remarkable deviation from the above mentioned attributes can be observed, as is demonstrated in Fig. 4. Specifically, in a highly confined domain, the droplet can assume a prolate shape, although the similar parametric settings predicts an oblate shaped droplet in an unbound domain, as can be appreciated from Fig. 4, by comparing the trend of α = 0.8 with others. It is important to mention in this context

Typical patterns of the temporal evolution of a droplet in an electric field are shown in Fig. 5, for different extents of the confinement ratio α, using the physical properties of system B specified in Table I. We have already shown the influence of the confinement in deciding the steady state shape of the droplet. Here we primarily focus on unveiling the underlying features of the temporal evolution of the droplet, leading towards such a steady state configuration. The temporal evolution of the deformation parameter follows an exponential pattern, as depicted in Eq. (48) and the same is demonstrated in Fig. 5. Using the property values of system B in Table I, the droplet temporally evolves towards the oblate configuration, for an unbound condition (α = 0), as depicted in Fig. 5(a). By progressively bringing the confining walls towards the drop surface, the similar evolution trend towards the oblate configuration persists. However, this occurs with reduced magnitude, as is shown in Fig. 5(a), for α = 0,0.3,0.4,0.5. In this way, one can approach a critical confinement ratio (αcr ≈ 0.5529, as we have found, for the property values chosen), for which the magnitude of the drop deformation parameter is zero. With further reduction of the space between the drop surface and the confining wall, the temporal evolution of the droplet leads towards attaining prolate configurations at the steady state. Figure 5(b) demonstrates those D − t trends, leading towards the prolate-shaped configuration, of the droplet under the strong influence of domain confinement (α = 0.6,0.7,0.8,0.9). However, it needs to be emphasized at this juncture that a close scrutiny of Fig. 5 reveals that the deformation at early times appears to be somewhat independent of the confinement ratio. We must appreciate that for the confinement effect to be realizable, an interaction between the incipient flow field with the confining wall needs to be taking place. This condition could be realized when the interface gets sufficiently charged, which essentially leads towards the development of the intrinsic flow field due to the considerable amount of mismatch of the tangential electrical stresses across the drop surface. The development of these

(a)

(b)

α = 0.6 α = 0.7 α = 0.8 α = 0.9

α =0 α = 0.3 α = 0.4 α = 0.5

t FIG. 5. (Color online) Temporal evolution of the deformation parameter (D) for different confinement ratio (α) using properties of system B from Table I. 053020-10

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . .

PHYSICAL REVIEW E 89, 053020 (2014)

(a)

(b)

α = 0.6 α = 0.7 α = 0.8 α = 0.9

α =0 α = 0.3 α = 0.4 α = 0.5 ss

t FIG. 6. (Color online) Temporal evolution of normalized deformation parameter (D/Dss ) for different confinement ratio (α) using properties of system B from Table I.

events, on the other hand, takes a certain time to occur, owing to the finiteness of the charge-relaxation time scale. For understanding the temporal relaxation of the drop shape, in Fig. 5 we primarily focus on the D − t patterns and the influence of domain confinement in deciding their approach towards the steady state configuration. However, those temporal relaxations exhibit certain subtle yet interesting additional features, as depicted through the D/Dss − t characteristics in Fig. 6, for the similar parametric settings used for constructing Fig. 5. The prominent fact that one can recognize is that there may exits nonmonotony in the temporal relaxation behavior of the droplet. This is reminiscent of the shape transition during temporal evolution, in addition to the shape transition in the final steady state configuration, induced by the domain confinement. A close scrutiny of Eq. (48), used for constructing Fig. 6, reveals that the temporal relaxation of the droplet actually follows a double exponential pattern in contrast to the single exponential pattern obtained by Esmaeeli and Sharifi [compare Eqs. (48) to (52)]. Without any loss of generality, thus, we can proclaim that it is this double exponential nature that endorses such nonmonotony, as shown in Fig. 6. It is now imperative to unveil the underlying physics leading towards such behaviors. From Eq. (48), one can appreciate that the double exponential nature of the D/Dss − t pattern is primarily induced by the involvement of the exp(−t/τ1 ) terms in addition to the exp(−t/τ2 ) term. While the contribution of the latter is dictated by the viscosity contrast and the domain confinement [consider the form of τ2 in Eq. (51)], the contribution of the former, on the other hand, primarily stems from the consideration of the finiteness of the charge relaxation time. For the infinitesimal charge relaxation time scale (τ1 → 0), one can easily recover the notion of Esmaeeli and Sharifi [Eq. (52)], for an unbound domain (α = 0). It is worth mentioning that some contemporary studies [14,32,35] have emphasized the contribution of a finite charge-relaxation time scale on such a nonmonotonic shape relaxation of a droplet in an electric field. We have already discussed, in the mentioned studies, the genesis of such behavior: initially the droplet starts deforming due to the polarization stress owing to the nonexistence of charges at the interface; later on when the interface gets charged, the deformation dynamics is governed by the combined confluence of the polarization stress and the

Coulombic force. Besides agreeing with their notion, here we have shown the influence of confinement in deciding the drop shape relaxation behavior, as is shown in Fig. 6. Essentially the hydrodynamic interaction between the intrinsic flow field with the confining wall becomes important. From Fig. 6(a), it is evident that for α < αcr the tendency to exhibit nonmonotony enhances with confinement, however, without any noticeable alteration in the time to reach the steady state. For α > αcr , on the other hand, the influence of confinement is more prominent in altering the time to reach the steady state configuration, in comparison to endorsing nonmonotony, as can be appreciated from Fig. 6(b). From Fig. 6(b), additionally, we note that the D/Dss − t pattern for α = 0.6 behaves somewhat differently than the other D/Dss − t patterns in the α > αcr regime. This can be attributed to the combined influence of the finite charge relaxation time scale and the confinement effect in determining the rate at which the droplet gets deformed. However, at higher values of the confinement ratio, the deformation rate is primarily determined by the confinement effect. For a critical assessment of the underlying physics, it is imperative to look into the behavior of the relaxation time scales τ1 and τ2 under different situations, as is endeavored in the subsequent section. D. Relaxation time scale

It has already been emphasized in the present article so far that the transient dynamics is governed by the time scales τ1 [represented by Eq. (21)] and τ2 [represented by Eq. (51)]. While τ1 is reminiscent of the consideration of the finite charge relaxation time scale, τ2 is intrinsic to the present situation under consideration. Even if the charge accumulation process on the drop surface occurs over an infinitesimal span of time (i.e., τ1 → 0), the shape evolution is governed by the intrinsic time scale τ2 . It is worth mentioning that for an unbound domain [imposing α → 0 in Eq. (51)], τ2 depends solely on λ, in the form of a rational polynomial (quadratic in the numerator and linear in the denominator). This is in line with the existing notions [11]. However, if confinement effects are introduced the relationship between τ2 , α, and λ may not be a trivial one, as can be appreciated form Eq. (51) and its supporting documentation in the Appendix. To illustrate this, in Fig. 7(a), we plot τ2 as a function of α for different λ. From Fig. 7(a) it is evident that for higher values of λ, there occurs a marked rise

053020-11

MANDAL, CHAUDHURY, AND CHAKRABORTY

PHYSICAL REVIEW E 89, 053020 (2014)

(a)

(b)

τ 2 τ1

τ2 λ = 10 λ =1 λ = 0.1

α

α

FIG. 7. (Color online) (a) Variation of τ2 with confinement ratio (α) for different viscosity ratio (λ) and (b) variation of τ2 /τ1 with confinement ratio (α) using properties of system B from Table I.

in τ2 . Another interesting fact represented by Fig. 7(a) is that for the increase in drop viscosity compared to the background liquid (higher values of λ), τ2 increases. Though this increment is significant for low to moderate confinement ratios, for high confinement ratios (α > 0.8) there is not much difference in relaxation time due to the change in λ, even when λ is varied over three orders of magnitude. Therefore, we can proclaim that for strong confinements, the intrinsic shape-relaxation time scale τ2 is mostly governed by the confinement effect. In Fig. 7(b), we also compare the relative importance of τ2 and τ1 , as eventually it is their combined influence that determines the shape relaxation of the droplet. It is evident from Fig. 7(b) that at the lower extent of confinement τ2 < τ1 . However, when the wall of the confining domain is brought very near to the drop surface, an enormous rise in τ2 well surpasses the effect of τ1 . The fact τ2 τ1 , in the limit α → 1, shows that at the higher extent of confinement, the shape relaxation is primarily governed by τ2 and its higher magnitude signifies the sluggishness of the droplet to attain the steady state configuration.

where the forms of τ2 and Dss remain unchanged, as shown in Eqs. (51) and (49), respectively. Now, starting from Eq. (51), bearing the definition of τ2 , it is straightforward to calculate (19λ + 16) (2λ + 3) τ2 (λ,α = 0) = . (54) 40 (λ + 1) Equation (54) is in line with the finding of Esmaeeli and Sharifi [11] and follows a similar trend with that reported in Haywood et al. [12] for an unbounded domain. In an unbounded domain, accordingly, the viscous effect alone determines the drop shape-relaxation behavior. In fact, at a higher value of λ, the drop relaxes more sluggishly. Interestingly, we have also noted that the confinement effect may hold the ability to enhance this sluggishness [note the remarkable rise in τ2 at higher values of α in Fig. 7(a)]. Based on this notion, we can say that confinement has an ability to produce additional viscous damping. In the above perspective, it would be interesting to draw an equivalence between the effective viscous effect in an

E. Intrinsic shape-relaxation time scale τ2 and the concept of equivalent viscosity

Confined α = 0.5

Although we have discussed the relaxation time scales τ1 and τ2 , it is nevertheless, important to appreciate that the genesis of τ2 is intrinsic to the present paradigm under consideration. If we consider that the time scale of the charge relaxation is infinitesimal (τ1 → 0), the drop shape dynamics is then solely governed by the time scale τ2 . It is, therefore, imperative to endeavor a critical scrutiny of τ2 in deciding the temporal behavior of the droplet. Towards this end, we consider the situations with infinitesimal charge-relaxation time scales (τ1 → 0), so as to eliminate any additional contribution of the charge accumulation process over the drop surface. The subsequent discussions hereinafter will be made in reference to the condition τ1 → 0, unless otherwise specified. With this notion, the description of the shape relaxation [Eq. (48)], assumes the form lim D(t) = Dss [1 − e−t/τ2 ],

τ1 →0

(53)

Unbounded α

=0

ss

λeq t

λ = 10 λ =1

λ = 0.1

α FIG. 8. (Color online) Representation of equivalent viscosity ratio (λeq ) as a function confinement ratio (α), as obtained from Eq. (55) for three different droplet viscosities in a confined domain. The inset shows a comparison of the equivalence, following Eq. (53).

053020-12

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . .

unconfined domain with the combined effect of viscosity and confinement in an bounded domain to produce identical drop relaxation characteristics. For this sake, we first consider Eq. (54) that is applicable for an unbounded drop. Then, we consider the case of a confined droplet with viscosity ratio λ by invoking the definition given in Eq. (51). Considering the equivalence between Eqs. (54) and (51), we obtain an equivalent viscosity ratio in a hypothetical unbounded situation that mimics the confinement effect under bounded conditions

λeq =



 −89 +

2

−26



+ 625

 76,

(55)

PHYSICAL REVIEW E 89, 053020 (2014)

 with = 2[19λ2 f1 (α) + λf4 (α) + f5 (α)]/[λf1 (α) + f3 (α)] (the coefficients f1 − f5 are given in the Appendix). If we know the extent of confinement (α), we can easily find out an equivalent viscosity ratio (λeq ), from Eq. (55), for which the drop shape-relaxation time scale (τ2 ) is the same. In Fig. 8 we represent the λeq − α relation as obtained from Eq. (55). From the figure we can see that at much lower confinement (for α ∼ 0.2 or less), λeq is not a strong function of α. The dependency is significant in moderate to higher confinement regimes (for α ∼ 0.3 to ∼0.9). To match the relaxation time scale at a much higher confinement ratio (for α ∼ 1), one has to increase λeq to a very large extent. From Fig. 8, therefore, one can appreciate a notable fact: confinement holds the ability to produce an

t

0

t=0.1

t=1.5

t=1

t=10

0.2

t=0.1

t=1

t=1.5

t=10

t=0.1

t=1

t=1.5

t=10

t=0.1

t=2.5

t=5

t=50

t=10

t=20

t=100

α 0.4

0.7

0.8

t=0.1

FIG. 9. (Color online) Evolution of streamline pattern in time for different confinement ratios (α), using the property value of system B from Table I, albeit, with consideration of instantaneous interface charging (Sa = 0). 053020-13

MANDAL, CHAUDHURY, AND CHAKRABORTY

PHYSICAL REVIEW E 89, 053020 (2014)

enormous effective viscous effect which eventually reflects in the alteration of the drop shape-relaxation time scale.

F. Distribution of the flow field

Besides the above mentioned considerations, the intrinsic hydrodynamic conditions also represent some remarkable features of the confinement effect, as demonstrated in Fig. 9. The figure is a representative case for the situations with infinitesimal charge relaxation time scale, i.e., τ1 → 0. Prior to scrutinizing the implication of the charge relaxation process, it is essential to unveil the influence of confinement on the intrinsic hydrodynamics. The figure clearly shows the confinement-induced alteration in the transient evolution of

streamline patterns. The enhanced crowding of streamlines, inside and outside of the droplet, is a notable feature of the confinement [20]. In a confined condition (α > 0), the relative distances between the streamlines are very small near the drop surface, signifying strong velocity gradients in the vicinity, though the same symmetric nature of the streamlines is retained in the confined domain. In Fig. 9, if we follow the α direction, we can see that along with enhanced crowding, streamlines are getting more closed. Moreover, if we follow the t direction, we can note that inside the drop, the streamlines are getting closed in due course. This tendency occurs at much earlier t, when α is high. This may be attributed to the fact that confinement has an effect of producing strong vortices and thereby produces crowded and closed streamlines. Such

t

0

t=1

t=0.1

t=0.5

t=0.1

t=0.5

t=1

t=0.1

t=0.5

t=1.1

t=0.1

t=1

t=0.1

t=10

t=10

0.2

t=10

α 0.4

t=10

0.7

t=1.5

t=10

0.8

t=35

t=100

FIG. 10. (Color online) Evolution of streamline pattern in time for different confinement ratios (α) considering finite charge relaxation time (Sa/Oh2 = 13.6) for system B in Table I. 053020-14

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . .

confinement-induced crowding and closing of streamlines can also be found out during the deformation of highly confined droplets in shear flows (see Ref. [45] and the references therein for details). With the aforementioned observations and the corresponding arguments in the background, we now demonstrate the influence of the charge-relaxation process on the intrinsic hydrodynamics of a confined droplet in an electric field, as is shown in Fig. 10. This figure is constructed for system B (see Table I), and corresponds to Sa/Oh2 = 13.6. By comparing Figs. 9 and 10 it is noteworthy that the effect of charge relaxation plays a dominant role when the confining wall is not very close to the drop surface. Initially the streamlines penetrate the drop surface, and the tendency is more for the case depicted by Fig. 10, as compared to the case represented by Fig. 9. In due course of time, streamlines outside and inside the droplet exhibit the tendency to form separated recirculation zones. If the charge relaxation process is considered to occur instantaneously, as in Fig. 9, the formation of separated recirculation zones is facilitated. In the postrelaxation regime, on the other hand, the recirculations primarily occur inside and outside the droplet almost separately; the penetration of streamlines is restricted within a very small region around the drop surface. It is, nevertheless, worth mentioning that for a highly confined domain, the temporal evolution of the streamline pattern is almost similar for the two cases represented by Figs 9 and 10 (compare the case α = 0.8 from Figs. 9 and 10). Thus, the charge relaxation process seems to contribute very little in the development of the inherent hydrodynamics in a highly confined environment. IV. CONCLUSION

In the present study, we provide an analytical estimation of the transient dynamics of confined liquid drops having finite charge relaxation time under the influence of a uniform dc electric field. The analysis brings out an explicit

PHYSICAL REVIEW E 89, 053020 (2014)

relationship between the flow field, the electric field, and the deformation time history as a function of the relative extent of confinement, viscosity, and electrical properties of the drop and the background liquid. The deformation parameter, in the limit t → ∞, is found to be the same as that obtained from steady state analysis of the electrohydrodynamics of confined drops [20]. Our results reveal that the droplet shape relaxation is governed by two important time scales, namely, the nondimensional charge-relaxation time scale and the intrinsic shape-relaxation time scale; domain confinement affects both of them. We uncover an explicit relationship for the intrinsic shape-relaxation time scale as a function of the relative extent of the confinement and viscosity ratio. The relationship shows that the combined effect of the confinement and viscosity holds the ability to elevate the intrinsic shape-relaxation time scale in a somewhat dramatic fashion. To realize the same intrinsic shape relaxation in an unconfined domain, one would, in effect, necessitate abnormally elevated viscosity ratios. The elevation in the intrinsic shape-relaxation time scale of a drop, in a highly confined environment, can even well surpass the contribution of the charge-relaxation time scale, although the interface charging process occurs over a considerable duration in real time. Thus, fluidic confinement provides an effective tool towards manipulating the time-relaxation properties of liquid drops under transient conditions, in line with the specific requirements commensurate with modern-day microfluidic applications, which may otherwise be realized in unbounded domains albeit with unphysically large viscosity ratios. APPENDIX

Equation (47) contain the functions g1 − g4 of the following form: g1n g2n g3n g4n g1 = , g2 = , g3 = and g4 = , (A1) g1d g2d g3d g4d where

g1n = Ca[(−(α − 1)2 (48α + 36λ + 72αλ + 108α 2 λ − 81α 3 λ − 270α 4 λ − 81α 5 λ + 108α 6 λ + 72α 7 λ + 36α 8 λ + 72α 2 − 204α 3 − 480α 4 − 294α 5 − 108α 6 − 72α 7 − 36α 8 + 24))p2 + (α + 1)2 (40α + 20λ + 40αλ + 60α 2 λ − 45α 3 λ − 150α 4 λ − 45α 5 λ + 60α 6 λ + 40α 7 λ + 20α 8 λ + 60α 2 + 30α 3 − 30α 5 − 60α 6 − 40α 7 − 20α 8 + 20)p5 ], g1d = (152α 10 − 950α 7 + 1596α 5 − 950α 3 + 152)λ2 + (−304α 10 + 50α 7 − 252α 5 + 150α 3 + 356)λ + 152α 10 + 900α 7 − 1344α 5 + 800α 3 + 192, g2n = [Ca((−36α 10 + 225α 7 − 378α 5 + 225α 3 − 36)λ + 36α 10 + 150α 7 − 462α 5 + 300α 3 − 24)p3 + ((20α 10 − 125α 7 + 210α 5 − 125α 3 + 20)λ − 20α 10 + 50α 7 − 50α 3 + 20)p6 ], g2d = [(152α 10 − 950α 7 + 1596α 5 − 950α 3 + 152)λ2 + (−304α 10 + 50α 7 − 252α 5 + 150α 3 + 356)λ + 152α 10 + 900α 7 − 1344α 5 + 800α 3 + 192], g3n = [(− 20(α 2 − 2α + 1)(4α 8 + 8α 7 + 12α 6 − 9α 5 − 30α 4 − 9α 3 + 12α 2 + 8α + 4))λ − 20(α 2 − 2α + 1)(−4α 8 − 8α 7 − 12α 6 − 6α 5 + 6α 3 + 12α 2 + 8α + 4)], g3d = [(76α 10 − 475α 7 + 798α 5 − 475α 3 + 76)λ2 + (−152α 10 + 25α 7 − 126α 5 + 75α 3 + 178)λ + 76α 10 + 450α 7 − 672α 5 + 400α 3 + 96], g4n = [Ca((−36α 10 + 225α 7 − 378α 5 + 225α 3 − 36)λ + 36α 10 + 150α 7 − 462α 5 + 300α 3 − 24)p1 + ((20α 10 − 125α 7 + 210α 5 − 125α 3 + 20)λ − 20α 10 + 50α 7 − 50α 3 + 20)p4 ], g4d = 2[(76α 10 − 475α 7 + 798α 5 − 475α 3 + 76)λ2 + (−152α 10 + 25α 7 − 126α 5 + 75α 3 + 178)λ + 76α 10 + 450α 7 − 672α 5 + 400α 3 + 96]. 053020-15

(A2)

MANDAL, CHAUDHURY, AND CHAKRABORTY

PHYSICAL REVIEW E 89, 053020 (2014)

The above expressions contains p1 − p6 of the following form: p1 =

9(S − R)L2R , 2(R + 2)2

 S (S − R)(1 − α 3 )  , p2 = −p1 1 + (S + 2)

 S (S − R)(1 − α 3 )  9R2 (R 2 − 2S + 1) , p4 = , (S + 2) 2(R + 2)2 9R2 S (S − R) [2R(2 + α 3 ) + 2(α 3 − 1) + 4S(1 − α 3 )2 ], and p5 = 2(R + 2)2 (S + 2) p3 = p1

p6 =

(A3)

9R2 S2 (S − R)2 [(2 + α 3 )2 + (α 3 − 1)2 − 2S(1 − α 3 )2 ]. (R + 2)2 (S + 2)2

Equations (48), (50), and (51) contain 1 − 3 and f1 − f5 of the following form: g2 2g2 g3 g1 τ1 + + , g4 g4 (g3 τ1 − 2) g4 (g3 τ1 − 1) g2 g3 τ1 , 2 = g4 (2 − g3 τ1 ) g3 g1 τ1 . 3 = g4 (1 − g3 τ1 ) 1 = 1 +

(A4)

and f1 = 4α 10 − 25α 7 + 42α 5 − 25α 3 + 4, f2 = −6α 10 − 25α 7 + 77α 5 − 50α 3 + 4, f3 = −4α 10 + 10α 7 − 10α 3 + 4, f4 = −152α f5 = 76α

10

10

(A5)

+ 25α − 126α + 75α + 178, 7

5

3

+ 450α 7 − 672α 5 + 400α 3 + 96.

[1] W. A. Macky, Proc. R. Soc. A 133, 565 (1931). [2] C. G. Garton and Z. Krasucki, Proc. R. Soc. A 280, 211 (1964). [3] K. J. Ptasinski and P. J. A. M. Kerkhof, Sep. Sci. Technol. 27, 995 (1992). [4] Y. Wu and R. L. Clark, J. Biomater. Sci., Polym. Ed. 19, 573 (2008). [5] O. A. Basaran, AIChE J. 48, 1842 (2002). [6] X. Chen, J. Cheng, and X. Yin, Chinese Science Bulletin 48, 1055 (2003). [7] A. Banerjee, E. Kreit, Y. Liu, J. Heikenfeld, and I. Papautsky, Lab on a Chip 12, 758 (2012). [8] K. Ikemoto, I. Sakata, and T. Sakai, Sci. Rep. 2, 289 (2012). [9] M. R. Abdelaal and M. A. Jog, Int. J. Heat Mass Transf. 55, 251 (2012). [10] A. Esmaeeli and P. Sharifi, J. Electrost. 69, 504 (2011). [11] A. Esmaeeli and P. Sharifi, Phys. Rev. E 84, 036308 (2011). [12] R. J. Haywood, M. Renksizbulut, and G. D. Raithby, AIChE J. 37, 1305 (1991). [13] O. Vizika and D. A. Saville, J. Fluid Mech. 239, 1 (1992). [14] J. Zhang, J. D. Zahn, and H. Lin, Phys. Rev. E 87, 043008 (2013). [15] G. Taylor, Proc. R. Soc. A 291, 159 (1966). [16] J. R. Melcher and G. I. Taylor, Annu. Rev. Fluid Mech. 1, 111 (1969). [17] D. A. Saville, Annu. Rev. Fluid Mech. 29, 27 (1997).

[18] S. Torza, R. G. Cox, and S. G. Mason, Proc. R. Soc. A 269, 295 (1971). [19] J. Hua, L. K. Lim, and C.-H. Wang, Phys. Fluids 20, 113302 (2008). [20] A. Esmaeeli and A. Behjatian, Phys. Rev. E 86, 036310 (2012). [21] J. Q. Feng, Proc. R. Soc. A 455, 2245 (1999). [22] J. Q. Feng and T. C. Scott, J. Fluid Mech. 311, 289 (2006). [23] J.-W. Ha and S.-M. Yang, J. Fluid Mech. 405, 131 (2000). [24] E. Lac and G. M. Homsy, J. Fluid Mech. 590, 239 (2007). [25] Y. Lin, P. Skjetne, and A. Carlson, Int. J. Multiphase Flow 45, 1 (2012). [26] J. M. L´opez-Herrera, S. Popinet, and M. A. Herrada, J. Comput. Phys. 230, 1939 (2011). [27] C. E. Rosenkilde, Proc. R. Soc. A 312, 473 (1969). [28] L. Capretto, W. Cheng, M. Hill, and X. Zhang, Top. Curr. Chem. 304, 27 (2011). [29] S. Moriya, K. Adachi, and T. Kotaka, Langmuir 2, 155 (1986). [30] T. Nishiwaki, K. Adachi, and T. Kotaka, Langmuir 4, 170 (1988). [31] C. Sozou, Proc. R. Soc. A 334, 343 (1973). [32] J. A. Lanauze, L. M. Walker, and A. S. Khair, Phys. Fluids 25, 112101 (2013).

053020-16

TRANSIENT DYNAMICS OF CONFINED LIQUID DROPS . . . [33] G. Supeene, C. R. Koch, and S. Bhattacharjee, J. Colloid Interface Sci. 318, 463 (2008). [34] A. Fern´andez, G. Tryggvason, J. Che, and S. L. Ceccio, Phys. Fluids 17, 093302 (2005). [35] X. Xu, Ph.D. thesis, University of California, Santa Barbara, 2007. [36] P. F. Salipante and P. M. Vlahovska, Phys. Fluids 22, 112110 (2010). [37] H. Sato, N. Kaji, T. Mochizuki, and Y. H. Mori, Phys. Fluids 18, 127101 (2006). [38] J. D. Wehking, L. Chew, and R. Kumar, Appl. Phys. Lett. 103, 054101 (2013).

PHYSICAL REVIEW E 89, 053020 (2014) [39] A. Behjatian and A. Esmaeeli, Int. J. Multiphase Flow 48, 71 (2013). [40] M. N. Reddy and A. Esmaeeli, Int. J. Multiphase Flow 35, 1051 (2009). [41] A. Esmaeeli and M. N. Reddy, Int. J. Multiphase Flow 37, 1331 (2011). [42] X. Xu and G. M. Homsy, J. Fluid Mech. 564, 395 (2006). [43] J.-W. Ha and S.-M. Yang, J. Rheol. 44, 235 (2000). [44] L. G. Leal, Advanced Transport Phenomena (Cambridge University Press, Cambridge, England, 2007). [45] V. Sibillo, G. Pasquariello, M. Simeone, V. Cristini, and S. Guido, Phys. Rev. Lett. 97, 054502 (2006).

053020-17