GZ曲線分段線性化對船舶橫搖運動響應的影響研究

柴 威 楊龍霞 章澤文 吳述慶

(武漢理工大學交通學院1) 武漢 430063) (中國船舶重工集團公司第七〇四研究所2) 上海 200031)(中國船級社浙江分社3) 寧波 315100) (武漢第二船舶設計研究所4) 武漢 430064)

0 引 言

國際海事組織(IMO)正在制定和完善第二代船舶穩性規范.在該規范中，癱船穩性、參數橫搖、純穩性喪失、騎浪/橫甩和過度加速度這五種波浪中的穩性失效模式被列為船舶設計中需要核準的危險狀態[1].癱船模型即失去動力的船舶，在橫向風浪作用下產生大幅橫搖運動的情況.對于癱船模型，在忽略各運動之間的耦合作用后，船舶橫搖運動可以由單自由度的運動方程進行描述[2].在橫搖運動微分方程中，船舶大幅橫搖運動的非線性特性以及外載荷與橫搖運動響應的隨機性可以得到充分體現.

對于單自由度的橫搖運動方程，國內外學者的研究主要集中在非線性動力學領域、隨機動力學領域以及船舶穩性領域.研究船舶橫搖運動非線性動力學特性的方法和技術手段主要有Melnikov方程[3]、Lypnov指數[4]和安全池理論等[5].在隨機動力學領域，當橫搖運動的激勵為白噪聲或過濾白噪聲過程時，船舶橫搖運動響應是一個馬爾科夫過程，其概率特性由相應的Fokker-Planck方程控制.Fokker-Planck方程的解析解僅限于線性系統和極少數非線性隨機系統，常用的數值或近似求解方法主要有路徑積分法[6-7]、隨機平均法等[8].

在船舶穩性領域，國內外學者提出和采用了分段線性化方法對橫搖運動微分方程進行簡化與求解[9-13].該方法將橫搖非線性阻尼進行隨機線性化，同時對GZ曲線進行分段線性化近似，從而得到了橫搖運動微分方程的解析解和船舶傾覆的概率.然而，對于該方法的可靠性目前還缺乏有效的評估，同時GZ曲線分段線性化近似對船舶橫搖運動隨機響應的影響目前也缺乏系統的研究.

文中基于路徑積分法，通過求解橫搖運動微分方程對應的Fokker-Planck方程，獲得了船舶橫搖運動的隨機響應，根據路徑積分法計算得到的橫搖運動平均穿越率(mean upcrossing rate)和譜分析方法，對三種GZ曲線分段線性化近似方案進行了評估，系統地研究了分段線性化方法對隨機橫搖運動響應的影響.

1 數學模型

1.1 癱船模型橫搖運動微分方程

隨機橫向風浪激勵下，單自由度船舶大幅橫搖運動微分方程為

ΔGZ(θ(t))=Mwind(t)+Mwave(t)

(1)

式中：θ(t)為船舶橫搖角；I44為橫搖慣性矩；A44為橫搖附加慣性矩；B44為橫搖線性阻尼系數；B44c為橫搖非線性阻尼系數；Δ為排水量；GZ為復原力臂曲線；Mwind和Mwave分別為風和浪的激勵.

文中主要研究目標為GZ曲線分段線性化近似對隨機橫搖運動響應的影響，因此暫不考慮外部激勵的有色性質.將方程式(1)右邊的隨機激勵簡化為白噪聲激勵，在不影響結論準確性的基礎上，可極大地簡化計算工作量.將簡化后的方程式(1)兩側除以I44+A44，橫搖運動微分方程為

(2)

式中：b44和b44c為相對阻尼系數；N(t)為標準白噪聲過程；σ為隨機激勵系數.采用曲線擬合船舶非線性復原力臂，式(2)中的復原力項為

(3)

1.2 分段線性化方法

GZ曲線的分段線性化目前沒有統一的標準，其簡單示意圖可見圖1.分段線性三角形的高一般與GZ曲線的最大值保持一致，三角形底長度通常取為θv(穩性消失角)，分段線性三角形的面積與GZ曲線圍成的面積近似相等.在引入分段線性化后，式(2)中的相對復原力項c(θ)為

圖1 GZ曲線分段線性化示意圖

(4)

式中：θm為分段線性三角形的頂點所對應的橫坐標值；k1和k2的值可根據頂點坐標以及θm和θv的值求得.

計及GZ曲線分段線性化近似后，船舶橫搖運動方程為

(5)

運用隨機線性化方法對式(5)中的非線性阻尼項進行處理后，可得到該運動方程的解析解，從而可以計算大幅橫搖運動引起的船舶傾覆概率.需指出的是，本文僅研究GZ曲線分段線性化對船舶橫搖運動響應的影響，在后續研究中將保留橫搖運動方程中的非線性阻尼項.

2 路徑積分法

將隨機微分方程式(2)轉化為如下二維狀態方程組.

將式(6)寫成如下伊藤隨機微分方程.

dx=m(x,t)dt+q(t)dt

(7)

當動態系統的外部激勵為白噪聲或過濾白噪聲過程時，橫搖運動響應x(t)=[x1(t),x2(t)]T是一個馬爾科夫過程.對于馬爾可夫過程，其轉移概率密度p(x,t|x′,t′)滿足如下Fokker-Planck方程，為

(8)

式中：x為時間t′對應的狀態向量，且t′

p(x,t|x′,t)=

(9)

(10)

根據動態系統的馬爾科夫特性，當時間t′的概率密度p(x′,t′)已知時，動態系統在時間t的概率密度為

路徑積分法的基本思想是在時間和空間分別對狀態方程式(7)進行離散，以路徑代替數值積分，即利用狀態矢量x的馬爾科夫特性，通過連接短時概率密度函數p(x,t|x′,t′)形成全局轉移概率密度函數，進而得到橫搖運動響應的聯合概率密度函數.因此，當狀態矢量x的初始概率分布p(x(0),t0)已知時，根據式(9)和(10)確定的短時轉移概率密度，可由如下迭代算法計算動態系統在時刻t的概率密度分布，為

p(x,t)=

p(x(2),t2|x(1),t1)dx(1)…dx(n-1)

(12)

式中：x=x(n)=x(tn)；t=tn=t0+(n-1)·τ；x(0)=x(t0)；x(1)=x(t1)；p(x(1),t1|x(0),t0)和p(x(2),t2|x(1),t1)等短時轉移概率密度函數僅與時間間隔τ有關，其值均為p(x,t|x′,t′).

假定橫搖運動是平穩隨機過程，則橫搖響應的概率密度函數為

(13)

另外，在運用路徑積分法時，初始概率可設定為二維高斯概率密度，相應的均值取值為零，方差則通過簡單的數值積分求解微分方程組式(6)進行估算.經過GZ曲線分段線性化近似后的橫搖運動微分方程式(5)對應的橫搖運動響應概率密度函數同樣可通過路徑積分法進行求解，在相關計算中僅需對漂移系數做出對應的修正.

3 計算與分析

3.1 CEHIPAR2792樣船

圖2 CEHIPAR2792樣船型線圖

圖3 GZ曲線及三種線性近似化方案

3.2 白噪聲激勵下船舶橫搖運動響應

選取隨機激勵系數σ=0.007，運用路徑積分法計算得到了船舶橫搖運動響應概率密度函數的穩態解.在白噪聲激勵下，橫搖角和橫搖角速度的聯合概率分布見圖4.將聯合概率密度函數對橫搖角速度進行積分，可以得到橫搖角的概率密度函數.為評估路徑積分法的計算精度與效率，采用蒙特卡洛仿真法(Monte Carlo simulation,MCS)計算船舶橫搖角的概率分布對上述計算結果進行驗證.蒙特卡洛仿真法通過數值求解橫搖運動微分方程(2)獲得橫搖運動響應時歷，進而對橫搖運動響應的概率特性進行統計分析.以上兩種方法獲得的橫搖角概率密度函數見圖5.

圖4 白噪聲激勵下船舶橫搖運動響應聯合概率密度函數

圖5 白噪聲激勵下船舶橫搖角概率密度函數

針對數值計算所得橫搖運動響應時歷，引入高斯分布(Gaussian distribution)對橫搖角的概率密度函數進行擬合.由圖5可知：當橫搖運動幅度較小時，橫搖角的概率密度函數服從高斯分布；當橫搖角絕對值超過25°后，橫搖運動響應呈現出非常明顯的非線性特性，其概率密度函數不再服從高斯分布.路徑積分法的優勢在于其計算精度高，且計算耗時較短，該方法對船舶橫搖運動響應的概率特性，尤其是對應較低概率水平的船舶橫搖響應極值都能準確預報.蒙特卡洛仿真法預測對應較低概率水平的橫搖運動響應時需大量計算樣本，使得計算效率大為降低.此外，當概率密度函數低于10-5時基于蒙特卡洛仿真法的橫搖運動響應預報將遭受很大程度的不確定性.因此，通過上述驗證和對比分析可知，路徑積分法可以作為有效的方法和技術手段對隨機橫搖運動響應進行預報.

3.3 GZ曲線分段線性方案

為系統地研究GZ曲線分段線性化對橫搖運動響應造成的影響，見圖3，本文提出了三種線性近似化方案.

方案2分段線性三角形的底和高與方案1相同，為確保分段線性化后船舶的初穩性高與樣船一致，k1的值等于原始GZ曲線在坐標原點處切線的斜率.對于CEHIPAR 2792樣船，該方案仍低估了GZ曲線所圍成的面積.

方案3既保障GZ曲線分段線性化近似后船舶的初穩性高不變，也保證了分段線性三角形的面積與GZ曲線所圍成的面積相等，該方案在文獻[10]中有詳細描述.

3.4 GZ曲線分段線性化對橫搖運動響應的影響

通過數值求解橫搖運動微分方程式(2)和(5)得到的橫搖運動響應時歷，可求得船舶橫搖運動的響應譜.除譜分析技術外，平均穿越率也可作為評估GZ曲線分段線性化近似對船舶橫搖運動響應造成影響的重要指標.平均穿越率v+(ξ)即橫搖運動在單位時間內穿越閾值ξ的次數.在路徑積分法求得橫搖運動響應的聯合概率分布函數后，v+(ξ)為

(14)

在隨機激勵系數σ=0.007的白噪聲激勵下，船舶橫搖運動的響應譜和平均穿越率見圖6.

圖6 白噪聲激勵下船舶橫搖運動的譜密度和平均穿越率

由于船舶橫搖阻尼系數較小，船舶橫搖運動的響應譜均為窄帶譜，且橫搖運動在單位時間內的過零次數v+(0)≈2π/ωp(ωp為響應譜的譜峰角頻率).對于原始動態系統(2)，由于非線性阻尼和非線性復原力的影響，ωp稍低于船舶橫搖自搖頻率ω0.

三種方案均對應動態系統(5)，其中方案1由于低估了船舶初穩性高且分段線性三角形的面積明顯小于GZ曲線所包圍的面積，響應譜的峰值顯著增大且譜峰角頻率ωp明顯小于ω0.采用方案2和3，則橫搖運動響應譜的譜峰角頻率ωp非常接近ω0且兩種方案對應響應譜的值在ω≥0.27 rad/s時大致相等.但方案3對應動態系統的復原力項較強，因此當0.15 rad/s≤ω<0.27 rad/s，橫搖運動響應譜的值比方案2低.此外，由圖6b)中各系統對應的平均穿越率可知，方案1所對應動態系統的橫搖運動響應與原始動態系統對應的橫搖運動響應有非常明顯的差別.

3.5 GZ曲線分段線性化對橫搖響應極值的影響

為評估以上三種分段線性化近似方案對船舶橫搖響應極值造成的影響，將每一種方案對應動態系統的橫搖運動平均穿越率分別與原始動態系統的橫搖運動平均穿越率進行對比，其結果見圖7.此外，為驗證路徑積分法的準確性，運用蒙特卡洛仿真法預測的橫搖運動平均穿越率也在上述圖中呈現.結果表明，路徑積分法精度高，運用該方法計算得到的平均穿越率可以作為橫搖運動響應極值的評估標準.

圖7 白噪聲激勵下原始動態系統和方案1～3對應的橫搖運動平均穿越率

由圖7a)可知，當橫搖角大于10°時，方案1對應動態系統的橫搖運動平均穿越率明顯高于原始動態系統對應的橫搖運動平均穿越率.結合圖6中對方案1的相關分析與評估，該方案對于橫搖運動響應的影響十分顯著，不推薦對CEHIPAR 2792樣船采用該方案進行GZ曲線線性化近似.在圖7b)中當橫搖角大于25°時，方案2對應動態系統的橫搖運動平均穿越率高于原始動態系統的橫搖運動平均穿越率.結合方案1和方案2可知，分段線性化近似造成復原力項變弱，將對橫搖運動的極值響應造成顯著影響.由圖6b)和圖7c)可知，在上述三種方案中，方案3的表現最優秀，但該方案在一定程度上增強了橫搖運動方程中的復原力項，使得橫搖運動響應極值較原始動態系統稍偏小.

4 結 論

1)路徑積分法能準確預報船舶橫搖運動的隨機響應，該方法具有精度高、計算時間短等優勢，可作為隨機橫搖運動響應預報，尤其是響應極值預報的高效數值方法.

2)通過對橫搖運動響應譜和平均穿越率進行分析與對比，結果表明GZ曲線分段線性化將改變橫搖運動方程的復原力項，對橫搖運動的極值響應造成顯著影響.

3)對于CEHIPAR 2792樣船，GZ曲線分段線性化方案1和方案2將造成橫搖運動響應極值顯著增大，方案3則引起橫搖運動的響應極值降低.方案3在三種方案中表現最為優秀，不推薦對該樣船使用方案1進行GZ曲線線性化近似.在實際的工程運用中，需要對更多船型進行相關研究，探討以上分段線性化方案的適用性.