Thermal and electrohydrodynamic plumes: a compartive study

THERMAL AND ELECTROHYDRODYNAMIC PLUMES: A COMPARATIVE STUDY P. A. V´azquez, A. T. P´erez and A. Castellanos Dpto. Electr´onica y Electromagnetismo, Facultad de F´ısica. Avda. Reina Mercedes s/n. 41012 Sevilla. Spain. This paper deals with self similar thermal and electrohydrodynamic (EHD) plumes. The former arises from hot lines or points, whereas the latter arises when sharp metallic contours submerged in non conducting liquids support high electrostatic potential, resulting in charge injection. Although the motive force is buoyancy in one case and Coulomb force in the other, it is shown that the solution for EHD plumes is the same as for thermal plumes in the limit of large Prandtl numbers. We present the analysis of axisymmetric plumes for large values of Prandtl number, and this analysis is subsequently applied to EHD plumes. The validity of the approximations for EHD plumes is discussed in the light of experimental data. I. INTRODUCTION Thermal plumes arising from line or point sources of heat have attracted the attention of the investigators for a number of years. These flow structures appear in many practical problems and its description have been important from the industrial point of view. When the fluid of interest is a viscous oil, the case of high Prandtl number deserves special consideration. Also plumes in high Prandtl number fluids are relevant in the context of convection in the Earth’s mantle.[1, 2] The boundary layer equations for natural convection, i. e., the equations for thermal plumes, are known to possess self similar solutions, that reduce the partial differential equations to a set of ordinary differential equations. These self similar solutions are also the starting point for the stability analysis of the laminar plumes.[3, 4] Although the problem of two-dimensional plumes for very large Prandtl numbers have been studied long time ago,[5] and this analysis was refined later,[6] the authors are not aware of the existence of a similar analysis for the axisymmetric plumes, except for some numerical investigations.[7] The first aim of this paper is to present such analysis. We demonstrate that the velocity of the fluid at the central line of the plume diverges as ln √Pr, provided the total heat flux is finite. On the other hand electrohydrodynamic (EHD) plumes have been observed experimentally.[8, 9, 11] When high voltage is applied to electrodes with sharp edges or points immersed in dielectric liquids, there is injection of ions from the electrode into the liquid.[12] The Coulomb force acting upon the injected space charge plays the role of buoyancy. There is a close analogy between EHD and thermal plumes. As a matter of fact the equations for EHD plumes, under a number of assumptions, are those of the thermal ones for Pr → ∞.In this context the second aim of this paper is to discuss the analogies and differences between EHD and thermal plumes, and to establish the validity of the approximations involved in the EHD case. As self similar solution of two dimensional EHD plumes is already discussed in,[11] we focus on the axisymmetric case in this paper. II. EQUATIONS FOR THERMAL PLUMES Natural convection resulting from a point source of heat can be considered as an axisymmetric laminar, steady flow. The governing continuity, momentum and energy equations simplify with the Boussinesq approximation and boundary layer assumptions to the following form:[13] ∂ru ∂x +∂rv ∂r = 0 (1) u∂u ∂x +v∂u ∂r =ν r ∂ ∂r(r∂u ∂r ) + gβ(T−Ta) (2) u∂T ∂x +v∂T ∂r =χ r ∂ ∂r (r∂T ∂r ) (3) Here xis the coordinate along the flow, rthe coordinate perpendicular to the flow, uand vare, respectively, the velocities in the xand radial directions, βthe coefficient of volumetric thermal expansion, νthe viscosity, T temperature, χthe thermal diffusivity and Tais the temperature of the ambient fluid. If the ambient fluid is at 2 rest, there is not any pressure contribution in equation (2). Viscous dissipation and compressibility effects have been neglected in the energy equation. These equations are complemented with the corresponding boundary conditions. On one hand at the central line (r= 0) the radial velocity has to be zero, and, hence, uhas a maximum. On the other hand velocity and temperature differences decay to zero far from the plume. That is: v=∂u ∂r =∂T ∂r = 0 at r= 0,(4) u=T−Ta= 0 at r→ ∞ (5) An additional condition is obtained considering the heat flux across a horizontal plane. This flux is: Q= 2πρcpZ∞ 0 (T−Ta)ur dr (6) where cpis the specific heat of the fluid. In the absence of other sources of heat integration of equation (3) on an horizontal plane implies the constancy of Q. It is possible to define a stream function Ψ so that equation (1) is automatically satisfied: u=1 r ∂Ψ ∂r and v=−1 r ∂Ψ ∂x (7) and we have now two unknowns: Ψ and T. The equations (1-3) can be transformed into ordinary differential equations using the similarity method. As there is not a typical length in the axial direction the velocity and thermal profiles scale with the plume width. Mathematically this means that velocity and temperature profiles are the same at each xwhen expressed as a certain combination of xand r,η=η(x, r). In our case it is easy to show that η= (Gr)1/4r x(8) is an adequate choice. Introducing similarity functions fand θthat depend only on η: Ψ = νxf(η) (9) T−Ta=ν2 gβx3Gr θ(η) (10) equations (1-3) give: f′′′ η+f−1 ηf′ η′ +θ= 0 (11) (ηθ′)′+ Pr(fθ)′= 0 (12) which are now a set of ordinary differential equations. The associate boundary conditions are: f′ η′ (0) = f η(0) −f′(0) 2=θ′(0) = f η ′ (∞) = θ(∞) = 0 (13) 1 = Z+∞ 0 f′θ dη (14) The quantity Gr = gβx3(Ts−Ta)/ν2θ(0) is the Grashof number. It measures the relative strength of buoyancy and viscous forces (Tsis the temperature at the central line of the plume). Pr = ν/χ is the Prandtl number. Gr is related to Qby equation (6) which gives: Gr = gβQ 2πρν3cp x2(15) 3 In this way the velocity of the plume is given by: u=ν xGr1/2f′ η(16) and its thickness by: ∆ = x Gr1/4(17) Velocity and temperature profiles are obtained integrating numerically equations (11-12) using a shooting method. The integration is carried out as if it was an initial conditions problem, with f′/η(0) and θ(0) known. We use a Runge-Kutta method to obtain f′/η(∞) and Rf′θ dη. This procedure defines these last values as functions of f′/η(0) and θ(0). Then, equations f′/η(∞) = 0 and Rf′θ dη = 1 are solved with a Newton-Raphson method. The profiles so obtained are shown in Figures 1 and 2 for different Prandtl numbers. The velocity at the central line r= 0 is proportional to f′(0)/η. This value, along with θ(0) is listed in Table I for different Prandtl numbers. This table is consistent with other authors’ results.[7] III. THE CASE OF LARGE PRANDTL NUMBER We are specially interested in axisymmetric plumes when Pr → ∞,because it is in this limit that thermal and EHD plumes behave in the same way. In this limit, the thickness of the thermal boundary layer is zero, and this allows us to solve the problem without solving the energy equation. Equation (12) can be integrated to give: θ(η) = θ(0) exp −Pr Zη 0 f(t) tdt(18) Condition f/η −f′/2 = 0 at η= 0 is fulfilled only if f(0) = f′(0) = 0.Therefore, when η << 1 : f(η) = f′′(0)η2/2 + O(η3) (19) And we get for θ: θ(η) = θ(0) exp −Pr f′′ (0) 4η2(20) It is clear that θ(η) goes to 0 at a distance of order η=δt∼1/√P r and the value of θ(0) can be obtained from equation (14): θ(0) = P r/2.For Pr → ∞,buoyancy acts only along the axis, and is negligible for any η6= 0. We seek now for an additional boundary condition in the limit P r → ∞ that allows us to solve the problem using only equation (11) for η > 0. Integrating (11), using (19), neglecting terms of order ǫand taking into account that in the limit Pr → ∞ is Rǫ 0ηθ dη = 1/f′′ (0) we obtain: f′′(ǫ)−f′ η|η=ǫ=−1 f′/η |η=0 (21) No regular function for η > 0 can fulfill this expression unless f′/η goes to infinity at η= 0. Let us analyze closer the way in which f′/η diverges. For high enough Pr the velocity is constant in the thermal boundary layer, as the thermal layer δtis much smaller than the momentum boundary layer (see figure 3). Therefore it is possible to think of a distance ǫmuch smaller than 1 but yet greater than δtso as to assume the following velocity profile: f′(η) η=αif η < δt h(η) if δt< η < ǫ (22) By continuity is α=h(δt).Outside of the thermal layer is f′′′ + (f−1) f′ η′ = 0 (23) 4 But f≃f′′(0)η2/2,so f << 1 for η < ǫ, and equation (23) can be approximated by: f′′′ −f′ η′ = 0 (24) Integrating once we have: f′′ −(f′/η) = λ, and this gives an equation for h(η) : ηh′=λ(25) whose solution is h=λln η. Therefore: f′(η) η=λln δtif η < δt λln ηif δt< η < ǫ (26) Equation (21) gives now λ2=−1/ln δtso f′ η=p−ln δt(27) inside the thermal layer. Since δt≃1/√Pr is: f′ η≃√ln Pr (28) and the velocity at the axis of the plume diverges as √ln Pr. Table I gives the results of the numerical solution of equations (11-12) for different Prandtl numbers. The best fit to those values of f′/η is: f′ η|0= 0.38 + 0.67√ln Pr (29) in complete agreement with the analysis. The conclusion is that the transport of a finite heat flux requires an infinite velocity at the axis for Pr → ∞. Even though the velocity divergence is very weak, it is in contrast with the two–dimensional case. In that case,[5] the velocity is finite in the limit P r → ∞ and it is: u(x, y) = 1.42ν xGr1/2(30) The Grashof number is related to the heat flux per unit length in transverse direction by: Gr = Qgβx3 2√2ρcpν3!4/5 (31) and the velocity in the central plane is finite provided Qremains finite. IV. EHD PLUMES Insulating liquids stressed by intense electric fields are subjected to strong injection of charge from metallic sharp points or edges. The injected ions are acted upon by the electric field, and through collisions with the neutral molecules give rise to strong convective motions of liquid. The flow so developed is similar to a thermal plume, the Coulomb force playing the role of buoyancy. In the following we will analyze the similarities and differences between thermal and EHD plumes. 5 A. EHD equations The general steady state equations for a perfect insulator liquid with permittivity ǫ, kinematic viscosity νand density ρare:[14] ∇·E=q ǫ(32) ∇·j= 0 (33) ∇·u= 0 (34) (u·∇)u=−1 ρ∇p+ν∇2u+q ρE(35) here Eis the electric field, qthe charge density, jthe current density, uthe liquid velocity and pthe pressure. The current density is given by: j=KqE+qu−D∇q(36) where Kis the ion mobility and Dthe diffusion coefficient. Usually in most EHD problems diffusion is neglected.[10] B. Axisymmetric EHD plumes The two-dimensional EHD plumes have been considered in other papers,[8, 11] so we will restrict here to the axisymmetric case. This is the EHD flow occurring between a point and a plate electrode a distance dapart. Our concern is with the central region of this EHD flow: far from the injection source and from the stagnation point on the opposite electrode. In this region the flow is expected to have developed self similar profiles of velocity and charge because the drift term KEq in the current is much smaller than the convection one qu. We also consider that the electric field is constant and primarily xdirected, and that the effects of recirculating fluid are negligible. Similarity solutions are also possible for electric fields with a potential dependency in the flow direction.[15] With all these assumptions equations (32-35) in cylindrical coordinates (x, r) are: ∂ru ∂x +∂rv ∂r = 0 (37) u∂u ∂x +v∂u ∂r =ν r ∂ ∂r(r∂u ∂r ) + 1 ρqE (38) u∂q ∂x +v∂q ∂r = 0 (39) with boundary conditions similar to those of the thermal plume: v=∂u ∂r =∂q ∂r = 0 at r= 0,(40) u=q= 0 at r→ ∞ (41) where ris the distance from the axis of the plume. The electrical current Imust be constant at any section of the plume, then I= 2πR∞ 0qur dr is constant. Equation (39) deserves some comments. From (36) and (33) it is, ∇·j=K(E·∇)q+ (u·∇)q+Kq2 ǫ= 0.(42) Consequently, equation (39) is valid if the first term, drift, and the third term, Coulomb repulsion, are negligible compared to the second one, convection. The drift term has two contributions K(E·∇)q=KEx ∂q ∂x +KEr ∂q ∂r ∼KV d q x+KEr q a.(43) Here ais the radius of the charged zone, Vis the applied voltage and dthe distance between the point and the opposite electrode. The longitudinal contribution KEx∂q/∂x is negligible compared to u∂q/∂x, since in most EHD flows the observed velocities are usually an order of magnitude greater than the ion drift velocities.[8, 10] On the other 6 hand the radial part KErq/a is of the same order of magnitude than Coulomb repulsion, as it is readily deduced from estimating the radial field Eras the one produced by a cylinder of radius a:Er∼qa/2ǫ. In conclusion, equation (39) holds if Coulomb repulsion is negligible compared to convection. This condition can be estimated as: Kq2/ǫ (u·∇)q∼Kq2/ǫ uq/x ∼KIx πǫa2u2≪1,(44) where the electric current is I∼qπa2u. We discuss below the validity of this approximation in some reported experiments. With the following similarity variables and functions: η=IE 2πρν31/4r x1/2(45) Ψ = νxf(η) (46) q=I 2πνxg(η) (47) the resulting equations are: f′′′ η+f−1 ηf′ η′ +g= 0 (48) (fg)′= 0 (49) with boundary conditions: f′ η′ (0) = f η(0) −f′(0) 2=g′(0) = f η ′ (∞) = g(∞) = 0 (50) 1 = Z+∞ −∞ f′g dη (51) These equations for the EHD plume are the same as for the thermal plume in the limit Pr → ∞.But we know from section III that this limit does not exist, being the velocity divergent. So in order to have a meaningful physical result we need to reconsider our model. As the thermal solution shows, the velocity is finite, although large, for finite large Prandtl number. Therefore the first idea is to take into account the equivalent of the Prandtl number in the electrical problem. This number will be ν/D, being Dthe diffusion coefficient for the electric charge. The inclusion of diffusion would resolve mathematically the singularity.[15] But this is not physically sound. Typically D≃10−11m2/s, and the expansion of the charged layer due to diffusion would be D/u ≃10−11m. This value is to be compared to the radius of the injecting point, about 10 microns, that gives an order of magnitude for the radius of the charged zone. Moreover in the similarity region the drift term is negligible in comparison to the convective transport of charge (typically K∼10−9m2/V s and E∼2×106V/m, so KE ∼2×10−3m/s while u∼1m/s). What really happens in practical situations is that the plume enters into the similarity zone with a charged region of finite radius, in opposition to the assumption of a point source of charge. This is due, first, to the finite radius of the tip and, second, to the strong Coulomb repulsion near the injecting tip, that expands the charged region to a certain thickness that will be larger than the radius of curvature of the tip. Therefore in the EHD case the singularity is avoided introducing a finite radius for the charge density. The approximation will be better the thinner this radius, otherwise equation (33) has to be solved to determine the charge distribution. We have solved numerically equations (48-49) assuming that g=g0for η < a 0 for η > a (52) The velocity and the charge at the center line are, using the best fit to the numerical data: (f′ η)η=0 = 0.51 + 0.65p−ln a2, g0= 0.60a−2(53) 7 TABLE I: Velocity and temperature at the central line for different Prandtl numbers in axisymmetric plumes. As it is clear from this equation and (29) we may go from the thermal to the EHD plume just substituting Pr with a−2.The physical reason is the extent of the volume force. Outside the limit radius 1/√Pr there is not any buoyancy force in the thermal plume, neither is there any electric force outside the radius ain the EHD case. Velocity and thickness of the plume are: u=f′ η(0) IE 2πρν 1/2 (54) ∆ = 2πρν3 IE 1/4 x1/2(55) If we define an electrical Grashof number Gr = qEx3/ρν2g(0),the velocity and the thickness have the same expression as in the thermal case, equations (16) and (17). This electrical Grashof number is related to the electric current by: Gr = IEx2 2πρν3(56) an expression to be compared with its thermal counterpart, equation (15).We have summed up the results in table II. In order to provide the reader with a complete overview of thermal and EHD plumes, we have included in this table the results for the two–dimensional case.[11, 13] Some experiments have been carried out in point-plane as well as in blade to plate geometry in our group and elsewhere.[8, 9, 11, 16] In blade to plane geometry the reported measurements agree with the equations listed in table II. Results for silicone oil of viscosity 50 centistokes and density ρ= 930Kg/m3are:[11] E= 106V/m, J = 2 ×10−7A/m and x∼1cm. which gives u∼9cm/s and ∆ ∼4mm. The corresponding Reynolds number u∆/ν is 7.3. The plume is expected to be laminar, as observed. In axisymmetric cases stable jets have not yet been observed.[9, 17] We have[17] I= 4 ×10−8A, E = 3 ×105V/m and x∼30mm. Equations (54) and (55) give u∼1m/s and ∆ ∼1mm so the Reynolds number is here (νwas 27 ×10−6m2/s) close to 40. Similar values have been obtained in our laboratory. From this data it is clear that condition (44) is not fulfilled in these experiments. The role of Coulomb repulsion seems to be important in the determination of the charged region thickness, both near the tip where it expands very quickly and far from the electrodes. Whether there could be experimental conditions in which a steady laminar axisymmetric EHD plume is observable or not is a matter of investigation. V. CONCLUSION We have shown that a particular solution of boundary layer equations describes both EHD plumes and thermal plumes for very large Prandtl number. For axisymmetric plumes the velocity at the axis becomes infinity for Pr → ∞, but in any case this divergence is very weak, it goes as ln √Pr. This result is not directly useful for the EHD case. The analogue of Prandtl number in EHD is the ratio between the diffusion coefficient and viscosity. But diffusion is not an important effect in EHD. An approximate solution to EHD axisymmetric plumes can be found supposing that the charged region has a finite, although small, radius. The inverse of the charge layer cross section plays the role of an effective Prandtl number. The velocity of the EHD plume increases without limit as this thickness decreases. VI. ACKNOWLEDGMENTS This work was carried out with financial support from DGICYT (Spanish Government Agency) PB93-1182 and Fundacin Cmara. 8 Pr (f′/η)(0) θ(0) 0.7 0.938 0.482 10 1.40 5.63 1021.79 53.1 1032.13 520. 1042.41 5152. 9 TABLE II: Summary of formulae for EHD and thermal plumes. Qand Jdenote, in the plane case, heat flux per unit length and electric current per unit length, respectively. Plane plume Axisymmetric plume Thermal Electrical Thermal Electrical Grashof Gr = Qgβx3 2√2ρcpν3!4/5 Gr = JEx3 2√2ρν34/5 Gr = Qgβx2 2πρν3cp Gr = IEx2 2πρν3 Thickness ∆ = 4 Gr1/4 x∆ = 1 Gr1/4 x Velocity u=2ν xGr1/2f′(η)u=ν xGr1/2f′(η) η Temperature Charge T−Ta=ν2Grθ(η) gβx3q=ρν2Grg(η) Ex3T−Ta=ν2Grθ(η) gβx3q=ρν2Grg(η) Ex3 [1] W. J. Morgan, “Convection Plumes in the Lower Mantle,” Nature 230, 43 (1971). [2] J. A. Whitehead, “Fluid Models of Geological Hotspots,” Annual Review of Fluid Mechanics 20, 61 (1988). [3] L. Pera and B. Gebhart, “On the Stability of Laminar Plumes: Some Numerical Solutions and Experiments,” International Journal of Heat and Mass Transfer 14, 975 (1971). [4] B. Gebhart, Y. Jaluria, R. L. Mahajan, and B. Sammakia, Buoyancy-Induced Flows and Transport (Hemisphere Publishing Corporation, New York, 1988). [5] D. B. Spalding and R. G. Cruddace, “Theory of the Steady Laminar Buoyant Flow Above a Line Heat Source in a Fluid of Large Prandtl Number and Temperature–Dependent Viscosity,” International Journal of Heat and Mass Transfer 3, 55 (1961). [6] H. K. Kuiken and Z. Rotem, “Asymptotic solution for plume at very large and small Prandtl numbers,” Journal of Fluid Mechanics 45, 585 (1971). [7] K. V. Rao, B. F. Armaly, and T. S. Chen, “Mixed Convection Plumes Arising From a Thermal Point Source,” Transactions of the ASME 107, 720 (1985). [8] F. M. J. McCluskey and A. T. P´erez, “The Electrohydrodynamic Plume Between a Line Source of Ions and a Flat Plate,” IEEE Transactions on Electrical Insulation 27, 334 (1992). [9] B. Malraison, P. Atten and A. T. P´erez, “Panaches charg´es r´esultant de l’injection d’ions dans un liquide isolant par une lame ou une pointe plac´ee en face d’un plan,” Journal de Physique III France 4, 75 (1994). [10] A. Castellanos, “Coulomb-driven Convection in Electrohydrodynamics,” IEEE Transactions on Electrical Insulations 26, 1201 (1991). [11] A. T. P´erez, P. A. V´azquez, and A. Castellanos, “Dynamics and Linear Stability of Charged Jets in Dielectric Liquids,” IEEE Transactions on Industry Applications 31, 761 (1995). [12] M. Haidara and P. Atten, “Role of EHD Motion in the Electrical Conduction of Liquids in a Blade-Plane Geometry,” IEEE Transactions on Industry Applications 21, 709 (1985). [13] Y. Jaluria, Natural Convection Heat and Mass Transfer (Pergamon Press, Oxford, 1980). [14] P. Atten and R. Moreau, “Stabilit´e ´electrohydrodynamique des liquids isolants soumis ´a une injection unipolaire,” Journal de M´ecanique 11, 471 (1972). [15] A. Zhakin, “Electroconvective Jets in Liquid Dielectrics,” Mekhanica Zhigkosty i Gaza. Izvestia Academy Science USSR 6, 13 (1984). [16] T. Takashima, R. Hanaoka, I. R., and A. Ohtsubo, “I-V Characteristics and Liquid Motion in Needle-to-Plane and Razor Blade-to-Plane Configurations in Transformer Oil and Liquid Nitrogen,” IEEE Transactions on Electronic Instrumentation 23, 645 (1988). [17] B. Atten, P.; Malraison and M. Zahn, “Electrohydrodynamic Plumes in Point-Plane Geometry,” Proceedings 1993 IEEE ICDL 534 (1993).