球形儲罐考慮儲液晃動時的簡化動力學(xué)模型基本理論

呂 遠(yuǎn)， 孫建剛， 孫宗光， 崔利富， 王 振

(1.大連海事大學(xué) 道路與橋梁工程研究所，遼寧 大連 116026；2.大連民族大學(xué) 土木工程學(xué)院，遼寧 大連 116650)

球形儲罐是石油化工領(lǐng)域非常重要的存儲設(shè)備，其存儲介質(zhì)往往為生產(chǎn)過程中易燃易爆的原料、過程料和成品。此類結(jié)構(gòu)一旦遭遇強震，可能發(fā)生儲罐容器墜落倒塌引起爆炸，進而引發(fā)火災(zāi)等。因此研究球形儲罐在地震作用下的動力響應(yīng)及抗震、減震措施是很有必要的?！稑?gòu)筑物抗震設(shè)計規(guī)范》[1]對球形儲罐的抗震設(shè)計做出了相應(yīng)規(guī)定，其忽略儲液晃動效應(yīng)，采用等效附加質(zhì)量法將球罐簡化為一個單質(zhì)點體系進行抗震設(shè)計，但在地震動作用下球罐的振動為罐內(nèi)儲液與罐壁及支承結(jié)構(gòu)的液固藕聯(lián)振動，其受力過程較為復(fù)雜。事實上，有關(guān)儲罐類的液固耦合振動國內(nèi)外學(xué)者已經(jīng)進行了大量研究，但主要針對立式圓筒形儲罐，而對舉架式球形儲罐的研究相對較少。20世紀(jì)，學(xué)者們對球形儲罐的晃動問題進行了理論和實驗研究[1-7]。Ramaneyulu等[8]采用有限元數(shù)值仿真手段，對LPG球形儲罐進行了地震荷載作用下的可靠性評估。Papaspyrou等[9]采用速度勢理論針對50%儲液球形儲罐，推導(dǎo)出了一種用于計算儲液動液壓力等儲液晃動效應(yīng)的數(shù)學(xué)解析模型，并研究了不同外部激勵條件時的儲液晃動效應(yīng)。Lazaros等[10-11]采用速度勢理論推導(dǎo)出了球形儲罐和水平圓柱形儲罐線性晃動效應(yīng)的數(shù)學(xué)模型，并以50%儲液時為例進行了地震動響應(yīng)研究。肖志剛[12]利用有限元軟件對球形儲罐進行了考慮儲液晃動的地震動響應(yīng)研究，認(rèn)為罐內(nèi)儲液的晃動吸收了部分地震能量，有利于球罐結(jié)構(gòu)安全。戴鴻哲等[13]采用有限元數(shù)值模擬手段，對液化石油球罐進行了地震響應(yīng)研究，并分析了儲液晃動對球罐受力的影響，認(rèn)為球罐罐體在地震作用下是偏于安全的，其支承體系才是薄弱點。Seyyed等[14]采用線性勢流體理論及Mehler-Fock變換求解球罐儲液晃動問題，并研究了不同儲液高度對動液壓力的影響。

目前球形儲罐考慮儲液晃動的地震動響應(yīng)研究已經(jīng)取得了一些成果，但仍處于理論與實驗的研究階段，且理論求解過程十分復(fù)雜，不便于工程應(yīng)用。因此為了更加真實地反應(yīng)球形儲罐在地震動作用時的液固藕聯(lián)振動動力響應(yīng)，本文采用速度勢剛性理論，根據(jù)邊界條件推導(dǎo)出合理的勢函數(shù)，并進一步推導(dǎo)出球形儲罐在地震動作用時的動液壓力、儲液晃動波高、支承底部剪力及傾覆彎矩表達(dá)式，構(gòu)建了便于工程應(yīng)用的球形儲罐考慮儲液晃動簡化動力學(xué)模型，進行了地震動響應(yīng)研究并與有限元數(shù)值仿真模擬計算結(jié)果進行對比，輔證了理論推導(dǎo)的準(zhǔn)確性。

1 基本理論

1.1 基本假定

由于球形儲罐相對于大型立式儲罐來說體積較小，但壁厚較厚，且由于其球形形狀及主要存儲輕質(zhì)油，所以在液體的晃動過程中罐壁發(fā)生的變形極其小，對工程應(yīng)用來說達(dá)到了可忽略不計的程度。故將其假設(shè)為剛性罐壁，設(shè)x0是球罐罐體在地震動作用時的相對位移，其是時間相關(guān)變量，同時我們也知道在由于球罐支承的存在，罐體各部位的位移未必相同，因此x0也是空間相關(guān)變量，記作x0(R,θ,β,t)。球形儲罐在地震動作用中主要激發(fā)第一階振型，故主要以第一階振型作為研究對象，現(xiàn)利用有限元軟件，對球形儲罐進行模態(tài)分析，研究其第一階振型的振動形態(tài)。取某一工程實例，建立有限元數(shù)值仿真模型，并進行模態(tài)分析。

由圖1可知，對球罐罐體來說，第一階振型可近似分解為球罐罐體水平平動及罐體繞赤道面的旋轉(zhuǎn)，由于球罐支承系統(tǒng)的抗彎剛度比抗剪剛度大得多，屬于剪切型結(jié)構(gòu)，故第一階振型以水平平動為主，轉(zhuǎn)動十分微弱，在工程計算時為簡化計算可忽略轉(zhuǎn)動影響。Lazaros等在進行儲液晃動效應(yīng)對工業(yè)球罐抗震設(shè)計影響的研究時便做過類似假定，其忽略了罐體轉(zhuǎn)動影響，假設(shè)球罐罐體整體平動，在地震動作用時罐體各部位位移僅為時間相關(guān)變量，記作x0(t)。

(a) 第一階振型

(b) 罐體水平平動

(c) 罐體轉(zhuǎn)動圖1 球罐第一階振型Fig.1 Spherical tank first order mode

同時假定儲液為無旋、無粘、不可壓縮的理想流體。在地震作用下，根據(jù)速度勢剛性理論[15-16]，儲液地震動響應(yīng)可分為對流運動和剛性運動，儲液總速度勢可分為剛性速度勢和對流晃動速度勢，記為Φ(r,θ,β,t)=φr(r,θ,β,t)+φs(r,θ,β,t)，滿足球坐標(biāo)系下的Laplace方程：

(1)

球罐計算簡圖如圖2所示，地面運動xg(t)作用下，罐體會產(chǎn)生相對位移為x0(t)的振動。因為φr、φs均滿足方程(1)，故分別討論剛性速度勢和晃動速度勢。

圖2 球罐簡圖Fig.2 Spherical tank diagram

1.2 剛性速度勢

根據(jù)上述假定，剛性速度勢φr(r,θ,β,t)應(yīng)滿足Laplace方程(1)和邊界條件：

φr(β+2π)=φr(β)

(2)

(3)

(4)

(5)

1.3 晃動速度勢

同剛性速度勢，晃動速度勢φs(r,θ,β,t)也滿足Laplace方程(1)和邊界條件：

φs(β+2π)=φs(β)

(6)

(7)

(8)

m=0,1,2,3…

(9)

其中hv為表面波動方程，因此可得：

(10)

而又因為x=rsinθcosβ，所以存在

式(10)可以寫為：

(11)

將速度勢方程式(5)、(9)代入式(11)可得：

(12)

由于式(12)在自由液面S2上均滿足，根據(jù)疊加原理，可將兩邊同時乘以調(diào)和函數(shù)φ*(r,θ)，并在圓域S2內(nèi)積分，所以可以得到如下方程：

sinθφ*(r,θ)rdrdβ=

(13)

(14)

同時由邊界條件(8)可得：

(15)

式(15)在儲液與罐壁耦合面S1(r=R)上均滿足，同樣可根據(jù)疊加原理，將兩邊同時乘以調(diào)和函數(shù)φ*(r,θ)，并在圓域S1內(nèi)積分，所以可以得到如下方程：

cosθR2sinθdθdβ=0

(16)

將式(14)、(16)相加可得：

(17)

式(17)可以寫成矩陣的形式：

(18)

其中：

(19)

(20)

(21)

(22)

[f]=[fn(t)]m×1

(23)

(24)

其中

[M]m×m=[P][N1]T

(25)

[K]m×m=[P]([N2]T-[N3]T+[N4]T)

(26)

(27)

(28)

(29)

(30)

參照結(jié)構(gòu)動力學(xué)振型疊加原理[18]，可以將式(24)看作是耦合在一起的m階線性無阻尼運動控制方程，現(xiàn)在需要對其進行解耦，首先求出各階振型及對應(yīng)晃動頻率：

(31)

式中：ωn為n階晃動頻率，{ψn}為對應(yīng)的n階晃動振型；根據(jù)上述內(nèi)容可知，式(24)為m階矩陣方程，其中m=1,2,3,4，為截斷數(shù)，其取值不同時結(jié)果精度會產(chǎn)生差異；

因此式(24)可以轉(zhuǎn)化為m個非耦合的方程：

(32)

其中

(33)

可將式(32)寫為如下形式：

(34)

其中

(35)

式(34)為晃動分量無阻尼運動方程，有阻尼運動方程如下所示：

(36)

由于儲液晃動主要以第一階振型為主，所以只以i=1時為主要研究對象，記xc1(t)=xc(t)。

根據(jù)式(32)～(35)可得儲液晃動速度勢為：

(37)

所以總速度勢為：

(38)

有了總速度勢，可得得出重力場作用下自由液面的振動形態(tài)，及儲液作用在罐壁上的動液壓力：

(39)

(40)

將式(38)分別代入式(39)、(40)中可得：

(41)

(42)

則由儲液運動產(chǎn)生的水平方向基底剪力為：

(43)

可以寫為：

(44)

由儲液慣性力而產(chǎn)生的對支柱底部的傾覆彎矩也是球罐地震動響應(yīng)研究時的主要研究對象：

M1(t)=

(45)

其中：

(46)

(47)

(a) e值

(b) 晃動分量系數(shù)

(c) 剛性分量等效高度

(d) 晃動分量等效高度

(e) 晃動頻率圖3 不同截斷數(shù)m時各參數(shù)曲線Fig.3 The parameters curve of different truncated

由圖3可知，隨著截斷數(shù)m的增大，各參數(shù)曲線趨近于某一值。當(dāng)m≥8時各參數(shù)數(shù)值變化已十分微弱，基本滿足計算精度要求，因此取截斷數(shù)m=9。由于上述各參數(shù)也受儲罐半徑R的影響，因此研究不同儲罐半徑時，各參數(shù)的變化情況，計算結(jié)果如圖4所示。

(a) e值

(b) 晃動分量系數(shù)

(c) 剛性分量等效高度

(d) 晃動分量等效高度

(e) 晃動頻率圖4 不同儲液高度時各參數(shù)曲線Fig.4 The parameters curve of different height of the fluid

(48)

s=0.517 9+0.538 1cos(1.377x)-0.072 68sin(1.377x)-0.059 57cos(2.754x)-0.064 12sin(2.754x)

(49)

(50)

hc=(0.035 01x5-0.149 3x4+0.232x3-0.147x2+0.372x-0.002 425)2R+h

(51)

h0=(-0.013 76x5+0.050 26x4-0.066 86x3+0.030 77x2+0.274 8x-0.000 725 7)2R+h

(52)

表1 曲線擬合優(yōu)度R2Tab.1 Goodness of fit of curves R2

將本文推導(dǎo)的晃動分量系數(shù)與晃動頻率近似解析解與Lazaros[11]值作對比，如圖5所示，發(fā)現(xiàn)兩者十分接近，由此也輔證了本文推導(dǎo)過程的正確性。

(a) 晃動分量系數(shù)對比

(b) 晃動頻率對比圖5 數(shù)值對比Fig.5 Numerical comparison

2 簡化力學(xué)模型

上述推導(dǎo)中所求出的Q1,M1分別為儲液在地震動作用下而產(chǎn)生的作用于支承底部的剪力和彎矩。地震作用時，支柱底部總剪力和總彎矩還應(yīng)包含由球罐罐體、支承質(zhì)量、球罐配件等產(chǎn)生的剪力和彎矩?！稑?gòu)筑物抗震設(shè)計規(guī)范》(GB 50191—2012)[1]中將球罐簡化為一個單質(zhì)量體系:

meq=m1+m2+0.5m3+m4+m5

(53)

式中:m1為球殼質(zhì)量；m2為液體有效質(zhì)量；m3為支柱和拉桿質(zhì)量；m4為保溫層質(zhì)量；m5為球罐其他附件質(zhì)量，包括各開口，噴淋裝置，梯子和平臺等。

由于m4、m5一般附加于罐體，1/2支柱和拉桿質(zhì)量集中于支承頂部附加在罐體上，為便于工程設(shè)計，此時也可將m1、0.5m3、m4、m5整合為一個集中質(zhì)量，記作

ms=m1+0.5m3+m4+m5

(54)

可設(shè)地震動作用時集中質(zhì)量ms對應(yīng)的相對位移為罐體相對位移x0(t)。

進而可得由于其慣性作用而產(chǎn)生的基底剪力為：

(55)

總基底剪力為：

(56)

由罐體，支承體系及附件等產(chǎn)生的基底彎矩表達(dá)式為：

(57)

式中：h1為球殼集中質(zhì)量等效高度；h3為支承質(zhì)量等效高度；h4為保溫層質(zhì)量等效高度；h5…h(huán)n為各附件等效高度，根據(jù)不同球罐結(jié)構(gòu)形式，表示多種附件對應(yīng)的等效高度。

對球殼來說其對支柱底部產(chǎn)生的彎矩可表達(dá)為

m1h1=

(58)

求解積分，可得：

(59)

式中：R1為球罐內(nèi)壁半徑；R2為球罐外壁半徑；ρ1為罐體密度；R1-R2為球罐罐壁厚度，相對于R1、R2來說十分小，因此式(59)可以寫為：

(60)

(61)

總基底彎矩為：

(62)

則根據(jù)式(56)、(62)可以構(gòu)建球形儲罐考慮儲液晃動時的簡化動力學(xué)模型

圖6中kc,cc,k0,c0,ms分別為：

kc=mcω2,cc=2ξωmc

(63)

式中:k0,c0分別為支承結(jié)構(gòu)的剛度和阻尼，可根據(jù)規(guī)范[1]相關(guān)公式可得出；ms為儲罐和支承結(jié)構(gòu)質(zhì)量，其等效高度為h+R。由Hamilton原理，可推導(dǎo)相應(yīng)的運動方程：

圖6 球罐簡化力學(xué)模型Fig.6 Spherical tank simplified mechanical model δ(T-V)dt+δWdt=0

(64)

式中：T,V分別為系統(tǒng)的動能和勢能；W為非保守力做的功。

由式(64)可得運動方程為：

(65)

3 算例分析

選取某一1 000 m3液化石油氣罐作為算例，儲液高度為H=1.5R，忽略其內(nèi)壓影響，儲液密度為480 kg/m3，球罐直徑為12.3 m，球心距地面8 m，拉桿上部連接處距地面6 m，球殼與支承體系均采用雙線性強化模型，算例參數(shù)如表2所示。分別選擇四類場地中滿足規(guī)范[1]的五條天然波和兩條人工波對球形儲罐簡化力學(xué)模型和有限元數(shù)值仿真模型進行地震動響應(yīng)對比研究，四種場地波加速度反應(yīng)譜如圖7所示，調(diào)整加速度時程曲線峰值為0.2 g，計算結(jié)果,如表3～6所示。

表2 模型參數(shù)Tab.2 Model parameters

從表3～表6中數(shù)據(jù)可以看出，對不同場地地震動輸入，球罐地震動響應(yīng)差異較大，從(類場地到Ⅳ類場地，晃動波高、基底剪力及傾覆彎矩均值逐漸增大。變異系數(shù)反應(yīng)了各工況計算數(shù)值在同場地不同地震波條件下的離散程度，從表中數(shù)據(jù)可知，變異系數(shù)均較小，說明離散程度較小，選波較為合理。而均值差異率則反應(yīng)了本文推導(dǎo)的理論解與有限元解的差異大小，從表中可以看出，均值差異率均在10%以下，有限元解與理論解十分接近，且理論解均值均大于有限元解均值，計算結(jié)果偏于安全，由此則進一步佐證了本文所構(gòu)建的簡化力學(xué)模型的準(zhǔn)確性及可靠性。

(a) Ⅰ類場地

(b) Ⅱ類場地

(c) Ⅲ類場地

(d) Ⅳ類場地圖7 四類場地加速度反應(yīng)譜Fig.7 Four types of site acceleration response spectrum

表3 Ⅰ類場地地震動響應(yīng)對比Tab.3 ClassⅠsite ground motion response comparison

注：差異率=(理論解-有限元解)/理論解

表4 Ⅱ類場地地震動響應(yīng)對比Tab.4 ClassⅡsite ground motion response comparison

表5 Ⅲ類場地地震動響應(yīng)對比Tab.5 Class Ⅲ site ground motion response comparison

表6 Ⅳ類場地地震動響應(yīng)對比Tab.6 Class Ⅳ site ground motion response comparison

4 結(jié) 論

(1) 采用速度勢剛性理論，根據(jù)邊界條件推導(dǎo)出合理的勢函數(shù)，并進一步推導(dǎo)出球形儲罐在地震動作用時的動液壓力、儲液晃動波高、支承底部剪力及傾覆彎矩表達(dá)式，并分析了不同截斷數(shù)m及求球罐半徑R對各參量的影響，將本文推導(dǎo)的晃動分量系數(shù)和晃動頻率近似解析解與與Lazaros值作對比，兩者十分接近，由此也輔證了本文推導(dǎo)過程的正確性。

(2) 構(gòu)建了便于工程應(yīng)用的球形儲罐考慮儲液晃動簡化動力學(xué)模型，進行了地震動響應(yīng)研究并與有限元數(shù)值仿真模擬計算結(jié)果進行對比，有限元解與理論解十分接近，且理論解均值均大于有限元解均值，計算結(jié)果偏于安全，由此則進一步佐證了本文所構(gòu)建的簡化力學(xué)模型的準(zhǔn)確性及可靠性。