Evaporation of liquid droplets of nano- and micro-meter size as a function of molecular mass and intermolecular interactions: experiments and molecular dynamics simulations

Abstract

Transport of heat to the surface of a liquid is a limiting step in the evaporation of liquids into an inert gas. Molecular dynamics (MD) simulations of a two component Lennard-Jones (LJ) fluid revealed two modes of energy transport from a vapour to an interface of an evaporating droplet of liquid. Heat is transported according to the equation of temperature diffusion, far from the droplet of radius R. The heat flux, in this region, is proportional to temperature gradient and heat conductivity in the vapour. However at some distance from the interface, Aλ, (where λ is the mean free path in the gas), the temperature has a discontinuity and heat is transported ballistically i.e. by direct individual collisions of gas molecules with the interface. This ballistic transport reduces the heat flux (and consequently the mass flux) by the factor R/(R + Aλ) in comparison to the flux obtained from temperature diffusion. Thus it slows down the evaporation of droplets of sizes R ∼ Aλ and smaller (practically for sizes from 10³ nm down to 1 nm). We analyzed parameter A as a function of interactions between molecules and their masses. The rescaled parameter, A(k_BT_b/ε₁₁)^1/2, is a linear function of the ratio of the molecular mass of the liquid molecules to the molecular mass of the gas molecules, m₁/m₂ (for a series of chemically similar compounds). Here ε₁₁ is the interaction parameter between molecules in the liquid (proportional to the enthalpy of evaporation) and T_b is the temperature of the gas in the bulk. We tested the predictions of MD simulations in experiments performed on droplets of ethylene glycol, diethylene glycol, triethylene glycol and tetraethylene glycol. They were suspended in an electrodynamic trap and evaporated into dry nitrogen gas. A changes from ∼1 (for ethylene glycol) to approximately 10 (for tetraethylene glycol) and has the same dependence on molecular parameters as obtained for the LJ fluid in MD simulations. The value of x = A(k_BT_b/ε₁₁)^1/2 is of the order of 1 (for water x = 1.8, glycerol x = 1, ethylene glycol x = 0.4, tetraethylene glycol x = 2.1 evaporating into dry nitrogen at room temperature and for Lennard-Jones fluids x = 2 for m₁/m₂ = 1 and low temperature).

---

Introduction

An evaporation process occurs at the interface between liquid and gas in the region of nanometer size. This region is very hard to access in direct thermodynamic experiments and therefore profiles of thermodynamic parameters governing evaporation in the vicinity of the interface are not easily available from experiments.^1,2 Computer simulations, on the other hand, give full information about evaporation in nano-sized systems, but unfortunately cannot be easily scaled to macroscopic systems studied in experiments.¹ In this paper we combine molecular dynamics (MD) simulations in nano-meter size liquid–gas systems and experiments performed on micro-meter size droplets evaporating into inert gas to study evaporation. We obtain a simple formula for evolution of liquid droplet radius R(t) at time t, during evaporation as a function of thermodynamic and material parameters characterizing the system.

Liquids evaporate into inert gas at the expense of internal energy of the gas. The heat transfer from the gas to the liquid is the limiting step in the process of evaporation. During evaporation the bulk temperature of the gas, T_b, is higher than the temperature of the evaporating liquid, T_liq. Far from the interface heat is transported towards the interface via temperature diffusion. Close to the interface temperature has a discontinuity (at the scale of mean free path in the gas phase). This effect was first measured by Ward and coworkers^3–6 in water close to the triple point, where the mean free path was of the order of micro-meters. This effect was later observed in molecular dynamics simulations⁷ of Lennard-Jones molecules in small nanoscopic systems. Due to the temperature discontinuity heat is transported in a ballistic manner in this region. The ballistic transport modifies the standard formula for evaporation where only temperature diffusion is taken into account. In our previous paper⁸ we have established in experiments and MD simulations the following formula for liquid droplet radius R(t) as a function of time t during evaporation into inert gas:

where λ is the mean free path in the gas phase, Δh = h(T_b) − h(T_liq) is the total change in enthalpy upon evaporation, κ_v is the heat conductivity in the gas phase at T_liq and ρ_liq is the liquid density at T_liq. All parameters in this equation apart from A can be simulated for model systems (see the ESI†) or taken from the NIST database for any given gas and liquid. Eqn (1) follows from the conservation of energy: the left hand side of eqn (1) gives the enthalpy consumed by the liquid droplet per unit time and surface area and the right hand side gives the heat flux from the vapour to the liquid. In this equation the ballistic mechanism of energy transfer reduces the heat flux from the vapour to the liquid near the interface by a factor R/(R + Aλ) in comparison to heat flux obtained purely from the temperature diffusion mechanism. Parameter A is a fit parameter. In the present paper we establish the dependence of A on thermodynamic and molecular parameters of the system (temperature of the gas phase, interaction between liquid molecules, molecular mass of liquid molecules and their molecular size). These parameters are important for energy transfer in molecular collisions between gas molecules and the liquid interface. In particular energy transfer from vapour to liquid is strongly reduced in molecular collisions when the molecular mass of gas molecules is much smaller than the molecular mass of liquid molecules. A similar role is being played by molecular interactions in the liquid. The temperature of gas (proportional to the kinetic energy of gas molecules) acts in the opposite manner.

We performed computer simulations of a two-component Lennard-Jones system. Both components differed in mass and size. We performed MD simulations at a few temperatures in the gas phase. Additionally we re-analysed our previous experimental data⁸ for very small liquid droplets (i.e. comparable to Aλ) to obtain reliable values of A. For R ≫ Aλeqn (1) is not sensitive to the ballistic transport of energy, and therefore the experimental data should be fit separately for large and small R values. We get ΔT = T_b − T_liq from fits to eqn (1) at large R. Next we fix ΔT and fit A parameter at small R.

In the computer simulations component 2 played the role of inert gas and component 1 formed a liquid droplet. During evaporation component 1 also entered into the gas phase, but component 2 did not dissolve in the liquid phase. The temperature of the bulk gas was set between the critical temperatures of the components. The gas component was well above its critical point and the liquid phase was below its critical point. The liquid droplet maintained a constant temperature during evaporation, with T_liq being smaller than T_b. After a short initial transient state, the temperature profile T(r,t) was quasi-stationary i.e. it quickly adapted to the changes of the droplet radius R(t) at time t. Thus we had T(r,t) = T(r/R(t)) during the evaporation. A schematic drawing of our system is shown in Fig. 1.

Fig. 1 A schematic drawing of the simulation system. In the left panel we show the simulation box of size L with a liquid droplet of size R consisting of one Lennard-Jones type of molecule 1. The droplet is surrounded by inert gas (component 2) and vapour molecules of component 1 (same as in the liquid phase). The temperature of vapour molecules leaving a liquid droplet and gas molecules as a function of distance r from the centre of the droplet are different (T₁ and T₂). Far from the droplet the temperature of both components is the same and equal to T_b, the temperature of the bulk gas phase. The typical temperature profile of component 1 as a function of the distance from the centre of the droplet, r, is shown in the right panel.

The system studied using computer simulations mimicked the typical two-component systems studied by us experimentally i.e. evaporation of ethylene glycol, diethylene glycol, triethylene glycol and tetra-ethylene glycol into a dry nitrogen atmosphere. We combined our MD simulations and experimental studies of these four compounds to find A as a function of molecular parameters characterizing the molecules in the system.

Computer simulations

Molecular Dynamics (MD) simulations were performed using the Verlet “leapfrog” method^9,10 to solve the Newton equations of motion that describe the time evolution of a two-component system composed of low volatile (the first component indexed by 1) and high volatile (the second component indexed by 2) molecules. The molecules of the same kind (ii = 11 or 22) interact via the short range Lennard-Jones (LJ) based potential given in the following form:and the repulsive interaction between components 1 and 2 was given byand 0 for d/σ₁₂ > 2^1/6. Here i = 1, 2 is the component index and d is the inter-molecular distance. The A₀, B₀, C₀, and C₀′ constants were adjusted to make the potential differentiable in the whole space. The particles of component 2 (forming an inert gas phase) had always the same mass m₂ = 1, the size parameter σ₂₂ = 1 and the energy parameter ε₂₂ = 1/4. The temperature was always given in 4ε₂₂/k_B units where k_B is the Boltzmann constant. All numerical values presented further are expressed in these units. The particles of component 1 (forming a liquid phase) had always the same energy parameter ε₁₁ = 2. The masses and the size parameters depended on the simulation run. The mass m₁ ranged from 1 to 8 for σ₁₁ = 1 and from 3 to 20 for σ₁₁ = 2. The mixing rule was ε₁₂ = 1 and σ₁₂ = (σ₁₁ + σ₂₂)/2. The Verlet scheme¹⁰ was applied with the time step δt = 0.005σ₂₂(m₂/ε₂₂)^1/2.

Simulations of droplet evaporation were performed in the three dimensional sphere of radius R_b. The list of all simulation runs is given in Table S1 of the ESI.† A schematic picture of the system is shown in Fig. 1. The distance from the centre of the sphere is denoted by r. Initially, the sphere was filled with the gas of component 2 at density ρ₂ = ρ⁰₂ and the liquid droplet that is composed of molecules of component 1 was placed with the centre of the mass at r = 0. Two systems were considered: the first of R_b = 340.6, ρ⁰₂ = 0.015, σ₁₁ = 1.0 with the total number of particles N = 2.91 × 10⁶ and the second of R_b = 557.3, ρ⁰₂ = 0.0065, σ₁₁ = 2.0, N = 4.91 × 10⁶. The number of molecules of component 1 was always N₁ = 2.15 × 10⁵. The particles of component 2 were reflected elastically and specularly from the boundary of the sphere at r = R_b. The way the particles are reflected has no influence on the evaporation process. (R_b − R(t))/λ is many times larger than 1. As a consequence, the particle distribution function close to the liquid surface does not depend on the way the conditions at R_b are imposed. The boundary temperature T_b was imposed by appropriate scaling of velocities of molecules of component 2 for r > R_b − kσ₂₂ where k = 3 to 5 depending on ρ⁰₂. The droplet volume was negligible if compared to the total volume of the system. Far from the droplet surface the density, ρ₂, was constant and equal to the initial value, ρ⁰₂. To make simulation conditions closer to real ones the molecules of component 1 could go out as well as come back into the simulation sphere. For r > R_b the molecules of component 1 moved independently in stochastic gas (see Fig. 1) characterized by the temperature, T_b, and the self-diffusion constant,⁹D₂. The stochastic medium was simulated introducing the Langevin force into the Verlet scheme (here with no 11 interaction). The random force was taken from the Euler–Maruyama approximation for the stochastic term.¹¹ The approximation was used to generate random force to move molecules of component 1 at distances r > R_b. Also for r > R_b gas 2 was replaced by the Brownian gas. For r < R_b no random force was used and the time evolution was always performed using a classical MD method. Molecules of component 1 could cross (in both sides) the R_b border and spread all over the space as in real systems for which movement of particles 1 is not restricted (i.e. like for R_b → ∞).

The molecules of component 1 formed the droplet of radius R of the number liquid density, ρ_liq, and the liquid temperature, T_liq. The droplet radius, R, was defined via the total number density i.e. ρ(r = R) = ρ_liq/2. Based on additional simulations we estimated the critical temperature for the interactions given by eqn (2) to T_C ≈ 0.65ε_ii (ESI†). The pressure value as a function of r was calculated during simulations using the standard formula (ref. 9). In the simulations we considered only a quasi-stationary stage during which the droplet evaporated at constant temperature and density. The pressure maintained its value in the initial non-stationary regime (ref. 7). During the whole quasi-stationary stage (ref. 7, 8 and Supplementary Information to ref. 8) the pressure value was constant.

The liquid was under thermodynamic conditions well below its critical point while the gas of component 2 (ε₂₂ = 1/4) formed an inert gas well above its corresponding T_C. Such simulations correspond to the experimental case of evaporation of liquids (e.g. water) into the inert gas e.g. air. The values of κ_v in eqn (1) were determined from additional MD simulations (ESI†) applying the method of Muller-Plathe¹² and Δh was calculated from the simulation data. As in our previous paper,⁸λ was estimated from

The values of the self-diffusion coefficient of molecules of component 2, D₂, were determined from the Einstein formula⁹ by performing additional MD simulations for the gas of particles of component 2 at T₂ = T_b and ρ₂ = ρ⁰₂ (ESI†). The value of A was determined from the fit to solution of eqn (1) given by eqn (5).

Kinetic parameter A from computer simulations

The evaporation process observed in simulations proceeds via an initial strongly non-equilibrium process characterized by the fast decrease of R, rapid decrease of T_liq and increase of ρ_liq.

Thus initially the enthalpy needed for evaporation comes from the decrease of temperature of the droplet. This stage is very short. In the second quasi-stationary stage the temperature in the gas phase develops a profile adjusting to the change in the radius of the droplet R(t), T(r,t) = T(r/R(t)). At this stage of evaporation the temperature of the droplet, T_liq, was constant. During this process we observed (Fig. 2) that molecules of component 1 leaving the droplet–gas interface had markedly lower temperature than the temperature of the gas.

Fig. 2 Temperature profiles for the quasi-stationary stage of evaporation for two components (T₁ (yellow, solid line) and T₂ (blue, dots)) and density profile ρ₁ for the first component (red colour, solid, dots). The upper panel shows the profiles for a system of liquid molecules having 8 times larger mass than gas molecules and the lower panel is for the same molecular masses of components 1 and 2. The molecules of component 1 which evaporate from the surface of the liquid droplet have smaller temperature than the gas made of component 2. The temperature profile of the gas (component 2) has a discontinuity at the interface. The temperature of component 1 becomes equal to the temperature of component 2 at distance Aλ from the interface. Both systems have the same values of σ₁₁ = 1 and T_b = 1.15.

In order to determine A we solved eqn (1):

where B = κ_V(T_b − T_liq)/(Δhρ_liq). The dependence of A on m₁ is a consequence of a decrease in efficiency of energy transfer via collisions if m₁/m₂ ≠ 1. The effect is seen in Fig. 2 in which the temperature discontinuity at the gas–liquid interface is much higher (i.e. the energy transfer is slower) for m₁/m₂ = 8 than for m₁/m₂ = 1. The temperature jump for the gas temperature profile at the interface is given by:^7,8 This discontinuity increases with A. Thus by comparing the upper and lower panels in Fig. 2 we find that A is a function of molecular mass. In Fig. 3 we plot A as a function of molecular mass for different molecular sizes. We removed any temperature dependence of A in Fig. 3, since we observed that A is inversely proportional to the square root of gas bulk temperature as shown in Fig. 4. The final approximate formula for A as a function of molecular parameters and temperature is as follows (Fig. 3 and 4):where α₁ and α₀ (two dimensionless parameters of the order of 1) weakly depend on their molecular size (range of interactions and structure of the liquid). The main mode of energy transfer between gas and liquid is via direct collisions of gas molecules with the surface of a liquid and exchange of kinetic energy. The interaction parameter between gas and liquid ε₁₂ does not affect this transfer appreciably in rarified gas. We have verified this physical observation in computer simulations. In our simulations gas 2 is always a rarefied gas well above its critical temperature. We set (ref. 8) ε₁₂ = ε₂₂ and ε₂₂ = γε₁₁ where the interaction parameter between liquid molecules ε₁₁ = 1. We changed γ from 0.05 to 0.0005 without any observable changes in A. We have also verified that A does not depend on the strength of interactions between molecules of component 2 which forms a gas phase. Finally the dependence of A on molecular interaction between molecules of component 1 (forming a liquid) and their masses were tested in experiments performed on a series of similar chemical compounds: ethylene glycol, diethylene glycol, triethylene glycol and tetra-ethylene glycol evaporating into dry nitrogen gas.

Fig. 3 Parameter A as a function of the molecular mass of liquid divided by the mass of gas molecules, m₁/m₂. The results for two sizes of molecules are shown. The inclination decreases with the size of particle 1. The triangles-up (blue) – σ₁₁ = 1.0, T_b = 1.05; the squares (red) – σ₁₁ = 1.0, T_b = 1.15; the circles (orange) σ₁₁ = 2, T_b = 1.05; the triangles down (black) σ₁₁ = 2, T_b = 1.25. The dashed lines are obtained from fits to eqn (5) for a given σ₁₁. We find (eqn (7)) α₁ = 0.65 and α₀ = 1.42 for σ₁₁/σ₂₂ = 1.0 and α₁ = 0.42 and α₀ = 1.40 for σ₁₁/σ₂₂ = 2.0. Thus α₁ weakly decreases with an increase in the size of molecules forming a liquid, while α₀ is a constant and probably depends on molecular shape (as we show in the Experimental section).

Fig. 4 Parameter A versus temperature of the gas phase, T_b, for different masses and sizes of molecules: the squares (blue) – σ₁₁ = 1, m₁ = 1; the triangles up (red) – σ₁₁ = 2, m₁ = 3; the triangles down (black) – σ₁₁ = 2, m₁ = 8. The dotted lines are given by the fit to eqn (7).

Experimental determination of parameter A

We have used previously obtained data⁸ on diethylene glycol and triethylene glycol and performed additional experiments on ethylene glycol and tetraethylene glycol.

In the case of real molecules, the dependence of evaporation coefficient A on molecular parameters is not as straightforward as in computer simulations. In MD simulations we could independently change the masses, interaction parameters and sizes of spherical molecules and in this manner we find A as a function of each of these parameters separately. In real molecules their masses, sizes and interaction parameters are inter-dependent. Eqn (7) shows which of the molecular parameters most strongly influence the value of A. Here we further verify eqn (7) in experiments.

We estimate the influence of both molecular mass and interactions upon A in a series of glycols: ethylene glycol (EG; 99.9% SPECTRANAL, GC, Riedel-de Haën), diethylene glycol (DEG; 99.99% BioUltra, GC, Fluka), triethylene glycol (TEG; 99.96% BioUltra, anhydrous, GC, Fluka) and tetraethylene glycol (TTEG; 99.7% puriss., GC, Fluka). The molecules of these four compounds have similar chemical structures, but differ in mass and interaction parameters ε₁₁.

Liquids of all these compounds were introduced into the droplet injector with a sterile syringe. Next the droplets were injected into the electrodynamic trap. They were charged during injection and thus could be trapped by appropriate combination of AC and DC electrodynamic fields. The initial droplet radius was several μm. In the electro-dynamic trap single droplets levitated in a small volume set by the trapping geometry and parameters of the system. The evolution of the droplet radius R(t) in time was obtained via the recorded angle-resolved light scattering and analysis based on the Mie scattering theory. Details of the experimental set-up are given in our previous publication.⁸

In Fig. 5 we show a radius of a droplet of triethylene glycol (very small droplet) and of ethylene glycol (a slightly large droplet) as a function of time. The data were fit to solution of eqn (1) (see eqn (5)). There were two free parameters in the fit: the kinetic parameter, A, and the temperature difference between the vapour temperature far from the droplet, T_b, and the droplet temperature, T_liq. The second of the fitting parameters was important for all sizes of studied droplets. The parameter A was important only for small droplets i.e. at times when the droplet size decreased to a value comparable to Aλ (where the mean free path in nitrogen was λ = 66 nm). Therefore in order to find both T_b − T_liq and A reliably, a droplet, which during evaporation reached a size below ∼2 μm, was selected. The presence of impurities influences the results very severely^1,12,13 and therefore we used pure liquids. The ballistic mode of evaporation is utterly negligible for droplets larger than ∼6 μm. Thus, fitting the model to the radius evolution above 6 μm was used to find the T_b − T_liq parameter only (see the lower panel of Fig. 5). Extending the fit to further evaporation enables finding both T_b − T_liq and A reliably. T_b − T_liq should not change in the extended fit and serves as a check against the influence of impurities. In short, after the fit of T_b − T_liq we fixed its value and for longer times (and thus smaller droplets) fitted only A.

Fig. 5 Upper panel: The radius R(t) as a function of time for triethylene glycol droplets evaporating into the dry nitrogen atmosphere. The following parameters characterized the system: Δh = 80 kJ mole⁻¹, ρ_liq = 1121 kg m⁻³, κ_v = 0.0258 J m⁻¹ s⁻¹ K⁻¹, λ = 66 nm. The blue line is a fit with T_b − T_liq as the only fitting parameter and A = 0. The red line is a fit with suitably chosen A and T_b − T_liq (see eqn (5)). Lower panel: The radius R(t) as a function of time for ethylene glycol droplets evaporating into the dry nitrogen atmosphere. The following parameters characterized the system: Δh = 62 kJ mole⁻¹, ρ_liq = 1109 kg m⁻³, κ_v = 0.0258 J m⁻¹ s⁻¹ K⁻¹, λ = 66 nm. The red line is a fit with suitably chosen A and T_b − T_liq in eqn (5). The blue line shows the region of droplet sizes affected only by the temperature difference and not by the kinetic parameter, A.

In Fig. 6 we show the dependence of A on molecular mass, temperature and intermolecular interactions of liquid molecules. The interaction between molecules is proportional to the enthalpy of evaporation ε₁₁ ∼ Δh (Δh is given in the ESI† for all studied compounds). The dependence of A(k_BT/ε₁₁) is nearly linear as a function of m₁/m₂ in very good agreement with MD results presented in eqn (7). The value of A for TTEG (tetra-ethylene glycol) is loaded with relatively high uncertainty since a TTEG droplet is more prone to contamination (e.g. with water) due to its hour-long evaporation. The dependence of A upon temperature, as predicted by MD calculations, is hardly accessible in our experiments. For the 273–298 K temperature range studied experimentally the square root of temperature changes is ∼4% while the accuracy of A is at least several %.

Fig. 6 Parameter A as a function of molecular mass, m₁, divided by the mass of nitrogen, m₂. The interaction energy ε₁₁ = 2Δh/(sN_A), where Δh is the enthalpy of evaporation¹⁴ (per mole; see also the ESI†) at 273 K, s is the number of nearest neighbours in the liquid (we set s = 4 from computer simulations of ethylene glycol (EG)¹⁵) and N_A is the Avogadro number. The dashed line is the linear fit to eqn (7). From the fit to eqn (7) we find α₁ = 0.3 and α₀ = −0.2, to be compared with α₁ = 0.42 to 0.65 and α₀ = 1.40 obtained in simulations of LJ fluid. The Lennard-Jones liquid of spherical molecules studied in MD simulations (Fig. 3) is very much different from the linear chain molecules studied experimentally. Nonetheless eqn (7) works very well for both systems giving very similar values of α₁ (compare Fig. 6 with Fig. 3). However α₁ and α₀ are not universal and vary from system to system depending probably on subtle differences in the range of intermolecular interactions and the structure of the liquid. In general we do not have a physical interpretation of these parameters. For example for water we found A(k_BT_b/ε₁₁)^1/2 = 1.8(m₁/m₂ = 0.64) and for glycerol A(k_BT_b/ε₁₁)^1/2 = 1(m₁/m₂ = 3.3). More simulations and experiments on different compounds are needed to elucidate the origin of α₁ and α₀. Nonetheless, since both α₁ and α₀ are of the order of 1 eqn (7) seems to grasp most of the physical ingredients of ballistic transport of energy in the system and its influence on the evaporation flux.

Conclusions

Kinetic effects and ballistic transport of energy is being currently studied in computer simulations.^16–25 The main difficulty in the quantitative prediction of the evaporation flux is the appropriate estimate of the energy flux from the gas to the surface of the liquid in the close vicinity of the interface. In the ballistic region close to the interface gas molecules individually collide with molecules in the liquid, but the collisions involve all nearest neighbours of a given molecule in the liquid and become multiple particle processes that cannot be simply reduced to binary collision. The finite range of intermolecular-interactions couples many liquid molecules during collisions with a gas molecule. Therefore we cannot describe this effect as a pure binary collision with well-known energy transfer. The phenomenological formula for A (eqn (7)) contains information about such processes. The practical message from our study is that slowly evaporating liquids consisting of high molecular mass are affected by the ballistic mode of energy transfer even for large, micro-meter size droplets. In the nanometric size range of droplets, evaporation of all liquids is very strongly influenced by parameter A. Eqn (1) together with eqn (7), the main result of this work, should apply to evaporation of all liquids from LJ liquids to water and for all sizes of droplets. The details of molecular structure affect the evaporation process. In our equations interactions between liquid molecules are mainly included in the enthalpy of evaporation, Δh, (eqn (1)) and in parameter A (eqn (7)). The origin of two parameters (α₀ and α₁) in eqn (7) is not known. Most probably they depend on the range of interactions and structure of the liquid and tell us how well molecules in the liquid are coupled during collisions with gas molecules. Therefore eqn (7)i.e. dependence of A on molecular mass, can be used for a series of compounds having similar chemical structures and the corresponding structure of the liquid as we have shown for a series of Lennard-Jones liquids (Fig. 3) and a series of glycols (Fig. 6).

The temperature discontinuity, ΔT, across an interface crossed by a thermal flux, j, is related to the Kapitza resistance, R_k = ΔT/j. The temperature discontinuity at the liquid droplet–gas interface (eqn (6)) (first observed by Ward and coworkers^2–6 and explained by Bond and Struchtrup²⁶) allows us to quantify the Kapitza resistance using energy flux from eqn (1). We find the Kapitza resistance for our system: R_k = Aλ/κ_v. The typical values of R_k are 10⁻⁵ to 10⁻⁶ m² K W⁻¹; for water 1.3 × 10⁻⁵, glycerol 10⁻⁵, TEG 1.3 × 10⁻⁶, TTEG 2.5 × 10⁻⁶ m² K W⁻¹ and similarly for the Lennard-Jones liquid.

The ballistic mode of energy transfer discussed in this paper is important for various phenomena, when characteristic sizes in the system become comparable to the mean free path in the gas. For example this mode of energy transfer arises when gas bubbles are formed in liquid by e.g. laser heating²⁷ or for the flow of gases through nanoholes of sizes comparable to the mean free path.²⁸

In this paper we have not discussed the issue of nucleation and growth just after the nucleation. All studied systems (LJ or glycols) were under thermodynamic conditions far from the region where nucleation plays a role.²⁹

Acknowledgements

This work was supported by the National Science Centre, Poland, under grant 2014/13/B/ST3/04414.

Notes and references

1. R. Hołyst, M. Litniewski, D. Jakubczyk, K. Kolwas, M. Kolwas, K. Kowalski, S. Migacz, S. Palesa and M. Zientara, Rep. Prog. Phys., 2013, 76, 034601.
2. A. Persad and C. Ward, Chem. Rev., 2016, 116, 7727–7767.
3. G. Fang and C. A. Ward, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1999, 59, 417.
4. C. A. Ward and G. Fang, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1999, 59, 429.
5. G. Fang and C. A. Ward, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1999, 59, 441.
6. A. J. H. McGaughey and C. A. Ward, J. Appl. Phys., 2002, 91, 6406.
7. R. Hołyst and M. Litniewski, Phys. Rev. Lett., 2008, 100, 055701.
8. R. Hołyst, M. Litniewski, D. Jakubczyk, M. Zientara and M. Wozniak, Soft Matter, 2013, 9, 7766–7774.
9. M. P. Allen and D. Tildesley, Computer Simulations of Liquids, Oxford University Press, 1987.
10. L. Verlet, Phys. Rev., 1967, 159, 98–103.
11. P. E. Kloden and E. Klaten, Numerical Solutions of Stochastic Diferential Equations, Springer, Berlin, 1992.
12. F. Muller-Plathe, J. Chem. Phys., 1997, 106, 6082.
13. D. Jakubczyk, M. Kolwas, G. Derkachov, K. Kolwas and M. Zientara, Acta Phys. Pol., A, 2012, 122, 709–716.
14. J. B. Pedley, R. D. Naylor and R. B. Kirby, Thermochemical Data of Organic Compounds, Chapman and Hall, New York, 1986, pp. 1–792.
15. A. V. Gubskaya and P. G. Kusalik, J. Phys. Chem. A, 2004, 108, 7151–7164.
16. D. Albernaz, D. M. Do-Quang and G. Amberg, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2015, 91, 043012.
17. N. M. Kovalchuk, A. Trybala and V. M. Starov, Curr. Opin. Colloid Interface Sci., 2014, 19, 336.
18. T. Ishiyama, S. Fujikawa, T. Kurz and W. Lauterborn, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2013, 88, 042406.
19. M. Heinen, J. Vrabec and J. Fischer, J. Chem. Phys., 2016, 145, 081101.
20. A. Lotfi, J. Vrabec and J. Fischer, Int. J. Heat Mass Transfer, 2014, 73, 303–317.
21. J. Zhang, F. Mueller-Plathe, M. Yahia-Ouahmed and F. Leroy, J. Chem. Phys., 2013, 139, 134701.
22. S. Cheng, J. B. Lechman, S. J. Plimpton and G. S. Grest, J. Chem. Phys., 2011, 134, 224704.
23. O. Wilhelmsen, T. T. Trinh, S. Kjelstrup and D. Bedeaux, J. Phys. Chem. C, 2015, 119, 8160–8173.
24. R. Ledesma-Aguilar, D. Vella and J. M. Yeomans, Soft Matter, 2014, 41, 8267–8275.
25. R. Hołyst, M. Litniewski and D. Jakubczyk, Soft Matter, 2015, 11, 7201–7206.
26. M. Bond and H. Struchtrup, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2004, 70, 061605.
27. J. Lombard, T. Biben and S. Marabia, Nanoscale, 2016, 8, 14870–14876.
28. M. Savard, C. Tremblay-Darvea and G. Gervais, Phys. Rev. Lett., 2009, 103, 104502.
29. J. L. Katz and J. Ostermier, J. Chem. Phys., 1967, 47, 478–487.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c7sm00804j

---