Evaporation of a Liquid Droplet

Two idealized models for the preheat stage of liquid droplets are analyzed theoretically. These models contain the effects of transient heat conduction and evaporation. It is assumed that the droplet surface area decreases linearly with time. This assumption necessitates the solution of moving boundary problems. These models, however, do not consider gas-phase mass transport. In the finite-gradient model, the temperatures of both the droplet and surrounding hot gases vary spatially and temporally. In the zero-gradient model the gas temperature varies spatially and temporally but the droplet temperature varies only temporally, i.e., the droplet temperature is spatially uniform. Numerical examples, which require extensive calculations of confluent hypergeometric functions, are presented for typical values of the droplet latent heat and evaporation rate constant. The temperature profiles given by the finite-gradient and zero-gradient models agree to within 20 percent of each other for all cases examined.


Introduction
The ignition of a droplet of conventional fuel consists essentially of two stages [l] 1 .During the first or preheat stage heat flows from the hot surrounding gas to the droplet, causing the droplet temperature to rise and liquid fuel to evaporate from the droplet surface.In the second stage, ignition occurs in the gaseous mixture of fuel and oxidizer surrounding the droplet.The preheat stage is dominated by transient processes and provides the motivation for the subject of this paper.
Wise and Ablow [2], Parks, et al. [3], and Waldman [4] have previously analyzed the effects of transient droplet heating by neglecting the existence of internal circulation and by neglecting selected terms in the heat conduction equations.In addition they assume that the droplet surface regresses linearly with time.In a previous publication, Bennett and Kayser [5] included the effects of internal circulation but also neglected evaporation and regression of the droplet surface.We now include these effects in the present paper.
The neglect of internal circulation can in many instances lead to gross underestimates of the heat transfer rates within the droplet.Available experimental evidence indicates that in most cases the heat transfer rate inside the droplet is much faster than is possible with heat conduction alone.El Wakil, et al. [6] observed vigorous circulation within the droplet and showed that the droplet temperature is uniform even during the initial transient period.This indicates that the rate at which heat is transferred to the droplet is much slower than the rate of internal mixing.As pointed out by Law [7], assuming the droplet temperature to be uniform in a sense circumvents the difficult task of describing internal circulation but still includes the effects due to its presence; that is, uniform droplet temperature represents the limit of rapid internal circulation.
Two models for the transient preheat stage are presented.In the first model, called the finite-gradient model, the effects of transient heat conduction and evaporation with linear surface regression rate are analyzed.The droplet and gas temperatures vary spatially and temporally and the effects of internal circulation are neglected.In the second model, called the zero-gradient model, we include the effects of transient heat flow, evaporation, and internal circulation.The latter is included by assuming the droplet temperature to be temporally varying but spatially uniform as discussed above.The gas temperature is again both time and space dependent.
The results of the finite-gradient and zero-gradient models are applied to the preheat stage of droplets in hot gases.The temperature profiles are used to predict the time duration between introduction of the droplet into a hot atmosphere and its reaching a given temperature.In addition we compare the analogous times predicted by a particularly simple and accurate approximation to the exact expression derived for the zero-gradient model.
In the next section we describe the physical assumptions for the preheat stage common to both models.Section 3 contains the mathematical description of these assumptions.We derive in section 4 the temperature profiles of the droplet and gas predicted by the finite-gradient model while in section 5 we give the corresponding expressions for the zerogradient model.We present in the last section numerical examples for both models.The temperature profiles obtained are very similar for both models and are within 20 percent of each other in all cases examined.

Preheat Stage
Williams has examined the problem of the assumptions underlying droplet vaporization and combustion theories in some detail [8].In models for preheating, many researchers consider the droplet to be a sphere with temperature independent thermal properties and neglect diffusion effects [1].Others include thermal diffusion effects but later in their theoretical analyses neglect some terms in the diffusion equa-tions [3].In addition, most theoretical investigators neglect the influence of gravity, forced convection, and radiant heat transfer.Combustion chambers and experimental techniques such as freely falling droplets exist to eliminate gravitational effects.Proper selection of some experimental parameters such as ambient pressure and initial droplet size may minimize the other effects not included in the model for preheating.Several researchers have shown that except for the very heavy fuel oils, radiant heat transfer from the hot gas or from adjacent droplets is negligible [7,9].
The models examined in this paper contain only thermal diffusion effects.We assume that a spherical liquid droplet of initial radius r = r^o is inserted at time t = 0 into a hot atmosphere (e.g., air).The initial temperature of the droplet at time t = 0 is T^o-The hot gas is unbounded and initially at a constant temperature TQQ.For t > 0 the temperature at infinity is kept at T g() .
The inclusion of evaporation demands that the droplet radius change with time.In the present work a linear surface regression rate is assumed:  [10].It must be stressed at this point that the solutions which we will obtain will be applicable only in those time regimes for which eq (1) is an accurate representation of the actual surface regression rate.We denote the density, specific heat, and the thermal conductivity of the droplet respectively by da, Cd, and Kd and of the gas respectively by dg, C g and K g .We assume that these thermal properties remain spatially and temporally constant during preheat.In addition, the ambient pressure of the host gas does not change.Using these and the above assumptions we seek to compute the dependence of the temperature T(r, i) upon the radial distance r and time t.In particular, we shall compute the time required for a droplet to reach a given temperature.

Theoretical Analysis
We first consider the finite-gradient model.The distribution of temperature within the system is governed by the Fourier heat conduction equation for the case in which the isothermal surfaces are concentric spheres.Let T d (r, t), T g (r, t), a d = (KJddCd), and a| = (K g /dgCg) denote the temperatures and thermal diffusivities of the droplet and surrounding gas respectively.The governing equation for t ^ 0 and 0 ^ r < R (inside the droplet) reads -^-= a A-^-+ 1-w dt \ dr 2 r dr J and for t ^ 0 and r > R (outside the droplet) The quantity R(t) specifies the motion of the interface with time and is given in eq. ( 1).The specification of the boundary conditions completes the statement of the problem.The temperature T(r, t) has the form t T d (r, t) = T d0 for 0 < r < R and t < 0 (4) and T g (r, t) = T g0 for r > R and t < 0 (5) and it becomes for t > 0 The temperature is finite everywhere, namely T d (0, t) is finite and

r->oo
The symmetry condition that no heat flux exists at the center of the droplet is dT d (0, t) I = 0 (7) or dT d (0, t)/dr is of course taken to mean the partial derivative of Td(r, t) which respect to r evaluated at r = 0.The problem statement for the finite-gradient model is completed by two conditions at the interface r = R which match the temperature and heat flux

T d (R~, t) = T g (R + , t) (8) v dT d (R~, t)_ tr dT g (R + , t) dR(t)
where R ± = lim (R ± Q and L^ is the latent heat of vapori-£."*°zation per unit mass of the condensed material.From eq (1) we note that dR/dt goes off to infinity near the end of the droplet lifetime so that the difference between the thermal conductivities times the temperature gradients must become unbounded.More will be said about this point in subsequent sections.
We find it convenient at this point to introduce dimensionless quantities which in this paper will be denoted by Greek letters.They are defined as follows: (1) The time r is chosen to move the singular point of complete droplet evaporation to infinity, r = 2{ajk e ) In ( r do/R(t))-As the evaporation rate constant k e goes to zero r reduces to r 0 = ddt/rdo- (2) The distance coordinate Tj is chosen to remove the time dependence of the position of the boundary, 17 = r/R(t).Distance is thus measured in instantaneous droplet radii.
Qdi 7 ), T ) is tne fraction of the maximum possible temperature rise.In addition, we define a = a d /a g , /3 = K g /K d , e =
We now consider the zero-gradient model.As was mentioned in section 1, the internal circulation inside the droplet where a is the dimensionless Laplace transform variable, may be sufficient in many cases to maintain a spatially j^e Bromwich integral uniform temperature within the droplet.When such conditions prevail, the Fourier heat conduction equation for the droplet, eq (2), and the boundary condition (9) are replaced by a single heat transfer equation.The temperatures in the zero-gradient model are denoted by T dz (i) and T gz (r, t) for the droplet and gas respectively.For the case in which the droplet temperature is spatially uniform, T d (r, t) = T dz (t) = expresses the temperature in terms of the Laplace transform.T g (R , t) -T gz (R , t) and the heat flux into the droplet, The quantity y is chosen sufficiently large so that the integral Jo exp (-orr)|e(i7, r)\dT
The operator in square brackets is the confluent hypergeometric differential operator and eqs (34) and ( 35) are confluent hypergeometric differential equations.We obtain the general solutions to eqs (28) and (29) as (31) f°r tne not surrounding gas.The quantities M(a, b, z) and U(a, b, z) are two independent solutions to the confluent hypergeometric differential equation and are called confluent hypergeometric functions.Their most useful general proper-(32) ti es ma y be found in Chapter 13 of Ref. [11].
The constants A, B, C, and/) must be determined from the conditions (30)-(33).The condition that no heat flux exists at the center of the droplet requires that B = 0 since c7(cr/€, 3/ 2, €T7 2 /4) is not well-behaved at t) = 0. Similarly, the condition that the reduced gas temperature goes to unity at large distances requires C -0 since M(cr/e, 3/2, a 2 €T7 2 /4) is not well-behaved at infinity.The two remaining conditions at the interfacial boundary r) = 1 constitute two inhomogeneous simultaneous equations for the remaining constants A andD.To solve the subsidiary differential eqs (28) and (29), we Solving these, we obtain the Laplace transformed solutions make a change of variable.Letting £ = €T7 2 /4 gives for the reduced temperatures: The differential relations of M and U have been used in deriving these equations.These expressions have been derived in a form suitable for going to the limit of zero evaporation rate constant, € = 0.The correct limiting cases may be found in Ref. [6] and the details of the limiting process in appendix A.

8-*» ATTl J y-i8
To evaluate the Bromwich integrals it is first necessary to investigate the nature of the integrand exp (OT)<I>(T7, cr).NO distinction will be made in this discussion between the droplet and the gas as the integrands in both cases exhibit the same general properties.First of all, exp ((TT)<$>(r), cr) is a single-valued, analytic function of cr except at its singularities.It has a simple second order pole at cr = 0 and an infinity of real, negative simple first order poles at the zeroes of//(cr/€) given in eq (40).To show that H(cr/e) possesses no complex, imaginary or positive zeros is straight-fore ward.General methods for locating such zeros may be found, for example, in Carslaw and Jaeger [12] and Ince [13].
We now evaluate the Bromwich integral.We complete the contour of eq (41) by a large semicircular arc in the left-hand complex plane and take y to be an arbitrarily small positive number.The final reduced temperature distributions are then obtained by applying the Cauchy Residue Theorem and letting 8 and the radius of the large semicircular arc tend to infinity.The integral along the large semicircular arc vanishes and we find We pause to make several comments about the expressions obtained for 0<*(i7, r) and O^TJ, T).Both contain a term proportional to the reduced time r which is a consequence of the second order pole at cr = 0.As mentioned after eq (9), the choice of linear surface regression rate implies that the difference between the thermal conductivities times the temperature gradients evaluated at the interface goes to infinity at large times.By large times of course, we mean times close to the end of the droplet's lifetime.These are times for which the linear surface regression rate of eq (1) is no longer an accurate representation of the true physical situation.For \, e, and r of the orders of magnitude of interest here, the linear term in r may be considered an artifact of the model and makes only a very small contribution to the reduced temperatures.We remark also that eqs (42) and ( 43) are not suitable for going directly to the limit € = 0. Further discussion of this point may be found in appendix A.
In section 6 are presented plots of 0<j(l, r) versus r for typical values of the dimensionless quantities \ and e.In the next section we derive expressions for the reduced temperatures for the zero-gradient model.

Zero-Gradient Model
The analysis for the zero-gradient model is mathematically very similar to that for the finite-gradient model and will not be quite as detailed.Equations ( 20)-( 25) are the working equations for the zero-gradient model.The Laplace transformed system of equations is for the gas where I(a/e) = 2I/(o/e, 3/2, a 2 e/4) (52) + 3/3a 2 U(l + o/e, 5/2 , a 2 e/4) These expressions have been derived in a form suitable for going to the limit € = 0, which represents the nonevaporative problem solved in Ref. [6].That they in fact reduce to the correct limit is shown in appendix A.
All that remains to be done is to evaluate the Bromwich integral of eq (8) with <£>(7}, a) = O^TJ, cr) The same general comments made about the corresponding integrands in the finite-gradient model are applicable here also.We complete the contour by a large semicircular arc in the left-hand complex plane, apply the Cauchy Residue Theorem and let 8 and the radius of the semicircular arc tend to infinity.Following this sequence of steps were find e d ,(r) = e«(i + , T)
It is possible to obtain some very simple and accurate approximations to eqs (54) and (55).The thermal diffusivity of the droplet is several orders of magnitude smaller than that of the surrounding gas.This suggests that before applying the Laplace inversion formula we expand the subsidiary solution <& gz (17, or) into a power series in a = a (Jag.We obtain only the simplest approximation by letting a go to zero in eq (51): The inverse Laplace transform of this remarkably simple expression yields the following approximations for the reduced temperatures in the zero-gradient model.
for the droplet, and for the surrounding gas.Note that the term linear in r has disappeared to the order of approximation used.However, if *&gz(yi °") is expanded out to first order in a, then a relatively simple infinite series and the linear term in r show up in the approximation of Q gz (ri, T).Unfortunately, the expressions obtained by letting a go to zero in the finite-gradient model are much the same as the exact solutions.
In the next section, we present plots of ©^(T) versus r for typical values of A and e.In addition, the times predicted by ©d0?> T )> Qdzi 7 ), T ) and 0do( T ) for the droplet to reach a given temperature are compared and discussed.

Numerical Examples and Conclusions
In this section, we give some illustrative numerical examples for the predictions made by finite-gradient and zerogradient models.The input data for the calculations are the four dimensionless quantities a, /3, X and €.We take a = 0.057, p = 0.186 and X = 0.01, 0.055, 0.10 and € = 0.10, 0.55, 1.0, independently, so that there are nine cases.These values correspond roughly to those of medium weight fuel oil in hot air at atmospheric pressure.Typical values of the physical parameters K, d, and C which yield the chosen values of a and )3 are given in table   A; e (cm 2 /s) 6.5 X 10~5 3.575 X 10~4 6.5 X 10~4 of the finite-gradient and zero-gradient models for oil droplets in air for X equal to 0.01, 0.055, and 0.10, respectively.0 d (l~, r) and O^(r), the reduced surface temperatures, are computed from eqs (42) and (54).It should be noted that due to the presence of the linear term in r, neither 0 d (l~, r) nor ©^(r) reach steady state during the droplet lifetime.Even if that term were absent, 0^(17, T) would never become uniform throughout the droplet.We must stress at this point that the apparent large difference between the temperature profiles for fixed A. but different € is deceiving.Although shown on the same graph, the curves are plotted versus er/2 which is a function of e, i.e., each curve has its own distinct time scale.The results have been shown in this fashion for clarity only.If they were all plotted on the same scale, real time t for example, the entire group of curves would lie very close to one another.As can be seen from the figures, 0<j(l -, r) and QaJj) are quite similar.Their difference is greatest when €r/2 is small.This is to be expected.The thermal gradients inside the droplet are largest near the beginning of its lifetime and for such times the approximations of the zero-gradient model are more suspect.Although not plotted, the simple approximation Odo(^") to &dz{ T )-> obtained in the last section, is an extremely good one as long as er/2 is neither too large nor too small.The difference between Odzi r ) and ©do( T ) is less than about 15 percent in the worst case (X = 0.1 and € = 1.0) and much smaller in the other cases for the r's we have considered.Due to its simple analytic form, accuracy, and ease of mathematical computation, we suspect that ©do( T ) is potentially a very useful expression for analyses of the initial preheat period of droplets in hot atmospheres.
A quantity of interest in the design of combustion chambers is the time required for fuel droplet to reach a given temperature T s .Let us denote this time by t e and consider a fuel droplet with the same physical parameters given in table 1 and an initial radius r d o = 50 /xm.For illustrative purposes, let us assume that T SL = 550 K, the initial droplet temperature T d Q = 300 K, and the initial gas temperature T g0 = 1966.67K.The droplet reading T s corresponds to 0 S = 0.15.Using the data displayed graphically as figs, la, lb, and lc, we compute the times required for 0 to increase from zero to 0.15.The results are shown in figure 2 as a plot of t e versus e.It can be seen that the times, t e , predicted by the finite-gradient and zero-gradient models are within about 20 percent of each other.For comparison, the £ e 's predicted by 0do( T ) are also presented and are seen to be within 10-20 percent of those predicted by 0^^) as long as e < 0.55. Figure 2 also provides an extremely powerful cross check on the analytic and numerical results of this paper and those of Ref. [9].By solving a totally different problem, Bennett and Kayser [9] obtained the times t e for the case of e = 0.These are the points on the f-axis of figure 2. By extrapolating the results obtained here for nonzero € to e equal to zero, it can be seen that the results are in agreement and this indicates that this work and that of Ref. [9] are internally consistent.Figure 3 contains a plot of the volume fraction of the droplet that has evaporated away by the time t = t e versus €.For e = 0.1, an evaporation rate constant fairly characteristic of the initial preheat stage, the volume fraction evaporated is less than 5 percent.If € is fairly small, this perhaps justifies neglecting evaporation and surface regression rate in analyzing the initial droplet preheat stage.A fairly decent approximate procedure for predicting the times t e and the droplet size at t = t e is to neglect evaporation, compute the time required for the droplet to reach a given temperature, and then to use that time and a relation like eq (1) to calculate the volume fraction evaporated.

1
Figures in brackets indicatethe literature references at the end of this paper.

FIGURE 2 .
FIGURE 2. The times t e required for a droplet to reach obtained from FIGURE 1.The reduced temperature O d (l~, T) (solid lines) and O dz (r) 0«/i~, r) (solid lines and closed circles), &dz( T ) (dashed lines and crosses) (dashed lines) plotted as a function ofer/2 with a = 0.057, /? = 0.186 and e plotted as a function of efor various values of X. = 0.10, 0.55, and 1.0 for A. equal to (a) 0.01, (b) 0.055, and (c) 0.10.The physical properties used are those of table 1 and the initial droplet radius, rd0t is 50 \un.The initial droplet radius, r^o, is 50 fim.The temperature Ts corresponds to Q(te) = 0.15.
1.In table 2, we give the values of L d corresponding to A. = 0.01, 0.055, and 0.10 and k e corresponding to e = 0.10, 0.55, and 1.0.The temperature difference T g0 -T d o is taken to be 1666.67K (corresponding to T d0 = 300 K and T g0 = 1966.67K) and a d 2 and C d are assigned the values given them in table 1. Figures la, lb, and lc compare the numerical predictions

TABLE 1 .
Typical values of the thermal conductivity K, density d, specific heat C at constant volume, and thermal diffusivity a 2 which yield a = ad/a g =

TABLE 2 .
The values of the droplet latent heat of vaporization per gram, L d , and the evaporation rate constant, k e , corresponding to the values of the dimensionless parameters A. = Lj{Ca{T go ~ T do )} and e = k e /a a 2 used.