顧思琪，朱仁慶，吳梓鑫

(江蘇科技大學(xué)船舶與海洋工程學(xué)院，江蘇鎮(zhèn)江212003)

隨著海洋工程向深海發(fā)展，海洋結(jié)構(gòu)物所處的海洋環(huán)境更趨惡劣.在惡劣海浪作用下，可能會(huì)產(chǎn)生嚴(yán)重的波浪爬升、相對(duì)運(yùn)動(dòng)和上浪等流體強(qiáng)非線性沖擊現(xiàn)象，這些都會(huì)對(duì)海洋浮式結(jié)構(gòu)物產(chǎn)生很大的局部沖擊破壞載荷.因此，探索強(qiáng)非線性波浪與海洋浮式結(jié)構(gòu)物的相互作用、發(fā)展穩(wěn)定高效的數(shù)值預(yù)報(bào)方法，不僅對(duì)海洋水動(dòng)力學(xué)的發(fā)展具有重要的理論意義，而且對(duì)海洋工程結(jié)構(gòu)物的設(shè)計(jì)和防護(hù)具有重要的工程價(jià)值.

早期，人們對(duì)波浪與結(jié)構(gòu)物完全非線性相互作用的研究大都以試驗(yàn)為主.近十年來(lái)，經(jīng)過(guò)水動(dòng)力學(xué)研究者們的不懈努力，數(shù)值方法不斷獲得突破，取得了令人鼓舞的成果［1］.目前，主要基于兩種理論模型對(duì)其展開研究，即勢(shì)流模型和粘性流體模型.基于勢(shì)流理論的數(shù)值模擬中，最經(jīng)典的是1976年由 Longuet-Higgins和 Cokelet［2］提出的混合歐拉-拉格朗日方法(MEL)數(shù)值模型，它可用于處理完全非線性波面運(yùn)動(dòng).文獻(xiàn)［3］中基于勢(shì)流理論和高階邊界元模型，建立了無(wú)限水域中波物作用的完全非線性模型.文獻(xiàn)［4-6］中基于完全非線性的勢(shì)流理論(fully nonlinear potential theory，F(xiàn)NPT)，采用QALE-FEM方法模擬了方形數(shù)值水池中水波與坐底波形水壩的流固耦合響應(yīng)，并將其推廣成功實(shí)現(xiàn)了二維陡波與浮體的完全非線性流固耦合模擬.而基于粘流理論的數(shù)值水池，主要結(jié)合自由表面追蹤技術(shù)，求解連續(xù)方程、N-S方程.通過(guò)引入 VOF，MAC，SPH 和 LEVEL -SET 等［7-10］重構(gòu)自由面的數(shù)值方法，不但可以充分考查由于粘性作用導(dǎo)致的各種復(fù)雜流動(dòng)特征，而且可以處理結(jié)構(gòu)物在波浪激勵(lì)下的大幅運(yùn)動(dòng)，以及波浪翻卷、破碎、合并等流體完全非線性運(yùn)動(dòng)等問(wèn)題.因此，基于 N-S方程構(gòu)建數(shù)值波浪水池，研究粘流中波浪與結(jié)構(gòu)物的相互作用等問(wèn)題，已成為當(dāng)今國(guó)際上的一個(gè)重點(diǎn)發(fā)展方向.文獻(xiàn)［11］中利用MAC差分格式，對(duì)歐拉動(dòng)力學(xué)方程進(jìn)行離散求解，利用VOF法處理自由表面的演變，研究了二維完全非線性波在艙內(nèi)的運(yùn)動(dòng)特性.文獻(xiàn)［12］中用VOF法對(duì)FPSO的上浪問(wèn)題作了研究，其計(jì)算對(duì)固定的FPSO模型可以得到較好結(jié)果，但對(duì)運(yùn)動(dòng)的船和波浪之間的相互作用問(wèn)題則尚不能給出理想結(jié)果.文獻(xiàn)［13-14］中采用區(qū)域分解思想、遠(yuǎn)處流場(chǎng)和船體的運(yùn)動(dòng)通過(guò)勢(shì)流理論來(lái)計(jì)算，船體附近的流動(dòng)通過(guò)求解N-S方程來(lái)模擬，這種分開求解的方式具有一定的局限性，與真實(shí)情況中波浪與船體的相互作用有較大差別.

文中以ANSYS-14.5軟件的Workbench為計(jì)算平臺(tái)，結(jié)合二次開發(fā)功能，建立了具有造波和消波功能的二維數(shù)值波浪水池.文中就干舷高度不同的結(jié)構(gòu)物遭遇波浪的完全非線性現(xiàn)象以及波浪對(duì)結(jié)構(gòu)物運(yùn)動(dòng)響應(yīng)和上浪的影響分別進(jìn)行了數(shù)值模擬研究，并采用VOF法追蹤自由液面，精確地描述了波浪與結(jié)構(gòu)物相互作用問(wèn)題，得到了結(jié)構(gòu)物遭遇波浪后甲板上的上浪水位高度和結(jié)構(gòu)物的縱搖、垂蕩、橫搖等運(yùn)動(dòng)響應(yīng).

1 物理模型

1.1 流體運(yùn)動(dòng)控制方程

控制方程包括連續(xù)方程和不可壓縮粘性流體的N-S方程.

連續(xù)方程

對(duì)于不可壓縮流體，方程簡(jiǎn)化為

式中:v為流體速度，ρ為流體密度.

動(dòng)量方程

式中:ν為流體的運(yùn)動(dòng)粘性系數(shù)，F(xiàn)b為作用于流體的質(zhì)量力，p為流場(chǎng)壓強(qiáng).

1.2 自由液面追蹤

文中采用VOF法對(duì)自由液面進(jìn)行求解，由此確定自由液面的位置.自由液面的標(biāo)記采用流體體積函數(shù)法，體積分?jǐn)?shù)αq定義為單元內(nèi)第q相流體所占體積與該單元的體積比.若αq=1則表示該單元內(nèi)全部為第q相流體，若αq=0則表示該單元內(nèi)沒(méi)有第q相流體，若0＜αq＜1則表示該單元為交界面單元.αq滿足方程

1.3 波浪理論

對(duì)于stokes三階波，其波面方程和水質(zhì)點(diǎn)的速度分別為［15］:

波面方程

x方向速度分量

z方向速度分量

2 數(shù)值模型

2.1計(jì)算模型

文中所采用的二維數(shù)值波浪水池模型見(jiàn)圖1.波浪水池長(zhǎng)15m，高2 m，其中水深1.5 m，水池末端5m長(zhǎng)的區(qū)域?yàn)橄▍^(qū).坐標(biāo)系原點(diǎn)位于入射邊界的靜水面處，ox軸上方全部定義為空氣，密度1.225 kg/m3，ox軸下方全部定義為水，密度998.2kg/m3.

圖1 數(shù)值水池示意Fig.1 Schematic diagram of the numerical wave tank

圖1中將浮式結(jié)構(gòu)物簡(jiǎn)化成一方形結(jié)構(gòu)，長(zhǎng)1 m，吃水0.18 m，干舷0.04 m，船體中橫剖面位于x=5.5 m處.在數(shù)值計(jì)算中設(shè)置水位監(jiān)測(cè)點(diǎn)P1、P2、P3，分別位于x=5.0 m，5.1 m和5.2 m處.

流場(chǎng)網(wǎng)格和結(jié)構(gòu)網(wǎng)格全部采用四邊形結(jié)構(gòu)化網(wǎng)格劃分.網(wǎng)格的質(zhì)量直接影響計(jì)算結(jié)果的精度，文中經(jīng)過(guò)多次試算，確定了較合理的網(wǎng)格劃分方案:沿水池長(zhǎng)度方向網(wǎng)格尺寸為dx=λ/100，沿水池高度方向網(wǎng)格尺寸為dz=H/8，在自由液面處對(duì)網(wǎng)格進(jìn)行加密以提高計(jì)算精度，其高度為dz=H/30.為了適當(dāng)?shù)販p少網(wǎng)格以減少計(jì)算工作量，在遠(yuǎn)離自由液面處網(wǎng)格以一定比例向上下逐漸稀疏.

2.2 邊界條件及參數(shù)設(shè)置

水池頂部邊界設(shè)置為壓力入口條件(Pressure-inlet)，底部邊界為無(wú)滑移固壁條件(Wall)，右側(cè)邊界設(shè)置為壓力出口條件(Pressure-outlet)，左側(cè)邊界定義為速度入口條件(Velocity-inlet).文中基于邊界造波法模擬stokes三階波，由于軟件Fluent的模塊中沒(méi)有提供完整的波浪模擬環(huán)境，因此需要通過(guò)用戶自定義模塊(UDF)，將UDF文件嵌入共享庫(kù)中并與Fluent連接.

同時(shí)，在Fluent中采用非定常分離隱式求解器求解，壓力方程選用加權(quán)體積力格式(Body Force Weighted)，壓力速度耦合方式采用PISO算法，體積分?jǐn)?shù)法為幾何重構(gòu)(Geo-Reconstruct)，壓力參考值為一個(gè)標(biāo)準(zhǔn)大氣壓，重力加速度為9.8066 m/s2，時(shí)間步長(zhǎng)為0.005s.

2.3 消波設(shè)置

為防止水池末端波浪反射而干擾入射波場(chǎng)，在水池末端應(yīng)設(shè)置消波裝置.文中采用阻尼消波方法，即在特定的消波區(qū)域內(nèi)添加阻尼項(xiàng)以達(dá)到減弱或消除該區(qū)域內(nèi)的波動(dòng).具體實(shí)現(xiàn)方式為在水池后端設(shè)置阻尼消波區(qū)，通過(guò)在其動(dòng)量方程中添加自編UDF程序 DEFINE_SOURCE(momentum，c，t，dS，eqn)實(shí)現(xiàn)消波.文中采用設(shè)置速度入口造波，為了保證流動(dòng)的連續(xù)性，維持質(zhì)量守恒，在消波區(qū)內(nèi)只需對(duì)動(dòng)量方程中 z方向的分量加載 DEFINE_SOURCE((z_momentum，c，t，dS，eqn)，x方向不作處理.

消波區(qū)內(nèi)，動(dòng)量方程如式(9).

式中:ν為運(yùn)動(dòng)粘性系數(shù);C(x)為消波系數(shù)，其函數(shù)表達(dá)式為:

式中:L為波長(zhǎng);x0，xL分別為消波區(qū)左、右邊界的x坐標(biāo)值;α為阻尼經(jīng)驗(yàn)參數(shù);ρ為流體密度.

3 結(jié)果與分析

3.1 波浪數(shù)值水池

穩(wěn)定的波浪是波物相互作用數(shù)值模擬的前提，為了驗(yàn)證波浪數(shù)值水池的正確性，文中首先對(duì)波浪周期T=1.13s，理論波長(zhǎng)λ=2 m，理論波高H=0.16 m的三階stokes波進(jìn)行了數(shù)值模擬.圖2給出了數(shù)值波浪水池在x=1 m和2 m處波面(A)隨時(shí)間的變化歷程曲線.在波浪的傳播過(guò)程中，隨著傳播距離的增加，波高是有衰減的，這是由于流體粘性作用的影響導(dǎo)致波浪能量耗散的緣故，因此，圖2中x=2 m處波幅計(jì)算值與理論值偏差較x=1 m處大.總體而言，數(shù)值模擬的波形比較穩(wěn)定，波浪周期、波幅參數(shù)等與理論解吻合良好，波形具有“上尖下平”的非線性特征.

圖3為不同時(shí)刻波浪水池消波區(qū)內(nèi)的波形圖，圖中較為清晰地顯示了波浪的傳播在受到消波區(qū)內(nèi)附加阻尼的影響后波形逐漸衰減的過(guò)程.從圖中可觀察到，在接近水池出口處，消波區(qū)自由液面高度值已趨于0 m，從而驗(yàn)證了文中阻尼消波方法的可行性.

圖2 數(shù)值波形與理論解對(duì)比Fig.2 Comparison of numerical and theoretical wave forms

圖3 不同時(shí)刻消波區(qū)內(nèi)的波形Fig.3 Waves’shape in the damping domain at different time

3.2 波浪與浮體的相互作用

3.2.1 二維甲板上浪驗(yàn)證實(shí)例

為驗(yàn)證文中計(jì)算方法的可行性和有效性，下面將基于ANSYS Workbench計(jì)算平臺(tái)，對(duì)迎浪狀態(tài)下二維固定結(jié)構(gòu)物甲板上浪過(guò)程進(jìn)行數(shù)值模擬，并將計(jì)算結(jié)果與有關(guān)文獻(xiàn)中的試驗(yàn)值進(jìn)行比較，以驗(yàn)證數(shù)值方法的有效性.波物相互作用問(wèn)題是雙向流固耦合問(wèn)題，文中基于軟件ANSYS，采用迭代耦合的方法求解流固耦合問(wèn)題，其主要思路是流體方程和結(jié)構(gòu)方程按順序相互迭代求解，各自在每一步得到的結(jié)果提供給另一部分使用，直到耦合系統(tǒng)的解達(dá)到收斂，迭代停止［16］.

1)計(jì)算模型說(shuō)明

數(shù)值波浪水池和結(jié)構(gòu)模型如圖4，模型尺寸與文獻(xiàn)［16］中甲板上浪試驗(yàn)?zāi)Ｐ捅３忠恢?圖中水池長(zhǎng)13.5 m，其中2 m為消波區(qū)，水深1.035 m.結(jié)構(gòu)模型已簡(jiǎn)化為方形結(jié)構(gòu)物，底部為半徑0.08 m的圓弧，長(zhǎng)和型深分別為1.5 m和0.248 m，吃水為0.198 m，干舷為0.05 m.

圖4 計(jì)算域示意圖Fig.4 Diagram of computational domain

2)計(jì)算結(jié)果

在甲板上浪的數(shù)值模擬中，波浪要素與文獻(xiàn)［17］中試驗(yàn)條件保持一致，波長(zhǎng)和波高分別為2 m和0.16 m.為了與試驗(yàn)結(jié)果對(duì)比，在模擬甲板上浪過(guò)程中對(duì)結(jié)構(gòu)上WL1和WL2兩點(diǎn)的波形時(shí)歷進(jìn)行了監(jiān)測(cè)，其中點(diǎn)WL1位于甲板迎浪側(cè)的首端點(diǎn)，點(diǎn)WL2與WL1的距離為0.075 m.圖5為文中數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果的對(duì)比圖，由圖可知，數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果基本吻合，具體表現(xiàn)在結(jié)構(gòu)甲板上浪時(shí)間以及甲板監(jiān)測(cè)點(diǎn)上浪水位高度(Hs)均與試驗(yàn)結(jié)果基本吻合，從而驗(yàn)證了文中計(jì)算方法的可靠性.如圖5所示，第一個(gè)波浪約在7s時(shí)到達(dá)WL1處，第二個(gè)波浪約在8s到達(dá)，且在第一個(gè)波浪涌上甲板時(shí)，部分水體被船首阻擋反射，使得船首附近的自由液面被抬高.因此，在第二個(gè)波到達(dá)船首時(shí)，波高增大許多，超過(guò)第一個(gè)波高的50%，引起的上浪現(xiàn)象也更為嚴(yán)重.

圖5 FPSO模型監(jiān)測(cè)點(diǎn)處波形時(shí)歷Fig.5 Time-history of wave at probe WL1 and WL2

從數(shù)值模擬結(jié)果與試驗(yàn)值的對(duì)比可以看出，文中的計(jì)算方法是有效可行的，計(jì)算結(jié)果誤差較小，滿足計(jì)算要求，所以該方法適用于下文數(shù)值計(jì)算.

3.2.2 不同波浪要素下浮體與波浪的相互作用

具有波高、波長(zhǎng)等波浪要素不同的規(guī)則波，對(duì)浮體的作用力也不盡相同.在與波浪接觸時(shí)，結(jié)構(gòu)物不僅會(huì)受到波浪的沖擊載荷，而且會(huì)受到上浪水體的影響.除此之外，由于波浪的形態(tài)變化，結(jié)構(gòu)物水下體積也隨之發(fā)生改變，浮心位置也會(huì)發(fā)生變化，因此結(jié)構(gòu)物在波浪中的運(yùn)動(dòng)響應(yīng)是由多種因素共同作用產(chǎn)生的.文中針對(duì)二維問(wèn)題，只考慮浮體的橫搖、垂蕩和橫蕩運(yùn)動(dòng)，研究浮體中橫剖面模型在波浪作用下的載荷與運(yùn)動(dòng)響應(yīng).

圖6為干舷高度為0.04 m時(shí)，浮式結(jié)構(gòu)物在波高H分別為0.12，0.16 m的三階stokes波作用下的運(yùn)動(dòng)響應(yīng)時(shí)歷曲線和上浪水位高度值.從圖中可看出，在不同的波高下，結(jié)構(gòu)物的垂蕩與橫搖曲線呈現(xiàn)一定的周期性，且其周期與波浪周期大致相同.第一個(gè)波浪到達(dá)甲板并發(fā)生上浪現(xiàn)象的時(shí)間大致為6s，在經(jīng)歷了第一個(gè)上浪過(guò)程后，水體的質(zhì)量改變了結(jié)構(gòu)物的重心，結(jié)構(gòu)物的垂蕩時(shí)歷曲線有整體下移趨勢(shì).就上浪水位而言，對(duì)于干舷相同的結(jié)構(gòu)物，波高大的波浪更容易造成上浪現(xiàn)象.

圖6 兩種波高的波浪作用下浮體的運(yùn)動(dòng)響應(yīng)與上浪水位Fig.6 Comparison between the data of floating body with different wave heights

圖7是干舷高度為0.04m時(shí)，λ分別為2m和2.5 m的三階stokes波作用下的運(yùn)動(dòng)響應(yīng)時(shí)歷曲線和上浪水位高度(圖中變量含義見(jiàn)上文).由于不同波長(zhǎng)的影響，浮體的橫蕩值也會(huì)不同，波長(zhǎng)較大的波浪作用下浮體的橫蕩值也隨之增大.波長(zhǎng)分別為2.0 m和2.5 m的波浪，其周期分別為1.13s和1.27s，在這兩種波長(zhǎng)作用下的結(jié)構(gòu)物垂蕩和橫搖時(shí)歷曲線的周期也不盡相同.觀察上浪水位高度圖可以發(fā)現(xiàn)，波長(zhǎng)小的波浪上浪水位較大，即上浪量也隨著波陡的增加而增加.這是由于在波長(zhǎng)較大的波浪作用下結(jié)構(gòu)物的橫搖幅值也較大，如圖7c)所示，因此在隨后的上浪過(guò)程中，上浪的水體較容易隨著結(jié)構(gòu)物的擺動(dòng)流入水池中，產(chǎn)生上浪水位高度較小的現(xiàn)象.

圖7 兩種波長(zhǎng)的波浪作用下浮體的運(yùn)動(dòng)響應(yīng)與上浪水位Fig.7 Comparison between the data of floating body with different wavelengths

3.2.3 不同干舷高度的浮體與波浪的相互作用

足夠的干舷不僅可以保證船舶有一定的儲(chǔ)備浮力，也可以減少船舶上浪.干舷高度的不同，受到波浪沖擊時(shí)結(jié)構(gòu)物的受力也不相同.圖8為在相同波浪下，干舷分別為0.02 m和0.04 m的浮式結(jié)構(gòu)物運(yùn)動(dòng)時(shí)歷曲線和上浪水位高度圖(圖中變量含義見(jiàn)上文).與圖6，7對(duì)比可發(fā)現(xiàn)，相對(duì)于波浪要素的影響來(lái)說(shuō)，浮體干舷對(duì)結(jié)構(gòu)物的橫蕩值影響更大些，干舷較大的結(jié)構(gòu)物的漂移距離要遠(yuǎn)遠(yuǎn)大于干舷較小的結(jié)構(gòu)物.而對(duì)結(jié)構(gòu)物垂蕩和橫搖運(yùn)動(dòng)來(lái)說(shuō)，干舷較大的結(jié)構(gòu)物其運(yùn)動(dòng)幅值也會(huì)較大，但差值不會(huì)特別明顯.在結(jié)構(gòu)物發(fā)生第一個(gè)上浪過(guò)程之前，干舷高度不同的兩個(gè)結(jié)構(gòu)物的垂蕩和橫搖曲線幾乎重合，在結(jié)構(gòu)物發(fā)生上浪后，兩者才產(chǎn)生了變化.原因在于，在相同的波浪下，干舷小的結(jié)構(gòu)物更容易上浪，沖上甲板的水體不會(huì)立刻流回水域，而會(huì)沿著甲板流動(dòng)并持續(xù)一段時(shí)間，此時(shí)這部分水體也會(huì)影響結(jié)構(gòu)物的運(yùn)動(dòng)響應(yīng).

圖8 兩種干舷高度浮體的運(yùn)動(dòng)響應(yīng)與上浪水位Fig.8 Comparison between the data of floating body with different freeboards

圖9為波浪與浮體相互作用時(shí)的流體形態(tài)和結(jié)構(gòu)應(yīng)力云圖.從圖中可以觀察到，上浪現(xiàn)象是波浪和船首下沉運(yùn)動(dòng)兩種機(jī)制組合作用的結(jié)果.其過(guò)程是:水體沿船首升高，波浪高于船首涌到甲板上，水體沿甲板流動(dòng)，從甲板另一側(cè)流出，主要對(duì)背浪側(cè)的流體環(huán)境產(chǎn)生擾動(dòng)，嚴(yán)重時(shí)還伴有氣泡和漩渦產(chǎn)生.在文中的計(jì)算模型中，浮體材料為結(jié)構(gòu)鋼，彈性模量E=2.1×1011Pa，因此結(jié)構(gòu)幾乎沒(méi)有變形產(chǎn)生.但是在t=7.6s左右時(shí)，上浪部分水體流經(jīng)甲板時(shí)，結(jié)構(gòu)的底部與甲板的中部所受應(yīng)力較大，應(yīng)引起注意采取結(jié)構(gòu)的加強(qiáng)措施.

4 結(jié)論

基于ANSYS Workbench平臺(tái)的流固耦合模塊，利用UDF程序進(jìn)行二次開發(fā)實(shí)現(xiàn)速度入口造波，采用VOF方法追蹤氣液兩項(xiàng)自由面，在水池后方采用動(dòng)量方程中添加源項(xiàng)形成阻尼消波，成功建立二維數(shù)值波浪水池，實(shí)現(xiàn)了波浪與結(jié)構(gòu)物相互作用的完全非線性數(shù)值模擬.

1)將設(shè)置造波邊界法應(yīng)用于Fluent模塊，成功模擬了效果良好的三階stokes波.消波區(qū)域采用阻尼消波，避免了波浪到達(dá)右端固壁邊界后產(chǎn)生反射波.

圖9 波浪與自由運(yùn)動(dòng)浮體的相互作用Fig.9 Interaction between wave and freely floating structures

2)采用雙向流固耦合的方法，分別對(duì)不同干舷高度的結(jié)構(gòu)物在不同波浪要素的波浪作用下的載荷與響應(yīng)進(jìn)行了對(duì)比研究，得到了結(jié)構(gòu)物遭遇波浪后的運(yùn)動(dòng)響應(yīng)和監(jiān)測(cè)點(diǎn)的上浪水位高度，發(fā)現(xiàn)結(jié)構(gòu)物的運(yùn)動(dòng)響應(yīng)幅值和甲板上浪水位高度都會(huì)隨著波陡的增加而增加.同時(shí)，增加結(jié)構(gòu)物干舷的高度能在一定程度上減少甲板上浪，但會(huì)使其運(yùn)動(dòng)響應(yīng)幅值尤其是橫蕩幅值增加，在結(jié)構(gòu)設(shè)計(jì)時(shí)應(yīng)加以權(quán)衡.

3)展示了波浪與浮體相互作用時(shí)的流體形態(tài)和結(jié)構(gòu)應(yīng)力云圖，觀察到了波浪與結(jié)構(gòu)物相遇后波浪沿船首爬升、破碎、飛濺等完全非線性現(xiàn)象.并注意到上浪部分水體流經(jīng)結(jié)構(gòu)物甲板時(shí)，結(jié)構(gòu)底部與甲板的中部所受應(yīng)力較大，應(yīng)注意采取結(jié)構(gòu)的加強(qiáng)措施.

References)

［1］ Ma Qingwei.Advances in numerical simulation of nonlinear water waves［M］.Singapore:World Scientific Publishing Co Pte Ltd，2010:1 -39.

［2］ Longuet·Higgins M S，Cokelet E D.The deformation of steep surface waves on water:a numerical method of computation［J］.Proceedings of the Royal Society of London，Series A 350，1976(1660):1 -26.

［3］ Zhou B Z，Wu G X，Teng B.Fully nonlinear wave interaction with freely floating non-wall-sided structures［J］.Engineering Analysis with Boundary Elements，2015(50):117-132.

［4］ Ma Q W，Yan S.Quasi ALE finite element method for nonlinear water waves［J］.Journal of Computational Physics，2006，212:52 -72.

［5］ Yan S，Ma Q W.Numerical simulation of fully nonlinear interaction between steep waves and 2D floating bodies using the QALE-FEM method［J］.Journal of Computational Physics，2007，221(2):666 -692.

［6］ Yan S，Ma Q W.Numerical simulation of interaction between wind and 2D freak waves［J］.European Journal of Mechanics B/Fluids，2010，29(1):18 -31.

［7］ Hirt C W，Nichols B D.Volume of fluid(VOF)method for the dynamics of free boundaries［J］.Journal of Computational Physics，1981，39(1):201 -225.

［8 ］ Stanley Osher，James A Sethian.Fronts propagating with curvature-dependent speed:Algorithms based on Hamilton-Jacobi formulation ［J］.Journal of Computational Physics，1988(01):12 -49.

［9］ Monaghan J J.Smoothed particle hydrodynamics［J］.Reports on Progress in Physics，2005，68(8):543 -574.

［10］ Yabe T，Xiao F，Utsumi T.The constrained interpolation profile method for multiphase analysis［J］.Journal of Computational Physics，2001，69(2):556 -593.

［11］ 朱仁慶，王志東，楊松林.完全非線性波的數(shù)值模擬［J］.船舶力學(xué)，2002，6(5):14-18.Zhu Renqing，Wang Zhidong，Yang Songlin.Fully nonlinear wave simulations［J］.Journal of Ship Mechanics，2002，6(5):14 -18.(in Chinese)

［12］ Nielsen K B，Mayer S.Numerical prediction of green water incidents［J］.Ocean Engineering，2004，31(3 -4):363-399.

［13］ Kleefsman K M，Loots G E.The numerical simulation of green water loading including vessel motions and the incoming wave field［C］∥Proceedings of the 24th International Conference on Offshore Mechanics and Arctic Engineering.Halkidiki，Greece:OMAE，2005:981 -992.

［14］ 朱仁傳，繆國(guó)平，林兆偉.運(yùn)動(dòng)船體甲板上浪的三維數(shù)值模擬［J］.水動(dòng)力學(xué)研究與進(jìn)展A輯，2008，23(1):7-14.Zhu Renchuan，Miao Guoping，Lin Zhaowei.3-D numerical simulation of green water occurrence on oscillating shuip［J］.Journal of Hydrodinamics Ser A，2008，23(1):7-14.(in Chinese)

［15］ 陶建華.水波的數(shù)值模擬［M］.天津:天津大學(xué)出版社，2005:17-252.

［16］ 朱仁慶，李辰，顧思琪，等.彈性液艙內(nèi)液體晃蕩研究［J］.江蘇科技大學(xué)學(xué)報(bào):自然科學(xué)版，2013，27(3):214-218.Zhu Renqing，Li Chen，Gu Siqi，et al.A study on liquid sloshing in elastic tanks［J］.Journal of Jiangsu U-niversity of Science and Technology:Natural Science E-dition，2013，27(3):214 -218.(in Chinese)

［17］ Greco M，F(xiàn)altinsen O M，Landrini M.Green water loading on a deck structure［C］∥Proceedings of 16th International Workshop on Water Waves and Floating Bodies.Japan:IWWWFB，2001.