7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Numerical Simulation of Evaporating Droplets with Chemical Reactions using a
Volume of Fluid Method
P. Kellert, P.A. Nikrityuk B. Meyert and M. MtillerHagedorn&
t Department of Energy Process Engineering and Chemical Engineering, TU Bergakademie Freiberg,
09596 Freiberg, Germany
CIC Virtuhcon, TU Bergakademie Freiberg, Freiberg, 09596 Freiberg, Germany
AIR LIQUIDE Forschung und Entwicklung GmbH
Peter.Keller@iec.tufreiberg.de and Petr.Nikrityuk@vtc.tufreiberg.de
Keywords: Evaporation, Combustion, VOF
March 24, 2010
Abstract
The main aim of this work is the detailed numerical investigation of the unsteady behaviour of chemical reacting
nheptane droplets in a forced convective environment under different temperature conditions. The numerical model
includes chemical reactions occurring on the deformable nheptane droplet surface, fluid dynamics of the droplets,
and mass and heat transfer due to the chemical reactions. In particular, based on the two and threedimensional
numerical simulations, the impact of heat transfer on the dynamics of evaporation and combustion of liquid hydro
carbons (nheptane as model substance for fuel oil) is investigated. Following the work done recently by Schlottke
and Weigand (2008) the Volume of Fluid (VOF) approach is used to describe the droplet behavior in gas flow.
Chemical reactions on the surface of the fuel droplets are taken into account by introducing source terms in the
mass conservation equations for nheptane, oxygen, carbon dioxide and steam. The evaporation source is included
in the VOF and nheptane equations. The model and the code have been validated against published numerical
and experimental data. In particular, the comparisons were performed for the time dependent decrease of droplet
diameter caused by evaporation of a nheptane droplet in a quiescent hot atmosphere. All validation test cases showed
acceptable agreement.
Nomenclature
Roman symbols
A preexponential factor (m3 mol 1 s 1)
cp heat capacity (m2 s 2 K 1)
C concentration (mol mn3)
d diameter (m)
D diffusion coefficient (m2 s 1)
err error ()
Ea activation energy (kg m2 m1ol 1 s2)
F, surface tension force (kg m 2 s 2)
g gravitational constant (m s1)
Ah, evaporation enthalpy (nm2 s 2)
H formation enthalpy (m2 s 2)
k reaction rate coefficient (m3 1ol 1 s1)
M molar mass (kg mol 1)
n surface normal ()
p pressure (kg m 1 s 2)
r reaction rate (kg m 3 s 1)
R universal gas constant (kg m2 mol 1 s2 K 1)
S source term (different units)
t time (s)
T temperature (K)
U velocity (m s 1)
Sv diffusion volumina (m3 kmol1)
X molar fraction ()
Y mass fraction ()
Greek symbols
a liquid volume fraction ()
K curvature (rn 1)
A heat conductivity (kg m s3 K 1)
p dynamic viscosity (kg m 1 s 1)
p volume fraction ()
p density (kg m 3)
a surface tension (kg s 2)
Subscripts
boil boiling
c critical
exp experimental
gp gas phase
i, j, k species index
1 liquid
m mass
max maximum
mid intermediate
M mixture
Nom data by Nomura
q heat
rel relative
st static
sim simulated
V volume
chem chemical
evap evaporation
sat saturation
Introduction
The gasification of liquid hydrocarbons is of great im
portance for chemical and energy industries. Basically,
gasification involves chemical reactions between two
and more phases, e.g. between liquids and/or solids
and gas, where one of the phases can be considered
as dispersed and the other one as continuum. In gasi
fication of liquid fuels, droplet evaporation is the key
phenomenon by achieving a better physical understan
ding of the whole processes. Practical relevance of fuel
droplets evaporation and vapor combustion has led to
numerous theoretical investigations of this phenomenon.
During the past decade igniliik.ii progress has been
reached in understanding droplet behavior in convective
environment, e.g. see works of Dwyer (1989), Dwyer et
al. (2000) and Zhang (2004).
In earlier numerical works devoted to the modeling
of droplet heating and evaporation, the droplets were
considered as unmovable spherical axisymmetric bod
ies placed in gas flow (see Dwyer (1989)). Rapid deve
lopment of powerful computers during the last ten years
allowed to perform simulations on threedimensional ar
rays of chemical reacting nheptane droplets at interme
diate Reynolds numbers (see Dwyer et al. (2000)). The
shape of the reaction zone was significantly influenced
by the position of droplets. The rate of vaporization in a
droplet array is strongly dependent on the distribution of
the single droplets.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Recently, ,ignilik .il success has been reported by the
modeling of multiphase chemical processes including
chemically reacting bubble swarms with deformable in
terfaces in a gasliquid reaction system, see e.g. Schlott
ke and Weigand (2008), Koynov et al. (2005) or Radl et
al. (2008). But basically, only mass transfer was consi
dered by the modeling of chemical reacting droplets.
Motivated by this facts, this work is devoted to the
numerical investigation of the unsteady behavior of che
mical reacting nheptane droplets in a forced convective
environment under different temperature conditions.
Model Formulation
In this work the VOF method (Hirt and Nichols (1981))
is utilized to describe the behaviour of chemical reac
ting deforming droplets in incompressible flow. In order
to simplify the conservation equations stated below, the
following basic assumptions have been carried out:
Marangoni effect on the droplet surface is neg
lected,
Radiative heat transfer has not been taken into ac
count,
Local thermodynamical equilibrium is reached on
the droplet interface,
Soret and Dufour effects are neglected.
In compliance with these assumptions the conservation
equations describing the behavior of chemical reacting
droplets in a gas flow have the following form:
Continuity equation
S+ V (pU) 0 (1)
at
Momentum equation
d(pU)
a + V (pU) U
at
V / t [V + (U)T]
+pg Vp F, (2)
where F, is the surface tension force F, = 7 aKi.
Heat transfer equation
d(pcT) +V.(pcpUT) V (AVT)+ q,chem q,evap*
at
(3)
Following Schlottke and Weigand (2008) the VOF
method is used to describe the behavior of liquid fuel.
For more information about this method the reader is re
ferred to Hirt and Nichols (1981). The liquid phase is
tracked by solving the following equation:
a+ v (aU)
att
SV,evap
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
with volumetric evaporation source Sv,evap described be
low. Diffusion effects of the liquid are neglected because
of their small order compared to diffusion effects of the
considered gases.
It should be noted that in contrast to the work of
Schlottke and Weigand (2008) eq. (1) does not have the
source term characterizing the generation or vanishing
of a new volume due to evaporation or condensation. In
this work the effect of phase change on the continuity
equation is included via density, which is not a constant
in interface cells.
To consider chemical reactions in the gas phase species
conservation equations are implemented for C7H16, 02,
HO0, C02 and N2 to describe the combustion of n
heptane
C7H16 + 1102 * 7CO + 8 H20.
The equation for inert species N2 has the form
) + V (pYN2 U) V (pDN, VYN2) (6)
at
and the one for the vaporized and combusted species
C7H16 the following one
8(pYc, H16)
PC + V.(PYCH16U)=
at
V (pDcH16VYc7H16)
+ Sm,evap + Sm,chem,C7Hi6. (7)
For all other species i = 02, H20, CO0 the species
equation reads
a(p) + V (pYiU) = V (pDiVYi) + Sm,chem,i. (8)
at
It should be noted here that all species conservation
equations are valid only for the gas phase including
the interface. Thus, diffusion coefficients of chemical
species were set to zero inside liquid phase. This con
dition implies the droplet surface to be nonpermeable,
which is equivalent to the following boundary condition
on the droplet surface:
DiVY = 0.
Mixture Properties, Heat and Mass Sources
To avoid mixing of liquid and gas mixture, diffusion of
gas into liquid is anticipated by
Di 0 {D 1 (10)
D i 1
DIAM ,a < 1
with diffusion coefficient Dim of species i into mixture
M. This coefficient is calculated according
(FAIRBANK/WILKE)
D i 1 ik
Ek=l,ki oDk
with binary diffusion coefficient Dik obtained by
(SHETTLER/GIDDINGS/FULLER)
Di= 1.013. 10 Tl 1 + k 1o5
pt [( ti) + (E V) ]
(12)
with diffusion volumina E v. To avoid overshoot of
diffusion coefficient in the interface region an artificial
minimum value for diffusion of liquid into gas mixture
10 M2
is assumed (i.e. DiM = 1010l ")
To determine mixture properties like mixture molar
mass MM, density p or heat capacity c, thermodynamic
relations are used e.g.
YC07HI + Y02 + YH o + Yco2 + YN2 + Y1
M YC H0+Y Yo2 Yco2 YH20 YN2
YC7H16+ Y 2 H2+
MC7H16 + M02 + MC H2 MN2
p = ap + (1 a)pp
CP = YC7H16CPC7H16 +Yo2C02o YH2("
+Yco02pco02 + YN2CPN2 + CYCp
with mass fraction of liquid Yi, density of gas phase pgp
and
YCH16 + Yo2 + YH2O + YC2 + YN2 + Y 1. (13)
Liquid and gas properties like evaporation enthalpy Ah,
or saturation pressure Psat and molar fraction Xsat are
calculated using formulas given by VDIWrirmeatlas
(2006)
Ah,
In Psat
In
Pc
A( ) +C +C (T )_
A T T 
+c(1
T\
T \3 (
^ +0 *
Tr 6]
Tr
with tabulated coefficients A, ..., E.
Chemical kinetics are assumed to follow Arrhenius law
(see e.g. Warnatz et al. (2006) for formulas). The equi
librium constant for reaction (5) is calculated according
ATb .exp( E)
k = 4T ex RT
with A 500 b 0 and Ea 66067.7 "9 as
mol!s' mol s2
given values by simplest CHEMKIN input file.
The source terms in equations (3), (4) and (7) are speci
fied using relations described e.g. in Turns (2000),
Schlottke and Weigand (2008), Zhang (2004) or Fra
ckowiak et al. (2010)
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
first order scheme and fluxes are interpolated using lim
ited linear TVD scheme (Total Variation Diminishing).
Validation and Results
DCTH16 PgpY HI
0 VYC7H1. i
Y07H16 P1
A I
VT flK,
D07H16 PgpVYCyHi6 r
1 YCY7H1
A
 VT fn,
Aht
 hvSm,evap,
,k 0.25 01.5
kC7H16 02
n
E ,Sm,chem,i
i1
with surface normal gradient n = V~ surface curva
ture K, formation enthalpies H, for all species i, concen
trations CC7H16 YCYH16 M H16 and Co = Yo02 M
and specialized gradients VYC H16 and VT.
If cells contain the interface of liquid and gas phase,
these gradients are calculated by setting saturation va
lues respectively boiling temperature at that location and
retain corresponding values in neighboring cells. It
should also be noticed that these source terms for evapo
ration are limited to values above zero to avoid conden
sation effects that are not discussed here.
Numerics
The open source code OpenFOAM was used to imple
ment and solve the equations above. In each time step
the properties of fluids are corrected with temperature
calculated at the end of the last time step followed
by adjusting evaporation and chemistry sources. In a
next step the volume of fluid equation (4) and species
transport for each species (equations (6), (7), (8)) are
solved. In a last step pressure correction (continuity
equation (1)), momentum equation (2) and temperature
equation (3) are discretized. Finally the next time step
size is calculated using iterative time advancement
controlled by Courant Friedrichs Lewy number.
In general, discretization of gradients is processed using
Gaussian linear interpolation central differencing. For
the gradients of liquid volume fraction a, concentration
gradient VYc7H16 and temperature gradient VT a
surface normal gradient discretization scheme is used.
Time differencing is discretized using Eulerian implicit
Because it is almost impossible to validate combustion
fiK or gasification of entire oil burners, simplified geome
tries are used and validation of single droplet evapora
(17) tion and combustion is shown.
In a first step analytical solutions are validated with the
implemented solver under following assumptions:
SV,evap
Sm,evap
Sq,evap
Sm,chem,C7 H16
Sq,chem
dd2 8DAB 1
dt P In 
dt P1 \1
with intermediate gas density p, diffusion coefficient
DAB of vapor into gas phase, surface mass fraction of
vapor YA,s and vapor mass fraction in sufficient far dis
tance from surface YA,oo. For analytical solution the sur
face temperature Ts is supposed to be the same as the
temperature of surrounding air T,. The surface mass
fraction of vapor is assumed to be at saturation mass
fraction YA,S = sat and YA,oo 0. The analytical
solution then yields
with
d (t)= :, Kt
K 8PDAB In 1
Pi 1
YA,mo
YAS)
The values for used properties are given in table 1 lea
ding to results given infigures 1 and 2.
Gas phase is threaten as ideal gas mixture,
Diffusion coefficients, heat conductivities and ca
pacities and densities of gas and liquid as well as
saturation pressure of vapor remain constant,
Transition to supercritical fluids is not considered,
Quiescient atmosphere is assumed and therefore no
momentum and continuity equation are solved,
Mixture properties are calculated at temperatures
next to the temperature of the liquid phase be
cause of higher heat capacities and conductivities
of the liquid (in comparison: Turns (2000) uses
T = (Tboil + Typ)/2 as intermediate temperature).
The considered environment is a water droplet with a
diameter below 100 pm surrounded by air with tempe
ratures in first case below and in second case above the
boiling temperature of water (Tboil 373 K). The dif
ferential equation of droplet evaporation with tempera
tures below the boiling temperature reads according d2
law (see Turns (2000))
InK
Table 1: Used properties for cases 1 till 3: T < Tboil
It can be seen in figure 1 that the results of simulation
almost perfectly match the analytical graph except of
some fluctuation caused by production of vapor at the
interface which is not included in analytics.
d 104 m
Ts To 363 K
Yat at Ts 0.58
Pi 1000
p at Ysat 0.7124
DAB 7.6.108 M
2
K 3.1 107 M
td 0.032
but with
d2(t)= ,Kt
S8A (Cp (T 
At
(25)
Tboil)
+1
0 0005 001 0015 002 0025 003
time in s
Figure 1: Case 1: constant intermediate density p at
interface.
2e09
0
0 0005 001 0015 002 0025 003 0035
time in s
Figure 2: Case 2/3: variable density p at interface with
different mesh sizes and time steps.
with intermediate heat conductivity of gas A, gas tem
perature To, set equal to critical temperature Tc and in
termediate heat capacity of gas mixture cp.
The values used for that case are given in table 2 and
results are given in figures 3 and 4.
d 104 m
Ts 373 K
To 647 K
Pi 1000 9
A 0.036 k'm
cP 1011 2176
Ah, 2265000 m
K 3.1108
td 0.32 s
Table 2: Used properties for cases 4 and 5: T > Tboil
Again analytical solution of differential equation can be
approximated best by simplification of properties. For
that purpose temperature was hold constant for calcula
ting surface temperature gradient VT in case 4 (figure
3). In case 5 (figure 4) results of real behaviour can be
seen. Obviously, temperature at the interface is rising
and therefore temperature gradient is increasing. Evapo
ration acceleration leads to decreasing slope of graph as
seen in figure 4.
In a next step the validation case of Zhang (2004)
was reconstructed to include inflow influence on evapo
ration. In contrast to the former settings, nitrogen was
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
In figure 2, however, the graphs of simulation show
grid depending behaviour. While in the case with coarse
grid the "equilibrium" mixture density is reached later
because of larger cell size and therefore smaller evapo
ration sources, the slopes of curves in both cases are
approximately equal. Coincidentally, coarse grid graph
is matching the lower boundary of analytical solution
while the fine grid one exceeds that solution. However,
the better approximation of real behaviour would be the
slowlier case.
Considering droplet combustion, where temperature of
surrounding gas is higher than boiling temperature of
liquid, another analytical solution is proposed. The dif
ferential equation then reads (Turns (2000))
dd2 8kgp (Cpg(Too Tboil) )
dt *. Ah, (24)
and analytical solution leads again to
S,,, in 1 ,in Ts=363 K
S I. Ts =const 
S_
2e09
005 01 015
time in s
02 025 03
Figure 3: Case 4: constant surface temperature.
2e09
0 002 004 006 008 01 012 014 016 018
time in s
Figure 4: Case 5: variable surface temperature.
injected at the inlet of the geometry. To show the ef
fect of different boundary conditions two different cases
were considered. Based on the experimental results
given by Nomura et al. (1996) two inlet temperatures
TN2, 471 K and TN,2 745 K were chosen to
evaporate one nheptane droplet of size d 690 pm at
pressure p 1 atm. To match the inflow parameters the
nitrogen inlet velocity was set to UN2 0.1 leading
to a Reynolds number of Re z 3 for the used geometry
settings. With grid size z 175000 cells the results show
acceptable agreement compared with the experiments of
Nomura listed in table 3.
As shown infigure 5 andfigure 6 the results match quite
well with experimental data. Even if the maximum va
lues of evaporation constants Ksim,max are overpredicted,
the intermediate values Ksim,mid are determined with ac
ceptable deviations.
The last considered case deals with the combination of
evaporation, convection and combustion. For that rea
son inspired by the investigations of Dwyer et al. (2000)
a two and a three dimensional droplet array is simulated.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Texp,Nom 471 K 745 K
Kexp,Nom 0.117 mm 0.390 mm2
2 2
Ksim,max 0.152 mm 0.559 mm
2 2
Ksim,mid 0.129 mm 0.387 m
errrel,max 0.299 0.43
errrel,mid 0.1026 0.0077
Table 3: Comparison simulation * Nomura's experi
ments
1 02
simulation with T=471 K
Kexp Nom=0 117mm /s experimental ....
1 simulation with T=745K ..........
*....... K p =0 390mm/s experimental
098
096
a 094
092
09 ".
088
086
0 005 01 015 02 025 03 035 04 045 05
t/d02 In s/mm2
Figure 5: Simulation results with TN2,1 471 K and
TN2,2 745 K compared with experimental data
1 06
1 04
1 02
1
098
 096
094
092
09
088
086
do=690tm, T=745K 
:... ......... .. ............ ....
0 005 01 015 02
t/d02 In s/mm2
025 03 035
Figure 6: Maximum and intermediate evaporation con
stants Ksim,max and Ksim,mid for TN2 = 745 K
In the 3D case the droplet positions are fixed in space by
setting the velocities to zero if a tolerance value for the
volume fraction of liquid is met (i.e. a > 0.01). Inflow
temperature was set to T = 1000 K with the effect of
enhancing evaporation at droplet surface but only slow
temperature increase. The results of simulation are sum
marized in figures 7 till 10 which show isosurfaces of
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
T, C02, C7H16 and 02, respectively. It can be seen
that at inflow value T z 1500 K sudden ignition occurs
and temperature increases towards adiabatic flame tem
perature (T z 2270 K). Only in some sections the tem
perature rises beyond this value T > 2800 K which is
caused by transient mixing and combustion mechanism
but also by approximated liquid and gas properties.
Finally, the 2D version of droplet array was considered.
Figure 9: 3D Droplet Array: isosurfaces of C7H16.
T: 400 700 1000 1300 1600 1900 2200 2300
Figure 7: 3D Droplet Array: isosurfaces of T.
Figure 10: 3D Droplet Array: isosurfaces of 02.
tivity is used (see VDIWarmeatlas (2006)):
YC02: 0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04 0.045 0.05 0.055
Figure 8: 3D Droplet Array: isosurfaces of C02.
With inflow temperature T2 = 1000 K and operating
pressure P = 20 bar that case is oriented towards the
model of Dwyer et al. (2000). Additionally, the droplets
are assumed to move free with flow and gravitional vec
tor is aligned with positive x coordinate. For this case
another thermodynamical relation for thermal conduc
AC7H16
with A
B
C
PD
and E
A +BT + CT2 + DT3 + 4 (26)
0.2149,
2.998 104,
3 107,
1.03 1010
1.2. 1013
Because the droplets considered in 3D case suffered in
finite heating above critical temperature, that relation,
which has shown better results, was used in the 2D array
version in order to test this modification and will be in
vestigated in future 3D simulations. Ty,, ... 11 and 12
reflect the nonuniform heating from interface to inner
droplet area within a time interval of At = 710 /s.
Nevertheless, depending on the volume fraction of gas
inside the cells, the temperature might also arise above
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
time=0.00211 s
720
10C
I C
600
500
400
Figure 11: 2D Droplet Array: contours of a with con
tour lines of T.
alpha
09
07
06
05
04
003
02
0 1
0
alpha: 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Figure 14: 2D Droplet Array: contours of a at t =
2.11 ms.
time.O.00211 s
YC7H16. 0 0.1 02 0.3 0.4 0.5 0.6 0.7
Figure 15: 2D Droplet Array: contours of C7H16 at
t 2.11 ms.
time..0.0493 s
2E 05 T
4E05 000
890
6E05 700
600
8E005 5000
0.00015 0.0002 0.00025 0.0003
X
Figure 12: 2D Droplet Array: contours of a with con
tour lines of T.
the boiling and critical temperature. This is one of the
disadvantages of VOF methods: as consequence of mo
ving droplets the interface is smeared over several cells
and hence diffusion occurs also into inner liquid cells.
With mixed properties, temperature rise of droplet is
faster than it should. This can be easily justified by using
the fact that heat conductivity of the gas is rising and the
one of the liquid is decreasing with increasing tempera
ture.
time=O s
alpha: 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1
Figure 13: 2D Droplet Array: contours of a at initial
ization.
alpha: 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 08 0.9 1
Figure 16: 2D Droplet Array: contours of a at t
4.93 ms.
time.0.00493 s
YC7H16: 0 0.1 02 0.3 0.4 0.5 0.6 0.7
Figure 17: 2D Droplet Array: contours of C7H16 at
t= 4.93 ms.
One of the advantages of VOF methods can also be
observed: because of the fact that droplets are not fixed
by given mesh, they have the ability to affiliate to larger
droplets and hence influence the subsequent flow, mi
xing and thus combustion. The time varying droplet
movement and affiliation are illustrated infigures 13, 14
and 16, the time varying formation of nheptane in fi
gures 15 and 17.
time0.00029 s
time0.001 .
time0.00424 s
1600 o /
S1200 60
S52000 s
Figure 18: 2D Droplet Array: contours of YC H,1 with
contour lines of T next to outlet.
The final remark to the results is related to the adiabatic
flame temperature. Considering complete combustion
of nheptane with hot air the adiabatic flame tempera
ture should yield T = 2274 K according Turns (2000).
Although the field has shown a quasisteady solution af
ter 0.00424 s and the distribution of nheptane and oxy
gen is not changing over several time steps, the tempe
rature is still rising. This can be explained by the simple
fact that density decrease with rising temperature is not
considered (see assumptions in section "validation and
ic',ilii ) and hence chemistry source is not affected by
concentration change. In that case only rising tempera
ture influences chemical reaction rates leading to higher
temperatures even if the gas reaches steady composition.
So the last time step with ,igiik .ill change in composi
tion as seen infigure 18 was taken to compare the values.
It can be seen that the temperature field has its maxi
mum quite close to the adiabatic temperature with the
small limitation that this value is obtained with incorrect
density values.
Discussions
In further studies more complex reaction mechanisms
will be included to get better approximations of tempe
rature and chemical product formation. This would al
low not only to evaporate more than one species but also
to consider pyrolysis.
In order to determine high pressure and high tempe
rature influence on fluid properties, real gas and liquid
behaviour will be investigated.
Conclusions
A mathematical model for the prediction of the unsteady
behaviour of chemical reacting nheptane droplets in
a forced convective environment has been developed.
A combination of more complex thermodynamical re
lations was used to describe a more detailed model
containing chemical reactions and further refinement of
properties. The model and the code have been validated
against published numerical and experimental data. In
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
particular, it was shown that temperature has a signifi
cant influence onto evaporation and hence interface dy
namics. Taking properties as constants leads to results
comparable with analytics but real behaviour cannot be
reconstructed. That effects could also be observed by in
vestigations concerning flame temperature and thermal
conductivities. The more values are kept constant du
ring the simulations the larger will be the deviation from
reality. The next step should be the implementation of
a model with all properties being temperature and pres
sure dependent leading to a system of relations which
describes thermodynamics exactly.
Acknowledgements
P.K. appreciates the German Federal Ministry of Eco
nomics and Technology (Project ID 0327113B), the
Saxon Ministry of Science and the Fine Arts (Project
ID 12372) by means of the European Regional Develop
ment Fund (ERDF) and Air Liquide for collaboration in
COORAMENT project and financial support.
P.A.N. appreciates the financial support of the Go
vernment of Saxony and the Federal Ministry of Edu
cation and Science of Federal Republic of Germany as a
part of CIC VIRTUHCON.
The authors are thankful to P. Seifert and K. Uebel for
their helpful comments.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
References
Dwyer, H.A., "Calculations of droplet dynamics in
high temperature environments", Progress in Energy and
Combustion Science 15, p. 131158, 1989
Dwyer, H.A.; Stapf, P.; Maly, R., "Unsteady vaporisa
tion and ignition of a threedimensional droplet array",
Combustion and Flame 121, p. 181194, 2000
Schlottke, J.; Weigand, B., "Direct numerical simula
tion of evaporating droplets", Journal of Computational
Physics 227 p. 52155237, 2008
Hirt, C.W.; Nichols, B.D., "Volume of fluid (VOF)
method for the dynamics of free boundaries", Journal
of Computational Physics 39 (1), 1981
OpenCFD, "OpenFOAM : open source CFD",
http://www.openfoam.com/, 20042010
Turns, S.R., "An Introduction to Combustion Con
cepts and Applications", McGrawHill Higher Educa
tion, 2000
VDIGesellschaft Verfahrenstechnik und Chemie
ingenieurwesen, "VDIWairmeatlas", Springer Berlin
Heidelberg, 2006
Warnatz, J.; Maas, U.; Dibble, R.W., "Combustion",
Springer Berlin Heidelberg, 2006
Zhang, H.; Gogos G., "Numerical research on a vapo
rizing fuel droplet in a forced convective environment",
International Journal of Multiphase Flow 30 p. 181198,
2004
Frackowiak, B.; Lavergne, G.; Tropea, C.; Strzelecki,
A., "Numerical analysis of the interactions between
evaporating droplets in a monodisperse stream", Inter
national Journal of Heat and Mass Transfer Vol. 53 p.
13921401, 2010
Koynov, A.; Khinast, J. G.; Tryggvason, G., "Mass
transfer and chemical reactions in bubble swarms with
dynamic interfaces", AIChE Journal Vol. 51 p. 2786
2800, 2005
Radl, S.; Koynov, A.; Khinast, J. G.; Tryggvason, G.,
"DNSbased prediction of the selectivity of fast multi
phase reactions: hydrogenation of nitroarenes", Chemi
cal Engineering Science Vol. 63 p. 32793291, 2008
Nomura, H.; Ujiie, Y.; Rath, H. J.; Sato, J.; Kono, M.,
"Experimental study on high pressure droplet evapora
tion using microgravity conditions", 26th Symposium
(Int.) on Combustion p. 12671273, 1996
