規(guī)則波中船舶騎浪運動條件及失穩(wěn)特性分析

李焱，王麗元,3，封培元，劉震，劉崢，唐友剛

(1.天津大學(xué) 水利工程仿真與安全國家重點實驗室，天津 300350；2.天津大學(xué) 天津大學(xué)建筑工程學(xué)院，天津 300350；3.江蘇科技大學(xué) 海洋裝備研究院，江蘇 鎮(zhèn)江 212003；4.上海市船舶工程重點實驗室，上海 200011；5.中國船舶及海洋工程設(shè)計研究院，上海 200011)

船舶騎浪運動是船舶第二代完整穩(wěn)性中的5種失效模式之一[1]，其非線性動力學(xué)特性最為顯著。船舶騎浪運動失穩(wěn)模式，是指船舶在特定的波浪中航行時，船舶推進控制系統(tǒng)失去對船的控制，并使船舶隨著波浪前進的現(xiàn)象。即船舶在隨浪條件下航行時會被前進的波浪所捕獲，并與波浪一起運動，使船舶的航速與波速相同。通常使船舶發(fā)生騎浪運動的波長約為船長的1.0～3.0倍，波浪足夠陡峭，且船速與波速相當，船速約為波速的75%[2]。

陶醉等[3-4]建立隨浪條件下船舶的單自由度縱蕩運動方程，通過將船舶騎浪運動與單擺運動同宿分岔現(xiàn)象進行類比，研究船舶縱蕩非線性系統(tǒng)中所存在的同宿分岔現(xiàn)象，并預(yù)報船舶是否發(fā)生騎浪運動。儲紀龍等[5]建立船舶騎浪運動方程，將其簡化為哈密頓系統(tǒng)，求出其對應(yīng)的異宿軌道方程，計算船舶發(fā)生騎浪運動的螺旋槳臨界轉(zhuǎn)速。顧民等[6]以內(nèi)傾船型為例，開展隨浪和尾斜浪中的船舶騎浪/橫甩試驗，分析了內(nèi)傾船型在不同波浪條件下的運動特性。

Kan等[7]通過模型試驗和數(shù)值模擬結(jié)果相對比，給出一種獲得臨界航速和臨界波高的簡便方法，判斷船舶是否發(fā)生騎浪運動：當船體的速度低于第一臨界速度時，船舶在縱蕩方向上作往復(fù)的周期運動，并不會發(fā)生騎浪運動；在第二臨界速度下船舶被波浪捕獲，與波浪之間的相對速度為零，發(fā)生騎浪運動。Spyrou等[8]研究在大波陡隨浪條件下船舶縱蕩運動的強非線性特性，研究表明在船舶騎浪條件的邊緣，船舶會發(fā)生大幅縱蕩運動。

本文基于非線性動力學(xué)理論，研究船舶騎浪運動特性。建立縱浪條件下船舶縱蕩運動方程，分別采用梅林科夫方法和龍格庫塔方法對該方程進行求解。通過解析分析方法，得到發(fā)生騎浪的臨界值，給出梅林科夫函數(shù)閾值圖。在數(shù)值計算中，分別繪制船舶騎浪時和未騎浪時的時歷曲線和運動相圖，并分析其非線性特性。以螺旋槳轉(zhuǎn)速為分岔參數(shù)畫出船舶騎浪運動分岔圖，分析了船舶騎浪的失穩(wěn)條件。

1 縱蕩運動方程建立

船舶在波浪上航行的坐標系如圖1所示。O(x,y,z)為大地固定坐標系，O′(x′,y′,z′)定義為隨船坐標系，其中以船舶的重心G為隨船坐標系坐標原點，O″(x″,y″,z″)為波浪坐標系，其原點固定于波浪的一個波谷處，并以波速c在大地坐標系中移動。

圖1 船舶運動坐標系

船舶在規(guī)則波波浪激勵下縱蕩方向的經(jīng)典的運動方程為[9]：

(1)

式中：m為船舶的質(zhì)量；mx為縱蕩附連水質(zhì)量；ξG為大地坐標系下船舶重心的縱坐標；R(u)為船舶的縱向阻力，其表達式為：

R(u)=r1u+r2u2+r3u3

式中：r1、r2、r3為阻力系數(shù)；u為縱蕩瞬時速度。

T(u,n)為縱蕩方向螺旋槳推力，其表達式為：

式中：tp表示為螺旋槳推力減額系數(shù)；ρ為海水密度；n為螺旋槳轉(zhuǎn)速；Dp為螺旋槳槳葉直徑；KT表示為螺旋槳推力系數(shù)，其表達式為：

KT(u,n)=k1(u,n)+k2(u,n)+k3(u,n)

式中：k1、k2、k3為螺旋槳性能參數(shù)，代表螺旋槳敞水特性曲線的回歸系數(shù)；JP為螺旋槳進速系數(shù)，表達式為：

JP=u(1-ωp)/(nDP)

式中ωp為螺旋槳性能參數(shù)中的伴流系數(shù)。

Xw(ξG,t)為船舶在縱蕩方向上受到的波浪力，按切片理論計算[10]：

(2)

式中；g表示重力加速度；A代表波浪幅值；k為波數(shù)；χ為船舶航向角；S(x)為船舶橫剖面浸水面積；d(x)為各站對應(yīng)吃水；B(x)為各站船寬；dx為站距;α為考慮繞射效應(yīng)的波浪力修正系數(shù)，按照下式進行計算：

式中：Cb為船體方形系數(shù)；Cm為船體中橫剖面系數(shù)。

考慮規(guī)則波條件計算船體受到的波浪力，波浪力隨時間正弦變化。由于各時刻波浪力大小與船-波相對位置有關(guān)，因此將波浪力寫作波浪力幅值乘以船舶與波浪之間相對位置的正弦函數(shù)的形式：

(3)

2 梅林科夫方法

梅林科夫(Melnikov)方法始于20世紀60年代，通過度量龐加萊(Poincare)映射的雙曲不動點的穩(wěn)定與不穩(wěn)定流形之間的距離來確定非線性動力系統(tǒng)是否存在橫截同(異)宿點，從而判斷動力系統(tǒng)是否發(fā)生運動的突變，在確定性系統(tǒng)和隨機系統(tǒng)中都有大量應(yīng)用[11-12]。

首先假設(shè)動力系統(tǒng)的參數(shù)：

x=f(x)+εg(x,t;μ)

(4)

其中：f:R2→R2是系統(tǒng)對應(yīng)的哈密頓矢量場，εg:R2×R×R→R2代表作用在系統(tǒng)上的一個小擾動，μ是系統(tǒng)參數(shù)。此后計算系統(tǒng)的梅林可夫積分函數(shù)：

(5)

式中xh(t)表示動力系統(tǒng)的異宿軌道方程。

在梅林可夫函數(shù)表達式中隱含時間t，可通過求解梅林可夫函數(shù)的零點來找到船舶縱蕩運動系統(tǒng)失穩(wěn)時的參數(shù)閾值，從而得到船舶發(fā)生騎浪失穩(wěn)的條件。

(6)

T(c;n)-R(c)

(7)

其中：

(8)

其中：

引入小參數(shù)ε，其中0<ε?1，則式(8)轉(zhuǎn)化為：

(9)

當小參數(shù)ε=0時，式(9)退化為沒有外激勵作用，沒有阻尼的哈密頓系統(tǒng)：

(10)

(11)

將哈密頓系統(tǒng)式(11)中的右端項取不定積分，并由第1項減去第2項，則得到哈密頓系統(tǒng)(11)對應(yīng)的Hamilton量H(y,z)：

(12)

當Hamilton量H(±π,0)=1時，整個哈密頓系統(tǒng)存在著連接(π,0)和(-π,0)這2個鞍點的2條異宿軌道，在異宿軌道之間，存在著幾條圍繞著(0,0)中心的同宿軌道，從同宿軌道變換成的第1條的異宿軌道則表示此縱蕩運動系統(tǒng)中的船舶是否發(fā)生騎浪失穩(wěn)運動的分界線。其中系統(tǒng)的異宿軌道方程為：

z=±2cos(y/2)

(13)

圖2中，虛線為本系統(tǒng)的異宿軌道線，虛線內(nèi)部的實線則表示系統(tǒng)的同宿軌道方程。

圖2 異宿軌道示意

省略式(9)右端項中的小參數(shù)，并將其沿著異宿軌道方程進行積分，得到縱蕩運動方程的梅林可夫函數(shù)為：

(14)

(15)

改變函數(shù)的積分路線，使梅林可夫函數(shù)沿著式(13)所示的異宿軌道進行積分，并使積分后的結(jié)果M=0，則：

(16)

(17)

對式(17)進行積分，得到I1=4，I2=π，I3=8/3，代入式(16)得：

(18)

由此，根據(jù)式(18)可求出使等式成立的臨界螺旋槳轉(zhuǎn)速n。再根據(jù)船舶在波浪中前行時的阻力和推力相等的關(guān)系式，求出船舶在波浪中行駛時的螺旋槳臨界船速，便確定了船舶騎浪失穩(wěn)的條件。

3 某內(nèi)傾船算例

以某內(nèi)傾船為例，分別采用解析方法和數(shù)值方法計算內(nèi)傾船的縱蕩運動響應(yīng)，并判斷其是否發(fā)生騎浪運動，給出內(nèi)傾船發(fā)生騎浪運動的臨界值條件，并通過2種方法結(jié)果的比較，驗證解析方法的正確性。內(nèi)傾船垂線間長為112 m，吃水為3.865 m，排水體積為2 789.72 m3，水線寬為13.399 m，水線長為111.997 m，濕表面積為1 518.31 m2，浮心縱向位置距舯后-2.456 m，重心垂向位置距基線5.785 m，中橫剖面系數(shù)為0.794 6。內(nèi)傾船的縱蕩運動方程中的推力系數(shù)和阻力系數(shù)由模型試驗的數(shù)據(jù)給出。

3.1 臨界轉(zhuǎn)速與臨界船速計算結(jié)果

選取波浪條件為λ/Lpp=1.5，航向角χ=0°?？紤]波幅為漸變量，采用解析方法求解船舶發(fā)生騎浪運動時螺旋槳臨界轉(zhuǎn)速與波幅之間的關(guān)系，采用梅林科夫方法計算得到的臨界螺旋槳轉(zhuǎn)速閾值曲線如圖3所示。由圖3可以看出當波長保持不變、波幅增加時，梅林科夫方法得出的螺旋槳臨界轉(zhuǎn)速逐漸減低。即在波幅較大時，螺旋槳臨界轉(zhuǎn)速較小，船舶易于發(fā)生騎浪現(xiàn)象。相反，當波幅較小時，螺旋槳臨界轉(zhuǎn)速則更高，騎浪現(xiàn)象更難發(fā)生。

圖3 螺旋槳臨界轉(zhuǎn)速隨波幅變化曲線

根據(jù)船舶的推力和阻力相等條件，可以進一步求出發(fā)生騎浪的臨界船速。船舶臨界船速隨波幅變化的關(guān)系圖如圖4所示。由圖4可以看出，船舶臨界船速隨波幅變化趨勢與船舶螺旋槳臨界轉(zhuǎn)速隨波幅的變化趨勢相同。即當波幅值增加時，梅林科夫方法得出船舶臨界航速逐漸減低。這一結(jié)果說明，在波幅較大時，船舶臨界航速較小，船舶就會易于發(fā)生騎浪現(xiàn)象。相反，當波幅較小時，船舶臨界航速則更高，不易發(fā)生騎浪。

圖4 船舶臨界船速隨波幅變化曲線

3.2 數(shù)值模擬計算結(jié)果

采用龍格庫塔方法對船舶縱蕩運動進行數(shù)值求解，計算不同螺旋槳轉(zhuǎn)速下船舶縱蕩運動，通過對比縱蕩運動的數(shù)值解形式，判斷螺旋槳的臨界轉(zhuǎn)速。

選取波浪條件為：λ/Lpp=1.5，航向角χ=0°。H/λ=1/10。在此波浪條件下，由梅林科夫方法得到的螺旋槳臨界轉(zhuǎn)速為3.38 r/s，對應(yīng)臨界航速為10.38 m/s，采用龍格庫塔方法得出的螺旋槳臨界轉(zhuǎn)速為3.5 r/s，對應(yīng)的臨界船速為10.69 m/s。將2種方法得出的結(jié)果進行對比，得出2種算法預(yù)報的螺旋槳臨界轉(zhuǎn)速誤差為3.43%，臨界船速誤差為2.90%，均在5%以內(nèi)，驗證了解析方法的準確性。

當螺旋槳轉(zhuǎn)速為3.3 r/s時，船舶與波浪間的相對速度如圖5(a)所示，與波浪的相對位移圖5(b)所示。從2個時間歷程曲線圖中可以看出，當螺旋槳轉(zhuǎn)速為3.3 r/s時，船舶的縱蕩速度相對于波浪速度是周期變化的，同時船舶與波浪的相對位移隨時間是逐步增大的。此時船由螺旋槳進行推進，并未被波浪所捕獲，所以在螺旋槳轉(zhuǎn)速為3.3 r/s時，船舶未發(fā)生騎浪現(xiàn)象。

圖5 轉(zhuǎn)速3.3 r/s時船舶縱蕩運動響應(yīng)

當螺旋槳轉(zhuǎn)速增加到3.5 r/s時，船與波的相對速度及相對位移如圖6(a)和圖6(b)所示。由2個時間歷程曲線圖中可以看出，當螺旋槳轉(zhuǎn)速為3.5 r/s時，船舶在波浪中運動約200 s后，船舶與波浪之間的相對速度變?yōu)榱?，同時船舶與波浪之間的相對位移也保持不變。從這2個時歷曲線圖中可以看出，當螺旋槳轉(zhuǎn)速為3.5 r/s時，船舶在運動一段時間后會以波浪速度前進，即船舶被波浪所捕獲。此時船的螺旋槳推力系統(tǒng)失去對船舶的控制，發(fā)生了船舶騎浪運動現(xiàn)象。

圖6 轉(zhuǎn)速3.3 r/s時船舶縱蕩運動響應(yīng)

圖7 船舶運動相圖

3.3 船舶騎浪運動分岔特性

取波浪條件為λ/Lpp=1.5，航向角χ=0°，H/λ=1/10。給定不同的螺旋槳轉(zhuǎn)速n，求船舶縱蕩運動的時間歷程。船舶相對速度隨螺旋槳轉(zhuǎn)速變化的分岔圖如圖8所示。從圖8中看出，此波浪條件下，當螺旋槳轉(zhuǎn)速低于3.4 r/s時，系統(tǒng)發(fā)生的是穩(wěn)定的周期運動，在分岔圖的表現(xiàn)上表示為一個點，即一個圓圈。當螺旋槳轉(zhuǎn)速為3.5 r/s時。分岔圖上不能畫出點，此處用點集區(qū)域示意船舶騎浪的非穩(wěn)定區(qū)域的范圍，不再是穩(wěn)定的往復(fù)周期性運動，此時船舶與波浪之間的相對速度為0。所以此條件下船舶縱蕩運動喪失穩(wěn)定性，即船舶發(fā)生騎浪失穩(wěn)運動。

圖8 船舶相對速度隨螺旋槳轉(zhuǎn)速變化分岔圖

4 結(jié)論

1)基于非線性動力學(xué)方法得到的騎浪條件和騎浪失穩(wěn)域，與數(shù)值方法得到的結(jié)論基本吻合，表明非線性動力學(xué)方法分析船舶騎浪運動的可行性。

2)當船舶在規(guī)則波中未發(fā)生騎浪運動時，船速相對于波速呈周期性變化的規(guī)律；當船舶發(fā)生騎浪運動時，船舶被波浪捕獲，相對速度為零。

3)當波幅較大，螺旋槳臨界轉(zhuǎn)速較小時，船舶容易發(fā)生騎浪運動，并且易發(fā)生騎浪失穩(wěn)導(dǎo)致橫甩。