基于CFD方法的船舶運(yùn)動(dòng)響應(yīng)仿真

韓翔希+馮志強(qiáng)+邱昂+符妃+范少濤

摘 要：本研究采用CFD方法建立三維波浪水池，并在所模擬的數(shù)值波浪水池中計(jì)算Wigley船型頂浪航行在規(guī)則波時(shí)的運(yùn)動(dòng)響應(yīng)。采用RANS方程求解非定常的不可壓縮粘性流體，并采用VOF方法對(duì)自由面進(jìn)行動(dòng)態(tài)模擬。通過(guò)編寫(xiě)UDF設(shè)置邊界入口速度、波高及在波浪水池尾部1～2倍波長(zhǎng)區(qū)域消波。采用動(dòng)網(wǎng)格技術(shù)結(jié)合六自由度求解器模擬船舶的縱搖與垂蕩運(yùn)動(dòng)。將采用CFD方法模擬得到的Wigley船運(yùn)動(dòng)響應(yīng)與通過(guò)邊界元法所得到的勢(shì)流結(jié)果進(jìn)行對(duì)比，發(fā)現(xiàn)兩者吻合良好。本文所探索的方法可有效計(jì)算船舶在規(guī)則波中的運(yùn)動(dòng)響應(yīng)，便捷高效，易于測(cè)量和控制，并能細(xì)致描述船體周圍的流場(chǎng)情況。同時(shí)，也可為模擬分析其他海洋結(jié)構(gòu)物在波浪中的運(yùn)動(dòng)響應(yīng)具有參考價(jià)值，豐富與擴(kuò)展數(shù)值波浪水池的應(yīng)用。

關(guān)鍵詞：數(shù)值波浪水池；船舶運(yùn)動(dòng)響應(yīng)；VOF；六自由度求解器；動(dòng)網(wǎng)格

中圖分類號(hào)：U661.3 文獻(xiàn)標(biāo)識(shí)碼：A

Simulation of Ship Motion Response Based on CFD

HAN Xiangxi1，2， FENG Zhiqiang1，2， QIU Ang3， FU Fei1，2 FAN Shaotao4

（ 1. Guangxi Engineering Technology Research Center of Marine Digital Design and Advanced Manufacturing， Qinzhou 535099；

2. Qinzhou University， Qinzhou 535099； 3. Guangdong Shipping Science Research Institute， Guangzhou 510000；

4. China Energy Engineering Group Guangdong Electric Power Design Institute Co.， Ltd. Guangzhou 510663 ）

Abstract： This paper establishes a three-dimensional （3D） numerical wave tank based on the theory of viscous flow to simulate the unsteady motion response of a Wigley ship navigating in regular heading waves. The governing equations， Reynolds Averaged Navier-Stokes and continuity equations are discretized by finite volume method， A six-degrees-of-freedom （6 DOF） solver is employed to predict the motions of the ship， and volume of fluid （VOF） method is adopted to capture the nonlinear free surface by writing user-defined functions （UDF） to set inlet boundary. The outgoing waves are dissipated inside an artificial damping zone located at the rear part （about 1-2 wave lengths） of the wave tank. The numerical simulation results show good agreement with the theoretical results. This research can be used to further analyze and predict hydrodynamic performance of the ship and marine floating structures in waves and help to extend the applications of the numerical wave tank.

Key words： Numerical wave tank； Ship motion response； Volume of fluid； Six DOF Solver； Dynamic mesh

1 引言

近年來(lái)，由于CFD技術(shù)的發(fā)展，基于CFD方法的數(shù)值造波技術(shù)逐漸成為研究熱點(diǎn)。Harlow等人提出的MAC方法和Hirt等提出的VOF方法，推動(dòng)了粘性自由表面數(shù)值計(jì)算技術(shù)的快速發(fā)展。此后，許多學(xué)者在粘性自由液面的數(shù)值水池上做了大量工作。Wang采用VOF方法模擬出二維數(shù)值波浪水池，并開(kāi)展波浪對(duì)近海平臺(tái)底部沖擊過(guò)程的研究；齊鵬和王永學(xué)在前人二維波浪水池的基礎(chǔ)上，采用VOF的方法，根據(jù)造波機(jī)的原理建立三維數(shù)值波浪水池；Xing-Kaeding等人基于RANS方法模擬計(jì)算船舶在波浪中的運(yùn)動(dòng)響應(yīng)，自由面同樣采用VOF方法進(jìn)行處理。

二維數(shù)值水槽的建立及分析已經(jīng)相當(dāng)成熟，但應(yīng)用到工程項(xiàng)目中帶有局限性，所以建立三維波浪數(shù)值模型分析計(jì)算波浪與結(jié)構(gòu)物相互作用具有特別重要的意義。另外，目前計(jì)算波浪作用在結(jié)構(gòu)物上的載荷大多通過(guò)線性勢(shì)流理論解決，但線性勢(shì)流理論有明顯的局限性，很難考慮粘性和非線性的影響，如果能夠模擬高階非線性波且結(jié)果可靠，則可用數(shù)值波浪水池模擬各種波況下、特別是極端波況下的波浪運(yùn)動(dòng)特性。本研究基于CFD建立三維數(shù)值波浪水池，自由液面采用VOF處理，通過(guò)編寫(xiě)UDF設(shè)置邊界入口速度及動(dòng)量方程中添加源項(xiàng)實(shí)現(xiàn)阻尼消波，采用動(dòng)網(wǎng)格技術(shù)結(jié)合六自由度求解器模擬船舶的縱搖與垂蕩運(yùn)動(dòng)，模擬了規(guī)則波頂浪中三維船體縱搖與垂蕩運(yùn)動(dòng)的計(jì)算，并與勢(shì)流理論邊界元法得到的結(jié)果相比較，兩者吻合很好。通過(guò)試驗(yàn)，本方法能準(zhǔn)確預(yù)測(cè)規(guī)則波頂浪中船舶的縱搖與垂蕩運(yùn)動(dòng)響應(yīng)，能細(xì)致描述船舶周圍的流場(chǎng)。本方法也可為模擬分析其他的海洋結(jié)構(gòu)物在波浪中的運(yùn)動(dòng)響應(yīng)具有參考價(jià)值，豐富與擴(kuò)展數(shù)值波浪水池的應(yīng)用。endprint

2 數(shù)學(xué)模型及數(shù)值方法

2.1 控制方程

本文采用兩相不可壓的RANS方程求解非定常不可壓粘性流體，RANS方程采用 SST k-ω湍流模型封閉，自由液面采用VOF法處理，消波區(qū)域采用在動(dòng)量方程中添加源項(xiàng)的辦法進(jìn)行消波。

連續(xù)方程： （1）

RANS方程： （2）

式中：U為速度矢量；流體密度定義為 ，其中體積分?jǐn)?shù)aq表示單元內(nèi)第q相流體占的體積與該單元的體積之比，本文采用氣液兩相模擬q=1，2，并且

；μ為相體積分?jǐn)?shù)平均的動(dòng)力粘性系數(shù)。

Menter提出了SST k-ω 兩方程模型，它在近壁面保留了原始k-ω的模型，在遠(yuǎn)離壁面的地方應(yīng)用了k-ω模型，其k方程和ω方程如下：

（3）

（4）

自由液面采用VOF法處理，aq滿足方程：

（5）

2.2 數(shù)值造波和消波方法

本文采用在邊界設(shè)置入口速度的方法造波，需要編寫(xiě)UDF控制入口速度以及波高。

根據(jù)線性波理論，小振幅波的速度場(chǎng)為：

（6）

（7）

小振幅波的自由液面高度

（8）

式中：U0為船舶航速；A是波幅；k為波數(shù)；ω是波浪圓頻率；ωe是遭遇頻率；x軸為波浪傳播方向。

為了避免入射波與反射波相互干擾，導(dǎo)致波浪出現(xiàn)破碎或不規(guī)則現(xiàn)象，采用消波的方法來(lái)減弱反射波。本文采用由Baker等人在二維邊界條件模式中首次提出的阻尼消波技術(shù)，從而可以實(shí)現(xiàn)在有限計(jì)算區(qū)域內(nèi)原本在無(wú)限計(jì)算域才能實(shí)現(xiàn)的消波效果，這就是所謂的數(shù)值阻尼消波原理。本文在計(jì)算域的尾部區(qū)域設(shè)置人工阻尼消波區(qū)，根據(jù)其他學(xué)者的經(jīng)驗(yàn)，消波區(qū)的長(zhǎng)度取波長(zhǎng)的1-2倍。

在阻尼消波區(qū)域內(nèi)，動(dòng)量方程為：

（9）

（10）

式中：μ（x）為在阻尼段起點(diǎn)為零的單調(diào)遞增函數(shù)（線性遞增、指數(shù)遞增等）。

3 模型建立

3.1 計(jì)算區(qū)域

用來(lái)計(jì)算規(guī)則波頂浪中Wigley船模縱搖與垂蕩運(yùn)動(dòng)的數(shù)值波浪水池尺度為：-1.5Lpp

3.2 邊界條件及FLUENT設(shè)置

數(shù)值波浪水池的邊界條件設(shè)置如下：水池前端為速度入口，水池后端為壓力出口，上邊界為壓力進(jìn)口，下邊界和兩個(gè)側(cè)面為對(duì)稱邊界。本文利用FLUENT軟件、使用VOF方法、結(jié)合SST k-ω模型、壓力耦合方程組的SIMPLE算法，對(duì)流項(xiàng)采用二階迎風(fēng)差分格式，求解非定常狀態(tài)下的紊流問(wèn)題。應(yīng)用DEFINE-PROFILE宏設(shè)置速度入口和VOF控制二相分布；應(yīng)用DEFINE SOURCE宏添加動(dòng)量源項(xiàng)實(shí)現(xiàn)數(shù)值波浪水池尾部區(qū)域的消波功能。采用動(dòng)網(wǎng)格技術(shù)結(jié)合六自由度求解器模擬船舶的縱搖與垂蕩運(yùn)動(dòng)，最終模擬Wigley船模在規(guī)則波中頂浪航行時(shí)的縱搖與垂蕩運(yùn)動(dòng)。

4 結(jié)果與分析

對(duì)航行于波浪中的船舶，一般情況下頂浪時(shí)的縱搖和垂蕩最嚴(yán)重，因此本文以模擬計(jì)算頂浪中船舶縱搖和垂蕩運(yùn)動(dòng)響應(yīng)作為所建立數(shù)值波浪水池的應(yīng)用例子。考慮到與試驗(yàn)結(jié)果進(jìn)行對(duì)比分析，在計(jì)算船模運(yùn)動(dòng)響應(yīng)時(shí)，入射波采用線性規(guī)則波。設(shè)波長(zhǎng)λ=1.55 m、4.5 m和6.36 m，波幅A=0.045 m，通過(guò)模擬計(jì)算分別得到船模的垂蕩和縱搖的時(shí)歷曲線，如圖3和圖4所示。

從圖3和圖4可以看出，在短波情況下（λ/L=0.516 7），其垂蕩的幅值相比波浪的波幅較小，縱搖角度同樣很小；在中等波長(zhǎng)情況下（λ/L=1.5），垂蕩運(yùn)動(dòng)幅值并未達(dá)到最大，但縱搖固有周期約等于波浪對(duì)船體的擾動(dòng)力周期，發(fā)生諧搖，縱搖幅值達(dá)到最大；在最大波長(zhǎng)情況下（λ/L=2.12），其垂蕩運(yùn)動(dòng)幅值比中等波長(zhǎng)情況大，縱搖運(yùn)動(dòng)幅值比中等波長(zhǎng)情況小，此時(shí)已經(jīng)處于隨波逐流狀態(tài)。

可以看出，對(duì)于垂蕩和縱搖運(yùn)動(dòng)而言，規(guī)則波對(duì)船體擾動(dòng)力的大小與波長(zhǎng)對(duì)船長(zhǎng)之比有很大關(guān)系：λ/L越小，垂蕩和縱搖越緩和，故大船在小的波浪上垂蕩和縱搖都不會(huì)很大；最大縱搖角度出現(xiàn)在中等波長(zhǎng)下，此時(shí)波長(zhǎng)略大于船長(zhǎng)，最大波長(zhǎng)的縱搖幅值介于小波長(zhǎng)和中等波長(zhǎng)之間；對(duì)于最大垂蕩幅值，筆者用勢(shì)流理論邊界元法得到在入射波波長(zhǎng)為20～28 m、波幅0.045 m時(shí)，船模垂蕩運(yùn)動(dòng)的幅值約等于波浪的幅值。

基于CFD方法的數(shù)值波浪水池可以處理復(fù)雜的運(yùn)動(dòng)響應(yīng)問(wèn)題，其計(jì)算結(jié)果考慮了粘性及非線性因素等的影響。為了較容易定量分析和比較計(jì)算結(jié)果，本文僅考慮線性因素，即只考慮船模的運(yùn)動(dòng)幅度是線性的，在這種情況下用最小二乘法以正余弦函數(shù)對(duì)運(yùn)動(dòng)幅度時(shí)歷進(jìn)行擬合，可以求出船體垂蕩和縱搖幅值。

規(guī)則波頂浪中船模的縱搖和垂蕩最嚴(yán)重，可以將垂蕩和縱搖運(yùn)動(dòng)表示成如下形式：

zG= za cos（ωet=ε3） （11）

θ= θa cos（ωet=ε3） （12）

式中：ZG和θ分別為船體的垂蕩和縱搖運(yùn)動(dòng)隨時(shí)間變化的函數(shù)；Za和θa分別為船體垂蕩運(yùn)動(dòng)幅值和縱搖運(yùn)動(dòng)幅值；ωe為遭遇頻率；ε3和ε5為相應(yīng)的相位。

在波幅A=0.045 m、頂浪海況下零航速時(shí)Wigley船模的垂蕩和縱搖運(yùn)動(dòng)響應(yīng)幅值曲線，如圖5和圖6所示，圖中同時(shí)給出了邊界元法得到的勢(shì)流理論結(jié)果。從圖5、圖6可知，CFD數(shù)值計(jì)算結(jié)果與勢(shì)流理論解吻合良好，CFD計(jì)算考慮了粘性影響，故CFD數(shù)值解總體上比勢(shì)流解要小。

圖7給出了波幅A=0.045 m和λ/L=0.516 7、λ/L=1.5、λ/L=2.12三種情況下，在第11 s時(shí)Wigley船模附近數(shù)值水池的波形圖。從圖中可以看出尾部消波區(qū)域自由面相對(duì)平整，說(shuō)明本文使用的消波方法效果很好。

在短波海況λ=1.55 m，A=0.045 m和Fr數(shù)分別為0.1、0.15、0.2、0.25、0.3時(shí)的運(yùn)動(dòng)響應(yīng)曲線，見(jiàn)圖8和圖9。由圖8可知，在Fr=0.1時(shí)垂蕩運(yùn)動(dòng)幅值最大；在Fr>0.1時(shí)船模垂蕩運(yùn)動(dòng)幅值隨著航速的增加而減緩。由圖9可知，縱搖運(yùn)動(dòng)的最大幅值出現(xiàn)在無(wú)航速情況下，并且縱搖運(yùn)動(dòng)幅值隨著航速的增加而減緩。

圖10給出了A=0.045 m、λ=1.55 m、Fr=0.2時(shí)，在一個(gè)周期內(nèi)（t=0 T， t=0.25 T， t=0.5 T， t=0.75 T） Wigley船模底部的壓力分布云圖。從圖中可以看出各時(shí)刻船模底部壓力值過(guò)渡緩和，說(shuō)明本文對(duì)Wigley船模附近網(wǎng)格劃分處理得當(dāng)。

圖11給出了同一個(gè)周期內(nèi)不同時(shí)刻自由面波形的變化以及Wigley船模的垂蕩和縱搖。從圖中可以清晰看到水線面的位置，整個(gè)船長(zhǎng)范圍內(nèi)波峰與波谷十分清晰且與圖10中對(duì)應(yīng)時(shí)刻船底壓力分布是一致的，即處于波峰的船體部分對(duì)應(yīng)的底部壓力值高，波谷部分則相反。

圖12給出了A=0.045 m、λ=1.55 m、不同F(xiàn)r數(shù)對(duì)應(yīng)的波形圖。從圖中可以看出，在航速較小情況下，數(shù)值水池波形完整，船尾基本上沒(méi)有出現(xiàn)興波；隨著航速加大，船行波漸漸明顯，由船行波與規(guī)則波的相互作用疊加形成的波形和原來(lái)規(guī)則波在船后尾流區(qū)域錯(cuò)分開(kāi)來(lái)，這些現(xiàn)象都與實(shí)際中觀察到的現(xiàn)象一致。造成這些現(xiàn)象的原因是船模航行時(shí)波流共同作用于彎曲的船體，沿船體表面的壓力分布不一樣，導(dǎo)致船體周圍的水面升高或下降，在重力和慣性作用下在船后形成實(shí)際的船波，并與規(guī)則波相互作用。當(dāng)Fr=0.25和0.3時(shí)，可以明顯看到船體尾部形成的凱爾文角，整個(gè)船波系基本集中在凱爾文角所限定的扇形面范圍內(nèi)。

綜上所述，本文CFD模擬計(jì)算結(jié)果與勢(shì)流理論邊界元法得到的結(jié)果吻合較好。一方面，相對(duì)于勢(shì)流理論的計(jì)算，CFD方法更能真實(shí)地反映模擬流場(chǎng)，計(jì)算結(jié)果更加精確；另一方面，CFD方法更加便捷高效，具有易測(cè)量和控制及能細(xì)致描述船體周圍的流場(chǎng)等情況優(yōu)點(diǎn)，因此采用CFD方法預(yù)報(bào)波浪中船模的性能比物理模型試驗(yàn)更具有吸引力。

5 結(jié)語(yǔ)

本研究基于CFD方法，建立了三維粘性數(shù)值波浪水池，并用其模擬了規(guī)則波頂浪中以多個(gè)航速前進(jìn)的Wigley船模的縱搖與垂蕩運(yùn)動(dòng)響應(yīng)。結(jié)果表明：本文的數(shù)值模擬方法能準(zhǔn)確地反映波浪與船模之間的干擾作用，并能給出較精確的垂蕩和縱搖運(yùn)動(dòng)幅值，表明了在數(shù)值波浪水池中計(jì)算船舶運(yùn)動(dòng)響應(yīng)的有效性；與物理模型實(shí)驗(yàn)相比，本文的方法更加便捷高效，具有易測(cè)量和控制及能細(xì)致描述船體周圍的流場(chǎng)情況等優(yōu)點(diǎn)，同時(shí)，該方法也為分析其他的海洋結(jié)構(gòu)物在波浪中的運(yùn)動(dòng)響應(yīng)具有參考價(jià)值，豐富與擴(kuò)展了數(shù)值波浪水池的應(yīng)用。

參考文獻(xiàn)

[1] 戴世強(qiáng).認(rèn)識(shí)水波，駕馭水波-評(píng)介陶建華的專著《水波的數(shù)值模

擬》[J].水動(dòng)力學(xué)研究與進(jìn)展，2007，22 （2）.

[2] Longuet-Higgins MS， Cokelet CD. The deformation of steep surface waves on

water， a numerical method of computation [C]. Proceedings of the Royal Society

of London 1976， A350.

[3] Baudic SF， Williams AN， Kareem A. A two-dimensional numerical wave

flume-Part 1： Nonlinear wave generation propagation and absorption [J]. Journal

of Offshore Mechanics and Arctic Engineering 2001， 123.

[4] Ryu S， Kim MH， Lynett PJ. Fully nonlinear wave-current interactions and

kinematics by a BEM-based numerical wave tank [J]. Computational Mechanics

2003， 32.

[5] Harlow FH， Welch JE. Numerical calculation of time-dependent viscous

incompressible flow of fluid with free surface [J]. Physics of Fluids 1965， 8.

[6] Hirt CW， Nichols BD. Volume of fluid （VOF） method for the dynamics of free

boundaries [J]. Journal of Computational Physics 1981， 3.

[7] Wang YX. Numerical wave channel with absorbed wave maker [J]. China

Ocean Engineering 1995， 9（2）.

[8] Wang YX. Wave slamming on offshore structure [A]. Proceedings of the 6th

International Off shore and Polar Engineering Conference： Vol 3 [C]. Los

Angeles， US， 1996. 1922196.

[9] 齊鵬，王永學(xué).三維數(shù)值波浪水池技術(shù)與應(yīng)用 [J]. 大連理工大學(xué)學(xué)報(bào)

2003， 43（6）.

[10] Xing-Kaeding Y， Jensen G， Hadzic I， Peric M. Simulation of flow-induced

ship motions in waves using a RANSE method [J]. Ship Technology Research

2004， 51.

[11] Menter FR. Zonal two equation k-w turbulence models for aerodynamic flows

[R]. AIAA-93-2906， 1993.

[12] Baker GR， Meiron DI， Orszag SA. Application of generalized vortex method to

nonlinear free-surface flows [C]. Proceedings of the 3rd International

Conference on Numerical Ship Hydrodynamics. Paris， France， 1981.

[13] 孫大鵬，李玉成.數(shù)值水槽內(nèi)的阻尼消波和波浪變形計(jì)算 [J].海洋

工程，2000，18 （2）.

[14] 高學(xué)平，曾廣冬，張亞.不規(guī)則波浪數(shù)值水槽的造波和阻尼消波 [J].

海洋學(xué)報(bào) 2002，24 （2）.

[15] 齊鵬，王永學(xué).三維數(shù)值波浪水池技術(shù)與應(yīng)用 [J].大連理工大學(xué)學(xué)

報(bào) 2003，43 （6）.endprint