, ,
(大連理工大學 海岸和近海工程國家重點實驗室,遼寧 大連 116024)
以往對液體晃蕩問題的研究均視液艙為剛性液艙[1-3],這與實際情況不相符。實際上液艙都具有一定的彈性。當液艙結構為彈性時,艙內液體與結構將會相互耦合作用,這將會急劇地改變晃蕩沖擊壓力與液體自由液面的非線性運動。尤其當晃蕩沖擊壓力的作用時間與彈性結構的濕周期相近或者處于同一量級時,較之剛性結構,彈性結構的動態響應與液體的壓力可能會急劇增大,這時就要研究其中的流固耦合和水彈性問題了[4-6]。本文采用ADINA-FSI模塊來研究縱蕩激勵下二維矩形彈性液艙內的液體晃蕩問題。
假設流體為不可壓縮的粘性流體,結構為線彈性材料。流體與結構在流固耦合接觸面上法向不直接脫離,只沿切向滑動,即流體與結構在接觸面上任意時刻法向位移都是相等。
1.2.1 結構運動基本方程
線彈性結構小變形時在直角坐標系中的運動方程為
(1)
應變率滿足
(2)

Fis——結構體積力分量;
ρs——結構質量密度;

1.2.2 流體域內N-S方程和連續方程
流體域內ALE描述二維不可壓縮粘性流體的N-S方程和連續方程為
(3)
式中:u,v——x和y方向的速度;
gx,gy——x和y方向的體積加速度;
p——流體壓力;
ρ——流體密度;
υ——流體運動粘滯系數。
流固耦合模型包括三部分:流體域、結構域和流固耦合接觸面,因此邊界條件由流固耦合接觸面和流體自由表面組成。
1.3.1 流固耦合接觸面上的邊界條件
1)運動學條件。在流固耦合接觸面上必須滿足位移協調性條件,用拉格朗日坐標系可以表示為
(4)

下劃線表示變量僅定義在流固耦合接觸面上。
由式(4)可推導出速度協調性條件。
固壁上的不可滑移和可滑移條件為
(5)

n——流固耦合接觸面的法向向量,方向向外。
2)動力學條件。
(6)

1.3.2 流體自由表面條件
自由液面是一個移動的邊界,流體需滿足運動學條件和動力學條件。
1)運動學條件。它要求組成自由表面的流體質點不能流出液面,或者說自由液面的法向速度要和流體的法向速度相同,這一條件對于不同的表面追蹤方法,會有不同的數學表達式。
2)動力學條件。必須滿足表面應力條件,若不計表面張力,則表面應力條件為
(7)
式中:un——垂直于自由表面的法向速度(對流體域來講,以外法向為正);
ut——切向速度;
p0——表面大氣壓強。

(8)
式中:Ff,Fs——流體與結構相應的有限元方程。
ADINA在求解流固耦合問題時有直接耦合求解和迭代耦合求解兩種方法。本文在研究液體晃蕩問題時采用直接耦合法,同時求解流體運動和結構的響應。直接耦合求解時采用Newton-Raphson求解耦合方程。
(9)

其中:λd,λτ——流體方程位移松弛因子和結構方程應力松弛因子。
求解過程概括如下:設初始解X0=Xt,對k=1,2,3…進行迭代,求t+Δt時刻的解Xt+Δt。
1)整理流體、結構和FSI界面方程組,得到耦合的系統方程組。
2)求解耦合系統的線性化系統方程組并更新得到的解。計算應力和位移殘量并與給定的迭代容差作比較。如果得到的解不收斂,且沒有達到FSI迭代的最大次數,則回到第一步繼續下一個迭代;否則,程序停止并顯示不收斂的信息。
3)保存并輸出流體和結構的結果。
考慮一個足夠長的液艙的橫向運動,將其作為一個二維問題處理,忽略沿長度方向的流動及液艙結構縱向變形的影響。矩形液艙縮比模型見圖1,尺寸為長L=1.0 m,高H=10 m,靜止水深為h=0.5 m,液艙壁厚δ=10 mm,結構的材料參數為ρ=7.8×103kg/m3,彈性模量E=210 GPa,泊松比v=0.3;流體采用室溫下的自來水,密度為ρ=1 000 kg/m3,粘度ν=1.005×10-3Pa·m。在ADINA-Structure中建立二維矩形結構模型,假定為平面應變問題,劃分為4節點2-D Solid單元,設定FSI流固耦合邊界條件;同時在ADINA-CFD中建立流體模型,劃分為4節點2-D流體單元,設定與結構對應的FSI流固耦合邊界條件和自由表面邊界條件,然后在ADINA-FSI中進行耦合計算。

a)模型尺寸 b)ADINA網格
圖1矩形液艙模型
縱蕩激勵下液艙內晃蕩液體固有頻率為ωn。當液艙所受激勵的頻率等于艙內晃蕩液體的固有頻率時,液體將會產生大幅晃蕩現象。對于給定幾何尺寸的液艙而言,晃蕩液體的固有頻率由液深決定,一般由Abramson線性理論公式進行估算。
(10)
式中:l——液艙長度;
n——模態數(n=0,1,2,…)。
將本文模型的相關尺寸代入式(10),求得液艙內液體晃蕩的固有頻率fn(fn=2π/ωn)。然后,采用ADINA來計算該模型液體晃蕩的固有頻率fn,將其與線性公式計算結果進行比較,基本一致,比較結果見表1。

表1 不同模態下液體晃蕩固有頻率
2.3.1 Faltinsen解析解
Faltinsen針對水平激勵下二維矩形液艙內的液體晃蕩,提出了一個基于勢流理論的線性解析解,廣泛地應用于各種數值模型的驗證。
二維矩形液艙內靜止水深為h,液艙長度為2a,遭受水平激勵,所遭受的水平周期激勵為
ue=-Acosωt
(11)
式中:ue——液艙激勵速度函數;
A——速度幅值,A=bω;
b——激勵位移幅值;
ω——激勵的圓頻率。
由速度勢函數φ,根據Faltinsen[6]線性解析解,可得到自由面位移η
(12)
(13)

Kn=
原點位于液艙中心、靜止水平面處。
2.3.2 ADINA數值解與解析解的比較
采用2.1模型,根據2.2分析得知:其一階固有頻率為f0=0.846 Hz。考慮兩種工況,分析x=a=0.5 m處自由液面的波高歷時曲線。
1)給定外界激勵頻率f=0.5f0,位移激勵振幅A=0.01 m,求解時間t=10 s。由圖2a)可知,ADINA計算所得數值解與解析解吻合較好。
2)給定共振頻率附近f=0.95f0,為避免大幅共振,將位移激勵振幅調小為A=0.000 4 m。由圖2b)可知,盡管此次位移激勵振幅僅為前一情況的4%,但是由于共振效應,當t=10 s時,其晃蕩波高幅值也達到了前者的幅值量級。
綜上所述,由1)、2)比較可知,采用ADINA模型計算所得數值解與Faltinsen解析解吻合較好,說明本文ADINA模型有效。

圖2 x=a處自由液面波高歷時曲線比較
為了研究艙壁彈性對液艙內液體晃蕩的影響,將2.1二維剛性液艙的壁厚減小為3和5 mm(液艙模型其它參數不變),并且將結果與10 mm進行比較。


表2 不同模態下液體晃蕩固有頻率
在位移激勵振幅A=0.002 m的情況下,給定5個不同的激勵頻率,分析彈性液艙內不同結構剛度、不同激勵頻率下晃蕩液體的自由液面運動情況。頻率分布為:f=0.80f0;f=0.90f0;f=0.95f0;f=1.00f0;f=1.10f0,所對應的自由液面的波高歷時曲線見圖3。

圖3 x=a處自由液面波高歷時曲線
由圖3可以看出,結構壁厚為δ=3 mm時,液體自由液面波高曲線圖中毛刺較多。這是由于液艙艙壁結構為彈性且彈性較大,在晃蕩過程中與液體有強烈的流固耦合作用,由于流體對艙壁的沖擊壓力,彈性艙壁產生了變形,影響了自由液面的波動。隨著艙壁剛度的提高,自由液面波高曲線圖中的毛刺情況有所好轉。這是由于結構剛度提高,艙壁變形較小,對自由液面影響較小。
另外發現,在共振頻率附近(f=0.95f0、f=1.00f0),自由液面波高曲線的毛刺較小,結構的彈性對晃蕩液體的波動影響較小;而其它非共振頻率時自由液面波高曲線的毛刺較多,結構的彈性對晃蕩液體的波動影響較大。


圖3f)為位移激勵振幅為A=0.002 m下3個不同剛度、5個不同頻率下液艙內晃蕩液體的自由液面波高曲線的最大值。可以看出,由于流固耦合作用,彈性液艙內液體的一階固有頻率要小于剛性液艙內液體的一階固有頻率,與表2采用ADINA計算所得的晃蕩液體的固有頻率的比較結果一致。因此,晃蕩液體的共振頻率已不能按照剛性液艙內的頻率來計算,在結構設計時應給予足夠重視。
1)通過與解析解的比較,基于有限元軟件ADINA建立的二維矩形液艙晃蕩模型的數值模擬是有效的。
2)彈性液艙內晃蕩液體與艙壁有強烈的流固耦合作用,由于流體對艙壁的沖擊壓力,彈性艙壁產生了變形,影響了自由液面的波動情況;隨著結構剛度的增大,影響逐漸減小。

4)由于流固耦合作用的影響,彈性液艙內液體的一階固有頻率要小于剛性液艙內液體的一階固有頻率。
[1] LIU Dongming, LIN Pengzhi. A numerical study of three-dimensional liquid sloshing in tanks[J].Journal of Computational Physics,2008(227):3921-3939.
[2] AKYILDIZA H, UNAL E. Experimental investigation of pressure distribution on a rectangular tank due to the liquid sloshing[J].Ocean Engineering,2005(32):1503-1516.
[3] NASAR T, SANNASIRAJ S A. Experimental study of liquid sloshing dynamics in a barge carrying tank[J].Fluid Dynamics Research,2008(40):427-458.
[4] LEE DONG Y. Study on sloshing in cargo tanks including hydro elastic effects[J]. Journal of Marine Scienceand Technology,1999(4):27-34.
[5] CHO J R, PARK S W, KIM H, et al. Hydroelastic analysis of insulation containment of LNG carrier by global-local approach[J].International Journal For Numerical methods in engineering,2008(76):749-774.
[6] ESWARAN M, SAHA U K, MAITY D. Effect of baffles on a partially filled cubic tank: Numerical simulation and experimental validation[J].Computers & Structures, 2009(87):198-205.
[7] WU G X, MA Q W, TAYLOR R E. Numerical simulation of sloshing waves in a 3D tank based on a finite element method[J].Applied Ocean Research,1998(20): 337-355.