A theory is presented which predicts the shock wave overpressure in a cryogenic explosion caused by heating the cryogen to the limit of thermodynamic metastability. It is determined that in a neighbourhood very close to the critical point the divergence of the bulk viscosity plays a very prominent role in the formation of a shock wave; and, in fact, is the necessary condition for its existence. A sample calculation is given for methane that yields an overpressure of the order of 500 Ib in "2.
Shock wave overpressure due to metastable phase transformation in single component cryogens A. H. Rausch and A. D. Levine
In a recent paper, we specified the triggering mechanism for a rapid phase change of a cryogen in thermal contact with a host liquid. The popular term in the natural gas utility industry for the shock wave that follows the phase change is 'LNG explosion on water'. It is now generally believed that a shock wave is formed at the interface between the cryogen and the host liquid and that the overpressure from the shock wave once in air and without a sustaining energy source, rapidly attenuates. It is the purpose of this paper to show that the conditions necessary for a shock wave do in fact exist and to determine a maximum bound on the overpressure resulting from the shock wave formation. It will be seen that the coefficient of bulk viscosity plays a central role in the theory. We first review the heat transfer mechanism between the surfaces of the adjacent liquids and conduct a short discussion of the law of corresponding states near the critical temperature of the cryogen, the region of unstable equilibrium. We follow this with a look at the behaviour of various transport parameters of interest in the neighbourhood of the critical point. The theory is developed using conventional hydrodynamic arguments for the shock wave formation. It is at this point that the necessity of examining closely the behaviour of the bulk viscosity becomes apparent. This we proceed to do and as a consequence show that the necessary condition for a shock wave to form exists. Utilizing some well-known resuits from statistical mechanics we then proceed to an estimate for the maximum overpressure in air immediately adjacent to the cryogen.
If the surface temperature is low enough (it will always be above the boiling point of the cold liquid) boiling will take place by the process of nucleation. If the temperature is above the critical temperature, Tc, film boiling will occur. At an intermediate point, however, a region of unstable equilibrium prevails. If the temperature of the surface lies close to the critical region, an accumulation of vapour will occur at the liquid interface. Yet boiling does not take place because the state of the cold substance does not sustain the vapour phase close to the critical temperature. We contend that a film of vapour accumulates at the surface without boiling, then discontinuities in the heat flux give rise to an accumulation of heat at the film. This, in turn, produces an innate instability which forces the liquidvapour boundary of the cold liquid to recede away from the interface in order to take up the excess heat. In essence, a shock wave occurs. It is likely that this model may explain the nature of instabilities that have been observed in such liquid-vapour contacts.1 In Fig. 1, we see a typical liquid-liquid contact. Let medium 1 be the host liquid, medium 2 the cold liquid. The primary mechanism for heat transfer is governed by the thermal conduction equation aT
pCv at = K V 2 T
(I,1)
Background We are concerned, in the first instance, with the problem of a liquid that is cold being spilled on another liquid, say a host liquid. The cold liquid will begin to evaporate. Then the mechanism for evaporation will depend on the surface temperature of the point of contact.
H e a t flow ~
~
--I nterfoc¢ Cryogen
5 c,p,
-.x2 Ts
AHR is with the American Gas Association, Arlington, Virginia, USA. ADL is with the Physics Department, University of West Virginia, Morgantown, West Virginia, USA. Received 20 September 1973.
C R Y O G E N I C S . M A R C H 1974
Fig.1
p, X
Plan of liquid--liquid contact
139
where O is density, Cv is specific heat at constant volume, and K is thermal conductivity. The boundary condition for continuity of heat flux is KilT,
(2)
"-~ = K2~T2 "-~
The surface temperature, Ts, will then be such as to satisfy this condition. The solution to the problem is well-known. 1, 2 Thus if we have a semi-infinite conducting medium, initially at a temperature To, with a plane surface, initially at a temperature Ts, the temperature at any given point is
=
+Ts
-
(3)
Pv-P---------~c-Pl-P----------c = (~---~)[~l Pc
(P + ap 2) (1 - pb) = pkT
1 e"~2d~
(4)
(11)
Pc
=
a
~-Pc
-
27b2
8a kTc
-
27b
(12)
If we expand (11) according to the Landau Liftshitz prescription we find
and
a2 = K/aCv
(5)
Now let Tlo and T2o be the initial temperatures of the liquids in Fig.1. Then, combining (2) and (3) we readily show that
Tlo + XT2o
T=
(6)
I+X where
A=4
3 C=2
B=6
(13)
Similar results can be obtained using other models (the Landau Liftshitz relation, as we have pointed out earlier, is model independent). The liquid-vapour boundary region can be obtained using a linear approximation to (8). Thus P - - "" 1 + m~
× = K2oq/K1 or2
(7)
If Ts is well below the critical point, or well above the critical point, the ordinary mechanisms of nucleation or film boiling will occur. One may study these by solving the hydrodynamic equations subject to the appropriate boundary conditions. An essential requirements for the success of such an approach is that the equation of state of the evaporating substance must be known. This will always be the case if we are not too close to the critical temperature. Little reliable information is available regarding the behaviour of substances close to the critical temperature. The law of corresponding states is widely accepted. A well-known model independent equation was proposed by Landau and Liftshitz 3 P-ec Pc
- At + B~3' + c3' 3
(8)
Pc, Pc, Tc are the critical pressure, density, and temperature. T-
~---
140
Pc
Ultimately, a realistic equation must be governed by knowing something about the mechanism of intermolecular interactions. Thus, a simple model which nevertheless provides us with a fundamental understanding of critical point behaviour is the van der Waals' model
y
erfO') = ~
(10)
The constants a, b are related to the critical point coordinates
where the error function is defined as
2/°
This is a crude example of a more general theory based on scaling hypothesis. 4 This model is valid only if T < Tc. Many important results may be derived from (8). Thus, if Pv, Pl are respectively the vapour and liquid densities,
Tc Tc
3'-
P--Pc Pc
(9)
(14)
Pc
where m =A + -4B 9
~(3-~C)
(15)
Now for the van der Waals' model m = 5.7
(16)
Then, when P = 0 T - - = 0.84 Tc
(17)
We do not believe that other models will yield substantially different results. Thus in the region 0.84 Tc < T < Tc
(18)
a configuration of unstable equilibrium involving liquidvapour transitions will occur.
CRYOGENICS. MARCH 1974
Host liquid r
~ Region I I superheated vapour ep, r
M oving boundary Region2 liquid
Shock wave theory applied to liquid surf'ace contacts We start with the hydrodynamic equations.
qp.r
Continuity X
rs
~.
ap
-+ bt
Fig.2 Planof liquid--vapour-liquidcontact
=0
(24)
Momentum transport
Recent investigations have shown that the Landau Liftshitz model does not agree with experiment, s Thus it would be better to let
P
[ ff+( v-3] 8t
(25) Here ~ is a tensor (traceless) with components
I P - P.____c= pltjla
(19)
Pc
dvj
avj
2
=__ +___
where a is a constant, and a is closer to 1/3 than to 1[2. An essential problem is that a knowledge of heat capacities and transport parameters is required. Although our knowledge i n this respect is incomplete at present, we quote some wellknown results, which we will use later, obtained from contemporary-scaling theory. 6 Specific heat
Cv ~- Cv01~ld
Thermal conductivity Viscosity (shear) Bulk viscosity
(26)
6i]~.~
and p ~ is the external potential, r/and X are coefficients of shear and bulk viscosity.
Energy transport
)
as p
p
~-+h
K ~- Kol~l c = ~" (K ~ T) + ~ (rt ~" A~+ X ~ ~ . ~ )
~ "-"constant X ~ Xol~l-b
Thermal conduction coefficient
~2 = ao 21~l-e
Here d -~ 0.1, c -~ 2/3, e ~ 2/3 and b is not known at present.
(20) (21)
At any rate, if such a boundary layer interposes itself at the surface of contact a discontinuity in the heat flux will result. This can be equalized by movement of the liquid-vapour boundary. As in Fig.2, we note that the host liquid is in contact with medium 1 (a vapour) which, in turn, is contigous with the cold liquid (medium 2). Then i f j is the particle flux, and l the heat of vaporization the boundary condition at surface
(a) is
(27)
Here e is the specific internal energy and h is the specific enthalpy. From the laws of thermodynamics we have P Tds = d e - ~ dp
(27a)
and 1
dh = Tds + - d p P
(27b)
P
h=e+P
(27c)
Heat flow
KH -VTrI"-~ - K] ~T, "-~ = lj
(22)
This will supply the energy needed to sustain the shock. The discontinuity in temperature gradients at suface (b) will be considered in the next section. The onset of such a shock will give rise to a pressure increment (Ap)2 2p
~- t/
These will now be considered.
CRYOGENICS. MARCH 1974
pT
=
+ x (~- v)2
(K
r) + '7
(28)
We will assume, as in Fig.2, that the liquid-vapour boundary is moving with velocity u. Then any of the variables can be represented by a series of step functions
(23) 0 (y)= 1
y>0
0 Cv) = 0
y < 0
(29a)
141
dO - - = 6 (v)
(29b)
dy
Subtract (36) from (35)
Then T (x, t) = T, (x, t) O (ut + x) + 1"2 (x, t) O (x - ut) p (x, t) = p, (x, t) 0 (ut - x) + P2 (x, t) 0 (x - ut)
(30)
--~- F~,X2 ~ -
Substitute (30) into (24)-(28) and, we obtain at the moving boundary the following relation s .
u-b-,,_~ • ~ - ~
+x2(~2-u-+)~ • v ~ - x , ( - v l - u - ~ . ~] (37) Continuity
Here
p, (~-u--~" ~=02 (~2-u-'~" ~'=j
(31)
where j is the mass flux. Now consider momentum transport
~. [p2~-~
-~
-~-p,~
-~1
u
=
h
-
(38)
Ts
where p is the chemical potential. Let us now combine (31) with (33)
(Pl --P2)
- ~" (,~2_~2 - ,~,~,) + ~ [e,2 - p,) ]2 =pl
-(X2~'~-X,~'~)]
=0
- ~).
=
~+p2
--
p,
--
(x2~
-.),.
(33)
~z
+/z2---V" P2
P2 [ ~ - u')" h+] [ ~ - ~')" ~n] - ~ . ('~2~2 " m) = p, [ ~
- y).
• n,
-
v2
xl (~" ~,) ] -
Equations 39 and 40 are the working relations for our shock wave analysis. As pointed out earlier, we expect to apply these to a substance close to the critical temperature. Then according to (20) we ought to ignore the shear viscosity terms. This finally yields, in place of (39)
[ (~+2-u-+)2+h212
vT2 - - K I
60, - , o 2 ) vT1
6o, - ps)
-~ (X2v'72
"1" T/2~2 ° ~ 2 - - U ~
--;kl V " ~1)>0
- 7 , 6~ • ~ -~')+ x 2 ~ - u-b~ • - X,(~ -u--~
(40)
Conditions necessary for a shock wave
j2=p, p 2(Pl _ _- P2) + - PIP2 [K2
n[
--i ,72 82" ~ - u - 3
-~,X,. ~-~)]
(34)
Transport o f energy
=h'~"
(39)
]-[
P,
This relation can be satisfied by conventional boundary conditions relating tangential components of velocity, and will not be considered further.
P2[(v-~--u--~)" ~]
+u,
~] [ ~ - ~ ) . ~]
, "m
> 0
Pz - P2
Similarly, if we combine (31) with (37)
Now postmultiply (32) by ~ a unit vector normal to
-n
(~2_,'3~-,7,L) • ~
+ P, P2 n •
" --" 'P2
- X,~" ~) = 0
P2 P l -- P2
h*] 2 _ p, [(7, - ~'). h*] 2 _ h*. (~_;$2
-,~,X,)__
+PI
P l -- P2
(32)
We can postmultiply this relation by ~ to obtain o2 [~
P2
• 71
(35)
and, in place of (40)
u-h ~ - 2? + ~
Finally, for heat flow
(41)
__P2 (~. 4 ) =
2
p2(-~2 -u-*) " -~ T2 s2 - p, T, s,(-~ -u-+) • -~ = k+" (K2 ~T2 - K, ~ T ~ )
142
(36)
+O, ---
Pl
(-V" ~ )
(42)
CRYOGENICS. MARCH 1974
It is clear that for a non-viscous medium the condition (41) reduces to
where
13= - -
z~o ___z"> 0
which means that the transition from state 1 to state 2 occurs via a series of thermodynamically accessible states. This is not true for the system we are considering since
Far away from the critical region, ~t is small and K is approximately a constant. Thus (48) becomes
Dt = Pl "~ Pc
P2 ~-- 0
Pl ~ Pv
P2 ~-- Pl
so that
(49)
Pc Tc Cv
(43)
zxp
Pc
K pCv (1 + t3rl~la- 1
V 2T
(50)
which is the ordinary equation for thermal conductivity. Note, that in obtaining (50) the minus sign was used in (48). Presumably the plus sign gives us an extraneous result. Close to the critical region we obtain, using (20)
Pl > P2 but
~ " ( K ~ T)~-c K(VT)2 I~1 Tc
Pl
and (48) yields
Ot = ariel---Za - ½
(51)
4(c:c,) --
~T
(52)
This occurs if T < Tc. Finally, substituting (52) into (56) ,
V " v = - - Dt p P
(44)
v.
U / j rl
(53)
Then referring back to (51) (note ?t2~ • ~ = 0) we see that the condition for a shock wave is not satisfied. Suppose, however, T > Tc, then (51) becomes
where ~P Dt p = ~ + ( ~ " ~)p
(45)
that is, it is the total time derivative of the density. Assume, for the moment, that the temperature is slightly below Tc. Then from (19), (44) becomes
(54)
~Tc
and (48) yields
Dt-
~ . V = + ariel" - I D t T Tc
ar~a_v "
[~TI
(46) and (48) yields
We now rewrite (28) as p Cv Dt r + p
K(~T) 2 V . (K~T) ~ - c - -
V'v-
ITI
(54a)
- ~-ff Dt p
which according to (41) allows for a shock wave if ~ is small enough. To estimate the size of ~, note from (41)
+ X ( ~ " v-~)2
(47) X~" v*>(Pl-P2)
Combine this with (19) and (46) and we obtain Since Pl -~ Pc
Dt
PCv Tc 2 (1 +/~ P I$1a - 1)
•V
2)~/2 172 1~12a-2
P2 = 0
• V>Pc
If ~"is the thickness of the boundary layer
{ 1 + ~[1 _ 4xa2r2 I/~1~-2 V " (K~ T)]
-
p-TC~Cv 2~'~ 0 + ~
CRYOGENICS. MARCH 1974
i~--T) T
_1}
(48)
AT
0.16Tc
143
Then
where fsv is the free energy of the saturated vapour at atmospheric pressure. But the last term in the brackets becomes (approximately)
0"16Tck 4 (~_~_c)
~.~
>Pc /2sv --/22 = 0
Another condition on the gradient is given (only approximately) by (22). Presumably one could combine these two relations to obtain estimates of ~"and ~. Unfortunately we do not have precision measurements on the transport coefficients close to the critical temperature. An interesting possibility is to define a Prandtl number for bulk viscosity
K This is approximately temperature independent. Then we would find
I~[-b-d/2-1VVc
0.16Tc - -
X
Cvo
>Pc
(02
(f_ Av) > 0
(59)
o~)
-
]2 iO12
]2 -
----
'
[=--
loge
If
z/
1
(60a)
where Zi is the internal partition function
Again, in the absence of reliable data on coefficients and exponents in (20) it is difficult to make a numerical estimate. We turn now to (42). If we assume that the transverse velocities are negligible, then from (31) _
/2 - 2022 022
Not it is quite difficult to obtain data on the free energy near the critical point. We can estimate the effect as follows. For an ideal gas the free energy is
XCv
p, _
so that finally (57) becomes
(55)
/922
Zi = X] gje"Effkr
(60b)
Here the sum runs over the internal molecular states and For a real gas, near the critical point, one may show that the density p in (60a) may be replaced by its value for an ideal gas, p/kT and a term ap/tm added. Since, by (12) a is 9/8 kTc/Pcin the van der Waals' model, (60a) becomes
g/is the degeneracy.
and, also taking
3,2
--+
V " v2 -~0
f= __
loge
m
\ ink/
(1)
(42) becomes P
]2
+
2 2 (.°22 - P22) - (/22 - / 2 2 )
2pip2
9 + i-g
(-~)~c
(56)
For the last term in (56) we may substitute the expression in (41) that is,
k2V" ~ - P2 --Pl /2+p2_p2 Pl P2
(p2-p2)~ O q - A ) ]
/2 =-[ 2p2 p2p----~2 + 2pI2 p22
9 kTc 16pcm
144
(62)
(57)
(58)
The right-hand side of (57) can be rewritten as (assuming a quasi equilibrium between the liquid and the saturated vapour along the wave front)
p:)
Pl T2 si2
(~sv- ~1) loge
3
?'~ Pc
m
f=~-Ts
- -
rn
After some rearrangement of terms this becomes > 0
where f i s the specific free energy
/SV -- f2
kTl loge
(Psv - p2 )
(F, - Fsv) = krc
0el - fSV) +
3/2 1
m
+ k Tsv loge Psv m (Tsv) s/z
+
Then (56) becomes
6o2-p~)
( 27rh2 )
k ( T s v - Tl) /i --.fsv -
60,
(61)
Then + - - V " v2 =0 P2
2Pt p22
]
+
-- --
Pl
+ (1 + ~sv) loge Psv( 1 + ~sv)"s/2
Pc - ( 1 +~l)loge Pt (1 +~1)'s/2
Pc 9 (,Osv- pl) ] 16
Pc
(63)
CRYOGENICS. MARCH 1974
where
Combining (66a) and (68) we finally have
X'=
1.15 X/l~ll] 2
/ 2 = 0 . 5 (~SV -- ~1) p 2 V c 2 [1 -
l( /22n\'m h 2- k) T
(64)
is a thermal wavelength. If we now use (I 4) as an approximation for the pressure
fsv -- fl --
To estimate ]2 take a limiting case where ~sv = - 0.16
kTc ) [(~sv-~')l°ge m (8e 3~pc
]2 = 0.5 pc 2 Vc2 [1 + 1.15 (0.4)] 2 [3.2 + 1.5 -- 9.3](69) where we have used (16) and the fact that
+ (1 +/Jsv) loge (1 + m~sv) (1 + ~sv)'S/2
3OcXc3
- ( 1 + ~s) loge (I +m~1)(1 +0.60
(~[~sv[
loge - 8e
+~i) -5/2
~1~,1)]
(65)
-
9.3
then j2 = 4.5 pc 2 Vc2
where (10)has been used. To first order in It[ (fl - f s v ) = Uc2 [ 0 - 3 3 ( ~ v - ~ l ) l o g e
(3 .-)
or
j = 2.2 Pc Vc Then since for methane
+o. o
Vc = 1500 ft s"1 andPc = 0.31 slugs ft "3 ] = 1035 slugs ft "2:1
(70)
(66) O ~ r e
The energy flux by (37) is
where Vc2 _ 3kTe
(67)
rCi
- --=] A dt
( -v+" u-c)2 + / a - ~ , ~ - " v" p
(71)
The net transport of energy is is used to represent a thermal velocity. After some rearrangvment of terms
(fI -- f ~ ) = 0.33 Vc2 ] loge
which by (37) is
L÷
(72)
A- /
l-
( m--~') A
÷
(~s,, - h )
1
(66a)
Using (59) and (10) we obtain
-
J[#l-/a2-Xl~"
~]
(73)
Using (56)
A
]
tq -/asv
(P2 - P 1) ] 2 _ P ll + p2 P2 Pl
(74)
By arguments analogous to those used to obtain (59) this becomes
j2 =
1.5 [t - 1.15x/|hl] 2 [1 + 1.15 v ' l h l ] 2 0 c 2 f f s v - A ) [x/l+ll + x/Ihl] 2
CRYOGENICS. MARCH 1974
(68)
,4-
/
fl - A v - p 2 - p l / 2 p2 P2
(75)
145
which by (59) is
[ 602 -
o~)2/2
L 2P12 p22
A
__,2]
P2 -- Pl p2 P2
_ ] 3 (pl 2 _ p22)
(76)
2 p l 2 p22 Now
(ae)2
(77)
A 2PaCa
where Pa is the density of air, and Ca the velocity of sound in air X/(2PaCa)
AP=]3/2
4tl (p22 -
pI2'~
I k2--~1:232-~ ]
(78)
]3/2 -
X/[2Ca (X/I~21 + X/[~xl Pa)]
(79)
Pc
=2.2Vc
=l.3]V c
PaCa(X/[~21 +X/I~I[)]
Conclusions and suggestionsfor further investigation
Ts = 0.84Tc
~[ j2pa ca(O'4)]pc Vc ~(Papc
C~aca)
(8o)
Now Pa = 2.7 x 10 -4 slugs ft "a
Ca = 1 100 ft s"1
j = 1 035 slugs ft "2 s"l
Vc = 1 500 ft s"1
pc = 7 x 10"2 slugs ft -3
AP= 538 lb in "2
(81)
(82)
Discussion of results We have derived a condition for the formation of shock waves when liquid surfaces come into contact and the surface temperature is near the critical temperature. We do not expect this phenomenon to occur at liquid-gas contacts because of the low thermal conductivity and density of gases. Thus, in (6) X will go to zero and the surface temperature will be essentially that of the liquid. Alternatively in liquid-solid contacts, if the host liquid is replaced by a metal, then the surface temperature is essentially that of the metal. A more complicated situation, governed by (6), occurs in the contact of liquids with non-metal surfaces. Presumably, if the temperature of the solid is right, the critical rggion would prevail. This is complicated, however, by the nature of the surface, and the possibility of running into mo~fying surface effects like adsorption. At any rate, under ordinary circumstances, we do not expect detectable critical surface effects in liquid-solid or liquid-gas contacts. Still, we have been unable to rule out the latter possibility completely, and the matter should be investigated more fully. An essential characteristic of the shock waves we are talking about is the dominant role played by the bulk viscosity. Thus, (41) would not be valid
146
Of greater importance, is the assumption of a pure substance. Evidently, the thermodynamically accessible paths for a mixture may be much greater than that for a pure substance. Presumably, if there is a mixture of loosely coupled .substances, the model which we have presented is still valid. The more volatile substance will boil off first, allowing the denser substance to reach a critical region. There is some evidence to support this hypothesis.
We have derived a condition for the onset of a shock wave in liquid-liquid contacts when the surface lies in the vicinity of the critical temperature, that is,
or
AP= /" 4[2/" Pc
if there were no term of this sort. Similarly, the assumption regarding energy transfer assumed in (22) is too naive, since the difference in energy levels is governed by the difference in chemical potentials, which is zero for a phase change. Moreover, by taking the chemical potential for an ideal gas, we find that the chemical potential decreases when one makes a transition from the saturated vapour to the superheated vapour. Evidently, the only reason why one obtains energy transfer is because of the extra squeeze introduced by the viscous force.
(83)
The model seems to agree fairly well with experimental results. It is, for the time being, limited by our lack of precise data on critical point transport parameters [(54a)]. Also, we have not realistically calculated the change in free energy in going from the saturated to the superheated vapour, though this is not an essential shortcoming. A greater difficulty is that of extending the model to mixtures. Thus, if the scaling hypothesis, which provides us with a convenient generalization of the law of corresponding states is essential in our argument, how will it work when the concentration of some component becomes a variable? Preliminary studies show that in the case of weakly coupled liquid mixtures our model still applies, but this point needs further study. Further attention should be directed to the internal partition function to see in what way a restructuring of molecular configurations can affect the relevant results. At any rate, we have studied a case where conventional boiling theory breaks down because of the branch point discontinuity in thermodynamic and hydrodynamic parameters that occurs at the critical point. Presumably, the justification for our model must rest on experimental verification.
References 1 Levine,A. D., Rausch, A. H. 2 3 4 5
6 7
Cryogenics13 (1973) 224
Drake, R. M., Eckert, E. R. Heat and Mass Transfer, 2nd edn
(McGraw Hill, New York, 1959) 93 Landau, L. D., Liftshitz, E. M. Statistical Physics (Addison Wesley, Reading, Mass, 1969) 268 Stanley, H. E. Introduction to Phase Transitions and Critical Phenomena (Oxford, New York, 1971) 175 Stanley, H. E. Introduction to Phase Transitions and Critical Phenomena (Oxford, New York, 1971) 172 Stanley, H. E. Introduction to Phase Transitions and Critical Phenomena (Oxford, New York, 1971) 203 Huang,K. Statistical Mechanics (John Wiley, New York, 1965) 198
CRYOGENICS . MARCH 1974