基于北極航線通航期冰情的多種船舶航速優化

閔德權，江可鑒，劉 茹

(大連海事大學 交通運輸工程學院，遼寧 大連 116026)

0 引 言

隨著全球氣候持續變暖，北極海冰的覆蓋面積逐漸縮小，北極航線的通航問題得到了更多的關注與研究。總體來看，北極航線可分成東北航道、西北航道以及穿極航線。北極航線東北航道的通航期從7月中下旬開始，在10月中下旬結束[1]。但相較傳統航線，其具有航線時間和航線成本上的優勢。

針對船舶的航速優化，需先對航線進行總體分析。ZHENG Binglei等[2]通過對北極航線上船舶航行的航速、干擾以及延遲等因素進行分析，建立以航線運營成本最小為目標的優化模型；關曉光等[3]在對北極航線的東北航線進行航次的收益和效率分析的基礎上，建立環繞模式集裝箱航線優化模型；J.H.NAM等[4]基于海冰和環境等數據，建立交互式的仿真系統，從而求解出最佳的北極航線。在上述研究的基礎上，可對船舶的航速做進一步分析。A.AFONIN等[5]利用包含的不同船舶有關信息的地理信息系統(GIS)，并結合冰情數據，統計處理冰區航行船舶的航速以識別航速減少的趨勢；FU Shanshan等[6]建立概率模型以預測船舶在北極水域的困阻概率，并使用貝葉斯信念網絡(bayesian belief networks)模型預測船舶的航速；李振福等[7]基于冰區氣象的條件，分析海況對船舶的干擾力和冰阻力等因素，建立以不同冰情航速優化模型。

基于上述研究，筆者以北極東北航線為基礎，建立通航期冰情條件下的混合整數雙目標船舶航速優化模型。使用線性加權和離散化的方法對模型進行處理，選擇MATLAB軟件對模型進行求解，以期得出不同船舶的航速優化方案。

1 航線通航期冰情分析

北極東北航線的可通航時間主要集中于7月到10月[8]，將上述4個月作為通航期，并以通航期的冰情數據為基礎對船舶航速進行優化分析。航線的冰情包括海冰厚度、海冰密集度等數據，目前可使用的數據庫有：美國國家冰雪數據中心(National Snow & Ice Data Center)數據庫、德國不來梅大學(The University of Bremen)的海冰數據庫等。

選擇德國不來梅大學海冰數據庫中2020年7月至10月數據，對航線經過海域的海冰密集度和海冰厚度數據進行分析。

1.1 通航期海冰密集度數據分析

不來梅大學海冰數據庫所提供的海冰密集度數據以HDF、netCDF等數據格式進行存儲[9]，可使用相應軟件進行讀取。

在航線的涉及海域中，渤海在冬季部分海域開始結冰，到次年的3月海冰完全消失[10]。巴倫支海西南部存在常年無冰區，且在夏季期間出現無冰情況[11]。挪威海和北海在夏季可通航，亦不做冰情數據分析。在現有海冰密集度數據中，白令海7至10月海冰密集度均值在1%以下，其海冰密集度僅作為計算參數進行使用。因此設定楚科奇海、東西伯利亞海、拉普捷夫海以及喀拉海作為冰區海域進行分析，其通航期月均海冰密集度數據如圖1。

從圖1中海域的海冰密集度數據可以看出，楚科奇海在2020年的4個月的月均冰密集度呈現出較高態勢，最高冰密集度接近60%。拉普捷夫海和喀拉海的月均冰密集度較低，均在10%以下。

圖1 冰區海域的海冰密集度

從總體來看，上述海域的冰密集度數據在8月或9月較低，但在7月或10月的波動較為劇烈。除喀拉海以外，其他3個海域在通航期4個月的變化趨勢基本一致。

1.2 通航期冰厚度數據分析

使用不來梅大學海冰數據庫的薄海冰厚度(the thickness of thin sea ice, SIT)庫所記錄的北極海冰數據[12]，其中部分數據存儲于內部資料庫。

選擇2020年7至10月的北極海冰月均厚度數據作為后期模型求解的參考數據，不再做具體海域的數值分析，其數值如表1。

表1 海冰的月均厚度

2 船舶阻力模型構建

在冰區航行的船舶不僅要考慮水的作用效應，還要考慮海冰對船舶的影響。因此，分別對船舶的靜水阻力和冰阻力進行模型構建。對阻力模型進行以下設定：①對靜水阻力和冰阻力區分分析；②不考慮船舶附屬體帶來的附體阻力。

2.1 靜水阻力模型

當船舶在無冰的靜水中航行時，其所受的阻力(即裸船體阻力)f0的計算如式(1)[13]：

(1)

式中：C0為裸船體阻力的阻力系數；ρ為海水密度；v為船舶的航行速度；S為船舶的濕表面積。

船舶的濕表面積計算如式(2)[14]:

(2)

式中：Lw為船舶的水線長度；d為船舶的吃水；B為船舶的型寬；Cb為方形系數。

依據文獻[13]，船舶在無冰的靜水中航行所受到的阻力f0可分為摩擦阻力f1、粘壓阻力f2以及興波阻力f3：

f0=f1+f2+f3

(3)

考慮到船體表面的粗糙度，計算船舶摩擦阻力時需額外加入粗糙度補貼系數CΔ[13]。

CΔ可由船長L進行選取，其具體值如表2[7]，設定表中船長不計邊界值。

表2 粗糙度補貼系數選取

2.1.1 摩擦阻力

摩擦阻力計算如式(4)：

(4)

依據文獻[13]，船舶的摩擦阻力系數可以使用1957年國際船模實驗池會議所提出的計算公式：

(5)

式中：Re為雷諾數；γ為海水的運動粘性系數。

(6)

2.1.2 粘壓阻力

依據文獻[13]，粘壓阻力系數可使用巴甫米爾提出的近似公式進行計算。粘壓阻力f2和粘壓阻力系數C2的計算如式(7)、式(8)：

(7)

(8)

式中：Am為船舶中橫剖面的水線下面積，后簡稱為中橫剖面積。Am可通過式(9)進行計算，Cm為中橫剖面系數。

Am=CmBd

(9)

Lr為船舶的去流段長度，為了避免嚴重的船舶的航行漩渦效應，選取去流段長度為[15]：

(10)

2.1.3 興波阻力

通過對Michell積分進行簡化，可對興波阻力進行計算[16]。具體公式如式(11)：

(11)

式中：λ=secθ，θ為凱爾文波系(Kelvin wave system)相對于x軸的角；ω1與ω2與船舶表面函數的縱向梯度有關。

簡化后的興波阻力公式也難以直接計算出船舶受到的興波阻力。而文獻[17]提供的測試程序可進行簡化測算，該程序通過計算Michell積分來估計Wigley船型的波阻力，且提供的阻力數據可作為計算興波阻力的基礎。

首先，在上述測試程序中輸入相應的船舶參數，運行程序得到多航速對應的阻力數據。其次，對結果數據進行MATLAB的傅里葉級數擬合，得到多參數的興波阻力和船舶航速的擬合函數，其擬合函數如式(12)：

(12)

式中：a0、aj、bj以及z為對應的擬合參數。

2.2 冰阻力模型

設定船舶在通航期航行所受到冰阻力為浮冰(pack ice)阻力fe，可由式(13)進行計算[18]。

fe=CeρegBhv2Wn

(13)

式中：ρe為海冰密度；h為海冰厚度；W為海冰密集度；Ce浮冰阻力系數；n為海冰密集度的冪次數；g為重力加速度。依據文獻[18]，使用n=2對于計算具有極小的誤差。

依據文獻[18]，加拿大海洋技術研究所經過多年試驗得出lnCe和lnFr呈線性關系，Ce和Fr的關系如式(14)。

Ce=4.4Fr-0.8267

(14)

式中：Fr為浮冰的弗勞德數，計算如式(15)：

(15)

3 航速優化模型構建

3.1 問題概述

在傳統航線的運營中，船舶通常有固定的掛靠港口[19]，且在港口間的貨運量相對確定。船公司可以根據自身的經營狀況，對航線上運營的多種船舶進行合理調配，以期得到最大的經濟效益。

在北極航線中，船舶一般只在通航期進行航行，難以在航線內進行多航段和多船組合的調配。因此，在現有船舶條件下，建立北極航線單船單航次的雙目標航速優化模型，在船舶的港口使費、固定成本及燃油成本在內的運營成本最低，以及船舶的航行時間和停泊時間在內的運營時間達到最低的條件下，并求解船舶的最佳航速。依據第2節中的船舶阻力公式，可對船舶克服阻力進行航行的做功進行轉化求出船舶的燃油成本，進而確定船舶運營成本。同時對航線進行航段劃分，使不同船舶的航速在不同航段達到最優。對船舶的航行過程做出3點假設：①航線上僅配置單個船舶，且船舶在航線的每個航段勻速航行；②不考慮航線極端氣象狀況，選取理想航行狀態。③不使用額外破冰船進行輔助航行。

3.2 規劃模型

第1個目標函數如式(16)，表示船舶完成航線航行的運營成本最小化。式(16)由船舶的港口使費、船舶航行的固定成本、船舶的燃油成本3個部分構成。

(16)

第2個目標函數如式(17)，表示船舶完成航線航行的運營時間最小化。式(17)由船舶航行時間、船舶停泊時間兩部分組成。

(17)

(18)

使用燃油的燃燒熱值對船舶做功進行燃油量的換算，從而求出燃油成本。船舶的燃機和傳動機構復雜，難以計算其實際效率，因此不計船舶的能量損耗。

(19)

式中：DL為輕質燃油單位價格；DH為重質燃油單位價格；ε1為輕質燃油的燃燒熱值；ε2為重質燃油的燃燒熱值；vki為船舶k在i航段的航速。

(20)

船舶于不同航段的航行時間如式(21)，分為3項：第1項為冰區航段的船舶航行時間；第2項為ECA航段的船舶航行時間；第3項為普通航段的船舶航行時間。

(21)

(22)

航段的航速限制范圍如式(23)：

(23)

式(24)表示航線只能使用單個船舶k。

(24)

4 算例的求解與分析

4.1 算法設計

建立的航速優化模型為雙目標模型，兩個目標函數均為線性函數，且存在離散化的參數條件。求解上述模型可使用線性加權法或啟發式算法。啟發式算法的求解速度較快，但求解結果通常為近似最優解。因此，針對雙目標航速優化模型，使用線性加權法進行離散化求解。求解軟件為MATLAB，求解步驟為：

1)使用線性加權法進行重要度賦權，將多目標問題轉化為單目標混合整數規劃問題。設置相應權值?1和?2，構造新目標函數Y3：

minY3=?1Y1+?2Y2

(25)

2)根據多船舶和多航段條件進行目標問題拆分，使決策變量xk失效，目標問題從而離散化為若干個一般的線性規劃子問題。

3)使用MATLAB求解上述若干子問題，得出對應的最優解。

4.2 船舶信息

選擇4艘冰級運輸船作為算例，船舶具體類型分別為：雜貨運輸船、集裝箱運輸船、原油運輸船、液化天然氣運輸船。上述船舶基于俄羅斯船級社標準的冰級分別為：Arc 4、Arc 4、Arc 7、Arc 7。船舶具體數據如表3。

表3 船舶數據

基于現有海冰密集度和海冰厚度數據，設定中度冰情(冰密集度在30%～60%和冰厚度在30 cm內)時船舶應按照預定航速的80%進行航行[4]。有冰情數據但低于上述值時，設為輕度冰情。船舶航速的80%為冰區航段中度冰情時的安全航速范圍。

4.3 航線的航段劃分

算例航線選擇天津港為起始掛靠港，鹿特丹港為終止掛靠港。并以符拉迪沃斯托克、彼得羅巴甫洛夫斯克、普羅維杰尼亞、佩韋克、季克西、迪克森以及摩爾曼斯克的港口為基準港進行航線長度的測算。設定船舶物資儲備充足，可滿足船舶完成航行的需求，船舶不掛靠基準港。設定船舶在起始港裝載貨物，在終止港卸載貨物。整條航線按照航行海域可劃分11個航段，航線內有1個起始港、1個終止港、7個途經港。

對航線的航段劃分做出以下設定：①對部分基準港區間內的非冰區海域進行整體航段劃分;②可收集到冰情數據(海冰密集度、海冰厚度)的海域跨越多基準港區間時：基于離散化的方法，對上述涉及海域的航段進行分屬基準港區間的子航段劃分。同屬一個海域的子航段的冰情參數和海水參數取值一致。為避免基準港區間的干擾，選取同海域的第一個子航段的航速作為航速分析的基礎。航線的總體航段劃分如表4。

表4 航線的總體航段劃分

4.4 優化結果分析

針對船舶的燃油成本，重質燃油的參考價格為430美元/t，輕質燃油的參考價格為530美元/t。在起始港，選取包括船舶停泊和引航(移泊)費用、貨物的港務費在內的港口使費，終止港視具體規定進行選取。依據既定條件對模型進行求解，可得出算例船舶在多航段中的不同月的最優航速，航速單位為kn。

表5為4艘船舶航速的優化結果。從表5可以看出，在非冰區航段，4艘船舶的航速波動較小，其航速均值分別為14.88、14.46、13.39、16.12 kn。而在冰區航段，相較于非冰區航段其航速有所下降，船舶在中度冰情航段航速均值分別為：12.72、12.42、11.86、13.91 kn。

表5 船舶的最優航速

圖2為4艘船舶在非冰區航段和冰區航段通航期不同月份的航速對比，冰區航段出現中度冰情顯示有降幅航速的均值。圖3為航速優化后的單航次船舶運營成本和運營時間。

圖2 船舶航速對比

由圖2可知：在非冰區航段，船舶的航速無明顯變化。而在冰區的中度冰情航段，船舶的航速有一定的降幅。依據表5的數據，在中度冰情航段，4艘船舶的航速均值相較非冰區航段的降幅最大達到14.53%、14.12%、11.41%、13.67%。由圖3可知，通航期冰情變化對船舶的運營成本以及運營時間的影響較小。因此，筆者所構建的模型能夠較好的反映船舶在北極航線通航期的實際航行狀態。

圖3 單航次運營成本和運營時間

5 結 論

筆者對北極東北航線的通航期冰情數據進行分析，考慮船舶的阻力因素，建立以船舶單航次的運營成本和運營時間最低為目標的航速優化模型。其次，使用線性加權和離散化的方法處理所建立的模型。最后，使用MATLAB軟件對模型進行求解，得出多種船舶不同航段的最優航速。從航速結果可以看出，在非冰區航段的船舶航速變化較小，而在中度冰情的條件下船舶航速有一定的降幅。相較非冰區航段，4艘算例船舶在中度冰情航段的航速均值最大降幅達到14.53%、14.12%、11.41%、13.67%。