Evaporation of Nanodroplets on Heated ... - American Chemical Society

Article pubs.acs.org/Langmuir

Evaporation of Nanodroplets on Heated Substrates: A Molecular Dynamics Simulation Study Jianguo Zhang, Frédéric Leroy, and Florian Müller-Plathe* Eduard-Zintl-Institut für Anorganische und Physikalische Chemie and Center of Smart Interfaces, Technische Universität Darmstadt, Petersenstrasse 22, D-64287 Darmstadt, Germany ABSTRACT: Molecular dynamics simulations of Lennard-Jones particles have been performed to study the evaporation behavior of nanodroplets on heated substrates. The influence of the liquid−substrate interaction strength on the evaporation properties was addressed. Our results show that, during the temperature-raising evaporation, the gas is always hotter than the droplet. In contrast to the usual experimental conditions, the droplet sizes in our simulations are in the nanometer scale range and the substrates are ideally smooth and chemically homogeneous. As a result, no pinning was observed in our simulations for substrates denoted either hydrophilic (contact angle θ < 90°) or hydrophobic (contact angle θ > 90°). The evaporative mass flux is stronger with increasing hydrophilicity of the substrate because the heat transfer from the substrate to the droplet is more efficient for stronger attraction between the solid and the droplet. Evaporation and heat transfer to the gas phase occur preferentially in the vicinity of the three-phase contact line in the hydrophilic system. However, in the case of a hydrophobic substrate, there is no preferential location for mass and heat fluxes. During the whole evaporation process, no pure behavior according to either the constant-angle or the constant-radius model was found; both the contact angle and contact radius decrease for the droplets on hydrophilic and hydrophobic substrates alike. extensively. For example, Langmuir34 pointed out that the droplet evaporation was hindered by the presence of the solid substrate underneath the droplet. Picknett and Bexon7 developed equations to describe the evaporation of droplets deposited on substrates. Birdi et al.35 studied the evaporation of water droplets placed on smooth solid substrates and found that the evaporation rates are linear with time for contact angles θ < 90° (e.g., water on glass), whereas the evaporation rates behave nonlinearly for contact angles θ > 90° (e.g., water on Teflon). Rowan et al.26 found that the evaporation rates are proportional to the droplet height rather than spherical radius. Recently Erbil et al.11 reported a linear relationship between V2/3 (where V is the volume of the droplet) and the evaporation time t in their experimental study of the evaporation of n-nonane, 1-butanol, n-octane, and toluene drops on a poly(tetrafluoroethylene) (PTFE) surface. Similarly, for hydrophobic substrates, a linear relationship between m2/3 (with m the mass of the droplet) and time was observed by Golovko et al.14 and Xu et al.18 Such a nonlinear evolution of the residual droplet volume was also found by Nguyen et al.29 through both theoretical and experimental analysis. In contrast, a linear relationship between the droplet mass and the evaporation time was reported by Hu et al.12 and Golovko et al.14 on hydrophilic substrates. Note that although the surfaces on which experiments are carried out are claimed to be very carefully prepared, they always exhibit surface (chemical and/or

I. INTRODUCTION The evaporation of a droplet is a spontaneous phenomenon which occurs when the surrounding vapor pressure is not saturated. Owing to its role in many industrial applications, such as cooling, printing, washing, and coating, the study of the evaporation of droplets has attracted a lot of attention. The topic has recently become even more popular because of its new applications in technology and biology, such as inkjet printing,1,2 micro/nanofabrication,3 DNA stretching and DNA mapping,4,5 and the manufacturing of novel optical and electronic materials.6 A significant amount of experimental7−22 and theoretical12,15,23−34 work has been carried out to understand the relationship between the mechanisms of droplet evaporation (e.g., how does the contact line or contact angle change, and what factors will influence the evaporation rate?) and its applications. Maxwell proposed the first theory of the evaporation of a free droplet surrounded by gas about 100 years ago,23 in which he assumed that the evaporation process is controlled by the mass transfer. Since then, theories of droplet evaporation have been developed and refined. For instance, on the basis of the dynamics of the contact line and contact angle, three different evaporation models have been proposed: (i) in the constant-contact-radius model,7−9,11,13,14,17−19,21,25,26,33,35 the three-phase contact line (TPCL) is pinned and the contact angle decreases during the evaporation (“pinning”), (ii) in the constant-contact-angle model,7,9,13,14,17−20,35 the contact angle remains constant over the evaporation while the TPCL moves (“depinning” or “shrinking”), and (iii) in the mixed mode,7,21 both the contact angle and the position of the contact line change. Besides, the evaporation rate has also been investigated © XXXX American Chemical Society

Received: April 30, 2013 Revised: June 27, 2013

A

dx.doi.org/10.1021/la401655h | Langmuir XXXX, XXX, XXX−XXX

Langmuir

Article

physical) inhomogeneities,36 which are considered to contribute to the droplet contact line pinning during evaporation. Therefore, the evaporation of sessile droplets is generally observed to take place in two stages of the TPCL dynamics,29 namely, pinning and depinning. Orejon et al.37 stressed that the dynamics of the contact line is dictated by the competition between pinning forces and depinning forces. The pinning forces arise from the contact line being anchored to the substrate due to surface heterogeneities, whereas the depinning forces are the result of the deviation of the equilibrium droplet profile due to evaporation. The pinning phenomena were also recently observed by Shin et al.21 and Grandas et al.38 In the former study of water droplet evaporation on an octadecyltricholorosilane (OTS) substrate, the pinning of the contact line was observed and claimed to probably be the result of the imperfections of the substrate surface. In the latter, the pinning of water on PTFE was also observed and explained by heterogeneities on the substrate. Anderson et al.39 suggested that the base diameter recedes while the contact angle remains constant (i.e., contact-angle model) during evaporation on ideal surfaces. However, the perfect quality of “ideal surface” cannot easily be achieved by the current experimental methods. Therefore, the evaporation mechanisms of droplets on ideal surfaces remain unclear. The importance of heat transfer in the evaporation process has been reported recently.10,40 Feng and Ward10 found through a series of experimental works that the temperature distribution across the interface between the liquid and the vapor phase of water was not continuous; the vapor temperature is higher than the liquid temperature. These results cannot be explained by the commonly accepted theories of evaporation, in which it is usually considered that the mass transfer rather than the heat transfer and energy balance dominates evaporation. Hence, the heat transfer is usually ignored in modeling the evaporation process, such as in the classical Hertz−Knudsen or Schrage model of evaporation.28 These simplified models can still be used to predict the evaporation rate, but will be questionable when applied to small (such as nanometer-sized) droplets. Actually, as reported in the simulation work by Holyst and Litniewski,32 the evaporation of a nanoscale droplet is limited by heat transfer and energy balance, rather than by mass transfer. They also found that whether the temperature at the liquid−vapor interface is continuous or discontinuous depends on the density ratio of liquid to vapor; the temperature is continuous when the value of the ratio is on the order of 10 and discontinuous otherwise. The local heat and mass fluxes in the vicinity of the threephase contact line area dominate the evaporation process, since most of the mass flux and heat flux occur in this region when the contact angle is less than 90°. On one hand, this nonuniform distribution of the evaporation flux along the droplet surface has been theoretically studied by Deegan et al.25,33 and Nguyen et al.29 or computationally studied by Hu et al.12 through a finite element method. On the other hand, Stephan and co-workers showed experimentally41 that about 50% of the heat was transferred within the microregion at the three-phase contact line. They42 recently showed that the local heat fluxes in the three-phase contact line area are about 6 times larger than the mean input heat fluxes at the foil. Moreover, when studying the evaporation of a sessile droplet, the wettability of the substrate, which is closely linked to its surface properties, such as the surface roughness and surface energy, should also be taken into account. The wetting behavior

actually plays an important role in evaporation. As reported by Sobac and Brutin,43 the dynamics and kinetics of the evaporation of water and simple organic molecular droplets are completely different if the droplets are deposited on the substrates with the surface properties modified from hydrophilic to hydrophobic. The evaporation time increases by 58% when the contact angle increases from 63° to 135°. Hence, the investigation of the processes in the vicinity of the three-phase contact line and under the influence of the wettability properties of the substrates appears to be crucial to understand the evaporation mechanism at the microscopic scale. Although significant progress has been achieved in the understanding of the evaporation of a macroscopic droplet, several problems of interest remain to be investigated in detail. A question often considered is how the evaporation of a nanometer-sized droplet will compare to the evaporation of macroscopic droplets (>100 μm in diameter). Usually, the evaporation laws for macroscopic droplets cannot be directly used to study nanometer-sized droplets, since specific dominant effects may have to be considered.14 For instance, as mentioned above, Holyst and Litniewski,32 reported that the evaporation of a nanoscale droplet is limited by heat transfer and energy balance, rather than by mass transfer. The latter mechanism is considered to dominate the evaporation of macroscopic droplets. In the case of small droplets, the modified Young’s equation44 cos(θ) = cos(θ∞) − τ/γLVr shows that size effects may have to be taken into account. In this equation, the line tension τ refers to the excess free energy per unit length of the three-phase contact line, γLV is the liquid−vapor interfacial tension, and θ and θ∞ are the apparent contact angle and contact angle of a macroscopic droplet (the spherical radius r is infinite). If the line tension τ is on the order of 10−11 J/m,36 it will affect the contact angle when the droplet is in the range of a few nanometers in diameter. To verify the applicability of the models for macroscopic droplets and to describe microscopic droplets, a few studies on subnanoliter (∼100 μm in diameter),14,17picoliter (∼50 μm in diameter),14,19 and femtoliter (∼1 μm in diameter)20 droplets have been performed recently by using advanced techniques. For instance, Furuta et al.17 found that, due to the effect of line tension, subnanoliter (∼100 μm in diameter) droplets have an opposite trend in contact angle change during evaporation compared with microliter droplets when deposited on either hydrophobic or hydrophilic substrates. Arcamone et al.20 reported that macroscopic models remain valid to describe droplets having a volume in the femtoliter range. These authors also suggested that it would be interesting to explore even smaller droplet volumes (i.e., diameters in the range of nanometers) where macroscopic models are expected to reach their limits. However, nanometer-sized droplets for which the line tension is expected to play a significant role are still inaccessible by experiments, but can be studied at the molecular level by molecular dynamics (MD) simulations. Indeed, such a methodology has proved to be a very useful tool to study evaporation properties of nanodroplets. To date a few simulation works have been performed. For instance, Long et al.45 studied the evaporation rate of an argon droplet surrounded by its vapor at subcritical conditions. Consolini et al.46 reported the time evolution of droplet shapes during the evaporation under both subcritical and supercritical conditions. Matsumoto47 studied the transport properties occurring at the liquid−gas interface using the MD method. Yang and Pan48 studied a thin water film using MD. Therefore, most of these B


Langmuir

Article

These three values of ε tune the substrate from weakly to strongly adsorbing, which, for compactness, we call in the remainder of the paper hydrophobic and hydrophilic, even though we are dealing with a Lennard-Jones fluid rather than water. From our results, the equilibrium contact angle is 114° for E4, i.e., typical for a hydrophobic interface, 88° for E6, namely, slightly hydrophilic, and 75° for E7, meaning hydrophilic, at a temperature of 0.67. In our simulations, the NVT ensemble was used, where the Berendsen thermostat49 was chosen to control temperature with a coupling time of 0.25τ. Note that the thermostat was only applied to the substrate atoms to mimic the heating process. Periodic boundary conditions in the three directions of space were used in all simulations. The simulations were done using LAMMPS.50 The time step and the cutoff distance were set to 0.0025τ (5 fs) and 4.4σ (1.5 nm), respectively. In the initial configurations, there is no gas phase, as all the fluid particles are placed on a lattice at a density close to the liquid density. The simulations were run on each surface of E4, E6, and E7 at a temperature T* = 0.67. Each system relaxed to its equilibrium configuration where a droplet surrounded by a vapor phase was formed. The temperatures of both the droplet and the gas phase were equal to the target value, and the droplet sizes no longer changed. Therefore, steady droplets surrounded by their vapor were formed. This typically took about 7500τ. The droplets were then heated by suddenly increasing the temperature of the substrate to T* = 0.83, yielding droplet evaporation. The whole evaporation was classified into two phases: the temperature-raising evaporation process, in which evaporation occurs, while the temperature of the droplet increases, and the subsequent constant-temperature evaporation process, where the droplet temperature has reached the target. The evaporation rate in phase 2 has been increased by suddenly increasing the box to twice its size in the z direction after 40000τ. At least 200 million time steps (50000τ) were performed for each simulation to collect enough data for further analysis of the evaporation properties. The sizes of the droplets and their variation upon evaporation were obtained by means of a cluster analysis. Two atoms were considered to be part of the same cluster if the distance between them was smaller than 1.5σ.

works studied evaporation of liquids in vacuum or within a gaseous environment. There is a lack of simulation work where the evaporation mechanism of surface-resident sessile droplets is studied. Such works could actually be very interesting and challenging, since when studying the evaporation of a droplet in contact with a surface, the complexity increases significantly due to the heat transfer from the substrate to the droplet and the appearance of the three-phase contact line. In this study, the evaporation behavior of nanodroplets deposited on an ideal (i.e., flat, completely smooth, and chemically homogeneous) substrate was investigated by MD simulations. The effect of the wetting properties of the substrate on the evaporation was studied. The following scientific questions have been addressed: Which evaporation model (dynamics of the contact line and contact angle) can be used to describe the evaporation behavior of a nanometer-sized droplet on an ideal substrate? Which evaporation scaling law do nanometer-sized droplets follow? Do most of the heat and mass fluxes still occur at the threephase contact line for evaporation on hydrophilic substrates when the droplet is reduced to a few nanometers in diameter? How do the local heat and mass fluxes behave along the droplet surface for the case of a droplet on hydrophobic substrates, which is rarely studied in both experiments and theories? This paper is organized as follows. In section II, we present the details about the model and methods used in our simulations. In section III, we show our main results, including the evaporation rate, temperature evolution of the droplet and gas phase, contact angles and contact radius, and local heat flux and mass flux as well as the evolution of the droplet shape. Finally, in section IV, we present the conclusions of this study.

II. MODEL AND METHODS MD simulations were carried out to study the evaporation of a droplet on heated substrates. There were three phases in the systems, namely, the liquid droplet surrounded by its gas and the solid substrate. Lennard-Jones atoms with the parameters of argon were used to model the fluid phases, because this atomic fluid has extensively been studied both experimentally and theoretically. To simplify the system, Lennard-Jones atoms were employed to model the solid substrate too. The total number of argon atoms for the liquid and gas phases was 27 000. These atoms were initially distributed on a cubic lattice with a density close to the liquid density. The solid substrate was modeled as an fcc lattice containing eight layers of 5000 atoms each. It was placed at the bottom of the simulation box underneath the lattice containing the liquid atoms. The simulation box size was 32.4 nm in the x and y directions and 50 nm in the z direction. The 12−6 Lennard-Jones potential (U(rij)) was used to model the interactions between the atoms in the systems:

⎡⎛ ⎞12 ⎛ ⎞6 ⎤ σ σ U (rij) = 4ε⎢⎢⎜⎜ ⎟⎟ − ⎜⎜ ⎟⎟ ⎥⎥ rij ⎠ rij ⎠ ⎝ ⎝ ⎦ ⎣

III. RESULTS AND DISCUSSION A. Density and Interface Thickness. The mass-density distribution along the droplet’s radial direction has been computed for each system under the assumption that the droplet was a cap of a perfect sphere, and the position of the geometric center of the sphere was derived from fitting the droplet isodensity contour by a circle. As an example, we show this density distribution for the E4 system at T* = 0.83 (Figure 1). Larger fluctuations around the small distances can be explained by the small statistical volume. Lower fluctuation is observed at large distances. The mass density decreases sharply

(1)

where rij is the distance between atoms i and j and σ and ε are the collision diameter and potential well, respectively. In our simulations, σ = 0.3405 nm and ε = 0.992 kJ/mol are the parameters for argon and σ = 0.4085 nm and ε = 9.92 kJ/mol are used to describe the solid atoms. In the following text, the parameters σ, ε, and m (= 6.63 × 10−26 kg) for argon are used to derive the reduced units of the other quantities. The unit of ε/kB = 119.8 K is used for temperature. The time unit is τ = (mσ2/ε)1/2. For the cross-interactions between the fluid and the solid atoms, a value of σ* = 1.1 was used. Three different values of ε were considered for the solid−fluid interaction to study the effect of the substrate wetting properties on the evaporation process. The systems were named E4, E6, and E7 in reference to the values of 0.4, 0.6, and 0.7 for the corresponding cross-interaction parameter ε*.

Figure 1. Mass-density distribution of the fluid particles with respect to the distance to the center of a spherical shape approximating the droplet for the E4 system at T* = 0.83. The empty circles refer to the simulation data. The solid line represents the fit to these data. C


Langmuir

Article

Figure 2. (a) Time evolution of the liquid mass densities in the E4, E6, and E7 systems during evaporation on a substrate at T* = 0.83. The horizontal dashed line shows c(t) = 1/e, and the vertical dashed lines show the corresponding time points (marked as a, b, and c). (b) Time evolution of the gas−liquid interface thickness ds of the E4, E6, and E7 systems during evaporation at T* = 0.83.

Figure 3. Time evolution of the temperature of the (a) gas phase and (b) liquid phase for E4 (red line), E6 (green line), and E7 (blue line) at a substrate temperature T* = 0.83. The black line is the substrate temperature.

Figure 2b shows the liquid−vapor interface thickness ds with respect to the simulation time. All the curves have an asymptotic behavior and converge to a value of 3.49 for E4,, about 2.99 for E6, and 2.87 for E7. The increase in thickness is driven by the increase of the droplet temperature. Theoretically, these three systems should have the same interface thickness after they reach the same temperature. The difference in the thickness of the three systems from our simulations may come from the fluctuations of the droplet shape in the evaporation as well as from the error on the fit presented in Figure 1. The difference (around 0.5σ) between E4 and E7 is minor, and all these values of thickness are sensibly within the same range when compared with the simulated values of 4.7 and 7.4 obtained in the simulations of a thin film by Wu and Pan52 for temperatures of 0.92 (110 K) and 1.04 (125 K), respectively. B. Temperature Evolution of the Droplet and Vapor. To study the heat transfer, the time dependence of temperature in the vapor and in the droplet was computed. It is shown in Figure 3a that the vapor temperature for E6 and E7 increases quickly from 0.67 to the target temperature of 0.83, and a relatively slower increase is observed for E4. In contrast, it takes a longer time for the droplets to reach equilibrium (Figure 3b). One reason is the larger contact area between the gas phase and substrate (15000σ2) compared to the droplet (2700σ2). More importantly, the molecules with higher kinetic energy evaporate with higher probability from the liquid−vapor interfacial region; thus, the gas-phase temperature reaches the target value faster than the droplet temperature. From Figure 3b, it can also be seen that the droplet temperature of E7 increases faster than the temperature in E6 and E4 and the droplet temperature in E6 converges faster than that in E4. Finally, they all converge to the temperature of the heated surface at 25000τ for E4, 12000τ for E6, and 7500τ for E7. This means that hydrophilic substrates (E6 and E7) transfer heat more

at a distance of approximately 17.6, indicating the liquid−vapor interface area. Four quantities, liquid mass density (ρl), vapor mass density (ρg), spherical radius (r0) of the droplet, and liquid−vapor interface thickness (ds), were derived through fitting the curve to a hyperbolic tangent equation: ρ=

⎡ 2(r − r0) ⎤ 1 1 (ρl + ρg ) − (ρl + ρg ) tanh⎢ ⎥ 2 2 ds ⎦ ⎣

(2)

After the calculation of these quantities for each time point in the evaporation, they are used to study the time evolution properties of evaporation. The mass-density curves for E6 and E7 are very similar to that of E4 and are not shown. The droplet densities were obtained by fitting eq 2 mentioned above. They are shown for the three systems in Figure 2a. The quick decrease of the mass-density curves at the beginning of the simulation is followed by a slower asymptotic decay. Finally, the three curves converge to the same value of ρ* = 0.775, which is very close to the experimental value51 of ρ* = 0.789 at T* = 0.83. The decrease is the fastest for system E7, followed by systems E6 and E4. This can be quantitatively described by the time evolution of the normalized density c(t) = (ρ(t) − ρ(∞))/(ρ(0) − ρ(∞)), where ρ(0), ρ(t), and ρ(∞) are the mass densities at time points 0, t, and ∞, respectively. Here, ρ(∞) is the average value over the last 5000τ of the curves in Figure 2a. Relaxation times are obtained as the times when c(t) = 1/e, which are 900τ, 1230τ, and 4000τ for E7, E6, and E4, respectively (marked as a, b, and c in Figure 2a). The evaporation rates become faster in the order E4 < E6 < E7. This arises from the difference in the heating rate of the droplet temperatures, which is shown in Figure 3b and will be discussed below in more detail. After the droplet temperature of the three systems has been fully thermalized to T* = 0.83, the mass densities of all the droplets converge to the same value. D


Langmuir

Article

liquid−vapor interface; because the evaporation happens in this area, heat is transferred from the droplet to the vapor by evaporating molecules. From the energy balance point of view,54 this heat would lead to a decrease of the interface temperature. However, this locally lower temperature area is obviously not observed in our simulation (Figure 4). The local temperatures of the interface areas (marked as vertical dashed blue lines) are not conspicuously lower than that of the bulk of the droplet. Possible reasons are that a small lowering of the temperature disappears in the thermal fluctuations or that the thermal conductivity is high enough to quickly replenish the lost heat. However, it still can be seen clearly in Figure 4b−d that the interface temperatures are lower than those of the respective gas phases. C. Droplet Size and Mass Flux. The droplet size was calculated by cluster analysis, as mentioned in section II. Two molecules were considered in the same cluster if their distance was less than 1.5σ. The initial total number of argon molecules in the box is 27 000; after equilibrium at T* = 0.67, a steady droplet was formed and the number of molecules inside the droplet was about 24 000 for E4 and 24 200 for E6 and E7. When a sudden increase in the temperature of the substrate is applied, the droplet size decreases quickly for all three systems as shown in Figure 5a. The decrease in E7 is faster than that in E6, and both are faster than that in E4; this can be explained by the heat transfer from the substrate to the droplet being more efficient on more hydrophilic substrates. As a result, the droplet temperature appears higher on the more hydrophilic substrate system at the beginning of the simulation (Figure 3b). Higher temperature and quicker heat transfer cause a faster mass loss.

efficiently than the hydrophobic substrate (E4 system). This can be qualitatively explained by the fact that, on one hand, the droplet has a larger contact area for E6 and E7 than for E4 and that, on the other hand, having stronger liquid−solid interactions results in higher interfacial thermal conductance. As found by Shenogina et al.,53 the interfacial conductance is directly proportional to the work of adhesion (or proportional to 1 + cos θ, where θ is the contact angle) and increases with increasing hydrophilicity. Therefore, the efficiency of heat transfer to the droplet follows the order E7 > E6 > E4. As was mentioned above, at the beginning of the simulation, the temperature increase in the droplets is slower than that in the vapor, suggesting that the average vapor temperature is larger than that of the corresponding droplet for all three systems for a certain period of time. Therefore, we assume that there is a temperature gap between the vapor and liquid phases before the droplet temperature is completely thermalized to the substrate temperature. To confirm this assumption, the temperature distribution along the radial direction of the droplet was computed. In Figure 4, five curves corresponding to the temperature distributions at 0, 2500τ, 5000τ, 10000τ, and 25000τ in the

Figure 4. Temperature distributions along the radial direction of the droplet in the E4 system for five different time points, (a) 0τ, (b) 2500τ, (c) 5000τ, (d) 10000τ, (e) 25000τ, at a substrate temperature T* = 0.83. The horizontal dashed red lines represent the mean values of the droplet and vapor temperatures, and the vertical dashed line represents the position of the liquid−vapor interface.

droplet of system E4 are presented. The temperature is shown as a function of the distance to the geometric center of the droplet. The larger temperature fluctuations in the vapor phase are explained by fewer degrees of freedom per unit volume in this phase. The average temperatures for both the droplet and vapor are initially 0.67, as can be observed from the curve at 0τ in Figure 4a. The average temperatures (as indicated by the horizontal dashed red lines in Figure 4b−d) of the droplet at 2500τ, 5000τ, and 10000τ are about T* = 0.753, 0.784, and 0.812, respectively, and are lower than the corresponding gasphase temperature of 0.83. The temperature gaps between the droplet and the vapor are thus 0.077, 0.046, and 0.018, respectively, becoming smaller during evaporation. At 25000τ in Figure 4e, the droplet temperature reaches the value of the substrate temperature of 0.83. Hence, both the droplet and the vapor have the same temperature, and the temperature gap no longer exists. Similar temperature gaps were also found for the E6 and E7 systems. Note that the temperature difference between the vapor and the liquid should have a contribution to the droplet evaporation. It has been reported by Fang and Ward10 that evaporation is limited by heat transfer due to the temperature gap between the droplet and its vapor. Theoretically, the temperature should be discontinuous at the

Figure 5. (a) Number of atoms inside the droplet as a function of time for E4, E6, and E7. The vertical dashed line shows the time point for the enlargement of the box size. The inset represents the variations of N2/3 as a function of time for E4, E6, and E7. N is the number of atoms inside the droplet. The blue dashed lines represent linear relationships. (b) Mass fluxes per unit liquid−vapor interface area per unit time in E4, E6, and E7 with respect to time at a substrate temperature T* = 0.83 for the first 40000τ. All the curves were smoothed by the adjacent-averaging method. (c) Mass fluxes of E4, E6, and E7 with respect to time at a substrate temperature T* = 0.83 starting after the gas-phase volume was suddenly doubled (40000τ). The dashed lines are fits to the exponential function y = y0 + AeBx. (d) Temperature dependence of evaporation (labeled as “Evap”) and condensation (labeled as “Condens”) rates (in atoms leaving/joining the droplet per 10τ) for the E4, E6, and E7 systems. All curves were obtained by means of the fits to the exponential function y = y0 + AeBx. E


Langmuir

Article

temperature of 0.67 and the gas pressure would be far from saturation at the higher temperature T* = 0.83 when the gas temperature went through a quick increase at the beginning of the simulation (cf. Figure 3a). As the gas pressure saturates, evaporation slows. As a result, a sharp peak appears in the mass flux curves. The peaks are ordered in height as E7 > E6 > E4, i.e., in the same order as the temperature increase (Figure 3b). Beyond their maxima, the curves decrease in the same order. Finally, it can be seen in Figure 5b that the three mass flux curves seem to converge asymptotically to zero. Figure 5c shows the mass flux for the evaporation after the box volume was doubled. The data were sampled with a time interval of 50τ (∼100 ps). There are large fluctuations in all curves. Negative values mean that the droplet size increases temporarily due to condensation being faster than evaporation. The appearance of condensation next to evaporation indicates that the systems approach equilibrium. To make the trend of the mass flux variation clearer, an exponential function (y = y0 + AeBx) was fitted to the data in all cases. Such an analytical shape gave satisfactory results. An almost constant mass flux is first observed for all systems, and the mass fluxes are ordered following E7 > E6 > E4. The mass fluxes of E7 and E6 are more than twice as large as that of E4. This indicates that, although in all systems the temperature is the same (0.83), and although E4 has a larger liquid−vapor interface area than E6 and E7, the mass flux from the hydrophilic systems is still more than twice faster than that from the hydrophobic system. The reason is the more efficient heat transfer from the hydrophilic substrates. The substrate wettability thus has a strong effect on the total mass flux. Similar phenomena have been shown by experiments18,21 recently. Shin et al.21 reported that sessile droplets with large contact angles (θ > 90° on hydrophobic substrates) have a longer lifetime (duration of evaporation) than those with smaller contact angles (θ < 90° on hydrophilic substrates). Xu et al.18 found that the evaporation rate decreases when the hydrophobicity of the surface is increased in the study of the evaporation of water droplets on superhydrophobic substrates. Therefore, our results are consistent with these experimental findings. Note that, toward the end of evaporation when the droplets have almost vanished, all three mass fluxes increase amid fluctuations which are also increasing. A number of reasons for a change of rate and mechanism for the evaporation of the last traces of liquid can be speculated about, but we have not analyzed this in detail. The number of molecules in the droplet corresponding to the starting time point of mass flux increase in Figure 5c is about 3000 for E4 and around 2000 for E6 and E7. Evaporation always takes place simultaneously with condensation. Net evaporation occurs when the evaporation rate is higher than the condensation rate and stops when they are equal. Both evaporation and condensation rates, which are defined as the number of atoms going out from or coming into the droplet per unit time, respectively, have been calculated. The curves in Figure 5d were obtained by fitting the data to an exponential function (y = y0 + AeBx) to reduce the large noise. All the curves increase first, indicating that both the evaporation and condensation rates increase as the droplet temperatures increase (Figure 3b). As expected, the evaporation rate is always larger than the condensation rate before they reach equilibrium. It can been seen that the difference between the evaporation rate and condensation rate first increases and then decreases, which explains part of the maxima in the mass flux (Figure 5b). Finally, the curves of evaporation and condensation of the same

As can also be seen in Figure 5a, the fast decay in the three curves is followed by a slower decay. After the droplet temperature reaches equilibrium (at a time point of 25000τ for E4, 12000τ for E6, and 7500τ for E7), the droplet size decreases very slowly to a plateau value, because the vapor pressure gets closer to saturation during evaporation. To achieve complete evaporation at a constant droplet temperature, the box length was suddenly doubled in the z direction at 40000τ. As can be seen in Figure 5a, the droplet size decreases quickly again from the time point of 40000τ and finally goes to zero, indicating that the droplet completely disappears. The shrinking speed of the three curves follows the same order as the increase in droplet temperature. During this process, the droplet temperatures in E4, E6, and E7 are the same (T* = 0.83), implying that the more hydrophilic substrate induces a faster mass loss by droplet evaporation. Besides, as already mentioned in the Introduction, for contact angles θ < 90° or in a constant-contact-radius evaporation model, the mass linearly depends on time, m ∝ t, as was found in experiments8,14,26 and finite element method (FEM) calculations.12 Another linear relationship where m2/3 ∝ t was observed for contact angles θ > 90° or in systems following a constant-contact-angle evaporation model.11,14,18 It can be noted that, in recent research on femtoliter sessile droplets of glycerol on a hydrophilic substrate,20 the constant-contact-angle model and the m2/3 ∝ t relationship for the macroscopic droplets were verified. However, there is clearly no evidence of m ∝ t in Figure 5a for E4, E6, and E7. To test the m2/3 ∝ t propotionality law, the curve of the variations of N2/3 with respect to time for the last period of evaporation (starting from a time of 40000τ) was plotted and is shown as the inset of Figure 5a. Here, N is the number of atoms inside the droplet and is thus proportional to the total mass m of the droplet. As indicated by the blue dashed line in the inset of Figure 5a, the N2/3 curves for E4, E6, and E7 deviate from the N2/3 ∝ t law as time increases. However, it can be observed that E4, E6, and E7 follow the trend defined by N2/3 ∝ t in the earlier steps of evaporation (until approximately 46000τ). A possible reason for the deviation to N2/3 ∝ t as time increases is the fixed volume of the systems in our simulations. Since the vapor pressure continuously increases during evaporation, we expect that evaporation slows, causing its reduction compared with that of opened systems. Moreover, although the numbers of argon atoms are initially the same in E4, E6, and E7 and the droplet volumes are also very close after equilibration at T* = 0.67, the liquid−vapor interface areas are different due to the different substrate wettabilities. The mass loss rate is also related to the interface area: the larger the liquid−vapor interface, the faster the evaporation expected. To study the effect of substrate wetting on the mass loss, the influence of the liquid−vapor interface area must be removed. To that end, another quantity, namely, the mass flux per unit liquid−vapor interface area, is defined:

fm =

Δm SsΔt

(3)

where Δm is the mass loss, Δt the time interval, and Ss the liquid−vapor interfacial area which is calculated from the fitted spherical radius. The mass flux (Figure 5b) quickly increases for all the systems at the beginning of the simulation when the substrate temperature is suddenly increased from T* = 0.67 to T* = 0.83. The reason for this quick increase in mass flux is that the gas was initially in equilibrium with the liquid droplet at a F


Langmuir

Article

increase of the droplet temperature as well as its evaporation rate. To some extent, it can be concluded that the degree of deformation is related to the evaporation rate (shown in Figure 5d): faster evaporation will cause larger deformations in the droplet shape. After 40000τ, as the box size was doubled in the z direction, faster evaporation was triggered, and the sphericity decreased quickly to less than 0.5 for E6 and E7 and slowly at first and later quickly for E4. It shows larger fluctuations at the end for all the curves. Again, the quicker decrease in E6 and E7 than that in E4 is caused by the higher mass flux, which is more than 2 times larger than that in E4. The small value of the sphericity at the end of all curves indicates that the droplet cannot really be considered as a spherical cap anymore. It could be more like a cylinder or a thin film. The larger fluctuation at the end could be caused by the smaller droplet size. In Figure 7a−c, it can be seen that both the contact angle and the contact radius of the three systems show fluctuations,

system almost converge to the same value of 77 for E4, 79 for E6, and 87 for E7, meaning the net evaporation almost vanishes. It can also be concluded that a stronger interaction strength results in a higher evaporation rate, and a higher evaporation rate goes along with a higher condensation rate. D. Sphericity, Contact Angle, and Contact Radius. To further study the effect of the interaction strength between the liquid and the substrate on evaporation, the contact angle has also been calculated. It should be noted that the droplet’s shape is not completely spherical. It can be slightly or even strongly (at high temperature) deformed during evaporation. To reduce the effect of droplet deformation on the contact angle, the mean value obtained from averaging over several contact angles coming from different positions along the TPCL is considered as the contact angle. In this study, we follow a method used by Leroy and Müller-Plathe:55 Out of the droplet, a vertical slice of width 1.2σ through its center of mass has been cut. A total of 180 slices for each droplet were obtained by rotating the slice 1° at a time around the droplet’s vertical axis. The mass-density distribution in each slice was computed, and then the isodensity contour for each slice was obtained by storing those points where the mass density was equal to half the liquid bulk density. The contour points were subsequently fitted to a circle, from which the contact radius and droplet spherical radius were derived and used further to obtain the contact angle. The final contact angle was obtained by averaging over the values from all slice orientations of the same droplet. As one example shown in Figure 6a, the 180 contact angles were calculated for the droplet in E7 at the time t = 10000τ.

Figure 7. Time evolution of the contact angles and the contact radii for (a) E4, (b) E6, and (c) E7 during the evaporation at a substrate temperature T* = 0.83. The pink and purple vertical dashed lines represent the time points where the droplet temperature reaches thermal equilibrium and when the simulation box size is doubled, respectively. (d) Correlation of contact angle variation and contact radius variation in systems E4, E6, and E7 from their equilibrium values at T = 0.67, i.e., before heating. The pink and purple vertical dashed lines represent the time points where the droplet temperature reaches thermal equilibrium and when the simulation box size is doubled, respectively.

Figure 6. (a) Contact angle along the TPCL of a droplet from E7 at a time point of 10000τ. (b) Time evolution of the sphericity of the droplets in E4, E6, and E7.

The curve shows large fluctuations around a value of 64.6°, and these fluctuations indicate deformation in the droplet shape. To quantify the deformation of the droplet during evaporation, the droplet sphericity S has been calculated via S=

SR smaller SR larger

because the droplet shape during evaporation is not perfectly spherical. Especially in the constant-temperature evaporation phase, larger fluctuations were found, meaning larger deformation in the droplet shape. To clearly see the underlying trends of the curves, some smoothing has been done (dashed lines in Figure 7). The whole evaporation was divided into three regions, separated by two vertical dashed lines. The first region from the beginning of the simulation to the time point marked by the pink vertical line is the temperatureraising evaporation process. In this region, both the contact angle and contact radius decrease for all systems. No pinning of the contact radius or constant angle was found. Hence, there is a mixed-mode behavior (cf. Introduction). After the liquid droplets are completely thermalized to the higher temperature of 0.83, the contact angle decreases from 113.20° to 109.10° for E4, from 87.76° to 81.02° for E6, and from 75.01° to 64.41° for E7 (Table 1). The contact angle decrease in the first region is

(4)

SRsmaller and SRlarger are the separate averages of the 90 smaller spherical radii and of the 90 larger spherical radii, respectively. The total 180 spherical radii come from the fit of the 180 isodensity contours corresponding to the 180 slices of one droplet in the calculation of the contact angle process as mentioned above. If S is equal to 1, the droplet is perfectly spherical; otherwise it is deformed. As shown in Figure 6b, for the first 40000τ, the sphericities of the E4, E6, and E7 systems decrease slightly from about 0.95 to about 0.92 for E4, to 0.91 for E6, and to 0.89 for E7, indicating that the droplet shape in all three systems was almost maintained during the first part of evaporation. The slight decrease in sphericity comes from the G


Langmuir

Article

Table 1. Variations of the Contact Angle and Contact Radius of E4, E6, and E7 in the Temperature-Raising Process system

time period (τ)

contact angle (deg)(T* = 0.67)

contact angle (deg)(T* = 0.83)

E4 E6 E7

0−25000 0−12000 0−7500

113.20 87.76 75.01

109.10 81.02 64.41

Table 2. Variations of the Contact Angle and Contact Radius of E4, E6, and E7 in the Slow Evaporation Process (Corresponding to the Second Region in Figure 7) time period (τ)

contact angle variation (Δθ, deg)

contact radius variation (Δrc, σ)

E4 E6 E7

25000−40000 12000−40000 7500−40000

−0.37 −0.76 −0.80

−0.25 −0.38 −0.64

−4.10 −6.74 −10.60

−2.76 −3.85 −3.89

ideal nature (flat, completely smooth, and chemically homogeneous) of the surface used in our simulations. As was already mentioned in the Introduction, the pinning phenomena during droplet evaporation are usually assumed to be the result of heterogeneities of the substrates.21,37,38 This assumption is supported by our simulation results. These observed variations of the contact angle and contact radius on hydrophilic substrates differ from the experimental observations for macroscopic droplets.11,26,35 They also differ from the behavior of microscopic droplets with subnanoliter (∼100 μm in diameter),17 picoliter (∼50 μm in diameter),14,19 and femtoliter (∼1 μm in diameter)20 volume, where the constant-contactangle model was always observed. In this context, Furuta et al.17 underlined the role played by line tension effects on the evaporation of small droplets. These authors found that the contact angle of droplets with a volume of 0.8 nL almost remains constant during evaporation on heated hydrophilic substrates. In contrast, they observed that the contact angle slightly decreases when the droplets evaporate on heated hydrophobic substrates. Instead, larger droplets with a volume of 2 μL followed a trend opposite that of the previously mentioned subnanoliter droplets. These authors considered that these phenomena were induced by line tension. Our simulation results confirm the trend17 of the contact angle of a subnanoliter droplet to decrease on a hydrophobic substrate. Apart from the possibility of treating two physically very different regimes (an 0.8 nL droplet has a diameter of 112 μm), our results could also indicate that the line tension is not large enough to keep the contact angle constant when the droplet is placed on hydrophilic substrates. At the same time, our results, which follow a mixed evaporation model in E4, E6, and E7, do not support the idea that on ideal surfaces the base diameter recedes while the contact angle remains constant (i.e., constantcontact-angle model) as proposed by Anderson et al.39 E. Local Mass Flux and Heat Flux. It is experimentally well-known that most of the heat and mass transfer during the evaporation happens in the vicinity of the TPCL.41 Note that the droplets studied in the experiments are usually significantly larger than those in our simulation, often of micrometer or millimeter size. Therefore, it takes a longer time for the heat to be transferred from the substrate to the top of the droplet than to the TPCL, as the thickness of the droplet is smaller near the TPCL than at its center. Therefore, most evaporation occurs at the TPCL. The evaporation of macro- or microscopic droplets is usually described by taking into account the existence of prewetting films close to the apparent TPCL whose thickness can be on the order of several tens to hundreds of nanometers. In previous theoretical investigations, the evaporation of droplets was atomistically described by a two-dimensional evaporating meniscus57 or the most evaporative region near the apparent three-phase contact line was numerically modeled by a thin film.58,59 In contrast, no such wetting films around the TPCL are found in our simulations of nanometer-sized droplets. However, one may wonder whether the heat and mass transfer during evaporation still occurs mainly in the vicinity of the

therefore ordered as E7 > E6 > E4. The same is true for the contact radius. Theoretically, a temperature increase can reduce the contact angle, as shown in previous simulations of static droplets:56 Contact angles of droplets on substrates of all wettabilities decrease with temperature. In our simulation, however, the droplet size was dynamically changing. Therefore, not only the temperature increase, but also evaporation should contribute to the contact angle decrease, as discussed below. In the second region, constant-temperature evaporation, both the contact angle and contact radius decrease only slowly, since the vapor pressure is very close to saturation after the droplet temperature reaches 0.83. As shown in Table 2, both the contact angle and contact radius change very little compared with those in the first region. Evaporation almost stops in this region for all three systems.

system

contact angle variation (Δθ, deg) contact radius variation (Δrc, σ)

As for the third region in Figure 7a−c, which is beyond 40000τ when the box size is doubled, the evaporation process was also in a constant-temperature regime. During this period, the contact radius in all systems decreases almost linearly and the contact angle decreases simultaneously. Again, no pinning of the contact radius or constant angle appears. In all systems, the contact angle decreases slowly at first, followed by a period of faster decrease until the droplet disappears. At droplet disappearance, the contact angle is about 20° in all cases. The trend of decrease is shown more clearly in Figure 7d, which plots the contact angle variation (Δθ = θ(t) − θ(0)) as a function of the contact radius variation (Δrc = rc(t) − rc(0)). Each of the curves is divided into three parts, separated by two vertical dashed lines, which delimit the same three regions as in Figure 7a−c. For all systems, there are two different regimes in region 3. Initially, the slopes are small compared to those of regions 1 and 2, indicating less contact angle change per contact area change. This is the behavior closest to a limiting scenario of “droplet shrinking at constant contact angle”. It reverses toward the end of evaporation, where the slopes become very high. Here the system is closer to the other extreme, “evaporation of a droplet of constant radius with shrinking contact angle”. Note, however, that neither extreme behavior is fully realized and that, at all times on all surfaces, a simultaneous reduction of both the contact angle and contact radius is observed (mixed model). The curves end when the droplets become too small to define a meaningful contact radius. The fact that there is no pinning of the droplets in E4, E6, and E7 during the evaporation is expected to arise from the H


Langmuir

Article

Figure 8. (a) Definition of the region of the three-phase contact line. The threshold distances between liquid atoms and the substrate are 1.5σ and 3σ for case A and case B, respectively. Mass flux ratio for the (b) E4 system, (c) E6 system, and (d) E7 system as a function of time. All the curves were smoothed. “Evap” and “Condens” mean evaporation and condensation, respectively. (e) Reduced temperature difference (ΔT = Tlayer − Tdroplet) between the average temperature of the “absorbed” layer with a thickness of 3σ inside the droplet and the average temperature of the rest of the droplet as a function of time. (f) Reduced gas density along the z direction for E4, E6, and E7 at 1500τ.

TPCL. To explore this problem, the local mass flux ratio KN and heat flux ratio KE have been calculated: KN =

RTN RN

RTN = NL /SL

KE =

RTE RE

RTE = E L /SL

RN = NT/ST R E = E T / ST

phase atoms have no potential energy. SL and ST are the same as in eq 5. For condensing atoms, similar definitions are made. In this case, the parameters NL, NT, EL, and ET in eqs 5 and 6 are the number of particles condensed onto the TPCL region per unit of time, the total number of particles condensed onto the whole droplet per unit of time, the kinetic energy carried by atoms condensed onto the TPCL region, and the kinetic energy carried by all atoms condensed onto the droplet per unit time, respectively. The surface areas SL and ST are the same as for evaporation. The values of KN and KE can be considered as the relative transfer abilities in the three-phase contact line region in comparison with that of the whole liquid−vapor interface. Figure 8b−d displays the mass flux ratio before box expansion (t < 40000τ). The KN curves for both choices A and B in all three systems decrease first and then converge to plateau values. The decrease of KN arises from the reduction in the temperature difference between the average temperature of the “absorbed” layer, which is the closest layer to the substrate with a thickness of 3σ, and the average temperature of the rest of the droplet. At the beginning of the simulation, when the temperature of the substrate is suddenly increased, the liquid atoms closest to the substrate first absorb heat and gain higher temperature than the atoms in the rest of the droplet. As shown in Figure 8e, the temperature differences (ΔT = Tlayer − Tdroplet)

(5)

(6)

where NL, NT, SL, and ST are the number of particles evaporated from the TPCL region per unit of time, the total number of particles evaporated from the whole droplet per unit of time, the surface area of the TPCL, and the total surface area, respectively. The surface areas are obtained from fitting spherical caps to the droplets as described above. Here, an atom is counted as being evaporated from the three-phase contact line region if, in one time step, it was within the droplet with a distance to the substrate less than a defined TPCL thickness and if it was in the gas phase in the following time step. Two different values of thickness were used to describe the three-phase contact line (Figure 8a), 1.5σ denoted as case A and 3σ denoted as case B. In eq 6, EL and ET are the kinetic energy carried by atoms evaporated from the TPCL region and the kinetic energy carried by all atoms evaporating from the droplet per unit time, respectively. Here, we assume that gasI


Langmuir

Article

Figure 9. Heat flux ratio for the (a) E4 system, (b) E6 system, and (c) E7 system during evaporation. For the definition of cases A and B, see Figure 7. “Evap” and “Condens” refer to evaporation and condensation, respectively.

resistance, therefore higher evaporation can be expected. For the E6 and E7 systems, KN following definition A is quite a bit larger than that according to definition B, suggesting strongly that the evaporation rate is nonuniform along the liquid−vapor interface. Moreover, the KN values of condensations for both cases A and B of the three systems have also been calculated and are shown in Figure 8 (labeled as “Condens”). From these curves, the same conclusion can be drawn; i.e., condensation also happens preferentially in the TPCL area for systems E6 and E7 and is almost uniform over the liquid−vapor interface for system E4. This could result from the higher gas density close to the TPCL. In Figure 8f, the gas density along z is plotted. At both ends of the curves (due to the periodic boundary conditions, the gas contacts both surfaces of the substrate), higher values are found for E6 and E7, especially for the end with a droplet (z < 10). The gas density close to the surface is higher for E6 and E7 because of the stronger attraction between the substrate and gas, as confirmed by the higher gas density values at the end where z > 140. The fact that the gas density is higher at the surface with the droplet (z < 10) than at that without a droplet (z > 140), which is especially noticeable for system E4, indicates the occurrence of droplet evaporation. As a result, a higher gas density around the TPCL can be expected, and the higher gas density leads to more condensation events in this area. The local heat flux ratios for both evaporation and condensation are shown in Figure 9. Their time evolutions are similar to the behavior of the mass fluxes. All curves decrease first and then converge to plateaus, which for evaporation are 1.05 and 1.0 for E4, 2.0 and 1.5 for E6, and 2.7 and 1.7 for E7 for definitions A and B, respectively. Hence, also the evaporative heat flux is enhanced in the three-phase contact line area for hydrophilic substrates (E6 and E7); however, it is almost homogeneous over the liquid−vapor interface for the hydrophobic case (E4). Again, the heat fluxes in condensation show the same trend as in evaporation.

for E4, E6, and E7 are over 0.012 at the beginning of evaporation. The evaporation occurring in the region of the three-phase contact line is therefore faster than that elsewhere. The values of KN being universally above 1 in our simulation confirms this assumption. Due to the continuous heat transfer from the substrate to the droplet, the temperature difference between the absorbed layer of liquid and the bulk of the droplet decreases, as shown in Figure 8e, gradually to zero. Thus, the relative importance of evaporation at the three-phase contact line decreases gradually. When the temperature difference disappears, KN converges to a characteristic value. For case A, KN for evaporation converges to values of 1.15, 2.2, and 2.8 for E4, E6, and E7, respectively. For case B, the KN values are 1.05, 1.6, and 1.8 for E4, E6, and E7, respectively. Hence, we know that, for E4, KN is closest to 1 for both definitions of the contact line region, which means that, for the hydrophobic substrate, evaporation near the TPCL is almost the same as that of other surface areas. However, the almost uniform mass flux along the droplet surface featured on the hydrophobic substrate (E4) in our simulations is different from the prediction of the theoretical analysis of Nguyen et al.29 By means of an analytical approach, these authors reported that the diffusive flux along the surface on a hydrophobic substrate decreases from the droplet apex to the edge, where it vanishes. However, the analytical approach does not account for the local evolution of the phenomenon which is accessible to our MD simulation. Therefore, experimental and theoretical work is needed to draw clear conclusions. For E6 and E7, however, evaporation preferentially occurs in the vicinity of the threephase contact line, which has also been observed for larger droplets (macroscopic droplets with contact angle θ < 90°) from experiments41 and finite element calculations12 as well as theoretical analysis.25,29 This could be explained by the thickness (or height) difference in different areas of the droplet; for the hydrophilic substrate, the thickness in the TPCL region is smaller than in the other area, since the contact angle is less than 90°. Smaller thickness results in less thermal J


Langmuir

Article

contact angle decreases much more rapidly than the droplet radius. Both mass flux and heat flux were locally resolved. They strongly indicate enhanced mass and heat flux near the threephase contact line. On hydrophilic substrates the contact line region contributes disproportionately to both fluxes, which is in good agreement with experimental results for macroscopic droplets. For hydrophobic substrates, however, which are rarely investigated in experiments, the evaporation rate is almost uniform over the droplet surface. Condensation, which is also present during the evaporation process, follows the same trends as evaporation; i.e., gas molecules preferentially condense into the droplet in the three-phase contact line region for droplets on a hydrophilic substrate, but show no preference on the hydrophobic substrates.

Therefore, the substrate properties have a very important effect on the local distribution of mass flux and heat flux.

IV. CONCLUSIONS We present an MD study on the evaporation of nanosized droplets from heated substrates. In contrast to the usual experimental conditions, where the droplets are >1 μm in diameter and where the droplets are placed on possibly heterogeneous substrates, the droplets in our simulations are a few nanometers in diameter and the substrates are ideally smooth and chemically homogeneous. To study the effect of substrate wettability on the evaporation, three different surface interactions were used by changing the Lennard-Jones parameter ε* from 0.4 to 0.7 to tune the substrate between hydrophobic and hydrophilic. The evaporation behavior has been investigated during both heating of the droplet and at constant droplet temperature. The temperature profiles of the droplet and gas show that there is a temperature gap between the gas phase and the droplet in the temperature-raising evaporation process; the gas temperature is always larger than the droplet temperature. The temperature increase is much faster in the gas phase than that in the droplet for both hydrophobic and hydrophilic substrates. The mass transfer has also been studied in our simulations. The mass flux from the droplet has a sharp increase at the beginning of the temperature-raising evaporation process due to the quick increase in the droplet temperature and is then followed by a decrease caused by the gas pressure approaching saturation. In the subsequent quick evaporation induced by the sudden pressure drop in the gas phase at constant temperature, the wettability of the substrate has a strong effect on the mass flux. The more hydrophilic the substrate, the faster the evaporation, since the substrate-to-droplet heat transfer is more efficient when the interaction between the solid and droplet becomes stronger. During the initial stage of the evaporations induced by the sudden pressure drop, the time dependence of the residual mass of nanosized droplets was found to follow the m2/3 ∝ t scaling law, which is experimentally or theoretically observed for larger droplets (>1 μm in diameter). The droplet shape is maintained during most of the evaporation time. However, this shape is strongly deformed and the droplets are no longer spherical at the late stage of the evaporation. Therefore, the droplets should not be considered as spherical close to the end of the evaporation. The contact angle and contact radius variations over time have also been computed. Different from the usual two-stage (pinning and depinning) evaporation observed in experiments, only one evaporation mode and no pinning of the contact line was found in our simulations. Such a difference between the experimental and simulation results is the consequence of the ideal nature of the substrate used in our simulations. The contact angle and contact radius always decrease simultaneously: for hydrophilic and hydrophobic substrates alike and for both temperature-raising and constant-temperature evaporation processes. The evaporation always follows the so-called mixed model with, however, different phases. The dynamics of the contact angle and contact line are different from the experimental observations on macroscopic droplets (>100 μm in diameter) and on microscopic droplets (1−100 μm in diameter), for which the constant-contact-angle model or/and constant-contact-radius model are usually reported. In our simulations, larger droplets show a faster decrease of the radius than the contact angle, whereas below a certain droplet size the

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS Funding from the Excellence Cluster 259 Center of Smart Interfaces at Technische Universität Darmstadt and the Collaborative Research Center Transregio 75 of the Deutsche Forschungsgemeinschaft is gratefully acknowledged. We also thank Prof. Peter Stephan and Dr. Tatiana GambaryanRoisman for very helpful discussions.

■

REFERENCES

(1) Sirringhaus, H.; Kawase, T.; Friend, R.; Shimoda, T.; Inbasekaran, M.; Wu, W.; Woo, E. High-Resolution Inkjet Printing of All-Polymer Transistor Circuits. Science 2000, 290 (5499), 2123−2126. (2) Lim, J. A.; Lee, W. H.; Lee, H. S.; Lee, J. H.; Park, Y. D.; Cho, K. Self-Organization of Inkjet Printed Triisopropylsilylethynyl Pentacene via EvaporationInduced Flows in a Drying Droplet. Adv. Funct. Mater. 2008, 18 (2), 229−234. (3) Bigioni, T. P.; Lin, X. M.; Nguyen, T. T.; Corwin, E. I.; Witten, T. A.; Jaeger, H. M. Kinetically Driven Self-Assembly of Highly Ordered Nanoparticle Monolayers. Nat. Mater. 2006, 5 (4), 265−270. (4) Jing, J.; Reed, J.; Huang, J.; Hu, X.; Clarke, V.; Edington, J.; Housman, D.; Anantharaman, T. S.; Huff, E. J.; Mishra, B. Automated High Resolution Optical Mapping Using Arrayed, Fluid-Fixed DNA Molecules. Proc. Natl. Acad. Sci. U.S.A. 1998, 95 (14), 8046−8051. (5) Hsieh, C. C.; Li, L.; Larson, R. G. Modeling Hydrodynamic Interaction in Brownian Dynamics: Simulations of Extensional Flows of Dilute Solutions of DNA and Polystyrene. J. Non-Newtonian Fluid Mech. 2003, 113 (2), 147−191. (6) Norris, D. J.; Arlinghaus, E. G.; Meng, L.; Heiny, R.; Scriven, L. Opaline Photonic Crystals: How Does Self-Assembly Work? Adv. Mater. 2004, 16 (16), 1393−1399. (7) Picknett, R.; Bexon, R. The Evaporation of Sessile or Pendant Drops in Still Air. J. Colloid Interface Sci. 1977, 61 (2), 336−350. (8) Birdi, K.; Vu, D.; Winter, A. A Study of the Evaporation Rates of Small Water Drops Placed on a Solid Surface. J. Phys. Chem. 1989, 93 (9), 3702−3703. (9) Bourges-Monnier, C.; Shanahan, M. Influence of Evaporation on Contact Angle. Langmuir 1995, 11 (7), 2820−2829. (10) Fang, G.; Ward, C. A. Temperature Measured Close to the Interface of an Evaporating Liquid. Phys. Rev. E 1999, 59 (1), 417− 428. (11) Erbil, H. Y.; McHale, G.; Newton, M. Drop Evaporation on Solid Surfaces: Constant Contact Angle Mode. Langmuir 2002, 18 (7), 2636−2641. K


Langmuir

Article

(36) Quéré, D. Surface Wetting: Model Droplets. Nat. Mater. 2004, 3 (2), 79−80. (37) Orejon, D.; Sefiane, K.; Shanahan, M. E. Stick−Slip of Evaporating Droplets: Substrate Hydrophobicity and Nanoparticle Concentration. Langmuir 2011, 27 (21), 12834−12843. (38) Grandas, L.; Reynard, C.; Santini, R.; Tadrist, L. Experimental Study of the Evaporation of a Sessile Drop on a Heat Wall. Wetting Influence. Int. J. Therm. Sci. 2005, 44, 137. (39) Anderson, D.; Davis, S. The Spreading of Volatile Liquid Droplets on Heated Surfaces. Phys. Fluids 1995, 7, 248. (40) Ward, C.; Stanga, D. Interfacial Conditions during Evaporation or Condensation of Water. Phys. Rev. E 2001, 64 (5), 051509. (41) Stephan, P.; Busse, C. Analysis of the Heat Transfer Coefficient of Grooved Heat Pipe Evaporator Walls. Int. J. Heat Mass Transfer 1992, 35 (2), 383−391. (42) Ibrahem, K.; Abd Rabbo, M.; Gambaryan-Roisman, T.; Stephan, P. Experimental Investigation of Evaporative Heat Transfer Characteristics at the 3-Phase Contact Line. Exp. Therm. Fluid Sci. 2010, 34 (8), 1036−1041. (43) Sobac, B.; Brutin, D. Triple-Line behavior and Wettability Controlled by Nanocoated Substrates: Influence on Sessile Drop Evaporation. Langmuir 2011, 27 (24), 14999−15007. (44) Boruvka, L.; Neumann, A. Generalization of the Classical Theory of Capillarity. J. Chem. Phys. 1977, 66, 5464. (45) Long, L. N.; Micci, M. M.; Wong, B. C. Molecular Dynamics Simulations of Droplet Evaporation. Comput. Phys. Commun. 1996, 96 (2), 167−172. (46) Consolini, L.; Aggarwal, S. K.; Murad, S. A Molecular Dynamics Simulation of Droplet Evaporation. Int. J. Heat Mass Transfer 2003, 46 (17), 3179−3188. (47) Matsumoto, M. Molecular Dynamics Simulation of Interphase Transport at Liquid Surfaces. Fluid Phase equilib. 1996, 125 (1), 195− 203. (48) Yang, T.; Pan, C. Molecular Dynamics Simulation of a Thin Water Layer Evaporation and Evaporation Coefficient. Int. J. Heat Mass Transfer 2005, 48 (17), 3516−3526. (49) Berendsen, H. J.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. Molecular Dynamics with Coupling to an External Bath. J. Chem. Phys. 1984, 81, 3684. (50) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117 (1), 1−19. (51) Daubert, T. E.; Danner, R. P.; Sibul, H.; Stebbins, C. Physical and Thermodynamic Properties of Pure Chemicals: Data Compilation; Hemisphere Publishing Corp.: New York, 1989. (52) Wu, Y.; Pan, C. Molecular Dynamics Simulation of Thin Film Evaporation of Lennard-Jones Liquid. Nanoscale Microscale Thermophys. Eng. 2006, 10 (2), 157−170. (53) Shenogina, N.; Godawat, R.; Keblinski, P.; Garde, S. How Wetting and Adhesion Affect Thermal Conductance of a Range of Hydrophobic to Hydrophilic Aqueous Interfaces. Phys. Rev. Lett. 2009, 102 (15), 156101. (54) Babin, V.; Holyst, R. Evaporation of a Sub-Micrometer Droplet. J. Phys. Chem. B 2005, 109 (22), 11367−11372. (55) Leroy, F.; Müller-Plathe, F. Solid-Liquid Surface Free Energy of Lennard-Jones Liquid on Smooth and Rough Surfaces Computed by Molecular Dynamics Using the Phantom-Wall Method. J. Chem. Phys. 2010, 133, 044110. (56) Shi, B.; Dhir, V. K. Molecular Dynamics Simulation of the Contact Angle of Liquids on Solid Surfaces. J. Chem. Phys. 2009, 130, 034705. (57) Freund, J. B. The Atomic Detail of an Evaporating Meniscus. Phys. Fluids 2005, 17, 022104. (58) Ma, H.; Cheng, P.; Borgmeyer, B.; Wang, Y. Fluid Flow and Heat Transfer in the Evaporating Thin Film Region. Microfluid. Nanofluid. 2008, 4 (3), 237−243. (59) Kunkelmann, C.; Ibrahem, K.; Schweizer, N.; Herbert, S.; Stephan, P.; Gambaryan-Roisman, T. The Effect of Three-Phase Contact Line Speed on Local Evaporative Heat Transfer: Experimental

(12) Hu, H.; Larson, R. G. Evaporation of a Sessile Droplet on a Substrate. J. Phys. Chem. B 2002, 106 (6), 1334−1344. (13) Soolaman, D. M.; Yu, H. Z. Water Microdroplets on Molecularly Tailored Surfaces: Correlation between Wetting Hysteresis and Evaporation Mode Switching. J. Phys. Chem. B 2005, 109 (38), 17967−17973. (14) Golovko, D. S.; Butt, H.-J. r.; Bonaccurso, E. Transition in the Evaporation Kinetics of Water Microdrops on Hydrophilic Surfaces. Langmuir 2008, 25 (1), 75−78. (15) Ristenpart, W.; Kim, P.; Domingues, C.; Wan, J.; Stone, H. Influence of Substrate Conductivity on Circulation Reversal in Evaporating Drops. Phys. Rev. Lett. 2007, 99 (23), 234502. (16) Choi, C.-H.; Kim, C.-J. C. Droplet Evaporation of Pure Water and Protein Solution on Nanostructured Superhydrophobic Surfaces of Varying Heights. Langmuir 2009, 25 (13), 7561−7567. (17) Furuta, T.; Sakai, M.; Isobe, T.; Nakajima, A. Evaporation Behavior of Microliter- and Sub-Nanoliter-Scale Water Droplets on Two Different Fluoroalkylsilane Coatings. Langmuir 2009, 25 (20), 11998−12001. (18) Xu, W.; Leeladhar, R.; Kang, Y. T.; Choi, C.-H. Evaporation Kinetics of Sessile Water Droplets on Micropillared Superhydrophobic Surfaces. Langmuir 2013, 29 (20), 6032−6041. (19) Lim, T.; Han, S.; Chung, J.; Chung, J. T.; Ko, S.; Grigoropoulos, C. P. Experimental Study on Spreading and Evaporation of Inkjet Printed Pico-Liter Droplet on a Heated Substrate. Int. J. Heat Mass Transfer 2009, 52 (1), 431−441. (20) Arcamone, J.; Dujardin, E.; Rius, G.; Pérez-Murano, F.; Ondarçuhu, T. Evaporation of Femtoliter Sessile Droplets Monitored with Nanomechanical Mass Sensors. J. Phys. Chem. B 2007, 111 (45), 13020−13027. (21) Shin, D. H.; Lee, S. H.; Jung, J.-Y.; Yoo, J. Y. Evaporating Characteristics of Sessile Droplet on Hydrophobic and Hydrophilic Surfaces. Microelectron. Eng. 2009, 86 (4), 1350−1353. (22) Cazabat, A.-M.; Guéna, G. Evaporation of Macroscopic Sessile Droplets. Soft Matter 2010, 6 (12), 2591−2612. (23) Maxwell, J. C. The Scientific Papers of James Clerk Maxwell; University Press: Cambridge, U.K., 1890; Vol. 2. (24) Schrage, R. W. A Theoretical Study of Interphase Mass Transfer; Columbia University Press: New York, 1953. (25) Deegan, R. D.; Bakajin, O.; Dupont, T. F.; Huber, G.; Nagel, S. R.; Witten, T. A. Capillary Flow as the Cause of Ring Stains from Dried Liquid Drops. Nature 1997, 389 (6653), 827−829. (26) Rowan, S.; Newton, M.; McHale, G. Evaporation of Microdroplets and the Wetting of Solid Surfaces. J. Phys. Chem. 1995, 99 (35), 13268−13271. (27) Cachile, M.; Benichou, O.; Cazabat, A. Evaporating Droplets of Completely Wetting Liquids. Langmuir 2002, 18 (21), 7985−7990. (28) Bond, M.; Struchtrup, H. Mean Evaporation and Condensation Coefficients Based on Energy Dependent Condensation Probability. Phys. Rev. E 2004, 70 (6), 061605. (29) Nguyen, T. A.; Nguyen, A. V.; Hampton, M. A.; Xu, Z. P.; Huang, L.; Rudolph, V. Theoretical and Experimental Analysis of Droplet Evaporation on Solid Surfaces. Chem. Eng. Sci. 2012, 69 (1), 522−529. (30) Girard, F.; Antoni, M.; Faure, S.; Steinchen, A. Evaporation and Marangoni Driven Convection in Small Heated Water Droplets. Langmuir 2006, 22 (26), 11085−11091. (31) Hu, H.; Larson, R. G. Analysis of the Microfluid Flow in an Evaporating Sessile Droplet. Langmuir 2005, 21 (9), 3963−3971. (32) Hołyst, R.; Litniewski, M. Heat Transfer at the Nanoscale: Evaporation of Nanodroplets. Phys. Rev. Lett. 2008, 100 (5), 55701. (33) 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 (1), 756. (34) Langmuir, I. The Evaporation of Small Spheres. Phys. Rev. 1918, 12 (5), 368. (35) Birdi, K.; Vu, D. Wettability and the Evaporation Rates of Fluids from Solid Surfaces. J. Adhes. Sci. Technol. 1993, 7 (6), 485−493. L


Langmuir

Article

and Numerical Investigations. Int. J. Heat Mass Transfer 2012, 55 (7), 1896−1904.

M


Evaporation of Nanodroplets on Heated ... - American Chemical Society

Recommend Documents