趙江艷， 劉 文， 張明陽， 徐強強， 謝 平,4

(1. 武漢大學 水資源與水電工程科學國家重點實驗室，武漢 430072； 2. 武漢理工大學 內河航運技術湖北省重點實驗室，武漢 430063； 3.國家水運安全工程技術研究中心，武漢 430063； 4. 國家領土主權與海洋權益協同創(chuàng)新中心，武漢 430072)

長江屬于季節(jié)性大型河流，受季節(jié)性因素和通航樞紐下泄流量等因素影響，航道水位變化明顯，尤其是長江中游航道，年最高水位與年最低水位的多年平均值之差高達13 m，對船舶通航安全的影響顯著。[1]通過對長江航道的事故進行統計分析發(fā)現，2013年和2014年枯水期發(fā)生的水上事故均占當年事故總數的50%以上，其中長江中游航道是長江干線的礙航瓶頸，有近20處礙航淺灘及10多個重點淺水道，通航條件十分不穩(wěn)定，極易發(fā)生擱淺事故，從而造成礙航局面。[2]因此，對水位變化下的長江中游船舶擱淺風險進行分析具有實際應用價值。

船舶擱淺風險分析已在水上交通工程領域引起廣泛關注，但不同通航環(huán)境下的船舶擱淺風險的定量計算依然是亟待解決的問題。通常采用的研究方法有故障樹方法、貝葉斯方法及概率模型等。[2-8]AMROZOWIEZ[6]將船舶擱淺事故分為動力性擱淺和漂移性擱淺，分別建立故障樹，確定失效概率；KITE-POWELL等[7]根據統計學中的貝葉斯方法提出船舶擱淺的風險計算模型，可有效預測船舶擱淺風險，定量計算相應的直接經濟損失和間接經濟損失；ARSHAM等[8]利用船舶交通和航道復雜度建立評估船舶擱淺概率的模型等。從以上研究成果中可看出，大多數研究是從整體的角度或從船舶的擱淺現象出發(fā)，闡述有關船舶擱淺分析的總體看法。但是，船舶擱淺分析中的影響因素與研究區(qū)域的交通條件、自然條件緊密相關，因此無法推廣到所有地區(qū)。此外，由于交通環(huán)境和自然環(huán)境會隨時間變化，因此船舶的擱淺分析應是一個適應新的風險因子的動態(tài)過程。對此，本文將船舶擱淺問題視為首次穿越問題，根據隨機振動的首次穿越失效理論建立首次穿越失效(船舶擱淺)模型。首次穿越問題是指一個受隨機參激或隨機外激作用的動力學系統，其響應在相空間上發(fā)生隨機擴散。當系統響應第一次離開某個安全域或允許域時，就稱系統發(fā)生首次穿越。首次穿越在工程[9]、物理[10]、生物[11]、醫(yī)學[12]及經濟[13]等學科中廣泛存在。QUY等[14]將首次穿越理論引入到船舶擱淺的概率模型中，初步驗證該理論在解決船舶擱淺問題上的適用性，但在具體的建模過程中僅考慮航速變化對船舶吃水的影響，沒有解決航道水位變化對船舶擱淺風險影響的問題。

對此，本文重點研究長江中游航道水位變化對船舶通航安全的影響。為保證分析結果的可靠性，運用隨機振動的首次穿越失效理論構建船舶擱淺概率計算模型，定量分析不同水位下船舶的擱淺概率；同時，給出船舶擱淺概率分別為0和1時對應的水位臨界值，用以驗證實際的設計最低通航水位是否符合通航安全的要求。

1 問題描述與建模

船舶擱淺事故是一種受多種因素影響的水上交通事故，在分析船舶擱淺風險的影響要素時，除了不同時期的水位變化之外，人為因素、船舶因素等不受水位變化影響的因子也會對船舶的擱淺風險造成一定的影響。若進行全面考慮，則會增加求解問題的難度，甚至無法求解。因此，本文在船舶擱淺風險的分析中主要考慮水位變化對船舶擱淺的影響，參考波浪載荷下船舶運動響應研究中的譜分析方法的基本思路分析船舶的垂向運動，基于首次穿越理論構建船舶擱淺風險計算模型。

1.1 水位變化對船舶擱淺的影響

擱淺事故表明，航道水深不足是導致船舶發(fā)生擱淺的直接原因。[15]由于航道深度在一定時期內是不變的，因此航道水位的變化也可間接描述航道水深的變化，從而建立水位與船舶擱淺風險之間的關系。富余水深、航道水深及船舶吃水之間的關系見圖1，其中：H為航道水深；W為當時當地實際水位；d為圖示水深；D為船舶吃水；Δh為富余水深。富余水深主要包括船舶航行下沉量和觸底安全富余量2部分，其主要作用是防止船舶觸底，保證船舶安全航行。有關系式

H=W+d

(1)

允許船舶的最大吃水

Dmax=H-Δh=W+d-Δh

(2)

一般而言，增大富余水深Δh的取值會提高船舶的通航安全，同時會降低航道的通航能力。因此，為平衡通航安全與通航能力之間的關系，須對富余水深Δh進行合理取值。《內河通航標準》[16]根據航道等級規(guī)定了相應的富余水深。

由式(2)可計算任意水位下允許船舶的最大吃水，當實際吃水大于該最大吃水時，船舶發(fā)生擱淺。由生產實踐可知，即使船舶的載貨量一定，船舶的吃水也會隨著通航環(huán)境的變化而變化。這是因為船舶航行時不僅會受波浪的影響產生升沉運動的響應位移，還會受周圍水流的影響出現下沉現象。船舶航行下沉量(即船舶動吃水)是指船舶航行中的吃水較靜止時的增加量，本文采用實際航運中推薦的公式，即

(3)

式(3)中：CB為方形系數；v為船速，kn。

因此，計算船舶擱淺概率的前提是確定當時當地的水位及船舶的動態(tài)吃水。對于載貨量一定的單艘次船舶來說，水位越高，船舶的擱淺風險越小。

1.2 水位變化下船舶擱淺風險計算模型

1.2.1建模的基本假設

為避免模型自身的復雜性帶來不確定性風險，參考波浪載荷下船舶運動響應研究的譜分析方法[17-20]，確定建模的2個基本假設：

(1)波浪表面是一個各態(tài)歷經的零均值高斯隨機過程；

(2)船舶為一個線性系統。

假設(1)與隨機波浪理論的思想一致，認為在任意較短的時間內，波浪符合零均值的平穩(wěn)正態(tài)隨機過程，其統計特性[21-23]可用功率譜密度來確定，即

(4)

式(4)中：S(ωi)為波浪功率譜密度函數，表示一個隨機波浪所包含的所有頻率及其對應的能量。

在波浪載荷的作用下，船舶會產生位移、切力、彎矩及應力等各種響應。[24]若將波浪的作用記為Q(t)，則有

(5)

式(5)中：mQ為Q(t)的均值；An和Bn為互不相關的隨機變量。An和Bn有均值為零、方差相等的特征，即mAn=mBn=0，DAn=DBn=Dn。

波浪的方差為

(6)

(7)

因此，當Δω→0時，可得

(8)

式(7)和式(8)中：SQ(ω)為波浪譜。

對于單自由度系統，假定受一簡諧干擾力作用，則振動方程為

(9)

式(9)中：q(t)為廣義位移；ω0為系統固有頻率；ξ為黏性阻尼系數；Q0為干擾力幅值；ω為干擾力頻率。

該系統的響應為

q(t)=q0eiωt=H(ω)Q0eiωt=H(ω)Q(t)

(10)

根據假設(2)，將垂向運動的船舶視為單自由度的線性系統，將波浪作為系統的輸入量，則船舶垂向運動的響應即為系統的一種輸出，其值由波浪載荷和響應函數共同決定。響應函數的物理意義是各種不同頻率的單位波幅余弦波產生的船舶響應。若定義H(ω)的模|H(ω)|為響應函數，則有

|H(ω)|=q/Q0

(11)

利用傅里葉變換可知

(12)

因此，作為船體—波浪線性系統的隨機輸出，波浪載荷在短時間內也服從零均值正態(tài)分布，其幅值遵循Rayleigh分布特性。[25-26]

1.2.2船舶垂向運動建模

建立坐標系O-xyz見圖2，當將船體作為剛體看待時，船舶的搖蕩運動是沿3個坐標軸方向的位移和轉動，具有6個自由度。在不考慮舵力作用的情況下，船舶的六自由度運動可分解為縱向和橫向2組耦合運動，其中：縱向運動由升沉、縱搖和縱蕩組成；橫向運動由橫蕩、橫搖和艏搖組成。

縱向耦合運動與橫向耦合運動之間的影響很小，且自由狀態(tài)下船舶的縱蕩運動對船舶升沉運動位移的分析影響不大。因此，在求解船舶在波浪中的升沉運動位移時，不用考慮船舶的縱蕩運動，只需研究船舶升沉和縱搖的耦合運動。

船舶升沉和縱搖的耦合運動的微分方程為

(13)

(14)

(15)

式(15)中：

水動力系數可通過計算近似得到。由此可得

(16)

式(16)中：Az=E-1M；Bz=E-1N；隨機波浪W為有色噪聲。

若選擇船舶升沉運動的位移和縱搖傾角為測量狀態(tài)，則測量方程為

Y=CzX+V

(17)

本文利用波浪擾動成形濾波器計算船舶受到的波浪干擾，其核心思想是將高斯分布的白噪聲通過一個波浪擾動成形濾波器G(s)，從而得到作為隨機波浪擾動仿真過程的有色噪聲。[17]若將成形濾波器視為船舶運動系統的一部分，則可將白噪聲信號視為新系統的擾動。新系統的組成見圖3。

圖3中：GZ(s)為波浪對船舶縱向擾動力的成形濾波器；GM(s)為波浪對船舶縱向擾動力矩的成形濾波器。通常選取成形濾波器

(18)

式(18)中：s為Laplace算子；G0為增益；ω0為中心頻率；ξ為阻尼系數。wZ(t)和wM(t)均為白噪聲信號，其均值為0，方差為1；Z(t)為船舶升沉運動干擾力；M(t)為縱搖干擾力矩；q(t)為最終求得的船舶垂向響應位移。

因此，當波浪的白噪聲信號wZ(t)和wM(t)已知時，可通過上述成形濾波器得到船舶升沉運動干擾力Z(t)和縱搖干擾力矩M(t)，隨后結合船舶運動的狀態(tài)方程求出船舶垂向運動響應位移q(t)。

可運用Simulink工具箱對上述新系統進行仿真，從而得到船舶垂向運動響應位移q(t)。仿真流程圖見圖4。

1.2.3船舶擱淺風險量化建模

由“1.2.1”節(jié)可知，船舶在波浪激勵下產生響應，由于激勵是服從高斯分布的，運動響應是線性的，最后得到的船舶運動也符合高斯分布。因此，若運用船舶擱淺的數學關系式D>H，即船舶吃水大于航道水深，則船舶擱淺問題可看成首次穿越失效問題：將航道水深表示成隨時間變化的閾值，當船舶垂向的位移響應首次穿越該閾值時，認為船舶發(fā)生擱淺。首次穿越概率(即船舶的擱淺概率)可通過相關的隨機理論計算得到。

通過分析可知，船舶在時域內穿越的失效閾值可表示為

φ(t)=W(t)+h-D

(19)

式(19)中：W為當時當地實際水位；h為航道圖示水深；D為船舶吃水。

圖5為船舶擱淺(首次穿越失效)示意，若在船舶航行期間(即[0，T]內)，其升沉運動響應位移q(t)穿過φ(t)，則船舶首次穿越失敗，即發(fā)生擱淺事故。

令β=φ(t)，則船舶響應的期望穿閾率[27]可表示為

(20)

式(20)中：m0和m2分別為船舶運動響應譜的零階矩及二階矩。m0和m2的計算式可表示為

(21)

(22)

式(21)和式(22)中：Sq(ω)為描述在頻率域ω上響應x(t)的均值分布的功率譜密度函數，即船舶運動響應譜。船舶運動的響應譜實質上是一個功率譜密度函數，描述信號功率隨頻率的分布狀況。

根據“1.2.2”節(jié)可求出船舶垂向運動的響應位移q(t)。將q(t)視為一個隨時間變化的隨機信號。以Ts對隨機信號進行等間隔抽樣，得到隨機序列q(nTs)，記為X(n)；隨后通過直接法、間接法、Bartlett法或Welch法計算得到其譜密度函數。

若應用首次穿越理論，則船舶首次穿越失效的概率(即船舶的擱淺概率)可表示為

P(β|T)=1-exp(-vβT)

(23)

2 實例計算與結果分析

選取一艘滿載排水量為2 858 t的船舶作為研究對象，運用本文提出的方法計算其擱淺風險。該船的具體參數為：船長75 m；垂線間長71.33 m；型寬13.4 m；型深5.4 m；設計吃水4.2 m；方形系數0.704；水線面系數0.980；中橫剖面系數0.980；滿載排水量2 858 t。

根據水文和氣象條件，將內河船舶航行區(qū)域劃分為A，B，C等3級(見表1)。長江中游主要屬于B級航區(qū)，根據該航區(qū)的波浪要素及船舶尺度確定

表1 我國內河各級航區(qū)波浪要素

波浪成形濾波器的參數見表2。

表2 波浪成形濾波器的參數

將上述船舶尺度、圖示水深(3.2 m)及船舶與波浪的遭遇頻率(π/5)等數據作為輸入，結合MATLAB/Simulink仿真平臺(如圖4所示)，便可求出波浪干擾力、波浪干擾力矩及船舶升沉運動的響應位移。設定仿真時間為60 s，波浪對船舶升沉運動的干擾力和干擾力矩分別見圖6及圖7。從圖6和圖7中可看出：船舶受到的波浪載荷符合正態(tài)分布，其幅值符合Rayleigh分布。因此，當將船舶視為一個線性系統時，船舶的運動響應過程是各態(tài)歷經的平穩(wěn)隨機過程，運動幅值也符合Rayleigh分布(見圖8)。

運用基于MATLAB的譜估計程序實現船舶升沉運動響應位移到船舶響應譜的轉化，最終得到船舶升沉運動的響應譜見圖9。隨后運用首次穿越理論計算得到船舶擱淺概率與水位的關系見圖10。

從圖10中可看出：隨著航道水位下降，船舶擱

淺風險逐漸增大；隨著航道水位升高，船舶擱淺風險逐漸減小。船速越大，船舶擱淺概率為0對應的水位臨界值越大；船速越小，船舶擱淺概率為0對應的水位臨界值越小。具體臨界值見表3。

表3 不同船速v下的水位臨界值

通過分析可知：船速越大，船舶的運動響應越顯著，船舶升沉運動的響應幅值越大；同時，船速越大，船舶航行下沉量越大，故其對通航水位的要求越高。因此，利用上述方法不僅可確定某一水位下船舶的擱淺概率，還能由船舶擱淺風險的容忍度確定航道的最低通航水位，從而驗證實際的設計最低通航水位是否符合要求。

此外，由于當時當地的實際水位與圖示水深之和即為航道水深，故該方法也可為航道部門確定航道維護水深提供一定的參考。

3 結束語

本文提出一種基于首次穿越失效理論的船舶擱淺風險的定量計算方法，建立水位變化下的首次穿越失效(船舶擱淺)模型。同時，以長江中游航道的典型船舶為例，定量計算其在淺水航段的擱淺風險，最終得到船舶擱淺概率與航道水位的對應關系。從結果上看：在同一船速下，航道水位越低，船舶的擱淺風險越大；船速越大，船舶擱淺概率為0和1對應的水位臨界值越大，與實際情況一致。研究結果表明，該方法不僅能確定某一水位下船舶的擱淺概率，還可確定船舶擱淺概率為0和1時對應的水位臨界值。此外，利用該方法還能求得任一可接受擱淺概率下的最低水位值。因此，該方法可用于確定船舶擱淺風險容忍度下的最低通航水位，為最低通航水位的設計提供新的思路，為航道部門確定航道維護水深和航運經營人對船舶進行配載提供參考。

[1] 劉婷.水位變化對長江中游航運物流成本的影響機理及模型研究[D]. 武漢: 武漢理工大學, 2014.

[2] 張笛. 枯水期長江通航風險評價和預測方法研究[D]. 武漢: 武漢理工大學, 2011.

[3] 陳婷婷, 施朝健, 雷琴. 船舶擱淺模糊故障樹的灰色模糊關聯分析[J]. 中國航海, 2015, 38(3):83-88.

[4] 吳金龍, 劉大剛, 范中洲,等. 大風浪條件下渤海海區(qū)重點船舶風險評估模型[J]. 中國航海, 2012, 35(2):89-92.

[5] MAZAHERI A, MONTEWKA J, KUJALA P. Modeling the Risk of Ship Grounding—A Literature Review from a Risk Management Perspective[J]. WMU Journal of Maritime Affairs, 2014, 13(2): 269-297.

[6] AMROZOWIEZ M D. The Need for a Probabilistic Risk Assessment of the Oil Tanker Industry and a Qualitative Assessment of Oil Tanker Groundings[Z]. Massachusetts Institute of Technology, 1996.

[7] KITE-POWELL H L, JIN D L, JEBSEN J L, et al. Investigation of Potential Risk Factors For Groundings of Commercial Vessels In U S Ports[J]. International Journal of Offshore & Polar Engineering, 1999, 9(1)：16-21.

[8] ARSHAM M, JAKUB M, PENTTI K, et al. Assessing Grounding Frequency Using Ship Traffic and Waterway Complexity[J]. Journal of Navigation, 2015, 68(1):89-106.

[9] ZHANG L Y, SUN L M, SHANG Z Q. Real-Time Reliability Assessment Based on Acceleration Monitoring for Bridge[J]. Science China Technological Sciences, 2016,56(8): 1-11.

[10] SHARMA K, KUMAR N. First-Passage Time: Lattice Versus Continuum[J]. Physical Review E, 2012, 86(3): 032104.

[11] HU Z, CHENG L, BERNE B J. First-Passage Time Distribution in Stochastic Processes with Moving and Static Absorbing Boundaries with Application to Biological Rupture Experiments[J]. Journal of Chemical Physics, 2010, 133(3): 034105.

[12] SIROVICH L, KNIGHT B. Spiking Neurons and the First Passage Problem[J]. Neural Computation, 2011, 23(7): 1675-1703.

[13] BAKSHI G, PANAYOTOV G. First-Passage Probability, Jump Models, and Intra-Horizon Risk[J]. Journal of Financial Economics, 2010, 95(1): 20-40.

[14] QUY N M, VRIJLING J K, GELDER P, et al. On the Assessment of Ship Grounding Risk in Restricted Channels [C]. 8th International Conference on Marine Sciences and Technologies, 2006: 294-299.

[15] 任亞磊, 牟軍敏, 李亞軍, 等. 進港航道中船舶擱淺概率的蒙特卡洛仿真[J]. 船舶力學, 2014, 18(5): 532-539.

[16] 交通運輸部.內河通航標準：GB 50139—2014[S]. 北京: 中國計劃出版社, 2014.

[17] 彭秀艷, 張文穎, 何志勇. 基于譜估計方法的船舶縱向運動建模預報研究[J]. 自動化技術與應用, 2008, 27(11): 1-4.

[18] DINIS M S, PIERSON J W J. On the Motion of Ships in Confused Seas[Z]. Routledge, 1953,22(23):8144-8154.

[19] 李成成. 船舶在波浪中運動時域數值模擬研究[D]. 上海：上海交通大學, 2010.

[20] 杭棟棟. 基于船舶運動響應數據的海浪估計方法研究[D]. 哈爾濱：哈爾濱工程大學, 2013.

[21] 徐德倫, 于定勇. 隨機海浪理論[M]. 北京：高等教育出版社, 2001:111-115.

[22] 俞聿修. 隨機波浪及其工程應用[M]. 大連: 大連理工大學出版社, 2011:60-81.

[23] KIMURA A. Statistical Properties of Random Wave Groups[Z]Coastal Engineering,1980:2955-2973.

[24] 王學亮. 大型起重船在海浪中的運動響應研究[D]. 天津: 天津大學, 2003.

[25] 吳小平. 基于切片理論的波浪載荷直接計算[J]. 船舶與海洋工程, 2010(4):21-25.

[26] PEREZ T, BLANKE M. Ship Roll Damping Control[J]. Annual Reviews in Control, 2012, 36(1):129-147.

[27] QUY N M, VRIJLING J K, VAN GELDER P. Risk-and Simulation-Based Optimization of Channel Depths: Entrance Channel of Cam Pha Coal Port [J]. Simulation, 2008, 84(1): 41-45.