Compurers & Sfrucrures, Vol. 8. pp. 31-39.
Pergamon Press 1978.
Printed in Greal Britain
DYNAMIC RESPONSE OF A CONTAINMENT VESSEL TO FLUID PRESSURE PULSES FRANK L. DIMAGGIot and HANS H. BLEICHt Department of Civil Engineering and Engineering Mechanics, Columbia University, NY 10027,U.S.A. and JOHN M. MCCORMICKZ Weidlinger Associates, Consulting Engineers, NY 10022,U.S.A. (Received November 1976)
Abstract-Two methods, one of which is a very simple approximation, are proposed for the dynamic analysis of the response of the wall of a nuclear containment vessel to the fluid pressure exerted on it when the relief valve discharge piping is cleared.
1. INTRODUCTION
2. MlRMUL4TION OF PROBLEM
The water-filled cylindrical steel shell, shown in Fig. 1, stiffened by steel T-beams as in Fig. 2, is acted upon by a cylindrical pressure wave whose time history is shown in Fig. 3 (with magnitude normalized to 1 psi), and whose variation with depth is illustrated in Fig. 4. This paper is an analysis of the dynamic response of this model of a portion of a nuclear containment vessel to an excitation approximating the pressure exerted on the vessel wall when the relief valve discharge piping is cleared. In the Appendix, an approximate analysis of the same problem is presented. This analysis is useful for preliminary design and as a check on the computer results of the more refined analysis.
tProfessor. SAssociate Partner.
The radial displacement of the shell, w(z, t), assumed positive outward is governed by (Ref. [l])
In eqn (1): z is the axial coordinate measured from the top, as in Fig. 1; h is the shell thickness; rj is the radius of the shell wall, as in Fig. 1; E is Young’s modulus for steel; I is the second moment of area, per unit of shell circumference, of the combination of shell and stiffeners shown in Fig. 2; pi(.z, t) is the incident pressure on the shell wall, as shown in Figs. 3 and 4, caused by relief valve clearing; p(r, z, t) is the fluid pressure induced by the shell wall motion; pS is the mass density of steel; and t denotes time. Assuming linear compressibility, the pressure p(r, z, t)
Symmetry
I-,?
Woter
Shell wall with longitudinal
stiffener?
Fig. 1. Model of containment vessel analyzed. 31
32
FRANK
L. DIMAGGIO et al. in which u(r. z. f) is the radial comnonent of the fluid velocity u. r is the radial coordinate as in Fig. 1. The boundary conditions to be satisfied are:
-
w(l, f) = 0
$y,f)=O.
(5) (6)
Equations (5) and (6) assume that the shell bottom is fixed. 1,as shown in Fig. 1, is the depth of the shell wall.
‘e $(OJ)=O
Fig. 2. Stiffenerconfiguration.
(7)
F: 0 25
Fig. 3. Time history of forcing function, p,(z,t).
~(O,f)=O. Equations (7) and (8) assume the shell top to be free. l((r,, 2, t) = 0 (0 S z c I,).
(9)
Equation (9) assumes the pedestal (rigid cylindrical surface of Fig. 1) to be rigid, so that the radial fluid velocity is zero there. 15is the depth of the flat rigid bottom from the free water surface, as shown in Fig. 1. u(r, L, t) = 0 (r, S r g r2).
(10)
Equation (10) states that the flat portion of the vessel’s bottom (see Fig. 1) is rigid, requiring that the axial component u(r, z, t) of the fluid velocity u vanish there.
13 B
Fig. 4. Variationof forcingfunction p,(z,t) with depth. in the contained water, and its (vector) velocity u(r, z, t), satisfy the equations (Ref. [2]).
[U sin (Y+ v cos Q]~=,+(~_,,~~~~ = 0 rz G r G rs. (11) Equation (11) states that the sloping bottom of the containment vessel (see Fig. 1) is rigid, so that the component of fluid velocity normal to it is zero. p(r, 0, t) = 0 (r, C r S r3).
in which: p is the density of water; c is the velocity of sound in water; and V is the (vector) gradient operator, and the dot denotes scalar multiplication. The continuity requirement, that the radial fluid velocity at the shell surface be equal to that of the shell, may be written as (4)
(12)
Equation (12) states that the fluid pressure on the free surface is zero. Initial rest conditions require that w(r, 0) = 0
!$(z,O)=O PV, z,O)= 0
(13)
(134 (14)
Dynamic response of a containment vessel to fluid pressure pulses
u(r, 2, 0) = 0.
33
(15) responding to eqns (2) and (3) are of the form
The solution of eqns (1-15) yields w(z, t), p(r, z, t) and u(r, 2, t).t From these, the stresses in the shell and stiffeners may be obtained as follows: The maximum circumferential stresses a&, t) are, referring to Figs. 1 and 2, 6vEIw”
(16)
but the bending components are small. Neglecting shell curvature. the axial stresses in the shell and stiffeners are (17)
in which 1, and q2 are distances from the neutral axis to the extreme shell and stiffener fibers, respectively (see Fig. 2).
3. NUMERICAL SCEEME-m
where i is an index in the radial direction and j is an index in the vertical direction. The superscript n or n + 1 denotes the time variable. Equations (20~(22) make use of central differences. The error in this approximation is of order (As)’ where As is the appropriate grid spacing. The error is of order As wherever central differences cannot be used, e.g. at boundaries and at points where the grid spacing changes, i.e. points along r = r, or z = 1. The finite difference equations for the shell, replacing eqn (l), are of the form
I
x [ ( WW,” + $ we”- jj (pii + p,“Od)
IuFmRENcEs
3
The basic set of equations to be solved is eqns (l)-(H). In W,“+’ =
the fluid we specify a mesh of points at which the variables p, u and u may be defined. The finite difference molecule is V
UP u
(19)
V
and the distribution of these molecules that was used in one of the calculations reported here is shown in Fig. 5. The finite difference equations for the fluid cor-
w,”
+
At,+=ih”+’
(23) (24
where (Yis an index along the shell and (WV),” is a fourth-derivative operator whose form depends on the value of the index a as shown in Fig. 6. For example, the fixed end is denoted by a = 1. At the next point, a = 2, (WV),” is given by h4( WZV)2”= Wq”-; Wj”+ ; W2” 23 22 21 20 19 I8 17 16 15 14 13 12 II IO 9 8 7 6 5 4 3
12
(24a)
a=N
II IO 9 8 7 6 5 4 3 2
Fii. 5. Coarse distribution of finite difference molecules used.
:
I
*o=I
.i’iiiiiip:‘:i:i:I:~:.:.~:~::::::~~~~:~:~:~~::~ :.:.:,:,:,:,:,:,:.:. _....,.,.,,,,,.,.,.,.,.,
tNote that eqns (3), (4), (12), (13) and (14) constrain the shell displacement at the fluid surface to be zero, i.e. w(0, t) = 0. CAS Vol. 8, No. I-C
CF.11
Grid spacing N=23
Grid spacing N=l2 Fig. 6. Grid indices on shell surface.
34
FRANK
while at all interior points except a = N - 1, (WV).” given by
L. DtMAooroet al. is
h4(WIV),” = w:& - 4w:_, t 6w,” - 4w:+, t w:+z. (24b) At (Y= N - 1, and Q = N, use is made of image points w:,, and WC+* which are determined from the requirement that the end of the shell be free at the fluid surface. The values used at these points are wh+, =
30~~” - 18wk_,t ~wE;_~ 14
(24c)
wPf+z =
30~~” - 32wf_, t 9w;;-, I
(24d)
Now ( WIV),” at the last two points on the shell is given by eqn (24b) which requires eqn (24~) for the case (I = N - 1 and requires both eqn (24~) and eqn (24d) for the case a = N. The boundary conditions for the fluid, replacing eqns (9)-(12), in finite difference form, are
where AV is the volume of the shaded element, AM,A, and AE are the surface areas of sides B, D and E and uq is linearly interpolated between the nodal values u, and u2: uq = (3u, t u&/4
(28b)
Note that side C makes no contribution to the total flux since the normal velocity vanishes there. The value p for the shaded element is obtained from eqn (28a) upon dividing both sides by AV. The process is repeated for all elements bounded by the slanted surface. In finite difference form, the initial conditions of eqns (13)-(E) are: p$=O
(29)
u&=0
(30)
f&=0
(31)
wo=o * kO=O 04 .
(32) (33)
4. NUMERICAL RESULlTFOR‘EXACT” RF.SPONSEOFsHEu,ANDFUJlD
(26)
0
In this section, the numerical scheme outlined in Sec(27) tion 3 was used to obtain the response of a stiffened shell and fluid with the following geometrical and material properties
l& = 0
where the index j=lisatz=O,theindexi=/3isatz=L,forr,crS r,. At the shell-fluid interface, the finite difference form of the continuity equation, eqn (4), is 3”ZU” = 7.a
(28)
where y is an index in the fluid for r = r+ The boundary condition on the slanted interface is satisfied by assigning the velocity normal to that surface a value of zero during the volume integration of eqn (2) after the divergence theorem had been used to transform the right hand side into a surface integral over the normal velocity. To illustrate, consider the shaded element shown in Fig. 7 bounded by sides B, C, D and E. In numerical form, the application of the divergence theorem to eqn (2) for this element leads to
0AV = -pc2(II&
- u3& - u&)
1
UC?
Side E
Normal veloaty equals zero on side C
Fig. 7. Detail of slanted boundary calculation.
h = 1.44in. 1=20ft2in. d = 5 ft 6 in. r3= 42 ft 10: in. r2=24ftOin. r, = 15 ft 2 in. L=30ftllin. p = 1.94lb sec*/fP ps= 14.9lb sec2/ft4 E = 30 x lo6 lb/in? c = 4790 ftlsec A = 40 in. I = 5.4 in.’
The moment of inertia I per unit of circumferential length recorded above was obtained using stiffeners and (28a) neglecting shell curvature between stiffeners. Calculations were made with the loading shown in Fig. 3[3] for the following structural configurations: (1) the shell alone (without the fluid interaction) to determine its vibration characteristics and to provide a check solution for part of the full problem. (2) the shell with fluid interaction but with a horizontal rigid bottom at the level of the base of the shell. These solutions were compared with the approximate analysis given in the Appendix and also serve as check solutions for the full problem. (3) the she11with fluid interaction and with the slanted rigid surface adjacent to the shell as shown in Fig. 1. The last configuration was used in several calculations to determine the proper parameters of computation (time step, spacing of mesh in the radial and vertical directions in the fluid, spacing of mesh points along the shell, etc.). The results from two of these calculations are reported here. They are: (1) Case A where the shell is divided into 11 equal segments of 22 in. each (i.e. 12 points on the shell). The
Dynamicresponse of a containmentvessel to fluid pressure pulses fluid mesh is as follows: (a) Ar = 23.56in. r, G r s r, (b) Ar = 30.20 in. r, < r d r3 (c) Az=22.OOin. Osr~l (d) AZ = 17.20in. I < z d L. The calculation was run for 4650 time steps of 0.0002 set, covering the entire 0.93 set for which the pressure loading history is specified. (2) Case B where the shell is divided into 22 equal segments of 11 in. each (i.e. 23 points on the shell). The fluid mesh is as in Case A except that AZ= ll.OOin. O
35
show that the maximum displacement occurs at point 6 (on the third negative peak) and has a value of 1.91 times the static deflection. Figure 9 exhibits the computational error in Case A in calculating the shell displacement at the fluid surface, which should be zero by eqn (F.l). Figure 10 presents the results for point 11 of Case B which is identical to point 6 of Case A. Again, calculations show that the maximum displacement occurs at point 11 (on the fifth negative peak) and has a value of about 2.03 times the static deflection. Comparison of corresponding curves from Case A and Case B [such as W(6) from Case A with W(ll) from Case B] shows excellent correlation of the response with regard to both peak values and frequency content. (2) Moment time histories at the base of the cantilevered shell. Figures 11 and 12 show the values for the entire loading computed by the formula, with error of order h*, EZ[ W(3) - 8 W(2)]/(2h*)
for both Case A and Case B.
0.00634in., which is the maximum static displacement of the shell under the positive peak value of the specified loading, and they are all plotted to the same scale. The pressure loading is also shown as a convenient and meaningful way of giving the time scale and showing correlations between the forcing function and the response. The first 2 curves (Figs. 8 and 9) present the full displacement time histories of two points on the shell for Case A. Point 1 is the fixed end of the cantilever shell and displacement is always zero there. Point 6 is the midpoint of the shell, while point 12 is at the tip of the shell (at the level of the top of the fluid). Calculations Fig. IO. Shell displacemdnt at middle of shell (point 11, case B) vs time.
I .6 1.2
2.4
0.6
I.6
0.4 i$
0
.s
z
1.2
.E 3
0.6
“9 .c
0
-04 -06 -1.2
$
-1.6
-0.6 -1.2
-2.oJ
Fig. 8. Shell displacementat middleof shell (point 6, case A) vs
-1.6 -2.4J
Fig. 11. Moment at base of shell, case A, vs time. 1.6 2.5 ,
t
-1.61 -2.0 J
Fig. 9. Computational error in computing shell displacement at fluid surface (point 12, case A).
-2.5j
Fig. 12. Moment at base of shell, case B, vs time.
36
FRANK
L. DIMAGGIO et al.
(3) Shear time histories at the base of the cantilevered shell. Figures 13 and 14 show the values for the entire
1.2
loading computed using the formula, in this case with error of order h,
06
i
3EI[ W(3) - 4 W(2)]/(2h3). (4) Depection patterns. Figures 15 and 16 show the deflection pattern along the shell at the times when W(6) from Case A and W(ll) from Case B reach peak values. These curves indicate the deformed shape of the shell during vibration, with the curve for Case B showing greater detail but only for the bottom half of the shell. Note the computational error exhibited by the non-zero displacement at the fluid surface (Fig. 15). (5) Moment diagrams. Figures 17 and 18 show the moment diagram for Case A and Case B at the time when the moment at the base is near its maximum value. (6) Shear diagmm. Figures 19 and 20 show the shear
Fig. 16. Shell displacement, case B, vs distance from base of shell (to cl. shell only) at time of peak displacement. 25OQ-
2000~ e $
1500.
.-d
60
.c
IOOO-
i 500
O-
-500-
Fig. 17. Moment at base of shell, case A, vs distance from base of shell. at time of maximum base moment.
-60
Fig.
2500
1
.50il\, CL. L-2
1.2
0.4
1
1
.s
Fig. IS. Shell displacement, case A, vs distance from base of shell, at time of peak displacement.
Fig. 18. Moment at base of shell, case B, vs distance from base of shell (to cl. shell only) at time of maximum base moment.
diagrams for Case A and Case B at the time when the shear at the base is near its maximum value. (7) Fluid pressure time histories. Figure 21 shows the time history of the induced fluid pressure adjacent to the shell at a distance of 220in. from the surface of the water (i.e. 22 in. up from the cantilevered end of the shell or opposite point W(2) of Case A). Figure 22 shows the sum of this time history and the forcing pressure (Fig. 3). This is the pressure in the fluid at this point independent of the hydrostatic pressure. If the tank were empty, the response curves of Figs. g-14 would look almost identical to the applied pressure shown in dotted lines on those figures. This is so because
31
Dynamic response of a containment vessel to fluid pressure puises
60. $
50.
= 1’
40. 30.
0. -IO-
/
---%
L-Z
Fig. 19. Shear at base of shell, case A, vs distance from base of shell, at time of maximum base shear.
-2.0 J
Fig. 22. Total (applied plus induced) fluid pressure at w(2), case A.
computer results, a very simple approximate method is also exhibited. The effects of cavitation and structural damping are not included. Since both of these are energy dissipating mechanisms, the response predicted using this procedure will be upper (or safe) bounds.
1. W. Fliigge, Static und Dynamik der Schalen. Edwards Brothers, Ann Arbor, Mich. (1943). 2. H. Lamb, Hydrodynamics. Dover Publications, New York (1932). 3. Mark II Containment Vessel Dynamic Forcing Function Information Report. Prepared by General Electric Company, Sargent and Lundy, NEDO-21061(Sept. 1975). L-Z
Fig. 20. Shear at base of shell, case B, vs distance from base of shell (to cl. shell only) at time of maximum base shear.
AFTKJWX-ANAWROXlMATg IWl’HOD In this appendix, an approximate procedure is used to obtain a one degree of freedom approximation to the shell wall response assuming the fluid to be incompressible. 1. Analysis Consider a simplified portion of the containment vessel shown in Fig. Al, which does not include the sloping portion of the configuration shown in Fig. 1, nor the rigid interior cylindrical surface representing the pedestal. Let
and assume
Fig. 21. Induced fluid pressure at w(2), case A.
the period of the empty shell (-T,=O.O16sec as obtained in Appendix A) is so small compared to the period of the excitation (0.1 set) that the shell would respond essentially statically. Thus the presence of the fluid results in both dynamic amplification and a phase shift between excitation and response, as exhibited in these time histories.
in which dots above the generalized coordinate q denote differentiation with respect to time. The fluid potential function of eqn (A3) satisfies Laplace’s
5. CONCL.USIONS
A method is presented, utilizing a digital computer, for determining the dynamic response of the wall of a nuclear containment vessel to fluid pressures exerted on the wall when the relief valve discharge piping is cleared. For preliminary design purposes, and to check the
I
1.
r3
_I
i
Fig. Al. Simplified shell and approximate response.
38
FRANK L.
DIMAGGIO et al.
equation governing the motion of an incompressible fluid
The kinetic energy of the shell mass is given by (A4)
(Al@
JA
L
which, using eqn (A2), becomes
the continuity condition
T, =; (2?rrJhp,)(i*.
(Al9)
In terms of the generalized shell mass, T, would be
and the rigid bottom condition
-- a4
Ir_o=o
al
T, = ; p$.
(A54
but the assumed she11 displacement does not satisfy the requirements that the shell bottom and top should have no radial displacement, i.e. w(0, f) f 0, w(l, t) # 0.
(A6)
(A20)
Therefore, from eqns (A19) and (A20). p = Zdhp,.
(-421)
The strain energy of the shell, which, because of the assumed displacement of eqn (A2) does not include bending, is
The induced fluid pressure is obtained from the potential as
LW
ad
P=P,
(A7) which using eqn (A2) becomes
which, using eqn (A3), gives
p=
f-6
(
W3)
> &
3
From eqn (A8), it is seen that the fluid pressure on the surface is not zero: Pk
1, t) f 0.
In terms of the generalized stiffness k, this strain energy would be
(A9)
The radial and axial components of the fluid velocity are obtained from the potential function as
From eqns (A23) and (A24), it follows that k = 2whEl
n(r, [, t) = - 2
u(r, l, t) = - 5
r3
(AlO)
(All)
The generalized force Qi corresponding to the incident fluid pressure pi is given by Q,(t)
which, using eqn (A3), give u,Lq r3 “_-2fd -
h
(A12)
.
=
I
A P&(P)
dA
LW
in which A is the shell surface and 1G_ is the spatial part of the assumed radial shell displacement. From eqn (A2)
(Al3)
*=1
(A27)
Qi = ZarJp,
L428)
so that Using eqns (A12) and (A13), the kinetic energy of the fluid in which the incident fluid pressure has been assumed constant over the shell depth, without the linear decrease to zero near the top shown in Fig. 4. Lagrange’s equations of motion for the shell may be written as
is obtained as
(P+m.)ii+kq=Q,. In terms of the virtual mass mu of fluid, the kinetic energy of the fluid would be written as Tf = ; m,cj*.
(‘416)
W9)
The square of circular frequency of the fluid filled shell is seen to be k (#Jz=p+m,
(A30)
while that of the empty shell is
From eqns (A15) and (A16), it follows that (Al7)
k wo== -. P
(A3l)
Dynamic response of a containment vessel to fluid pressure pulses The corresponding periods are
39
1.0
l-,211
(A=)
0
and
If the maximum incident fluid pressure p, is taken as unity, the static shell displacement it would cause would be wsr= 4.1 =
1 i;
(A34)
-1.5I Fig. A2. Solution of eqn (A35).
2.
Numerical
results
Assuming
lb Qi = 528Op,, pi in 2 /I = l&in. 1=2oft r, = 42 ft lb se? p=lHft p, =&=
lb sec2 14.97
Equations (A31)-(A34)yield o = 93 rad/sec o0 = 400rad/sec
The equation of motion, eqn (A29) becomes 4 + 85OOq
Equations (A17), (A21), (A25) and (A28) yield m. =
172,000!b$
lb se? a=%OOft k=
156x 107’b ft
7’= 0.063set T,,= 0.016sec.
= O.O28p,.
In Fig. A2, the computed response q/qs, of the shell to the incident pressure of Fig. 3 (superposed in dashed curves) is exhibited as a function of time. The dynamic amplilication factor is seen to be about 1.4. Considering the grossness of the approximation, this result, though not conservative when compared with the value of approximately 2 obtained from the exact analysis, is a useful measure, expecially in the preliminary design stage, of the magnitude of the dynamic amplification to be expected.