TS10052 Revised-Paper
TS10052 Revised-Paper
TS10052 Revised-Paper
net/publication/236270772
CITATIONS READS
30 785
2 authors:
Some of the authors of this publication are also working on these related projects:
All content following this page was uploaded by Chaouki Habchi on 30 May 2014.
Abstract
proposed. Compared to previously published models, it has two new features. First, an
expression of the Stephan velocity is proposed which ensures gas mass conservation. In
addition, the evaporation rate of each species is obtained by the integration of the exact
equation of species mass fraction. Second, the heat flux due to species diffusion is taken into
account in addition to the classical conduction heat flux between the gas and the liquid
above two features is presented for high and low pressure conditions, for which a real and a
perfect fluid equation of state (EOS) has been used, respectively. Free convection is also taken
into account using the Grashof number in the Kulmala-Vesala correlations [1] for the
Sherwood and Nusselt numbers. The model is compared with very accurate experimental data
which were recently obtained by Chauveau et al. [2] at atmospheric pressure and temperature
ranges of 473-973K for n-heptane and 548-623K for n-decane droplets of 400 µm initial size.
A very good agreement with the experimental data including micro-gravity conditions has
been obtained. Indeed, the results have confirmed that the free convection process plays a
significant role in the evaporation rate of liquid droplets under earth gravity and quiescent
conditions. This shows the relevance of the new features of the model. The numerical results
1
Corresponding author: e-mail Chawki.HABCHI@ifp.fr, Tel: +33147526151, Fax : +33147527068
have also shown that real fluid EOS is not necessary at atmospheric pressure for the
temperature range given above. In addition, the numerical results of the new model are also
compared with the experimental data of Birouk [3] for two-component droplets of n-heptane
and n-decane with different compositions of the liquid mixture. Finally, the non-ideality of the
mixture is shown to become significant at high ambient pressures and especially at low
complex diffusion.
1 Introduction
during the past decades. Indeed, this phenomenon is of interest to many fields of engineering
fuel sprays in automotive engines or gas turbines is a prominent example. The influence of the
volatility between the numerous components of practical fuels on pollutant emissions still
received a considerable attention in the past, few experimental accurate data are available for
the validation of droplet evaporation models even at atmospheric pressure. Most isolated
droplet evaporation experiments have been conducted with the drop suspended on a support
fiber to avoid the experimental difficulties for free-falling droplets [4-7]. The support fiber
which has a relatively large diameter (more or less 150µm) increases the drop evaporation rate
due to the heat transfer through the fiber to the drop even at low ambient temperatures. Yang
& Wong [8] have studied the effects of heat conduction into the droplets through the fiber and
the liquid-phase absorption of the radiation from the furnace wall on the evaporation rate of
droplets at micro-gravity condition. They found that without considering these effects, there is
a large discrepancy between their theoretical results and the experimental data of Nomura et
al. [4]. More recently, Chauveau et al. [2] studied the effects of fiber diameter and gas
diameter of 14µm is used to reduce the effect of heat conduction on the droplets evaporation
rate. To consider the effect of free convection, Chauveau et al. [2] applied both normal and
Droplets vaporization models from the literature such as [9, 10] have compared their
numerical results mainly on the experiments with large sizes of fibers [4-7]. Consequently,
these models have overestimated usually the vaporization rate in order to be close to the
measurements. Moreover, most previous models neglect the heat flux between the gas and the
droplet due to the enthalpy transport by species diffusion during the evaporation process. In
addition, the mass conservation is not rigorously respected when Fick's law is used to model
the diffusion velocity in the mass conservation equation of the gas mixture. Thus, a new
model that could have a better estimation of droplets vaporization is needed. Furthermore, a
real gas equation of state (EOS) have to be used in order to take into account non-ideality
effects for pressure and/or temperature supercritical conditions [11, 12]. In addition, the case
of droplets vaporization under quiescent conditions (i.e. without forced convective velocity)
will be considered. Indeed, the effect of natural convection on the evaporation rate may
natural convection, normal and micro-gravity effects have been performed by some authors
[4, 13-16]. The natural convection effect can be presented in term of the Grashof number. In
fact, the Grashof number plays the same role in natural convection that the Reynolds number
plays in forced convection [17]. To consider the effect of natural convection, Ranz and
Marshall [18] introduced the Grashof number in the definition of Sherwood and Nusselt
numbers. In the present study, the most recent correlation of Kulmala & Vesala [1] is used for
the evaluation of Sherwood and Nusselt numbers with consideration of natural convection
effect.
evaporation model is proposed which takes into account a modified Stephan velocity and
enthalpy diffusion of gaseous species to compute the evaporation rate. In Section 3, the
numerical results of a single-component isolated drop are compared with the experimental
data of Chauveau et al. [2] and Yuen and Chen [19]. Results for the vaporization of droplets
of droplets diameter. The effects of high pressure and fluid ideality on vaporization are also
investigated in Section 3.2 using the real and ideal gas EOS. The experimental data of
Nomura et al. [4] at high pressure conditions have not been used in this paper because of the
uncertainty of the experiments mentioned above. To our knowledge, few experimental data
exist for the evaporation of multi-component droplets. In section 3.3, the numerical results are
compared with the experimental data of Birouk [3] in the case of a two-component droplets
vaporization.
The droplet evaporation models are generally based on the conservation equations of
mass, momentum and energy. In the basic conservation equation of species ( see Equation
(2.1)), the total convective velocity includes a first component u, usually called Stefan
velocity and a second component Vi, called the diffusion velocity of species i. This last may
result from the pressure gradients (i.e. pressure diffusion), the temperature gradients (i.e. Soret
effect), the external force fields (i.e. forced diffusion) and the concentration gradients (i.e.
ordinary diffusion). The mass conservation equation of gas mixture is conserved by diffusive
mass flux which is commonly approximated by Fick's law. In this case, the effects of pressure
diffusion and Soret effects are usually assumed to be small [20-22], although, sometimes the
effects of Soret may be significant [11, 23]. A general expression for the determination of the
diffusion velocity is proposed by Williams [24] which is very complex and costly, especially
for multi-species gas. A simplified expression for the diffusion velocity in multi-species gas
has been already presented by Hirschfelder et al. [25] which neglects the pressure gradients
and volume forces. Using Hirschfelder's law reduces the complexity of the exact resolution of
Vi; but, it does not recover the mass conservation equation of the gas mixture. To overcome
this problem, an additional term is suggested by some authors [20, 26] to ensure gas mass
conservation equation. Similarly, a modified Stephan velocity which ensures the mass
conservation of gas mixture is derived below for the droplets vaporization modeling
considering Hirschfelder's law. The model uses spherical droplet hypothesis with no
interaction between droplets. The effects of radiation, Soret and Dufour are neglected and
there is no chemical reaction in the gaseous environment. The one-third rule is used for the
properties of the gaseous mixture in the film region around the droplet [27]. Gas phase quasi-
steadiness and isobaric assumptions are applied to the model. The discrete fluid method is
used for multi-component droplet which resolves the continuity equation for all droplet
components. The droplet temperature is assumed to be uniform during all the process of
evaporation (i.e. the infinite thermal conductivity assumption is invoked). The last hypothesis
gives an average temperature in the droplet which is different from the real surface
temperature. Some methods have been presented to evaluate the temperature profile in
droplets [22, 28]. Nevertheless, temperature uniformity is a reasonable assumption, at least for
small droplets for which the heating period is short. Liquid species mass fraction in liquid
phase is time dependent; but there is no species diffusion in the liquid phase. The components
mass fraction in the liquid phase assumed to be constant in all parts of the droplet (i.e. there is
no liquid component mass fraction profile inside the droplet). However, there are some
valuable works which consider the species diffusion in the liquid phase [29, 30] and
consequently the profile of liquid components mass fraction inside the liquid phase.
The resulting governing equations for the model consist of the two phase flow equations
for the gas and liquid phases along with the thermodynamic equilibrium condition at the
liquid-gas interface.
∂ρ g Yi
+ ∇. ( ρ g Yi ( u + Vi ) ) = 0 (2.1)
∂t
where u and g are the velocity and density of the gas mixture and Yi and Vi are the mass
fraction and diffusion velocity of species i respectively. Summation of Equation (2.1) over all
∂ρ g
+ ∇. ( ρ g u ) = 0 (2.2)
∂t
The internal energy conservation equation for the gas phase is:
∂ρ g I
+ ∇. ( ρ g uI ) = −∇.J + σ : ∇u − P∇.u (2.3)
∂t
where σ : ∇u represents the heat production due to the shear stress, σ , and P∇.u is the
compressibility effect. These two terms are neglected in the present study. J represents the
heat flux vector. It includes the contributions due to heat conduction and enthalpy diffusion
J = −λg ∇T + ρ g hi (YV
i i) (2.4)
i
where g is the thermal conduction coefficient of the gas mixture in the film region (Figure 1),
temperature of the droplet (Ts=Td). The energy conservation equation for the two phase
system consisting of the droplet and the surrounding gas mixture gives the change of liquid
dTd
md C pf = Ql (2.5)
dt
where md is the droplet mass, Cpf is the specific heat of liquid at constant pressure, Td is the
temperature of the droplet and Ql is the heat penetrating into the liquid phase. According to
Figure 1, the energy balance at the liquid surface could be written as:
Ql = Qg − mLv , g (2.6)
where Qg is the heat flux from the gas to the liquid (Qg = J.A), and A is the droplet surface
area. m and Lv , g are the total evaporation mass flow rate and latent heat of the evaporated
species mixture, respectively. Since the real gas EOS is applied to gas mixture, Lv , g is a
function of (P,T). However, Lv , g is only a function of temperature in case of ideal gas EOS
assumption.
The evaporation of droplet consists of resolving two phase flow equations as well as
equilibrium is based on the assumption that at the liquid-gas interface, the chemical potential
for liquid and gas phases are equal for each species i. At low pressure conditions, ideal gas
hypothesis with the ideal gas EOS seems to be a good assumption. In this condition, the ideal
solution could be said to follow the Raoult's law [31]. However, at high pressure conditions
and particularly in supercritical conditions, the gaseous species do not follow the ideal gas
behavior and the effects on some parameters like transport and physical properties variation,
gas solubility increase into the liquid phase and thermodynamics’ non-ideality should be
In this study, both ideal and real gases EOS are evaluated. They are used especially to
determine the mass fraction of gas species at droplets surface as well as latent heat of
Based on the previous works, some authors [9, 10] have considered the gaseous
boundary layer around the droplet to evaluate heat and mass fluxes. These models give the
instantaneous droplet vaporization rate from the integration of the quasi-steady species
balance around the droplet. By definition [26], the diffusion velocity in Equation (2.1) has to
i i =0
YV (2.7)
i =1
N
Yi = 1 (2.8)
i =1
Considering Equations (2.1) and (2.8), there are N unknowns (Yi) with N+1 equations. Thus,
the system is over determined. However, if all species of Equation (2.1) are added and the
identity of Equation (2.7) is used, the mass conservation equation (Equation (2.2)) is
recovered:
∂ρ g N
+ ∇. ( ρ g u ) = −∇. ρ g i i =0
YV (2.9)
∂t i =1
Solving the above equation with the general diffusion velocity is too complex and most codes
Vi X i = − Dig ∇ ( X i ) (2.10)
where X is the mole fraction of species i and Dig the diffusion coefficient of species i into the
gas mixture. In the present study, the above Hirschfelder's law will be used in the form of
mass fraction to be compatible with the model. Using the following relation between the mole
W
X i = Yi (2.11)
Wi
one obtains the Hirschfelder's law in the form of mass fraction of species i as:
N
Xk
ViYi = − Dig 1 − X i + Yi ∇ ( Yi ) (2.12)
k =1 Yk
k ≠i
Using the Hirschfelder's law, the summation of Equation (2.1) over the N gaseous species
∂ρ g N N
Xk
+ ∇. ( ρ g u ) = ∇ ρ g Dig 1 − X i + Yi ∇ ( Yi ) (2.13)
∂t i =1 k =1 Yk
k ≠i
where the right hand side of the above equation is not zero, unless the diffusion coefficients
Dig are assumed to have the same value for all the species (this is the case of the KIVA code
[34] for instance). Thus, the global mass conservation may not be ensured. One method to
implement Hirschfelder's law with different diffusion coefficients while maintaining global
mass conservation is to add a correction velocity uc [26] to the absolute velocity in the species
Equation (2.1):
∂ρ g Yi
( ) Xk
N
+ ∇. ρ g Yi ( u + u c ) = ∇ ρ g Dig 1 − X i + Yi ∇ (Yi ) (2.14)
∂t k =1 Yk
k ≠i
If all species equations are summed, the mass conservation equation must be recovered:
∂ρ g N N
Xk
+ ∇. ( ρ g u ) = ∇ ρ g Dig 1 − X i + Yi ∇ (Yi ) − ρ g u c = 0 (2.15)
∂t i =1 k =1 Yk
k ≠i
Thus,
N N
Xk
uc = Dig 1 − X i + Yi ∇ (Yi ) (2.16)
i =1 k =1 Yk
k ≠i
Assuming quasi-steady evaporation and integrating Equation (2.15) for liquid species leads to
Nl
m= mi = 4π rd2 ρ g u (2.17)
i =1
where rd is the droplet radius. Integrating Equation (2.14) for liquid species between the
∂Yi
Yi = Yi ∞ , u∞ = u∞c = 0 , =0 (2.18)
∂r r =∞
Dig N
X ks ∂Yi
mi = 4π rd2 ρ g Yi s u + u c − 1 − X is + Yi s (2.19)
Yi s k =1 Yks ∂r
k ≠i
r = rs
mi = mYi s + 2π rd ρ g ×
N
X ks Nl N
X as (2.20)
i( i i ) Shk (Yks − Yk∞ )
∞
Dig 1 − X is + Yi s s
Sh Y s
− Y − Yi
s
Dkg 1 − X s
k + Yk
s
s
k =1 Yk k =1 a =1 Ya
k ≠i a≠k
In the above equation, Nl is the total number of liquid components, Yis and Yi are the mass
fraction of component i at the droplet surface and at infinity respectively. The dimensionless
Sherwood number is defined as the mass fraction gradient at the droplet surface over the
To consider the effect of natural convection due to volume forces such as gravity, the
( ( ))
1/ 2
Shi = 2.0009 + 0.514 max Rei , max ( Gri , 0 )
1/2
Sci1/3 (2.22)
with
where, Vrel is the relative velocity between the gas and droplet, ig, ig and Dig are the dynamic
viscosity, partial density and binary diffusion coefficient of species i in the gas mixture of the
film region, g0 is the gravity acceleration; T and ν i∞ are the gas mixture temperature,
kinematic viscosity of species i at infinity; Td is the droplet temperature and Gri, Rei and Sci
are the Grashof, Reynolds and Schmidt numbers of species i respectively. Hence, the
Sherwood number has different value for different components. The physical parameters in
the film region (with index g) are evaluated at the reference temperature [27]:
Tref = Td + Ar (T ∞ − Td ) (2.24)
where Ar is the averaging parameter. For the one-third rule: Ar=1/3. Equation (2.22) is the
Equation (2.20) over all liquid components, gives the total mass evaporation rate:
Nl Nl Nl
N
X ks N
X as
2π rd ρ g Dig 1 − X is + Yi s Shi (Yi s − Yi ∞ ) − Yi s . Dkg 1 − X ks + Yks Shk (Yks − Yk∞ )
i =1 k =1 Yks i =1 k =1 a =1 Yas
k ≠i a ≠k (2.25)
m= Nl
1− Yi s
i =1
Consequently, from the above Equation (2.25) and Equation (2.17), the Stephan velocity is
given as:
Nl Nl Nl
N
X ks N
X as
Dig 1 − X is + Yi s Shi (Yi s − Yi ∞ ) − Yi s . Dkg 1 − X ks + Yks Shk (Yks − Yk∞ )
i =1 k =1 Yks i =1 k =1 a =1 Yas
u=
k ≠i a≠k
(2.26)
Nl
2rd 1 − Yi s
i =1
The total mass evaporation rate could also be obtained from the combination of internal
energy conservation equation of gas mixture (Equation (2.3)) with Equation (2.17). This
During the evaporation process of a drop, the internal energy of the surrounding gas
changes simultaneously to its composition. The heat flux in the gas phase comprises three
contributions. The first one is the thermal conduction flux usually modeled by a Fourier law;
the second contribution is associated with compositional changes resulting from species
diffusion; and the third contribution is the Dufour effect which is related to concentration
gradients [35, 36]. Although the Dufour effect [11, 37] and enthalpy diffusion [21] have been
shown to be significant, most of the droplet vaporization models in the literature [9, 10]
assumed the heat flux only due to the first contribution thermal conduction. Droplet
evaporation models which do not have an appropriate estimation for the heat flux, can lead to
rate. It is shown that Navier-Stokes solvers which incorporate enthalpy diffusion can provide
much more accurate results than without considering this term [20]. Similarly, the Dufour
effect has been neglected in this study. In this case, the heat flux from the gas to the drop is
then assumed to be the sum of the two first contributions where the (YiVi) term in Equation
(2.4) is given by Equation (2.12). Hence, the Equation (2.4) could be rewritten as follows:
Nl N
X ks
Qg = Qconduction − Qenthalpy = 2π rd λg Nu (T ∞ − Td ) − ρ g Dig 1 − X is + Yi s hi Shi (Yi s − Yi ∞ ) (2.27)
i =1 k =1 Yks
k ≠i
Similarly to the Sherwood number, the following correlation is used for the Nusselt
number[1]:
( ( ))
1/ 2
Nu = 2.0009 + 0.514 max Re, max ( Gr , 0 )
1/ 2
Pr1/ 3 (2.29)
C pg µ g
Pr = (2.30)
λg
where, Cpg is the specific heat of the gas mixture in the film region at constant pressure.
An iterative method is needed for the implicit procedure of calculating the mass flow
rate, droplet temperature and composition of the liquid and gas mixtures at the liquid-gas
interface and also in the liquid phase. In this work, a Newton iterative method has been used.
(Equation (5.5) for real gas model and Equation (5.16) for ideal gas model) and with the
initial temperature and composition of the mixture, one can determine the species mass
fraction in the liquid and gas phases. The physical and transport properties of the liquid and
the gas can be obtained from the estimation techniques and mixing rule as recommended by
Reid et al. [31] (see appendix A). From Equations (2.25) and (2.27) one can obtain the total
evaporation mass flow rate and the heat flux from the gas phase to the droplet. The energy
balance at the liquid-gas interface (Equation (2.6)) gives the new droplet temperature. The
variation of droplet radius is evaluated using the droplet mass mdn . The new droplet radius
1/ 3
mdn − m.∆t
rdn +1 = (2.31)
4
πρl
3
where t is the time step and l is the droplet density which is assumed to be constant in each
time step. However, the density change of the drop with the temperature is taken into account
In the case of a single-component isolated droplet, the expression for the mass flow rate
N
X ks
m = 2π rd ρ g Dvg 1 − X vs + Yvs Sh (Yvs − Yv∞ ) (2.32)
k =1 Yks
In addition, in case of binary mixture (for instance: n-heptane (1) and Nitrogen (2)) one has:
N
Xk X X W X
Dvg 1 − X v + Yv = D1g 1 − X 1 + Y1 2 = 2 2 D12 1 − X 1 2 = D12 (2.33)
k =1 Yk Y2 X2 W Y2
From Equation (2.27), the following expression for the heat flux can be obtained:
Qg = 2π rd λg Nu (T ∞ − Td ) − mhg (2.35)
The multi-component droplet evaporation model presented above has been implemented in
IFP-C3D code [38] already including the Abramzon and Sirignano droplet evaporation model
While there are some works that take into account the convective flow around the droplet
[39-42], we apply this model to a single-component isolated droplet to be able to compare its
results with the most realistic vaporization experiment of an isolated single droplet obtained
recently by Chauveau et al. [2]. The results of the new model are also compared with the
center of a large cube with the dimension of 1.5×1.5×1.5 m3 in order to avoid numerical
border effects. A uniform mesh of 3×3×3 cells is used. There is no air flow in the cube and
the calculation is performed at atmospheric pressure with several ambient temperatures. Table
The droplet vaporization phenomenon consists of two main periods; the heating-up
period when heat transfer increases the droplet temperature, and the vaporization period,
when the heat transfer from the gas to the droplet is dedicated to vaporize the droplet. In the
experiment of Chauveau et al. [2], the droplet was initially at room temperature. However, the
droplet surface regression is recorded while it is already heated during the placement of the
droplet inside the heating chamber (i.e furnace). Thus, the exact initial temperature of the
droplet T0 is unknown at the beginning of the measurement procedure. Figure 2 illustrates the
variations of squared dimensionless droplet diameters with time for several initial droplet
temperatures. According to this figure, it is shown that T0 is proportional to the ambient gas
unaffected by the heater. However, T0 reaches a temperature up to 340 K when the ambient
identical for all curves depicted in Figure 2 with the different T0, and is very close to the
experimental results. These numerical results may also indicate that the temperature
uniformity inside the droplet may be considered as an acceptable assumption especially for
high T0 values.
The results of the single-component droplet evaporation are presented in Figure 3 for
n-heptane and in Figure 4 for n-decane droplet for the different gas temperatures given in
Table 1 at normal gravity condition. As discussed in the Introduction, the differences between
the experimental data of Nomura et al. [4] and Chauveau et al. [2] (Figure 3(a) and (b)) are
due to the heat transfer through the fiber to the droplet in the experiments of Nomura et al.
[4]. According to Figure 3 and Figure 4, results of the new model are much more closer to the
experimental results of Chauveau et al. [2] than the results of the AS-1989 model (Abramzon
and Sirignano [9]). In fact, the new Stephan velocity (Equation (2.26)) and the enthalpy
diffusion in Equation (2.4) seem to be two crucial features for the evaporation rate and heat
flux from the gas to the liquid. Also, the results show that at atmospheric pressure and
regardless of the evaporation models, both ideal and real gas assumptions present a similar
behavior. Indeed, the temperature variation (as can be seen in Figure 3 for n-heptane and in
Figure 4 for n-decane droplets) leads to a similar behavior for both real and ideal gas models.
Because of the low computational time of the ideal gas model compared to the real gas model,
it is a good approximation for the thermodynamic equilibrium condition at low pressures for
the range of gas temperature used in this work. However, the deviation from the ideality may
At the micro-gravity condition, the natural convection effects are eliminated. Figure 5
ambient temperature for micro-gravity condition. As can be seen in this Figure, the
vaporization rate K has a quasi-linear evolution with gas temperature. In addition, the new
model has closer numerical results to the experiments than the AS-1989 model. The numerical
results of the calculation for both the normal and micro-gravity conditions are compared in
Figure 6 for two different gas temperatures. At the normal gravity condition, the droplet is
deformed during its lifetime due to the buoyancy force and non-spherical evaporation. The
buoyancy force increases the evaporation rate for both ambient temperatures. The differences
between the two conditions are more obvious in the vaporization period rather than the
This temperature represents the average temperature of the droplet which increases from the
initial value up to a constant. Indeed, as for the surface drop regression curves, two different
periods could be distinguished: the heating period, when the temperature of the droplet
increases, and the evaporation period, when the droplet temperature has a constant value
defined as the wet-bulb temperature. According to Figure 7, the droplet temperature is always
less than the saturation temperature of the liquid even at super-critical conditions (Figure 7
(a)). The wet-bulb temperatures obtained at different gas temperatures and at atmospheric
pressure are shown in Figure 8 for n-heptane droplet. These numerical results are compared
with the experimental data of Yuen and Chen [19] for normal gravity condition. They have a
similar behavior as the wet-bulb temperatures obtained from the new model are very close to
the experimental data. Indeed, the maximum deviation from the experimental results is less
than 2%.
Figure 9(a) shows the ratio of conduction heat flux to the enthalpy diffusion heat flux
for n-heptane droplet in the temperature range of 473 K to 973 K. In this Figure, the enthalpy
diffusion heat flux represents 55 to 75% of the conduction heat flux for the mentioned
temperatures. One can obviously conclude that the contribution of the enthalpy heat flux on
the total gas heat flux could not be neglected in a droplet evaporation model.
Figure 9(b) and (c) shows the heat fluxes in the gas and liquid phases. The heat fluxes
due to conduction and enthalpy diffusion of species have not the same behavior in the
heating-up period. Indeed, the heat flux due to conduction reduces as the droplet temperature
increases; whilst, the enthalpy diffusion of species increases as the surface mass fraction of
liquid increases. Since the temperature of the droplet reaches a constant value (at the
beginning of the vaporization period), there is no more heat flux inside the droplet. In this
situation, the heat flux from the gas to the liquid, which is dedicated to vaporize the droplet,
begins to reduce gradually due to the reduction of both conduction heat flux and heat flux due
to the enthalpy diffusion of species. According to Figure 9(b) and (c), in the vaporization
period, the conduction, the enthalpy diffusion and consequently the total gas heat flux
The instantaneous vaporization rate for n-heptane and n-decane droplets at different
ambient temperatures and at atmospheric pressure is illustrated in Figure 10. The vaporization
rate first increases during the heating-up period, when the temperature of the droplet (Figure
7) as well as the mass fraction of liquid at the droplet surface (Figure 9 (b) and (c)) increase,
and then decreases (since the mass fraction and temperature gradients are constant) due to the
Finally, a sensitivity study on the effect of the averaging parameter Ar in Equation (2.24)
on the evaporation rate has been conducted. Two values of Ar (Ar=1/2 and Ar=1/3) have been
tested using n-heptane droplets. The numerical results are shown in Figure 11 (a) and (b) in
terms of evaporation rates and in Figure 11 (c) and (d) in terms of droplet temperatures. For
the high gas temperature case, the best evaporation rate is obtained when the one-third rule
Ar=1/3 is used for determination of properties of gaseous mixture in the film region around
the droplet. On the other hand, at low ambient temperature, the difference between the
ideality. The effects of ambient pressure on n-heptane droplet vaporization are investigated at
sub/super-critical conditions using both ideal and real gas models. The results are presented in
Figure 12 for two ambient temperatures of 473 K and 973 K and for different pressures. It
appears that the ideal gas model is more sensitive to the pressure variations than the real gas
model especially at very high ambient pressures. Comparison of Figure 12 (a) and (b) shows
that the real gas model tends to vaporize the droplet more rapidly. In addition, one can
observe that at super-critical pressure conditions (more than 5 MPa), there is a large
difference between the real and ideal gas models prediction. In this situation the non-ideality
of the fluid becomes significant. At low or sub-critical pressure conditions (e.g., 2 MPa), the
in Figure 12 (c) and (d) for sub and super-critical ambient temperatures respectively. To
satisfy the equilibrium conditions in the deal gas model, the droplet temperature at
supercritical ambient pressures tend to be higher than in the real gas model. This behavior has
also been discussed by Givler and Abraham [43]. As shown in Figure 12 (d), at super-critical
pressures, the temperature of the droplet reaches the critical temperature of the liquid in the
ideal gas model. In fact, this limitation is applied to the model to avoid the temperature of the
droplet to be higher than the critical temperature of the liquid. The temperature limitation
approach is different for the real gas model. As illustrated in Figure 13, liquid temperature for
the real gas model is limited by the dew point line in the phase envelope diagram. The phase
envelope diagram depends on the composition of the mixture which changes at each time
step. In other words, there is not a unique phase envelope diagram for the entire droplet life
time. In Figure 13, phase envelope diagrams have been drawn for four different cases (cases
in Table 2). The mole fractions of each component and the wet-bulb temperatures (given in
Table 2) are obtained at the end of droplet life time. According to Figure 13, the dew point
and consequently the wet-bulb temperature of the droplet increases as the pressure increases
(from 2 to 10 MPa). Still, increasing the pressure from 10 to 15 MPa causes a reduction in the
dew point and subsequently the wet-bulb temperature at 15MPa. This phenomenon is in
accordance with the droplet temperature curves presented in Figure 12 (d). This numerical
study leads to the conclusion that the pressure and the temperature are the parameters that
affect the ideality of the fluid. At high pressure and/or low temperature conditions, different
behaviors of the droplet vaporization are observed with real and ideal gas models.
The model is also compared with the experimental data of Birouk [3] obtained by
nitrogen at a temperature of 293 K and a pressure of 0.1 MPa. In the experiment of Birouk
[3], the droplet is suspended on a vertical quartz fiber with 0.2 mm diameter and the
temperature of the drop is initially of 293 K. The numerical results obtained for the stagnant
condition are presented in Figure 14. A classical distillation is obtained. In the numerical
results, the most volatile component (n-heptane) is first completely depleted followed then by
n-decane. Indeed, in Figure 14 the numerical results show that the evaporation of n-decane
starts later than in the experiment. To explain this behavior, we have introduced a virtual
diameter dn-decane which is proportional to the initial volume fraction of n-decane αn-decane in
the droplet. In Figure 14(a), the value of αn-decane = 0.3 leads to a virtual diameter for n-decane
single-component, at least when (d/d0)2 reaches the value of 0.45. Actually, n-decane begins
to evaporate in the same time as n-heptane (at the beginning of the evaporation process) with
less evaporation rate which causes a reduction of the virtual diameter dn-decane. Thus, the
Figure 14(a) as well as in Figure 14(b) for the numerical model. This reasoning leads to the
conclusion that the initial volume fraction of n-decane in droplets was probably overestimated
by Birouk [3] in their experiments. A probable scenario may be that a certain fraction of n-
heptane was already evaporated before the beginning of the measurement procedure of Birouk
[3]. Reconsidering Figure 14(a) with a new initial composition of liquid mixture (65% n-
heptane) in the numerical calculation shows closer behavior between the numerical results
and experimental data. In this case, the results of Birouk [3] shows faster evaporation of the
first component (n-heptane) which is due to the effects of heat conduction through the fiber.
As can be seen in Figure 14(c), droplet temperature decreases at the beginning of evaporation
process and then increases to the ambient temperature. In the first part of the temperature
curve (Part 1: which mostly corresponds to the evaporation of n-heptane), the differences
between the droplet temperature and ambient temperature cause the heat flux through the fiber
to be important. Since the second component (n-decane) starts to evaporate, there is not a
considerable temperature difference between the droplet and the ambient gas. Thus, the effect
of conduction heat flux through the fiber is negligible in part 2 (Figure 14(c)). That is why the
slope of the experimental and numerical curves in Figure 14(a) are identical.
4 Conclusion
new Stephan velocity which recovers the mass conservation equation when Hirschfelder's law
is used for the diffusion velocity. In addition, the heat flux is assumed to be the summation of
heat conduction and enthalpy diffusion of species. Enthalpy diffusion accounts for energy
changes associated with compositional changes resulting from species diffusion. The
• The numerical results of the new model are in very good agreement with the
• The numerical results have shown that the heat flux due to the enthalpy diffusion of
• Infinite thermal conductivity assumption inside the droplet which has been invoked in
this study gives acceptable results compared to the experimental results of Chauveau et al. [2].
• The results show that the gravity and consequently, natural convection affects the
• A sensitivity study on the choice of the averaging parameter has shown that more
satisfying results are obtained with the one-third rule (Ar = 1/3) at sub and super-critical
ambient temperatures.
• At low pressure conditions, the ideal gas model could be a good assumption for the
droplet vaporization in the range of gas temperatures used in this study. However, at high
pressure and/or low temperature conditions, there is a large difference between the real and
ideal gas models. This leads to the conclusion that at high pressure conditions, a real gas
model that could better represents the thermodynamic equilibrium condition is needed.
component droplet (n-heptane/n-decane) test case. The numerical results are compared with
the experimental results of Birouk [3]. It is shown that in both numerical model and
experiments, the n-decane (the less volatile component) evaporates with the same evaporation
rate. However, this could not be observed in the evaporation of n-heptane (the most volatile)
component. It may be due to the over-estimation of experimental initial volume fraction of the
most volatile component in the droplet and due to the heat flux through the fiber.
Acknowledgment
The authors wish to acknowledge the experimental data provided by C. Chauveau and co-
workers from ICARE-CNRS and LME, Polytech’Orléans, France. The discussion with Dr B.
Cuenot from CERFACS Toulouse and Prof. T. Poinsot from IMFT Toulouse and their advices
Some properties like g, g, Dg, hi and g are calculated at the reference temperature Tref
(Equation (2.24)) while others like and are obtained at T . Density of the gas mixture is
calculated from the ideal/real gas equation of state with the appropriate temperature. For the
binary diffusion coefficient Dg in the gas mixture, the correlation of Wilke and Lee [31] has
been used. Then, the diffusion coefficient of each species into the gas mixture is obtained
N
1 − Yi X αWα
Dig = = α ≠i
(5.1)
N Xα N Xα
α ≠i
W α ≠i
Dα i Dα i
Dynamic viscosity and thermal conduction coefficients ( g (W.m-1s-1) and g (Kg.m-1s-1)) are
obtained at Tref:
32 32
AT A1 = 252 × 10−5 AT A1 = 1457 × 10−9
λg = 1 ref
with , µg =
1 ref
with (5.2)
Tref + A2 A2 = 200 K Tref + A2 A2 = 110 K
is obtained from Equation (5.2) at the ambient temperature T . Specific heat of liquid at
constant pressure Cpf is calculated from the enthalpy gradient at the reference temperature Tref.
Specific enthalpy hi for each component is obtained from the tabulated thermodynamic data
At the liquid-gas interface, for every component i in the mixture, the condition of
f i v = fi l (5.3)
where f is the fugacity. The v and l exponents stand for vapor and liquid phases respectively.
composition of that mixture. Reid et al. [31] introduced the fugacity coefficient φi as follows:
fi
φi = (5.4)
XiP
where Xi is the mole fraction of component i and P is the total pressure. Assuming a non-ideal
solution, the surface mole fraction of the droplet can be determined using the thermodynamic
equilibrium condition for the partial vaporization at the given temperature and pressure as
[45]:
where the fugacity coefficients are evaluated from an equation of state. The two-parameter
cubic equation of state is often used for hydrocarbon mixtures. It can be expressed by
RT am
P= − 2 . (5.6)
V − bm V + kbmV + k ′bm2
Where
am P bm P
A∗ = , B∗ = (5.8)
R 2T 2 RT
In Equation (5.7), Z is the compressibility factor which represents the deviation from the ideal
gas model. One of the well-known cubic equations is the Peng-Robinson (PR) equation [32].
For this equation, k and k' take on the integer values of 2 and -1 respectively. The two
X i , g X j , g ( ai a j ) (1 − k ) ,
1/ 2
am = ij bm = X i , g bi (5.9)
i j i
In the above Equations, R is the universal gas constant, i, Tc,i and Pc,i are the acentric factor,
critical temperature and pressure of component i and kij is the binary interaction parameter
which is assumed to be zero for hydrocarbon mixtures [31]. From the Peng-Robinson
bi
ln φi = ( Z − 1) − ln ( Z − B∗ ) +
A ∗
b i
− δ i ln
2 Z + B∗ k + k 2 − 4k ′ ( ) (5.12)
bm b
B k − 4k ′ m
∗ 2
2 Z + B ∗ k − k 2 − 4k ′ ( )
where
bi Tc ,i / Pc ,i 2ai1/ 2
= , δi = X j ,l a1/j 2 (1 − kij ) (5.13)
bm X j , g Tc , j / Pc , j am j
j
The mole fraction of each component is obtained from Equation (5.5) and Equation (5.12)
with an iterative method. The latent heat of vaporization for each component, Lv,i(Td), can be
calculated from the fugacity at constant pressure and composition of mixture, by:
RT 2 ∂
Lv ,i = −
Wi ∂T
( ln φi, g − ln φi ,l ) (5.14)
P ,n i
where Wi is the molar weight of component i. The latent heat of vaporization, Lv,g used in
N
X i , g Lv ,i
Lv , g = i =1
N
(5.15)
X i,g
i =1
In the case of an ideal fluid mixture, the thermodynamic equilibrium condition at the liquid-
where Pv,i(Td) is the vapor pressure of species i at the droplet temperature. The mass fraction
N
X iWi
W
Yi = X i i , W= i =1
N
(5.17)
W
Xi
i =1
(c) T ∞ = 973K
Figure 9: Heat flux in gas and liquid phases along with the surface mass fraction of n-heptane droplet with real
gas EOS at normal gravity condition, ( P ∞ = 0.1MPa ).
(a) n-heptane (b) n-decane
Figure 10: Vaporization rate versus time for different gas temperatures at the normal gravity condition and
( P ∞ = 0.1MPa ). The droplet initial diameter is 400 µm.
(a) T ∞ = 473K (b) T ∞ = 973K
(c) New model : 70% n-heptane, 30% n-decane, Fitted: 65% n-heptane, 35% n-decane
Figure 14: Comparison of vaporization rate of droplet with real gas EOS and normal gravity condition with the
experiment of Birouk [3] and the model of Abramzon & Sirignano (AS-1989) [9] ( P ∞ = 0.1MPa ). Part 1:
Transient temperature during the mixture evaporation, Part 2: Quasi-steady temperature during n-decane
evaporation
Table 1: Droplet vaporization case studies
473 300
623 310
n-Heptane
823 320
973 340
548 300
n-Decane
623 310
Table 2: High gas pressures numerical studies. Effects of the pressure on the wet-bulb temperature at T = 973K
[2] CHAUVEAU C., HALTER F., LALONDE A. and GOKALP I., An experimental
study on the droplet vaporization: effects of heat conduction through the support fiber,
ILASS 2008, 2008.
[4] NOMURA H., UJIIE Y., RATH H.J, SATO J. and KONO M., Experimental study on
high-pressure droplet evaporation using microgravity conditions, Proceedings of the
Combustion Institute, 26, p.1267-1273, 1996.
[5] HIROYASU H., KADOTA T., SENDA T. and IMAMOTO T., Evaporation of a
Single Droplet at Elevated Pressures and Temperatures : Part 1, Experimental Study,
Transactions of the Japan Society of Mechanical Engineers, 40,339, pp.3147-3154,
1974.
[6] MORIN C., CHAUVEAU C. and GOKALP I., Droplet vaporization characteristics of
vegetable oil derived biofuels at high temperatures, Exp.Therm.Fluid Sci., 21,1-3,
pp.41-50, 2000.
[7] GHASSEMI H., BAEK S.W. and KHAN Q.S., Experimental study on binary droplet
evaporation at elevated pressures and temperatures, Combus.Sci and Tech., 178,6,
pp.1031-1053, 2006.
[8] YANG J.R., WONG S.C., On the discrepancies between theoretical and experimental
results for microgravity droplet evaporation, Int.J.Heat Mass Transfer, 44,23,
pp.4433-4443, 2001.
[9] ABRAMZON B., SIRIGNANO W.A., Droplet Vaporization model for spray
combustion calculations, Int.J.Heat Mass Transfer, 32,9, pp.1605-1618, 1989.
[10] HOHMANN S., RENZ U., Numerical simulation of fuel sprays at high ambient
pressure: the influence of real gas effects and gas solubility on droplet vaporisation,
Int.J.Heat Mass Transfer, 46,16, pp.3017-3028, 2003.
[11] HARSTAD K., BELLAN J., Isolated fluid oxygen drop behavior in fluid hydrogen at
rocket chamber pressures, Int.J.Heat Mass Transfer, 41,22, pp.3537-3550, 1998.
[12] HARSTAD K., MILLER R.S. and BELLAN J., Efficient high pressure state
equations, AIChE J., 43,6, pp.1605-1610, 1997.
[13] CHAUVEAU C., CHESNEAU X. and GOKALP I., High pressure vaporization and
burning of methanol droplets in reduced gravity, Advances in Space Research, 16,7,
pp.157-160, 1995.
[14] NOMURA H., UJIIE Y., Evaporation behavior of fuel droplet at high temperatures
and high pressures under microgravity conditions, Nippon Kikai Gakkai Ronbunshu,
B-hen, 61,591, pp.4137-4143, 1995.
[15] GOGOS G., SOH S. and POPE D.N., Effects of gravity and ambient pressure on
liquid fuel droplet evaporation, Int.J.Heat Mass Transfer, 46,2, pp.283-296, 2003.
[16] ZHANG H, RAGHAVAN V. and GOGOS G., Subcritical and supercritical droplet
evaporation within a zero-gravity environment: Low Weber number relative motion,
International Communications in Heat and Mass Transfer, 35,4, pp.385-394, 2008.
[17] INCROPERA F, P. and DEWITT D.P., Fundamentals of Heat and Mass Transfer,
New York, USA, 1996.
[18] RANZ W.E., MARSHALL W.R., Evaporation from drops part I, Chemical
Engineering Progress, 48,3, pp.141-146, 1952.
[19] YUEN M.C., CHEN L.W., On drag of Evaporating Liquid Droplets, Combus.Sci and
Tech., 14,4-6, pp.147-154, 1976.
[21] RA Y., REITZ R.D., A vaporization model for discrete multi-component fuel sprays,
Int.J.Multiphase Flow, 35,2, pp.101-117, 2009.
[22] SAZHIN S., Advanced models of fuel droplet heating and evaporation, Progress in
Energy and Combustion Science, 32,2, pp.162-214, 2006.
[25] HIRSCHFELDER J.O., CURTISS C.F. and BIRD R.B., Molecular Theory of Gases
and Liquids, New York, 1954.
[27] HUBBARD G.L., DENNY V.E. and MILLS A.F., Droplet evaporation: Effects of
transients and variable properties, Int.J.Heat Mass Transfer, 18,9, pp.1003-1008,
1975.
[28] SAZHIN S., KRISTYADI T., ABDELGHAFFAR W.A. and HEIKAL M.R., Models
for fuel droplet heating and evaporation: Comparative analysis, Fuel, 85,12-13,
pp.1613-1630, 2006.
[29] SAZHIN S.S., ELWARDANY A., KRUTITSKII P.A., CASTANET G., LEMOINE
F., SAZHINA E.M. and HEIKAL M.R., A simplified model for bi-component droplet
heating and evaporation, International Journal of Heat and Mass Transfer, 53,
pp.4495-4505, 2010.
[30] MAQUA C., CASTANET G. and LEMOINE F., Bicomponent droplets evaporation:
Temperature measurements and modelling, Fuel, 87, pp.2932-2942, 2008.
[31] REID R.C., PRAUSNITZ J.M. and POLING B.E., The properties of gases and
liquids, New York, 1988.
[34] AMSDEN A., O'ROURKE P. and BUTLER T. Kiva-II: A Computer Program for
Chemically reactive Flows with Sprays. LA-11560-MS, Los Alamos.,LA-11560-MS,
Los Alamos., 1989.
[35] RAMSHAW J.D., Fluid dynamics and energetics in ideal gas mixtures, Am.J.Phys.,
70, p.508, 2002.
[36] BIRD R.B., STEWART W.E. and LIGHTFOOT E.N., Transport Phenomena, New
York, 1960.
[37] HARSTAD K., BELLAN J., An all-pressure fluid drop model applied to a binary
mixture: Heptane in Nitrogen, International Journal of Multiphase Flow, 26,10,
pp.1675-1706, 2000.
[38] BOHBOT J., GILLET N. and BENKENIDA A., IFP-C3D: an Unstructured Parallel
Solver for Reactive Compressible Gas Flow with Spray, Oil & Gas Science and
Technology-Revue de l Institut Francais du Petrole, 64,3, pp.309-335, 2009.
[39] BIROUK M., GOKALP I., A new correlation for turbulent mass transfer from liquid
droplets, Int.J.Heat Mass Transfer, 45,1, pp.37-45, 2002.
[41] ACKERMAN M., WILLIAMS F.A., Simplified model for droplet combustion in a
slow convective flow, Combustion and Flame, 143,4, pp.599-612, 2005.
[43] GIVLER S.D., ABRAHAM J., Supercritical Droplet Vaporization and Combustion
Studies, Progress in Energy and Combustion Science, 22,1, pp.1-28, 1996.