大風浪環(huán)境下的船舶氣象航線優(yōu)化設計

崔康靖 鄭元洲 陳國成 胡衛(wèi)東

(武漢理工大學航運學院1) 武漢 430063) (內河航運技術湖北省重點實驗室2) 武漢 430063)(湄洲灣港引航站3) 莆田 351100)

0 引 言

隨著經濟全球化，海運量已占據全球商品貿易的80%以上[1]，而由氣象因素導致的海難事故時有發(fā)生．在氣象航線及其優(yōu)化方面國內外學者已有眾多研究．Park等[2]通過分階段來優(yōu)化油耗問題，在航速恒定下用A*算法求解轉向點位置，結合幾何規(guī)劃法重新分配各航路航速．Zaccone等[3]提出了一種基于三維動態(tài)規(guī)劃的船舶航線優(yōu)化方法，通過分配航路點尋求最優(yōu)航線．張進峰等[4]針對我國近海船舶進行了避臺航線優(yōu)化設計，基于動態(tài)規(guī)劃算法確定各航路點．在航線中船舶各航段的航速優(yōu)化方面，李明峰等[5]以浪高為約束條件構建氣象航線模型，引入懲罰函數(shù)應用于遺傳進化算法進行優(yōu)化．劉洋等[6-7]擬定危險風區(qū)為障航物，提出一種基于向量的多步滾動窗口優(yōu)化算法進行航線設計．

上述研究在模擬海上氣象時考慮的氣象條件都比較單一，應用算法在規(guī)避惡劣氣象區(qū)域時也不易設置合適的安全距離．針對大風浪問題，文中提出了一種結合人工勢場法和模擬退火算法的優(yōu)化方法．該方法在考慮風速和浪高的雙重氣象因素下對各航段航速影響，并以航程和航行時間為優(yōu)化目標，利用人工勢場法能較好的避讓大風浪區(qū)，并通過模擬退火算法替代原有的負梯度下降法對相應的航路點位置進行調整，增加可選航路點位置來優(yōu)化人工勢場法航線，以尋求最優(yōu)氣象航線．

1 船舶氣象航線設計的數(shù)學模型

1.1 航線模型定義

在考慮船舶航線設計時設定為整條航線由多個航路點組成，即將整條航線分為多個航路段，每段航路以恒向線方式表示，包含該段航路端點的航行速度、經緯度位置和航向角等信息．

設定航線中有n個航路點，用向量R表示．各個航路點的位置坐標為

R={(Lat1,Long1),…,(Lati,Longi),…,

(Latn,Longn)}

式中：Lati為第i個航路點坐標的緯度；Longi為第i個航路點坐標經度．

船舶在每段航路間將以恒向方式航行，航線中各航段的航向變量φ為

φ=[φ1,…,φi,…,φn]

(1)

式中：φi為第i段航路的航向z，由第i個和第i+1個航路點經緯度坐標決定，φi∈[0°,360°]，以正北方向為0°，順時針方向遞增．

V0為船舶在靜水中航行速度，而在實際航行中船舶容易受到風浪流等外界因素的影響而導致速度發(fā)生變化，因此用向量V表示各航路段航速，即

V=[V1,…,Vi,…,Vn]

(2)

式中：Vi為第i段航路的航速．

1.2 風浪對船舶速度的影響

船舶在海上航行時，受到風、浪、涌等外界環(huán)境因素影響，船舶行駛受到的阻力增加，在主機輸出功率不變時，船舶實際航速低于靜水時航速，這種現(xiàn)象稱為船舶失速．采用文獻[8]以最小二乘法為基礎迭代計算船舶失速公式.

Va=V0-(1.08h-0.126qh+2.77×10-3×

Fcosα)(1-2.33×10-7DV0)

(3)

式中：Va為船舶在風浪中實際速度，kn；V0為所給定的船舶靜水速度，kn;F為風速大小，m/s；D為船舶實際排水量，t;h為有義波高，m;q為船舶航向與浪向間的夾角，rad；α為船舶航向與風向之間的相對夾角，rad．

1.3 大風浪區(qū)的界定

在海上航行時，不同船舶所能應對的風浪情況并不相同．根據船舶的船型、噸位、結構強度、操作性能，以及裝載貨物等，風浪對船舶造成的影響會有很大差異．因此，在大風浪區(qū)的具體界定時，應當考慮該船的現(xiàn)實情況，只要某片海域的風浪環(huán)境超出該船的承受能力，便能將其定義為大風浪區(qū)域．船舶抗風浪等級的公式為[9]

(4)

h=2(M0/g+K1DLw1-K2W0Lw1-

(5)

式中：V10為船舶在改裝載狀態(tài)下可安全航行的最大風速；D為船舶排水量；lq為最小傾覆力臂；Af為船舶受風面積；K為穩(wěn)性衡準數(shù)；Zf為風作用力臂；h為船舶船體強度所能承受的最大浪高；M0為船舶的船中彎矩；g為重力加速度；K1為平水排水量力矩系數(shù)；Lw1為船舶水線長；K2為空船重量力矩系數(shù)；W0為空船重量；∑PiXi為所有載荷對船中的力矩；K3為波形排水量力矩修正系數(shù)；B為船寬．

2 動態(tài)航線規(guī)劃模型

2.1 目標函數(shù)

1) 航線航程 根據前文中航線模型的設定，航線總航程應為各個航路段相加之和，為

(6)

式中：L為總航線航程；Li為各分段航程．

在求航線航程時，需要根據坐標點的經緯度信息求實際航行距離，由前文介紹在各航路段中航線表示為恒向線，選用墨卡托算法計算航程，為

(7)

L=(Lat2-Lat1)×secφ

(8)

式中：Lat1,Long1分別為第一個點的緯、經度坐標；Lat2,Long2分別為第二個點的緯、經度坐標；φ為恒向線方向；L為兩點間航程；e為地球偏心率．

2) 航線航時 從起點到終點的總航行時間可由各航路段航行時間求和得到，為

(9)

式中：Ttotal為總的航行時間；Vi為第i段的實際航速．

2.2 人工勢場算法

文中在模型仿真中采用的是人工勢場算法，其原理為人通過引力場和斥力場共同作用使物體能夠沿梯度下降最快的方向運動，其中目標點對物體產生引力，引導物體朝向其運動．障礙物對物體產生斥力，避免物體與之發(fā)生碰撞．物體在路徑上每一點所受的合力等于該點所有斥力和引力之和．具體的算法模型如下.

常見的引力場：

(10)

式中：ξ為引力尺度因子；ρ(q,qgoal)為物體當前狀態(tài)與目標的距離．

引力即為引力場對距離的導數(shù)：

Fatt(q)=-Uatt(q)=ρ(qgoal-q)

(11)

斥力場：

(12)

式中：η為斥力尺度因子；ρ(q,qobs)為物體與障礙物之間的距離；ρ0為每個障礙物的影響半徑．

斥力為斥力場的導數(shù)：

(13)

引力和斥力尺度因子的取值對人工勢場算法航路規(guī)劃結果影響較大，如果取值不當，會導致航路不可行或航路點冗余[10]．取ξ=10且保持不變，η分別取5、8、10、12、15、18時人工勢場算法規(guī)劃結果見圖1．

圖1 η取不同值時對航線規(guī)劃的影響

由圖1可知，當η取5、8、10時，會出現(xiàn)航路途徑大風浪障礙區(qū)的情形，規(guī)劃航線不可取；而相比η=12，當η取15和18時，發(fā)現(xiàn)航路點有明顯冗余現(xiàn)象，因此本文參數(shù)選擇取ξ=10，η=12，此時人工勢場算法得到的規(guī)劃航路最佳．

2.3 模擬退火算法

經由人工勢場法得到航線并非是最優(yōu)航線，船舶駛于航線中途的航路點時，一方面可能會面臨引力和斥力的合力矢量為零陷入局部最小而終止路徑規(guī)劃．另一方面人工勢場法可能受引、斥力增益系數(shù)的影響造成迭代步伐過大，引起抖動現(xiàn)象，所以需采用模擬退火算法進行動態(tài)航線優(yōu)化．

在人工勢場法航線中加入模擬退火算法，可以通過增加隨機擾動值調整新航路點，避免陷入局部最小值，還可以搜索大風浪區(qū)附近的航路點位置，替代原有的負梯度下降受力控制的方法，避免人工勢場法中參數(shù)系數(shù)對航路點選取的影響，可有效緩解抖動現(xiàn)象．

每次經過調整之后的航路點組合便構成了一條新航線，可以利用目標函數(shù)值比較以及劣航線的接受準則進行處理．將新航線的目標函數(shù)值設為f(A)，當前最優(yōu)航線目標函數(shù)值設為f(B)．若f(A)

Metropolis準則判斷具體方法為

設當前航線為xi，新航線為xj，其目標函數(shù)值分別記為f(xi)和f(xj)，新航線被接受可以作為最優(yōu)航線的概率Pt為

(14)

在判斷過程中，系統(tǒng)將產生一個隨機數(shù)ε∈(0,1)，并與Pt進行比較，若Pt>ε，則將新航線設為當前最優(yōu)航線；否則，舍棄新航線，繼續(xù)使用當前航線作為最優(yōu)航線．從式(14)可見，溫度越低，接受概率Pt的值也就越小，即劣航線被接受的概率會越來越小，當退火到最小值時，劣航線基本被全部淘汰．

經過多次實驗，選用了以下參數(shù)：初始溫度t0=100；衰減系數(shù)α=0.7；終止溫度tf=10；迭代次數(shù)L=500．

模擬退火算法的具體流程為

步驟1初始參數(shù)設置 初始溫度t0，衰減系數(shù)α，終止溫度tf，迭代總次數(shù)L．

步驟2初始值設置 初始航線x0，當前航線xcurrent，新航線為xnew，最優(yōu)航線為xbest，初始目標值設為E0，當前航線航程為Ecurrent，新航線航程為Enew，最優(yōu)航線航程為Ebest，當前迭代次數(shù)l=0．

步驟3調整航路點 生成隨機數(shù)n和n個2維擾動變量ΔLat,ΔLon，隨機選擇當前航線上n個航路點(除始點和終點)，令Lat′=Lat+ΔLat,Lon′=Lon+ΔLon,得到新航線xnew，計算其目標值Enew，迭代次數(shù)l=l+1．

步驟4計算新航線xnew，與最優(yōu)航線xbest之間的航程差值ΔE，ΔE=Enew-Ecurrent，如果ΔE<0，那么令xbest=xcurrent=xnew，Ebest=Ecurrent=Enew，轉至步驟6．否則轉至步驟5．

步驟5令ΔE=Enew-Ecurrent，若exp(-ΔE/ti)>ε，則xcurrent=xnew，Ecurrent=Enew，其中ε為(0，1)內的隨機數(shù)；否則xcurrent和Ecurrent的值不變．

步驟6若l≥L，令ti+1=αti、l=0，轉至步驟7；否則返回步驟3．

步驟7若ti+1≤tf，算法終止；否則返回步驟3．

3 仿真結果及分析

3.1 氣象資料獲取和處理

研究采用歐洲中期天氣預報中心(european center for medium-range weather forecasts，ECMWF)提供的氣象數(shù)據．ECMWF主要提供中長期天氣預報如實時的天氣、海況信息來為遠洋航行船舶提供服務．由于本文在模型中主要考慮風浪對船舶的影響，所以采用ECMWF提供海平面10 m高處的風速與風向數(shù)據．氣象數(shù)據文件格式采用網絡通用數(shù)據格式NetCDF(Network Common Data Format)，數(shù)據更新間隔為6 h，數(shù)據精度的網格大小為1°×1°．利用Matlab可以將氣象數(shù)據可視化，更直觀顯示風浪分布情況[11]．圖2為2019年6月太平洋區(qū)域第8日12:00時刻的風速數(shù)據信息．圖3為2019年6月太平洋區(qū)域第8日中12:00時刻的浪高數(shù)據信息．

圖2 風速與風向數(shù)據信息

圖3 浪高數(shù)據信息

3.2 結果及分析

選用大海域太平洋作為航行水域，起點為日本橫濱港(35°N,141°E)，終點為美國舊金山港(37°N,123°W)．仿真選用的船舶為標準排水量23 740 t的集裝箱船S175，船舶靜水速度設為15 kn，浪高與風級合并考慮，僅將船舶與風浪的遭遇角作為海況考慮[12]．選用2019年5月30日和6月8日的大風浪數(shù)據在仿真中展示效果圖(見圖4)，航線為大圓航線．

圖4 航行途徑大風浪區(qū)示意圖

大圓航線為球面上兩點之間距離最短的航線，即在理想環(huán)境下船舶走大圓航線將是最優(yōu)選擇，但在現(xiàn)實中船舶并不容易以大圓航線航行，如考慮危險氣象環(huán)境時，航線必將有所調整．因此本文在進行航線規(guī)劃時將以大圓航線為基準并利用人工勢場算法進行動態(tài)調整，在經過人工勢場法避障之后的航線其航程里數(shù)并不是最優(yōu)解，因此再利用模擬退火算法對該航線進行優(yōu)化，優(yōu)化結果見圖5，幾種航線的航程與航時對比見表1．

圖5 優(yōu)化航線圖

表1 三種航線結果對比

由圖5和表1可知：單純的人工勢場算法所規(guī)劃航線可以很好的避開大風浪區(qū)域，且相較大圓航線船舶平均航速有所提高，但航行距離和航行時間有明顯增加．在人工-模擬法航線中則有較好的改進，雖總航程相較大圓航線有所增加，但總航時明顯降低而且平均航速也有明顯提高．

4 結 束 語

文中提出了一種結合人工勢場法和模擬退火算法的氣象航線優(yōu)化方法，并用Matlab進行模擬仿真，在考慮風速和浪高氣象因素下，調整風浪中各航段航速，優(yōu)化縮短航程和航時．該方法能根據船舶自身條件設置安全距離，有效避過大風浪區(qū)，人工-模擬退火算法對航線優(yōu)化效果良好，可以顯著提高船舶平均航速，縮短航行時間．

文中的氣象因素僅考慮了風速和浪高，為更準確模擬海上環(huán)境，后續(xù)研究中還應考慮云霧、暴雨等更多海上氣象要素的影響，并從經濟層面考慮，構建如油耗等模型以進一步完善優(yōu)化算法．