Abstract
The ambient conditions surrounding liquid droplets determine their growth or shrinkage. However, the precise fate of a liquid droplet expelled from a respiratory puff as dictated by its surroundings and the puff itself has not yet been fully quantified. From the view of airborne disease transmission, such as SARS-CoV-2, knowledge of such dependencies are critical. Here we employ direct numerical simulations (DNS) of a turbulent respiratory vapour puff and account for the mass and temperature exchange with respiratory droplets and aerosols. In particular, we investigate how droplets respond to different ambient temperatures and relative humidity (RH) by tracking their Lagrangian statistics. We reveal and quantify that in cold and humid environments, as there the respiratory puff is supersaturated, expelled droplets can first experience significant growth, and only later followed by shrinkage, in contrast to the monotonic shrinkage of droplets as expected from the classical view by William F. Wells (1934). Indeed, cold and humid environments diminish the ability of air to hold water vapour, thus causing the respiratory vapour puff to super-saturate. Consequently, the super-saturated vapour field drives the growth of droplets that are caught and transported within the humid puff. To analytically predict the likelihood for droplet growth, we propose a model for the axial RH based on the assumption of a quasi-stationary jet. Our model correctly predicts super-saturated RH conditions and is in good quantitative agreement with our DNS. Our results culminate in a temperature-RH map that can be employed as an indicator for droplet growth or shrinkage.
Significance Statement Influence of environmental conditions on airborne diseases transmission is an important issue, especially during the pandemic of COVID-19. Human-to-human transmission is mediated by the transport of virus-laden respiratory droplets. Here we investigate the problem from a fluid mechanics perspective by conducting numerical simulations to quantify the fate of respiratory droplets in a warm humid coughing puff under different ambient conditions. We reveal a non-intuitive regime with considerable growth of respiratory droplets, dominated by a super-saturated vapour field, preferentially occurring in cold and humid environments. We further propose a theoretical model that accurately predicts the condition for droplet growth. Our work should inform socializing policies and ventilation strategies for controlling indoor ambient conditions to mitigate dispersion of droplets from asymptomatic individuals.
The flow physics of respiratory droplet is crucial to understand the airborne transmission of respiratory diseases (1, 2), and recently due to the pandemic of COVID-19, the interest on this subject has been renewed (3–12). Tiny saliva and mucus droplets, mostly of the size of tens to hundreds of micrometers (13, 14), are expelled from the mouth at speaking, screaming, shouting, singing, coughing, sneezing, or even breathing. The expelled micro-sized droplets have long been recognized to be the carriers for viruses which are responsible for the viral transmission from one host to others (15).
Airborne disease transmission is a difficult subject because many intricate factors simultaneously affect the transmission and a complete understanding requires multi-disciplinary collaboration (16, 17). In addition to the understanding from epidemiology (18–20) and virology (21), fluid mechanics also offers essential insight into the mechanisms of airborne disease transmission (2). However, surprisingly little is known on the transport and fate of the respiratory droplets, once they have been expelled from the mouth. For example, how the air temperature and the relative humidity influences the actual fate of the droplets remains an open question. It is a complicated process involving interactions between turbulent eddies, vapour field, and temperature field with the respiratory droplets. Therefore, accurate quantification of these effects is necessary.
On the fate of respiratory droplets in fluid flow, classical work was done by William F. Wells in 1930s (22, 23). The picture that Wells developed, at that time in connection with the transmission of tuberculosis, considers only the settling and evaporation of the droplets. After droplets are expelled from the mouth, they settle to the ground by gravity, and, in this simple model, they evaporate with their surface area decreasing linearly with time, the so-called d2-law (24). Then for a given ambient temperature and relative humidity, one can estimate the lifetime of the respiratory droplet based on the d2-law.
However, in an actual respiratory event, the assumptions of settling and evaporation are insufficient to describe the full physical process. In reality, more complicated flows are possible: Abkarian et al. (25) demonstrated that speech can produce jet-like flow transport, and Bhagat et al. (11) illustrated that body or breathing plumes from a person can also spread droplet nuclei in enclosed spaces. Recently, a new paradigm has been suggested by Bourouiba et al. (1, 3, 26), which states that instead of the isolated droplet being considered, one should also take into account the importance of turbulent vapour puff. Inspired by this paradigm, using direct numerical simulations, Chong et al. (27) showed that droplets have O(100) longer lifetimes in a turbulent vapour puff at high ambient relative humidity than that predicted by the d2-law. Past numerical studies that consider different ambient conditions typically adopt mean flow modelling approaches (28–30). The role of ambient conditions requires large-scale parameterizations, (31–33) but run the risk of obscuring the underlying physics.
In this study, we aim at quantifying the effect of ambient temperature (θamb) and relative humidity (RHamb) on the respiratory droplets using direct numerical simulations of a turbulent puff. Normal coughing with an actual injection profile has been simulated with 5000 droplets being tracked by the Euler-Lagrangian approach. We reveal a regime in which the respiratory droplets considerably grow before their shrinkage due to evaporation. The droplets growth happens under specific combinations of ambient conditions, preferentially in cold and humid environment. To predict the criteria of having droplet growth, we further propose a theoretical model to estimate the local relative humidity experienced by the droplets, which agrees excellently to our simulation results. Whilst our model can help to formulate public health guidelines or indoor ventilation strategies in order to mitigate the spreading of respiratory droplets, the infectivity of virus-containing droplets and droplet nuclei, which remain after evaporation still requires further investigation, in particular on the virus dose required for infections and whether dried-out droplet nuclei may also be infectious. This issue remains controversial (15) and cannot be answered within a fluid dynamical study.
Results and Discussion
A. Effect of Ambient Temperature on Turbulent Vapour Puff and Droplet Counts
First, we highlight the differences in the local RH values for two different θamb (10°C and 30°C) at RHamb = 90%. Note that θamb and RHamb both constitute the control parameters of our study. The local RH field is shown in figure 1 together with the instantaneous droplet distributions and instantaneous count histogram. At θamb = 30°C, the local RH field is almost indistinguishable from RHamb, as shown in Fig. 1A. In contrast, at θamb = 10°C (Fig. 1B), the local RH field reaches super-saturation values larger than 100% (whitish colour shading). In the count histograms in Fig. 1B, most of the droplet diameters remain at least around 20 µm throughout the spread distance, whereas in Fig. 1A, droplet diameters can go below 20 µm. This means that droplets expelled into colder ambient surroundings shrink slower because of the super-saturated local RH field. On the other hand, droplets expelled into hotter ambient surroundings shrink much more rapidly because the local RH rapidly under-saturates.
To illustrate the temperature sensitivity of the droplet size distributions, in Fig. 2, the distributions across three θamb, 10°C, 20°C, and 30°C are shown at t = 0.6s (Fig. 2A) and 1.2s (Fig. 2B). In the figure, the initial droplet size distribution (solid line) is also included for comparison. In Fig. 2A, which corresponds to the cut-off time of the simulated puff, many smaller droplets ≲ 15um already exist at θamb = 20°C and 30°C as compared to θamb = 10°C. Interestingly, the histogram distribution of droplet diameters for θamb = 10°C is slightly higher and shifted to the right from the initial droplet size histogram, reflecting that a portion of the smaller droplets in fact grew. The reason for this increase is that due to the cold surrounding air, which can take up less vapour than the warm one, the puff gets super-saturated. This effect will be quantified in detail in Section B. At the later time of t = 1.2s (Fig. 2B), most of the droplet sizes for θamb = 20°C and 30°C are shifted to the left towards much smaller sizes. Remarkably, the size distribution for θamb = 10°C remains roughly similar to the initially prescribed size distribution, with some lower counts at sizes ≲ 15µm.
After having shown the distinct droplet shrinkage behavior under different ambient conditions, in the following section, we will further and in more detail quantify how the local fluid properties affect the shrinkage rate of the droplets.
B. Lagrangian Statistics of the Droplets
For the following Lagrangian statistics, we focus on the droplets with initial diameter around 15µm, given that for those we observe the most pronounced temperature effect, see the count histogram in Fig. 2B. First, we track the droplet surface area ∝d2 versus time as shown in Fig. 3A,B. Note that droplets are expelled from the mouth at different time instants until the respiratory event stops (0.6s for the coughing event considered here). In order to compare the droplet statistics, we have shifted the time frame by their respective ejecting times.
For θamb = 30°C (Fig. 3A), the normalized surface area of the droplets decreases monotonically with time (the nomenon of droplets condensation during the respiratory event ensemble average of the value is depicted by the solid curve), which indicates pure evaporation of the expelled droplets. The surface area decreases linearly with time with the averaged droplets surface area being halved after 1s. Although a linear decline of the surface area is consistent with the classical d2-law (24), one should note that the magnitude of the droplet shrinkage rate is still much smaller than that predicted by William F. Wells (22). Such a discrepancy is caused by the fact that the turbulent vapour puff engulfs small droplets and the puff contains higher RH than the ambient fluid (3, 27).
Surprisingly, for the low θamb of 10°C (Fig. 3B), the respiratory droplets evolve in a very different and distinct way: Instead of the pure evaporation during the time evolution of the droplets, there exists an initial stage (from the expelled time to about 0.4s) during which the droplets can grow with the averaged d2 increasing by 10%. This non-monotonic trend of the shrinkage markedly differs from the classical d2-law. It also contradicts the well-accepted view that the respiratory droplets should simply evaporate after expulsion. Such a phenomenon of droplets condensation during the respiratory event is the first time, to the best of our knowledge, to be quantified thanks to the Lagrangian statistics obtained from DNS, which has been neglected in the hitherto studies of droplet transmission in exhalation.
Another indicator that for 10°C the droplet size first increases is that the droplet temperature θdrop increases relative to the fluid temperature at the location of the droplet θfluid (Fig. 3D). As latent heat is released during droplet condensation for 10°C (Fig. 3D), the difference θdrop− θfluid increases during the early time. In contrast, θdrop− θfluid remains less than zero for θamb = 30°C (Fig. 3C), due to the latent heat needed for evaporation. Eventually, the droplet temperature saturates at a value lower than the fluid temperature. This is reasonable because the Jakob number, which is the ratio of the sensible heat to the latent heat due to the phase change, is much smaller than one (0.04 and 0.007 for θamb = 10°C and 30°C, respectively). Therefore, the evaporation cooling effect plays a significant role here.
To understand why the droplets can grow during the early time, a key quantity to examine is the local RH experienced by the droplet. In Fig. 3E,F, we compare the time evolution of the local RH at the droplet location for the two different θamb, showing quite different behaviours: For 30°C, the local RH decreases monotonically to RHamb which is 90% in this case. However, for 10°C, the droplets first experience super-saturation within 0.4s after expelling and the maximum RH reaches almost 115%. This duration of having supersaturation is consistent with the duration of droplet growth, which reflects that it is the local super-saturation of the surrounding air which makes the droplets grow.
In fact, our daily life experience gives us intuition about this intriguing droplet condensation phenomenon. In cold and humid weather, one may observe the “white mist” coming out from the mouth while breathing or speaking. The physical explanation is that: Warm air can contain more moisture than cold air. Therefore, the exhaled warm vapour puff becomes super-saturated when it enters the cold ambient fluid, and the droplet nuclei and dust in air favours vapour condensation. Such super-saturated conditions become more prominent if the temperature difference between the exhaled vapour and ambient fluid is larger. Here, we have clearly demonstrated the significance of super-saturation in the respiratory event that leads to the possibility of vastly different droplet dynamics.
C. Theoretical Model Predicting Relative Humidity Profiles and Super-saturation Criteria
Motivated by the strong effect of the local RH field on the droplets, here, we propose a simple model to calculate the local RH. Several assumptions first need to be made to justify this model. First, the axial scalar quantities (vapour mass fraction and temperatue) are the dominant fluid properties that determine droplet growth/evaporation, and secondly, the exhaled puff at the time of t = 0.2s exhibits jet-like properties such that the quasistationary mean properties of the puff admits self-similarity. Here, we employ the Antoine equation in order to compute the saturated vapour mass fraction for a given temperature field (34). In short, given the local values of vapour mass fraction and temperature along the jet centreline, the axial RH can be directly computed.
Following the theory for axisymmetric turbulent jets (35– 37), the axial scalar concentration distribution (say Φaxial) can be approximated by the relation, where Φinlet is the inlet scalar concentration value, and x/dinlet is the axial distance normalised by the inlet diameter. Since this model diverges as x/dinlet→ 0, we further impose the requirement that max[Φaxial/Φinlet]= 1. Both the exhaled vapour mass fraction and temperature are modelled using Eq. (1) along the axis normal to the inlet.
To assess the validity of our assumptions, we plot the resulting model of the axial RH profiles for different θamb in Fig. 4. The agreement with our DNS (lighter curves) is quite promising, which is consistent with our assumption that the puff retains self-similar jet-like characteristics at cut-off time. This agreement lends some confidence to our ability to model the axial RH of our flow at a specific time instant, while the cough will approximate to a turbulent puff or an interrupted jet at large t after the injection velocity becomes negligible (38).
Based on the robustness of this model, we further propose the use of the maximum local RH as a metric to indicate the likelihood of the growth for small droplets. To illustrate this point, in Fig. 5, we plot the axial RH map for various θamb versus distance at the fixed RHamb of 90%. The dashed line shows the crossover between super-saturated and undersaturated axial RH values. One interesting observation from this model is that with decreasing θamb, the overall axial distance that experiences super-saturation increases. In other words, it is more likely that with decreasing θamb (at RHamb = 90%) small droplets will encounter longer super-saturated distances and grow in size. This likelihood to grow is important in understanding how far the droplets can be advected by the puff in the initial stages of the cough, before evaporating and dispersing.
Finally, we extend our idea of modelling the axial RH to map out the maximum local RH values for a range of RHamb and θamb. This map is shown in Fig. 6. From the figure, we find that the maximum local RH values are dependent on both RHamb and θamb. Several interesting observations can also be made from the figure. For instance, at a given indoor ambient conditions of RHamb = 40-50% (within acceptable occupancy comfort levels), the maximum local RH is at most 100% as long as θamb > 20°C. However, once θamb < 20°C, the local RH tends to super-saturate (RH > 100%), implying that in cooler environments, smaller droplets within the humid puff will first undergo supersaturation-driven growth in the early stages of the cough before evaporating upon fall-out. Such conditions become particularly exacerbated in outdoor environment, for example, in stadiums or sports field during autumn or winter, when the temperatures tend to be colder. Similar conditions can also occur in crowded bus and train stations. Super-saturated conditions are more easily attainable (easily seen in exhaled mists) and therefore alters the respiratory droplet dynamics.
To substantiate our proposed model, in Fig. 6, we indicate on the plot 12 different control parameter combinations (RHamb and θamb) for which droplets with an initial size of ≈15µm grow (indicated by a circle symbol) or shrink (indicated by a cross symbol) in our DNS and compare with the predictions of the model. (The profiles of d2 versus time is provided in Fig. S1 of the Supplementary Material). By comparing the symbols from the DNS and the map from the model, a clear correlation between growth and super-saturated local RH fields can be seen. Thus we propose that this map can be used as an indicator to determine the droplet size distributions of a similar cough.
It is also worth noting that we did not consider the influence of background wind, which is an important aspect both in ventilated indoor and outdoor environments. We also stress that we make no claims to correlate growth or shrinkage of droplets with virus transmissions or infectivity. The subject of the viability of droplet sizes for viral media is still the topic of ongoing investigations in other communities and are not the focus of this work.
Conclusions
Through the Lagrangian statistics of the respiratory droplets, we have identified that when droplets are expelled from the mouth to the cold environment with high enough relative humidity, the droplet will first grow instead of immediately shrink, because the turbulent vapour puff becomes super-saturated. More importantly, we have provided a theoretical framework to accurately predict (i) how the relative humidity varies with distance in a respiratory puff, and (ii) the threshold of ambient temperature (θamb) and relative humidity (RHamb) for which the super-saturated vapour field can result. Our theoretical model can also be applied to other respiratory events whenever there is jet-like transport, such as speaking (25), which is extremely relevant for asymptomatic and presymptomatic spreading of the coronavirus. In the situation of continuous speaking, when the hot vapour is continuously injected into the cold environment, we expect that the super-saturation zone may remain in front of the mouth for a long period of time. Our finding implies that tiny droplets of the size of tens of micrometers can propagate further away when the weather is cold and humid due to the protection from the super-saturated vapour field. The result of the very distinct fate of respiratory droplets under different RHamb and θamb may be instructive for developing timely strategies in mitigating the COVID-19 pandemic in different seasons. For example, during the Winter time, where the maritime climate for many coastline European countries have much higher RHamb, the transmission in outdoor environments can presumably not be overlooked. However, we note that obviously also non-fluid-dynamical aspects affect the transmission and infectivity of viruses, such as from virology, epidemiology, and sociology, to name a few. These are not covered within our study and presumably interfere with our fluid dynamical results.
Materials and Methods
Governing equations
We consider an incompressible flow (∂ui/∂xi = 0) of gas phase, with both temperature and vapour concentrations coupled to the velocity field by employing the Boussinesq approximation. The governing equations read: For droplets, we apply the spherical point-particle model, and consider the conservation of momentum (Maxey-Riley equation (39)), energy, and mass as follows: The notations we used in equations are as follows: ui, ui,n, and ui,g,n are the velocity of gas, droplets, and gas at the location of droplets, respectively. θg, θn, and θg,n are in Kelvin and used to represent the temperature of gas, droplets, and gas at the location of droplets, respectively. c is the vapour mass fraction, rn the droplet radius, An the surface area of the droplets, Vn the volume of the droplets, mn the mass of the droplets. Also p denotes the reduced pressure. hm is the heat transfer coefficient. L is the latent heat of vaporisation of the liquid. ρl and cp,l are the density and specific heat capacity of the droplets assumed to consist of water. νair is the kinematic viscosity of air. kg is the thermal conductivity of gas which relates to the thermal diffusivity of gas Dg by kg Dg ρg≡ cp,g. ρg and cp,g are the density and specific heat capacity of gas, Note that ρg cp,g = (ρacp,a + ρv cp,v), with ρa and ρv being the densities of air and vapour and cp,a and cp,v being the specific heat capacities of air and vapour. cdrop and cfluid denote the vapour mass fractions of the droplet and the fluid at the location of the droplet. β is a dimensionless measure of the droplet density relative to the fluid density, and is defined as β ≡ 3ρg /(ρg + 2ρl) − 1. fd is the prefactor for the drag corrections defined as (40).
The estimations of hm and Sh drop for a single spherical droplet are given by the Ranz-Marshall correlations (41): where we have droplet Reynolds number To calculate the local relative humidity, the saturated vapour mass fraction csat,vap is determined by ideal gas law csat,vap = Psat/(ρg Rθg), where R is the specific gas constant of water vapour and Psat is the saturated vapour pressure, determined through Antoine’s relation as:
Simulation details
The computational domain size is 0.37 m (spanwise length) ×0.37 m (height) ×1.47 m (streamwise length) with corresponding grid points chosen as 512 ×512 ×2048. We model the mouth as a circular inlet centred at mid-height of the domain and the cough temporal profile we apply is a gamma distribution as where α = 60.9 s−1, such that the entire cough process lasts for about 0.6 s. Based on an experimental measurement (13), we seeded the respiratory event with droplets with initial diameters ranging between 10 µm and roughly 1000 µm. We employ a similar droplet size distribution as (13, 26) with 5000 droplets. The droplets are randomly positioned at the inlet area and evenly injected in time with the velocity matching the local inlet velocity. Droplets are assumed not to collide or coalesce since the estimates of the droplet volume fractions give O(10−6). These droplets also do not set the relative humidity of the puff, which at inflow is assumed to be 100%.
Data Availability
Data is available upon request.
ACKNOWLEDGMENTS
The authors thank J.G.M. Kuerten for suggesting the idea to study droplet condensation. The authors also thank N. Hori for the insightful discussions. This work was funded by the Netherlands Organisation for Health Research and Development (ZonMW), project number 10430012010022: “Measuring, understanding & reducing respiratory droplet spreading”, the ERC Advanced Grant DDD, Number 740479 and by several NWO grants. The funders have no role in study design, data collection and analysis or decision to publish. The simulations were performed on the national e-infrastructure of SURFsara, a subsidiary of SURF cooperation, the collaborative ICT organization for Dutch education and research, the Irish Centre for High-End Computing (ICHEC) and Irene at Très Grand Centre de calcul du CEA (TGCC) under PRACE project 2019215098.
Footnotes
C.S.N., K.L.C., R.V., and D.L. designed research; C.S.N., K.L.C., R.Y., M.L., R.V., and D.L. performed research; C.S.N., K.L.C., R.Y., and M.L. analyzed data; C.S.N., K.L.C., R.Y., M.L., R.V., and D.L. wrote the paper.
The authors declare no competing interest.