何廣華, 陳麗敏, 王佳東

(哈爾濱工業(yè)大學(xué)(威海) 船舶與海洋工程學(xué)院，山東 威海 264209)

船舶在波浪中運(yùn)動(dòng)的強(qiáng)非線(xiàn)性時(shí)域模擬

何廣華, 陳麗敏, 王佳東

(哈爾濱工業(yè)大學(xué)(威海) 船舶與海洋工程學(xué)院，山東 威海 264209)

為準(zhǔn)確預(yù)報(bào)在非線(xiàn)性波浪中航行船舶的大幅運(yùn)動(dòng)，基于CIP方法建立了波浪中船舶強(qiáng)非線(xiàn)性運(yùn)動(dòng)響應(yīng)的時(shí)域分析模型.該模型基于粘性流理論，對(duì)N-S方程進(jìn)行求解，可以考慮船舶的運(yùn)動(dòng)響應(yīng)和自由面變形過(guò)程中的粘性效應(yīng)和強(qiáng)非線(xiàn)性效應(yīng).首先，通過(guò)模擬船舶的強(qiáng)制垂蕩運(yùn)動(dòng)，對(duì)該數(shù)值模型的收斂性展開(kāi)了研究；然后，基于該模型計(jì)算了不同入射波幅和不同入射波長(zhǎng)對(duì)S175船模的垂蕩和縱搖運(yùn)動(dòng)響應(yīng)的影響；最后，模擬了S175船模在大幅波浪中的強(qiáng)非線(xiàn)性運(yùn)動(dòng)響應(yīng)現(xiàn)象.模擬結(jié)果表明：基于自編程建立的船舶強(qiáng)非線(xiàn)性運(yùn)動(dòng)響應(yīng)模型具有較好的收斂性和精確性，與其他數(shù)值計(jì)算結(jié)果吻合較好；尤其擅長(zhǎng)模擬水花飛濺、甲板上浪、底部砰擊等強(qiáng)非線(xiàn)性問(wèn)題.

CIP；強(qiáng)非線(xiàn)性；運(yùn)動(dòng)響應(yīng)；大幅度運(yùn)動(dòng)；時(shí)域模擬

船舶耐波性的準(zhǔn)確預(yù)報(bào)是船舶設(shè)計(jì)的關(guān)鍵要素之一.船舶在波浪中航行時(shí)，其產(chǎn)生的垂蕩和縱搖運(yùn)動(dòng)響應(yīng)會(huì)明顯增加船舶阻力；船舶的大幅運(yùn)動(dòng)會(huì)造成螺旋槳出水，推進(jìn)效率降低；同時(shí)還會(huì)影響船上設(shè)備的正常運(yùn)行，嚴(yán)重時(shí)導(dǎo)致設(shè)備的損壞，甚至船舶會(huì)因此短暫或徹底失去控制，發(fā)生傾覆.關(guān)于船舶耐波性的研究，可分為物理模型試驗(yàn)和數(shù)值模擬研究.物理模型試驗(yàn)最大的優(yōu)點(diǎn)是可靠性高，但其具有成本高、周期長(zhǎng)、難以進(jìn)行系統(tǒng)性研究等缺點(diǎn).此外，對(duì)于甲板上浪、砰擊、水花飛濺等強(qiáng)非線(xiàn)性現(xiàn)象，實(shí)驗(yàn)條件和測(cè)量精度仍不成熟.因此目前關(guān)于強(qiáng)非線(xiàn)性船舶運(yùn)動(dòng)響應(yīng)物理模型試驗(yàn)的報(bào)道仍較為罕見(jiàn)；數(shù)值模擬具有經(jīng)濟(jì)性、快捷性、靈活性等優(yōu)點(diǎn)，傳統(tǒng)的船舶耐波性計(jì)算方法有切片理論法[1-2]，EUT (enhanced unified theory)法[3]，Rankine源法[4-10]等.這些計(jì)算工具大都是基于勢(shì)流理論求解的，忽略了流體的粘性作用，無(wú)法考慮船體大幅運(yùn)動(dòng)、波浪破碎、水花飛濺、甲板上浪等強(qiáng)非線(xiàn)性現(xiàn)象，而這些強(qiáng)非線(xiàn)性效應(yīng)對(duì)船舶的運(yùn)動(dòng)響應(yīng)有著顯著的影響.相對(duì)于勢(shì)流方法，基于粘性流理論的CFD (computational fluid dynamics)，可充分考慮船舶航行過(guò)程中的粘性效應(yīng)和強(qiáng)非線(xiàn)性效應(yīng).近年來(lái)，得到了國(guó)內(nèi)外廣大學(xué)者的關(guān)注，并取得了若干研究成果.如：沈志榮等[11]基于開(kāi)源代碼OpenFOAM開(kāi)發(fā)的數(shù)值模型計(jì)算分析了Wigley Ⅲ 船模在迎浪中的運(yùn)動(dòng)響應(yīng).石博文等[12]基于粘性流理論建立了三維數(shù)值水池，分析了DTMB5512船模在不規(guī)則波浪中頂浪航行時(shí)的垂蕩與縱搖運(yùn)動(dòng)時(shí)間歷程曲線(xiàn).

CIP (constrained interpolation profile)方法由Takewaki等[13]于1985年首先提出；Yabe等[14-15]將該方法推廣到物理學(xué)、電磁學(xué)和激光學(xué)等學(xué)科.Hu等[16]于2004年首次將CIP方法應(yīng)用到船舶與海洋工程領(lǐng)域，用來(lái)計(jì)算極端非線(xiàn)性自由液面的流動(dòng).隨后，He等[17-20]將其應(yīng)用到處理潰壩問(wèn)題、波浪與建筑物相互作用、船體耐波性等海洋工程和船舶工程中的強(qiáng)非線(xiàn)性問(wèn)題.CIP法利用三次多項(xiàng)式進(jìn)行插值，并同時(shí)對(duì)物理量和其空間導(dǎo)數(shù)進(jìn)行求解，保證了數(shù)值求解在時(shí)間、空間上的三階精度.其具有精度高、穩(wěn)定性好、數(shù)值耗散小等優(yōu)點(diǎn).本文采用自編程方式，基于CIP法建立了船舶耐波性強(qiáng)非線(xiàn)性分析模型；并基于該數(shù)值模型對(duì)波浪中航行的S175船模的運(yùn)動(dòng)響應(yīng)展開(kāi)了研究，分析了不同入射波高和不同入射波長(zhǎng)對(duì)S175船模的垂蕩和縱搖運(yùn)動(dòng)的影響，計(jì)算結(jié)果與其他數(shù)值方法進(jìn)行了對(duì)比驗(yàn)證.最后，計(jì)算了S175船模在大浪中的大幅度運(yùn)動(dòng)響應(yīng)，結(jié)果表明，該計(jì)算模型能較好地模擬船體大幅運(yùn)動(dòng)、甲板上浪、底部砰擊、水花飛濺等強(qiáng)非線(xiàn)性現(xiàn)象.

1 理論與數(shù)值模型

本研究采用自編程技術(shù)，基于CIP方法建立船舶耐波性分析的粘性流數(shù)值模型.文獻(xiàn)[17]采用CIP數(shù)值模型對(duì)Wigley船的強(qiáng)制垂蕩運(yùn)動(dòng)和縱搖運(yùn)動(dòng)進(jìn)行了線(xiàn)性研究，計(jì)算結(jié)果與實(shí)驗(yàn)值吻合較好.限于篇幅，本文不再詳細(xì)、全面地?cái)⑹鲈揅IP模型的基礎(chǔ)理論和數(shù)值技術(shù)；關(guān)于該部分，二維模型具體可參考文獻(xiàn)[16,20]，三維模型可參考文獻(xiàn)[17].

1.1 控制方程

模型以連續(xù)性方程、N-S方程為基本控制方程.不考慮溫度變化，假設(shè)流體是不可壓縮的，可以得到以下方程：

(1)

(2)

式中：ρ為密度；ui(i=1, 2, 3)為速度；μ為粘性系數(shù);Sij=?ui/?xj+?uj/?xi；fi為質(zhì)量力，如：重力等. 其中式(2)可用分裂法，將其可分成：

(3)

(4)

(5)

對(duì)流項(xiàng)，式(3)采用CIP法計(jì)算；擴(kuò)散項(xiàng)，式(4)采用中心差分法；壓力-速度耦合項(xiàng)，式(5)通過(guò)引入連續(xù)性方程，可得以下泊松方程:

(6)

式(6)為壓力泊松方程，同時(shí)適用于固體、液體、氣體, 可采用SOR(超松弛迭代法)來(lái)迭代求解.

1.2CIP法

CIP法由文獻(xiàn)[13]首先提出；文獻(xiàn)[16]將其用于船舶與海洋工程領(lǐng)域的強(qiáng)非線(xiàn)性問(wèn)題的求解.隨后文獻(xiàn)[17-20]將其在海洋工程和船舶工程中進(jìn)一步推廣.

圖1給出了一維CIP方法的原理圖.為降低數(shù)值耗散，CIP法通過(guò)補(bǔ)償網(wǎng)格間的內(nèi)部信息來(lái)提高數(shù)值精度，基本思想為：對(duì)于對(duì)流項(xiàng)變量f，不僅要計(jì)算輸運(yùn)方程的值，而且要計(jì)算輸運(yùn)方程的空間梯度(圖1(d))，即gi=?f/?xi，其中f可為ρ、ui和p.f的輸運(yùn)方程可以寫(xiě)成:

(7)

將式 (7)對(duì)空間坐標(biāo)進(jìn)行離散，可以得到gi的輸運(yùn)方程：

(8)

式(8)的計(jì)算可以分為對(duì)流項(xiàng)和非對(duì)流項(xiàng)兩步，對(duì)于式(7)和式(8)的對(duì)流項(xiàng)計(jì)算，可以用以下半拉格朗日方法求得(如圖2所示)：

圖1 一維CIP法原理示意

式中，多項(xiàng)式的系數(shù)ai，bi，ci，di，可通過(guò)強(qiáng)加在xi-1和xi兩個(gè)網(wǎng)格點(diǎn)處的物理量fn的值以及其梯度值gn共4個(gè)量來(lái)確定.

圖2 CIP法的半拉格朗日時(shí)間步進(jìn)法

1.3 界面捕捉技術(shù)

為實(shí)現(xiàn)自由面的精確重構(gòu)，本數(shù)值模型選擇高精度THINC法[21]來(lái)重構(gòu)自由面.其基本思路是將流體的密度函數(shù)近似為雙曲正切函數(shù)，然后通過(guò)半拉格朗日法積分得到流量，詳細(xì)論述請(qǐng)參考文獻(xiàn)[21].流體的密度函數(shù)的插值函數(shù)如下：

1)當(dāng)ui+1/2≥0時(shí),

2)當(dāng)ui+1/2≤0時(shí),

式中：F(x) 為插值函數(shù)；u為流速；g為流量； 未知參數(shù)α、β、γ、δ可通過(guò)等量關(guān)系和經(jīng)驗(yàn)公式得到.圖3給出了THINC法的示意圖(ui+1/2≤0).

圖3 THINC法插值函數(shù)示意

2 數(shù)值模擬與結(jié)果分析

在計(jì)算波浪中船舶的運(yùn)動(dòng)響應(yīng)時(shí)，船體一般可看作6自由度運(yùn)動(dòng)的剛體，且需要定義兩個(gè)坐標(biāo)系統(tǒng)：一個(gè)是大地坐標(biāo)系，另一個(gè)是隨體坐標(biāo)系.本模型中采用的計(jì)算域隨著船體航速移動(dòng)，船體的運(yùn)動(dòng)響應(yīng)也在該隨體坐標(biāo)系下求解.該坐標(biāo)系采用右手笛卡爾坐標(biāo)，如圖4所示. 原點(diǎn)置于在船中與靜水面的交點(diǎn)，x-方向以船尾指向船首為正，z-方向以豎直向上為正.

2.1 船體模型

本文采用ITTC-2010[22](international towing tank conference)算例對(duì)比中的標(biāo)準(zhǔn)船型——S175.S175是一艘長(zhǎng)為175 m的集裝箱船，其具體主尺度參數(shù)見(jiàn)表1, 船體表面的重構(gòu)采用虛擬粒子法，其網(wǎng)格(粒子)模型如圖5所示. 計(jì)算域?yàn)榉叫?，如圖6所示.計(jì)算域的大小x、y、z方向分別設(shè)置為：-2.7L～5.9L, -1.1L～1.1L, -0.8L～1.3L；其中L為船長(zhǎng)，已在表1中列出.采用非均勻網(wǎng)格，如：在x、y、z方向的網(wǎng)格數(shù)分別為221×81×81.為提高計(jì)算精度，在靠近船首附近采用精細(xì)網(wǎng)格，最小的網(wǎng)格尺度達(dá)Δx=0.010L；在靠近船體和自由水面處也采用較精細(xì)網(wǎng)格，最小的網(wǎng)格尺度為Δy=Δz=0.005L.

圖4 計(jì)算模型

表1 S175集裝箱船的主尺度參數(shù)

2.2 收斂性驗(yàn)證

通過(guò)模擬S175船模的強(qiáng)制垂蕩運(yùn)動(dòng)，對(duì)建立的數(shù)值模型展開(kāi)了收斂性研究. 計(jì)算域大小x、y、z方向分別設(shè)置為：-2.7L～5.9L, -1.1L～1.1L, -0.8L～1.3L. 為保證數(shù)值計(jì)算的穩(wěn)定性, 船舶從靜止?fàn)顟B(tài)逐漸加速到設(shè)定航速U. 鑒于篇幅，收斂性研究只展示了對(duì)數(shù)值模擬影響較大的關(guān)鍵參數(shù)的計(jì)算結(jié)果, 如：網(wǎng)格數(shù)量和時(shí)間步.

船舶強(qiáng)制垂蕩運(yùn)動(dòng)可以表示為

式中:j= 3表示垂蕩運(yùn)動(dòng)模態(tài)；為垂蕩幅值;ωe為垂蕩運(yùn)動(dòng)圓頻率. 將縱蕩力Fx和垂蕩力Fz以ρωe進(jìn)行無(wú)因次化, 其中，ρ為流體密度,為排水體積; 時(shí)間的無(wú)因次化因子為T(mén)e= 2π/ωe(Te為船舶在波浪中的遭遇周期或強(qiáng)制運(yùn)動(dòng)周期). 圖7～8為S175船模在航速Fn= 0.2處以頻率為kL=14、振幅為kα=0.28進(jìn)行強(qiáng)制垂蕩運(yùn)動(dòng)的數(shù)值模擬結(jié)果.

圖5 船體模型(粒子法)

圖6 計(jì)算域與網(wǎng)絡(luò)劃分

圖7給出了3種網(wǎng)格下無(wú)因次化垂蕩力的計(jì)算結(jié)果.3種網(wǎng)格數(shù)量分別約為Mesh-1 (1218181)、Mesh-2 (1818181) 和Mesh-3 (2218181). 相比于Mesh-1，Mesh-2與Mesh-3的計(jì)算結(jié)果更加接近, 為保證計(jì)算精度，以下的數(shù)值模擬采用Mesh-3 (2218181).

圖8給出了3種不同時(shí)間步長(zhǎng)對(duì)無(wú)因次化垂蕩力的影響，時(shí)間步長(zhǎng)分別取為T(mén)e/1 000,Te/2 000和Te/4 000, 計(jì)算結(jié)果并沒(méi)有較大差異,Te/2 000更接近于Te/4 000.綜合考慮計(jì)算的時(shí)間和精度, 在以下算例中時(shí)間步長(zhǎng)采用Te/2 000.

圖7 計(jì)算網(wǎng)格對(duì)垂蕩力的影響 (Fn =0.2,kL=14,kα=0.28)

Fig.7 Sensitivity of hydrodynamic force to domain grids for the S175 model

2.3 波浪中航行船舶的運(yùn)動(dòng)響應(yīng)

經(jīng)過(guò)收斂性研究, 以下數(shù)值模擬均采用Mesh-3網(wǎng)格和Te/2 000的時(shí)間步長(zhǎng). 本文對(duì)S175船模在波浪中以航速Fn=0.2迎浪航行的工況展開(kāi)了系統(tǒng)的研究，分析了不同入射波幅和不同入射波長(zhǎng)對(duì)船舶運(yùn)動(dòng)響應(yīng)的影響.鑒于目前關(guān)于強(qiáng)非線(xiàn)性船舶運(yùn)動(dòng)響應(yīng)物理模型試驗(yàn)仍較為罕見(jiàn)，ITTC于2010年以S175為標(biāo)準(zhǔn)船模進(jìn)行了多家數(shù)值模型計(jì)算結(jié)果的綜合比對(duì)，供后繼程序的驗(yàn)證與校驗(yàn)之用.在本文中每個(gè)波況均計(jì)算了12個(gè)運(yùn)動(dòng)周期Te, 選取穩(wěn)定后的最后2個(gè)周期的垂蕩和縱搖運(yùn)動(dòng)響應(yīng)時(shí)間歷程線(xiàn)，與ITTC-2010[22]中的兩種數(shù)值模型的結(jié)果展開(kāi)了對(duì)比.計(jì)算采用主頻為3.40 GHz的Intel I3-3240 CPU, 單核計(jì)算占用內(nèi)存為680 M, 以Te/2 000的時(shí)間步長(zhǎng)計(jì)算12個(gè)運(yùn)動(dòng)周期Te所需計(jì)算時(shí)間約為11 h(隨輸出數(shù)據(jù)類(lèi)型及數(shù)據(jù)量的大小有所變化，如：Tecplot流場(chǎng)數(shù)據(jù)等).

圖8 時(shí)間步長(zhǎng)對(duì)垂蕩力的影響 (Fn =0.2,kL=14,kα=0.28)

Fig.8 Sensitivity of hydrodynamic force to time step for the S175 model

在圖9～11中, 帶實(shí)心方點(diǎn)的曲線(xiàn)為本數(shù)值模型的計(jì)算結(jié)果, 三角點(diǎn)為NMRI (national maritime research institute)的數(shù)值模擬結(jié)果; ‘×’ 為NTUA (national technical university of athens) 的數(shù)值模擬結(jié)果. NMRI所開(kāi)發(fā)的數(shù)值模型是基于弱非線(xiàn)性理論，采用二維邊界元方法(2D-BEM)進(jìn)行求解.而NTUA開(kāi)發(fā)的數(shù)值模型是基于三維邊界元方法，可以考慮線(xiàn)性和弱非線(xiàn)性問(wèn)題.圖9～11給出了S175船在入射波為λ/L=1.0中，以Fn=0.2航行時(shí)的垂蕩和縱搖運(yùn)動(dòng)的響應(yīng)時(shí)間歷程線(xiàn); 圖9～11(a) 和 (b) 分別為無(wú)因次化的垂蕩運(yùn)動(dòng)和縱搖運(yùn)動(dòng)的計(jì)算結(jié)果,kA為無(wú)因次化入射波波幅.

從圖9～11可以看出：1) 3種數(shù)值方法計(jì)算的垂蕩和縱搖運(yùn)動(dòng)響應(yīng)結(jié)果總體吻合較好.2) 從圖9(a)～11(a)可以看出：入射波幅對(duì)船體垂蕩運(yùn)動(dòng)的影響.當(dāng)入射波幅較小時(shí) (kA=0.04, 圖9 (a))，船體運(yùn)動(dòng)的非線(xiàn)性現(xiàn)象不是太明顯，仍可以按線(xiàn)性理論進(jìn)行考慮; NTUA的線(xiàn)性結(jié)果與NMRI的結(jié)果吻合較好.隨著入射波幅增加到kA=0.08和kA=0.12，船體運(yùn)動(dòng)和自由面變形的非線(xiàn)性現(xiàn)象更加突出，此時(shí)必須采用弱非線(xiàn)性或強(qiáng)非線(xiàn)性理論才能更好地反映和模擬真實(shí)的結(jié)果.從圖10 (a) 和圖11 (a) 可以看出：在kA=0.08和kA=0.12時(shí)，本文的強(qiáng)非線(xiàn)性模型的計(jì)算結(jié)果與非線(xiàn)性理論的NMRI的結(jié)果更趨于一致.而NTUA的結(jié)果與其他兩非線(xiàn)性理論(NMRI和本數(shù)值方法)的計(jì)算結(jié)果相差偏大.3) 無(wú)論是垂蕩還是縱搖運(yùn)動(dòng)，總的來(lái)說(shuō)：本計(jì)算方法所得的運(yùn)動(dòng)響應(yīng)幅值要略小于其他兩種數(shù)值方法.其原因?yàn)椋词乖趉A=0.04時(shí)，本模型的計(jì)算結(jié)果也能看到波浪破碎等非線(xiàn)性現(xiàn)象，因此能量也會(huì)因?yàn)椴ɡ似扑榈榷兴鶞p小，從而造成船體運(yùn)動(dòng)幅值的減小.而基于邊界元方法 (NMRI和NTUA) 是無(wú)法模擬波浪破碎等強(qiáng)非線(xiàn)性現(xiàn)象，也不存在因?yàn)閺?qiáng)非線(xiàn)性因素造成的能量損失.因此船舶運(yùn)動(dòng)響應(yīng)的計(jì)算結(jié)果會(huì)較實(shí)際值偏大.

圖9 船體的垂蕩和縱搖運(yùn)動(dòng) (Fn =0.2,λ/L=1.0, kA=0.04)

圖10 船體的垂蕩和縱搖運(yùn)動(dòng) (Fn=0.2,λ/L=1.0,kA=0.08)

圖11 船體的垂蕩和縱搖運(yùn)動(dòng) (Fn=0.2,λ/L=1.0,kA=0.12)

圖12給出了S175船模在不同入射波長(zhǎng)下，以航速Fn=0.2航行時(shí)的運(yùn)動(dòng)響應(yīng)時(shí)間歷程線(xiàn).眾所周知，當(dāng)船體運(yùn)動(dòng)的固有周期(如：縱搖、橫搖等周期)與波浪周期產(chǎn)生共振時(shí)，船體將會(huì)激起較大的運(yùn)動(dòng)幅值.不同的運(yùn)動(dòng)模態(tài)(如：垂蕩、縱搖、橫搖)，其各自的固有周期是不同的.從圖12可以看出：隨著入射波長(zhǎng)的增長(zhǎng)，船體垂蕩運(yùn)動(dòng)在λ/L=1.2時(shí)產(chǎn)生較大的運(yùn)動(dòng)響應(yīng)振幅，而縱搖運(yùn)動(dòng)則在λ/L=1.4產(chǎn)生較大的運(yùn)動(dòng)響應(yīng).

圖12 3種不同入射波長(zhǎng)下，船體的垂蕩和縱搖運(yùn)動(dòng) (Fn=0.2,kA=0.04)

2.4 船舶強(qiáng)非線(xiàn)性運(yùn)動(dòng)與波形分析

圖13給出了S175船在λ/L=1.2，kA=0.16的波浪中，以航速Fn=0.2迎浪航行時(shí)的4個(gè)不同時(shí)間的運(yùn)動(dòng)截圖.圖14給出了無(wú)因次化的船體運(yùn)動(dòng)及其所受波浪載荷隨時(shí)間的變化歷程；其無(wú)因次化形式與上述結(jié)果一致，即：縱蕩運(yùn)動(dòng)為Surge/A, 垂蕩運(yùn)動(dòng) Heave/A, 縱搖運(yùn)動(dòng)Pitch/A；縱蕩力為Fx/ρAωe2, 垂蕩力Fz/ρAωe2, 縱搖力矩My/ρALωe2. 圖13 (a) 為船體處于接近平衡位置；圖13 (b) 為船首正俯沖下水，激起較大浪花，船首部受到較大的沖擊載荷；圖13 (c) 為船體發(fā)生甲板上浪，甲板受到上浪拍擊的時(shí)刻；圖13 (d) 為船首出水后繼續(xù)前行，水花飛濺和波浪破碎在船首附近發(fā)生.從圖13 (b)～ (d) 可以看出, 該船體運(yùn)動(dòng)過(guò)程是一個(gè)伴隨著艏部沖擊、甲板上浪和水花飛濺等復(fù)雜現(xiàn)象的強(qiáng)非線(xiàn)性過(guò)程.因此，相應(yīng)的船體運(yùn)動(dòng)響應(yīng)和所受的波浪載荷也將是一個(gè)強(qiáng)非線(xiàn)性過(guò)程，這一點(diǎn)可以從圖14中較易看出，尤其是圖14 (b) 中，由于沖擊載荷的存在，波浪力曲線(xiàn)的局部呈現(xiàn)為鋸齒形.

圖13 船舶強(qiáng)非線(xiàn)性運(yùn)動(dòng)一個(gè)周期內(nèi)的4個(gè)不同時(shí)刻的快照(Fn =0.2,λ/L=1.2, kA=0.16)

3 結(jié) 論

1)本文建立的波浪中航行船舶的強(qiáng)非線(xiàn)性運(yùn)動(dòng)響應(yīng)的數(shù)值模型，具有較好的收斂性和精確性.

2)本模型給出的船模垂蕩和縱搖運(yùn)動(dòng)響應(yīng)，與其他數(shù)值計(jì)算結(jié)果吻合較好；由于本模型考慮了粘性和強(qiáng)非線(xiàn)性效應(yīng)，隨著非線(xiàn)性現(xiàn)象的增強(qiáng)，本模型的垂蕩運(yùn)動(dòng)計(jì)算結(jié)果與非線(xiàn)性理論的NMRI的結(jié)果更趨于一致.

3)對(duì)S175船在大幅波浪中航行時(shí)的強(qiáng)非線(xiàn)性運(yùn)動(dòng)響應(yīng)和自由液面大變形(或破波)的成功模擬，表明本模型具有較好地模擬水花飛濺、甲板上浪、底部砰擊等強(qiáng)非線(xiàn)性問(wèn)題.

圖14 船舶強(qiáng)非線(xiàn)性運(yùn)動(dòng)響應(yīng)和船舶所受波浪力的時(shí)間歷程線(xiàn) (Fn=0.2,λ/L=1.2,kA=0.16)

[1] 戴仰山，沈進(jìn)威，宋競(jìng)正. 船舶波浪載荷[M]. 北京: 國(guó)防工業(yè)出版社，2007.

[2] 陳良偉，肖桃云. 艦船波浪中非線(xiàn)性運(yùn)動(dòng)響應(yīng)與砰擊研究[J]. 船海工程, 2011, 40(6): 31-34. DOI: 10.3963/J.ISSN.1671-7953.2011.06.009. CHEN Liangwei， XIAO Taoyun. Study on nonlinear motion response and slam of ship in waves[J]. Ship & Ocean Engineering, 2011, 40(6): 31-34. DOI: 10.3963/J.ISSN.1671-7953.2011.06.009.

[3] KASHIWAGI M. Prediction of surge and its effect on added resistance by means of the enhanced unified theory [J]. Transactions of the West-Japan Society of Naval Architects, 1995, 89: 77-89.

[4] HE Guanghua, KASHIWAGI M. A time-domain higher-order boundary element method for 3D forward-speed radiation and diffraction problems [J]. Journal of Marine Science and Technology, 2014, 19(2): 228-244. DOI: 10.1007/s00773-013-0242-1.

[5] 陳京普, 朱德祥. 船舶在波浪中運(yùn)動(dòng)的非線(xiàn)性時(shí)域數(shù)值模擬[J]. 水動(dòng)力學(xué)研究與進(jìn)展(A輯), 2010, 25(6): 830-836.DOI: 10.3969/j.issn.1000-4874.2010.06.015. CHEN Jingpu, ZHU Dexiang. Numerical simulations of nonlinear ship motions in waves by a Rankine panel method[J]. Journal of Hydrodynamics, 2010, 25(6): 830-836.DOI: 10.3969/j.issn.1000-4874.2010.06.015.

[6] LIN W M, YUE D K P. Numerical solution for large amplitude ship motions in the time domain [C]// Proceedings of 18th Symposium on Naval Hydrodynamics. Washington DC: National Academy Press, 1991: 41-66.

[7] NAKOS D, SCLAVOUNOS P. Ship motions by a three-dimensional Rankine panel method [C]// Proceedings of 18th Symposium on Naval Hydrodynamics. Washington DC: National Academy Press, 1991: 21-39.

[8] KRING D , SCLAVOUNOS P. Numerical stability analysis for time-domain ship motion simulations [J]. Journal of Ship Research, 1995, 39(4): 313-320.

[9] HUANG Yifeng. Nonlinear ship motions by a Rankine panel method [D]. Massachusetts: Massachusetts Institute of Technology, 1997.

[10] KIM Y, KIM K H, KIM J H, et al. Time domain analysis of nonlinear motion responses and structural loads on ships and offshore structures: development of wish programs [J]. International Journal of Naval Architecture and Ocean Engineering, 2011, 3(1): 37-52.DOI: 10.2478/IJNAOE-2013-0044.

[11] 沈志榮, 葉海軒, 萬(wàn)德成. 船舶在迎浪中運(yùn)動(dòng)響應(yīng)和波浪增阻的RANS數(shù)值模擬[J]. 水動(dòng)力學(xué)研究與進(jìn)展(A輯), 2012, 27(6): 621-633.DOI: 10.3969/j.issn1000-4874.2012.06.001. SHEN Zhirong, YE Haixuan, WAN Decheng. Motion response and added resistance of ship in head waves based on RANS simulations[J]. Journal of Hydrodynamics, 2012, 27(6): 621-633.DOI: 10.3969/j.issn1000-4874.2012.06.001.

[12] 石博文, 劉正江, 吳明. 船模不規(guī)則波中頂浪運(yùn)動(dòng)數(shù)值模擬研究[J]. 船舶力學(xué), 2014, 18(8): 906-915. DOI: 10.3969/j.issn.1007-7294.2014.08.005. SHI Bowen, LIU Zhengjiang, WU Ming. Numerical simulation of ship motions in irregular head waves[J]. Journal of Ship Mechanics, 2014, 18(8): 906-915. DOI: 10.3969/j.issn.1007-7294.2014.08.005.

[13] TAKEWAKI H, NISHIGUTI A, YABE T. Cubic interpolated pseudo-particle method (CIP) for solving hyperbolic-type equations [J]. Journal of Computational Physics, 1985, 61(2): 261-268. DOI:10.1016/0021-9991(85)90085-3.

[14] YABE T, AOKI T. A universal solver for hyperbolic equations by cubic-polynomial interpolation I. One-dimensional solver [J]. Computer Physics Communications, 1991, 66(2/3): 219-232. DOI: 10.1016/0010-4655(91)90071-R.

[15] YABE T, ISHIKAWA T, WANG Peiyuan, et al. A universal solver for hyperbolic equations by cubic polynomial interpolation II. Two- and three-dimensional solvers [J]. Computer Physics Communications, 1991, 66(2/3): 233-242. DOI:10.1016/0010-4655(91)90072-S.

[16] HU Changhong, KASHIWAGI M. A CIP-based method for numerical simulations of violent free-surface flows [J]. Journal of Marine Science and Technology, 2004, 9(4): 143-157. DOI: 10.1007/s00773-004-0180-z.

[17] HE Guanghua, ISSHIKI T, KASHIWAGI M, et al. Prediction of Radiation Forces by Means of a CIP-based Cartesian Grid Method [C]// Proceedings of 21st International Offshore and Polar Engineering Conference. Maui, Hawaii: International Society of Offshore and Polar Engineers, 2011, 2: 633-638.

[18] ZHAO Xizeng, HU Changhong. Numerical and experimental study on a 2-D floating body under extreme wave conditions [J]. Applied Ocean Research, 2012, 35: 1-13. DOI: 10.1016/j.apor.2012.01.001.

[19] LIAO Kangping, HU Changhong. A coupled FDMFEM method for free surface flow interaction with thin elastic plate [J]. Journal of Marine Science and Technology, 2013, 18(1): 1-11. DOI: 10.1007/s00773-012-0191-0.

[20] HE Guanghua. A new adaptive Cartesian-grid CIP method for computation of violent free-surface flows [J]. Applied Ocean Research, 2013, 43: 234-243. DOI: 10.1016/j.apor.2013.10.004.

[21] XIAO Feng, HONMA Y, KONO T. A simple algebraic interface capturing scheme using hyperbolic tangent function [J]. International Journal for Numerical Methods in Fluids, 2005, 48(9): 1023-1040. DOI: 10.1002/fld.975.

[22] KIM Y. Comparative study on linear and nonlinear ship motion and loads [C]// Proceedings of ITTC Workshop on Seakeeping. Seoul: Seoul National University, 2010: 283-310.

(編輯 張 紅)

Strongly-nonlinear simulation of ship motions in head waves

HE Guanghua, CHEN Limin, WANG Jiadong

(School of Naval Architecture and Ocean Engineering, Harbin Institute of Technology in Weihai, Weihai 264209, Shandong， China)

To simulate the ship advancing in large waves with large-amplitude motion, a CIP-based strongly-nonlinear seakeeping analysis tool is developed. The effects of viscosity and strongly-nonlinearity are considered in the numerical model, which is based on the viscous flow theory. First, convergence studies are carried out by calculating the radiation problems of S175 hull with forward speed. Then, the ship response of S175 hull including the heave and pitch motions in waves with amplitude and length of waves varied is calculated. Finally, strongly-nonlinear ship responses in large head waves are calculated and captured successfully. It is confirmed that the developed seakeeping analysis tool has numerical stability and accuracy, and has the capability of simulating the strongly-nonlinear seakeeping problem, such as, the phenomena of wave breaking, slamming, and green water on deck.

CIP; strongly-nonlinear; ship motion response; large amplitude motion; time domain simulation

10.11918/j.issn.0367-6234.201508063

2015-08-21

國(guó)家自然科學(xué)基金 (51579058); 山東省自然科學(xué)基金 (ZR2014EEQ016)

何廣華(1980—)，男，教授，博士生導(dǎo)師

何廣華，ghhe@hitwh.edu.cn

U661.31

A

0367-6234(2017)04-0142-07