Modeling of conductive particle motion in viscous medium affected by an electric field considering particle-electrode interactions and microdischarge phenomenon

PHYSICS OF FLUIDS 28, 107102 (2016) Modeling of conductive particle motion in viscous medium affected by an electric field considering particle-electrode interactions and microdischarge phenomenon Ghiyam Eslami,1,a) Esmaeil Esmaeilzadeh,1,b) and Alberto T. Pérez2,c) 1Department of Mechanical Engineering, University of Tabriz, Tabriz, Iran 2Departamento de Electrónica y Electromagnetismo, Facultad de Física, Universidad de Sevilla, Avenida Reina Mercedes s/n, 41012 Sevilla, Spain (Received 11 February 2016; accepted 28 September 2016; published online 17 October 2016) Up and down motion of a spherical conductive particle in dielectric viscous fluid driven by a DC electric field between two parallel electrodes was investigated. A nonlinear differential equation, governing the particle dynamics, was derived, based on Newton’s second law of mechanics, and solved numerically. All the pertaining dimensionless groups were extracted. In contrast to similar previous works, hydrodynamic interaction between the particle and the electrodes, as well as image electric forces, has been taken into account. Furthermore, the influence of the microdischarge produced between the electrodes and the approaching particle on the particle dynamics has been included in the model. The model results were compared with experimental data available in the literature, as well as with some additional experimental data obtained through the present study showing very good agreement. The results indicate that the wall hydrodynamic effect and the dielectric liquid ionic conductivity are very dominant factors determining the particle trajectory. A lower bound is derived for the charge transferred to the particle while rebounding from an electrode. It is found that the time and length scales of the post-microdischarge motion of the particle can be as small as microsecond and micrometer, respectively. The model is able to predict the so called settling/dwelling time phenomenon for the first time. Published by AIP Publishing. [http://dx.doi.org/10.1063/1.4964683] I. INTRODUCTION Small conductive particles exposed to a DC electric field are encountered in various technical and industrial situations, such as electrostatic technologies, high voltage power systems, printing equipment, surface coating, separation of mixed granular materials,1and biological cell manipulation. A deep understanding of the particles dynamics is required for designing and maintaining such systems. Also, experimental progresses during past few years in the areas of a microelectromechanical systems (MEMS), nano-science, and microfluidics have intensified the interest in the control of motion and manipulation of particles by electrical forces.2Accurate prediction and control of conductive particles motion would improve the performance and efficiency of high voltage technologies. One of the common configurations involving tiny conductive particles under influence of a DC electrical field refers to the recycling system of plastic waste.3In this configuration, an electric force is applied to remove conductive impurities, like carbon particles, from liquefied plastic.4,5Also, it has been shown that the motion of a conductive particle in a typical parallel electrode system can be used for measuring the viscosity of fluids.6 a)Author to whom correspondence should be addressed. Electronic mail: [email protected]. Telephone: (+98) 3016044. b)Email: [email protected] c)Email: [email protected] 1070-6631/2016/28(10)/107102/29/$30.00 28, 107102-1 Published by AIP Publishing. 107102-2 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) Another interesting potential application of similar configurations is the heat transfer enhancement in fluid flow by the use of micron sized solid particles dispersed in the liquid bulk. Bologa et al.7studied the heat transfer enhancement caused by the electric field induced motion of solid particles in a gas-solid suspension. Added solid particles can improve the heat transfer rate between a fluid flow and a solid surface by advection, and disturbing and thinning the laminar sub-layer in the wall vicinity.8Dispersed particles would improve turbulence intensity of a flow up to 100% when their size is of the order of the turbulence length scale.9A theoretical study of local heat transfer enhancement through bombardment of a surface by spherical particles can be found in a paper by Murray.10 The electric force may help to overcome some difficulties faced by these technologies, like particle sedimentation and particle clogging. Augmentation of heat transfer rate in a dielectric fluid medium with the aim of electric force induced motion of added particles may become more important everywhere the gravitational force is weak, specifically in space transport systems. Also, heat transfer in colloidal suspensions with charged conductive particles is a field of interest.11 There are many papers on the subject of conductive particle motion induced by electric forces. Dascalescu and Mihailescu1proposed a mathematical model for motion of a conductive particle exposed to DC corona field of coaxial wire-cylinder electrode system. Asano et al.4conducted an analytical and experimental study on the acceleration and deceleration of a conductive particle in the DC electric field in a viscous fluid medium. They found a significant disagreement between their experimental data and the theoretical model prediction. They attributed this discrepancy to two factors. One is the hydrodynamic interaction between particle and electrode. The other is the fact that the particle obtains a charge smaller than expected after the particle-electrode contact.3,12 Recently, Jiang et al.13 introduced an improved model for predicting the trajectories of conductive particles in roll-type electrostatic separator. Dascalescu et al.14 conducted a numerical simulation of conductive particle behavior at the surface of a plate electrode affected by a DC corona field and employed the model to study the behavior of coal particles in fly-ash corona separators. Wang et al.15 fulfilled an experimental study about effects of conductive particles on the breakdown voltage of mineral oils and esters. Karunanayake and Hoshino16 accomplished a numerical simulation for predicting the electrostatic force acting on a conducting particle in the configurations of two parallel electrodes for the case in which the particle size is of the same order of the gap between electrodes. Besides the solid particle case, there is a large amount of studies carried out on the behavior of water droplet. Valuable experimental data are available from investigations of Jalaal et al.17 and Khorshidi et al.18 As is the case in the work of Asano,4two main simplifications are usually done by previous researchers. The first one concerns the electric force being exerted on the particle. It is often considered a constant electric force equal to E0Q0, where E0is the uniform electric field far from the particle and Q0is the charge acquired by the particle after its collision with the electrodes. The second simplification is to neglect the hydrodynamic interaction between the particle and the electrode surface and simply employ the well-known Stokes’ formula for the viscous drag (6πµRup). The main goal of this paper is to construct an accurate model for the quantitative prediction of the motion of a conductive particle in the gap between two parallel electrodes. We take into account all the pertinent effects, including the influence of electrode presence on both electrical and hydrodynamic forces, as well as the microdischarge that takes place between an approaching conductive particle and the counter electrode (the electrode with charges opposite to that of the particle). The electric force was calculated by an expression derived by Pérez19 in which the image forces (caused by image charges and dipoles) have been taken into account. Particle-wall hydrodynamic interaction was introduced to the problem using appropriate wall-effect functions to modify steady Stokes drag and the so called added mass forces acting on the particle. Also, the microdischarge effect on the particle motion has been included by accounting the changes of electric force (both in value and orientation) immediately after the phenomenon (before the particle-electrode collision). II. PROBLEM DESCRIPTION AND RELATED ISSUES Consider a conductive spherical particle with radius, R(=D/2), mass density, ρp, between two parallel electrodes (see Figure 1). We assume that the particle diameter is greater than 10 µm. This assumption allows us to neglect any adhesion force between the particles and the electrodes.20 The 107102-3 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 1. Schematic presentation of the problem. upper electrode has been connected to a high voltage supply with positive polarity while the lower one has been grounded. The electrodes are large enough to consider that the electric field is uniform everywhere except in a small region in the vicinity of the particle. This electric field has an average magnitude E0=ϕ/d, where ddenotes the distance between the two electrodes and ϕis the applied electric potential. In what follows, it is described how the particle motion is started and how it is continued by the microdischarge mechanism. A. Particle starts to move When the particle is in mechanical contact with the lower electrode, it acquires an electric charge Q0. Assuming that R/d≪1, the magnitude of Q0, is given by the well-known formula21,22 Q0=(π2 6)·4πεR2E0,(1) where εis the permittivity of the surrounding liquid. After charging, the particle is subjected to a total electric force, Fe0=0.832Q0E0,(2) which tends to lift it. When the sum of the buoyancy force Fband the electric one Fe0overcomes the gravitational force Fg, the particle detaches from the electrode and accelerates towards the upper electrode. After detachment, the particle behavior changes depending on the medium properties and charge exchange between the medium and the particle. Electrical charges of the particle would leak into the surrounding liquid due to the liquid electrical conductivity according to the formula23 Qp(t)=Q0e−t/τc,(3) where tis the time elapsed after the last detachment and τc=ε/σ is the charge relaxation time of the liquid (σrefers to the electrical conductivity of the liquid). From Eq. (3) it can be deduced that if the particle flight time, tflight, (the time interval between two successive particle rebound) becomes considerably lower than the charge relaxation time of the surrounding fluid, it is reasonable to ignore the charge leakage. Therefore, θc=tflight/τcis a dimensionless parameter that can affect the particle dynamics. Other dimensionless groups are presented in Section V C. B. Continuation of the motion through the microdischarge mechanism After detachment, while the conductive particle approaches to the upper electrode, the electric field between the particle and the electrode increases unboundedly until the surrounding dielectric liquid breaks down and a microdischarge takes place.20 During the microdischarge, first, the particle 107102-4 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) is neutralized and then it acquires the charge needed to be equipotential with the facing electrode. From the theoretical point of view, in the configuration under study, the charge transferred to the particle during the microdischarge event depends on its position at the discharging instant, ¯ Smd, and is obtained by Qth =4πεR2E02*,sinh2αmd ∞  n=1 cosh(nαmd) sinh2(nαmd)+-,(4) where αmd =cosh−1(1+¯ Smd)and ¯ Smd =smd/R(see Eq. (37) in Ref. 19). The value of the function within parentheses in Eq. (4) equals to π2/6 in the limit ¯ Smd →0, applying Qth →Q0(see Eq. (1)). However, this is an ideal theory, and in practice only a fraction of Qth might be transferred to the particle (see Subsection VII A for more discussion). Immediately after the microdischarge, the orientation and magnitude of the electric force Fe change, and it acts as a breaking repulsion force on the approaching particle. Nevertheless, the particle may keep on its way of motion until it completely losses its inertia and changes the orientation of motion without physical contact with the electrode. Hereafter, we refer to this part of the particle motion as “post-microdischarge motion” which is discussed in Subsection VII B. However, if the initial inertia of the particle is great enough (usually big particles), the particle could collide with the electrode. In that case, a proper contact model is required to predict the behavior of the particle during and after the collision. In both the cases (with and without mechanical collision with the facing electrode), the particle rebounds and the motion is repeated in the same way. There is a large amount of research dealing with microdischarge phenomena in different circumstances.24–26 It is rather difficult to exactly predict the distance between the particle and the electrode at which the breakdown happens (smd in Figure 1), because in addition to the electrical properties and purity of the interstitial fluid, it depends on the topographical details of the particle surface, as well. For example, even a tiny roughness on the particle and electrode surface would intensify the electric field significantly, since its size may be of the same order of magnitude as the gap between the particle and the electrode. However, as reported in Ref. 27, the microdischarge occurs at total electric field strengths (between the particle and the electrode) Emd, similar to the breakdown electric field of the interstitial fluid. III. EFFECTIVE FORCES ACTING ON THE PARTICLE Ignoring adhesion forces between the particle and the electrodes, assuming that no charge injection takes place into the surrounding liquid from the surface of high voltage electrode and that no Electrohydrodynamic (EHD) flow sets up between the electrodes, the total external force acting on the particle will be Fexternal =Fe+Ff+Fg+Fb,(5) where Fe,Ff,Fg, and Fbdenote the electric force, the fluid force, the gravitational force, and the bouncy force exerted on the particle, respectively. In addition to the microdischarge phenomenon, the other main contribution of the present work pertains to the modification of Feand Ffterms in Eq. (5). The hydrodynamic force Ffconsists of three distinguished terms, including steady drag, added mass, and history forces. Inclusion of the history force is extremely time consuming, and we ignore the last one in the present study. A more detailed discussion on the modification of the electric and hydrodynamic forces is presented in Secs. III A–III C. A. Electrostatic force In the parallel electrode system shown in Figure 1, the electrostatic force applied on a conductive sphere of radius R, having arbitrary charge Qp, immersed in a medium of dielectric constant ε, can be calculated using three distinct methods: •numerical solution of the Laplace equation 107102-5 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) •a Stokesian dynamic approach introduced recently by Drews et al.28 •the image charge method used by Pérez.19 The first method refers to the classical direct method and gives “exact” results, while the two other ones are approximate methods. The Stokesian dynamic approach is more accurate and also more complex than the image charge method. In this study we used the third method due to its simplicity, explicit expression, and less computational cost. It has been shown in Ref. 19 that for R/d≤0.1, its error is less than 3%. According to this method the total electrostatic force consists of four distinct terms: charge-external field, charge-image charge, dipole-image dipole, and charge-dipole interaction. By adding these four terms, Pérez obtained an explicit expression for the total electric force as follows:19 Fe(H)= 4πεR2E02f0(H−1)+QpE0f1(H−1)+Qp2 4πεR2f2(H−1),H≤d 2R −4πεR2E02f0(d R−H−1)+QpE0f1(d R−H−1)−Qp2 4πεR2f2(d R−H−1),H>d 2R , (6) where H=h/R.his the distance between the particle center and the lower electrode. E0denotes the uniform electric field between the two parallel electrodes in the absence of the particle, and f0,f1and f2are known functions of H(see Ref. 19 for details). Other symbols are as before. The particle electric charge Qp=Q(t)can be obtained from Eq. (3) for any given time. Using Eq. (6), a modification factor can be introduced for the electric force as λe=Fe E0Q0 = g1(H,Qp),H≤d 2R g2(H,Qp),H>d 2R ,(7a) where g1(H,Qp)and g2(H,Qp)read g1(H,Qp)=6 π2f0(H−1)+(Qp Q0)f1(H−1)+π2 6(Qp Q0)2 f2(H−1),(7b) g2(H,Qp)=−6 π2f0(d R−H−1)+(Qp Q0)f1(d R−H−1)−π2 6(Qp Q0)2 f2(d R−H−1).(7c) Some numerical values of the functions f0,f1and f2are given in Table I. The values in Table I illustrate the fact that these functions play an important role in the vicinity of the electrodes. This is caused by the electric field interaction between particle and the electrodes. This electric field interaction has been considered in the study of a conductive particle behavior in a roll-type corona-electrostatic separator by Dascalescue et al.29 using an approximate formula, Fe0=0.832Q0E01+0.21(H−1) H−0.25 ,(8) TABLE I. Numerical values of f0,f1, and f2functions in vicinity of an electrode. S=h/R−1f0(S)f1(S)f2(S) 0.0001 −220.0924 269.1032 −81.7464 0.001 −34.1924 43.0022 −12.9994 0.01 −5.4889 8.0143 −2.3355 0.1 −0.7177 2.0687 −0.4703 1−0.03 1.0685 −0.0671 1.5 −0.0100 1.0334 −0.0414 2−0.0047 1.0190 −0.0283 5−2.9004 ×10−41.0023 −0.0070 107102-6 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 2. Dimensionless electric force exerted on a detaching particle (a) and an approaching particle (b). while the same formula was used in Ref. 20 for modeling of a conductive particle behavior in the parallel electrode system. In the present work, we use a more accurate one, Eq. (7a), which is derived specifically for the parallel electrode system and is accurate for R/d≤0.1. It should be noted that the electric force applied on a detaching particle is completely different from the force on an approaching one. To illustrate this fact, we plot the force obtained from Eq. (7a) in Figure 2for different ratios of the particle radius to the inter-electrode gap. As seen in Figure 2(a), all detaching particles see the same electric force predicted in Eq. (2) regardless of the ratio R/d. Also, as expected, deviation of the electric force from E0Q0(Coulomb force) is more pronounced for bigger particles so that it reaches to nearly 25% of the gap between the two electrodes for R/d=0.09. Approaching particles, in contrast to the detaching ones, are attracted towards the electrode with an unbounded force, as shown in Figure 2(b). B. Hydrodynamic force exerted on the particle In order to determine the total force exerted on a moving particle by the surrounding fluid, it is needed to solve the continuity and Navier-Stokes equations with the appropriate boundary and initial conditions. However, when the Reynolds number (based on the particle velocity and diameter) is sufficiently low (Re ≪1), the nonlinear inertial terms can be ignored. In this way the Navier-Stokes equation is reduced to the unsteady Stokes flow equation. Landau and Lifshitz30 solved this equation for a spherical particle starting to move at t=0, in an infinite domain of fluid, using boundary conditions of no-slip on the particle surface ( u=0 at r =R) and free flow velocity in the far field ( u=  Uat r → ∞). They determined the total force applied on the particle by the surrounding fluid as31 Ff=9µυp 2R2(V−U)+1 2ρfυp dV dt −3 2ρfυp dU dt +9ρfυp 2Rνf π t  0 (dV dτ−dU dτ)dτ √t−τ,(9) where Vand Udenote the particle and fluid free velocity, respectively, υpis the particle volume, τ is a dummy time variable, and d/dt is the Lagrangian time derivative following the particle. For a stagnant medium (U=0), Eq. (9) simplifies to Ff=9µυp 2R2V+1 2ρfυp dV dt +9ρfυp 2Rνf π t  0 dV dτ dτ √t−τ.(10) The first term in Eq. (10) is the well-known Stokes drag, Fsd =9µυp 2R2V=24 Re 1 2ρfAV2=6πµRV,(11) 107102-7 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) and the second term is known as the conventional added (or virtual) mass force. This term is due to the particle acceleration in the fluid medium. Indeed, the added mass is an inertia added to the moving particle because it must accelerate (or decelerate) some volume of ambient fluid during its motion. And the last one is the history force (or fading memory force) that is often called the Basset force. It is caused by a temporary perturbation in the flow field created due to additional vortices at the solid particle surface which result from the particle acceleration. In other terms, the unsteady surface motion interacts with the developing boundary layer on the particle surface and results in a viscous force additional to the steady motion one. Further information about the concept and impact of the history force on the particle motion in various circumstances can be found in Refs. 31–37. Particle to fluid density ratio (ρp/ρf) is a key parameter dictating the history force effect on the particle motion. The smaller ρp/ρf, the larger the history force effect. We did not account for the history force in this study. C. Particle-electrode hydrodynamic interaction (wall effect) The particle motion in the vicinity of a solid wall is influenced by the wall. As a particle approaches a solid wall, the interstitial liquid is squeezed out from the gap between the particle and the wall, increasing the value of the hydrodynamic force on the particle.38 The wall changes the form and pattern of streamlines around the moving particle and subsequently affects all terms of the hydrodynamic force acting on the particle. In this section, we present the latest available modification functions of wall effect on the steady drag and added mass forces acting on a moving spherical particle. Note that the wall effect on the particle motion depends on the orientation of its motion with respect to the wall, as well as on the Reynolds number. Here, only motion perpendicular to the wall (towards or away from) is considered. To take into account the wall effect, a common practice is to introduce a correction factor λw, such that the steady Stokes drag force, Fsd =6πµRU, is modified as Fd=λwFsd.(12) 1. Wall effect on the steady drag for low Reynolds numbers (λS) For a spherical particle moving perpendicularly towards or away from a solid plane wall with constant velocity in a stagnant viscous medium having Reynolds number much less than unity, the wall correction function, λS, was derived independently by Maude39 and Brenner40 by solving the quasi-steady Stokes equation in a bi-spherical coordinate, λS=4 3sinh(α)∞  n=1 n(n+1) (2n−1)(2n+3) 2 sinh(2n+1)α+(2n+1)sinh 2α 4sinh2(n+1 2)α−(2n+1)2sinh2α−1,(13) where α=cosh−1(H)and the subscript Sin λSrefers to the Stokes flow. It is important to note that although λShas been derived disregarding the inertial terms in the Navier-Stokes equation, as Brenner40 discussed in his seminal work, its validity goes substantially above 0.5, the value normally cited as the upper limit of Stokes’ law for an unbounded fluid. In this study the instantaneous Reynolds number is considerably low near the electrodes but it grows up to 1 ∼2 in the middle of the interelectrode distance. The following function gives a very good approximation to the exact formula for the entire range of gap between the particle surface and the wall (S=h1/R−1 for the lower electrode and S=h2/R−1 for the upper one), λS=1+1 S+3 19S+4.(14) The wall effect correction function obtained from Equations (13) and (14), as well as experimental data of Adamczyk et al.,41 has been depicted in Figure 3(a). As it is seen, there is a very good agreement among three of them. Therefore, Eq. (14) can be used in calculations of the proposed model instead of the complicated formula Eq. (13) without losing accuracy. Figure 3(b) indicates 107102-8 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 3. Wall correction function. Comparison among results of Eq. (13), Eq. (14), and experimental data (a). Comparison between the nearest wall and wall-superposition approximations (b). that the particle-wall interaction starts to influence the particle motion at a distance equal to nearly ten times the particle diameter. This means that the drag force acting on a typical particle with D=1 mm submerged in a parallel electrode system with 10 mm inter-electrode distance will be continuously affected by the presence of electrodes via the wall effect on the whole gap. Hence, it seems that λSplays an important role in the up and down motion of conductive particle submerged in a viscous liquid. In order to take into account the simultaneous influence of both electrodes on the particle motion, one approximate approach is to superimpose the effect of both single electrodes, λsup =λ1+λ2−1,(15) in which λ1and λ2are the wall effect correction functions of the nearest and farthest electrodes, respectively. The unity has been subtracted from the right side of Eq. (15) in order to avoid the double counting of the drag force applying on the particle in an unbounded medium. Since the governing equations for the creeping flow (steady Stokes equations) are linear, it seems that the wall-superposition approximation is a reasonable method to take into account the simultaneous influence of both electrodes. The nearest wall and wall-superposition approximations have been compared in Figure 3(b) for R=1 mm, and d=10 mm. As it is clear, they give the same magnitudes in the vicinity of the wall. The wall-superposition approximation predicts slightly greater values in the central region. 2. Wall effect on the added mass force (λam) As it has been said, the added mass force is due to the acceleration of the fluid caused by the particle motion. Obviously, when a particle is accelerated in an unbounded medium, the surrounding fluid can displace freely in every direction. However, the existence of a wall in the vicinity of the particle would restrict the displacement of the fluid, causing an increase of the added mass force. In spite of this consideration, it has been shown that, in the case of a spherical particle normally approaching to a solid plane wall, the added mass enhancement is small. For example, in an inviscid flow, the wall effect imposes a factor of (1+0.375/H)to the unbounded one, corresponding to a maximum enhancement of 37.5%.33 Most recently, Simcik and Ruzicka42 conducted a 3D numerical simulation of full Navier-Stokes equations and obtained an accurate correlation for the added mass force acting on a spherical particle near a solid wall as Fam =(1 2+0.1785 H)ρfυdV dt ,(16) 107102-9 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) which gives a wall modification factor as λam =1+0.357 H.(17) We use Eq. (17) as the wall modification factor of the added mass force in our calculations everywhere is needed. It is worth to note that according to Eq. (17) the maximum value of λam is equal to 1.357, which corresponds to the moment of contact between the particle and the wall (i.e., H=1). Therefore, based on the analysis of Simcik and Ruzicka,42 the added mass force will increased just up to 37% which is not much in comparison with the wall effect on the steady drag (see Eq. (13) or Eq. (14)). So it seems that the steady drag force modified by the wall correction function (λS) is completely dominant near the wall and wall effects on the added mass forces can be reasonably neglected. However, we did include it in our calculations in order to give a quantitative comparison. IV. PARTICLE–WALL COLLISION AND REBOUNDING PROCESSES Particle–wall as well as particle–particle collisions have an effective impact on liquid–solid flows because they influence inter-phase transportation phenomena.43 Particularly, particle–wall collisions would affect the heat transfer rate between a submerged surface and its surrounding medium by thinning the thermal and viscous boundary sublayers. The change in particle velocity after rebounding defines the restitution coefficient as kr=Vr Vi ,(18) where Viand Vrdenote the particle velocity before and after a rebound, respectively. It should be noted that the impact velocity in a viscous medium is defined in a somewhat different way from that of a dry collision. In the viscous case, contrary to the dry collision case, the impact velocity is not the particle velocity in the moment just before the real physical contact. What is taken as the impact velocity is the maximum particle velocity before the rebound. In a gaseous environment, the maximum particle velocity is the same as the particle velocity at the real contact. This is not true for a viscous collision due to a sudden and sharp decrease in the particle velocity right prior to its collision with the wall.44,45 Such an enormous reduction in the approaching particle velocity is originated from the particle-wall interaction, which is brought about by a lubrication effect due to the interstitial liquid and can be quantified by the correction function, λS(for the Stokes flow), represented in Eq. (13). The restitution coefficient is often expressed as a function of the ratio of the particle inertia to viscous forces, namely, the Stokes number, St =1 9 ρp ρf Rei,(19) where Rei=ρfViD/µfis the particle Reynolds number based on the impact velocity Vi. Almost all the experimental studies indicate that there is a critical Stokes number, Stc≈10, below which no rebound occurs.44–47 In the current study, the Stokes number is quite smaller than Stc≈10. Therefore, the restitution coefficient is not considered. V. MODELING APPROACH AND EQUATION OF MOTION In this section, first, we describe a model through which the microdischarge characteristics are treated. Then the particle equation of motion is presented. Also, we talk about how the model can be extended to intermediate Reynolds numbers. Finally, dimensionless form of the governing equations is obtained, and all the relevant dimensionless groups are introduced and discussed. A. Modeling of the electrostatic force In order to take into account the microdischarge, an additional parameter, denoted by smd, is introduced (see Figure 1). This parameter is the distance between the surface of the particle and the 107102-16 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) The minimum charge that the particle needs to be able to rebound (to lift) from the lower electrode (Qp,min), can be determined by Fe≥Fg+Fb, which gives g1(Hmin,Qp)E0Q0≥(mp− mf)g, where Hmin is the particle center position at the instant of rebound (the instant that velocity of particle switch sign). This inequality can be converted to an equality by replacing Qp with Qp,min, i.e., g1(Hmin,Qp,min)E0Q0=(mp−mf)g. Using Eqs. (27) and (28), it is simplified to g1(Hmin,Qp,min)=Πg. After substituting g1(Hmin,Qp,min)from Eq. (7b), it results in a second order polynomial equation as π2 6*, Qp,min Q0+- 2 f2(¯ Smin)+*, Qp,min Q0+-f1(¯ Smin)+6 π2f0(¯ Smin)−Πg=0,(31a) where Smin =Hmin −1. A similar analysis is applicable for the particle bouncing (dropping) from the upper electrode. In this case the necessary condition is Fe≥ −(Fg+Fb). Indeed, for possibility of the particle rebound from the upper electrode, the attractive electric force (negative force) should be necessarily weaker than the particle apparent weight. Such a constrain results in π2 6(¯ Qp,min Q0)2 f2(¯ Smin)+(¯ Qp,min Q0)f1(¯ Smin)+6 π2f0(¯ Smin)+Πg=0,(31b) where ¯ Qp,min stands for the minimum charge required for the particle rebound from the upper electrode and ¯ Smin denotes the dimensionless minimum gap between the particle and the electrode at the instant of rebound. In Figure 10, the acceptable roots of Eqs. (31a) and (31b) are depicted as a function of ¯ Smin for different values of Πg. As it is observed, the minimum charge required for the particle while bouncing from the upper electrode might be considerably smaller than that of the lower electrode. Also, as expected, for a buoyant particle (Πg=0), there is no difference between Qp,min and ¯ Qp,min. But, as Πgincreases their difference grows significantly. For a given condition, a weighty particle may rebound from the upper electrode surface with less charge in comparison with a lighter one. However, it should be noted that such conclusion is correct only when the rebound takes place at the same ¯ Smin for the two compared particles. Another remarkable result from Figure 10 relates to the effect of ¯ Smin on the minimum required charge. As seen irrespective of both the electrode location (lower or upper one) and Πgvalue, lesser minimum separation at the instant of rebound requires more charge for the rebounding particle. FIG. 10. The minimum required charge as a function of the minimum particle-electrode gap at the instant of rebound for different Πg. 107102-17 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 11. Three distinct zones determining the possibility of the particle rebound from the upper and lower electrodes. In Figure 11, the electric-to-apparent weight ratio is plotted versus ΠQ=Qp/Q0for different values of the minimum particle-electrode gap. As a general rule, three distinct zones are detectable. In thezoneg1(Hmin,Qp)/Πg>1,theparticlereboundfromboththeupperandlowerelectrodeispossible. In the zone g1(Hmin,Qp)/Πg≤ −1, the particle can rebound neither from the upper electrode nor from the lower one. In the middle zone −1< g1(Hmin,Qp)/Πg≤1, the particle can separate from the upper electrode and drop low, but it cannot lift from the surface of the lower electrode. B. Post-microdischarge motion In order to reveal the effect of the microdischarge on the particle dynamics, first we studied the length of the post-microdischarge motion, ∆hpmdm, for input parameters of the case 1 in Table II. The particle trajectories were derived by solving Eq. (24) with the initial conditions: h(0)=R+smin and h′(0)=0 for different values of smd and E0in the range 1 µm≤smd ≤100 µm and 0.6 MV/m≤E0≤1.5 MV/m. Eq. (24) cannot be employed to predict the motion of a particle initially in full contact with the electrode, because as it is clear from Equation (13) or (14), the wall effect correction function is singular at S=0, representing no possible contact. In other words, since S→0 implies λS→ ∞, it is essential to use a minimum distance for the particle-wall gap, say smin. This strategy has been used by many researchers, for example, Drews et al.28 used smin =0.0247R(0.7µm) in their modeling approach based on the assumption that the microdischarge takes place when the particle-electrode gap reduces to 0.7µm (i.e., smd =0.7µm). Indeed, they assumed that the particle motion is immediately reversed after the microdischarge. We evaluated this assumption and found out that it is correct when the microdischarge occurs in a point very close (in the order of one micron) to the next electrode (see Figure 13). However, as discussed in Subsection II B, an earlier microdischarge can take place for many reasons, particularly due to asperities on the particle and electrode surfaces. In that case, the “minimum distance” is not exactly the same as the “microdischarge point,” and the particle can keep its direct motion towards the facing electrode. For such cases, the “minimum distance” should be determined by the post-microdischarge motion. Variable time steps were utilized for solving the equation of motion. But in order to capture precisely the details of the particle trajectory and velocity, a constant time step size of 1 µs was used in the region close to the microdischarge position. Independency of the results was carefully checked against the time step size. Residual of a typical solution is depicted in Figure 12. It is in the order of 10−18 that indicates an excellent accuracy for the obtained numerical solution. The results are shown in Figure 13 for the third microdischarge using logarithmic scale for x and yaxes. As seen, both smd and E0can influence the length of the post-microdischarge motion. A stronger electric field can drive the particle closer to the opposite electrode resulting in smaller smin. However, the wall hydrodynamic effect acts on the opposite side and prevents the particle from 107102-18 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 12. Residual of a typical numerical solution in the present study. more advancing. Indeed, the competition between λeand λSdetermines the position of the reverse point in the particle trajectory. It is clear from Figure 13 that ∆hpmdm decreases significantly with decreasing smd, and ∆hpmdm ≈0 for smd ≤1µm which conforms to the assumption of Drews et al.28 Generally, Figure 13 gives smin =f(smd), where fis a function between the microdischarge point position, smd, and the possible minimum particle-electrode separation, smin. This function can be used as the initial position of the particle in using the proposed model as h|t=0=R+smin where,smin =f(smd).(32) (Also, instead of R+smin, one would use an arbitrary initial position for the particle. But, in that case, it should be noted that the time evolution of the predicted trajectory is under the influence of the particle’s initial position until the first microdischarge occurs. Thus, the experimental trajectory should not be compared with the first cycle of the predicted trajectory because the experimental data pertain to a consistent periodic bouncing motion of the particle.) We studied the trajectory of a glassy carbon particle for three different values of smd: 100 µm, 10 µm, and 1 µm. The results were derived by solving Eq. (24) with the initial conditions: h(0)=R+smin and h′(0)=0, where smin was determined by smin =f(smd)using Figure 13. All other parameters of the model were the same as in case 1 (see Table II). The particle trajectory, velocity, and acceleration are depicted in Figures 14(a)–14(c), respectively. It is seen that the location of the microdischarge point can significantly influence the temporal evolution of the particle trajectory. Table IV shows some relevant quantitatives. Focusing on the case smd =100 µm, as an instance, it is seen that right at tmd =101.970 ms, the particle reaches the position which has been previously specified for the microdischarge phenomenon (point C in Figure 4). The particle velocity is almost Vmd =47.14 mm/s, and it suddenly decreases to nearly zero in an interval of ∆tpmdm =583 µs. During this interval, the particle experiences a huge negative FIG. 13. The length of the post-microdischarge motion, ∆hpmdm, as a function of smd and E0. Other parameters are the same as case 1 in Table II. 107102-19 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 14. Effect of the microdischarge phenomenon on the particle dynamics. (a) The particle trajectory; (b) temporal evolution of the particle velocity; (c) dimensionless acceleration of the particle for smd =1µm ( ¯ Smd =10−3) immediately after the microdischarge near the upper electrode. Other input parameters are the same as in case 1. acceleration and can move forward just 10 µm (∆hpmdm =10 µm). Consequently, its distance from the front electrode deceases from 100 µm to almost 90 µm. Therefore, the proposed model predicts no physical contact between the particle and the electrode provided that the microdischarge happens at smd =100 µm. This is in agreement with the experimental observation reported in Ref. 23. Overall behavior of the particle for the other two cases is similar to this case. But there is significant change in quantitative values of the pertaining parameters, particularly in λeand ∆hpmdm. Comparison between values of these two parameters demonstrates that lesser smd results in lesser ∆hpmdm. 107102-20 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) TABLE IV. Some parameters that characterize the particle dynamic behavior after the microdischarge near the upper electrode. smd (µm)¯ Smd λetmd (ms)Vmd (mm/s)∆tpmdm (µs)∆hpmdm (µm) 100 10−13.24 101.970 47.14 583 10 10 10−215.06 123.666 24.53 121 0.9 1 10−383.87 142.516 13.78 16 0.05 This means that in a real case it is difficult to realize whether the particle rebound occurs after a real contact with the electrode or not. As seen in Figure 14(b), in the discharging point (point C in Figure 4), there is a sharp change in the particle velocity. This is also clear from Table IV. As an instance, for smd =1µm, the particle velocity decreases from Vmd =13.78 mm/s to V≈0 during an interval time of ∆tpmdm =16 µs. Such a very rapid change in the velocity causes a large discontinuity in the particle acceleration as observed in Figure 14(c). This large discontinuity in the particle acceleration would cause the history force to be important in the rebounding trajectory. (Remember that the history force contains the acceleration term.) An outcome of this section is the demonstration of the fact that the time and length scales of the particle post-microdischarge dynamic are, respectively, microsecond and micrometer. This is a very important result because such small time and length scales would influence the amount of charge obtained by the particle through the post-microdischarge motion. Indeed, since the microdischarge prevents the particle from a “mechanical” contact with the facing electrode, the charge exchange occurs throughout a thin layer of the interstitial liquid. In this case, the time scale of the charge exchange event might be comparable with that of the post-microdischarge motion. Based on such hypotheses, in addition to the particle inertia, even the chemical properties of the dielectric liquid plasma and roughness of the particle and electrode surfaces would play a role in the process of charge exchange during the microdischarge event. Further theoretical and experimental investigations are needed to shed more light on the details of this issue. An appropriate discussion can be found in Ref. 51. Due to very small time and length scales, a laser interferometry method55 is required for capturing the details of the post-microdischarge motion. Also, it should be noted that after the microdischarge, the viscosity of interstitial liquid may change due to tiny bubbles which are produced because of vaporization of water content in the liquid.24 Also any possible sudden change in the particle temperature can influence the viscous drag force as explained in Ref. 53. We did not account for such effects in our study. C. Comparison with experimental data and discussion In order to use the proposed model, it is required to determine Qmd. Here, first we explain a simple way to solve this issue. For applying this approach, at least three experimental data points, representing the particle position-time, are needed. These data points should be quite close to each other to save the accuracy. They would be extracted through image processing of the particle motion in the middle of the interelectrode distance. Assume that the three data points are (ha,0),(hm,tm), and (hb,tb), where, hm≤d/2 and 0 <tm<tb. Using these data, we have Va=(hm−ha)/tm,Vm= (hb−hm)/(tb−tm), and V′ m=(Vm−Va)/tm. By applying Eq. (24) for upward motion of the particle at h=hm, after some simplification, it reads (mp+λam,m 1 2ρfυp)V′ m=λe,mE0Q0−(mp−ρfυp)g−λS,m(6πµR)Vm.(33) The parameters λam,m,λS,mand Q0can be obtained, respectively, from Eqs. (17), (13), and (1) at h=hm. By equating the dimensionless electric force λe,m, from Eqs. (33) and (7b), we have 6 π2f0(Hm−1)+(Qm Q0)f1(Hm−1)+π2 6(Qm Q0)2 f2(Hm−1) 107102-21 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) =1 E0Q0(mp+λam,m 1 2ρfυp)V′ m+(mp−ρfυp)g+λS,m(6πµR)Vm,(34) where Qmstands for Qat h=hm. Now, the value of ΠQ=Qmd/Q0can be easily obtained from Eq. (34). Since Qm=Qmd exp(−tm/τc), therefore Qmd =ΠQQ0exp(tm/τc).(35) The above mentioned procedure is used for obtaining Qmd from Eq. (35) whenever the proposed model is applied for the prediction of the particle trajectory. Figure 15 shows a comparison between the predicted results of the model and the experimental data for the trajectory and velocity of a glassy carbon particle of 2 mm diameter and of 1500 kg/m3 density submerged in silicon oil of 100 cSt viscosity. The data labeled as “experimental” velocity in this and the following figures have not been measured directly. They have been calculated by first order forward step differentiation using the experimental data of the particle trajectory. All the required parameters are given in Table II under the column of case 1. These values are those reported in Ref. 3. The dielectric constant and conductivity of the working liquid have not been given in the mentioned reference, so we assumed the sounded values of ε=2.5, σ =1.3E−12 S/m which results in almost τc=17 s for the charge relaxation time. The three experimental data points (4.26 mm, 88.18 ms), (5.00 mm, 96.2 ms), and (5.83, 103.88 ms) were selected from the particle trajectory graph given in Ref. 3to be used for obtaining Qmd from Eq. (33). Thus, these data give (ha,0)=(4.26 mm,0 ms),(hm,tm)=(5.00 mm,8.02 ms)and (hb,tb)=(5.83 mm,15.81 ms). Using these data, it is resulted that Qmd =0.99Q0. This is in contradiction with Qmd ≈0.4Q0−0.6Q0 FIG. 15. Numerical and experimental data of the glassy carbon particle motion (case 1); temporal evolution of the particle trajectory (a) and the particle velocity (b). Numerical results were depicted by lines and obtained by using smd =2µm (¯ Smd =2×10−3). 107102-22 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) which has been reported in Ref. 3. In fact, simplification of the particle equation of motion particularly, elimination of the wall effect, has caused such a significant error. As seen in Figure 15, neglecting the wall effect (λS=1) causes a very large error in predicting the temporal evolution of both the trajectory and the velocity of the particle. In contrast, when the wall effect is taken into account, the model results show good agreement with the experimental data. In particular, the upward motion is predicted very well when an approximate value of 2 µm is considered for the microdischarge point ( ¯ Smd =2×10−3), which is a reasonable value. In contrast, the downward motion is not well predicted by the model. As reported in Ref. 3, the glassy carbon particle settles on the upper electrode for a while (tdwell ≈35 ms) before detaching from it. The origin of this behavior is not clear yet. However, as seen in Figure 15(a) even after shifting the rebound trajectory from t=135.7 ms to t=135.7 ms +tdwell, again considerable discrepancy exists between the numerical prediction and the experimental data. Possible reasons are as follows: The most probable reason for such disagreement is the charge transfer between the glassy carbon particle and the upper electrode. Probably less charge is transferred to the particle during the “collision” with the upper electrode. Unfortunately, the charge exchange phenomenon in the configuration under study has not been well understood up to now, and as we mentioned earlier in Subsection VII A, there is a contradiction in the available literature in this regard and complementary investigation is required to find more details. Here, a natural question comes to mind: Is there any special amount of charge, Qmd,2, for which the proposed model predicts correctly the downward motion of the particle? (Qmd,2denotes for the particle initial charge at the instant of particle rebound.) The answer is positive. Using try and error technique, we found that the downward motion can be correctly predicted assuming Qmd,2=0.81Qth as seen in Figure 16. This approach has been applied previously in some models. For example, Bernanek et al.54 used the exchanged charge value as a free parameter to adopt their model predictions with the experimental data. Figure 16 indicates that the reduction of the particle charge in downward motion can significantly retard the particle downward trajectory. We shall refer to Figure 8in order to find out the reason for this considerable retardation in the downward motion of the particle. As seen in Figure 8, when the particle is very close to the facing electrode (some microns), even a little reduction in the particle charge can result in large decrease in the repulsive electric force. Even attractive electric force can be achieved by extra reduction of the particle charge which can stick the particle to the electrode and prevent from falling. It seems that settling/dwelling time originates from the particle charging phenomenon. In Figure 17, the numerical results of the present model are compared with the experimental data for a steel particle of 2 mm diameter (case 2 in Table II). The electric potential is 10 kV (in contrast to 8 kV in the case 1). Other parameters are the same as in case 2 and can be seen in Table II. Two series of experimental data are seen in the figure. One of them has been adopted from Ref. 3(square symbols) and the other pertains to the present work (circle symbols). Three FIG. 16. Comparison between the experimental and predicted trajectory in which the dimensionless parameter ΠQhas been considered as free parameter for the upper electrode. 107102-23 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 17. Numerical and experimental data of the steel particle motion (case 2). Temporal evolution of the particle trajectory (a) and the particle velocity (b). The legend depicted in figure (b) holds for figure (a), as well. All the numerical results have been obtained using smd ≈1µm ( ¯ Smd ≈1×10−3). For better prediction of the experimental data, the model is run with ΠQ=0.85 for the downward motion of the submerged particle. different frames showing the upward motion of the particle pertaining to the current study are shown in Figure 18. Present experiment has been carried out with the same condition as Ref. 3. However, as seen in Figure 17, there is a little difference between the two series of experimental data. Most probably, this is due to a little difference in the viscosity of silicon oil (100 cS) which can occurr due to a limited difference in temperature of the two laboratory. Such estimation is confirmed by the model because the model results in a good prediction when it is run with a 5% lesser viscosity (see dashed-dotted line in Figures 17(a) and 17(b)). Generally, the results in Figure 17 are similar to the previous case. The position-time function is incorrectly predicted without taking the wall hydrodynamic effect into account (see Figure 17(a)). Very good agreement between the model numerical results and the experimental data for the upward motion is observed implying Qmd =Qth (ΠQ=1). However, a little faster downward motion is predicted for Qmd =Qth. Again, and similar to the case 1, the particle-electrode charge exchange is responsible for such discrepancy. This minor defect is improved by applying Qmd =0.85Qth for the downward motion (V<0). Also, in order to show the effect of the particle apparent weight, we plot the trajectory and velocity of an equivalent buoyant particle in Figure 17 (red dashed line). For this end, the solution code was run with −→ Fg+−→ Fb=0 (see Eq. (21)). In this case, there is no difference between the upward and downward motion, and the only driving force is the electric force. Therefore, a comparison between E0Q0and (mp−mf)gcan produce a scale to determine the possible difference between the 107102-24 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) FIG. 18. Three different frames of the steel particle motion. The particle diameter and the gap between the two parallel electrodes are 2 mm and 10 mm, respectively. The parallel electrode system and the particle submerged in silicon oil (100 cS). upward and downward motion of a typical particle. The apparent weight to the electric force ratio is 2(ρp−ρf)gR/π2εE02which is exactly the same as Πg(see Eq. (28)). In Figure 19, period ratio of the downward-to-upwardhalf-cycles,Td/Tu,isdepictedasafunctionofΠg.ItisclearthatasΠg→1,more difference appears between the upward and downward half-cycles periods in the particle trajectory. The model predictions are in good agreement with the experimental data as seen in Figure 15. FIG. 19. The period ratio of the downward-to-upward half-cycles as a function of apparent weight to the electric force ratio. 107102-25 Eslami, Esmaeilzadeh, and Pérez Phys. Fluids 28, 107102 (2016) In spite of very good prediction shown in Figures 15–17, it is important to note that the history force has been ignored in the model. It is known that the history force makes a delay in the accelerated motion of a particle.31 Also, higher acceleration causes longer delay time. Although, oftentimes the history force can be ignored without anxiety about the accuracy of the results, it is important to note that in predicting the particle rebound trajectory, it might be important mainly due to its old wake.55 This has been confirmed by some studies. For example, Gondret et al.33 conducted an experimental study on the bouncing motion of solid spheres onto a solid plate. They compared some experimental data with numerical results and concluded that it is essential to take into account the history force in prediction of the bouncing trajectory for Reynolds number up to 103. Rostamei et al.56 extended Gondret’s study in a similar way to higher Reynolds number and found out similar result. They remarked that in prediction of the particle motion, from the starting point of motion to the wall impact, the history force effect is small (about 1%–4.3%) and can be ignored in comparison with the other hydrodynamic forces. But it significantly influences the bouncing trajectory of the particle after the first collision, even for the Reynolds numbers as much as 5000. Based on these findings, it deserves to investigate the history force effect on the particle trajectory in the configuration under study. We will focus on this issue in our next work. D. Effect of ionic conductivity of the working liquid It is worthy of particular attention that the conductivity of the surrounding liquid would strongly affect the particle dynamic behavior, since it determines the charge relaxation time and, subsequently, the Coulomb force. Therefore, it is very important to have a sufficiently accurate value of the liquid conductivity for the correct prediction of the particle velocity and trajectory. On the other hand, the conductivity of dielectric liquids would be affected by several agents including chemical impurities, water content, temperature, strengths of the applied electric field, and even the surface characteristics of the measuring electrodes.22,57 For instance, a very small amount of impurities would increase the conductivity of a dielectric liquid over hundred times. So, it is desirable to measure the liquid conductivity in the same conditions of the test. For this reason, the ionic conductivity was measured by two different methods including IRLAB system (model LDTRP-2) and directly using the two parallel electrode system of our test cell. The measuring procedure in our test cell was as follows. First, the plastic box and the electrodes were washed with cleaning liquid and dried properly. Then, an adequate amount of the sample liquid was poured into the plastic box. After adjusting the gap between electrodes, a DC voltage of a given value ϕwas applied on the upper electrode and the consequent electric current, I, was recorded by the multimeter after almost 2 minutes (steady state). The conductivity was calculated by σ=(A/d) (ϕ/I)in which Ais surface area of the electrode and other parameters are as before. There was no significant difference between the values given by IRLAB and the one measured by our system for low values of the applied voltage. However, the electric current showed deviations from the Ohmic behavior for voltages high enough. For the computations, we used the value of the liquid conductivity obtained using the parallel electrode system. A series of experiments were conducted in order to examine the effect of fluid ionic conductivity on the dynamic behavior of the particle. All parameters pertaining to these experiments have been given in Table II (case 3 A, B, and C). We used a small amount of butyl-alcohol as an additive to increase silicon oil conductivity. The results are depicted in Figure 20. First of all, it is observed that there is a very good agreement between the experimental data and the results obtained by the model. Such agreement demonstrates that the model has much predictive value. Note that in Figure 20(a) for each time t, the reported corresponding experimental data, h, have been derived by an arithmetic averaging of several data points (usually more than five). As seen in Figure 20(a), when the surrounding liquid is pure silicon oil with a conductivity of σ=2.63E−12 S/m (case 3A), the particle keeps its upward motion until it collides with the upper electrode and then rebounds downwards. In this case the charge relaxation time is τc=8.4 s, which is sufficiently larger than the upward motion time (θc=0.125). So the charge leakage does not modify the electric force during the particle motion. However, by adding 0.25 wt. % butyl-alcohol into the pure sample oil, its conductivity increases nearly five times (σ=5.22E−11 S/m) which leads to τc=1.7 s (case 3-B).