Computer Simulations of Evaporation of Pinned Sessile Droplets

Computer Simulations of Evaporation of Pinned Sessile Droplets: Influence of Kinetic Effects. Sergey Semenov†, Victor M. Starov*†, Ramon G. ... v...
0 downloads 4 Views 3MB Size
Download PDF
Article pubs.acs.org/Langmuir

Computer Simulations of Evaporation of Pinned Sessile Droplets: Influence of Kinetic Effects Sergey Semenov,† Victor M. Starov,*,† Ramon G. Rubio,‡ and Manuel G. Velarde§ †

Department of Chemical Engineering, Loughborough University, LE11 3TU Loughborough, United Kingdom Department of Química Física I, and §Instituto Pluridisciplinar, Universidad Complutense, 28040 Madrid, Spain

‡

ABSTRACT: The aim of the current work is to present results of computer simulations, which show the influence of kinetic effects on evaporation of pinned sessile water droplets of submicrometer size placed on a heat conductive substrate. The computer simulation model also takes into account the following phenomena: influence of curvature of the droplet’s surface on saturated vapor pressure above the surface (Kelvin’s equation), the effect of latent heat of vaporization, thermal Marangoni convection, and Stefan flow inside an air domain above the droplet. The suggested model combines both diffusive and kinetic models of evaporation. The obtained results allow the characteristic droplet sizes to be estimated, where each of the mentioned above phenomena becomes important or can be neglected. are equal locally, but the intensity of each flux is always less than that predicted by either pure diffusive or pure kinetic models. The latter means that overall process of evaporation is limited either by the process of molecules transition across the liquid−gas interface or by the diffusion process into ambient gas. Computer simulations are performed using the software COMSOL Multiphysics. The dependences of total molar evaporation flux, Jc, on the radius of the contact line, L, and the value of contact angle, θ, of pinned droplets are studied below using a quasi-steady state approximation.

1. INTRODUCTION The evaporation of sessile liquid droplets plays a significant role in industrial applications such as spray cooling,1,2 inkjet printing,3 tissue engineering,4 printing of microelectromechanical systems (MEMS),5 surface modification,6,7 various coating processes,8,9 as well as biological applications.10 It is a reason why a number of theoretical and experimental investigations have been focused on investigations of this phenomenon.11−19 Investigation of evaporation of microdroplets can help understanding the influence of Derjaguin’s (disjoining/ conjoining) pressure acting in a vicinity of the apparent three-phase contact line.20−22 The aim of the presented computer simulations is to show how the evaporation of pinned sessile submicrometer size droplets of water on a solid surface differs from the evaporation of bigger millimeter size droplets. The obtained results prove the importance of kinetic effects, whose influence becomes more pronounced for submicrometer droplets. The model used below includes both diffusive and kinetic regimes of evaporation. Our model differs from a purely diffusive model, because Hertz-Knudsen-Langmuir equation23,24 is used as a boundary condition at the liquid−gas interface instead of a saturated vapor condition. The adopted model differs from a purely kinetic model of evaporation, because it also includes the vapor diffusion into the surrounding gas similar to the diffusion model. In this way the evaporation rate is controlled by both rates of vapor diffusion into ambient gas and molecules transfer across the liquid−gas interface. As a result the vapor concentration at the liquid−gas interface falls in between its saturated value and its value in ambient gas. This intermediate value of vapor concentration at the liquid−gas interface drives both transition of molecules from liquid to gas (kinetic flux) and vapor diffusion into ambient gas (diffusive flux). According to the mass conservation law these two fluxes © 2012 American Chemical Society

2. PROBLEM STATEMENT Below only small droplets (less than 1 mm) are under consideration, so that the influence of gravity is neglected. It is assumed that an axisymmetric sessile droplet forms a sharp three-phase contact line and maintains a spherical-cap shape of the liquid−gas interface due to the action of liquid−gas interfacial tension. The geometry of the problem is schematically presented in Figure 1. Below computer simulations and a quantitative comparison of the models with and without convective terms in transport

Figure 1. Axisymmetric sessile droplet on a solid substrate. θ and L are the contact angle and radius of the droplet base. Received: February 28, 2012 Published: October 9, 2012 15203

dx.doi.org/10.1021/la303916u | Langmuir 2012, 28, 15203−15211

Langmuir

Article

where I is the identity tensor. Continuity equation for fluids is also used

equations of the gas phase are undertaken. However, it can be expected that due to a small size of evaporating droplets under consideration (radius of the contact line, L, is less than 1 mm), the diffusion of vapor in the gas phase dominates the convection. The latter is confirmed by small values of both thermal and diffusive Peclet numbers: Peκ = Lu/κ < 0.05; PeD = Lu/D < 0.04, where L is the radius of the contact line, u ( 10−5 m (triangles in Figure 5). In this case the evaporation rate, Jc, is lower than that for the isothermal model and higher than that if the latent heat is included but without Marangoni convection (squares in Figure 5). For water droplets of size L < 10−5 m, the influence of the Marangoni convection is negligible and the evaporation rate, Jc, coincides with the one for a model which includes only the latent heat of vaporization (squares in Figure 5). The effect of Stefan flow in surrounding gas slightly changes the evaporation rate in the present model (circles in Figure 5) and makes it lower due to an appearance of an outward convective heat flux in the gas above the droplet. This reduces the heat flux from the ambient environment to the droplet’s surface through the gas phase. Thus the temperature of the droplet’s surface becomes lower, which reduces the evaporation rate. In our particular case this effect appeared to be much weaker than the effects of the latent heat of vaporization. The Stefan flow effect is also weaker than the effect of thermal Marangoni convection for L > 10−4 m but a bit stronger for L < 10−4 m. However, in any case the influence of the Stefan flow is small and can be neglected. It is interesting to note that the influence of the thermal effects on the kinetics of evaporation is less than 5% (according to Figure 5). Note, in our calculations we used copper as a solid support, that is, of a relatively high thermal conductivity. The latter conclusion is not valid in the case of lower thermal conductivity of the solid support. Thus, additional computer simulations show that the influence of thermal effects is more pronounced for a Teflon substrate: the evaporation flux is reduced by 19% on Teflon as opposed to 5% on copper (see Figure 6). 4.3. Influence of the Contact Angle. Everywhere above, pinned droplets with a fixed contact angle equal to 90° were used. Below we consider the dependence of the evaporation rate on the contact angle only for the isothermal model (no heat fluxes) and without the Stefan flow effect because as it was shown above it can be neglected. Such dependence for the nonisothermal model was discussed in ref 25, where nonisothermal diffusion limited evaporation was studied. As it was mentioned above, for an isothermal, pure kinetic regime of evaporation (no Kelvin’s equation), the droplet’s evaporation rate, Jc,i, should be proportional to the area, S, of the liquid−gas interface, Jc,i ∼ S. For a spherical cap shape of the droplet: S = 2πL2(1 − cos θ)/sin2 θ; that is, Jc,i ∼ (1 − cos θ)/sin2 θ, or

As already mentioned, there is a lower limit of L, where the present model can be applied. This limit is determined by the range of surface forces action, which is usually around 0.1 μm.21 Below this limit the whole droplet is in the range of surface forces action and there is no spherical part any more even on the droplet top (microdroplets according to ref 21). The evaporation process in the latter case should be substantially different from that considered above. 4.2. Influence of Thermal Effects. Computer simulations were also performed including both Kelvin’s equation and kinetic effects in the case when the thermal effects were taken into account. The latter was made to show the influence of the latent heat of vaporization, Marangoni convection and Stefan flow on the evaporation process. Droplet’s evaporation rates, Jc, were normalized using those, Jc,i (circles in Figure 3), from the isothermal model. Results are presented in Figure 5.

Figure 5. Influence of latent heat of vaporization, Marangoni convection, and Stefan flow on the evaporation rate in the case when kinetic effects and Kelvin’s equation are included into the model. Jc and Jc,i are the total molar flux of evaporation and the one in the isothermal case, respectively.

Figure 5 shows that latent heat of vaporization reduces the evaporation flux as compared to the isothermal case (that is Jc/ Jc,i < 1) in all cases considered. The reason is a temperature decrease at the droplet’s surface due to heat consumption during the evaporation process. This reduces the value of the saturated vapor pressure at the droplet’s surface and, subsequently, the rate of vapor diffusion into the ambient gas. The relative reduction of the evaporation rate (caused by the latent heat of vaporization) reaches the maximum for droplets with L ∼ 10−5 m. The latter size according to Figures 3 and 4 is in the range of diffusion limited evaporation. For smaller droplets, when the kinetic effects come into play, the influence of the latent heat on evaporation rate, Jc, decreases. The reason for that is that vapor at the outer droplet surface according to the kinetic model of evaporation is not saturated and therefore its pressure is less influenced by the local temperature but more influenced by the relative humidity of the ambient air. If we exclude kinetic effects and Kelvin’s equation from both nonisothermal (Jc) and isothermal (Jc,i) models, and include only latent heat of vaporization, then the ratio Jc/Jc,i becomes

Jc,i /Jc,i (θ = 90°) = (1 − cos θ )/sin 2 θ

(13)

On the other hand for a pure diffusive isothermal regime of evaporation this dependence is obtained by Picknett and Bexon31 15208

dx.doi.org/10.1021/la303916u | Langmuir 2012, 28, 15203−15211

Langmuir

Article

model which includes vapor diffusion, kinetic effects, and Kelvin’s equation. For comparison we selected 2 droplet sizes only: L = 10−7 and 10−6 m. The present model for the radius of the contact line L = 10−6 m results in a dependence on the contact angle close to that one for the diffusive model of evaporation. However, in the case of L = 10−7 m, the latter dependency is getting close to the kinetic model of evaporation. Note, in Figure 7, that dependencies on the contact angle are presented but not the absolute values of the evaporation rate. It is important in Figure 7 that the scale values of Jc,i(θ = 90°) are different for all curves and equal to the corresponding evaporation rate at 90°. That is why all curves go through the point (90°; 1). Let us give a qualitative explanation of dependencies in Figure 7. The total evaporation rate, Jc,i, is obtained by the integration of the local evaporation rate, jc,i, over the whole liquid−gas interface: Jc,i = ∫ Γ jc,idA (dA is an area element). Thus, the change of the local evaporation rate, jc,i, and the change of the area of the liquid−gas interface, both influence the value of the total evaporation rate, Jc,i. In a pure kinetic model of evaporation the local evaporation rate, jc,i, does not depend on the size or shape of the droplet, because it is limited by the local rate of molecules transfer across the liquid−gas interface. Therefore, in the kinetic model of evaporation the local evaporation flux, jc,i, is constant and uniformly distributed along the liquid−gas interface. Thus, the total evaporation rate is directly proportional to the area of the liquid−gas interface. The increase of the contact angle, θ, (while the contact line radius, L, is fixed) leads to the increase of the area of the liquid−gas interface and, as a consequence, to the increase of the normalized evaporation rate, Jc,i/Jc,i(θ = 90°); see kinetic model in Figure 7. For a diffusion limited evaporation, not only the area of the liquid−gas interface but also the value of the local evaporation rate, jc,i, depends on the contact angle. The local evaporation rate, jc,i, in this model is determined by the solution of the diffusion equation for vapor in surrounding gas (Laplace equation for vapor concentration). In fact, the local evaporation rate, jc,i, in the diffusive model of evaporation depends on the size of the droplet, which in turn depends on the contact angle, θ. Increase of the contact angle leads to the increase of the area of the liquid−gas interface and to the decrease of the local evaporation rate. The rate of growth of this area with contact angle is higher than the rate of decrease of the local flux, jc,i. As a result, the total evaporation rate, Jc,i, increases with the contact angle, θ, (according to the eq 14) but not that fast as for the pure kinetic model, in which the local evaporation rate is constant. 4.4. Experimental Data. To the best of our knowledge there are no experimental measurements of the contact angle, θ, and the contact line radius, L, for submicrometer (L < 1 μm) sessile droplets. However, some experimental data from the literature32,33 give measurements of θ and L for bigger (5 μm < L < 25 μm for Golovko et al.32 and 10 μm < L < 30 μm for Liu et al.33) evaporating sessile droplets of water attached to a cantilever. Comparing those measurements with the theory of nonisothermal diffusion limited evaporation, developed for sessile droplets with contact angle hysteresis,19 shows reasonably good agreement (see Figure 8). The authors of refs 32 and 33 did not present error bars for their measurements; therefore, it is possible that deviations of experimental curves from the theoretical ones lie within the experimental error bars.

Figure 6. Computer simulations: influence of a substrate material (heat conductivity of the substrate) on the evaporation rate of the droplet. Jc is the total molar flux of evaporation, when the model comprises Kelvin’s equation, kinetic effects, the effect of latent heat of vaporization, and thermal Marangoni convection. Jc,i is the total molar flux of evaporation in the isothermal case, when effects of droplet cooling and subsequently thermal Marangoni convection are switched off. Jc,i /Jc,i (θ = π /2) ⎧(0.6366θ + 0.09591θ 2 − 0.06144θ 3) θ < π /18 ⎪ ⎪ /sin θ , =⎨ ⎪(0.00008957 + 0.6333θ + 0.116θ 2 θ > π /18 ⎪ 3 4 ⎩ − 0.08878θ + 0.01033θ )/sin θ ,

(14)

where θ is in radians. In Figure 7 we present the dependence of droplet’s evaporation rate on contact angle for two above-mentioned models (eqs 13 and 14) as well as for our present isothermal

Figure 7. Dependence of normalized evaporation rate on contact angle, θ, according to diffusive (eq 14) and kinetic (eq 13) models. Symbols are calculated values according to the isothermal model of evaporation. Note, scalings Jc,i(θ = 90°) are taken from Figure 3 and are different for all plots and equal to the corresponding evaporation rate at 90°. 15209

dx.doi.org/10.1021/la303916u | Langmuir 2012, 28, 15203−15211

Langmuir

Article

Figure 8. Comparison of the experimental data on evaporation of sessile microdroplets of water32,33 with the earlier developed theory for a diffusion limited evaporation in case of contact angle hysteresis:19 (a) contact angle on nondimensional time τ̃19 and (b) nondimensional contact line radius S on nondimensional time τ.̅ 19

convection in water droplets is negligible for droplets of size L < 10−5 m. For the system considered above, Stefan flow effect appeared to be weaker than the effect of thermal Marangoni convection for L > 10−4 m but stronger for L < 10−4 m. However, in all cases its influence is small and can be neglected. According to Figure 5 the influence of latent heat of vaporization on the kinetics of evaporation is less than 5%. However, our calculations were performed for the case of copper substrate, that is, a relatively high thermal conductivity. This conclusion could not be valid in the case of a lower thermal conductivity of the solid support. The presented model can be applied for evaporation of any other pure simple liquids, not water only.

Comparison of experimental and theoretical curves in Figure 8 supports our conclusion that the evaporation rate for droplets of size above 1 μm is not affected by kinetic effects. Evaporation of considered droplets32,33 is still governed by the vapor diffusion into ambient atmosphere.

5. CONCLUSIONS Computer simulations of the evaporation of small sessile droplets of water are performed. The present model combines diffusive and kinetic models of evaporation. The effect of latent heat of vaporization, thermal Marangoni convection, and Stefan flow in the surrounding gas was investigated for a particular system: a water droplet on a heat conductive substrate (copper) in air at standard fixed conditions. Results of modeling allow the characteristic droplet sizes to be estimated when each of the mentioned above phenomena become important or can be neglected. The presented model is valid only for a droplet size bigger than the radius of surface forces action, which is around 10−7 m = 0.1 μm. That is, the data in Figure 3 for a droplet size smaller than 10−7 m are presented only to show the trend. Figure 3 shows that if the radius of the droplet base is bigger than 10−7 m then (i) deviation of the saturated vapor pressure caused by the droplet curvature (Kelvin’s equation) can be neglected, (ii) a deviation from the pure diffusion model of evaporation can be neglected for the droplet size bigger than 10−6 m, and (iii) this deviation becomes noticeable only if the droplet size is less than 10−6 m. This deviation is cause by an increasing influence of the kinetic effects at the liquid−gas interface (Hertz−Knudsen− Langmuir equation), and this theory should be applied together with the diffusion equation of vapor in the air if the droplet size is less than 10−6 m. The latter conclusions show that a consideration of evaporation of microdroplets completely covered by the surface forces action (that is less than 10−7 m) should include both deviation of the saturated vapor pressure caused by the droplet curvature and the kinetic effects. The latent heat of vaporization results in a temperature decrease at the surface of the droplet. Due to that, the evaporation rate is reduced. This effect is more pronounced in the case of diffusion limited evaporation (L > 10−5 m), when vapor pressure at the droplet’s surface is saturated and determined by local temperature. The effect of Marangoni

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel: +44(0)1509222508. Fax: +44(0)1509223923. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This research was supported by the European Union under Grant MULTIFLOW, FP7-ITN-2008-214919. The work of R.G.R. was supported in part by the Spanish Ministerio de Ciencia e Innovación through Grant FIS2009-14008-C02-01 and by ESA through Project MAP-AO-00-052. Both V.M.S. and R.G.R. recognise support from the European Space Agency (PASTA project).

■

REFERENCES

(1) Sodtke, C.; Stephan, P. Spray cooling on micro structured surfaces. Int. J. Heat Mass Transfer 2007, 50, 4089−4097. (2) Cheng, W.-L.; Han, F.-Y.; Liu, Q.-N.; Zhao, R.; Fan, H.-L. Experimental and theoretical investigation of surface temperature nonuniformity of spray cooling. Energy 2011, 36, 249−257. (3) Peng, D.; Luhai, L.; Wen, Z.; Xian, L.; Xuwei, H. Study on the printing performance of coated paper inkjet ink. Adv. Mater. Res. 2011, 174, 358−361. (4) Campbell, P. G.; Weiss, L. E. Tissue engineering with the aid of inkjet printers. Expert Opin. Biol. Ther. 2007, 7 (8), 1123−1127. (5) Fuller, S. B.; Wilhelm, E. J.; Jacobson, J. M. Ink-jet printed nanoparticle microelectromechanical systems. J. Microelectromech. Syst. 2002, 11, 54−60.

15210

dx.doi.org/10.1021/la303916u | Langmuir 2012, 28, 15203−15211

Langmuir

Article

(6) Haschke, T.; Wiechert, W.; Graf, K.; Bonaccurso, E.; Li, G.; Suttmeier, F. T. Evaporation of solvent microdrops on polymer substrates: from well controlled experiments to mathematical models and back. Nanoscale Microscale Thermophys. Eng. 2007, 11, 31−41. (7) Pericet-Camara, R.; Bonaccurso, E.; Graf, K. Microstructuring of polystyrene surfaces with nonsolvent sessile droplets. ChemPhysChem 2008, 9, 1738−1746. (8) Karlsson, S.; Rasmuson, A.; Björn, I. N.; Schantz, S. Characterization and mathematical modelling of single fluidised particle coating. Powder Technol. 2011, 207, 245−256. (9) Ehrhard, P.; Davis, S. H. Non-isothermal spreading of liquid drops on horizontal plates. J. Fluid Mech. 1991, 229, 365−388. (10) Kim, J. H.; Shi, W.-X.; Larson, R. G. Methods of stretching DNA molecules using flow fields. Langmuir 2007, 23, 755−764. (11) Deegan, R. D.; Bakajin, O.; Dupont, T. F.; Huber, G.; Nagel, S. R.; Witten, T. A. Contact line deposits in an evaporating drop. Phys. Rev. E 2000, 62, 756−765. (12) Guena, G.; Poulard, C.; Voue, M.; Coninck, J. D.; Cazabat, A. M. Evaporation of sessile liquid droplets. Colloids Surf. A 2006, 291, 191−196. (13) Girard, F.; Antoni, M.; Sefiane, K. On the effect of Marangoni flow on evaporation rates of heated water drops. Langmuir 2008, 24, 9207−9210. (14) Hu, H.; Larson, R. G. Marangoni effect reverses coffee-ring depositions. J. Phys. Chem. B 2006, 110, 7090−7094. (15) Sefiane, K.; Wilson, S. K.; David, S.; Dunn, G. J.; Duffy, B. R. On the effect of the atmosphere on the evaporation of sessile droplets of water. Phys. Fluids 2009, 21, 062101. (16) Ristenpart, W. D.; Kim, P. G.; Domingues, C.; Wan, J.; Stone, H. A. Influence of substrate conductivity on circulation reversal in evaporating drops. Phys. Rev. Lett. 2007, 99, 234502. (17) Bhardwaj, R.; Fang, X.; Attinger, D. Pattern formation during the evaporation of a colloidal nanoliter drop: a numerical and experimental study. New J. Phys. 2009, 11, 075020. (18) David, S.; Sefiane, K.; Tadrist, L. Experimental investigation of the effect of thermal properties of the substrate in the wetting and evaporation of sessile drops. Colloids Surf. A: Physicochem. Eng. Aspects 2007, 298, 108−114. (19) Semenov, S.; Starov, V. M.; Rubio, R. G.; Agogo, H.; Velarde, M. G. Evaporation of sessile water droplets: Universal behaviour in presence of contact angle hysteresis. Colloids Surf., A 2011, 391, 135− 144. (20) Moosman, S.; Homsy, G. M. Evaporating menisci of wetting fluids. J. Colloid Interface Sci. 1980, 73, 212−223. (21) Starov, V.; Velarde, M.; Radke, C. Dynamics of wetting and spreading. Surfactant Sciences Series; Taylor&Frances: London, 2007; Vol. 138. (22) Ajaev, V. S.; Gambaryan-Roisman, T.; Stephan, P. Static and dynamic contact angles of evaporating liquids on heated surfaces. J. Colloid Interface Sci. 2010, 342, 550−558. (23) Kryukov, A. P.; Levashov, V. Yu.; Sazhin, S. S. Evaporation of diesel fuel droplets: kinetic versus hydrodynamic models. Int. J. Heat Mass Transfer 2004, 47, 2541−2549. (24) Sazhin, S. S.; Shishkova, I. N.; Kryukov, A. P.; Levashov, V. Yu.; Heikal, M. R. Evaporation of droplets into a background gas: kinetic modelling. Int. J. Heat Mass Transfer 2007, 50, 2675−2691. (25) Semenov, S.; Starov, V. M.; Rubio, R. G.; Velarde, M. G. Instantaneous distribution of fluxes in the course of evaporation of sessile liquid droplets: computer simulations. Colloids Surf., A 2010, 372, 127−134. (26) Yi, S.; Yue, C. Y.; Hsieh, J.-H.; Fong, L.; Lahiri, S. K. Effects of oxidation and plasma cleaning on the adhesion strength of molding compounds to copper leadframes. J. Adhesion Sci. Technol. 1999, 13, 789−804. (27) Krahl, R.; Adamov, M.; Avilés, M. L.; Bänsch, E. “A model for two phase flow with evaporation”. Weierstraß-Institut für Angewandte Analysis und Stochastik (WIAS), Preprint ISSN 0946-8633.

(28) Fuller, E. N.; Schettler, P. D.; Giddings, J. C. A new method for prediction of binary gas-phase diffusion coefficients. Ind. Eng. Chem. 1966, 58 (No 5), 19−27. (29) Galvin, K. P. A conceptually simple derivation of the Kelvin equation. Chem. Eng. Sci. 2005, 60, 4659−4660. (30) Bligh, P. H.; Haywood, R. Latent heat - its meaning and measurement. Eur. J. Phys. 1986, 7, 245−251. (31) Picknett, R. G.; Bexon, R. The evaporation of sessile or pendant drops in still air. J. Colloid Interface Sci. 1977, 61, 336−350. (32) Golovko, D. S.; Butt, H.-J.; Bonaccurso, E. Transition in the evaporation kinetics of water microdrops on hydrophilic surfaces. Langmuir 2009, 25, 75−78. (33) Liu, C.; Bonaccurso, E. Microcantilever sensors for monitoring the evaporation of microdrops of pure liquids and mixtures. Rev. Sci. Instrum. 2010, 81, 013702.

15211

dx.doi.org/10.1021/la303916u | Langmuir 2012, 28, 15203−15211