key: cord-275303-8mj8gjv9
authors: Chaudhuri, Swetaprovo; Basu, Saptarshi; Kabi, Prasenjit; Unni, Vishnu R.; Saha, Abhishek
title: Modeling the role of respiratory droplets in Covid-19 type pandemics
date: 2020-06-01
journal: Phys Fluids (1994)
DOI: 10.1063/5.0015984
sha: 
doc_id: 275303
cord_uid: 8mj8gjv9

In this paper, we develop a first principles model that connects respiratory droplet physics with the evolution of a pandemic such as the ongoing Covid-19. The model has two parts. First, we model the growth rate of the infected population based on a reaction mechanism. The advantage of modeling the pandemic using the reaction mechanism is that the rate constants have sound physical interpretation. The infection rate constant is derived using collision rate theory and shown to be a function of the respiratory droplet lifetime. In the second part, we have emulated the respiratory droplets responsible for disease transmission as salt solution droplets and computed their evaporation time, accounting for droplet cooling, heat and mass transfer, and finally, crystallization of the dissolved salt. The model output favourably compares with the experimentally obtained evaporation characteristics of levitated droplets of pure water and salt solution, respectively, ensuring fidelity of the model. The droplet evaporation/desiccation time is, indeed, dependent on ambient temperature and is also a strong function of relative humidity. The multi-scale model thus developed and the firm theoretical underpinning that connects the two scales—macro-scale pandemic dynamics and micro-scale droplet physics—thus could emerge as a powerful tool in elucidating the role of environmental factors on infection spread through respiratory droplets.

It has been well established that the SARS-CoV-2 virus responsible for the Covid-19 pandemic transmits via respiratory droplets that are exhaled during talking, coughing, or sneezing. 1 Each act of expiration corresponds to different droplet sizes and myriad trajectories for the droplets embedded in the corresponding jets. Wells 2, 3 was the first to investigate the role of respiratory droplets in respiratory disease transmission. Expelled respiratory droplets from an average human being contain dissolved salt with a mass fraction of about 0.01 as well as various proteins and pathogens in varying concentrations. 4, 5 In this paper, to model the outbreaks, we extensively use the evaporation and settling dynamics of NaClwater droplets as a surrogate model of the infectious droplets. Stilianakis and Drossinos 6, 7 included respiratory droplets in their epidemiological models. However, they neglected the droplet evaporation dynamics and assumed that characteristic post-evaporation droplet diameters are half of the pre-evaporation droplet diameters based on Nicas et al. 8 In the context of the present Covid-19 pandemic, while the role of droplet nuclei and corresponding "aerosol transmission" route are not clear, 1 it is widely accepted that respiratory droplets are definitely a dominant vector in transmitting the SARS-CoV-2 virus. This merits a detailed investigation of the evaporation dynamics of respiratory droplets and development of a pandemic model that is explicitly dependent on the respiratory droplet characteristics. As such, the evaporation mechanism of respiratory droplets are laced with complexities stemming from droplet aerodynamics, initial droplet cooling, heat transfer, mass transfer of the solvent and solute, respectively, and finally, crystallization of the solute-a phenomenon known as efflorescence. All these Physics of Fluids ARTICLE scitation.org/journal/phf are strongly affected by ambient conditions in which the droplet evaporates. These urgently necessitate a model based on first principles, which connects the detailed evaporation dynamics of respiratory droplets with the pandemic evolution equations. In this paper, a model for the infection rate constant based on collision theory incorporates the evaporation physics of respiratory droplets, ab initio. The droplet evaporation model thus developed is first validated with new experimental results obtained from droplets observed to evaporate in an acoustic levitator. While very interesting insights can be obtained from sessile droplet evaporation, [9] [10] [11] [12] [13] after an expiratory event, the floating droplet evaporates in the absence of surface contact. Thus, the levitated droplets are similar to the droplets in atmosphere [14] [15] [16] compared to their sessile counterpart. Furthermore, the desiccation dynamics necessitates a contact-less environment for the droplet. Alongside a droplet evaporation model, a chemical kinetics based reaction mechanism model is developed with final rate equations similar to that yielded by the SIR (Susceptible, Infectious, Recovered) model. 17 In general, the resemblance of the equations modeling kinetics to those of population dynamics is well known. However, the rigorous framework (analytical as well as computational) of chemical reaction mechanisms that can at present handle few thousands of species and tens of thousands of elementary reactions seems particularly attractive. 18 This could be utilized toward adding further granularity in the pandemic model, if required large mechanisms can be reduced systematically with mechanism reduction techniques. 19 Furthermore, it can be integrated into advectiondiffusion-reaction equations, and their moments could be solved using appropriate moment-closure methods. 20, 21 However, for any reaction mechanism, the key inputs are the parameters for the reaction rate constant. In our case, one rate constant is shown to be a strong function of the droplet lifetime. Therefore, next, the droplet lifetime is evaluated over a wide range of conditions relevant to the ongoing Covid-19 pandemic, and the growth rate exponents (eigenvalues) are presented. The results do not suggest that factors not considered in this paper play a secondary role in determining the outbreak spread. Rather, this paper aims to establish the mathematical connection between the pandemic and the respiratory droplet dynamics using a well defined framework rooted in physical sciences. This paper is arranged as follows: first, we provide details of the experiments used to obtain the evaporation characteristics of the water and salt solution droplets. This is followed by the reaction mechanism model that yields the equations for the growth rate and the infection rate constant of the outbreaks. This infection rate constant provides the connection and motivation for modeling the droplet evaporation time scales. Next, to evaluate the rate constant, detailed modeling of the droplet evaporation is presented. This is followed by results and discussions. Finally, we summarize the approach and findings in Sec. VII.

The experiments with isolated evaporating droplets were conducted in a contact-less environment of an ultrasonic levitator (tec5) to discount boundary effects, generally present in suspended, pendant, or sessile droplet setups. 22, 23 The experimental setup with the diagnostics is shown in Fig. 1 . A droplet was generated and positioned near one of the stable nodes of the levitator by using a micropipette. The levitated droplet was allowed to evaporate in the ambient condition of the laboratory at 30 ○ C and at about 50% relative humidity (RH). The transient dynamics of evaporation and precipitation of the evaporating droplet was captured with the shadowgraphy technique using a combination of a CCD camera (NR3S1, IDT Vision) fitted with a Navitar zoom lens assembly (6.5× lens and 5× extension tube) and a backlit-illumination by a cold LED light source (SL150, Karl Storz).

A set of ten images at a burst speed of 30 fps is acquired every 2 s for the entire duration of the droplet lifetime. The spatial resolution of the images was ≈1 μm/pixel. The temporal evolution of the diameter of the evaporating droplet was extracted from the images using the "Analyze Particles" plugin in ImageJ (open source platform for image processing). The final precipitate was carefully collected on carbon tape and observed in the dark-field mode under a reflecting microscope (Olympus BX-51). A range of initial droplet diameters varying from 300 μm to 1000 μm were investigated in experiments.

In this section, we model the infection spread rate using the collision theory of reaction rates, well known in chemical kinetics. 18 The connection between droplets and the outbreak will be established later. In this model, we adopt the following nomenclature: P represents a Covid-19 positive person infecting a healthy person(s) susceptible to infection. The healthy person is denoted by H (who is initially Covid-19 negative), and R represents a person who has recovered from Covid-19 infection and hence assumed to be immune from further infection, while X represents a person who dies due to Covid-19 infection. We consider one-dimensional head on collisions, and the schematic of a collision volume is shown in Fig. 2 . Here, one healthy person denoted by H with the effective diameter σH is approached by a Covid-19 positive person P of the same effective diameter with an average relative velocity ⃗ VDH. σH can be considered as the diameter of the hemispherical volume of air that is drawn by H during each act of inhalation, which comes out to be approximately 0.124 m.

It is widely believed that Covid-19 spreads by respiratory droplets 24 resulting from breathing, coughing, sneezing, or talking. Thus, we assume that a volume in front of P is surrounded by a cloud of infectious droplets exhaled by P. The droplet cloud is denoted by D, and the maximum cloud diameter is given by σD. Clearly, σD should be determined by the smaller of evaporation or settling time of the droplets ejected by P, the horizontal component of the velocity with which the droplets traverse, as well as the dispersion characteristics. In each such cloud, we assume that there are numerous droplets containing the active Covid-19 virus. The velocity of this droplet cloud relative to H is given by ⃗ VDH. In such a scenario, we assume that in a unit volume, there are nP infected persons and nH healthy persons. For a collision to be possible, the maximum separation distance between the centers of D (the droplet cloud) and H is given by

The collision volume-the volume of the cylinder within which a collision between the droplet cloud of P and air collection volume of H should lie for the collision to occur in a unit time-is given by πσ 2 DH VDH. Thus, the number of collisions between H and the droplet cloud D of P, per unit time per unit volume, that will trigger infections, is given by

where nP and nH represent the number of P and H, respectively. Now, given that each collision between P (basically, its droplet cloud D) and H results in conversion of the healthy individual to the infected individual, we can write

Now, we can define [P] = nP/n total and [H] = nH/n total , whereas n total is the total number of people those who are capable of transmitting the infection, as well as accepting the infection per unit volume, in that given volume. This implies

where k = n total πσ 2 DH VDH.

Here, ω is the reaction rate. Furthermore, if we assume that the mortality rate is about 3% for the ongoing Covid-19 pandemic, we can convert the kinetics of infection spread to a complete reaction mechanism given by the following:

It is to be recognized that H does not become P immediately on contact with the droplet cloud. The virus must proliferate for a finite time after contact to render a person infectious. A person who has just come in contact with the virus and does not have the capability to infect others yet is denoted by P * . k 1 , k 2 , and k 3 are the rate constants of reactions [R1], [R2], and [R3], respectively. All rate constants must have dimensions of [T] −1 (inverse of time). Clearly, k 1 > k 3 for the rapid outbreak to occur. It is to be recognized that this framework implies that k 1 , the rate constant of the second order elementary reaction [R1] resulting from collisions between the droplet cloud from an infectious individual and healthy individual, is purely controlled by physical effects. The rate constants k 2 and k 3 of the other two first order elementary reactions [R2] and [R3] are essentially decay rates emerging from the time by which the respective concentrations reach e −1 levels of the initial concentration for the respective reactions. Thus, k 2 and k 3 are purely determined by interaction between the virus and the human body. We know that the approximate recovery time from the Covid-19 disease is about 14 days. Thus, we can assume k 3 = 1/14 day −1 . We also assume the latency period (not incubation period) to be 1 day; hence, k 2 = 1 day −1 . Given the importance of k 1 in determining the outbreak characteristics, we will refer to k 1 as the infection rate constant. The major contribution of this work is imparting a rigorous physical interpretation to k 1 and calculating it.

Using Eq. (4), we can write the system of ODEs for d[P]/dt and d[P * ]/dt as

In this paper, we are interested in modeling the initial phases of the outbreaks where [H] ≫ [P]. Hence, we can safely assume [H] ≈ [H] 0 , i.e., the concentration of healthy people remains approximately constant during the early phase of the outbreak and is equal to the initial concentration, which is very close to unity at t = 0, i.e., at the onset of the outbreak. The time of the beginning of the outbreak denoted by t = 0 for a particular location can be assumed to be the day when the number of Covid-19 positive persons equaled 10.

[P] 0 is [P] at t = 0. Then, [P] can be solved as an eigenvalue problem and is given by

C 1 and C 2 are constants to be determined from the eigenvectors and the initial conditions [P] 0 and [P * ] 0 . λ 1,2 are the eigenvalues.

These can be termed growth parameters and are given by

By Eq. (5), k 1 = n total πσ 2 DH VDH. As mentioned before, k 2 = 1 day −1 and k 3 = 1/14 day −1 , which yields λ 1,2 = −0.5357 ± √ 0.2156 + k 1 . If k 2 → ∞, i.e., a healthy person becomes infectious immediately on contact with an infectious person,

Clearly, this model does not yet account for the preventive measures such as "social distancing," "quarantining" after contact tracing, and population wide usage of masks. We will call this "social enforcement." However, it can be included by accounting for the time variation in [H] . Social enforcement measures reduce the concentration of healthy, susceptible individuals from [H 0 ] to [HSE] where the concentration of healthy population susceptible to infection after implementing strict social distancing (at time t = tSE) [HSE] < [H 0 ]. In the case of social enforcement, [P] will be given by

Here, [P] = [P]SE at t = tSE and λ 1,SE , λ 2,SE are the eigenvalues from Eq. (6) with [H] = [HSE]. k 1 , the infection rate constant, remains to be completely determined. It is to be recognized that two of the key inputs of k 1 are σDH and VDH since k 1 ∝ VDHσ 2 DH by Eq. (5). As already mentioned, σH is the diameter of the hemisphere from which breathable air is inhaled. σD is the diameter of the droplet cloud. The aerodynamics of the respiratory droplets needs to be analyzed to evaluate these quantities.

The droplets ejected during respiratory events, such as sneezing and coughing, co-follow the volume of air exhaled during the event. Studies have confirmed that due to entrainment, the exhaled air volume grows in diameter, while its kinetic energy decays with time. Specifically, Bourouiba et al. 16 showed that initially, for a short duration, the droplets evolve inside a turbulent jet, while in later stages, the jet transitions to a puff. Recognizing that the ejected droplets during the respiratory event is surrounded by this dynamically evolving air volume and that the motion of the droplets will be strongly coupled due to the aerodynamic drag, we first model the surrounding air in two parts using the analytical results of the turbulent jet and puff, respectively. The axial location, axial velocity, and radial spread of a transient turbulent jet and puff can be expressed, respectively, as 25, 26 

and

where subscripts j and pf denote the jet and puff, respectively. R 0 and U 0 are the radius and axial velocities at a distance x 0 . K is a characteristic constant for the turbulent jet and is reported to be 0.457. 25 At the inception of the respiratory event (t = 0), the jet is assumed to have a velocity Uj ,0 = 10 m/s and a radius Rj ,0 = 14 mm-the average radius of human mouth. The characteristic constants for a puff are a ≈ 2.25 and m = (xp ,0 a)/(3Rp ,0 ). 26 Since the continuous ejection of air from the mouth lasts only for the duration of a single respiratory event, the jet behavior persists only for this period and beyond which the puff behavior is observed. The average duration of such events is roughly 1 s. 27 Hence, the velocity and the radial spread of the air surrounding the exhaled droplets will be

The horizontal displacement (Xp) of the exhaled droplet and its instantaneous velocity (Up) due to the drag can be solved with 14

where Rs is instantaneous radius of the droplet, ρv and ρ l are gas phase and liquid phase densities, μg is gas phase dynamic viscosity, and CD is the drag coefficient, which can be taken as 24/Rep for the gas phase Reynolds number, Rep = (2ρv|Ug − Up|Rs)/μg < 30. 14 As it will be stated later, Rep for the respiratory droplets were found to be mostly less than 0.1. By solving Eqs. (10)-(13) over the droplet lifetime, τ, the axial distance traveled by the droplets, Xp, can be evaluated. The average velocity of the droplet cloud relative to P is VD,P = Xp/τ. The diameter of the droplet cloud ejected by P can be approximated as twice the radial spread of the exhaled air, σD = 2Rg(t = τ). It is to be recognized that while the above equations are analytically tractable, given the complexities of the associated turbulent jet/puff, a detailed description of the motion of the droplets necessitates time resolved Computational Fluid Dynamics (CFD) simulations in three dimensions. This has been recently reported in Ref. 28 , which simulated dispersion of water droplets using a fully coupled Eulerian-Lagrangian technique including the wind effects. In this paper, we worked with salt solution droplets, accounting for salt crystallization, but did not include wind effects to retain analytical tractability. Nevertheless, the results presented in Subsection VI B are qualitatively consistent with the CFD results.

Due to evaporation or settling, the droplet is present only for a short time τ after it has been ejected. Therefore, the steady state k 1 can be defined as

Just like in collision theory, not all molecules are energetic enough to effect reactions; in our case, the droplet cloud is not always present. The last fraction (τ/tc) is the probability that the droplet cloud with the average diameter σD is present. tc is the average time period between two vigorous expiratory events. VDH = (VD,P + VP) + VH.

We can assume VP = VH. It is thus apparent that τ appears in σDH, VDH, and in the last fraction in Eq. (14), thereby emerging as a critical parameter of the entire pandemic dynamics. Hence, τ merits a detailed physical understanding. Given the composition of the respiratory droplets, modeling τ is highly non-trivial and is taken up in Sec. V.

It is well documented in the literature that an average human exhales droplets (consisting of water, salt, proteins, and virus/bacteria) in the range of 1 μm-2000 μm. 5, 29, 30 In this section, we offer a detailed exposition of the evaporation dynamics of such droplets as ejected during the course of breathing, talking, sneezing, or coughing.

The small droplets (<2 μm-3 μm) have a very short evaporation timescale. This implies that these droplets evaporate quickly (<1 s) after being ejected. However, the same conclusion does not hold for slightly larger droplets ejected in the form of cloud (>5 μm). These droplets exhibit longer evaporation time, leading to increased chances of transmission of the droplet laden viruses. In particular, when inhaled, these droplets enable quick and effective transport of the virus directly to the lungs airways causing a higher probability of infection. In general, the smaller droplets (<30 μm) have low Stokes number, thereby allowing them to float in ambient air without the propensity to settle down. For larger droplets (>100 μm), the settling timescale is very small (∼0.5 s). In effect, based on the diameter of the exhaled droplets, there are three distinct possibilities:

• Small droplets (<5 μm) evaporate within a fraction of second. • Large droplets (>100 μm) settle within a small time frame (<0.5 s), limiting the radius of infection. • Intermediate droplets (∼30 μm) show the highest probability of infection due to a slightly longer evaporation lifetime and low Stokes number.

In this work, we particularly focus our attention to the modeling of droplets over a large range of diameters from 1 μm to 100 μm. Based on the available literature, we assume that the droplets exhaled during breathing are at an initial temperature of 30 ○ C. 31 The ambient condition, however, vary strongly with geographical and seasonal changes, etc. Hence, in the following, we conduct a parametric study to determine the droplet lifetime across a large variation of temperature and relative humidity conditions. The droplet evaporation physics is complicated by the presence of non-volatile salts (predominantly NaCl) as present in our saliva. 4 We would also look into simultaneous desiccation of the solvent and crystallization of such salts Subsections V A and V B. Once exhaled and encountering ambience, the droplet will evaporate as it undergoes simultaneous heat and mass transfer.

For the modeling purpose, the exhaled droplets are assumed to evaporate in a quiescent environment at a fixed ambient temperature and relative humidity. In reality, during coughing, talking, or sneezing, the droplets are exhaled in a turbulent jet/puff. 16 However, as shown in Eqs. (10) and (11), the puff rapidly decelerates due to entrainment and lack of sustained momentum source, rendering the average VD,P to be less than 1% of the initial velocity. Furthermore, since the Prandtl number, defined as ratio of kinematic viscosity and thermal diffusivity, is approximately unity (Pr = ν/α ≈ 0.71) for air, we can safely assume that the temperature and relative humidity that the droplets in the puff experience are on average very close to that of the ambient. At the initial stages, the puff will indeed be slightly affected by buoyancy, which will influence droplet cooling and evaporation dynamics. Quantifying these effects accurately, merit separate studies, see for e.g., Ref. 32 for buoyant clouds. In a higher dimensional model, these could be incorporated. Nonetheless, the evaporation rate of the droplet is driven by the transport of water vapor from the droplet surface to the ambient far field. Assuming the quasi-steady state condition, the evaporation mass flux can be written asṁ

Here,ṁ 1 is the rate of change of the droplet water mass due to evaporation, Rs is the instantaneous droplet radius, ρv is the density of water vapor, Dv is the binary diffusivity of water vapor in air, and αg is the thermal diffusivity of surrounding air. BM = (Y 1,s − Y 1,∞ )/(1 − Y 1,s ) and B T = Cp, l (Ts − T∞)/h fg are the Spalding mass transfer and heat transfer numbers, respectively. Here, Y 1 is the mass fraction of water vapor, while subscripts s and ∞ denote the location at the droplet surface and at the far field, respectively. The numerical subscripts 1, 2, and 3 will denote water, air, and salt, respectively. Cp, l and h fg are the specific heat and specific latent heat of vaporization of the droplet liquid. For the pure water droplet, the vapor at the droplet surface can be assumed to be at the saturated state. However, as indicated earlier, the exhaled droplets during talking, coughing, or sneezing are not necessarily pure water; rather, they contain plethora of dissolved substances. 5 The existence of these dissolved non-volatile substances, henceforth denoted as solute, significantly affects the evaporation of these droplets by suppressing the vapor pressure at the droplet surface. The modified vapor pressure at the droplet surface for binary solution can be expressed by Raoult's Law, Pvap(Ts, χ 1,s ) = χ 1,s Psat(Ts), where χ 1,s is the mole-fraction of the evaporating solvent (here water) at the droplet surface in the liquid phase 14 and χ 1,s = 1 − χ 3,s . The far field vapor concentration, on the other hand, is related to the relative humidity of the ambient. Considering the effects of Raoult's law and relative humidity, the vapor concentrations at the droplet surface and at the far field can be expressed as 

M 1 and M 2 denote the molecular weights of water and air, respectively. For evaporation, the droplet requires latent heat, which is provided by the droplet's internal energy and surrounding ambient. It has been verified that the thermal gradient in the liquid phase is rather small. Therefore, neglecting the internal thermal gradients, the energy balance is given by

where Ts is instantaneous droplet temperature, m = (4/3)πρ l R 3 s and As = 4πR 2 s are the instantaneous mass and surface area of the droplet, ρ l and e l are the density and specific internal energy of the binary mixture of salt (if present) and water, and kg is the conductivity of gas surrounding the droplet. ∂T ∂r |s is the thermal gradient at the droplet surface and can be approximated as (Ts − T∞)/Rs, which is identical to convective heat transfer for a sphere with a Nusselt number of 2. As such, including aerodynamic effects, the Nusselt number is given by Nu = 2 + 0.6Re 0.5 p Pr 0.3 . The droplet Reynolds number, Rep, was observed to be mostly less than 0.1, and as such, the aerodynamic enhancement of the Nusselt number, i.e., the second term in the right-hand side, is ignored.

Evaporative loss of water leads to an increase in the salt concentration in the droplet with time. As shown before, Pvap(Ts, χ 1,s ) is a function of the salt concentration in the droplet, which thus must be modeled using the species balance equation, as shown in the following equation:

Here, Y 3 is the dissolved solute (salt) mass fraction in the droplet. m 3,out , which represents the rate at which solute (salt) mass leaves the solution due to crystallization, is modeled below. Clearly, Eq. (18) shows that as water leaves the droplet, Y 3 increases. When Y 3 is sufficiently large such that the supersaturation ratio S = Y 3 /Y 3,c exceeds unity, crystallization begins. Here, we use Y 3,c = 0.393 based on the efflorescent concentration of 648 g/l reported for NaCl-water droplets in Ref. 33 . The growth rate of the crystal could be modeled using a simplified rate equation from 34,35

Here, l is the crystal length. Following Ref. 35 , for NaCl, we find the constant Ccr = 1.14 × 10 4 m/s, the activation energy Ea = 58 180 J/mol, and constant gcr = 1. Using this, the rate of change of the crystal mass, which equalsṁ 3,out , is given by 35

We note that while crystallization process could involve complex kinetics of solute, solvent, and ions; a well-studied 35 single-step crystallization kinetics has been used here for tractability. It will be shown that this model is able to predict the experimentally studied droplet lifetime reasonably well. The governing equations [Eqs. (15)- (20) ] manifest that several physical mechanisms are coupled during the evaporation process. For Ts ,0 > T∞, the droplet undergoes rapid cooling from its initial value. The droplet temperature, however, should eventually reach a steady state limit (wet bulb). This limit is such that the droplet surface temperature will be lower than the ambient, implying a positive temperature gradient or heat input. The heat subsequently transferred from the ambient to the droplet surface after attaining the wet-bulb limit is used completely for evaporating the drop without any change in sensible enthalpy. For a droplet with pure water, i.e., no dissolved non-volatile content, the mole-fraction of the solvent at the surface remains constant at 100%, and at the limit of steady state, the droplet evaporation can be written in terms of the well-known D 2 law, 14,18

where

However, for a droplet with the binary solution, the evaporation becomes strongly dependent on the solvent (or solute) molefraction, which reduces (or increases) with evaporative mass loss. The transient analysis, thus, becomes critically important in determining the evolution of the droplet surface temperature and instantaneous droplet size. During evaporation, the mole-fraction of the solute increases and attains a critical super-saturation limit, which triggers precipitation. The precipitation and accompanied crystallization dynamics, essentially, reduce the solute mass dissolved the in liquid phase, leading to a momentary decrease in its mole-fraction. This, in turn, increases the evaporation rate as mandated by Raoult's law, which subsequently increases the solute concentration. These competing mechanisms control evaporation at the latter stages of the droplet lifetime. At a certain point, due to continuous evaporation, the liquid mass completely depletes and evaporation stops. The droplet after complete desiccation consists only of salt crystals, probably encapsulating the viruses and rendering them inactive. If the SARS-CoV-2 virus would remain active within the salt crystal, also known as droplet nucleus or aerosol, Covid-19 could spread by aerosol transmission in addition to that by droplets. In this paper, we focus on infection spread exclusive by respiratory droplets since the role of aerosols is not clear for transmission of Covid-19.

To validate the model, few targeted experiments were conducted to observe isolated levitated droplets evaporating in a fixed ambient condition. Particularly, the droplets with (1% w/w) NaCl solution vaporized to shrink to 30% of its initial diameter during the first stage of evaporation, as shown in Fig. 3 . Hereafter, a plateaulike stage is approached due to increased solute accumulation near the droplet's surface, which inhibits the diameter from shrinking rapidly. However, as shown in Fig. 3 , shrinkage does occur (until Ds/Ds ,0 ≈ 0.2) as the droplet undergoes a sol-gel transformation. The final shape of the precipitate is better observed from the micrographs presented in Fig. 3 . Figure 3 shows the final precipitate morphology for the desiccated droplets. The precipitates display a cuboid shaped crystalline formation, which is consistent with the structure of the NaCl crystal. The size and crystallite structure does show some variation, which could be linked with the initial size of the droplet. Only precipitates from larger droplets could be collected since smaller sized precipitates tend to de-stabilize and fly-off the levitator post-desiccation. While the precipitate from larger sized droplets tend to yield larger and less number of crystals, smaller droplets seem to degenerate into even smaller crystallites. However, this work does not investigate the dynamics of morphological changes of crystallization in levitated droplets. Figure 3 also displays the comparison between results obtained from experiments and modeling. Experiments were performed with both pure water droplets as well as with droplets of 1% salt solutions. Experiments have been described in Sec. II. For the pure water cases shown in the left panel, simulation results follow the experiments rather closely. In the pure water case, classical D 2 law behavior could be observed. For the salt water droplets, a deviation from the D 2 law behavior occurs, and droplet evaporation is slowed. This is governed by Raoult's law where the reduced vapor pressure Pvap(Ts, χ 1,s ) on the droplet surface results from the increasing salt concentration with time. The evaporation rate approaches zero at about Ds/Ds ,0 for 0.3 (experiments) and 0.25 (simulations), respectively. However, the salt concentration attained at this stage exceeds the supersaturation S ≥ 1 required for onset of crystallization. Thus, the salt crystallizes, reducing its concentration and increasing Pvap(Ts, χ 1,s ) such that evaporation and water mass loss can proceed until nearly all the water has evaporated and only a piece of solid crystal as shown in Fig. 3 is left. It can be observed from Fig. 3 that in all cases, the final evaporation time is predicted within 15% of the experimental values. This suggests that despite the model being devoid of complexities (associated with inhomogeneities of temperature and solute mass fraction within the droplet and simple one step reaction to model the crystallization kinetics), it demonstrates reasonably good predictive capability. It is prudent to mention again that although we have done the analysis for the single isolated droplet, in reality, coughing or sneezing involves a whole gamut of droplet sizes in the form of a cloud.

Humans expel respiratory droplets while sneezing, coughing, or talking loudly. Such droplets have a size range from 5 μm to 2000 μm, 36 while the dispersion could depend on the severity of the action. For example, while talking, an average human being will expel ∼600 droplets in a size range of 25 μm-50 μm, but this number goes upto ∼800 in the case of sneezing. For any given act (sneezing, coughing, or talking), the highest number of droplets fall in the range between 25 μm and 50 μm, while the smaller or larger droplets (than the above) are comparatively fewer in number. Nevertheless, the expelled volume of air contains a very small fraction of liquid droplets. To illustrate this point, the droplets are assumed to be in a uniform dispersion. The total volume of air expelled by a human being is estimated to be Vg ∼ 0.0005 m 3 . The total volume of liquid for a given mean droplet size is simply (NVD,s ,0 ), where N is the total number of droplets of size Ds ,0 and VD,s ,0 is the volume of such a droplet. The total volume occupied by droplets of different sizes in a given act of coughing or sneezing is ∑ i NiVD s,0 ,i. Thus, the volume fraction (ϕ) of liquid droplets during a given act is ∑ i NiVD s,0 ,i/Vg . Using size distribution, reported in Ref. 36 , one can show that ϕ for sneezing, coughing, coughing with covered mouth, and talking loudly is 59 ppm, 549 ppm, 361 ppm, and 263

scitation.org/journal/phf ppm, respectively. These numbers indicate that respiratory spray is rather sparse, implying that collective evaporation of droplet clusters may not be significant. This also justifies the modeling based on an isolated, contact free droplet.

Next, we set out to use this model to predict the droplet lifetime characteristics over a wide range of ambient conditions. The droplet evaporation time tevap(Ds ,0 ) is calculated from the analysis presented in Secs. V A and V B. Indeed, the droplet evaporation competes with gravitational settling. 2 The settling time t settle (Ds ,0 ) is calculated by accounting for the decreasing diameter using the equation for the Stokes settling velocity,

The settling time is estimated as that time by which the droplet gets out of the radius from which breathable air is collected-already defined as σH/2 in Sec. III. Mathematically, t settle is obtained by the following equation:

Clearly, for any condition, while tevap monotonically increases with Ds ,0 , t settle monotonically decreases with Ds ,0 . In view of this, it is necessary to estimate the maximum time an exhaled droplet can remain within the collection volume without being evaporated or settled. Such a time can be estimated by defining a characteristic droplet lifetime τ, where τ = min{ tevap | tevap ≥ t settle ∀D s,0 }. (25) τ is essentially the time where the two curves tevap, tevap as a function of Ds ,0 intersect and represents the maximum time a liquid droplet of any size can exist before it is removed either by evaporation or gravity. For a given ambience specified by the ordered pair (T∞, RH∞), Ds ,0 corresponding to τ can be defined as Dcrit. Droplets with Ds ,0 > Dcrit settle due to gravity, while Ds ,0 ≤ Dcrit evaporate, earlier than the droplets with Ds ,0 = Dcrit. While droplets with Ds ,0 ≠ Dcrit can certainly transmit the disease, those with Dcrit establishes the boundaries in terms of the lifetime, cloud diameter, and maximum distance traversed. Dcrit is dependent on ambient conditions, i.e., temperature and relative humidity. The distribution of Dcrit over a wide range of relevant ambient conditions is shown in Fig. 4(a) . Interestingly, at high T∞ and low RH∞, where the evaporation rate is very fast, even a large droplet rapidly shrinks before it can settle. By the same argument at low T∞ and high RH∞, a relatively smaller droplet cannot evaporate quickly; therefore, tevap = t settle is attained for smaller droplet sizes. This explains why we observe large Dcrit at high T∞, low RH∞, and small Dcrit at low T∞ and high RH∞. The evaporation time of the droplet of diameter Dcrit, which has been established as the characteristic lifetime τ of the droplet set, is shown in Fig. 4(b) . Despite the different initial sizes, we find that τ is minimum at high T∞ and low RH∞ conditions, whereas it is maximum at low T∞ and high RH∞. At the same time, longer lifetime, i.e., large τ, allows the droplet cloud to travel a longer distance axially (Xp) and disperse radially (σD). Thus, XP and σD are 4(d) . This also shows how the minimum required "social distance" represented by Xp is not constant but depends on ambient conditions. Combining these results, it can be concluded that the size, lifetime, distance traveled, and the radial dispersion of the longest surviving droplet is not constant and is a strong function of ambient conditions. In particular, low temperature and high RH enhance the droplet lifetime significantly. Relative humidity strongly affects the droplet lifetime compared to temperature. An increase in droplet lifetime also implies that such droplets stay in the ambient for longer periods and hence travel longer distances as reflected by Xp. This implies that such droplets can lead to higher infection propensities. The critical size, droplet lifetime, distance traveled, or size of the droplet cloud for any practical condition are readily obtainable from Figs. 4(a)-4(d), respectively. In Fig. 5 , we look into the evolution of the normalized mass and temperature of the droplets for two cases named as case A and case B, also identified with black dots in Fig. 4 . For case A, (T∞, RH∞) = (8, 55), while for case B, (T∞, RH∞) = (28, 77). The unit of T∞ is ○ C and that of RH is %. These conditions have been specifically chosen to loosely represent the spring weather in North America and South East Asia, respectively, for the early days of the Covid-19 pandemic at both these locations. Figure 5 clearly explains why τA > τB. Indeed, higher RH at case B implies that the temporary hiatus in evaporation due to reduced vapor pressure is reached at a higher water mass load of the droplet than at the case A condition. In both cases, this occurs at about 3 s. However, the higher temperature in B results in faster crystallization kinetics due to the Arrhenius nature of the equation given by Eq. (19) , which causes an eventual faster crystallization rate than at A. Figure 5 clearly shows that although the knee in the solvent depletion profiles are attained at the same time, it is the crystallization and simultaneous desiccation dynamics that governs the eventual difference in lifetime of the droplets at two conditions. It remains to be seen whether this result holds for a detailed crystallization reaction mechanism. The temperature evolution plots shown in the right panel of Fig. 5 also reveals how the droplet initially exhaled at 30 ○ C rapidly cools to the corresponding wet-bulb temperature to subsequently allow heat transfer into the droplet leading to evaporation. However, as the salt concentration reduces due to crystallization, the temperature rises subsequently above the corresponding wet-bulb limits.

With the droplet lifetime available over a wide range of conditions, the corresponding infection rate constant and eigenvalues given by Eqs. (5) and (6) could be evaluated. Just to recapitulate, τ determines the infection rate constant k 1 by Eq. (14) . In turn, k 1 affects the exponents of the time dependent infection equation [Eq. (7)], the growth parameters-eigenvalues λ 1 , λ 2 through Eq. (8) . The contours of λ 1 , λ 2 , and k 1 as a function of T∞ and RH∞ are shown in Figs. 6(a)-6(c), respectively. The direct correspondence between τ and k 1 is immediately apparent upon comparing the respective Figs. 4(b) and 6(c). The infection rate constant is highest at low T∞ and high RH∞ where the droplet evaporation is slowed due to the slow mass loss rate and enhanced crystallization time. On the other hand, faster droplet evaporation leads to small infection rate constant values at high T∞ and low RH∞. The temperature and relative humidity dependency of the eigenvalues λ 1 and λ 2 are shown in Figs. 6(a) and 6(b), respectively. The direct correspondence of Figs. 6(a) and 6(b) with Fig. 6 (c) and Fig. 4 are established through λ 1,2 = −0.54 ± √ 0.22 + k 1 given by Eq. (8) . It should, however, be noted that due to the inherent negative sign of λ 2 , its influence on determining the growth rate of the infected population is rather limited. It is λ 1 that primarily drives the growth of the infected population as apparent from Eq. (7). From Fig. 6 (a), we observe that for a fixed T∞, λ 1 increases with RH∞, while for a fixed RH∞, λ 1 decreases with T∞. Furthermore, we observe that the iso-λ 1 contour lines bend and converge at RH∞ > 75%. This means that for RH∞ > 75%, large λ 1 > 0.4 is expected over a wider range of temperatures between 5 ○ C < T∞ < 20 ○ C. This is a manifestation of the greatly reduced evaporation potential-the difference between water vapor concentration on the droplet surface and in the ambient at high RH∞ conditions. This is further reflected in the dramatic difference in the rate ratio-the ratio of the cumulative number of positive cases on a particular day to the cumulative number of positive cases seven days before: NP/NP ,0 in Fig. 6(d) . This figure has been arrived at by assuming the local population density to be 10 000 km −2 . Furthermore, we calculate tc in Eq. (14) as tc = 3600 × 24/Nexp. Nexp is the number of infecting expiratory events per person per day and is assumed to be 3 based on the coughing frequency of 0-16 in normal subjects. 37 We find that purely based on ambient conditions, implying all other factors have been held constant, the rate ratio can be different by an order of magnitude between (T∞, RH∞) = (5, 75) vs (35, 20) . Practically, such a contrast might be less apparent in real data in which other important factors such as population density, social enforcement, travel patterns, and susceptible supply 38 exert significant influence.

Respiratory flow ejected by human beings consists of a polydisperse collection of droplets. In this paper, we have presented 

ARTICLE scitation.org/journal/phf a model for the early phases of a Covid-19 like pandemic based on the aerodynamics and evaporation characteristics of respiratory droplets. The model and its inter-dependencies on the different physical principles/sub-models are summarized in Fig. 7 . To our knowledge, this is the first model that utilizes the structure of a chemical reaction mechanism to connect the pandemic evolution equations with respiratory droplet lifetime by first principles modeling of the reaction rate constant. However, it must be recognized that the model assumes conditions where transmission occurs solely due to inhalation of infected respiratory droplets alongside many other simplifying assumptions. The evolution of the droplets is characterized by a complex interaction of aerodynamics, evaporation thermodynamics, and crystallization kinetics. As such, after being ejected, smaller droplets attain the wet-bulb temperature corresponding to the local ambience and begin to evaporate. However, due to the presence of dissolved salt, the evaporation stops when the size of the droplet reaches about 20%-30% of the initial diameter, but now, the droplet salt concentration has increased to levels that trigger onset of crystallization. Of course, these processes compete with settling-the process by which larger droplets fall away before they can evaporate. The smaller of the two, complete evaporation time and settling time, thus dictates the droplet lifetime τ. The infection rate constant derived using collision theory of reaction rates is shown to be a function of the respiratory droplet lifetime (τ), where τ is sensitive to ambient conditions. While the infection rate constant in reality is dependent on numerous parameters, the present approach allows us to compute its exclusive dependence on ambient conditions through respiratory droplet modeling. We find that the respiratory droplets exclusively contribute to the infection growth parameters and infection growth rate, which decrease with ambient temperature and increase with relative humidity. As such, the model could be used for providing fundamental insights into the role of respiratory droplets in Covid-19 type viral disease spread. Furthermore, the model could be used, with extreme caution and in cognizance of its limitations, toward estimating the risk potential of infection spread by droplet transmission for specific ambient conditions of interest from purely physics based calculations.

Modes of transmission of virus causing Covid-19: Implications for IPC precaution recommendations: Scientific brief

On air-borne infection: Study II. Droplets and droplet nuclei

Airborne Contagion and Air Hygiene. An Ecological Study of Droplet Infections

Natural Ventilation for Infection Control in Health-Care Settings (World Health Organization

How far droplets can move in indoor environments-Revisiting the wells evaporation-falling curve

Dynamics of infectious disease transmission by inhalable respiratory droplets

Spatial dynamics of airborne infectious diseases

Toward understanding the risk of secondary airborne infection: Emission of respirable pathogens

Capillary flow as the cause of ring stains from dried liquid drops

On the lifetimes of evaporating droplets with related initial and receding contact angles

Retention and infiltration of bacteria on a plant leaf driven by surface water evaporation

Precision stacking of nanoparticle laden sessile droplets to control solute deposit morphology

Insights into vapor-mediated interactions in a nanocolloidal droplet system: Evaporation dynamics and affects on self-assembly topologies on macro-to microscales

Fluid Dynamics and Transport of Droplet and Sprays

Evaporating droplets in turbulence studied with statistically stationary homogeneous direct numerical simulation

Violent expiratory events: On coughing and sneezing

A contribution to the mathematical theory of epidemics

Combustion Physics

A directed relation graph method for mechanism reduction

Turbulent Combustion

Modeling and Simulation of Turbulent Combustion

Effects of acoustic-streaming-induced flow in evaporating nanofluid droplets

Experimental analysis of thermo-physical processes in acoustically levitated heated droplets

Unsteady turbulent round jets and vortex motion

Environmental Fluid Mechanics

Characterizations of particle size distribution of the droplets exhaled by sneeze

On coughing and airborne droplet transmission to humans

The numbers and the sites of origin of the droplets expelled during expiratory activities

Exhaled droplets due to talking and coughing

Validation of the exhaled breath temperature measure: Reference values in healthy subjects

Laboratory simulations show diabatic heating drives cumulus-cloud evolution and entrainment

Drying kinetics of salt solution droplets: Water evaporation rates and crystallization

Evaporation with sodium chloride crystallization in a capillary tube

Salt transport and crystallization in porous limestone: Neutron-xray imaging and poromechanical modeling

The size and the duration of air-carriage of respiratory droplets and droplet-nuclei

Coughing frequency in patients with persistent cough: Assessment using a 24 hour ambulatory recorder

Susceptible supply limits the role of climate in the early SARS-CoV-2 pandemic

The authors thank Professor Chung K. Law 

The data that support the findings of this study are available from the corresponding authors upon reasonable request.