Applied Energy 64 (1999) 317ą321
www.elsevier.com/locate/apenergy
Simulation of vapour explosions
Nail Suleiman Khabeev*
Department of Mathematics, University of Bahrain, Bahrain PO Box 32038, Bahrain
Abstract
Vapour explosions, also called thermal detonations, occur when a liqui�ed gas comes into
contact with water. During the contact of the two liquids with considerably di�erent tem-
peratures, intensive boiling of one of them takes place which is accompanied by an explosive
increase in pressure. A similar phenomenon also arises in the cooling systems of nuclear-
power stations, when, as a result of some accident, the heated particles of the nuclear fuel settle
in the cold water. This leads to the explosive boiling of the liquid and to a rapid increase of the
pressure. This paper presents analytical and numerical modelling of the behaviour of small-
scale single drops of liqui�ed gas in water. # 1999 Elsevier Science Ltd. All rights reserved.
Keywords: Explosive boiling; Liqui�ed gas; Sea water; Vapor explosion
1. Introduction
Plans to exploit oil and gas �elds in the Gulf region usually include options to
transport LNG (liqui�ed natural gas) to markets in Europe and America. The pos-
sibility of sinking, ship-ship collision or grounding needs to be considered in the
development plans. Such accidents can of course occur near towns and other den-
sely-populated areas. If a LNG ship is involved in an accident, such that large
quantities of liqui�ed gas escape, a large explosion can occur [1].
LNG ships are designed and operated to very high standards, because of the dan-
gers involved in handling LNG. Nuclear power stations, similarly, are designed and
operated to very high standards. Nevertheless, accidents do occur in nuclear-power
stations [2]. It does not take much imagination to visualise an accident along the Gulf
coast involving a LNG ship. The potential consequences of such an accident need to
be known.
* Corresponding author. Tel.: [phone]; fax: [phone].
E-mail address: [e-mail] (N.S. Khabeev).
0306-2619/99/$ - see front matter # 1999 Elsevier Science Ltd. All rights reserved.
PII: S0306-2619(99)00052-5
318 N.S. Khabeev / Applied Energy 64 (1999) 317ą321
2. Formulation of the problem
Let us consider the behaviour of a spherical drop of liqui�ed gas, surrounded by a
spherical vapour layer in the liquid. The processes occuring within the drop, bubble
and around the bubble in the surrounding liquid are assumed to be spherically sym-
metric, so all parameters depend on time t and the radial Eulerian coordinate r only.
We assume that the ``cold'' drop and surrounding ``hot'' liquid are incompressible
and ideal, and the vapour can be described by the equation of state for the perfect gas.
We will suppose also pressure uniformity ensues within the bubble, which prevails
when the size of the bubble is much less than the length of a sound wave in the gas
[3]. Within the framework of the assumptions made, the continuity, heat and state
equations for the gas phase will have the form
@ 1 @ Wgr2
g
g
� � 0
@t r2 @r
@Tg @Tg 1 @ @Tg dPg
1ą
Cg � Wg � lgr2 �
g
@t @r r2 @r @r dt
Pg � RgTg; btą4r4atą
g
where Cg; lg; Rg; Wg are the speci�c heat of the gas at constant pressure, thermal
conductivity, gas constant, and velocity of the gas; btą is the drop radius, atą is the
radius of the bubble. Here the subscripts d, g and ` refer to the parameters of the
drop, gas, and surrounding liquid, respectively.
The corresponding equations for the liquid drop and surrounding liquid will have the
form [4]
@Td 1 @ @Td
Cd � ldr2
d
@t r2 @r @r
Wd � 0; � const; 04r4btą
d
2ą
@Te @Te 1 @ @Te
Ce � We � ler2
e
@t @r r2 @r @r
Wer2 � Weaa2; � const; atą4r
e
Here Wea is the velocity of the liquid at the bubble surface �Wea � Weaą�. The
boundary conditions for Eqs. (1) and (2) at the centre of the drop, on the moving
boundaries btą and atą and at in�nity are [4]
@Td
r � 0 : � 0
@r
8
: :
>
> � Wgb btą � btą
g d
>
<
Tgb � Tdb � Tsg Pg
r � btą :
>
> @Tg @Td
>
:
lg ld � l
@r @r
N.S. Khabeev / Applied Energy 64 (1999) 317ą321 319
8
@Tg @Te
>
>
lg � le
>
<
@r @r
Tga � Tea
r � atą : 3ą
:
>
>
Wga � a � Wea
>
:
Pea � Pg 2 =a
r ! 1Te ! T1 � const; P ! P1 � const
Here ` is the latent heat of evapouration of the drop, is the of phase transition
rate per unit interfacial surface ( >0 for evapouration). TsgPą is the saturation
temperature; is the surface tension.
Being in the saturated state at the interface, the vapour obeys the Clapeyron-
Clausius equation
dTsg Tsg gs
� 1 4ą
dPg ` gs d
The Rayleigh equation for bubble oscillations in incompressible Żuid has a form [5]
dWea 3 Pg P1 2 =a
a � W2 � 5ą
ea
dt 2
e
When the pressure uniformity condition in the gas is satis�ed, there is an integral
of the heat equation for the gas which can be obtained by substitution from the
continuity equation to the heat transfer equation and using the equation of state [4]
a
dPg 3 @Tg a 3 Pg
� lgTąr2 r2Wg ; 6ą
dt a3 b3 @r a3 b3 b
b
frąja � faą fbą
b
The continuity equation for the gas with the use of the pressure uniformity con-
dition and the boundary conditions, yields the velocity pro�le in the gas:
b2 1 @Tg r r3 b3 dPg
Wg � Wgb � lgTąr2 7ą
r2 r2Pg @r 3 Pgr2 dt
b
So, the system of basic equations consists now of three non-linear equations of
convectional heat conductivity in di�erent space regions (in the drop, gas, and
liquid), four ordinary di�erential equations [Rayleigh equation (5), Clapeyroną
Clausius equation (4), equation for pressure (6), and equation for radius of the drop
(3), and the boundary conditions (3)].
The considered problem is characterized by considerable variations of the vapour
thermophysical parameters along the vapour layer. For example, for the system, sea
320 N.S. Khabeev / Applied Energy 64 (1999) 317ą321
water (300 K)ąliqui�ed hydrogen (Ts � 22 K for P=1 atm) or sea water (300 K) and
liqui�ed helium (Ts=4.2 K for p=1 atm), the thermal di�usivity coe�cient Dg �
lg= cg changes by a factor of 40.
g
The high-temperature gradient in the thin vapour-layer Te=Ts >> 1ą leads to the
strong variation of vapour density � Pg=RgTg and thermal conductivity along
g
the vapour layer. That is why it is very important to take into account the depen-
dence of the thermal conductivity on the temperature in the heat transfer equation
for the gas phase. The other thermophysical parameters, namely ld; ;
d
Cd; Cg; `; le; and Ce can be considered as constants.
e
3. Discussion of the results
The problem was solved numerically using a �nite-di�erence technique by dividing
the whole system into spherical layers inside the drop and bubble and outside the
bubble. There are two moving boundaries btą and atą. We used new variables for
``freezing'' the moving boundaries, namely
r
r 2 o; btąą : � ; 2 0; 1ą
btą
8ą
r b
r 2 btą; 1ą : � ; 2 0; 1ą
a b
So the boundaries now are not moving and �xed at the points
r � 0 $ � 0; r � btą $ � 1; � 0
r � atą $ � 1; r !1 $ ! 1
With allowance made for the �nite thermal-conductivity of the liquid, the bound-
ary condition at in�nity can be applied to the last layer of the liquid.
Calculations were made for di�ent sizes (0.1!1.0 mm) of drops. The surrounding
liquid was water at 293 K. The initial pressure in the system was 1 atmosphere. We used
di�erent liqui�ed gases at the saturation temperature. Thermophysical parameters
were taken from reference [6].
Calculations made for di�erent initial temperature-distributions in the vapour indicate
that the initial temperature distribution has an e�ect only on the initial time-interval. That
is why we used as the initial condition, a linear temperature distribution in the vapour.
Calculations show that the variations of the saturation temperature at the drop's
surface and of the interface bubble surface-temperature are very small, i.e.,
TsPą TsPoą Tea Teao
<< 1; << 1
TsPoą Teao
Here the subscript 0 refers to the parameters as at the initial state.
So, it is possible to consider these temperatures as constants, and the surrounding
liquid as a thermostat. Hence it is possible to simplify the problem; i.e. do not solve the
heat equations inside the drop and surrounding liquid, but solve them in the gas only.
N.S. Khabeev / Applied Energy 64 (1999) 317ą321 321
Analysis of the evapouration process shows that, because << , the char-
g d
acteristic time of the drop's radius change due to evapouration is much greater than
the characteristic time of the heat processes. That is why, at the time intervals,
comparable with the characteristic time of the heat processes the drop's radius can
be considered as constant. This was con�rmed by calculation.
The calculations have shown that there are two stages of the process. First or
dynamic stage is characterized by intensive pressure-oscillations in the gas during
the time interval b2=Dg. At this stage, the pressure in the bubble can be con-
0
siderably greater than initial atmospheric-pressure. The second or thermal stage is
characterized by the monotonic growth of the vapour bubble; the pressure in the
vapour is settled and equal to the external pressure Pg � P1; the temperature pro�le
approaching the quasi-stationary temperature-distribution. The initial value of the
vapour-layer thickness has a signi�cant e�ect on the process at the dynamic stage.
The calculations have shown that a decrease of the initial vapour-layer thickness
leads to an increase of the pressure-oscillations amplitude and the oscillations attenua-
tion time. The decrease of the initial vapour-layer size leads also to a sharp increase of
:
the evapouration rate and the velocity a of the bubble surface. This is connected with
the growth of the temperature gradient in the vapour phase due to the decrease of the
vapour-layer thickness. The calculations have shown that, at the �rst dynamic stage of
the process, there is a slight change of the initial linear temperature-distribution in the
gas phase, but a strong change of the velocity �eld. At the second thermal stage of the
process, the velocity pro�les approach quasi-steady con�gurations.
4. Conclusion
A mathematical model describing the dynamics of a bubble containing the drop of
liqui�ed gas, with allowances for the non-linear non-steady interphase heat and
mass transfers, has been proposed. General regularities of the behaviour of such a
``two-phase'' bubble, depending on the inital conditions and thermophysical prop-
erties of the phases have been studied. The proposed model can be useful for the
analysis and understanding of the complicated phenomenon of a vapour explosion.
References
[1] Banko� SG. Vapour explosions: a critical review, In: Proc. of the 6th int. heat-transfer conf., vol. 6,
Toronto, Canada. Toronto: National Research Council of Canada, 1978. p. 355ą66.
[2] Cronenberg AW. Recent developments in the understanding of energetic molten fuelącoolant inter-
actions. Nuclear Safety 1980;21(3):319ą37.
[5] Feng ZC, Leal LG. Non-linear bubble dynamics. Ann Rev Fluid Mech 1997;29:201ą43.
[3] Nigmatulin RI, Khabeev NS, Nagiev FB. Dynamics, heat and mass transfer of vapour-gas bubbles
in liquid. Int J Heat Mass Trans 1981;24(6):1033ą44.
[4] Nigmatulin RI. Dynamics of multiąphase media. Washington, DC: Hemisphere, 1990.
[6] Vargaftik NB. Handbook of thermophysical properties of gases and liquids. (in Russian) Nauka:
Moscow, 1972.