Revisiting $D2$-law for the evaporation of dilute droplets

The prediction of droplet evaporation in turbulent flows is crucial in science and applications, e.g., in internal combustion engines¹ or respiratory flows.² The outbreak of the COVID-19 pandemic is highly increasing the scientific relevance of this topic:^3–5 Virus-laden droplets dispersed into the environment by an infected person spread the disease,⁶ but their dispersion range, although at the base of health guidelines, is still scientifically debated.^2–5 For practical estimates, the prediction of droplet dispersion is still based on the work of Wells,⁷ which is founded on the $D2$-law to compute droplet evaporation rate. The $D2$-law, proposed in the seminal works of Langmuir⁸ and others,^9,10 states that the surface of an evaporating droplet decreases linearly with time, at a rate fixed by the ambient properties. Although deviations from the classical $D2$-law have been observed in some peculiar conditions, e.g., for droplet sizes of the order of the mean vapor's free path,¹¹ trans-critical evaporation of nano/micro-droplets¹² and multi-component droplets,^13,14 a linear decrease in the droplet surface is usually observed for single-component micro/millimetric droplets. Nonetheless, recent studies on respiratory events^3,4 found that droplet lifetime may increase even 150 times as compared to Wells's estimate; although a linear decrease in the droplet surface could still be observed, the evaporation rate is found to be much lower than that predicted by the $D2$-law. Considering that the temperature of an evaporating and isolated droplet tends to an asymptotic temperature lower than the environmental one, we propose a revision of the $D2$-law for single component, millimetric/micrometric droplets. This effect, caused by the balance between heat flux and latent enthalpy, is not accounted for in the classical formulation of the $D2$-law, which assumes a fast mixing of the gas around the droplet. We show that, by including this effect in a revised $D2$-law, the droplet evaporation time can be accurately estimated in dilute conditions. To this purpose, the revised and the classical $D2$-laws are compared to reference data obtained from direct numerical simulations (DNSs) under a hybrid Eulerian–Lagrangian framework and the point-droplet approximation.^15–17 The temporal evolution of droplets smaller than the smallest flow length scale can be Lagrangianly described using the so-called point-droplet approximation.^15,17 Indeed, droplet temperature and radius equations read

where T_d is the droplet temperature, r_d the radius, and T_m the temperature of the carrier mixture evaluated at the droplet center. The variable $τd=2ρlrd2/(9ρ¯ν¯)$ is the droplet relaxation time, whereas the Schmidt and Prandtl number are $Sc¯=ν¯/Dg,v$ and $Pr¯=ν¯ρ¯cp¯/k¯$⁠, with $ν¯$ the kinematic viscosity, $ρ¯$ the density, $c¯p$ the isobaric specific heat capacity, and $k¯$ the thermal conductivity of the gas. The mass diffusivity is $Dg,v$⁠, whereas c_l and ρ_l are the specific heat capacity and density of the liquid, and $ΔHv$ is the latent heat of vaporization (over-bars refer to Eulerian quantities evaluated at droplet centers). In the model, a fast conductivity is implicitly assumed, i.e., the liquid thermal conductivity is assumed to be much higher than that of the carrier gas (low Biot number¹⁸) such that the temperature inside a droplet is uniform but time-varying. In analogy with $(Tm−Td)$⁠, which is the forcing term for convective heat transfer, the term H_m drives the mass transfer rate:

with $Yv,m(χv,m)$ the vapor mass fraction evaluated at the droplet center and $Yv,d(χv,d)$ the mass fraction of a saturated gas–vapor mixture evaluated at the temperature of the droplet, T_d. The mass fractions depend on the vapor molar fractions which, in turn, are related to the vapor pressure as follows:

where p_a is the ambient pressure, $psv(T)$ the saturated vapor pressure evaluated at T, and RH_m the relative humidity of the moist carrier evaluated at the point-droplet center. The saturated vapor pressure, $psv(T)$⁠, depends on temperature and pressure and can be calculated by using the Clausius–Clapeyron equation. Finally, the Nusselt and Sherwood numbers are correlated with the droplet Reynolds number via the Frössling correlations, see, e.g., Wang et al.¹⁹ for details.

Based on Eqs. (1) and (2), for an isolated droplet in an environment at temperature T_m = T_a and relative humidity RH_m = RH_a, the classical formulation of $D2$-law⁸ can be derived by assuming a fast mixing of the droplet atmosphere with the environment. Under this hypothesis, the quantities regulating the evaporation of the droplet are set by the bulk thermodynamic state of the ambient (⁠$ρ¯=ρa, ν¯=νa,S¯c=Sca$⁠, and $P¯r=Pra$⁠). With the droplet temperature fixed to T_d = T_a, Eq. (2) can be exactly integrated as

where $rd,0$ is the initial droplet radius and K is a function of the bulk properties of the ambient T_a, p_a, and RH_a:

Then, the droplet evaporation time can be estimated as $td,e=rd,02/K$⁠. Equation (5) is historically referred to as $D2$-law since it predicts a quadratic temporal evolution of the droplet diameter. It assumes that the droplet temperature is fixed and equal to the environmental one, $Td≃Ta$⁠; in turn, this condition implicitly assumes that the convective heat transfer term in the droplet temperature equation, Eq. (1), dominates the latent heat term. Nonetheless, in real cases, this condition is valid only for a short time, whose duration scales with the droplet relaxation time as sketched in the inset of Fig. 1. After a time $t≃6τd,0$⁠, the droplet temperature sets to a constant value, $Td,s$⁠, which is lower than the environmental temperature, T_a. The ratio between T_a and $Td,s$⁠, determined by the balance between convective heat transfer and latent heat absorption, depends on the ambient temperature itself and relative humidity, as shown in Fig. 1. The droplet temperature is similar to the ambient one only for $t≪τd,0$⁠, whereas, for the major part of the vaporization process, it keeps closer to the asymptotic temperature, $Td,s$⁠. Since the prediction of the droplet lifetime using the $D2$-law is highly sensitive to the droplet temperature, the latter may induce a significant overestimation of the evaporation rate.^3,19 Hence, we propose to revise the $D2$-law by assuming a constant droplet temperature equal to the droplet asymptotic temperature and not to the ambient one: $Td≃Td,s$⁠. Under these conditions, Eqs. (1) and (3) become

which can be easily solved for $Td,s$⁠. The latter can be employed to compute a revised decay constant, K_r,

leading to a revised evaporation time, $td,e=rd,02/Kr$⁠. Figure 2 shows the temporal evolution of the square diameter of a water droplet in a quiescent environment at $Ta=50 °C$ and RH_a = 0, providing the numerical solution of Eq. (2) and the predictions by the classical and revised $D2$-laws. In the initial evaporation phase, $t≪τd,0$⁠, since $Td≃Ta$⁠, the vaporization rate predicted by the classical $D2$-law is close to the reference one [numerical solution of Eq. (2)]. Nonetheless, for $t≫τd,0$⁠, the behavior predicted using the revised $D2$-law guarantees a much better agreement with the reference data. Figure 3 shows the value of the decay constants K and K_r for water droplets vs ambient temperature at different relative humidities. The relative difference between K_r and K increases with T_a and reduces with RH_a. It is worth remarking that this difference spans over nearly one order of magnitude, leading to a significant mismatch in the evaporation rates estimated by using the bulk ambient temperature or the droplet asymptotic temperature.

To test the proposed revision of the $D2$-law, we consider a turbulent evaporating jet-spray. The well-established framework of the point-droplet equations is used to predict the temporal evolution of droplet temperature and radius [Eqs. (1) and (2)] together with the Lagrangian equations of motion.^15,19,20 The droplet equations are coupled with the low-Mach number asymptotic expansion of the Navier–Stokes equations, which are directly solved to reproduce the dynamics of the carrier flow. A fully developed turbulent jet, generated from a companion DNS of a turbulent pipe flow²⁰ at $Re0=2U0 R0/ν=6000$⁠, is injected into an open, quiescent environment (⁠$U0=9.3 m/s$ and $R0=4.9×10−3 m$ the pipe bulk velocity and radius). Both the jet and the environment gases consist of dry air at temperature $Ta=20 °C$⁠. The computational domain is a cylinder extending for $2π×22R0×70R0$ in the azimuthal, θ, radial, r, and axial, z, directions and is discretized by a staggered grid of $Nθ×Nr×Nz=128×223×640$ nodes. A convective boundary condition is prescribed at the outlet section, and a traction-free boundary condition is imposed on the domain side boundary to allow entrainment. A monodisperse population of water droplets of initial radius $rd,0=4 μm$ and temperature $Td,0=Ta=20 °C$ is randomly distributed over the inflow section at each time step. Four simulations have been performed changing the liquid mass fraction at the inflow: $ψ0=ml/mg=0.0007; 0.0028; 0.0058; 0.01$⁠, with m_l and m_g the mass of the liquid and gas injected per unit time. The latter correspond to $8.7×10−7; 3.5×10−6; 7.0×10−6; 1.2×10−5$ in volume fraction, respectively. Two additional simulations were conducted suppressing the numerical integration of Eqs. (1) and (2): In the former, the classical formulation of the $D2$-law, Eq. (5), is used to evolve the droplet radius; in the latter, the revised $D2$-law is used, Eq. (10). In all cases, the position and velocity of droplets are numerically integrated according to the standard equation of point-droplet.^19,20 For a detailed description of the methodology and tests, the reader is referred to Refs. 19–21. It is worth remarking that droplet number density, mass, and volume fractions are related parameters to characterize the spray dilution, governing the droplet collision/coalescence rate, mutual interaction via emitted vapor, and droplet momentum back-reaction. Since the mass fraction is the parameter that mainly controls the vapor concentration and determines the balance between heat and latent enthalpy fluxes during evaporation, we selected this parameter to characterize the spray dilution in the present study. Figure 4 provides the contours of the normalized mean liquid mass fraction and mean droplet diameter in the spray for the simulations based on the classical $D2$-law, the revised one, and on the full model for the most dilute case, $ψ0=0.0007$⁠. The average is computed considering a time range of $100t0$ (200 independent snapshots) after the establishment of a statistically steady regime,²⁰ with $t0=R0/U0$ the jet advection timescale. Both the distributions of the mean mass fraction and diameter computed using the revised $D2$-law model are in excellent agreement with the reference data. Besides, the classical $D2$-law results in a much shorter spray evaporation length, due to the overestimation of the vaporization rate. The evaporation length of a turbulent jet-spray, z_e, can be defined as the distance from the inflow where the 99% of the injected liquid mass has evaporated. This distance corresponds to a mass fraction level of $ψ/ψ0=10−2$ and, hence, to a mean droplet radius $rd,99/rd,0≃(ψ/ψ0)1/3≃10−2/3$⁠. By neglecting possible two-way coupling effects, an estimate of the evaporation length can be obtained by considering the self-similar behavior of the mean jet centerline velocity, $Uz,c$ and supposing that the droplet mean axial velocity is approximately $uz,d≃Uz,c$⁠,

with $B≃6$ a universal constant, z₀ the so-called jet virtual origin,²⁰ and z_d the droplet position along the jet axis. Then, by assuming a vanishing virtual origin, $z0/R≃0$⁠, and considering that, according to Eqs. (5) and (10), the time required for a droplet to reduce its radius from $rd,0$ to $rd,99$ is $td,e=(rd,02−rd,992)/K$⁠, the integration of Eq. (11) along z leads to

where K is either the decay rate defined according to Eq. (6) or K_r defined in Eq. (9). By using Eq. (12), the resulting evaporation lengths are $ze/R0≃28$ and $ze/R0≃46$⁠, computed by using K and K_r, respectively. These values are consistent with the axial locations where 99% of the injected liquid mass has evaporated according to Fig. 4 (⁠$ψ/ψ0≃10−2$⁠) showing the superior estimate of the proposed model. To better quantify the evaporation length, the upper panel of Fig. 5 provides the mean liquid mass fraction computed along the jet centerline. In the near-field, up to 10 jet radii away from the inlet, all the curves are close to each other. In this region, the actual droplet temperature is close to the environmental one, $Td≃Ta$⁠; being the droplet and ambient temperature similar, the actual evaporation rate is close to the one imposed by the environmental conditions. Moving downstream, the droplet temperature decreases and the classical $D2$-law underestimates the mean liquid mass fraction, since droplets evaporate at a constant rate driven by the ambient temperature. On the other hand, the simulation based on the revised $D2$-law better approximates the actual droplet temperature, providing a more accurate estimate of the mean evaporation rate and axial evolution of the mean mass fraction. For the latter, the agreement with all the reference simulations is excellent, except for the far-field of the higher mass fraction cases, $ψ0>0.0028$⁠. In non-dilute cases, the mutual droplet interactions, via the emitted vapor, cause a deviation of the local evaporation rate from the one predicted by the $D2$-laws. Similar results are obtained for the radial profiles of the mean mass fraction provided in the lower panel of Fig. 5 and computed at $z/R0=10$ and $z/R0=25$⁠. At $z/R0=10$ all curves are still similar, although the simulation employing the classical $D2$-law shows differences with respect to reference data. Moving downstream, the underestimation of the liquid mass fraction resulting from the classical $D2$-law is more apparent, whereas the revised $D2$-law still provides an excellent estimate.

Considering Lagrangian statistics, the top panel of Fig. 6 shows the mean droplet square diameter, $dd2/dd,02$ vs the flight time, $td,f/t0$⁠, calculated since the injection. In all cases, the evolution of $dd2/dd,02$ computed from full DNSs is compared with the estimates from the classical and revised $D2$-laws. The predictions obtained from the proposed revised model appear much more accurate with respect to the classical formulation. For mass fraction equal, or below, to $ψ0=0.0028$ the agreement is excellent, whereas the classical $D2$-law predicts an evaporation time shorter than the half with respect to the reference DNS. Increasing the injected mass fraction, the mutual interactions of droplets slow down the evaporation. From the figure, we conclude that the proposed revised $D2$-law optimally describe droplet evaporation in dilute conditions up to mass fractions of the order of $ψ0∼10−3$⁠. The effect of the mutual interactions of the evaporating droplets is better highlighted in the lower panels of Fig. 6, which provides the joint probability density functions (JPDFs) of $dd2/dd,02$ vs $td,f/t0$ for the lowest and highest mass fraction cases, $ψ0=0.0007$ and $ψ0=0.01$⁠. It is worth noting that, in the most dilute conditions, the JPDF is superimposed on the straight line determined by the revised $D2$-law: Each droplet evaporates as an isolated one, which is perfectly described by the proposed model. On the contrary, for the highest mass fraction, the JPDF spreads over a wider range and the evaporation times are longer than that predicted by the proposed model (around 35% in mean). This indicates that droplets exhibit different Lagrangian histories and, as expected, the proposed model cannot account for this complex interaction, typical of non-dilute conditions.

To highlight the generality of the present study, we employed our model against the dataset from recent independent works by Ng et al.⁴ and Chong et al.³ on respiratory droplets evaporation. We focus on their numerical simulation on droplets expelled in a turbulent respiratory puff (cough) at $Ta=30 °C$ and $RHa=90%$⁠. The authors observed a linear decrease in the droplet surface but found an evaporation rate much lower than that estimated by the classical $D2$-law. In Fig. 7, we compare the mean droplet square diameter, $dd2/dd,02$⁠, vs mean droplet flight time obtained by Ng et al.⁴ with that predicted by using the classical and the revised $D2$-laws. Our revised formulation accurately reproduces the temporal evolution of mean square droplet diameter, whereas the classical prediction strongly overestimates the vaporization rate, leading to shorter evaporation times. It is worth noting that the additional cases presented in Ng et al.⁴ pertain to conditions where the emitted gaseous jet strongly differs from the environmental ones, so the present model cannot be directly applied.

Summarizing the applicability limits of the present model, the droplet size should be smaller than the smallest flow length scales (point-droplet approximation) and well above the non-continuum length scales, i.e., nano-droplets.¹¹ Hence, we can state the model is applicable for droplet size below millimeters and above a micrometer. In addition, beyond a dilute mass fraction (⁠$ψ0<10−2$⁠), the basic assumptions behind the presented revision of the $D2$-law require the thermodynamic ambient conditions to be sufficiently homogeneous and similar to that of the gaseous jet. These assumptions are also at the base of the classical $D2$-law.

To conclude, the present study demonstrates that reformulating the $D2$-law by using the asymptotic temperature of isolated, evaporating droplets instead of the ambient one allows a superior description of the evaporation of dilute, millimetric/micrometric, and single component droplets in turbulent flows with respect to the classical $D2$-law. An excellent agreement of the predictions by our revised model against reference DNS is found for jet-sprays up to a liquid mass fraction of the order of $ψ0∼10−3$⁠, whereas the agreement is still acceptable up to $ψ0∼10−2$⁠. The authors believe that the proposed revision of the $D2$-law will contribute to improve practical estimates of droplet evaporation times. In the context of respiratory droplet dispersion, we expect that present findings could be significant for a revision of the Wells theory,⁷ which is based on the classical $D2$-law.

J.W. is grateful to the China Scholarship Council (CSC) for supporting this research (No. 201806250023). We acknowledge PRACE for providing computational resources via the COVID-DROPLETS grant.