地震激勵下大型儲罐內液體晃動響應分析

梁俊杰，王瓊瑤，林國勉，徐 凱

（五邑大學智能制造學部，廣東 江門 529000）

0 引言

液體晃動是在有限空間內受到外部激勵后具有自由液面的流體的運動[1]。由于其強大的非線性和隨機性，它一直是許多研究的焦點，特別是在航空航天、陸海運輸工程和液體儲存領域。以儲油罐為例，儲油罐通常被視為剛性的，一般用于儲存石油或其他化學液體產品，這類儲液通常易燃或有毒。儲罐中的儲存液體晃動會給儲罐壁帶來額外的力，并沖擊儲罐頂部，這可能會導致儲罐結構破壞和無法正常使用。嚴重情況下，儲存的液體會泄漏[2]，這可能會污染周圍環境，甚至引發火災。為了減少這種情況，一種方法是降低儲罐的液體填充率，并為儲罐中的存儲液體保留足夠的振動空間，以有效地減少儲罐頂部受到的沖擊載荷，但這將導致儲罐的利用率不足，對于較大的儲罐經濟性較差。相比之下，在儲罐中安裝擋板[3]、浮板[4]和彈性氣囊[5]等阻尼裝置可以減少流體沖擊載荷。其中，擋板由于其安裝簡單、抑制液體晃動的效率高而被公認為最常用的防晃裝置。因此，研究擋板參數對液體晃動固有特性的影響是一個非常有意義的課題。

1 液體晃動研究方法的發展

對于充液容器，由其滿足不同理論下建立的控制方程，通常可將液體晃動的研究方法分為兩類：一類是基于勢流理論的分析方法。Abramson[6]對各種規則幾何形狀（包括矩形、環形、圓柱體等）儲罐中的液體晃動進行了系統研究，應用線性勢流理論分析了橫向微晃動下儲罐中儲存的液體。Hasheminejad[7]對二維半充液剛性橢圓儲罐進行了晃動分析，并利用連續保角技術和勢流理論推導了一套計算帶有垂直擋板的橢圓儲罐固有振動特性的公式。

一般而言，分析方法適用于結構簡單且晃動幅度較小的情況，而數值方法適用于不規則儲罐形狀且晃動幅度較大的情況。此外，數值方法也可用于研究具有各種阻尼裝置（如擋板、浮板等）的儲罐中的液體晃動。根據控制方程的不同離散化方法，數值方法可分為有限元法（FEM）、有限差分法（FDM）、有限體積法（FVM），邊界元法（BEM），光滑粒子流體力學（SPH）。

Biswal 和Bhattacharyya[8]使用FEM 評估了液體和儲罐-擋板系統的耦合響應，研究了不同參數（如復合擋板、層壓方案）對液體填充復合儲罐中晃動頻率和耦合振動頻率的影響。ünal 等[9]應用FDM 研究了封閉的二維部分填充矩形儲罐中的液體晃動，并比較了T 形擋板和無擋板儲罐中的晃動情況，結果表明，當擋板的高度大于液位的80%時，擋板在壓力和波浪阻尼方面完全有效。Chen Yang[10]采用FVM 方法和Lattice Boltzmann 方法相結合的方法對不同容器設計下的液體晃動行為進行了快速評估，結果表明，LBM 模型的計算時間開銷僅為FVM 模型時間開銷的7.0%左右，且不影響數值精度。Cho 等[11]應用BEM 獲得了多孔水平擋板晃蕩的解析解，考慮了矩形儲罐中心和兩側壁的兩個擋板位置，以確定不同的孔隙率、長度和浸沒深度。結果表明，與安裝在儲罐中心的水平多孔擋板相比，安裝在兩個儲罐壁上的水平多孔隔板顯著抑制了劇烈的共振晃動響應。McNamara[12]提出了二維顯式不可壓縮SPH 模型，用于研究安裝在調諧液體阻尼器內的矩形儲罐中的液體晃動響應。發現SPH 模型與TLD 儲罐尺寸范圍內的實驗數據非常吻合，可進一步用于研究矩形TLD的非線性響應。

除了上述用于研究液體晃動現象的分析和數值方法外，實驗方法是一種直接的方法。Chia Ren Chu 等[13]使用LES 模型研究了矩形水箱中的晃動，并通過添加多個底部擋板測試了其抑制晃動的能力。實驗結果表明，多個擋板對水動力的抑制效果優于單個擋板。多個擋板可以改變容器的固有頻率。防晃效果將隨著擋板數量和高度的增加而增加。當擋板的相對高度大于0.75 時，擋板對晃蕩的影響將減小。T.Nasar和Sannasiraj[14]分別在罐壁的l/3 和21/3 處安裝了孔隙率分別為15%、20.2% 和25.2%的擋板，測試了它們在規則束激勵下的晃動動力響應，并討論了這些擋板方案的有效性和駁船效應。

本文建立了一個基于FVM 的數值模型來研究直立圓柱形儲罐中的液體晃動現象。由于二維模型獲得的模擬結果的精度與三維模型相當[15]，本文采用了二維模型，該模型模擬了直立圓柱形儲罐的最大橫截面。根據對自由液面、動態晃蕩力和施加在罐壁上關鍵點的沖擊應力的評估，提出了幾種擋板配置，以研究其對地震激勵下液體晃蕩的防晃蕩效果。

2 理論背景

2.1 控制方程

基于OpenFoam 進行液體晃動模擬。OpenFoam 是基于FVM 的開源軟件。由于封閉罐中晃動液體的速度較低，可采用層流模型。忽略液體的壓縮，封閉罐中的控制方程為：

式中：?為哈密頓量；u→為速度；ρ為液體密度；p為壓力；為部分黏性應力張量，（式中：μ為動態黏度，S為平均應變率張量，S=，且I=δi）j；fb為身體力。數值模擬采用了氣相和液相兩相流模型。

基于流體體積（VOF）方法捕獲兩相之間的界面，該方法需要定義指數函數，并需要了解界面是單獨被某一流體占據還是同時被兩種流體占據，這需要求解相分數方程。

式中：α為每個單元的液相分數。α=0表示網格完全被氣體占據；α=1 表示網格完全被液體占據；而0 ＜α＜1 表示細胞同時被氣相和液相占據，即網格位于氣相和液相的界面。

每個網格區間的密度、速度等物理性質由網格的液體體積分數計算，其方程為：

式中：αc為每個細胞液體的體積分數；ρl為液體的密度；ρa為氣體的密度；v為運動黏度；vl為液體的運動黏度；va為空氣的運動黏度。

2.2 OpenFoam 及其求解器interFoam

使用OpenFoam 中的兩相流求解器“interFoam”與動態網格技術一起解決儲罐中液體的晃動情況。模擬實驗首先會將罐壁設置為移動邊界，將地震激勵下的移動理想化。對于臨界自由曲面，求解器基于VOF 理論方法中的網格重建和方程執行變形捕捉。壓力-速度耦合問題通過PIMPLE 算法解決，該算法結合了PISO（壓力隱式算子分裂）和SIMPLE（壓力關聯方程的半隱式方法）。

2.3 模型創建

為研究地震運動作用下儲罐中液體的晃動以及各種阻尼擋板配置的防晃效果，本文采用一個實際尺寸的工業容器進行了模擬實驗。尺寸如圖1、表1所示。

圖1 立式儲罐的最大橫截面圖

為了分析儲罐在真實地震運動下的晃動，本次模擬選擇了kocaeli 波[16]。如圖2 所示，該地震數據選自Yarimca（YPT060）平臺，該長周期地面運動根據地面加速度進行縮放。最大加速度為1.32 m/s2。為了確保模擬符合實際情況，選擇晃動液體作為原油，密度為857 kg/m3，運動黏度為0.000 032 m2/s。

圖2 地震加速度變化情況

如圖1 所示，在左側罐壁和罐底設置了4 個監測點。監測點P1至P3用于監測地震運動下罐壁壓力隨時間的變化，P1至P3相對于罐底的高度分別為9.0、6.0 和1.0 m。P4是波高監測點，用于分析左側罐壁附近自由面高程的時程。為了研究環形擋板配置對液體晃蕩響應的防晃蕩效應，在模擬中考慮了具有不同相對寬度（擋板寬度與儲罐半徑之比，w/R）和相對安裝高度（擋板高度與靜態時液體高度之比，hb/hl）的環形擋板。選擇了4 種不同的相對寬度，即0.21、0.32、0.49和0.72，相對安裝高度選擇為0.39、0.50、0.61、0.72、0.83 和0.94。在不同的擋板配置下，比較了晃動響應，如監測點處的壓力、施加在罐壁上的動態晃動力、自由表面位移的評估，以研究其防晃動效果，并最終找到最合適的擋板配置。

3 采用不同擋板方案的儲罐中油液晃動特性的數值分析

3.1 自由液面變化

通過保持一個參數不變ua 改變另一個參數的方式研究了擋板的相對寬度和相對安裝高度對自由液面變化的影響。將帶不同擋板的儲罐與不帶擋板的儲罐的仿真結果進行了比較。

圖3 顯示了具有不同相對擋板寬度和相對安裝高度的折流板儲罐和空罐的最大波高。當擋板相對寬度超過0.32 時，隨著擋板相對安裝高度的增加，最大波高顯著降低，而當擋板相對高度小于0.32 時（擋板相對安裝深度從0.39增加到0.83），最大波高也隨之降低。然而，隨著擋板寬度的進一步增加，最大波高顯著增加，這是因為擋板寬度太短，無法覆蓋液體晃動的波峰，從而導致波峰和擋板邊緣發生碰撞。

圖3 波高峰值變化

為了更具體和直觀地表達每種擋板配置的晃動抑制效果，引入了波高抑制率的概念，可用式（6）表達。

式中：η0為空罐的最大波高；η為安裝折流板后的最大波高；δ為波高抑制率。

圖4 為各種擋板配置下的波高抑制率。對于相對寬度為0.72的擋板，當擋板相對安裝高度在0.39～0.94范圍內時，波高抑制率為5%～10%，而對于相對寬度0.21、0.32 和0.49 的擋板，相應的波高抑制率分別為0.5%～6.4%、2%～6.5%和3.5%～7.5%。從圖中可以看出，除了擋板寬度小于0.32的情況外，隨著擋板安裝高度的增加，波高抑制率增加。這意味著，當立式圓柱形儲罐中擋板的相對寬度大于0.32，擋板越靠近液體自由表面，擋板的防晃蕩效果更加明顯。此外，從圖中可以看出，當擋板相對安裝高度從0.83 增加到0.94 時，相應的波浪抑制率呈現緩慢增長。因此，折流板直立圓柱形儲罐的最佳相對安裝高度在0.72～0.83之間。

圖4 波高抑制率

3.2 壁面上的壓力變化

圖5 為各種擋板配置下施加在左壁上壓力的最大值。可以看出，對于不同相對寬度的擋板，隨著相對安裝高度的增加，施加在左側罐壁上的力的最大值增加，尤其是當擋板接近液體自由表面時。這是由于擋板的添加確實限制了晃動幅度，但同時限制了動能轉化為重力勢能，從而導致動態晃動力的增加。這意味著折流板具有顯著的阻尼效應。此外，對于具有給定相對寬度的擋板，相對安裝高度越大，壓力減少所需要的時間越短，這意味著擋板的阻尼效應對于安裝高度較高的擋板更明顯。其原因是液體和擋板之間存在相互作用。相對地，擋板的阻尼效應能使壓力更快地達到穩定狀態。

圖5 壓力峰值變化

3.3 關鍵點的沖擊應力變化

監測點P1 設置在相對高度為0.11（相對于液體高度hl）的較低位置，因為擋板的最低相對安裝位置為0.39，這使其始終位于擋板下方。圖6 顯示了各種擋板配置下P1 的峰值壓力。從圖中可以看出，除非其接近自由表面，擋板相對安裝高度的變化只對壓力峰值有輕微影響。相對地，擋板相對寬度對壓力峰值有很大影響。隨著擋板寬度的增加，壓力峰值急劇增加。如圖7 所示，地震激勵在30 s時停止，儲罐中的液體晃動可視為自由振動，這可以證明各種擋板的阻尼效果。圖7 顯示了地震激勵停止后P1處的壓力峰值，當擋板相對寬度處于較低水平時，P1 處的峰值壓力隨著擋板安裝高度的增加呈線性下降。當擋板寬度足夠大時，擋板安裝高度對P1處壓力峰值的影響是有限的。

圖6 P1的峰值應力

圖7 地震激勵結束后P1的峰值應力

與P1 不同，監測點P3 被設置在相對高度為1.0（相對于液體高度hl）的位置，這使得它總是在擋板上方。圖8 為不同擋板配置下P3 處的峰值壓力值，結果表明擋板相對寬度和安裝高度對P3 處的壓力峰值有顯著影響。當擋板安裝高度從0.83 增加到0.94 時，P3 處的峰值壓力值反而在一定程度上增加，這是因為存在擋板，與圖3所示的峰值自由表面評估的情況一致。

圖8 P3的峰值應力

圖9 為地震激勵停止后P3 處的壓力峰值，因為擋板始終低于監測點P3，無論擋板高度如何，擋板寬度對P3處的峰值壓力變化有輕微影響。相對地，擋板的相對安裝高度對該點的峰值壓力有顯著影響，隨著擋板安裝高度的增加，P3處的峰值顯著降低。

圖9 地震激勵結束后P3的峰值應力

圖10 顯示了各種擋板配置下P2 處的峰值壓力值。監控點P2 不同于點P1 和P3，它被設置在相對高度為0.67（相對于液體高度hl）的位置。根據擋板的安裝位置，監控點并不總是在擋板上方。從圖中可以看出，隨著安裝高度的增加，P2 處的峰值壓力幾乎以線性方式降低，除了擋板從監控點下方的位置轉移到監控點上方的位置以外，即相對安裝高度從0.61 增加到0.72 時，P2 處的峰值壓力顯著地提高了。圖11 顯示了地震激勵停止后擋板配置對P2 峰值壓力值的影響，與P1 和P3 的情況類似，峰值壓力值隨著擋板安裝高度的增加而減小，然而，當擋板相對安裝高度達到0.72 時，隨著擋板安裝高度的進一步增加，峰值壓力值幾乎保持恒定。

圖10 P2的峰值應力

圖11 地震激勵結束后P3的峰值應力

4 結束語

本文模擬分析了Kocaeli 地震波激勵下立式儲罐在最大垂直截面處的晃動行為。使用基于FVM 數值方法的OpenFoam 平臺建立2D 模型，采用VOF 方法和動態網格技術捕捉了液-氣兩相界面。實驗中儲罐模型和擋板都被假設為剛性的。分析表明，擋板的相對寬度和安裝高度對液體晃動有顯著影響。

（1）當擋板接近液體自由表面時，只要擋板的相對寬度大于0.32，擋板的防晃效果就非常明顯。

（2）擋板的添加使得液體在晃動過程的初始階段更頻繁地振蕩。之后，擋板的阻尼效應使液體更快地達到穩定狀態。

（3）隨著相對安裝高度的增加，施加在左罐壁上的力的最大值增加，特別是當擋板接近液體自由表面時。這是因為添加擋板確實限制了晃動幅度。

（4）隨著安裝高度的增加，監測點的峰值壓力幾乎呈線性下降。并且位于擋板上方的監測點的峰值壓力明顯大于擋板下方的監測點。

（5）在地震激勵下，折流板立式圓筒形儲罐的最佳相對安裝高度為0.72～0.83。