基于VOF法的液艙晃蕩數(shù)值模擬及載荷計(jì)算

劉楨兵

(斯特拉斯克萊德大學(xué)，格拉斯哥G40LZ)

舶液艙晃蕩作為一種常見(jiàn)的流體運(yùn)動(dòng)現(xiàn)象，通常發(fā)生在部分充滿液體的液艙中［1］。當(dāng)液艙與船體運(yùn)動(dòng)發(fā)生耦合時(shí)，將使船體運(yùn)動(dòng)加劇甚至傾覆，同時(shí)液體晃蕩是非常復(fù)雜的流體運(yùn)動(dòng)現(xiàn)象，呈現(xiàn)出強(qiáng)的非線性和隨機(jī)性。

1 VOF法

流體體積法最早出現(xiàn)在DeBar(1974)的文章中，Hirt＆Nichols［2］形成了它完整的理論體系和實(shí)現(xiàn)方法，現(xiàn)在的VOF方法都是在Hirt＆Nichols的基礎(chǔ)上發(fā)展起來(lái)的。VOF法的基本原理是通過(guò)研究網(wǎng)格單元中流體和網(wǎng)格體積比函數(shù)來(lái)確定自由面，追蹤流體的變化，而非追蹤自由液面上質(zhì)點(diǎn)的運(yùn)動(dòng)。VOF法使得流體的質(zhì)量得到守恒，并可跟蹤自由面上發(fā)生翻轉(zhuǎn)、吞并、飛濺等復(fù)雜的自由面現(xiàn)象。VOF法簡(jiǎn)單而有效，得到了廣泛使用，在實(shí)際工程中得到了令人滿意的計(jì)算結(jié)果。

2 建模及網(wǎng)格類型

本文所有的建模和網(wǎng)格劃分都是在Gambit軟件上完成的。Gambit可以生成結(jié)構(gòu)化網(wǎng)格和非結(jié)構(gòu)化網(wǎng)格，雖然非結(jié)構(gòu)化網(wǎng)格對(duì)復(fù)雜幾何形狀的物理模型更有吸引力，但由于其在處理邊界層還存在著不足，而且在同等精度下非結(jié)構(gòu)網(wǎng)格的網(wǎng)格數(shù)是結(jié)構(gòu)化網(wǎng)格的好幾倍。由于使用全非結(jié)構(gòu)化網(wǎng)格來(lái)處理時(shí)對(duì)計(jì)算硬件計(jì)算能力要求較高，本文的計(jì)算模型采用結(jié)構(gòu)化網(wǎng)格。另外網(wǎng)格數(shù)量與時(shí)間步長(zhǎng)的選擇對(duì)數(shù)值計(jì)算是有影響的。網(wǎng)格數(shù)量太多時(shí)要求時(shí)間步長(zhǎng)很小，這樣會(huì)大大增加計(jì)算花費(fèi)的時(shí)間;網(wǎng)格數(shù)量太少時(shí)計(jì)算花費(fèi)的時(shí)間縮短許多，但是求解所得的結(jié)果不夠精確。因此網(wǎng)格數(shù)量應(yīng)當(dāng)適中，一方面可以節(jié)省計(jì)算花費(fèi)的時(shí)間，讓計(jì)算可以有效率地運(yùn)行;另一方面所得的計(jì)算結(jié)果能夠達(dá)到理論研究所需要的精度。本文對(duì)離散的RANS方程都是采用分離求解器求解(segregated solve)，壓力速度耦合以PISO法計(jì)算;本文所有的計(jì)算模型都帶有自由液面，壓力項(xiàng)離散使用體積力加權(quán)方法(body force weighted)，用gauss-seidel方法迭代求解，以代數(shù)多重網(wǎng)格技術(shù)加速迭代收斂速度。

3 程序流程及數(shù)值計(jì)算模型

程序流程及數(shù)值計(jì)算模型見(jiàn)圖1、2。

4 計(jì)算實(shí)例

4.1 算例1:低液深大幅晃蕩下的非線性和隨機(jī)性

在滿載狀況下的晃蕩固然要引起我們的重視，但是在空載出港或到港情況下，往往液艙內(nèi)液體深度很低，并且考慮到壓載艙內(nèi)不同水位的壓載水的晃蕩影響，低液深的大幅晃蕩對(duì)于船舶的穩(wěn)性和操縱性會(huì)帶來(lái)很大程度的影響，這種現(xiàn)象不容忽視。

4.1.1 計(jì)算模型及問(wèn)題描述

液艙模型如圖3所示。液艙做強(qiáng)迫橫蕩運(yùn)動(dòng)，橫向激勵(lì)為x=Asin(ωt)，其中橫蕩周期為1.3 s，頻率為0.77 Hz，橫蕩幅值為0.06 m，P1為壓力監(jiān)測(cè)點(diǎn)。表1給出了幾種計(jì)算工況。

圖3 幾何模型(單位:m)

表1 幾種工況

4.1.2 計(jì)算結(jié)果及分析

1)圖4為自由液面在不同時(shí)刻的運(yùn)動(dòng)情況，液面變化可以說(shuō)明以下幾點(diǎn):

①自由液面發(fā)生沖頂?shù)乃俣葮O為迅速，當(dāng)晃蕩比較充分時(shí)波的破碎、翻卷現(xiàn)象比較明顯，液體晃動(dòng)比較激烈。

②晃動(dòng)比較充分時(shí)，可以看出雖然液艙做周期性運(yùn)動(dòng)，但自由表面波的變化卻并不是完全符合周期性的變化，表現(xiàn)了明顯的非線性特征。

③大幅晃蕩時(shí)自由表面的運(yùn)動(dòng)存在太多的隨機(jī)性，因而壓力的測(cè)量可能存在很多的不確定性，影響因素較多。

④晃蕩載荷還具有離散性和飽和性［3］。離散性說(shuō)明拍擊力不呈現(xiàn)任何規(guī)律性，需要用概率分布來(lái)描述:飽和性則表示拍擊力的大小并不與外界激勵(lì)成正比關(guān)系，存在飽和現(xiàn)象，在非規(guī)則的激勵(lì)下表現(xiàn)得更為明顯。

圖4 自由表面的變化

2)圖5～6為不同工況的數(shù)值計(jì)算結(jié)果。

可以看出晃蕩壓力曲線表現(xiàn)了明顯的非線性特征，此外各工況之間的壓力曲線在不同程度上存在一定的差異，分析可以得出以下結(jié)論:

①由于網(wǎng)格的疏密以及時(shí)間步長(zhǎng)的不同，不同工況的計(jì)算結(jié)果均會(huì)出現(xiàn)一些壓力突變值，分析原因可能是時(shí)間步長(zhǎng)相對(duì)網(wǎng)格間距太大或太小而產(chǎn)生的截?cái)嗾`差和舍入誤差。

②由于晃蕩發(fā)生沖頂時(shí)，自由表面發(fā)生的波的破碎、翻卷等現(xiàn)象是短暫而迅速的，進(jìn)而對(duì)壁面所造成的沖擊也是快速的。因此，不同的時(shí)間步長(zhǎng)下所捕捉的壁面壓力的大小可能會(huì)有所不同。

③由于在監(jiān)測(cè)壁面壓力時(shí)采用的是“貼片法”，即在壁面監(jiān)測(cè)點(diǎn)選取一單元作為代替該點(diǎn)進(jìn)行監(jiān)測(cè)，所得出的壓力是該單元面的平均壓力值。因此在發(fā)生壁面沖擊時(shí)，單元?jiǎng)澐值拇笮?duì)壓力的捕捉也有一定影響。

④當(dāng)晃蕩為微幅晃蕩時(shí)，數(shù)值計(jì)算的結(jié)果受網(wǎng)格的疏密、時(shí)間步長(zhǎng)的變化等影響不大，而當(dāng)晃蕩比較激烈(大幅晃蕩，運(yùn)動(dòng)周期靠近共振周期時(shí)，液體晃蕩更加劇烈［4］)，晃動(dòng)液體的沖頂現(xiàn)象比較嚴(yán)重，進(jìn)而自由表面的波的破碎、翻卷等現(xiàn)象也尤為明顯。

4.2 算例2:高液深微幅晃蕩下同一主尺度下的矩形液艙與菱形液艙晃蕩對(duì)比

縱觀各類液貨船，出于不同用途以工程實(shí)際應(yīng)用，液艙形式從單一的矩形衍變成菱形、圓形等，尤其是隨著LNG船舶的日益盛行，菱形艙室逐漸成為晃蕩研究的一個(gè)重點(diǎn)。

4.2.1 計(jì)算模型及問(wèn)題描述

根據(jù)文獻(xiàn)［5］給定液艙模型，如圖7，其中矩形和菱形的主尺度是相同的，所不同的是菱形在頂部和底部的折角有所改變。液艙受強(qiáng)迫橫蕩運(yùn)動(dòng)，橫向激勵(lì)為x=Asin(ωt)，其中橫蕩周期為1.404 s，頻率為0.71 Hz，該液深時(shí)的諧振頻率為0.69 Hz，橫蕩幅值為0.015 m。表2給出艙室結(jié)構(gòu)劃分情況，其中菱形1的折角處尺寸為0.1 m×0.1 m，菱形2的折角處尺寸為0.1 m×0.2 m。

圖8 矩形液艙與菱形液艙模型(單位:m)

表2 艙室結(jié)構(gòu)劃分情況

在以下的壓力圖中分別以rhombus和Rhombus命名菱形1和菱形2，以示區(qū)別。

4.2.2 計(jì)算結(jié)果及分析

1)圖9為10 s與20 s時(shí)刻自由表面的位置。在圖中可以看出，在相同的外界激勵(lì)下，菱形液艙在晃蕩過(guò)程中自由液面擺動(dòng)幅度比較大，矩形液艙則相對(duì)緩和一點(diǎn)。

2)圖10則給出了壁面點(diǎn)P1、P2的壓力的時(shí)間歷程。

圖9 不同時(shí)刻自由表面的形狀

圖10 3個(gè)模型在P1和P2點(diǎn)壁面壓力的時(shí)間歷程

3)圖11給出了2個(gè)菱形液艙在壁面點(diǎn)P1和P2處的壓力對(duì)比

圖11 菱形1與菱形2在P1和P2處壁面壓力的比較

比較圖中數(shù)據(jù)可以得出:

①矩形液艙的壁面所受壓力介于2種菱形液艙所受壓力值之間，說(shuō)明2種液艙形式的孰優(yōu)孰劣視具體情況而定。

②菱形液艙壁面壓力的時(shí)間歷程曲線存在很多壓力突變點(diǎn)，而矩形液艙壁面壓力的時(shí)間歷程曲線則表現(xiàn)了很好的周期性。

以上分析表明，采用矩形液艙與菱形液艙的防晃性能優(yōu)劣視具體情況而定，但在一定折角尺寸下，菱形艙室抑制晃蕩作用更明顯。在某個(gè)尺寸范圍外，菱形艙室的折角處面板相當(dāng)于一塊擋板會(huì)對(duì)運(yùn)動(dòng)的來(lái)流施加一個(gè)反沖擊力，故在各點(diǎn)處的壓力值都有所增大;而當(dāng)進(jìn)入某個(gè)尺寸范圍內(nèi)，折角處則起了抑制晃蕩的作用。在LNG船液艙形式上之所以選擇菱形艙室考慮的因素不僅僅是壓力，在某個(gè)壓力差額范圍內(nèi)，菱形艙室可能更有其他方面的優(yōu)勢(shì)，比如折角處較矩形的直角更適合船體結(jié)構(gòu)，或者從避免應(yīng)力集中上更勝一籌，4個(gè)折角處亦可設(shè)為底邊艙、頂邊艙或者壓載水艙。菱形艙室作為一種新型艙室，具有很大的發(fā)展前途，須不斷改善。

5 結(jié)束語(yǔ)

分別從小液深大幅晃蕩下的非線性和隨機(jī)性驗(yàn)證以及大液深微幅晃蕩下2種液艙形式優(yōu)劣比較展開(kāi)討論，應(yīng)用matlab強(qiáng)大的圖表生成功能以及tecplot的自由液面繪制功能進(jìn)行數(shù)據(jù)比較討論，結(jié)果清晰準(zhǔn)確。

［1］ZHU Ren-qing，WU You-sheng，Incecic Atilla.Numerical simulation of liquid sloshing［J］.Shipbuilding of China，2004(3):21－22.

［2］Hirt C W，Nichols B D.Volume of Fluid(VOF)method for the dynamics of free boundaries［J］.J Comput.Phys，1981，39:201－225.

［3］ 劉澤民.關(guān)于船艙中液體晃蕩問(wèn)題［J］.哈爾濱船舶工程學(xué)院學(xué)報(bào)，1992，13(3):241－249.

［4］ 祁江濤.基于VOF法的液艙晃蕩數(shù)值模擬［D］.無(wú)錫:702研究所，2007:32－33.

［5］ 侯玲.液艙晃蕩與彈性防晃結(jié)構(gòu)的相互耦合作用研究［D］.鎮(zhèn)江:江蘇科技大學(xué)，2009:47－48.