- Email: [email protected]
A simplified model for real gas expansion between two reservoirs connected by a thin tube
Pergamon
Chemical Engineerin 0 Science, Vol. 51, No. 2, pp. 295-308, 1996 Copyright © 1995 Elsevier Science Ltd Printed in Great Britain. All fights reserved 0009 2509/96 $9.50 + 0.00
0009-2509(95)00256-1
A S I M P L I F I E D M O D E L FOR REAL GAS E X P A N S I O N B E T W E E N TWO RESERVOIRS C O N N E C T E D BY A THIN TUBE S. CHARTON ~, V. BLET~ and J. P. CORRIOU t* ~LSGC-CNRS-ENSIC, BP451, 54001 Nancy Cedex, France; ~CEA, BP12, 91680 Bruy~res-le-Chfitel, France (First received 27 December 1994; accepted in revised form 19 July 1995) Abstract--A simplified model is derived for simulating gas blowdown through a thin and long tube connecting a high-pressure gas-filledreservoir to a vacuum vessel. During the depressurization process, the flow is assumed to be quasi-steady and approximated as one-dimensional.The transient compressibleflow can therefore be solved analytically. The model accounts for gas nonideality and heat transfer with the ambient in both reservoirs, and for wall friction in the pipe. An experimental study has been carried out using helium and deuterium in order to validate this simplifiedapproach. The discharge of the pressure vessel being quite rapid (a few seconds long), the measurements reliability and the required adjustments are discussed. The model predictions are shown to be in good agreement with our experimental results.
INTRODUCTION The depressurization of large pressure vessels occurs in process plants either accidentally or voluntarily. Consequently, and for the sake of security, the system behavior must be understood. However, the related literature is not abundant. Concerning the vessels themselves, pressure-relief devices must be sized, with respect to the highest pressure admitted in the vessel, in order to ensure a sufficient exiting flow rate in case, for instance, of runaway reactions. Furthermore, blowdown time should be evaluated through a given orifice (Leung, 1986). Internal phenomena are also to be described, and particularly the thermal history: the wall temperature may indeed reach the ductile-brittle transition point of the metal. Assuming that the flow is steady and using a constant heat transfer coefficient between the vessel wall and the fluid, Xia et al. (1993) derived a simplified model to predict the temperature and the pressure evolutions during the discharge. A more complete computer model has been developed by Haque et al. (1992), which simulates the tank pressure and temperatures histories, as well as the transient flow rate, during a rapid depressurization process. On the other hand, Levenspiel (1977) studied the atmospheric discharge of an infinite reservoir through a pipe. Focusing on the phenomena induced by the pipe, he determined, by means of thermodynamic considerations, the adiabatic, frictional and steady flow rate of a perfect gas for various tube size and reservoir pressures. This paper is related to the simultaneous occurrence of these two extreme situations: the transfer of a real gas from a high-pressure vessel to a low-pressure one through a thin tube. It is then important to describe correctly the discharge of the container, as
*Corresponding author.
well as the simultaneous filling of the receiver, in terms of mass and temperature variations, since they both define the boundary conditions of the transient flow occurring in the pipe. A rigorous investigation, concerned with dynamic flow-patterns and local temperature fields, would imply the resolution of the mass, energy and momentum balances of the compressible fluid in the whole system. Many schemes are available for compressible pipe flow calculations (Fletcher, 1991) although not ensuring the correctness of the calculated results because of the system complexity. Most of these algorithms are suitable when inflow and outflow boundary conditions are specified. However, these latter are often difficult to define and the computational task would become tricky and time consuming in our configuration. Indeed, because of the abrupt constriction at the tube inlet, the use of a fine grid would be essential, at least locally. Furthermore, another difficulty arises due to the significant pressure gradient prevailing between the two enclosures: pressure waves are generated, which propagate in both directions of the tube, and reflect on the vessels' walls until a steady flow regime is achieved. The numerical difficulties encountered when handling such discontinuities can be overcome by means of the characteristics theory. Issa and Spalding (1972) derived a hybrid scheme which combines the use of a rectangular grid with the use of characteristics. It is therefore particulary efficient for shock-tube simulations. The object of this paper is to present a simplified way of modelling this specific transfer, while keeping a sufficient accuracy for chemical engineering purposes. In order to reach an adequate compromise between calculation accuracy and computation efficiency, the choice of the model assumptions is of prime importance, in order to account only for the most significant effects and to neglect secondary phenomena.
295
296
S. CHARTON e t
The derived model is inspired from the analytical procedure described by Levenspiel (1977) insofar as a steady, adiabatic and frictional flow of a perfect gas is considered at each time step in the tube. Contrary to this previously studied case, the flow system is such that temporal variations of pressure and temperature take place in both the discharging and downstream reservoirs. Hence, rigorous transient mass and energy balances, accounting for real gas behavior and heat transfer with the ambient, are implemented in the model for the reservoirs treatment. The Soave cubic equation of state is used to estimate compressibility factors and departure functions. An experimental study of the process has been carried out in order to test the model reliability. Two gases of close mole weight were chosen: helium and deuterium, respectively monoatomic and diatomic. Initial pressures ranging from 10 to 45 MPa were investigated for both a long (1.1 m) and a short (0.3 m) tube. As fast variations are involved in the studied process, compared to the dynamics of the sensors, measurement distortions are encountered. For this reason, prior to any confrontation of the model predictions with experimental data, the raw signals have been treated to account for the sensors characteristics. This experimental investigation validates the transfer-time prediction and reveals, in each case, a good agreement between predicted and measured profiles in the vessels. 2. T H E O R E T I C A L D E V E L O P M E N T
The model concerns the transient transfer of a real gas from a high-pressure reservoir to a low-pressure reservoir through a thin and long tube (Fig. 1). Assumptions made in the development of the model are listed in the following: 1. In the vessels:
(a) For the purpose of the simplified approach, temperature and pressure are assumed to be spatially uniform. So, homogeneous conditions are prevailing in both vessels. (b) Moreover, the vessels are large with respect to the tube section. Consequently, the gas velocity is neglected and stagnant conditions are admitted. (c) The treatment of heat losses is simplified by considering a constant wall temperature, equal to the
al.
temperature of the surroundings. Indeed, the transfer rapidity prevents the heavy metallic-shell temperature from being significantly affected by the thermal variations, however important they are, of the small amount of gas involved, in accordance with the total heat-balance: (mCe) . . . . ImTmetal = (mCe)gasATgas. (d) Like in industrial vessels discharge, we assume that natural convection prevails in the pressure container (Haque et al., 1992), since it is driven by the density gradients generated by the depressurization. (e) In the downstream reservoir, to be consistent with the stagnant conditions specified, we will also assume, according to the Boussinesq approximation, that the gas motion results only from the density gradients, and hence that natural convection is dominant. These two last assumptions will be discussed later. (f) The global heat transfer coefficient is calculated in both vessels by eq. (1), valid for isothermal surfaces (Cess, 1973):
{
-N'-uU = a l R a
Nu
°'25
aERa °'33
for Ra <~ 109 for Ra > 109
where the global Nusselt number and the Rayleigh number are respectively defined by S~ =
Qt°talfluxL 2S(Tw-
T~)
__
h~L 2
and R a = Gr Pr =
gfl(Tw - T ~ ) p 2 L 3 Cplt It 2
2
L is the characteristic length of the exchange surface, Tw is the wall temperature, T~ the bulk gas temperature, al is set to 0.47 for a cylinder of diameter L and 0.49 for a sphere, and a2 = 0.1 for a cylinder of length L and for a sphere. The dependence of relation (1) on the Rayleigh number allows to account for the onset of turbulent convection (still generated by the density gradients) in the reservoirs. Since the gas properties depend on temperature, the reference temperature defined by Sparrow and Gregg (Cess, 1973) is adopted to estimate C v, #, 2, and p: T = T w - 0.38(Tw - T®).
(2)
Only the dilatation factor fl = (1/p)(dp/aT)v is calculated at the bulk temperature.
® RECEIVER
® ;
® i I
i
K
x]
'~11 Hish-Pretsu~
(1)
gl Vacuum Conditiom
Fig. 1. Schematic illustration of the process.
A simplified model for real gas expansion
297
1.40
"O 1.30
~J
1.20
! t..
.= 1.10
°m
1.00 L
E 0.90 0.90
1.00
1.10
1.20
1.30
1.40
Compressibility Factor - Experimental Values
Fig. 2. Comparison of compressibility factors (Z = PV/RT) in the range 102 ~
2. In the tube: (a) The flow is quasi-steady. Indeed, as it was mentioned previously, and will be illustrated later, at the beginning of the transfer, a rarefaction wave proceeds upstream in the tube, while a compression or shock wave proceeds downstream, thus ensuring that pressure conditions consistent with the frictional flow are prevailing between the tube openings. Under such choked conditions, the transient flow turns to a quasisteady one. Thanks to this major assumption, the numerical task is greatly simplified. (b) Unlike the assumption taken inside the vessels, the ideal gas law is assumed, thus allowing the derivation of an analytical procedure. (c) Because of the small tube section, radial variations of the fluid properties are negligible compared to the longitudinal ones. Therefore, the flow is treated as one-dimensional. (d) Finally, the flow is supposed to be adiabatic. This is reasonable, first of all, because the gas resi-
dence time in the tube is very short (around 3 ms). Furthermore, at the tube wall, the gas velocity is almost zero; its temperature at the wall is therefore equal to the stagnation temperature T o of the flowing gas. Hence, as mentioned by Levenspiel (1977), the constant wall temperature in the tube is better represented by the adiabatic flow assumption: d T °/dz = O. (e) To set the inlet conditions, we assume that the fluid expands isentropically at the pressure-vessel oririce. The significant effect is indeed the mechanical work of the tube itself and not the pressure drop induced by the orifice. Based upon these assumptions, mass and energy balances in the reservoirs yield:
• In the pressure vessel:
dpi
V1 - -
dt
= - p2u2[2
dU~=dt - p 2 u 2 ( h 2 + ~ )
(3)
+hwSw(Tw- T1)" (4)
• In the downstream reservoir: do. V4 - ~ = p3u3~ d t = psU3 h3 +
+ hwSw(rw- r4)
(5) (6)
where U refers to the internal energy of the gas contained in the vessel and h2, h3 to the specific enthalpy of the gas respectively at the tube inlet and outlet. They are estimated from the corresponding ideal gas
298
S. CHARTON et al.
properties by means of the departure functions. Pressure, temperature and density are related by the real gas law: P = p Z ( R / M ) T where Z is given by the cubic equation of Soave (1980): Z 3 - - Z 2 "~ [A* -- B*(1 + B*)]Z - A ' B * = 0
extremities and leads to + 1 In ( J/2{1 + [(7 - 1)/2] J/2}'~ 2
\M/32 {1 + [(y
(7) -
=
with A* = aP/R2T 2 and B* = bP/RT. • In the tube, the flow is described by the following mass, momentum and energy conservation equations d(pu) = 0 dx
(8)
d(pu u) dP 2Cf 2 dx + -~x + ---ff- pu = 0 d--x pu h + ~ -
P3 P2
---if-
(15)
L
J[2 N/~ + [(Y -- 1)/2] d/2 J43 + [(y 1)/2] j [ 2 -
(16)
The inlet boundary conditions have been previously defined by assuming an isentropic expansion of the gas leaving the container, that is
(9)
=0
1)/2] J/~}J
, P1 _ P2 ~7-=7 Pl P2
(10)
=
ht = 112 -{ us 2
and the ideal gas law R
P = p - ~ T.
(11)
In eq. (9), C: is the friction factor. It is commonly admitted that incompressible fluid dynamic correlations remain valid for subsonic flows of compressible gas (Shapiro, 1953). Churchill's correlation has been chosen for its accuracy in both laminar and turbulent flow regimes: C.r V ( 8 " ~ ~2 1 1 ~/2 -2- = [_\-~eJ + (A + B) 3/2
(12)
with
=
}
and (37,530~ 10 B = \----RT--e) " Re = puD/# is the flow Reynolds number, e is the pipe roughness. Introducing the dimensionless velocity ,A¢
U
For an ideal gas, enthalpy is only temperature dependent: dh = CvdT. It comes in terms of the Mach number: PI=P2
T , = T2(1 + ~ _ ~ ~ 2 ) .
~h
hO=hc+IC-2 2i U [
1 - d¢ 2 d~1¢ ~g2{ 1 + [ ( T - 1)/2] ~//a} #// 1 +(7+1).//2 -
(18)
Finally, eqs (15)-(18) and eq. (11) constitute a set of five equations involving six unknowns, Pa, Pa, T2, T3, ,g/2 and ~¢3, where subscripts refer to points indicated in Fig. 1. All but one can be deduced from P1 and Tt. Indeed, the determination of,//2 requires another set of boundary conditions since ~/¢2 depends both upon the upstream and downstream flow boundaries. The one-dimensional, adiabatic and frictional flow we are dealing with is known as Fanno flow (Shapiro, 1953). It is characterized by a constant stagnation enthalpy h ° = h + u=/2, and is usually illustrated in an enthalpy/entropy diagram by the set of constant flow rate curves. The upper branch of each Fanno curve refers to subsonic flows, while the lower one is only described by supersonic flows (Fig. 3).
eqs (8)-(11) yield
e
(17)
U
c
dP
7 - 1 j / 2 ) ~/(~-1) 1 + ----~-
2yCrdx D
"~"
(13)
de//
1 + [(~ - 1 ) / 2 ] ~ '~ ~
(14)
"
Equations (13) and (14) describe respectively the Mach number variations along the tube for the adiabatic and frictional flow, and the corresponding pressure variations with respect to ./¢. Their analytical integration can be derived between the two
C*
~s
Fig. 3. Fanno curves in the enthalpy/entropy diagram.
A simplified model for real gas expansion Since there are no divergent devices in the flow system of interest here, we are only concerned with subsonic flows: ~¢¢ ~< 1. Let F s represent the steady flow at the pipe inlet 2 at time t. As the gas progresses in the tube, its velocity is enhanced by frictional pressure drop and F S moves to the right on the upper branch of the Fanno line. If the tube is long enough, the sonic velocity associated to the critical point C* is reached: F S = 3 = C*, but cannot be passed beyond. Indeed, this transition from a subsonic to a supersonic flow regime, together with an entropy diminution, is actually a violation of the second law of thermodynamics. Instead, choking occurs and acts to decrease the flow rate until a steady-state solution again becomes possible with 3 = C*. From these thermodynamic considerations, it comes out that, regardless of the pressure ratio P1/P4, the highest velocity attainable is restricted to the speed of sound c = ~ ) T , and can only be reached at the end of the pipe. Thus, under choking conditions, the exhaust pressure is higher than the backpressure (the receiver pressure in this case). The missing set of boundary conditions is deduced from the foregoing discussion. Substituting eq. (17) in eq. (16) together with ~¢¢3 = 1 yields the gas pressure at the critical point C*: p*
/ 2 ~ g 2 { 1 + [(~' - 1)/2],/(2}~1+y)/.-~'~
~-1 = ~ /
y+l
(19)
Consequently, under choking conditions, the flow is such that P1
>(P3 =P*)>IP4
JCa = 1
Container
7 X
t
299
otherwise P1 > (P3 = P4)/> P* J//3 < 1.
Using these relations, the system of equations (15), (16) and (17) is solved iteratively in order to find J/¢2 at each time step. The other flow properties are then easily derived. 3. EXPERIMENTAL
3.1. Equipment The experimental apparatus (entirely made of stainless steel) is constituted by a set of two spherical instrumented reservoirs, connected by a thin tube (Fig. 4). One of them is filled at a relatively high pressure with the pure gas (container, C) and the second one is kept under vacuum (receiver, R). Their dimensions are respectively • container: internal diameter 0.062 m, thickness 0.03 m • receiver: internal diameter 0.197m, thickness 0.003 m. C is connected on one side to either a storage vessel or a vacuum pump, and is kept isolated from the tube by a computer-driven electro-valve on the opposite side. The tube is directly opened to R, the latter being separated by a manual valve from a precision gauge, on which the equilibrium pressure is read after transfer. The tube internal and external diameters are respectively 1 and 1.5 mm, with a roughness of 5 pro, typical of a drawn tube (Hodge et al., 1989). Absolute pressure transducers are used to give direct measurements in each reservoir. They respectively
Receiver
Precision gauge
~,e
aE
Vacuum
Fig. 4. Schematic illustration of the experimental apparatus.
300
S. CHARTON et al.
range from 0 to 4 and 0 to 60 MPa, with an accuracy o f l % of the full-scale length. Their output signals are linearized and amplified. The whole chain time constant is estimated as 0.1 s. Piezoelectric sensors connected to fast charge amplifiers are also used to quantify the rapidity of the pressure variations. Fast transient temperature measurements are more difficult to implement. The dynamic characteristic of a thermocouple is indeed strongly dependent on the heat transfer kinetics between the fluid and the junction. In recently published papers, time constants as low as 3.4 ms have been reached with microthermocouples (Beckman et al., 1993). In this study, we have used commercial chromel-alumel thermocouples. The 1 m m diameter isolation sheath is waisted to 0.5 m m over 2 cm from the external wires junction (in contact with the fluid). These features lead to a time constant of 0.4 s for gas-phase measurements. The output from the container sensor is sent through a fast amplifier and is converted afterwards. The receiver thermocopule signals are directly ampli-
fled and converted, providing a total time constant of the order of 0.8 s. Data are collected with an AD-2000 acquisition board. Several experimental conditions were investigated in order to test the model reliability: initial pressures ranging from 10 to 45 MPa, total tube length (from the container orifice to the receiver) of 0.3 and 1.1 m. Furthermore, since the ratio of the heat capacities appears to be an important parameter of the transonic flow involved, helium (7 = 1.67) and deuterium (y = 1.41) gases were chosen for their close mole weights. The experimental conditions are summarized in Table 1, where experiments I - V I I I refer to helium, and experiments I X - X V I to deuterium. 3.2. Data processin# For all the experiments, mass values in each reservoir result from calculations based on experimental data of temperature and pressure. The overall mass balance of the system which is deduced points out the influence of the sensors dynamics on the physical
Table 1. Experimental conditions investigated Gas nature I II III IV V VI VII VIII IX X XI XII XIII XIV XV XVI
He He He He He He He He D2 D2
D2 D2 D2 D2 D2 D2
Tube length L (m)
Container initial pressure (MPa)
Mean absolute pressure deviation (MPa)
Discharge duration ~90 (s)
1.1 1.1 1.1 1.1 0.3 0.3 0.3 0.3 1.1 1.1 1.1 1.1 0.3 0.3 0.3 0.3
10.0 22.0 31.9 45.0 9.9 22.0 31.8 43.1 9.6 20.4 31.0 42.4 10.7 20.1 31.0 40.8
0.452 0.162 0.954 0.208 0.397 0.482 0.844 0.849 0.236 0.410 0.421 0.512 0.967 0.963 1.350 0.926
2.30 2.15 2.00 1.95 1.25 1.15 1.10 1.05 2.55 2.45 2.40 2.30 1.40 1.35 1.32 1.30
18 16
i~
14
12 i0
~
\
6
2 0
0.5
1
1.5
2
2.5
3
3.5
4
t i m e (s) Fig. 5. Apparent mass fluctuation during experiment III.
4.5
A simplifiedmodel for real gas expansion meaning of the measurements. For instance, the apparent gas mass variation for the whole system, calculated for experiment III, is plotted in Fig. 5. From the experimental values of temperature and pressure, it would be concluded that the amount of gas enclosed in the hermetically sealed apparatus does not remain constant during the transfer. This apparent mass vari-
301
ation, noticeable in the first seconds of each experiment (and attaining sometimes 20% in the case of the short tube experiments), is of course inconsistent. It is attributed to the dynamics of the sensors. Obviously, the time scales of the responses of the sensors, especially thermocouples, are of the same order of magnitude as the time scales of the pressure
Pressure (MPa) 35 30
25
20
model 15
''~
",,~ •
10
'
~ •
corrected model ex )eriment
•
It--
0
0.5
1
1.5
2
2.5
3
3.5 4 time (s)
Fig. 6. Pressure evolution in the container during helium transfer (experiment III).
Pressure (MPa) 0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0.5
1
1.5
2
2.5
3
3.51m e (s)4
Fig. 7. Pressure evolution in the receiver during helium transfer (experiment III). CE$ 51-2-J
302
et
S. CHARTON
and temperature variations during the transfer. Consequently, the temporal characteristics of the sensors must be taken into account for a correct model validation. The time responses of the temperature and pressure measurement chains are represented by first-order transfer functions of gain unity, with respective time
al.
constants given earlier. U n d e r these conditions, model predictions can be corrected by convoluting the simulated profiles by the inverse of the sensor transfer function. Since a first amplifier is coupled to the container thermocouple, we found it more meaningful to proceed in the reverse order. The heat transfer between
Temperature (K) 313.15 293.15'
b"
2~.15
'.
253.15
•
:
model
':
233.15
D
c4~L~__~ measuremenls
-
-
r a w me~m~tn~t~
•
213.15
b.
9
•, •
193.15
O..•
.*"
"". 173.15
• "',D "- " ' .
.~.',t'" •
[] "''-,.
"DO
~h ...o-''"
[] .........
.
iC~.._- "
[] .'"~
•
-
153.15 133.15 0
0.5
1
1.5
2
2.5
3
3.5
4
time (s) Fig. 8. Temperature evolution in the container during helium transfer (experiment III).
Temperature (K)
,'°'1
" 423.15
model .
...... •
' ". 398.15
348.15 323.15
corrected model experiment
~b ah~w~lqmu
/
,-°.. ......
i. . . . . . . .
273.15 0
0.5
1
1.5
2
2.5
3
3.5
4
time (s) Fig. 9. Temperature evolution in the receiver during helium transfer (experiment III).
A simplified model for real gas expansion the fluid ( f ) and the thermocouple junction (j) is indeed expressed by
303
Equation (20) can be rewritten as a simple firstorder differential equation, introducing the time constant • j of the thermocouple:
(20)
Tj = - Ajhj_f(Tj - Tf) ( p V C v ) j dd-T
dT~ = dt
where h~_[ is the convection heat transfer coefficient at the junction surface.
( T ~ - Tj-) zj
(21)
Pressure (MPa) 35
20
\
.......
•
model corrected model expenment
•
10
0
.... 0
0.5
1
1.5
2
2.5
- ---3
" 3.5
4 t i m e (s)
Fig. 10. Pressure evolution in the container during deuterium transfer (experiment XI). Pressure (MPa)
0.6
0.5
o., .,
|
•
;- ~/•
0.3
model 0.2
:
;"
corrected model experiment
•
q
0.1
0 0
0.5
I
1.5
2
2.5
3
3.5
4 t i m e (s)
Fig. l l. Pressure evolution in the receiver during deuterium transfer (experiment XI).
304
S. CHARTON
The obtention of corrected experimental temperature T I is achieved by solving eq. (21) from the junction temperatures measurements in the container.
et al.
the temperature profiles in the pressure vessel, Figs 8 and 12, is characteristic of the competition between the two thermal phenomena occurring during the discharge:
4. RESULTS AND DISCUSSION
The results from each experiment have been compared with the corresponding predicted data from the model. The mean absolute deviation between the experimental and the corrected simulated pressure profiles has been estimated in each case (Table 1). Those small departures are mainly attributed to the transducer uncertainty which is about 0.6 MPa. Two representative comparisons, III and XI (referred to as A and B), are discussed in the following. The pressure variations with respect to time in both reservoirs are respectively illustrated in Figs 6, 10 and 7, 11. The raw evolutions predicted by the model have been corrected by introduction of the sensor dynamics as indicated earlier, and the resulting curves, labelled as corrected model, can thus be compared with the measured profiles. The same convention is used for the corresponding temperature histories, Figs 8, 12 and 9, 13, except for the container, where the raw thermocouple output signals have been corrected according to eq. (21). Therefore, in Figs 8 and 12, the simulated temperature profiles are compared with the real temperature evolutions in the container. After the valve opening, the gas rapidly expands from the container to the receiver. The pressure balancing goes on afterwards, together with the achievement of the thermal equilibrium. The shape of
• the rapid cooling, due to the important depressurization, strongly dominates in the first second, • after a while, the endothermic process becoming less significant, the gas is slowly warmed up by the surroundings. In the receiver, conversely, a temperature increase resulting from both compression and kinetic energy conversion into internal energy is first observed. Then, these exothermic processes are counterbalanced by the heat losses through the metallic shell (Figs 9 and 13). This qualitative behavior of the system appears to be correctly described by our simplified model. Furthermore, there is clearly a very good agreement between the predicted profiles and the measured ones in the container: • The departures on pressure are of the same order of magnitude as the transducer uncertainty (Figs 6 and 10); • the mean bulk temperatures predicted in this vessel are very close to the corrected experimental values (Figs 8 and 12), the latter being measured locally in the vessel.
T e m p e r a t u r e (K)
293.15 ~ 273.15 - - Q
• • - -
253.15
r-~
233.15
[]
• []
"
model
corrected measurements
raw measurements
"1% %.
213.15
~...
•
[]
193.15
°. o~I
6 .... b'-'-.
• 0...<2...~-~-"
:~
i1
•
173.15 153.15 133.15 0
0.5
1
1.5
2
2.5
3
3.5
4 t i m e (s)
Fig. 12. Temperature evolution in the container during deuterium transfer (experiment XI).
A simplified model for real gas expansion
305
Temperature (K) 423.15
398.15
373.15 ~
348,15
3
2
0
3
0.5
.
I
1
1.5
5
2
~
2.5
3
3.5
4
time (s) Fig. 13. Temperature evolution in the receiver during deuterium transfer (experiment XI). In the receiver, however, the model obviously overestimates the thermal variations: a sharp temperature increase is predicted in the first milliseconds, which is hardly compensated by the convective heattransfer. The measurements on the other hand reveal a smoother variation and lower temperature values (Figs 9 and 13). Too high predicted pressures follow from this temperature departure as shown in Figs 7 and 11. The weaker agreement observed in the receiver may be due to some inconsistencies in the assumptions made during the model development. Among them, the less reliable are the following: • spatially uniform temperatures • stagnant conditions • natural convection dominance. Indeed, because the receiver has finite dimensions, the sonic gas jet is likely to be reflected by the walls and then deviated by the continuous incoming gas flow. Therefore, nonstagnant conditions are prevailing in the enclosure and prevent the fluid kinetic energy from being entirely and instantaneously restituted in thermal energy. Sharp velocity fluctuations are taking place, thus generating a complex flow pattern in this enclosure, and strongly promoting forced convection because of the resulting buoyancy. Conversely, the velocity of the gas at the pipe entrance is considerably less significant than that of the exhausting sonic jet; therefore the assumptions of gas stagnation and natural convection dominance in the container are justified.
To take into account the influence of forced convection in the receiver, a different correlation for Nusselt number Nu = 0.023 Re °8 Pr °aa
(22)
has been introduced in the model instead of eq. (1). The Reynolds number has been given empirically, as there is no simple way to obtain a mean value of it. Two values of Re were used: 104 and 105 . When Re increases, the temperature decrease is stronger and the curves obtained with this model (Figs 14 and 15) are situated closer to the experimental curves and well under the natural convection curve. Unfortunately, the main problem with the model including forced convection is that the Re number cannot be predicted in the frame of a simplified approach; consequently this model cannot be considered as totally predictive. The particularly important discrepancy between raw experimental and corrected simulated temperatures may also be due to the significant time constant (around 0.8 s) of the corresponding measurement chain. The precision of the correction brought by the convolution process in this case may be insufficient. Lastly, the temperature field in the downstream reservoir is obviously very complex and strongly dependent on time (through the relative positions of the incident and reflected gas jets). It must be kept in mind that in this work we compare a mean value of the temperature (uniform temperature predicted by the simplified model) to a local in situ measurement. For the reasons evoked above, this procedure is appropriated in the discharging vessel, but the comparison
306
S. CHARTON et al.
Pressure (MPa) 0.7
0.6
0.5
'°
°°r
~
0.4
°t
..o~o
.dlo
, °, 4~°'°
j°;
0.3
0.2
f
--
."1
0,j
-.
•
natural convection forced convection - Re=le4 forced convection - Re=le5 ex ~eriment
0 0
0.5
1
1.5
2
2.5
3
3.5 4 t i m e (s)
Fig. 14. Pressure evolution in the receiver during helium transfer (experiment III)--influence of forced convection heat transfer.
Temperature (K) 473.15
448.15 423.15
natural convection forced convection - Re=le4 forced convection - Re=le5 experiment
......
39& 15
•
373.15
323°|5
f ~
•
!1,
/6
It
/,;,~
41.o
"
.... a ....... !; ...... o ...... o
o./,
lw
•
•
•
, - - I ~ a m
oo •
•
273.15
0
0.5
l
1.5
2
2.5
3
3.5
4 t i m e (s)
Fig. 15. Temperature evolution in the receiver during helium transfer (experiment III)--influence of forced convection heat transfer.
is n o t so representative in the receiver. However, in the receiver, the pressure evolution is correctly described by the model, a n d the transfer d u r a t i o n as well as the final pressure are accurately predicted (Figs 2 a n d 7).
Since the model simulates perfectly the p h e n o m e n a occurring in the container, the predicted discharge d u r a t i o n Z9o, defined as the time required to decrease the initially stored mass of gas of 90%, h a s been plotted with respect to the initial pressure, in o r d e r to
A simplified m o d e l for real gas e x p a n s i o n
307
Discharge Duration "c~(s) 3 • •
He- L=l.lm D2 - L = I . I m
2.5 -0 ................................... ....... 0
---..~.
""A
1.5 P. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
/
-0
,.o
.......A-......................... ............ .A
0.5
t,
He - L = 0 . 3 m
0
D2 - L = 0 . 3 m
0
10
30 40 50 Container Initial Pressure (MPa)
20
Fig. 16. V a r i a t i o n s of the d i s c h a r g e d u r a t i o n Zgo with the initial s t o r a g e pressure.
illustrate the properties of the frictional flow (Fig. 16). For a given gas and a given tube length, %0 decreases as the storage pressure increases. The latter is indeed the driving force of the discharge. It appears in Fig. 16 that %0 is also strongly dependent on the tube length and the gas molarity: (a) According to the Fanno flow characteristics, the transfer time is nearly doubled when the tube length is increased from 0.3 to 1.1 m (all other things being equal). The longer tube inducing a higher frictional pressure drop and hence a higher acceleration of the subsonic flow, the thermodynamic consistency requires the inlet gas flow rate to be reduced in a higher proportion for this long tube than for the shorter one. (b) For a given initial pressure and a fixed tube length, the gas transfer between the two reservoirs is shown to be slower in the case of deuterium than when helium is concerned. Keeping in mind that most of the gas is blown out under sonic conditions, the gas velocity at the pipe outlet is equal to the sonic velocity c = x / - ~ ( R / M ) T . Thus, for a given temperature, the deuterium over helium velocity ratio is given by the square root of their ?-ratio. Figure 16 reflects that the mean value of the helium over deuterium discharge duration is effectively close to the reciprocal of this ratio:
~
= 0.92
while
= 0.89.
a thin tube. The flow is assumed to be one-dimensional and quasi-steady. Real gas effects and heat transfer phenomena are neglected in the tube but taken into account in the reservoirs. Finally, providing sensible assumptions, the fast transient flow of the compressible fluid can be solved with minimal calculation efforts, simply by setting the system geometry and the fluid properties. Besides the apparent difficulty of the simulation task, the experimental investigation of the process appeared to be critical. Owing to the fact that the major amount of gas is blown in less than half a second, the measurements accuracy is strongly dependent on the dynamic characteristics of the sensors involved. The actual pressure and temperature evolutions can be estimated anyway by taking into account the transfer function of the measurement chains. The convolution product of the predicted data by the inverse of this function leads to a great improvement in the meaning of the comparison between simulations results and experimental data. The model predictions have been compared with several experimental data. The qualitative behavior of the process appears to be very well described and a very close quantitative agreement is reached in the discharging vessel. A weaker agreement is reached in the downstream reservoir probably because of the stagnant bulk approach derived in the model. However, the receiver conditions have a feeble influence on both the container discharge and the choked flow in the pipe which are correctly described.
5. C O N C L U S I O N
A simplified model has been proposed for simulating gas transfer between two enclosures through
NOTATION
c C/
velocity of sound, m/s friction factor
S. CHARTONet al.
308 Cp, C~
D # h hw L M P R
S t
T u
U V x Z
specific heat, J/(kg K) tube diameter, m gravity, m/s 2 specific enthalpy, J/kg heat transfer coefficient, W/(m 2 K) tube length, m mole weight, kg/mol Mach n u m b e r absolute pressure, Pa universal ideal gas constant, J/(K mol) area of the heat exchange surface, m 2 time, s absolute temperature, K velocity, m/s internal energy, J volume, m 3 abscissa compressibility factor
Greek letters isobaric dilatation coefficient, K - 1 ratio of specific heats 7 tube roughness, m 8 2 thermal conductivity, W / m /z viscosity, kg/(m s) density, kg/m a P characteristic time tube section, m 2 fZ Subscripts f fluid j junction w wall oo bulk Exponents 0 stagnation property
S
choking condition steady condition REFERENCES
Beckman, P., Roy, R. P., Whitfield, K. and Hasan, A., 1993, A fast-response microthermocouple. Rev. Sci. lnstrum. 64, 2947-2951. Cess, R. D., 1973, Free-convection boundary-layer heat transfer, in Handbook of Heat Transfer (Edited by W. M. Rohsenow and J. P. Harnett), Section 6, p. 613. McGrawHill, New York. Fletcher, C. A. J., 1991, Computational Techniques for Fluid Dynamics, Vol. 2: Specific Techniques for Different Flow Cate#ories, 2nd Edition. Springer, Berlin. Haque, M. A., Richardson, S. M. and Saville, G., 1992, Blowdown of pressure vessels. I. Computer model. Trans. Instn Chem. Engrs 70, 3-9. Hodge, B. K., Taylor, R. P. and Coleman, H. W., 1989, Predicting turbulent rough-wall skin friction and heat transfer, in Encyclopedia of Fluid Mechanics, Vol. 8: Aerodynamics and Compressible Flows (Edited by N. P. Cheremisinof), Chap. 13. Gulf Publishing Company, Houston. Issa, R. I. and Spalding, D. B., 1972, Unsteady one-dimensional compressible frictional flow with heat transfer. J. Mech. Engng Sci. 14, 365-369. Leung, J. C., 1986, Simplified vent sizing equations for emergency relief requirements in reactors and storage vessels. A.I.Ch.E.J. 32, 1622-1634. Levenspiel, O., 1977, The discharge of gases from a reservoir through a pipe. A.I.Ch.E.J. 23, 402-403. Reid, C., Prausnitz, J. M. and Poling, B. E., 1987, The Properties of Gases and Liquids, 4th Edition. McGrawHill, New York. Shapiro, A. H., 1953, The Dynamics and Thermodynamics of Compressible Fluid Flow, Vol. I. The Ronald Press Company, New York. Soave, G., 1980, Rigorous and simplified procedures for determining the pure-component parameters in the Redlich-Kwong-Soave equation of state. Chem. En#ng Sci. 35, 1725-1734. Xia, J. L., Smith, B. L. and Yadigaroglu, G., 1993, A simplified model for depressurization of gas-filled pressure vessels. Int. Comm. Heat Mass Transfer 20, 653-664.
Chemical Engineerin 0 Science, Vol. 51, No. 2, pp. 295-308, 1996 Copyright © 1995 Elsevier Science Ltd Printed in Great Britain. All fights reserved 0009 2509/96 $9.50 + 0.00
0009-2509(95)00256-1
A S I M P L I F I E D M O D E L FOR REAL GAS E X P A N S I O N B E T W E E N TWO RESERVOIRS C O N N E C T E D BY A THIN TUBE S. CHARTON ~, V. BLET~ and J. P. CORRIOU t* ~LSGC-CNRS-ENSIC, BP451, 54001 Nancy Cedex, France; ~CEA, BP12, 91680 Bruy~res-le-Chfitel, France (First received 27 December 1994; accepted in revised form 19 July 1995) Abstract--A simplified model is derived for simulating gas blowdown through a thin and long tube connecting a high-pressure gas-filledreservoir to a vacuum vessel. During the depressurization process, the flow is assumed to be quasi-steady and approximated as one-dimensional.The transient compressibleflow can therefore be solved analytically. The model accounts for gas nonideality and heat transfer with the ambient in both reservoirs, and for wall friction in the pipe. An experimental study has been carried out using helium and deuterium in order to validate this simplifiedapproach. The discharge of the pressure vessel being quite rapid (a few seconds long), the measurements reliability and the required adjustments are discussed. The model predictions are shown to be in good agreement with our experimental results.
INTRODUCTION The depressurization of large pressure vessels occurs in process plants either accidentally or voluntarily. Consequently, and for the sake of security, the system behavior must be understood. However, the related literature is not abundant. Concerning the vessels themselves, pressure-relief devices must be sized, with respect to the highest pressure admitted in the vessel, in order to ensure a sufficient exiting flow rate in case, for instance, of runaway reactions. Furthermore, blowdown time should be evaluated through a given orifice (Leung, 1986). Internal phenomena are also to be described, and particularly the thermal history: the wall temperature may indeed reach the ductile-brittle transition point of the metal. Assuming that the flow is steady and using a constant heat transfer coefficient between the vessel wall and the fluid, Xia et al. (1993) derived a simplified model to predict the temperature and the pressure evolutions during the discharge. A more complete computer model has been developed by Haque et al. (1992), which simulates the tank pressure and temperatures histories, as well as the transient flow rate, during a rapid depressurization process. On the other hand, Levenspiel (1977) studied the atmospheric discharge of an infinite reservoir through a pipe. Focusing on the phenomena induced by the pipe, he determined, by means of thermodynamic considerations, the adiabatic, frictional and steady flow rate of a perfect gas for various tube size and reservoir pressures. This paper is related to the simultaneous occurrence of these two extreme situations: the transfer of a real gas from a high-pressure vessel to a low-pressure one through a thin tube. It is then important to describe correctly the discharge of the container, as
*Corresponding author.
well as the simultaneous filling of the receiver, in terms of mass and temperature variations, since they both define the boundary conditions of the transient flow occurring in the pipe. A rigorous investigation, concerned with dynamic flow-patterns and local temperature fields, would imply the resolution of the mass, energy and momentum balances of the compressible fluid in the whole system. Many schemes are available for compressible pipe flow calculations (Fletcher, 1991) although not ensuring the correctness of the calculated results because of the system complexity. Most of these algorithms are suitable when inflow and outflow boundary conditions are specified. However, these latter are often difficult to define and the computational task would become tricky and time consuming in our configuration. Indeed, because of the abrupt constriction at the tube inlet, the use of a fine grid would be essential, at least locally. Furthermore, another difficulty arises due to the significant pressure gradient prevailing between the two enclosures: pressure waves are generated, which propagate in both directions of the tube, and reflect on the vessels' walls until a steady flow regime is achieved. The numerical difficulties encountered when handling such discontinuities can be overcome by means of the characteristics theory. Issa and Spalding (1972) derived a hybrid scheme which combines the use of a rectangular grid with the use of characteristics. It is therefore particulary efficient for shock-tube simulations. The object of this paper is to present a simplified way of modelling this specific transfer, while keeping a sufficient accuracy for chemical engineering purposes. In order to reach an adequate compromise between calculation accuracy and computation efficiency, the choice of the model assumptions is of prime importance, in order to account only for the most significant effects and to neglect secondary phenomena.
295
296
S. CHARTON e t
The derived model is inspired from the analytical procedure described by Levenspiel (1977) insofar as a steady, adiabatic and frictional flow of a perfect gas is considered at each time step in the tube. Contrary to this previously studied case, the flow system is such that temporal variations of pressure and temperature take place in both the discharging and downstream reservoirs. Hence, rigorous transient mass and energy balances, accounting for real gas behavior and heat transfer with the ambient, are implemented in the model for the reservoirs treatment. The Soave cubic equation of state is used to estimate compressibility factors and departure functions. An experimental study of the process has been carried out in order to test the model reliability. Two gases of close mole weight were chosen: helium and deuterium, respectively monoatomic and diatomic. Initial pressures ranging from 10 to 45 MPa were investigated for both a long (1.1 m) and a short (0.3 m) tube. As fast variations are involved in the studied process, compared to the dynamics of the sensors, measurement distortions are encountered. For this reason, prior to any confrontation of the model predictions with experimental data, the raw signals have been treated to account for the sensors characteristics. This experimental investigation validates the transfer-time prediction and reveals, in each case, a good agreement between predicted and measured profiles in the vessels. 2. T H E O R E T I C A L D E V E L O P M E N T
The model concerns the transient transfer of a real gas from a high-pressure reservoir to a low-pressure reservoir through a thin and long tube (Fig. 1). Assumptions made in the development of the model are listed in the following: 1. In the vessels:
(a) For the purpose of the simplified approach, temperature and pressure are assumed to be spatially uniform. So, homogeneous conditions are prevailing in both vessels. (b) Moreover, the vessels are large with respect to the tube section. Consequently, the gas velocity is neglected and stagnant conditions are admitted. (c) The treatment of heat losses is simplified by considering a constant wall temperature, equal to the
al.
temperature of the surroundings. Indeed, the transfer rapidity prevents the heavy metallic-shell temperature from being significantly affected by the thermal variations, however important they are, of the small amount of gas involved, in accordance with the total heat-balance: (mCe) . . . . ImTmetal = (mCe)gasATgas. (d) Like in industrial vessels discharge, we assume that natural convection prevails in the pressure container (Haque et al., 1992), since it is driven by the density gradients generated by the depressurization. (e) In the downstream reservoir, to be consistent with the stagnant conditions specified, we will also assume, according to the Boussinesq approximation, that the gas motion results only from the density gradients, and hence that natural convection is dominant. These two last assumptions will be discussed later. (f) The global heat transfer coefficient is calculated in both vessels by eq. (1), valid for isothermal surfaces (Cess, 1973):
{
-N'-uU = a l R a
Nu
°'25
aERa °'33
for Ra <~ 109 for Ra > 109
where the global Nusselt number and the Rayleigh number are respectively defined by S~ =
Qt°talfluxL 2S(Tw-
T~)
__
h~L 2
and R a = Gr Pr =
gfl(Tw - T ~ ) p 2 L 3 Cplt It 2
2
L is the characteristic length of the exchange surface, Tw is the wall temperature, T~ the bulk gas temperature, al is set to 0.47 for a cylinder of diameter L and 0.49 for a sphere, and a2 = 0.1 for a cylinder of length L and for a sphere. The dependence of relation (1) on the Rayleigh number allows to account for the onset of turbulent convection (still generated by the density gradients) in the reservoirs. Since the gas properties depend on temperature, the reference temperature defined by Sparrow and Gregg (Cess, 1973) is adopted to estimate C v, #, 2, and p: T = T w - 0.38(Tw - T®).
(2)
Only the dilatation factor fl = (1/p)(dp/aT)v is calculated at the bulk temperature.
® RECEIVER
® ;
® i I
i
K
x]
'~11 Hish-Pretsu~
(1)
gl Vacuum Conditiom
Fig. 1. Schematic illustration of the process.
A simplified model for real gas expansion
297
1.40
"O 1.30
~J
1.20
! t..
.= 1.10
°m
1.00 L
E 0.90 0.90
1.00
1.10
1.20
1.30
1.40
Compressibility Factor - Experimental Values
Fig. 2. Comparison of compressibility factors (Z = PV/RT) in the range 102 ~
2. In the tube: (a) The flow is quasi-steady. Indeed, as it was mentioned previously, and will be illustrated later, at the beginning of the transfer, a rarefaction wave proceeds upstream in the tube, while a compression or shock wave proceeds downstream, thus ensuring that pressure conditions consistent with the frictional flow are prevailing between the tube openings. Under such choked conditions, the transient flow turns to a quasisteady one. Thanks to this major assumption, the numerical task is greatly simplified. (b) Unlike the assumption taken inside the vessels, the ideal gas law is assumed, thus allowing the derivation of an analytical procedure. (c) Because of the small tube section, radial variations of the fluid properties are negligible compared to the longitudinal ones. Therefore, the flow is treated as one-dimensional. (d) Finally, the flow is supposed to be adiabatic. This is reasonable, first of all, because the gas resi-
dence time in the tube is very short (around 3 ms). Furthermore, at the tube wall, the gas velocity is almost zero; its temperature at the wall is therefore equal to the stagnation temperature T o of the flowing gas. Hence, as mentioned by Levenspiel (1977), the constant wall temperature in the tube is better represented by the adiabatic flow assumption: d T °/dz = O. (e) To set the inlet conditions, we assume that the fluid expands isentropically at the pressure-vessel oririce. The significant effect is indeed the mechanical work of the tube itself and not the pressure drop induced by the orifice. Based upon these assumptions, mass and energy balances in the reservoirs yield:
• In the pressure vessel:
dpi
V1 - -
dt
= - p2u2[2
dU~=dt - p 2 u 2 ( h 2 + ~ )
(3)
+hwSw(Tw- T1)" (4)
• In the downstream reservoir: do. V4 - ~ = p3u3~ d t = psU3 h3 +
+ hwSw(rw- r4)
(5) (6)
where U refers to the internal energy of the gas contained in the vessel and h2, h3 to the specific enthalpy of the gas respectively at the tube inlet and outlet. They are estimated from the corresponding ideal gas
298
S. CHARTON et al.
properties by means of the departure functions. Pressure, temperature and density are related by the real gas law: P = p Z ( R / M ) T where Z is given by the cubic equation of Soave (1980): Z 3 - - Z 2 "~ [A* -- B*(1 + B*)]Z - A ' B * = 0
extremities and leads to + 1 In ( J/2{1 + [(7 - 1)/2] J/2}'~ 2
\M/32 {1 + [(y
(7) -
=
with A* = aP/R2T 2 and B* = bP/RT. • In the tube, the flow is described by the following mass, momentum and energy conservation equations d(pu) = 0 dx
(8)
d(pu u) dP 2Cf 2 dx + -~x + ---ff- pu = 0 d--x pu h + ~ -
P3 P2
---if-
(15)
L
J[2 N/~ + [(Y -- 1)/2] d/2 J43 + [(y 1)/2] j [ 2 -
(16)
The inlet boundary conditions have been previously defined by assuming an isentropic expansion of the gas leaving the container, that is
(9)
=0
1)/2] J/~}J
, P1 _ P2 ~7-=7 Pl P2
(10)
=
ht = 112 -{ us 2
and the ideal gas law R
P = p - ~ T.
(11)
In eq. (9), C: is the friction factor. It is commonly admitted that incompressible fluid dynamic correlations remain valid for subsonic flows of compressible gas (Shapiro, 1953). Churchill's correlation has been chosen for its accuracy in both laminar and turbulent flow regimes: C.r V ( 8 " ~ ~2 1 1 ~/2 -2- = [_\-~eJ + (A + B) 3/2
(12)
with
=
}
and (37,530~ 10 B = \----RT--e) " Re = puD/# is the flow Reynolds number, e is the pipe roughness. Introducing the dimensionless velocity ,A¢
U
For an ideal gas, enthalpy is only temperature dependent: dh = CvdT. It comes in terms of the Mach number: PI=P2
T , = T2(1 + ~ _ ~ ~ 2 ) .
~h
hO=hc+IC-2 2i U [
1 - d¢ 2 d~1¢ ~g2{ 1 + [ ( T - 1)/2] ~//a} #// 1 +(7+1).//2 -
(18)
Finally, eqs (15)-(18) and eq. (11) constitute a set of five equations involving six unknowns, Pa, Pa, T2, T3, ,g/2 and ~¢3, where subscripts refer to points indicated in Fig. 1. All but one can be deduced from P1 and Tt. Indeed, the determination of,//2 requires another set of boundary conditions since ~/¢2 depends both upon the upstream and downstream flow boundaries. The one-dimensional, adiabatic and frictional flow we are dealing with is known as Fanno flow (Shapiro, 1953). It is characterized by a constant stagnation enthalpy h ° = h + u=/2, and is usually illustrated in an enthalpy/entropy diagram by the set of constant flow rate curves. The upper branch of each Fanno curve refers to subsonic flows, while the lower one is only described by supersonic flows (Fig. 3).
eqs (8)-(11) yield
e
(17)
U
c
dP
7 - 1 j / 2 ) ~/(~-1) 1 + ----~-
2yCrdx D
"~"
(13)
de//
1 + [(~ - 1 ) / 2 ] ~ '~ ~
(14)
"
Equations (13) and (14) describe respectively the Mach number variations along the tube for the adiabatic and frictional flow, and the corresponding pressure variations with respect to ./¢. Their analytical integration can be derived between the two
C*
~s
Fig. 3. Fanno curves in the enthalpy/entropy diagram.
A simplified model for real gas expansion Since there are no divergent devices in the flow system of interest here, we are only concerned with subsonic flows: ~¢¢ ~< 1. Let F s represent the steady flow at the pipe inlet 2 at time t. As the gas progresses in the tube, its velocity is enhanced by frictional pressure drop and F S moves to the right on the upper branch of the Fanno line. If the tube is long enough, the sonic velocity associated to the critical point C* is reached: F S = 3 = C*, but cannot be passed beyond. Indeed, this transition from a subsonic to a supersonic flow regime, together with an entropy diminution, is actually a violation of the second law of thermodynamics. Instead, choking occurs and acts to decrease the flow rate until a steady-state solution again becomes possible with 3 = C*. From these thermodynamic considerations, it comes out that, regardless of the pressure ratio P1/P4, the highest velocity attainable is restricted to the speed of sound c = ~ ) T , and can only be reached at the end of the pipe. Thus, under choking conditions, the exhaust pressure is higher than the backpressure (the receiver pressure in this case). The missing set of boundary conditions is deduced from the foregoing discussion. Substituting eq. (17) in eq. (16) together with ~¢¢3 = 1 yields the gas pressure at the critical point C*: p*
/ 2 ~ g 2 { 1 + [(~' - 1)/2],/(2}~1+y)/.-~'~
~-1 = ~ /
y+l
(19)
Consequently, under choking conditions, the flow is such that P1
>(P3 =P*)>IP4
JCa = 1
Container
7 X
t
299
otherwise P1 > (P3 = P4)/> P* J//3 < 1.
Using these relations, the system of equations (15), (16) and (17) is solved iteratively in order to find J/¢2 at each time step. The other flow properties are then easily derived. 3. EXPERIMENTAL
3.1. Equipment The experimental apparatus (entirely made of stainless steel) is constituted by a set of two spherical instrumented reservoirs, connected by a thin tube (Fig. 4). One of them is filled at a relatively high pressure with the pure gas (container, C) and the second one is kept under vacuum (receiver, R). Their dimensions are respectively • container: internal diameter 0.062 m, thickness 0.03 m • receiver: internal diameter 0.197m, thickness 0.003 m. C is connected on one side to either a storage vessel or a vacuum pump, and is kept isolated from the tube by a computer-driven electro-valve on the opposite side. The tube is directly opened to R, the latter being separated by a manual valve from a precision gauge, on which the equilibrium pressure is read after transfer. The tube internal and external diameters are respectively 1 and 1.5 mm, with a roughness of 5 pro, typical of a drawn tube (Hodge et al., 1989). Absolute pressure transducers are used to give direct measurements in each reservoir. They respectively
Receiver
Precision gauge
~,e
aE
Vacuum
Fig. 4. Schematic illustration of the experimental apparatus.
300
S. CHARTON et al.
range from 0 to 4 and 0 to 60 MPa, with an accuracy o f l % of the full-scale length. Their output signals are linearized and amplified. The whole chain time constant is estimated as 0.1 s. Piezoelectric sensors connected to fast charge amplifiers are also used to quantify the rapidity of the pressure variations. Fast transient temperature measurements are more difficult to implement. The dynamic characteristic of a thermocouple is indeed strongly dependent on the heat transfer kinetics between the fluid and the junction. In recently published papers, time constants as low as 3.4 ms have been reached with microthermocouples (Beckman et al., 1993). In this study, we have used commercial chromel-alumel thermocouples. The 1 m m diameter isolation sheath is waisted to 0.5 m m over 2 cm from the external wires junction (in contact with the fluid). These features lead to a time constant of 0.4 s for gas-phase measurements. The output from the container sensor is sent through a fast amplifier and is converted afterwards. The receiver thermocopule signals are directly ampli-
fled and converted, providing a total time constant of the order of 0.8 s. Data are collected with an AD-2000 acquisition board. Several experimental conditions were investigated in order to test the model reliability: initial pressures ranging from 10 to 45 MPa, total tube length (from the container orifice to the receiver) of 0.3 and 1.1 m. Furthermore, since the ratio of the heat capacities appears to be an important parameter of the transonic flow involved, helium (7 = 1.67) and deuterium (y = 1.41) gases were chosen for their close mole weights. The experimental conditions are summarized in Table 1, where experiments I - V I I I refer to helium, and experiments I X - X V I to deuterium. 3.2. Data processin# For all the experiments, mass values in each reservoir result from calculations based on experimental data of temperature and pressure. The overall mass balance of the system which is deduced points out the influence of the sensors dynamics on the physical
Table 1. Experimental conditions investigated Gas nature I II III IV V VI VII VIII IX X XI XII XIII XIV XV XVI
He He He He He He He He D2 D2
D2 D2 D2 D2 D2 D2
Tube length L (m)
Container initial pressure (MPa)
Mean absolute pressure deviation (MPa)
Discharge duration ~90 (s)
1.1 1.1 1.1 1.1 0.3 0.3 0.3 0.3 1.1 1.1 1.1 1.1 0.3 0.3 0.3 0.3
10.0 22.0 31.9 45.0 9.9 22.0 31.8 43.1 9.6 20.4 31.0 42.4 10.7 20.1 31.0 40.8
0.452 0.162 0.954 0.208 0.397 0.482 0.844 0.849 0.236 0.410 0.421 0.512 0.967 0.963 1.350 0.926
2.30 2.15 2.00 1.95 1.25 1.15 1.10 1.05 2.55 2.45 2.40 2.30 1.40 1.35 1.32 1.30
18 16
i~
14
12 i0
~
\
6
2 0
0.5
1
1.5
2
2.5
3
3.5
4
t i m e (s) Fig. 5. Apparent mass fluctuation during experiment III.
4.5
A simplifiedmodel for real gas expansion meaning of the measurements. For instance, the apparent gas mass variation for the whole system, calculated for experiment III, is plotted in Fig. 5. From the experimental values of temperature and pressure, it would be concluded that the amount of gas enclosed in the hermetically sealed apparatus does not remain constant during the transfer. This apparent mass vari-
301
ation, noticeable in the first seconds of each experiment (and attaining sometimes 20% in the case of the short tube experiments), is of course inconsistent. It is attributed to the dynamics of the sensors. Obviously, the time scales of the responses of the sensors, especially thermocouples, are of the same order of magnitude as the time scales of the pressure
Pressure (MPa) 35 30
25
20
model 15
''~
",,~ •
10
'
~ •
corrected model ex )eriment
•
It--
0
0.5
1
1.5
2
2.5
3
3.5 4 time (s)
Fig. 6. Pressure evolution in the container during helium transfer (experiment III).
Pressure (MPa) 0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0.5
1
1.5
2
2.5
3
3.51m e (s)4
Fig. 7. Pressure evolution in the receiver during helium transfer (experiment III). CE$ 51-2-J
302
et
S. CHARTON
and temperature variations during the transfer. Consequently, the temporal characteristics of the sensors must be taken into account for a correct model validation. The time responses of the temperature and pressure measurement chains are represented by first-order transfer functions of gain unity, with respective time
al.
constants given earlier. U n d e r these conditions, model predictions can be corrected by convoluting the simulated profiles by the inverse of the sensor transfer function. Since a first amplifier is coupled to the container thermocouple, we found it more meaningful to proceed in the reverse order. The heat transfer between
Temperature (K) 313.15 293.15'
b"
2~.15
'.
253.15
•
:
model
':
233.15
D
c4~L~__~ measuremenls
-
-
r a w me~m~tn~t~
•
213.15
b.
9
•, •
193.15
O..•
.*"
"". 173.15
• "',D "- " ' .
.~.',t'" •
[] "''-,.
"DO
~h ...o-''"
[] .........
.
iC~.._- "
[] .'"~
•
-
153.15 133.15 0
0.5
1
1.5
2
2.5
3
3.5
4
time (s) Fig. 8. Temperature evolution in the container during helium transfer (experiment III).
Temperature (K)
,'°'1
" 423.15
model .
...... •
' ". 398.15
348.15 323.15
corrected model experiment
~b ah~w~lqmu
/
,-°.. ......
i. . . . . . . .
273.15 0
0.5
1
1.5
2
2.5
3
3.5
4
time (s) Fig. 9. Temperature evolution in the receiver during helium transfer (experiment III).
A simplified model for real gas expansion the fluid ( f ) and the thermocouple junction (j) is indeed expressed by
303
Equation (20) can be rewritten as a simple firstorder differential equation, introducing the time constant • j of the thermocouple:
(20)
Tj = - Ajhj_f(Tj - Tf) ( p V C v ) j dd-T
dT~ = dt
where h~_[ is the convection heat transfer coefficient at the junction surface.
( T ~ - Tj-) zj
(21)
Pressure (MPa) 35
20
\
.......
•
model corrected model expenment
•
10
0
.... 0
0.5
1
1.5
2
2.5
- ---3
" 3.5
4 t i m e (s)
Fig. 10. Pressure evolution in the container during deuterium transfer (experiment XI). Pressure (MPa)
0.6
0.5
o., .,
|
•
;- ~/•
0.3
model 0.2
:
;"
corrected model experiment
•
q
0.1
0 0
0.5
I
1.5
2
2.5
3
3.5
4 t i m e (s)
Fig. l l. Pressure evolution in the receiver during deuterium transfer (experiment XI).
304
S. CHARTON
The obtention of corrected experimental temperature T I is achieved by solving eq. (21) from the junction temperatures measurements in the container.
et al.
the temperature profiles in the pressure vessel, Figs 8 and 12, is characteristic of the competition between the two thermal phenomena occurring during the discharge:
4. RESULTS AND DISCUSSION
The results from each experiment have been compared with the corresponding predicted data from the model. The mean absolute deviation between the experimental and the corrected simulated pressure profiles has been estimated in each case (Table 1). Those small departures are mainly attributed to the transducer uncertainty which is about 0.6 MPa. Two representative comparisons, III and XI (referred to as A and B), are discussed in the following. The pressure variations with respect to time in both reservoirs are respectively illustrated in Figs 6, 10 and 7, 11. The raw evolutions predicted by the model have been corrected by introduction of the sensor dynamics as indicated earlier, and the resulting curves, labelled as corrected model, can thus be compared with the measured profiles. The same convention is used for the corresponding temperature histories, Figs 8, 12 and 9, 13, except for the container, where the raw thermocouple output signals have been corrected according to eq. (21). Therefore, in Figs 8 and 12, the simulated temperature profiles are compared with the real temperature evolutions in the container. After the valve opening, the gas rapidly expands from the container to the receiver. The pressure balancing goes on afterwards, together with the achievement of the thermal equilibrium. The shape of
• the rapid cooling, due to the important depressurization, strongly dominates in the first second, • after a while, the endothermic process becoming less significant, the gas is slowly warmed up by the surroundings. In the receiver, conversely, a temperature increase resulting from both compression and kinetic energy conversion into internal energy is first observed. Then, these exothermic processes are counterbalanced by the heat losses through the metallic shell (Figs 9 and 13). This qualitative behavior of the system appears to be correctly described by our simplified model. Furthermore, there is clearly a very good agreement between the predicted profiles and the measured ones in the container: • The departures on pressure are of the same order of magnitude as the transducer uncertainty (Figs 6 and 10); • the mean bulk temperatures predicted in this vessel are very close to the corrected experimental values (Figs 8 and 12), the latter being measured locally in the vessel.
T e m p e r a t u r e (K)
293.15 ~ 273.15 - - Q
• • - -
253.15
r-~
233.15
[]
• []
"
model
corrected measurements
raw measurements
"1% %.
213.15
~...
•
[]
193.15
°. o~I
6 .... b'-'-.
• 0...<2...~-~-"
:~
i1
•
173.15 153.15 133.15 0
0.5
1
1.5
2
2.5
3
3.5
4 t i m e (s)
Fig. 12. Temperature evolution in the container during deuterium transfer (experiment XI).
A simplified model for real gas expansion
305
Temperature (K) 423.15
398.15
373.15 ~
348,15
3
2
0
3
0.5
.
I
1
1.5
5
2
~
2.5
3
3.5
4
time (s) Fig. 13. Temperature evolution in the receiver during deuterium transfer (experiment XI). In the receiver, however, the model obviously overestimates the thermal variations: a sharp temperature increase is predicted in the first milliseconds, which is hardly compensated by the convective heattransfer. The measurements on the other hand reveal a smoother variation and lower temperature values (Figs 9 and 13). Too high predicted pressures follow from this temperature departure as shown in Figs 7 and 11. The weaker agreement observed in the receiver may be due to some inconsistencies in the assumptions made during the model development. Among them, the less reliable are the following: • spatially uniform temperatures • stagnant conditions • natural convection dominance. Indeed, because the receiver has finite dimensions, the sonic gas jet is likely to be reflected by the walls and then deviated by the continuous incoming gas flow. Therefore, nonstagnant conditions are prevailing in the enclosure and prevent the fluid kinetic energy from being entirely and instantaneously restituted in thermal energy. Sharp velocity fluctuations are taking place, thus generating a complex flow pattern in this enclosure, and strongly promoting forced convection because of the resulting buoyancy. Conversely, the velocity of the gas at the pipe entrance is considerably less significant than that of the exhausting sonic jet; therefore the assumptions of gas stagnation and natural convection dominance in the container are justified.
To take into account the influence of forced convection in the receiver, a different correlation for Nusselt number Nu = 0.023 Re °8 Pr °aa
(22)
has been introduced in the model instead of eq. (1). The Reynolds number has been given empirically, as there is no simple way to obtain a mean value of it. Two values of Re were used: 104 and 105 . When Re increases, the temperature decrease is stronger and the curves obtained with this model (Figs 14 and 15) are situated closer to the experimental curves and well under the natural convection curve. Unfortunately, the main problem with the model including forced convection is that the Re number cannot be predicted in the frame of a simplified approach; consequently this model cannot be considered as totally predictive. The particularly important discrepancy between raw experimental and corrected simulated temperatures may also be due to the significant time constant (around 0.8 s) of the corresponding measurement chain. The precision of the correction brought by the convolution process in this case may be insufficient. Lastly, the temperature field in the downstream reservoir is obviously very complex and strongly dependent on time (through the relative positions of the incident and reflected gas jets). It must be kept in mind that in this work we compare a mean value of the temperature (uniform temperature predicted by the simplified model) to a local in situ measurement. For the reasons evoked above, this procedure is appropriated in the discharging vessel, but the comparison
306
S. CHARTON et al.
Pressure (MPa) 0.7
0.6
0.5
'°
°°r
~
0.4
°t
..o~o
.dlo
, °, 4~°'°
j°;
0.3
0.2
f
--
."1
0,j
-.
•
natural convection forced convection - Re=le4 forced convection - Re=le5 ex ~eriment
0 0
0.5
1
1.5
2
2.5
3
3.5 4 t i m e (s)
Fig. 14. Pressure evolution in the receiver during helium transfer (experiment III)--influence of forced convection heat transfer.
Temperature (K) 473.15
448.15 423.15
natural convection forced convection - Re=le4 forced convection - Re=le5 experiment
......
39& 15
•
373.15
323°|5
f ~
•
!1,
/6
It
/,;,~
41.o
"
.... a ....... !; ...... o ...... o
o./,
lw
•
•
•
, - - I ~ a m
oo •
•
273.15
0
0.5
l
1.5
2
2.5
3
3.5
4 t i m e (s)
Fig. 15. Temperature evolution in the receiver during helium transfer (experiment III)--influence of forced convection heat transfer.
is n o t so representative in the receiver. However, in the receiver, the pressure evolution is correctly described by the model, a n d the transfer d u r a t i o n as well as the final pressure are accurately predicted (Figs 2 a n d 7).
Since the model simulates perfectly the p h e n o m e n a occurring in the container, the predicted discharge d u r a t i o n Z9o, defined as the time required to decrease the initially stored mass of gas of 90%, h a s been plotted with respect to the initial pressure, in o r d e r to
A simplified m o d e l for real gas e x p a n s i o n
307
Discharge Duration "c~(s) 3 • •
He- L=l.lm D2 - L = I . I m
2.5 -0 ................................... ....... 0
---..~.
""A
1.5 P. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
/
-0
,.o
.......A-......................... ............ .A
0.5
t,
He - L = 0 . 3 m
0
D2 - L = 0 . 3 m
0
10
30 40 50 Container Initial Pressure (MPa)
20
Fig. 16. V a r i a t i o n s of the d i s c h a r g e d u r a t i o n Zgo with the initial s t o r a g e pressure.
illustrate the properties of the frictional flow (Fig. 16). For a given gas and a given tube length, %0 decreases as the storage pressure increases. The latter is indeed the driving force of the discharge. It appears in Fig. 16 that %0 is also strongly dependent on the tube length and the gas molarity: (a) According to the Fanno flow characteristics, the transfer time is nearly doubled when the tube length is increased from 0.3 to 1.1 m (all other things being equal). The longer tube inducing a higher frictional pressure drop and hence a higher acceleration of the subsonic flow, the thermodynamic consistency requires the inlet gas flow rate to be reduced in a higher proportion for this long tube than for the shorter one. (b) For a given initial pressure and a fixed tube length, the gas transfer between the two reservoirs is shown to be slower in the case of deuterium than when helium is concerned. Keeping in mind that most of the gas is blown out under sonic conditions, the gas velocity at the pipe outlet is equal to the sonic velocity c = x / - ~ ( R / M ) T . Thus, for a given temperature, the deuterium over helium velocity ratio is given by the square root of their ?-ratio. Figure 16 reflects that the mean value of the helium over deuterium discharge duration is effectively close to the reciprocal of this ratio:
~
= 0.92
while
= 0.89.
a thin tube. The flow is assumed to be one-dimensional and quasi-steady. Real gas effects and heat transfer phenomena are neglected in the tube but taken into account in the reservoirs. Finally, providing sensible assumptions, the fast transient flow of the compressible fluid can be solved with minimal calculation efforts, simply by setting the system geometry and the fluid properties. Besides the apparent difficulty of the simulation task, the experimental investigation of the process appeared to be critical. Owing to the fact that the major amount of gas is blown in less than half a second, the measurements accuracy is strongly dependent on the dynamic characteristics of the sensors involved. The actual pressure and temperature evolutions can be estimated anyway by taking into account the transfer function of the measurement chains. The convolution product of the predicted data by the inverse of this function leads to a great improvement in the meaning of the comparison between simulations results and experimental data. The model predictions have been compared with several experimental data. The qualitative behavior of the process appears to be very well described and a very close quantitative agreement is reached in the discharging vessel. A weaker agreement is reached in the downstream reservoir probably because of the stagnant bulk approach derived in the model. However, the receiver conditions have a feeble influence on both the container discharge and the choked flow in the pipe which are correctly described.
5. C O N C L U S I O N
A simplified model has been proposed for simulating gas transfer between two enclosures through
NOTATION
c C/
velocity of sound, m/s friction factor
S. CHARTONet al.
308 Cp, C~
D # h hw L M P R
S t
T u
U V x Z
specific heat, J/(kg K) tube diameter, m gravity, m/s 2 specific enthalpy, J/kg heat transfer coefficient, W/(m 2 K) tube length, m mole weight, kg/mol Mach n u m b e r absolute pressure, Pa universal ideal gas constant, J/(K mol) area of the heat exchange surface, m 2 time, s absolute temperature, K velocity, m/s internal energy, J volume, m 3 abscissa compressibility factor
Greek letters isobaric dilatation coefficient, K - 1 ratio of specific heats 7 tube roughness, m 8 2 thermal conductivity, W / m /z viscosity, kg/(m s) density, kg/m a P characteristic time tube section, m 2 fZ Subscripts f fluid j junction w wall oo bulk Exponents 0 stagnation property
S
choking condition steady condition REFERENCES
Beckman, P., Roy, R. P., Whitfield, K. and Hasan, A., 1993, A fast-response microthermocouple. Rev. Sci. lnstrum. 64, 2947-2951. Cess, R. D., 1973, Free-convection boundary-layer heat transfer, in Handbook of Heat Transfer (Edited by W. M. Rohsenow and J. P. Harnett), Section 6, p. 613. McGrawHill, New York. Fletcher, C. A. J., 1991, Computational Techniques for Fluid Dynamics, Vol. 2: Specific Techniques for Different Flow Cate#ories, 2nd Edition. Springer, Berlin. Haque, M. A., Richardson, S. M. and Saville, G., 1992, Blowdown of pressure vessels. I. Computer model. Trans. Instn Chem. Engrs 70, 3-9. Hodge, B. K., Taylor, R. P. and Coleman, H. W., 1989, Predicting turbulent rough-wall skin friction and heat transfer, in Encyclopedia of Fluid Mechanics, Vol. 8: Aerodynamics and Compressible Flows (Edited by N. P. Cheremisinof), Chap. 13. Gulf Publishing Company, Houston. Issa, R. I. and Spalding, D. B., 1972, Unsteady one-dimensional compressible frictional flow with heat transfer. J. Mech. Engng Sci. 14, 365-369. Leung, J. C., 1986, Simplified vent sizing equations for emergency relief requirements in reactors and storage vessels. A.I.Ch.E.J. 32, 1622-1634. Levenspiel, O., 1977, The discharge of gases from a reservoir through a pipe. A.I.Ch.E.J. 23, 402-403. Reid, C., Prausnitz, J. M. and Poling, B. E., 1987, The Properties of Gases and Liquids, 4th Edition. McGrawHill, New York. Shapiro, A. H., 1953, The Dynamics and Thermodynamics of Compressible Fluid Flow, Vol. I. The Ronald Press Company, New York. Soave, G., 1980, Rigorous and simplified procedures for determining the pure-component parameters in the Redlich-Kwong-Soave equation of state. Chem. En#ng Sci. 35, 1725-1734. Xia, J. L., Smith, B. L. and Yadigaroglu, G., 1993, A simplified model for depressurization of gas-filled pressure vessels. Int. Comm. Heat Mass Transfer 20, 653-664.