Transient failure analysis of liquid-filled shells Part I: Theory

Nuclear Engineering and Design 117 (1989) 107-140 North-Holland, Amsterdam

107

TRANSIENT FAILURE ANALYSIS OF LIQUID-FILLED SHELLS PART

I: T H E O R Y

Wing Kam LIU

* * *

and Rasim Aziz URAS

* * *

Northwestern University, Department of Mechanical Engineering, Evanston, Illinois 60208, USA Received 14 March 1989 In this paper, the governing equations which consider dynamic fluid-structure interaction, modal coupling in both axial and circumferential directions, and dynamic buckling are derived. The various pressure components acting on the shell wall due to a seismic event are also analyzed. The matrix equation of motion for liquid-filled shells is obtained through a Galerkin/Finite Element discretization procedure. The modal coupling among the various combinations of axial and circumferential modes are identified with a particular reference to the fluid-structure system under seismic excitation. Finally, the equations for the dynamic stability analysis of liquid-filled shells are presented.

1. I n t r o d u c t i o n

Although several important analytical and computational methods have been developed and analysis/ design procedures are now available for storage tanks, nevertheless, damage recorded in recent earthquakes has revealed a much more complex behavior than implied by design assumptions and some indicated important gaps in our understanding of the mechanisms involved in damage and failure. It is not clear whether this damage can be explained by phenomena associated with simply horizontal, vertical and rocking ground motions, or whether it requires a consideration of nonlinear fluid-structure interaction, or phenomena such as lift-off. Until recently, it was believed that only beam type (cos 0) mode was important in the analysis of the vibrational behavior of anchored storage tanks under seismic excitation. It has been shown by Shih and Babcock [15] in their buckling test of two plastic (Mylar) tank models that buckling is predominantly influenced by the stresses resulting from the lowest response mode and the higher order shell modes seem to have only a secondary role. However, shaking table experiments with aluminum tank models conducted by Clough [2], Niwa [13] and Clough and Niwa [3] showed that cos(n0)-type modes were significantly excited by earthquake-type motions, especially the cos(30) and cos(40) out-of-round shell modes. Clough and Babcock explained that the initial imperfections of the tank geometry are probably the reason for the appearance of out-of-round modes. These out-of-round modes have also been observed in the vibration tests of the full-scale thin shell water tanks conducted by Housner and Haroun [6]. In mathematical modelling, the task of identifying the loading mechanism onto an anchored liquid storage tanks which leads to these kinds of damage is also quite intriguing. Two mechanisms have been postulated. (1) Axial membrane stresses caused by vertical seismic excitation; this could lead to substantial periodic axial stresses acting on the shell and hence, the possibility of buckling. This type of dynamic buckling analysis has been studied by Kana and Craig [8], and Chiba and Tani [1]. In the Kana and Craig analysis, an uncoupled set of Mathieu's equations were obtained from the fluid-structure governing * This research is supported by National Science Foundation Grant No. CES-8614957. ** Professor. *** Graduate Student.

0029-5493/89/$03.50 © Elsevier Science Publishers B.V.

108

W..K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

equations and only hoop stresses were considered in the stability analysis. Whereas, in the Chiba and Tani analysis, the effect of axial modal coupling and all components of the membrane stresses are included, consequently, a set of coupled Mathieu's equations in the axial direction is obtained. Unfortunately, their studies are aimed at high frequency range (250-900 Hz), and hence, this type of stability analysis is not applicable to explaining the buckling of anchored liquid-filled tanks under seismic excitations. (2) Moments caused by high lateral pressures exerted by the fluid onto the upper part of the tank; this could lead to substantial compressive stresses at the bottom of one side of the tank, and hence, the possibility of buckling. This buckling mechanism has been studied by Liu and Lam [9], Liu and Lam [10] and Rammerstorfer et al. [14]. However, these analyses are limited to static buckling analyses. In order to truly understand the damage mechanisms of anchored liquid-filled tanks due to seismic excitations, a buckling analysis which considers dynamic fluid-structure interaction and modal coupling in both the axial and the circumferential directions is needed. The emphasis should also be on determining which phenomena are crucial to explaining the observed damage of tanks. From the experimental observations and previous studies, several interesting questions arose in developing design and computational methods for tanks in view of the above-mentioned damage. (1) Which of the loads (horizontal, vertical and rocking) that the tank sustains during seismic excitation is responsible for buckling damage? (2) Does sloshing in the broad tank contribute significantly to the mechanisms which cause buckling damage? (3) Is the behavior of the enclosed fluid of importance in correctly identifying the buckling loads and modes? (4) For failure conditions in which the buckling modes are not exactly known from the limited shaking table experiments, can the significant bifurcation solutions be identified (dynamic instability)? (5) Is buckling predominantly influenced by the stresses resulting from the lowest response mode (i.e. beam type (cos 0)mode)? If this is true, how can the cos(nO)-type modes, especially the cos(30) and cos(40) modes for a tall tank, be excited? The main objectives of these papers (Part I and Part II) are to develop the governing finite element equations for the dynamic buckling analysis of this type of fluid-structure interaction systems and to provide answers to the above questions. These papers should shed some light on the understanding of the buckling damage mechanisms of anchored storage tanks under seismic loads. This also promises to be the first comprehensive study in which fluid-structure interaction and the dynamic buckling analysis of anchored liquid-filled tanks are made. In the next section, the fluid equations are analyzed to identify the various pressure components acting on the shell wall. In Section 3, a variational formulation is introduced for fluid-structure interaction problems. The governing equations for a thin cylindrical shell are obtained in Section 4. Section 5 is devoted to the Galerkin/Finite Element discretization. Modal coupling for the resulting matrix equations is studied in Section 6. The dynamic stability equations are derived in Section 7. In Part II, applications to liquid storage tanks under seismic excitation are presented. Finally, some key answers to pertinent questions on buckling failure of anchored tanks are discussed.

2. The derivation of the total fluid pressure

The governing equation of a compressible inviscid fluid is the wave equation expressed in terms of a velocity potential, `#: 1 02,#

g7"2`# C2 i~/2'

(2.1)

109

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

Lineariz~dSloshingat the Fluid Surface:~

Cylindrical S~e'try C~nditionT

z a

I..~-~

~

Fluid-Tankwall Condition: GoverningEquation: 0..~= h(O,z,t) Or

~~~.~...~

~nO

~w,(~t)

=

n=-

~ ~ . . ~

I

v2,o7

(tankvibration)

a25

~

'L H

oo

+ gh(t)COSO

+ z ci(Ocos0

(horizontal pound excitation)

(rockingmotion)

0 c g gh

velocitypotential speedof sound gravitational~u~elea'afion ~tal grounddisplacement q angleof rotation w tankdisplacement W platedispl~ement

L

= f(r, 0,t) = ~W.(r,t) J nO

- r qCt)cos0 Fig. 1. Fluid equations and boundary conditions.

~h~e vibration) (rockingmotion)

where in cylindrical coordinates the Laplacian operator is given as V2

.

l a(ra . . r ar

~r

) .

1 a2 + r 2 ~0 2

a2 + --

~z 2"

(2.2)

r, 0 and z denote the radial, circumferential and axial coordinates, respectively, and z is measured from the bottom of the tank; c is the velocity of sound and t denotes time. In order to establish a well-posed boundary value problem, proper boundary conditions must be specified (refer to fig. 1): (i) Free surface boundary condition

0% + aq, at 2

g'~'Z = 0

at z = H,

(2.3)

where H is the height of the fluid inside the tank, and g is the gravitational acceleration. Eq. (2.3) represents the linearized sloshing behavior of the fluid surface. (ii) Cylindrical symmetry condition q, is finite or ~-~ r = 0

at r = 0,

(2.4)

(iii) Fluid-tank interface normal velocity compatibility

aep =h(O, z, t) ar

at r = R ,

where R is the radius of the fluid-structure interface.

(2.5)

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

110

(iv) Fluid-base plate velocity compatibility

aq,

O-"z-= f ( r ' 0, t)

at z = 0.

(2.6)

Moreover, all the initial conditions are assumed to be homogeneous. In seismic analysis, it is a common practice to assume harmonic behavior in the time domain, which can also be obtained directly from the application of the method of separation of variables. Hence, d~(r, 0, z, t ) = ~ ( r , 0, z, o~) e -j~t

(2.7)

where j = v/Z-1 and a superposed " - " denotes the frequency dependence for the functions. After substituting eq. (2.7) into the differential equation (2.1), it is transformed to the frequency domain (O~) 1 O r'~-r

r Or

1 O2t~ +----+--=

r 2 302

02~

OJ2 -

OZ2

C2

q~.

(2.8)

Equation (2.8) will be solved by employing a standard separation of variables technique, viz,

~(r, O, z, o~)= R ( r , w ) Z ( z ) O ( O ) .

(2.9)

The solution in the circumferential direction can be found as a Fourier series with the index n and the constant qn: O ( 0 ) = q, e j"°

n ~ (- ~,

+o¢).

(2.10)

Consequently, the governing equation, eq. (2.1), and the boundary conditions, eqs. (2.3)-(2.6), can be rewritten as

n 2 ~o2 --'~ +-~-+

1 d (rdR,) rR. dr~ ~

1 d2Z~

- 2 =0, Zn -dz

(2.11)

- o~2~,~+ gao-~ = 0

at z = H,

(2.12)

--=0

at r = 0,

(2.13)

at r = R,

(2.14)

at z = 0,

(2.15)

ar

Or

=L(z,

Oz =L(,, and

q,,(r, z, oo)= R , ( r , ¢o)Z~(z),

(2.16a)

where 1

r2,~-.

~n(r, z, to) = ~ Jo dp(r, O, z, o~) e-in°de.

(2.16b)

The frequency-dependent Fourier coefficients for h(0, z, t) and f(r, O, t) are denoted by h,(z, 00) and

f~(r, o~), respectively; and their expansions are given as +oo

h(O, z, t) =

E

h.(z, oo) e -J~'e j"°,

(2.17)

W.K. Liu, R.A. Uras/ Transientfailureanalysisof liquid-filledshellsI

111

and

-4-oo

f ( r , 0, t) =

~., L ( r , ,o) e -j~' e in'. tls

(2.18)

--CO

The last term in eq. (2.11) is a function of z only and therefore, it can be further separated as 1 dZZn Z, dz 2 = Yn,

(2.19)

The form of the solutions of eq. (2.19) depends upon the choice of the constant "t,. Three cases are possible: (a) v,, = O:

Zn( z ) =An z + Bn

(2.20)

and the radial equation becomes r ~rr t - - ~ r J

r2

c2 R n = 0

(2.21)

with the solution ~O

60

Rn(r, co)= EnJn(cr) + F n ~ ( c r ) ,

(2.22)

where Jn and Yn are the Bessel functions of the first and second kind, respectively; and according to eq. (2.13) F. --- 0 must hold since Y,(r) --* oo as r --* 0. (b) ~'n = -/-t2, /"t2 > 0:

Zn(z ) = Cn cos/~z + Dn sin/~z,

(2.23)

and eq. (2.9) becomes the Bessel's equation

/

(2 .24)

whose solution is given as

R~(r, ~) = Cjn(~r) + Hnrn(~r),

(2.25)

where I n and K n are the modified Bessel functions of the first and second kind, respectively; and ~2 __/~2 __ 602//C2.

(2.26)

Similarly, eq. (2.13) imposes H n - 0 since K n goes to infinity as r --, 0. (C) Vn=e 2, £ 2 > 0 :

Zn(z) = Pn cosh ez + Qn sinh ez,

(2.27)

and the modified Bessel's equation is

d ( r d R--~-r n ) + (~2r_ --;-)Rn=O, n2

(2.28)

W.K, Liu, R.A. Uras / Transientfailure analysis of liquid-filled shells I

112

whose solution is given as

R , ( r , to) = U,J,(~r) + VnYn(~r )

(2.29)

~2 = e2 + to2/cZ"

(2.30)

and

Due to the fact that Yn approaches infinity as r --* 0, eq. (2.13) dictates that V~= 0. These results can be summarized as for ~'n = 0:

a=Jn(-~r)[an(to)

H + bn(to)],

(2.31)

for Yn= _~2:

~2 = in(Ftr)[cn(to) cos(/~z) + dn(to ) sin(/~z)]

(2.32)

and for "in = e2:

~3 =jn(~r)[en(to) cosh(ez)+ gn(to)sinh(ez)],

(2.33)

where all integration constants are conveniently combined into six coefficients a,(00), b,(to) . . . . etc. The general solution for this problem can be obtained by superimposing all possible solutions, i.e., ~ n = ~ +~2 + ~

(2.34)

and applying the boundary conditions to obtain the unknown constants. However, this approach would give rise to a general eigenvalue problem which requires the determination of six coefficients through three boundary conditions eqs. (2.12), (2.14) and (2.15). Moreover, the three possible solutions, eqs. (2.31-33) would be coupled, and the response of the physical system to certain boundary conditions can not be easily identified. Instead, each ~k, k = 1, 2, 3 is determined by applying a particular set of boundary conditions such that the linear superposition of ~k would yield the general solution. The appropriate sets of boundary conditions are obtained by analyzing the natures of ~k, and exploiting the orthogonality relations of the sinusoidal and Bessel functions. This procedure is explained as follows: (a) Each ~ , k = 1, 2, 3, satisfies the free-surface condition, eq. (2.3), yielding

an [ -

to2

q_ g / H ] - to2bn =

0,

(2.35a)

G[ _to2 cos(ttH) _ #g sin(#H)] + d n [ _to2 sin(#H) + #g cos(#H)] = 0,

(2.35b)

en[--to 2 cosh(eH) + eg sinh(eH)] + g, [_to2 sinh(eH) + eg cosh(eH)] = 0 .

(2.35c)

(b) At the tank bottom, z = 0, ~ and ~ are chosen to satisfy ~z

= 0

for k - -

1, 2

(2.36a)

whereas

3z = f n ( r ' to).

(2.36b)

After applying the condition (2.36a) on eqs. (2.31) and (2.32) an(to ) and d,(to) are determined as a n( to ) = 0, d,(to) = 0, respectively.

(2.37) (2.38)

W.K. Lit~ R.A. Uras / Transient failure analysis of liquid-filled shells I

113

Substitution of eqs. (2.37) and (2.38) into eqs. (2.35a) and (2.35b) yields

b,(to) = 0,

(2.39)

and -to 2 cos(Pill) =

#ig sin(~tiH) for i = 1, 2 . . . .

(2.40)

respectively. As a consequence of eqs. (2.37) and (2.39) ~, vanishes. Before applying the nonhomogeneous boundary condition on ~3, eq. (2.36b), the eigenvalue ~ needs to be determined. (c) At the tank wall, r = R, ~k, k = 1, 2, 3 are assumed to satisfy the following conditions:

0K Or = O,

(2.41a)

-2

Oq~, v, w), Or~=n.~,z,

(2.41b)

and ar

-- 0.

(2.41c)

Since ~ vanishes, eq. (2.41a) is automatically satisfied, and eq. (2.41c) yields

J,'('~,R) = 0

i= 1, 2, ...,

(2.42)

where [iR are the zeros of eq. (2.42) and a superposed prime denotes the derivative with rcspect to the argument. (d) The coefficient c,(to) is evaluated using eq. (2.41b) together with eqs. (2.32), (2.38) and the following orthogonality property of the cosine-functions

foHcos(I.t,Z)COS(pjZ) d z = ½ ( H - ~ 2 2 Sin2(piH))8,j.

(2.43)

where 8u is the Kronecker delta. The coefficient is

2 fonh,(z, w) c.(o~) = [ H---d~gsin2(P~'H)]giI:(~iR )

cos(/~,z) dz.

(2.44)

(e) Employing eq. (2.36b) and the orthogonality condition for the Bessel functions R2

foRrJ. ( ~,r ) J. ( ~jr ) dr g.(w) is found

(1 n2) ]

(eiR) 2 J.2(~,R) 8ij ,

J'2(~,R)+

(2.45)

to be

2

gn( w ) = eiR2( 1 _ n2/( ~iR )2) j~( ~iR ) f:fn( r, w) Jn( ~ir )r dr.

(2.46)

This expression is substituted into eq. (2.35c) to obtain e,(o~) = -

- ,02 sk

(en) + es cosh(

H)

- w 2 cosh(eH) + eg sinh(eH)

g,(w).

(2.47)

W.K. Liu, R.A. Uras/ Transientfailureanalysisof liquid-filledshellsI

114

Finally, the general solution is given as ~ , ( r , z , to)=¢~2+~ 3,

(2.48)

where ~2 is due to the tank wall motion

oo

in(Ftir)

~ 2 = 2 i=lE [ n---~g sin2(IzgH)]~iI'(~i R )

C0S(/~iz)f0Hh,,(~, to) cos(/.ti~) d~

(2.49a)

and ~3 is due to the base plate motion 2

[sinh(eiz) - cosh(egZ)

"R i=l L

to: sinh( eiH ) - egg cosh( eiH ) to2 cosh(eiH ) -eig sinh(eiH)

J.(~,r)

×

fRL(r, to)S,(~ir)r dr. ao

(2.49b)

The fluid dynamic pressure, Pf, is the pressure induced by the fluid flow and is given by Pf = - oF-if 7 .

(2.50)

In the frequency domain it is rewritten as

ff f = JtoPF~,

(2.51a)

and additionally it can be expanded in Fourier series +oo

Pf=

E

(Pf,) eJ"°.

(2.51b)

Since ~ is given by eq. (2.48), the nth harmonic pressure component becomes

oo I,(~ir ) cos(#iz) H Pr"=2OF~=, [ H _ _ Jg sin2(l~iH)]~,i~(~i R) fo jtoh.(~, to) cos(/~/~) d~" _

o~2 sinh(eiH ) - eig cosh(e,H) egg sinh(e,H)

+ -2PF -~

sinh(e,z) - cosh( e,z ) to: cosh( e,H )

i=l

J,(~,r) ×

-

fo~tof,( r, to )Jn( ~ir )r dr.

(2.52)

To identify the various components of the fluid pressure induced by the prevailing boundary conditions at the tank wall and at the base plate, the proper forms of h(O, z, t) and f(r, O, t) will be discussed. A summary of different boundary excitations is presented in table 1. The tank wall boundary condition can be interpreted as ]~(0, z, /)--- f t ( 0 , z, t ) = f + f h + f r m , where f t is the total tank wall acceleration, and f , f h and firm are defined below:

(2.53)

W..K Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

115

Table 1 Boundary conditions for the fluid equation

structural vibration due to boundary excitation

rigid Ixxly morion due to

relative deformation horizontal¢xcitatino rocking motion

at the tank wall II(e,Z,t)

~ = ~ ~n(Z.t) eJno

Gh(t) cosO

Zq(t) case

I1=-~

at the base plate f(r,O,t)

W= ~g'.(r,t) d no

- r q(t) cos0

n==-~

w shell defonr~fion W plate defocmatino Gh horizontal ground acceleration q angle of tank rotation A supcl~os~ "." denotes the time derivative.

(a) The relative wall acceleration is expressed as +oo

fi(O, z, t)=

~'. fi,(z, t) e j"e, rl B

(2.54)

--O0

where w is the tank radial displacement. (b) The rigid body acceleration due to horizontal seismic excitation is given as

fih(O, Z, t) = Gh(t ) COS 0,

(2.55)

where Gh(t ) is the horizontal ground acceleration. It should be noted that only cos 0 modes are excited if the tank does not have any geometric imperfections. (c) When the tank (with the base plate) is subjected to a rocking motion, a rigid body acceleration is induced by the rigid rotation of the structure: firm(a, Z, t) = Z#(t) COS a,

(2.56)

where q(t) is the angle of rotation. Here, the absence of any geometric imperfection is assumed again. The base plate boundary condition can be treated similarly

f(r, a, t) = l~t(r, a, t ) - - w + l~rm,

(2.57)

where ff.t is the total base plate acceleration. In general it may consist of two components: (d) The relative plate acceleration is given by +oo

l~(r, O, t) =

E

l~?~(r, t) e j"°,

(2.58)

where W is the transverse plate deformation. (e) The rigid body motion component due to the rocking motion is l~rrm(r,

a,

t) -- -r~/(t) cos 0,

(2.59)

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

116

In the frequency domain the general form of the tank wall and the base plate boundary excitations can be rewritten as

-josh(O, z, (o)ffi -to 2 ~_, ~ , ( z , ~o) eJn0q- ah(0)) COS 0--Zo)2q(03) COS 0,

(2.60)

/'/~ --00

and +oo

-joaf(r, O, 60) = -(o 2 ~.,

W,(r, to)

eJn°q-rw2q(6o) cos

0.

(2.61)

The pressure induced by the fluid flow, Pf, has been introduced in eq. (2.52). After substituting eqs. (2.60) and (2.61) into eq. (2.52), Pf is decomposed into four components for a better understanding of the physics of the problem:

Pf ---- ?d q- ?b q- ?h -t- P r m ,

(2.62)

where Pd and Pb denote the shell and base plate vibrational pressures, Ph and Prm are the pressures due to horizontal ground excitation and rocking motion under the rigid tank assumption, respectively. The total pressure in the fluid is determined by superimposing the pressure due to the induced fluid motion, Pf, with the pressure due to vertical ground acceleration, Pv, and the hydrostatic pressure, PL: P ffi Pd + Pb + Ph + Prm + Pv + PL,

(2.63)

where

Pv = --PFGv(t)[ H - z],

(2.64)

eL ffi OFg[H-- zl

(2.65)

and

with Gv(t) and g being the vertical ground and the gravitational accelerations, respectively. The expressions for the pressure components are summarized in table 2. Since the pressure components listed in table 2 are functions of the frequency ,0, it requires a rather complicated iterative solution technique. For most liquids, such as water, the compressibility effect can be shown to be negligible. By analyzing the relative magnitude of the nondimensional parameters ~oH/c and ,oR/c (refer to eqs. (2.26) and (2.30)) estimates of the limits of this assumption can be established. Recall eq. (2.26) in the following form: (#iA) 2 = (p,A) 2 - M 2,

i = 1, 2 . . . . .

(2.66)

where M A = oJA/¢,

A ~ H

or R.

(2.67)

In order to neglect the effect of the compressibility, M A must be small compared to (#~A) in eq. (2.66). For example, if M A is chosen to be less than 10% of the lowest axial eigenvalue, (/~IA), then (~A) are within 0.5% of (/~IA) for all values of i. For water-filledtanks the fluid-structureinteraction frequencies typically range from 2 to 40 Hz. Employing this knowledge together with the speed of sound in the fluid the range for H and R can be determined. A similar analysis on eq. (2.30) which is rearranged as (~iA) 2 -- (eiA) 2 + M 2,

i = 1, 2 . . . . .

(2.68)

W..K. Lit~ R.A. Uras / Transientfailure analysis of liquid-filled shells I

Table 2 Fluid pressurecomponentsdue to variousboundaryconditions Expression

Pressure component due to Symbol

vibration of the shell

Pd

~

20~2 Fin(r,z,t) ein0

nffi-m H

Table 3 Fluid pressurecomponentsin the absenceof sloshing Pressure component due to

Expression

Symbol

¢ibfatioftof the shen

Pd

"n+ffi~."~1 2 Fin(r,z) O"e H *0 [. ~v.(~.O cm(~.,O d~;

vibrationOfthe plate

Pb

"n+~ ~1 ~Qm(r,z) e~s

* j fvn(~,m) cos(~ti~)d~ Pb

vibration of the plate

~ flffi-~

Qm(r,z,t)e}nO =

117

R

R

* J W.(~,t)J.(~) ~ d~ hodmnmlgroundrnofion

Ph

- ~. 2 H ~ i h ( m ) ~ c o s 0 l=~= i,ptH)

mofice

P~

~ 2¢o2H2~(m)cos0

P,

[.[piH(-I) i+l- l]Fil(r,z,t) R2Qil(r,z,t)J2(~'iR)] (l~iH)2 H2(~'IR) J - PF G,(t) [ H - z ]

PL

ppg[H-z]

rocking

vertical ground motion weight Of the fluid

where

Fin(r,z,t) = PF

In(~ir)cos(p.iz) e'Jet [I - ~--~Hsin2(~iH)](~iH) l'n(~iR)

....... ~2siah(~iH) " ~ig cosh(~iH)] Qin(r,z,0 =/smnt¢iz; - cosnt¢iz? / L ¢a2cosh(l[iH)- ¢ig sinh(EiH)J ,

horizontal ground motion

[rocking rnodon

oo

-| i+lF r z 2HGh(t)~ c o s 0

Ph

-

Prm

"l=~ 2H2q(t) cos0 .~[p,iH(-I)i+t- l]Fiffr,z) R2Qil(r.z)J2(£1R)"

motion J~ei~,htof the fluid where Fin(r,z) ffipF~ ver~cui Found

(plrl) tn~'R)

P,

PL

"pp G,(t) [ H - z ] pvF(H-z]

COS(~iZ)

prJn(¢ir) Qin(r,z.t)= [sinh(f.lz)- cosh(¢iz)tanh(¢iH)] (EIR~|" ~ ) J2(eiR)

ppJn(~'ir)e"j'u~

(ciRri - ~ )

J n2(E'iR)

where MA = w A / c ,

A ffi H or R

(2.69)

would yield restrictions on H and R. The values of H and R determined through this dimensional analysis are rather conservative since the pressure terms are relatively insensitive to the addition of M R and MH. When the compressibility is neglected, the following changes take place: ~,~Pi,

(2.70)

ei --' ei-

(2.71)

and

The sloshing frequencies constitute the lowest range of the Fluid-Structure Interaction (FSI) spectrum. For the fluid-structure vibration frequencies (2-40 Hz) the contribution from sloshing is insignificant. A closer investigation of the terms in eq. (2.52) reveals the fact that this contribution is mainly governed by the non-dimensional parameter -- oj2H/g.

(2.72)

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

118

In the range of the fluid depth encountered in practical applications the value for t3 ranges from 500 to 1000, and therefore, the effect of sloshing can be neglected. Under this assumption the roots of eq. (2.40), viz., (/~H) can be approximated by

~iH

2i-1 2

[

~r 1 +

1]

~"

2i-1 2 lr.

(2.73)

It should be pointed out that the natural frequencies depend upon the depth of the liquid in the tank. If the fluid depth is small, then the fluid-structure natural frequencies are considerably higher. Consequently, t3 generally remains within the limits mentioned earlier. With these approximations the pressure terms presented in table 2 can be simplified greatly. The corresponding expressions in the time domain are given in table 3. Finally, the total pressure acting upon the tank wall is evaluated by substituting r = R into the total pressure expression, eq. (2.63):

Psh(O, z, t)= P(r--R, O, z, t).

(2.74)

Similarly, the total pressure acting on the base plate is determined by considering

Ppl(r, O, t)---P(r, O, z = O, t).

(2.75)

3. Variational formulation

The principle of virtual work statement for a structure subjected to hydrodynamic loads due to an inviscid, incompressible fluid can be expressed as

8W+ST-frs.,t'.(a,x,t),,,dr-frs.,t'F(x,t),,,dr-fvs.,t'L(x),,,dr=O,

(3.1)

where r , n i, xi and u~ denote the fluid-structure interface boundary, the outward normal to the structural surface, the spatial coordinates and the structural displacement components, respectively. A superposed dot designates the temporal derivative. The pressures acting on the tank wall are tabulated in table 4, and PF --

Ph

"{-Pv + Prm"

(3.2)

The internal work of the structure is given as

8V= f~ 8ui.j %.ydJ~,

(3.3)

where f~ is the domain of the structure; ~-~jis the Cauchy stress tensor. The spatial derivatives are denoted by a subscripted "," and repeated indices indicate sum. The structural inertia is •T--

fop

a, dJ2,

where p is the structural mass density.

(3.4)

W.K. Liu, R.A. Uras / Transientfailure analysis of liquid-filled shells I

119

Table 4 Expressions of pressures acting upon the shell wall Pressure component due

to

~brafionofthcshell

Symbol

Pd

Expression

X

~ Fin(z) cosn0J Ga(~,t) cosO.ti~)d~

n=0 i=l ~zontal groundmotion

ph

rocking modon

P~

l_~lHGh(t)( ~ cos0 = LP.IH)

~

[~.iH('l)i+ | " 1]Fil(z). R2QiI(z)J2(~iR)] (v.iH)2 (eiR)Hz .]

vertical ~rtmnd motion ~vtfightof the fluid

P, PL

- pF Gv(t) [ H - z ]

pvg[H-z]

where

F~(z)= ~ c o s ( I J i z ) (P-iH)In(v.IR ) Qit(z)= 2[sinh(elz)-cosh(eiz)tanh(~iH)]

pFll(giR)


After expanding ,%, ui, Pd, PF, and n, into the zero-th and first order parts via a consistent linearization procedure, eq. (3.1) becomes

f.8.,,,[~,, +,%] d. + f: 8"`i[u,+ A~,]dS~- fr8"`i[I~d(x't)+ l:d(x,t ) ] .

i dF

-fr 8"`,P~(~, t),,dr- fr8,,, ?L(x)n,dr-fr ~'`,P~(~, ,) ,~n,dr - fr 8u, PL(X) A n i dF = O,

(3.5)

where the added fluid inertia and the surface normal are defined as

I~d(x, t)--- aPd(~, x, t)..

a~j

I('(x,

t)=- aPd(ii' x, t)

~iij

"`J' Ai~j,

(3.6a)

(3.6b)

and Ani =--AUk,k ni -- Aum,i rim,

(3.6c)

respectively. Note that the change in the surface traction components is neglected in eq. (3.5) as the pressure components are derived in terms of the initial configuration.

120

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

The linearized Cauchy stress tensor, Arij, can further be separated into two parts representing the material response and the initial stress effect denoted by superscripted "m" and "0", respectively:

Arij= A~i7 + A~i°.

(3.7a)

The components are defined as

Ari7 = C~tjkmAu~,,,

(3.7b)

AriO=sik ,jo AUk, m

(3.7C)

and where C~km is the material response tensor and 8~k is the Kronecker delta. After substituting eqs. (3.7) into eq. (3.5), the zeroth-order and the first-order equations can then be identified:

zeroth -order ~pSuiiJidI2-~SuiI~3d(x,t)nidF+~Sui,j

*ijdl2

= fr nu, eF x, t)n, dF + fr *u, PL(X)FIi dr.

(3.8)

first-order

f o'u, aU,da--fr*',Z¢"(x,t)n, - f S u , [PF(x, t) + PL(x)] Ani d F = 0.

(3.9)

"F

4. Reduction to the thin shell theory

Preliminaries: assumptions (refer to fig. 2): (i) The shell thickness is small compared to the radii of curvature. (ii) The deformations are small compared to the shell thickness.

~ / ~ ~ ~ u

,

1, U

~w

.~

J

T

112, V

ru3w

L

1

midsurfacc

Fig. 2. Shell geometry and coordinate system.

W.K. I_du,R.A. Uras / Transient failure analysis of liquid-filled shells I

121

(iii) Planes normal to the midsurface remain normal and unchanged in length after the deformation. This implies that the shear strains are negligible, and the normal strain in the radial direction can be omitted. (iv) The normal stress in the radial direction, Or, is neglected. The components of the displacements at a distance r from the midsurface are defined as l l l ( r , 0, Z, t ) = U(Z, O, t ) -

0w r-~Z (Z, 0, t),

(4.1a)

r 0w u2(r, 0, z, l) = v(z, O, t) - --~ O'--ff(z, 0, t),

(4.1b)

u3(r, O, z, t ) = w ( z , O, t).

(4.1c)

In the case of thin shells the plane stress assumption is introduced, so that the stress and the strain tensors reduce to

T = [otz, (io, ozo, or = O, OrO, arz]T

(4.2)

' = [%, 'o, 2"z#, ~'r, 2,,0 ----O, 2,rz = 0] T.

(4.3)

and

In cylindrical coordinates the domain and the boundary integrals become

fa d~2 = f h / 2 fr d F =

f2~r f L R dz d0 dr,

"0

"0

(4.4a)

~0

~ h~2~0

R dz dO,

(4.4b)

where R, h and L are the radius, the thickness and the length of the cylinder, respectively; H is the height of the fluid. After eqs. (4.1) are substituted into eq. (3.9) and integrated in the radial direction, the individual terms for the linearized shell variational equation are illustrated separately: (i) The structural inertia term can be divided into the translational and rational parts:

f p 8U i A~ i d[2 =

It +

(4.5a)

Ir,

where (4.5b)

I t = fAoh(Su Aii+ 8v A6 + 8w A ~ ) d A and I r=

--~

8w,z A f f , z + - ~

8w.oAf¢,o

dA.

(4.5c)

A denotes the area of the midsurface, and a subscripted "," denotes the partial derivative. In the subsequent analysis, the higher order rotational inertia, Ir, is neglected. (ii) Since the surface normal vector to the structure at the interface has only one non-zero component, i.e. n = [0, 0, 1] T

(4.6)

W.K. Liu, R.A. Uras / Transientfailure analysis of liquid-filled shells I

122

and since the dynamic pressure is a function of the radial acceleration, the added mass expression becomes

at'~(~, x, t)

~'?"(~, t ) -

a¢0

A#

(4.7)

and hence

dr = £ 8~

; 8u i I ~ a ( x , t ) n ,

aPd(ii', x, t)

za~ d r .

(4.8)

(iii) With the plane stress assumption, the material response term for a linear elastic material can be given in terms of the total strains, % Im----

£~U(i,j) citjkm AU(k,m) dO=- £~Cij citjkm A~'km dl?,

(4.9)

where ~, = e, - rK,,

~o = eo - rKo,

2c,o = 2e~0 - 2rK~o.

"e" denote the midsurface strains; K z and Ko are the changes of curvature of the midsurface and Kzo is the twist of the rnidsurface. The strain-displacement and the curvature-displacement relationships for a circular cylindrical shell is given by Flugge [5]. /A

1 eo = - ~ ( v o + w ) ,

let

"~

(4.10b)

1

2e,o = -~u,o + v,,,

(4.10c)

K, = w,~,

(4.11a)

and

Ko

= ~2(woo + w),

2K~o =

(4.11b) (4.11c)

2w,,o + V z - -~u,o ,

respectively. Through a standard reduction procedure applied in plane-stress problems, eq. (4.9) can be rewritten as Ira= - - - ~ 1

[seTcAe--r(SeT CAK+SKT

C A e ) + r 2 8KT C A K ]

d~,

(4.12a)

where e = [e,, eo, 2e, o] v,

(4.12b)

K = [ K , , K o, 2K, o] T,

(4.12c)

o]

and C=

1

o

0

(I

-

~)/2

123

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

E and v are the Young's modulus and the Poisson's ratio, respectively. After integration in the radial direction eq. (4.12) reduces to

h~eTCAe+-~KTCAK

Ira-- _ " " ~

1

dA.

(4.12d)

(iv) Assuming that the initial state of stress is governed by the membrane stresses and neglecting the higher order curvature terms, the initial stress term can be written as

fa 8ui'j ~o Aui, m d ~ = ~ 8~TN ° A~ dA,

(4.13a)

where

[

1 1 1 ]T = u,~ v,~ w~ -~u,e -~vo -~wo

(4.13b)

and

No = [ N~°l [N°I

N° I ] N°I J"

(4.13c)

I is a 3 × 3 identity matrix and the initial membrane stresses are defined as

N°(z, O, t) = haa(z, O, l),

a = z, g or (z/9).

(4.13d)

In deriving eqs. (4.13), ui, 1 = ui,~; ui,e = (1/R)ui. o and ui,3 -- ui,r ~ 0 have been employed. (v) The load correction term is (eq. (3.9))

In= fr 8ui [PF(X' t) + PL(X)] An i dF,

(4.14)

where In arises due to a change in the normal, after the structure has undergone a deformation, and (4.15a)

An 1 ~ --Aw,z ,

1

(4.15b)

An 2 = - - - ~ A w , o ,

1

1

(4.15c)

An 3 = Au,z + --~ Av,a + --~ A w .

The total pressure acting upon the shell wall can be redefined as

PF +

PL = PL(Z) +

Pv(z, t) + ea(z, t)

(4.16)

c o s O --- e w ,

where P1 stands for the 0-independent part of Ph and Prm" After substituting eqs. (4.15) and (4.16) into eq. (4.14) the load correction term becomes

In=

fr{- ~ u e W A w z - - -8O_w ~ e A w , a + S w e w( A u z + - ~

1 Av,o+-~1 Aw )} dF.

(4.17)

Although this expression appears to be non-symmetric, it is possible to symmetrize it through integration by parts of the following two terms: 8upWAw,~-----(~upWAw)~, - S u z, pWAw_BupW,zAw

(4.18a)

124

W.K. Liu, R.A. Uras / Transientfailure analysis of liquid-filled shells I

and 8v pW Aw,o = ( 8v Pw A w ) , o - S v

pW A w - S v

(4.18b)

pW,o A w .

Both the last terms in eq. (4.18a) and (4.18b) vanish since the total pressure is determined at the initial configuration. After the substitution of eqs. (4.18) into eq. (4.17) and rearranging terms yields I n =

fr[suz

P W A w + 8w e W a u , z ] d E

l fr[ 8v,o pW A w + Sw pW Av o] d E + ½fr 8w pW Z~w d I ' -

+ -~

I z - Ie,

(4.19)

where Iz and I0 are line integrals in z and 0 directions, respectively; and they are given by Io=

dO,

8u

(4.20)

2~r

I t = ~1 fo H 80 pw A w

dz.

(4.21)

0 I 0 vanishes if the load is conservative, i.e. the pressure at the free surface is zero PW(z=H,

O, t ) = 0 .

If either A w ( z = O, 0, t) = 0

(cantilevered shell at the bottom)

of the test function 8 u ( z = O, O, t) = O,

Iz drops out due to the singie-valuedness of all functions appearing in the equations 8v pW A w l o = o _ 8v pW A w l o = 2 ~ = O.

Finally, eq. (4.14) can be rewritten in matrix form I n = f#~e*rPe * dr,

(4.22)

where e* =

[(U z + - ~1v , o) 1, w •,

p _ _ [ 0i

1 1 1/R

(4.23a)

(4.23b)

(vi) For the zeroth order equation, the structural inertia and the fluid added inertia are the same as those given in eqs. (4.5) and (4.7), respectively; however, with the "'A" sign removed. The internal work term is

8Wint= f 8ui, j "rij dI2. •t o

(4.24)

W.K. Liu, R.A. Uras / Transientfailure analysis of liquid-filled shells I

125

Since the interface force acts only in the radial direction, the external work term can be rewritten as ~W ex'=

fJF ~u, pWn, d r =

f ~w pw d r . Jr

(4.25)

5. Galerkin/finite element discretization For the shell structure investigated in this paper, the displacement components can be expanded in terms of the mode shapes in the coordinate directions I

N

u(z, 0, t)= ~., ~_, ui,(t)Ui(z)Un(O),

(5.1a)

i=l n=O I N

v(z, 0, t)= ~, Y'. vi,(t)Vi(z)V,(O),

(5.1b)

i=1 n=0 1

N

w(z, 0, t)= ~., ~.. w,,(t)W~(z)W,(O),

(5.1c)

i=1 n=0

where I and N are the total numbers of the terms taken in the axial and circumferential directions, respectively. The discretization of the terms introduced in the previous section is performed in the same order below. Repeated indices denote sum unless otherwise stated. The matrices defined below are given in table 5. (i) The structural mass matrix:

i..) L o h 8d~,N~,No,,No,,,N~j , T • Adj,. dA

(5.2)

Consequently, the structural mass matrix is defined as

(5.3)

nijnm-~:J)phRNzT[fo2~trNoTNomdO]Nzjdz. (ii) The added mass matrix: From eq. (4.7) and table 4: q-oo

oo

I~d(x, t)= -- E E F,q(R)7iqCOS(l~qz)

cos(n0)

Aibin(t),

(5.4)

n~O q=l where

2PvI,(IXqR ) Fnq ( R ) ' ~ (Id,q n ) i ~ ( ~ q R ) '

v~=

fo"w~(~) cos(~,~) d~,

(5.5) (5.6)

and let ~,~ = f0~"W~(O) cos(nO) dO.

(5.7)

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells 1

126

Table 5 Explicit forms of the matrices employed in the discretization din = [ Bin Vin Win ]T a=0. b=n

Nab

[J~ 0

0

0 V~

0

0

Un

Blb

B.~b

0

0

0

1

0

~'CCn,o 0 ~ W ~

1 -

~" U.,o

0

I

0 "~,,

0

0

a=z. b=i

0

~¢n

0

.]

Wn

0

0

0

0

0

W~+W~.oo

0

0

0

0

0

R-'r

[- Ui.z 0

0

[

0

Ui

0 0 L 0

0

"

Vi 0 Vi,~ 0 0 Wi

F- 0 O WL=" 0 0 Wi Ui 0 0 0 Vi.z 0 L 0 0 Wi.z

" 2g-'~- Yg - I U "

Un 0 0 0 0 V~ 0 0 0 0 W~ 0

0 0 0

0 0 0

o

o

o U..o o

0

0

0

0

Vn, o

0

_ 0

0

0

0

0

W.,o.

0,

B:b

Ui0 0 ] 0 Vi 0 0 0 Wi

1

~-V,.0

0

0

~'n

0

[

o !

I

-ROiz 0 0 " 0 ' RViz 0 1 0 0 "RWiz ~" Ui 0 O' 0 Vi 0 -

0

0

Wi .

Euooj 0

Vi 0

0

0 Wi

S u b s t i t u t i n g eqs. (5.4) t h r o u g h (5.7) i n t o eq. (4.8) y i e l d s the a d d e d m a s s matrix ad

Fnq ( R )'yiq'Yjq6OnmA ,

Mijnm =

(5.8)

q=l

where

[o o] 0

A=

0 0 0

0

•

(5.9)

T h e total m a s s matrix is the s u m of the structural a n d the a d d e d m a s s matrices

m, jom ~- M,%m-~ nijnm" ad

(5.10)

(iii) T h e material stiffness matrix: T h e strain a n d curvature vectors, eqs. (4.12b) a n d (4.12c), respectively, are discretized as e = B~. B ) , d ,. ,

(5.11)

K=B~.B~di.

(5.12)

W..K. Liu, R.A. Uras / Transientfailure analysis of liquid-filledshells I

127

so that I m in eq. (4.12d) becomes Ira----

1

---TfAEh ~din(nzi1 T (ndn)1 TC(B~m)(Blzj)I

dA Adjm

Eh 3

+ 12(7--- v 2) fa 8d~ (B~)a" (B~n)TC(B~")(B~)

dA Adjm.

(5.13)

The material stiffness matrix is defined as D

ehR

KiJn'~. . . .1 - p2 +

jo

(B 1o) T c(B

Eh3R

120 Z ~,2)

.)dO

L 2 T 2~r 2 T 2 fo (Bzi) [fo (B~n) C(B°m) dO] (B2")

dz.

(5.14)

(iv) The geometric stiffness matrix: The discretization of eq. (4.13) leads to L

t T

2v"

t T

0

t

p

K~nm= fo (B~i) [fo (B~) N (B~,~)dOl(B~j)Rdz.

(5.15)

The initial membrane force matrix, N 0, is composed of a hydrostatic part, N s, and a dynamic part, N d, (due to vertical, horizontal and rocking ground motions):

N°=N~(z)

+ N d ( z , 0, t),

(5.16)

where

NS=

(5.17)

[ 03×3

03x3

t 0 ×3

p gR(n-z)i ×3

and

N d = No(z, t) + Nl(z, 0, t).

(5.18)

NO includes only the axial'-and' hoop forces induced by the vertical ground motion and is therefore independent of 0:

03×3]

No(z, t)= [ N°(z' t)I3x3 t

03x3

N°o(Z'

t)13x3

(5.19) :

The matrix, N1, includes the initial membrane forces induced by the horizontal and rocking motion. The resulting axial, hoop and torsional forces are denoted by N°, N° and N°l, respectively. Furthermore, the membrane forces in the coordinate directions can be separated from the torsional component: N 1 = ~l(0)N~0(z, t) + ~2(e)N~l(z, t),

(5.20)

where

= [ N°(z' NI°

L

t)13×3

03x3

03×3

N°(z'

]

t)13×3 '

(5.21a)

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

128 and

[ 03×3 UA(z, t)I3×3] N,I= LNOl(z' t)~3×3 03×3 ]

(5.21b)

The choices for i l l ( 0 ) and fi2(0) will be explained in the next section. After substituting eqs. (5.17)-(5.21) into eq. (5.15) and rearranging terms the overall geometric stiffness matrix is obtained (5.22)

g ~ n m __ Gs GO 4- K G 1 0 4- K G l x -- Kijnm --I-gijnm __ --iJnm -- --ijnm"

The physical significance of each term is identified below: (a) hydrostatic part

K,j,,~ =

Bz,

B~,) N (B~,,) dO B~j)R dz.

(5.23)

(b) vertical ground excitation part G0 L t T 2, t T t x~j.~=f0 (n~,)[f0 (n~.)No(n~.)

dO](azj)R

dz.

_ '

(5.24)

(c) horizontal ground excitation and rocking motion parts (i) due to the induced axial and circumferential membrane forces ijnm

(ii) due to the induced torsional membrane forces

,:,,,=

B~,

Bd,) Nll(Bd,.)u2(O ) dO B,j)R dz.

(5.26)

(v) The load correction matrix is composed of three parts, a hydrostatic pressure component, a pressure component due to vertical ground motion, and a third pressure component due to horizontal and rocking ground motions. They are designated by subscripts "Ps", "P0" and "PI", respectively. Eq. (4.19) yields the load correction matrix after the discretization (refer to table 5): X/Pj,m = Kijnm Ps + Kijnm PO + Kijnm, m

(5.27)

where

K,j,,. p~ = f0Lp v g ( H - z ) ( B ~ . ) v [ f 0 2 ,(Bo,) , . TP(B~m . ) dO ](B,~)Rdz, --,j,m

"Io _~(z, t)(B~*) v

[:o

(5.28a)

]

(5.28b)

cos O(B~)Tp(Bo *) dO (Bz~)R dz.

(5.28c)

2,~(Be,) . v P( Be,.) . dO ( B,~ )R dz,

and

Kij, _-

Pl(z, t)(Bzi )

(vi) From eq. (4.24) the internal force vector becomes

/,:'= f

dO+ f

* dO,

(5.29)

W.K Lib R.A. Uras / Transient failure analysis of liquid-filled shells I

129

where T* = [az, 0t0, Oz0] T.

The external force vector follows from exl. (4.25)

= fr , [eF(z, o, t) + PL(Z)]

dE,

(5.30)

where

a,~=[o,o,w~(z)W.(O)]T The governing matrix equations can be summarized as follows: the zeroth-ordcr equation is:

M,j.,,.d),. + ~i,,m=~,,,

(5.31)

and the first-order equation is:

Mijnm Ad)m + gijnm Adjm + gi~,m(t ) Adj~ + KPnm(t) Adjm=O.

(5.32)

The hydrostatic part of the fluid pressure acting upon the interface has a significant contribution to the overall response of the fluid-structure system. In order to account for this, the static part of the geometric stiffness and the load correction matrices are combined with the material stiffness to yield K i j n m = KiDnm ..{_Kijnm Gs + K i jPs nm.

(5.33)

Similarly, the overall geometric stiffness and load correction matrices can be cast into the following forms: o

ICijnm

= g--ijnm ° ° +go10 goll -- --ijnm + --ijnm

(5.34)

and x~.,. = --tJnm xr.o -. --tjnm" xr.~

(5.35)

6. Modal coupling of the fluid-stnmm~ ma~ix equations

To exploit the modal deformation patterns a Fourier series is taken in the circumferential direction whereas a standard beam deformation pattern is used in the axial direction. After ufflizing the symmetry conditions for cylindrical shell in the circumferential direction, the functions are given by (Donnell [4]).

V,(z) =,/,,(z).

~.(o) = ~ s ( , , o ) ,

(6.1a)

E ( z ) = ,,(~),

F,,(o) = sin(,,o),

(6.1b)

W,(z) = ¢,(~),

~ . ( 0 ) = cos(,0).

(6.1c)

The axial functions are chosen according to the boundary conditions: (no sum on i) cki(z) --AI~ cosh

h~z +A2i cos ~kiz +A3i sirth ~kiz + A 4 i sin big ,

(6.2a)

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

130

and

~i( Z ) ~--Ali sinh X i z - Azi sin )tiz

cosh 2t~z + A4i

+ A3i

COS

)kiZ ,

(6.2b)

where the coefficients A a i , o~ = 1, 4 and the axial eigenvalues, hi, are determined through the application of the boundary conditions. For a shell fixed on one end and free on the other the constants in eqs. (6.2) are determined as cosh X~L + cos X~L A3i = - A 4 i = -- sinh X~L + sin X,L

All = - A 2 i = 1,

and XiL are the roots of cosh ~,,L cos X,L = - 1 . As tabulated in table 4, the fluid pressure components induced by the horizontal and rocking seismic excitations contain cos 0 terms. Hence, the functions in eqs. (5.20) and (5.21) are identified as: fil (0) = cos 0,

(6.3a)

fi2(0) = sin 0.

(6.3b)

and

The Fourier expansion in the circumferential direction greatly simplifies the analysis due to the orthogonality properties of the sinusoidal functions. The matrices defined in eqs. (5.3), (5.8), (5.14) and (5.22) through (5.27) can be reduced to axial integrals, after using the integral identities given in Appendix A. (i) The structural mass matrix (refer to eq. (5.3)): s

rL

T

(6.4)

Mij,,,, = 6nm~Jo RhRN:iN:/ dz. (ii) The added mass matrix (refer to eq. (5.8)): oo ad m M~j,,,,8,m~r E F,q( R)YiqYjoA.

(6.5)

q=l

(iii) The material stiffness matrix (refer to eq. (5.14)): (6.6) where R 2 0

C 1=

Eh R ( 1 - v 2)

nvR 0

pR

0 1 -v

2

muR nm

0 1 - PmR

2 0

0

nm 0

m

vR

0 -

1 -v

--nR 2 0 1 -

2

VR2 0

0 n

0 1

(6.7)

131

W.K. Liu, R.A. Uras / Transientfailure analysis of liquid-filledshells I

and

R4

(1

1-v

Eh 3 12R3(I + v)

2

1 -v

1 -v ( 1 - n 2 ) ( 1 - m 2) 1-v

0

0

0

0

(1- n2)vR 2 C2-

-- m 2 ) v R

0

0

0

0

0

0

nm

nR

- nmR

8

8

4

mR

R2

8 -

0

nmR 4

(6.8)

mR 2

8

4

nR 2

nmR 2

4

2

(iv) The geometric stiffness matrix." (a) hydrostatic part (refer to eq. (5.23))

L

, T 1"03X3 03X3]

t

(6.9)

KijG:m=*.,.~rfo nmp~g(H-z)(Bzi ) [03X3 13×3 (B#)dz. (b) vertical ground excitation part (refer to eq. (5.24))

co = Bn,*rfot.(Bz,) , TiN°(z,03×3t)13x3 nmN°o(z,03×3t)/3x3 a]( Bzj)R , K,j,m dz.

(6.10)

(c) horizontal ground excitation and rocking motion: N° and N~, and N°zl include a cos 0 and a sin 0 term, respectively. Therefore, circumferential integrals yield two non-zero stiffness terms due to the orthogonality properties. Using the integral identities presented in Appendix A, the geometrical stiffness terms due to axial and hoop membrane forces, and torsional membrane force become

KGIO __ ~'(~n+l,m+ ~ --ijnra-- -2

"~[L(Rt'~T[NOz'(z' t)13x3

n--l,m]Jo , ' z i J [

03×3

l(nz'j)R dz

nmNO(z ' t)13x3

03x3

(6.11)

and

K,~

=

-~(8~+l,m--8,_l,m) fo Ndzl(z, tW - ' ' T [ 03x3 qr

L

0

,~,Mzi )

[ni3×3

-_ml3x3](B~)R dz,

I}3×3 J

(6.12)

respectively.

(o) The load correction matrix: When the circumferential integrals are carried out in eqs. (5.28), the components of the load correction matrices take the following form:

K,jn,~

0

( Bz* ) dz,

(6.13a)

m Kijn,nfSnmer

fo ,. .,

z, t)(Bzi )

0

m

(Bz~ ) dz,

(6.13b)

W.K. l_a'u,R.A. Uras / Transientfailure analysis of liquid-filledshells I

132

and K,j.,. - ~(8.+1,,. + 8 L P1 _ ' / r n - l , m ) f 0 PI(Z,

r 0

t)(Bz*)[0

0 m

Bz~)dz.

(6.13c)

Finally, for each circumferential mode number n for n - - 1 . . . . . N, the governing finite element equations can be rewritten as

(m)..(d). + (X)..(d). + (X*~-').,._l(d)._l + (X*°)..(d). + (X*l)...+l(d).+l = o (no sum on n),

(6.14)

where the i and j components of each of matrices are defined below: 31 × 3I total mass matrix:

(M).,.={M,~.,.+Mi~a.,.}

forn=m,

(6.15a)

the 31 x 31 stiffness matrix: Ps ( K ) n m - - { gp--,~nm-+gPs--,~n,.+ Kij,r~ }

for n = m ,

(6.15b)

the 31 × 31 time-dependent geometrical stiffness matrix due to vertical ground excitation GO + KO.,. PO } ( K * ° ) . m -- { Kiy.m

for n = m,

(6.15c)

the 31 x 31 time-dependent geometrical stiffness matrices due to horizontal ground excitation and rocking motion

"'IKOllijnm "+'Kijnm}P'

(/[" * ( - 1') nm ~--" I.~"XG10-'ijnm+"

for n = m - 1,

(6.15d)

( K . 1 ) . , ffi t[ x--,j.m olo + o u + Kij,m P1 } - s--q,m

for n = m + 1,

(6.15e)

orn=m+l.

(6.150

and the displacement vector

(d)m={dj,}

forn--m

orn=m-1

After considering eq. (6.14) for each n = 1. . . . . N, the resulting set of equations can be rearranged into the global matrix equation of motion:

Md'+ K d + K * ( t ) d = O ,

(6.16)

where 0

0

0

0

0

0

(M)2 2

0

0

0

0 0 0 0

0 0 0 0

0 0 0

0 (M).. 0 0

0 0

0 0 0 0

0

(M),~N

- (-M ) 1-1 Mffi

'

(6.17)

W.K. Lit~ KA. Uras / Transientfailure analysis of liquid-filled shells I

(g),,

0

0

0

0

0

0

(K)2 2 0 0 0

0 0 0

0 0 (g).. 0

0 0 0

0 0 0 0

0

0

0

0

(K)N N

K=

0 0 0 0

(6.18)

(K*°)n

(K*')I 2

0

0

0

0

(K*(-1))2,

(K*°)22

(K*')23

0

0 0

(g*(-')o.,_,

(g*°)~n

(K*')°.n+,

0

o 0

0 0 0

0

0

0

0

(K*(-'))N,N_~

(K*°)N~

0

K*(t)=

133

o

0

(6.19) and

d-- [ ( d ) l , (d)2 . . . . . (d)n . . . . . ( d ) N ] T.

(6.20)

The tri-diagonal form of K * ( t ) is a direct consequence of the circumferential modal coupling due to n = 1 translational excitation, as encountered in horizontal and rocking seismic excitations.

7. Governing equations for the dynamic stability analysis of liquid-filled shells The matrix equation of motion, eq. (6.16), includes the time-dependent coefficient matrix K * ( t ) which incorporates the effects of the geometric stiffness and load correction matrices. Physical systems modelled through equations of this type are referred to as parametrically excited systems. For stability analysis, eq. (6.16) is cast into a standard form through an orthogonality transformation

(d). = (Q).(u)n

(no sum on n),

(7.1)

where (Q). are the diagonal entries of the orthogonal matrix Q

Q=

(Q)I

o

o

o

o

o

o 0 o 0 o

(Oh 0 o 0 o

0

o 0 (QL 0 o

o 0 o

o 0 o 0 (Q).

o 0 o

o

(7.2)

By the application of this transformation the total mass matrix is normaliT.e,d to the identity matrix and the stiffness matrix is reduced to a diagonal matrix, A, of natural frequencies, ~2, for i -- 1 through I and n -- 1 through N. Consequently, eq. (6.16) becomes

ii + [A + G( t)]u=O,

(7.3)

W.K. Liu, R.A. Uras / Transient failure analysis of fiquid-filled shells I

134

where

-(ah~

0

0

0

(a h2

0

0 0 0 0

0 0 0 0

0 0

A=

(G°)n (G(-1))2 a

a(t) =

0 0 0 (a).. 0 0

0

( G1)12 (a°)22

0

0

0 0 0

0 0 0 0 (A)N N

0

(7.4)

0 (a~123

0

0

0

0

0

0

0

0

(a°)..

(a~).,.+~

o

0

(G(-I))N,N_ 1

(GO)NN

0

o

o

(a(-').,._~

0

0

0

0

0

0

(7.5)

and

u = [(u),, (U)z ..... (u) . . . . . . (U)N] T.

(7.6)

The submatriees in eqs. (7.4) and (7.5) are defined as: eigenvalue submatriees: ( Q ) . T ( K ) . . ( Q ) . -- ( A ) . .

(no sum on n),

(7.7)

geometric stiffness submatrices due to vertical ground excitation: (Q)T~(K*°)..(Q). = (Go)..

(no sum on n),

(7.8)

geometric stiffness submatrices due to horizontal ground excitation: ( Q ) .T( K . 1) . , . + I ( Q ) . + ~ = ( G 1).,.+1

(no sum on n),

(7.9)

(Q)T,(K*(-1)),,,_I(Q),_1 = ( v¢7(-D~J , , , - 1

(no sum on n).

(7.10)

Therefore, A and G(t) are N × N matrices and each submatrix is (31) × (31). As depicted in fig. 3 the block-diagonal terms in G(t) arise from the effect of vertical ground excitation, whereas the non-zero block-off-diagonal terms arise as a result of horizontal ground excitation and rocking motion. The modal coupling generated by the breathing (n = 0) and the translational (n = 1) loads is explained better by choosing only two axial modes for an arbitrary circumferential wave number n (refer to fig. 4). In this case, the transformation matrix consists of two 31 × 1 vectors, qa and q2 ( Q ) . = [(q,), (q2)] .-

(7.11)

The diagonal submatrices ( G ° ) . . correspond to the same circumferential mode number, n, and to axial mode numbers one and two. That is,

(G°)-.

1°°1 [ ao oo a°~ /

for a fixed n,

(7.12a)

135

W.K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

1

2

n

1

_ (G°) l l ~

2

< 1 2 ~1 ~ ) 2 2 ~ ) 2 3 ~

G(O=

~ n

(GI)12"~

~ 0

0

N

0

. •

•

0

0

0

0

0

0.~.

0

(G

0

andg

ground excitation N

~('l~N,N-~(G )l~J vertical

/

horizontal and rocking ground excitation

ground excitation

(G (.I~n,n-h (G)n.n 0 and (G)n.n.l 1 arc 31x31 axial coupling submatriccs for each n=l,..,N

Fig. 3. The form of the time-dependentgeometricstiffnessand load correctionmatrices. where ( G ° ) . n --- (q,)T (K *°).n (qj).

i, j = 1, 2,

(7.12b)

For a given circumferential mode n, the diagonal terms, G°l and GO, are due to the uncoupled axial modes, (1,1) and (2,2), of the fluid-structure system. The off-diagonal terms, G°2 and GO, are due to the coupling between the axial modes one and two as a result of vertical excitation. The off-diagonal blocks, (G(-l))n,,_ 1 and (G1)n.~+l, are brought about by the horizontal and rocking seismic motion and represent the coupling of the (n)th circumferential mode to (n - 1)th and (n + 1)th circumferential modes, respectively. They are given by (G(-'))n.n_ 1 =

<;'1

G2(l__l) G(_I)/ 22

for a fixed n,

(7.13a)

i, j----- 1, 2,

(7.13b)

Jn,n--1

with n=

(q,). ( t o * ( - " ) . . . - l ( q j ) . - i

and

( l)n'n+l

__[<1 <2] Gh

for a fixed n,

(7.14a)

n,n+l

where (G,~).,.+I = ( q i ) T . ( /(.1 ) . . n + l ( q j ) . + l

i, j = 1, 2.

(7.14b)

The diagonal terms, G~ -1), and G2t21), G~ and G~= arise from the above-mentioned circumferential coupling for each axial mode in consideration, (1,1) and (2,2). The cross-coupling of axial and circumferen-

136

W..K. Liu, R.A. Uras / Transient failure analysis of fiquid-filled shells I

[o o0]o02L

]

i_i21

~i22Jn,n+l

r~<-,> ,-,<-,1

~ ~

,G(-lk /'-'n '-'12 / )~n-I =/,..,(-t) ,.-,(a)I

t

modal coupling vertical ground excitation part

axial

GO 1 = (q~)nT (K *0)nn (ql)n 0 T *0

G 12 (ql)n ( K ) m (q2)n GO1 = (q2) Tn(K *0)nn (ql)n =

0 T *0 G22 = (q2) n ( K ) n n

(q2)n

cireurnferential

(1,1)

(n,n)

(1,2)

(n,n)

(2,1)

(n,n)

(2,2)

(n,n)

horizontal and rocking ground excitation parts 1

T K*I.

GH = (ql)n (

)n,n+l(ql)n÷l

(1,1)

(n,n+l)

Gl12 = (q0nT(K*l)n.n+l (q2)n+l

(1,2)

(n,n+l)

GII

(2,1)

(n,n+l)

(2,2)

(n,n+l)

(ql)n-I

(1,1)

(n,n-l)

12 = I"111n*, )n.n-I (qz)n-I (-1) .T.K*(-I). . G21 = (q2) a t )n.n-~ (ql)n-1 f,(-l) . xT,v*(-1),

(1,2)

(n,n-1)

(2,1)

(n,n-1)

(2,2)

(n,n-1)

1

T

*1

T

*1

= (q~n (K)n,n+l (ql)n+l

G 22 = (q:z)n (K)n,n+l (q2)n+l

G11( -=~ql)n I. ) .T...*(-1). 11~ )n.n-I G (-1) (- ~T:K*(-I),

~22 = tqVn ~

)n.n-1(qgnd

Fig. 4. Axial and circumferential modal coupling due to ground excitation.

tial modes axe illustrated by the off-diagonal terms, G~21) and G2~i-1), G]2 and G~I. The twelve axial and circumferential modal coupling due to breathing and translational loads axe also depicted in fig. 4. To study the dynamic stability characteristics of a liquid-filled shell, the stresses in a shell are determined from the zeroth-order equation. The dominant stresses are the membrane stresses in comparison to the bending stresses (Shih and Babcock [15]). Hence, the membrane theory yields a good approximation for the actual stress distribution away from the built-in end. The resulting membrane stresses are given as

N°(z, O, t) = NZ(t)Fz(z)[az + bz cos 0],

(7.15a)

N°(z, O, t) = N°(t)ro(z)[ao + bo cos 0], N~°(z, O, t) ffi NZ°(t)Fzs(z)b~o sin 0,

17.15b) (7.15c)

where N~(t) axe the time-dependent membrane force amplitudes; Fa(z) represent the axial distribution of membrane forces; a~ and ba account for the relative weights of vertical and horizontal ground excitations, respectively; and a = z or # or zO. It should be indicated that under seismic excitation the axial and torsional membrane forces induced in the n = 0 breathing mode axe very small. Nevertheless, the

W.K. Liu, R.A. Uras / Transientfailure analysis of liquid-filledshells I

137

coefficient a~ is included in eqs. (7.15) to account for any additional external axial and torsional forces, respectively. The membrane force amplitudes can be obtained by expanding N*(t) into a Fourier series S

N~(t) = Y'. [N~. cos(s~t) + N~ sin(soJt)],

a = z or 0 or zO,

(7.16)

s--1

where S is the truncation limit for the expansion. N~ and N~, are the Fourier coefficients for the cosine and sine series, respectively. For stability analysis, three non-dimensional parameters e~, e0 and e~e are defined for the axial, circumferential and torsional stresses, respectively:

e~ N~/N~,

Eh 2

z_

=

N:,-

R73(1 -

E# = N~l//Ner,

Ndr -

0.927

~=0--......o,,

N~--0.755 g

,

.~)

(7.17)

(he),

Z

(7.18)

(hE),

(7.19)

where the subscript "cr" refers to the respective static critical buckling forces (Yamald [16]). When eqs. (7.17) and (7.18) are substituted into eq. (6.10), the geometric stiffness matrix arisen from vertical ground excitation becomes S

Kij. moo _ e=a= 8.,.~rKfj.m ~ [N;= cos(so~t) + N~, sin(s60t)] S

+ eoao 8,,.~rK~., E [ NlS: cxis(sa~t) + Nf, sin(so~t)].

(7.20)

Similarly, the geometric stiffness submatrices arisen from horizontal and rocking seismic excitation are cast into the following forms, after utilizing eqs. (6.11) and (6.12) together with eqs. (7.17)-(7.19): s K G,.i.m I 0 ffi

z e=b:7~ ( 8.+ a,m + 8.- i,,. ) K i.#,,. E [Ni~=cos(so~t) + N~ sin(s~0t)] s~1 S

7r

+ eobo-~( 8.+ l,m + 8._ l,m ) K~.m ~'. [ Nl°~cos( sa~t) + N °, sin(so~t)],

(7.21)

and Kcax ¢r ,j.m=~.ob.o~(8.+l,

S

_8._l,m )Kijnm ,o E

[N(°s cos(stot)+ N~fl sin(soJt)],

(7.22)

sml

respectively; where the coefficients of the normalized time-independent matrices are given by

N;I foLr=(~)(B;')T

o3×31

,.

°

N:Or fLHml-,o(z)(B~)T[03×3 K~.m = No1 Jo [03,,3

:

Nff

X'Zj~m N(f

S0L F,o (=)(Bz,)# T [03x3

1ni3×3

0 3 × 3 "1

t"

/3×3] (I~zj)R dz,

-m13x3 ]" ,. 03×3

JtB=jJR d:.

(7.24) (7.25)

138

W . K . Liu, R . A . Uras / T r a n s i e n t f a i l u r e a n a l y s i s o f liquid-filled shells I

Since the load correction terms make no significant contribution to the overall response of the fluid-structure system, their time-dependent parts are neglected. Finally, using the transformation the submatrices of G(t), eqs. (7.8) through (7.10), are rewritten as (no sum on n): S z z z z 0 0 ( G ° ) . . =It Z { [ e:a:g,,,,Ni. + eoaog~nN;~] cos(so~t) + [ e:a:g,,,,N~s + eoaogn,,N~s ] sin(so~t)}, s=l

(7.26) S

(G(-1))n'n-1

---- "It

7 E s=l

1,+eobog~,,,-lN;,-e:ob~og,,,,,-1Ni,]

z {[Ezbzgn,n-1

N z

zO

zO

cos(so~t)

z : - o NoZs-ezoO~ogn,,,-, - :o N 2.] :°1 sin(so~t) } , +[e:bzg.,,,-,N~.+eooog,~,n-1

(7.27)

and S "ff

z

z

0

0

zO

(Ga)""+a = 7 E {[ezb~g,,,.,+lNis+eobog,,,,,+lNls+ezob:og,,,,,+lNis

zO

] cos(so~t)

s=J,

+ [e:b:gn,.+lN~s o o + e:ob:ogn,n+lN~. :o :o ] sin(s~t)}, : : + eobogn,.+lN~s

(7.28)

where

g~,. = (QIT,(K:),m(q),,

re=n-l,

n, n + l ,

(7.29a)

0 g.,.

m=n-1,

n, n + l ,

m=n-1,

n, n + l .

(7.29b) (7.29c)

= (q)T(KO).,. (Q),.

zO - -

g,,,,,- (QI~(K:°I,,,,,(Qlm

Hence, the Fourier expansion of the time-dependent geometric stiffness matrix becomes S

G ( t ) = E [G~ cos(s~0/) + G~ sin(s~0/)],

(7.30)

s=l

where G~ and G~ can be identified through eqs. (7.26) through (7.28). In summary, the governing equation of motion, eq. (7.3), takes the following form: S

ii+ Cil+ A u + E [G~ cos(sto/) + G~ sin(so~t)])u=O,

(7.31)

s=l

where C, a diagonal damping matrix, has been included for the completeness of the stability analysis. The i th and n th components of the C matrix are given by C,, = 2~',,,o,,,

(7.32)

where ~, is the damping coefficient associated with eigenfrequency ~ . This form of equations is commonly referred to as the coupled Hill's equation. The stability characteristics of this equation has been studied by many investigators using the notion of parametric resonance (Hsu [7], Nayfeh and Mook [12], among others). 8. Conclusions

This paper and the companion paper, Liu and Uras [11] are aimed at searching for fundamental understanding of the various possible failure mechanisms of anchored storage tanks under seismic excitation.

W..K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

139

The fluid-structure interaction problem is formulated using the concept of added mass. The contribution of each individual fluid pressure components are identified. A Galerkin/Finite Element discretization is applied to obtain the governing matrix equations. For a better understanding of the transient failure mechanisms, the shell membrane forces are decomposed into a hydrostatic part, and into a dynamic part due to vertical, horizontal and rocking ground motions. The modal coupling in the axial and circumferential directions due to different types of seismic loads is identified. The significance of the axial and circumferential cross-coupling in the failure of tanks is also revealed. For dynamic stability analysis, the matrix equations are cast into a set of coupled Hill's equations by employing an orthogonality transformation. The applications of the present formulation to anchored liquid storage tanks are presented in a companion paper (Liu and Uras [11]).

Appendix A. Circumferential integrals for n > 0

02~ cos(n0) cos(me) dO =

~rO.,.,

02'~ sin(n0) sin(me) dO =

~rO..,,

cos(n0)

d0---0,

f02~r cos O cos(n0) cos(m0) dO =

~[8n+1. m + 8~_1.m],

f02~r cos e sin(nO) sin(m0) d e =

"~[8~+l.m+8._l.m],

f2~ sin 0 cos(n0) sin(n0) dO = ~[8.+1, m - 8._1,m],

f02,, sin O sin(nO) c o s ( m e ) d e = - -~ [8n+1.,. - 8,_,,,,]. References [1] M. Chiba and J. Tani, Dynamic stability of liquid-filled cylindrical shells under vertical excitation, Part II: Theoretical results, Earthquake Engineering and Structural Dynamics 15 (1987) 37-51. [2] D.P. Clough, Experimental evaluation of seismic methods for broad cylindrical tanks, University of California Earthquake Engineering Research Center, Report No. UC/EERG 77-10 (May, 1977). [3] D.P. Clough and A. Niwa, Static tilt tests of a tall cylindrical liquid storage tank, University of California Earthquake Engineering Research Center, Report No. UCB/EERG 79-06 (Feb., 1979). [4] L.H. Dormell, Beams, Plates and Shells (McGraw-Hill, New York, 1976). [5] W. Flugge, Stresses in Shells (Springer Verlag, New York, 1973). [6] G.W. Housner and M.A. Haroun, Vibration tests of full-scale liquid storage tanks, Proceedings of 2nd U.S. International C o n f e r e n c e o n Earthquake Engineering, Stanford University, Aug., 1979. [7] C.S. Hsu, On the parametric excitation of a dynamic system having multiple degrees of freedom, Journal of Applied Mechanics 30 (1963) 367-372. [8] D.D. Kana and R.R. Craig, Jr., Parametric oscillations of a longitudinally excited cylindrical shell containing liquid, Journal of Spacecraft 5 (1968) 13-21.

140

W..K. Liu, R.A. Uras / Transient failure analysis of liquid-filled shells I

[9] W.K. Liu and D. Lain, Nonlinear analysis of hquid-filled tank, Journal of Engineering Mechanics 109 (1983) 1344-1357. [10] W.K. Liu and D. Lam, Numerical analysis of diamond buckles, Finite Elements in Analysis and Design 4 (1989) 291-302. [11] W.K. Liu and R.A. Uras, Transient failure analysis of liquid-filled shells. Part II: Apphcations, Nucl. Engrg. Des. 117 (1989) 141-157, next article in this issue. [12] A.H. Nayfeh and D.T. Mook, Nonlinear Oscillations (John Wiley & Sons, New York, 1979). [13] A. Niwa, Seismic behavior of tall liquid storage tanks, University of California Earthquake Engineering Research Center, Report No. U C / E E R G 78-04 (Feb., 1978). [14] F.G. Rammerstorfer, K. Scharf, F.D. Fischer and R. Seeber, Collapse of earthquake excited tanks, Res. Mechanica 25 (1988) 129-143. [15] C.F. Shih and C.D. Babcock, Buckling of oil storage tanks in SPPL tank farm during the 1979 Imperial Valley earthquake, Journal of Pressure Vessel Technology, ASME 109 (1987) 249-255. [16] N. Yamaki, Elastic Stabihty of Circular Cylindrical Shells (North-Holland, Amsterdam, 1984).