海上施工水域船舶航線規(guī)劃數(shù)學(xué)建模及求解

謝新連，劉 毅，何 傲

(大連海事大學(xué) 交通運(yùn)輸工程學(xué)院 , 遼寧 大連 116026)

0 引 言

水上施工會(huì)給通行船舶帶來許多限制和安全風(fēng)險(xiǎn)，海上施工水域往往也是交通事故的高發(fā)區(qū)域。當(dāng)施工環(huán)境比較復(fù)雜、通航限制較多時(shí)，船舶安全航線的設(shè)計(jì)就顯得格外重要，也具有難度。特殊水域上船舶航線自動(dòng)規(guī)劃的發(fā)展將很大程度決定現(xiàn)代無人駕駛船舶技術(shù)的發(fā)展水平。早在現(xiàn)代智能化無人船出現(xiàn)之前，路徑規(guī)劃已經(jīng)有大量的相關(guān)研究[1]

傳統(tǒng)的船舶航線規(guī)劃是指在靜態(tài)和動(dòng)態(tài)障礙物并存的環(huán)境中，尋找一條從已知起點(diǎn)到終點(diǎn)的滿足一定評(píng)價(jià)標(biāo)準(zhǔn)的航行路徑，使得船舶在航行過程中能夠安全可靠地避開所有的障礙物[2]。常見的航線規(guī)劃方法有動(dòng)態(tài)規(guī)劃法、最速下降法、最優(yōu)控制法、啟發(fā)式搜索法、神經(jīng)網(wǎng)絡(luò)法、模擬退火法、遺傳算法等。由于每種航線規(guī)劃算法都或多或少存在一些缺陷，所以實(shí)際應(yīng)用的方法大多數(shù)是基本算法的改進(jìn)算法，或者融合了數(shù)種基本航線規(guī)劃算法的算法[3]。孟祥杜[4]做了無人船路徑規(guī)劃算法研究，針對(duì)無人船工作水域具有不規(guī)則邊界且常常存在障礙物的特點(diǎn)，也介紹了幾種常用的建立環(huán)境模型的方法。

環(huán)境建模也是智能駕駛技術(shù)研究最基本、最關(guān)鍵的環(huán)節(jié)之一，在環(huán)境建模方面的研究，王飛等[5-7]采用柵格法進(jìn)行環(huán)境建模，提出了基于Maklink圖和遺傳算法的改航路徑規(guī)劃方法，并且利用遺傳算法對(duì)Maklink圖中的初始飛行器航線進(jìn)行了改航規(guī)劃；李明富等[8]在Maklink圖的基礎(chǔ)上采用變參數(shù)螢火蟲算法對(duì)路徑進(jìn)行優(yōu)化，實(shí)驗(yàn)表明改進(jìn)后的螢火蟲算法能夠較好的解決離散路徑規(guī)劃問題；孟祥杜[4]做了無人船路徑規(guī)劃算法研究，針對(duì)無人船工作水域具有不規(guī)則邊界且常常存在障礙物的特點(diǎn)，也介紹了幾種常用的建立環(huán)境模型的方法。

由以上的文獻(xiàn)可以發(fā)現(xiàn)，目前研究的問題主要集中在普通或正常環(huán)境下各種交通工具的路徑的規(guī)劃。然而在特殊環(huán)境下的路徑規(guī)劃研究則較為罕見。而且船舶的路徑規(guī)劃和航線設(shè)計(jì)有其自身的特殊性，船舶在水中的操作性能和陸上的交通工具或機(jī)器人也有較大差異，因此不能生搬硬套其它交通工具的路徑規(guī)劃的方法。船舶在海上航行時(shí)事故頻發(fā)的區(qū)域就是海上的施工區(qū)域，海上施工種類繁多，對(duì)其附近船舶的影響也各不相同，為了保證船舶能夠順利的通過施工水域，有必要對(duì)施工水域的船舶航線規(guī)劃進(jìn)行相關(guān)的研究。

1 問題描述及環(huán)境模型的構(gòu)建

船舶在海上航行時(shí)可能會(huì)經(jīng)過某些施工水域，并且在有的水域中同時(shí)包含多個(gè)水面和水下施工地點(diǎn)，當(dāng)在同一片施工水域中包含多個(gè)施工點(diǎn)時(shí)，船舶在其中的航行將變得極為復(fù)雜和危險(xiǎn)，假設(shè)船舶在某施工水域中的環(huán)境示意如圖1。

圖1中，S為航線規(guī)劃的起點(diǎn)，T為航線規(guī)劃的終點(diǎn)，其中的陰影部分表示海上的施工區(qū)域。

為了方便施工水域的航線規(guī)劃數(shù)學(xué)模型的建立，首先利用Maklink圖論的方法來處理施工水域的危險(xiǎn)區(qū)域，通過在施工水域中構(gòu)造Maklink鏈接圖來構(gòu)造施工水域中的可航行網(wǎng)絡(luò)圖，構(gòu)造Maklink航線網(wǎng)絡(luò)圖的具體步驟可以參考文獻(xiàn)[5, 8]。由于船舶始終是處于海平面高度，因此筆者建立數(shù)學(xué)模型時(shí)不考慮海拔高度的影響，將環(huán)境模型的構(gòu)建簡(jiǎn)化為二維平面。另外對(duì)施工水域中所有的不可航行區(qū)域都采用凸多邊形或多個(gè)凸多邊形的組合來表示。在模型化表示不可航行區(qū)域的基礎(chǔ)上進(jìn)一步構(gòu)造可航行的自由空間，海上施工水域的可航行自由空間是由鏈接線構(gòu)成的凸多邊形。在圖1的環(huán)境模型中構(gòu)造的可航行自由空間如圖2。

圖1 施工水域環(huán)境示意Fig. 1 Schematic diagram for construction waters

圖2 可航行自由空間Fig. 2 Nautical free space map

構(gòu)造可自由航行空間圖時(shí)主要需要遵循以下兩個(gè)基本準(zhǔn)則：①所構(gòu)造的鏈接線不能夠穿越施工水域中的危險(xiǎn)區(qū)域；②鏈接線的一端是施工水域的危險(xiǎn)區(qū)域所構(gòu)成的凸多邊形的頂點(diǎn)，另一端是另一個(gè)施工危險(xiǎn)區(qū)域構(gòu)成的凸多邊形的頂點(diǎn)或者是施工水域的自由空間的邊界點(diǎn)。

在施工水域的自由空間中，將所有鏈接線的中點(diǎn)兩兩相連所形成的圖就是需要構(gòu)造的可航行區(qū)域的Maklink航線網(wǎng)絡(luò)圖，如圖3。

圖3 施工水域的MaklinkFig. 3 Maklink diagram for construction waters

2 海上施工水域航線規(guī)劃的數(shù)學(xué)模型

將圖3中構(gòu)造的Maklink航線圖，進(jìn)一步表述為可求解的數(shù)學(xué)模型。假設(shè)由起點(diǎn)S到終點(diǎn)T的一條在施工水域中的航線是由點(diǎn)序列(S,p1,p2,…,pk,T)組成，數(shù)學(xué)建模時(shí)其目標(biāo)方程就是點(diǎn)序列所構(gòu)成的航線的距離最短。在點(diǎn)序列中，船舶在除了起點(diǎn)S和終點(diǎn)T以外的所有中間點(diǎn)處都有可能轉(zhuǎn)向，為此我們需要考慮轉(zhuǎn)向角對(duì)船舶航行的約束，船舶在一定的航行速度下，船舶在pi處的轉(zhuǎn)向角我們用θpi來表示，它不能超過其轉(zhuǎn)向角的極限角度，用φmax表示船舶在設(shè)計(jì)航速下的轉(zhuǎn)向極限角度。另一方面船舶的轉(zhuǎn)彎半徑也是航線設(shè)計(jì)中必須考慮的一個(gè)重要因素，船舶轉(zhuǎn)彎的轉(zhuǎn)彎半徑和船舶所處的航行速度有很大的關(guān)系，根據(jù)參考文獻(xiàn)[9-10]中指出船舶在沒有拖輪的輔助的情況下在10 kn的設(shè)計(jì)航速下船舶的最小轉(zhuǎn)彎半徑可按表1 的設(shè)計(jì)參數(shù)選取。

表1 航道轉(zhuǎn)彎半徑設(shè)計(jì)指標(biāo)Table 1 Channelturning radius design index

根據(jù)十一五科技攻關(guān)項(xiàng)目“沿海深水港航道選線及設(shè)計(jì)主要研究”中[11]，對(duì)10萬噸級(jí)集裝箱船、15萬噸集裝箱船和10萬噸級(jí)LNG船舶在不同轉(zhuǎn)向和航速情況下的轉(zhuǎn)彎進(jìn)行的模擬研究，得出了航道的轉(zhuǎn)彎半徑R也可按照公式(1)估算：

(1)

式中：R為航道設(shè)計(jì)的轉(zhuǎn)彎半徑,m；Vt為航道的設(shè)計(jì)行駛速度，kn;L為航道的通行船長,m；θ為船舶在某點(diǎn)需要轉(zhuǎn)向的角度,(°)。

研究發(fā)現(xiàn)航道的轉(zhuǎn)向角達(dá)到60°時(shí)，船舶的操縱出現(xiàn)操縱困難的現(xiàn)象，船舶操縱人員不易掌握船位、轉(zhuǎn)向角速率及在航道上的姿態(tài)；當(dāng)轉(zhuǎn)向角達(dá)到90°時(shí)，操縱極為困難，對(duì)操作人員的壓力極大，特別是在風(fēng)、浪、流等自然因素的共同作用下，過彎操作的失誤率明顯增加[12]。因此在進(jìn)行航線規(guī)劃時(shí)，得到的航線必須滿足船舶在設(shè)計(jì)航速下在各個(gè)轉(zhuǎn)向點(diǎn)處的轉(zhuǎn)向角度小于等于90°。

考慮這些因素后，建立船舶在海上施工水域中航線規(guī)劃的數(shù)學(xué)模型如式(2)～式(6)。

(2)

s.t.

?l(pi，pi+1)∩Sj=?,i=1,2,3,…,n,j=1,2,3,…,m

(3)

?θPj∈(1,…,k)≤φmax

(4)

RPj∈(1,.…,k)≥λL

(5)

D≥Dmin

(6)

模型的解釋：公式(2)表示目標(biāo)方程，其中(pi,pi+1)表示由pi和pi+1兩點(diǎn)構(gòu)成的航段，d(pi,pi+1)表示航段(pi,pi+1)間的距離，式中k表示路徑終點(diǎn)T前的最后一個(gè)節(jié)點(diǎn)。目標(biāo)方程由三部分組成，第一部分表示組成航線的中間航段的長度總和，后面兩部分表示連接起點(diǎn)和終點(diǎn)的航段長度。公式(3)表示構(gòu)成航線的任意航段與施工水域的危險(xiǎn)區(qū)域的交集為空集，其中l(wèi)(pi,pi+1)表示由pi和pi+1兩點(diǎn)連接形成的航段的點(diǎn)集；Sj表示施工水域中危險(xiǎn)區(qū)域的點(diǎn)集。公式(4)表示船舶在第pj節(jié)點(diǎn)處的轉(zhuǎn)向角必須小于航線設(shè)計(jì)航速下船舶的最大安全轉(zhuǎn)向角φmax，由前面的研究結(jié)果可知一般情況下φmax≤90°。公式(5)表示的是在pj節(jié)點(diǎn)處船舶進(jìn)行轉(zhuǎn)向時(shí)，航線必須滿足船舶最小轉(zhuǎn)彎半徑的要求，其中L表示船舶的長度,m；λ表示航道在不同設(shè)計(jì)航速下所取的最小轉(zhuǎn)彎半徑與船長的比值。公式(6)表示規(guī)劃航線與危險(xiǎn)區(qū)域的最小距離的約束；D表示航線某處與危險(xiǎn)區(qū)域的最小距離；Dmin表示船舶與危險(xiǎn)區(qū)域的最小安全距離,m。

3 Dijkstra算法求解航線規(guī)劃數(shù)學(xué)模型

對(duì)建立的數(shù)學(xué)模型的求解，采用兩階段優(yōu)化求解的方法：第一階段采用Dijkstra算法，在建立的Maklink可航行航線網(wǎng)上搜尋經(jīng)過網(wǎng)絡(luò)節(jié)點(diǎn)的最短航線；第二階段再利用蟻群算法對(duì)得到的初始解優(yōu)化，最終得到最優(yōu)的航線。

Dijkstra算法是用于解決非負(fù)權(quán)值網(wǎng)絡(luò)中尋找指定兩點(diǎn)間最短距離的最好的算法之一。Dijkstra算法的基本思想是采用標(biāo)號(hào)法，給每個(gè)節(jié)點(diǎn)預(yù)先設(shè)置兩個(gè)標(biāo)簽，其中一個(gè)用于標(biāo)記路長，另一個(gè)是用于標(biāo)記從起點(diǎn)到終點(diǎn)路徑的最后一條弧的起始點(diǎn)號(hào)。網(wǎng)絡(luò)節(jié)點(diǎn)的標(biāo)號(hào)分為兩類，一類是永久標(biāo)號(hào)，另一類是臨時(shí)標(biāo)號(hào)。當(dāng)?shù)恋趉步時(shí)，獲得永久標(biāo)號(hào)的點(diǎn)意味著已經(jīng)找到起點(diǎn)到該點(diǎn)的最短路徑，將獲得永久標(biāo)號(hào)的點(diǎn)放在Sk集合中，獲得永久標(biāo)號(hào)的點(diǎn)的路長d值和路徑標(biāo)號(hào)ξ值不再修改。

在利用Dijkstra算法求解前，需要得到Maklink圖中網(wǎng)路節(jié)點(diǎn)各點(diǎn)間的鄰接距離矩陣，航線網(wǎng)絡(luò)圖的鄰接距離矩陣為一個(gè)19×19階的矩陣，即：

在鄰接距離矩陣中，∞表示兩個(gè)節(jié)點(diǎn)不連通，其距離為無窮大，其中的數(shù)字25則表示p1和p3是連通的，并且兩點(diǎn)間的距離為25 km。利用Dijkstra算法求解得到的結(jié)果的總距離為221.573 0 km，其路徑依次為S-p8-p9-p2-p1-p12-p13-T。航線規(guī)劃的結(jié)果如圖4。

圖4 Dijkstra算法求解路徑Fig. 4 Diagram for Dijkstra algorithm solution

4 蟻群算法對(duì)初始航線的優(yōu)化

在構(gòu)造海上施工水域中的航線圖時(shí)，選擇的是可自由航行鏈接線的中點(diǎn)來構(gòu)造Maklink航線網(wǎng)絡(luò)的節(jié)點(diǎn)，但是根據(jù)模型約束條件式(6)的實(shí)際要求，只要與施工水域中的危險(xiǎn)區(qū)域的距離大于Dmin即可自由航行，因此最優(yōu)的航線不一定經(jīng)過航線的網(wǎng)絡(luò)的中點(diǎn)。即節(jié)點(diǎn)pi的位置是可以在自由空間的鏈接線上進(jìn)一步優(yōu)化的，因此下面就使用蟻群算法對(duì)初始航線進(jìn)行進(jìn)一步的優(yōu)化。

4.1 工作空間離散化處理

(7)

式中：pi(hi)為鏈接線上的任意一點(diǎn)；hi為一個(gè)0到1的隨機(jī)數(shù)，用它來表示當(dāng)前點(diǎn)所處在兩端點(diǎn)間的比例參數(shù)；N為某鏈路所經(jīng)過的節(jié)點(diǎn)的個(gè)數(shù)[13]。

在得到Dijkstra算法的解的基礎(chǔ)上，只要給出一組(h1,h2，…,hN)參數(shù)，就可以得到一組由起點(diǎn)到終點(diǎn)的新解。在采用蟻群算法優(yōu)化施工水域的船舶航線前，需要離散化工作空間，對(duì)鏈接線的劃分筆者采用固定長度劃分法。為了滿足模型約束條件距離危險(xiǎn)區(qū)域的最小距離為Dmin的要求，采用的劃分長度ξ取Dmin，每條鏈接線L將劃分成的分?jǐn)?shù)為

(8)

式中：Int()函數(shù)表示取整函數(shù)，當(dāng)Int(Li/Dmin)為奇數(shù)時(shí)，為了保證鏈接線的中點(diǎn)也是一個(gè)等分點(diǎn)，將劃分的等分?jǐn)?shù)加1，經(jīng)過離散處理后那么從鏈接線Li-1到其相鄰的鏈接線Li將會(huì)有Ni+1條路徑可以選擇。

4.2 螞蟻路徑的搜索和信息素的更新

每一只螞蟻在完成從起點(diǎn)到終點(diǎn)的路徑搜索后都會(huì)產(chǎn)生一組與之相對(duì)應(yīng)的路徑參數(shù)集合(h1，h2，…，hk)，在移動(dòng)的過程中，一只螞蟻在鏈接線Li上時(shí)，選擇下一個(gè)鏈接線Li+1上的節(jié)點(diǎn)j的方法為，依次計(jì)算當(dāng)前鏈接線節(jié)點(diǎn)i到下一條鏈接線節(jié)點(diǎn)j的選擇概率Pij，然后根據(jù)選擇概率Pij采用輪盤賭法找出下一個(gè)節(jié)點(diǎn)j，Pij的計(jì)算公式為[14]

(9)

式中：ηi,j為啟發(fā)值；τi,j為信息素，信息素的濃度越高，選擇概率Pij越大。當(dāng)螞蟻在移動(dòng)過程中選擇了某個(gè)節(jié)點(diǎn)，則在該節(jié)點(diǎn)處需要釋放信息素，則對(duì)該節(jié)點(diǎn)的信息素需要進(jìn)行更新的公式為

τi,j=(1-ρ1)×τi,j+τ0

(10)

式中：τ0為信息素的初始值；ρ1為[0,1]區(qū)間的可調(diào)參數(shù)，文中ρ1取0.1。當(dāng)所有的螞蟻都完成一次搜索，保存最短路徑的長度，同時(shí)再次更新每個(gè)點(diǎn)上的信息素。此次更新稱為信息素的揮發(fā)更新，信息素?fù)]發(fā)的公式為

τi,j=τi,j(1-ρ2)

(11)

式中：ρ2為信息素?fù)]發(fā)參數(shù)，文中ρ2的值取0.003。

4.3 蟻群算法對(duì)航線轉(zhuǎn)向角的優(yōu)化

圖5 轉(zhuǎn)向角優(yōu)化示意Fig. 5 Diagram for steering angle optimization

τi,j=λij×[(1-ρ1)×τi, j+τ0]

(12)

(13)

公式(12)為修正后的信息素更新公式，式中λij為新引進(jìn)的信息素更新公式的修正系數(shù)；公式(13)為信息素的修正系數(shù)的計(jì)算公式，其中αij為當(dāng)前位置和備選位置的連續(xù)與起點(diǎn)和終點(diǎn)位置連線的夾角，為了消除角度正負(fù)方向的影響，公式(13)中對(duì)夾角 取絕對(duì)值。另一方面因?yàn)榇暗淖畲蟀踩D(zhuǎn)向角度為90°，所以取π/2為參考基準(zhǔn)，當(dāng)αij絕對(duì)值大于90°時(shí)修正系數(shù)λij的值將小于1，將會(huì)使路徑所含有信息素濃度τi,j減少，從而達(dá)到迭代過程選中轉(zhuǎn)向角大于90°的概率不斷減小的目標(biāo)。同時(shí)當(dāng)夾角αij的絕對(duì)值越小時(shí)λij的值也將越大，當(dāng)αij的絕對(duì)值取最小值0時(shí)，λij取得最大值2。

4.4 蟻群算法優(yōu)化求解最終結(jié)果

在利用蟻群算法進(jìn)行優(yōu)化求解時(shí)，每次迭代產(chǎn)生10只螞蟻來搜索最優(yōu)路徑，在每一次迭代結(jié)束后都記錄下這10只螞蟻尋找到的最優(yōu)路徑和最優(yōu)路徑的長度，總共進(jìn)行了500次迭代。最終通過蟻群算法求解的航線總長度為205.43 km，比Dijkstra算法求解的路徑總長度221.573 0 km減少了16.139 3 km，較初始解縮短了7.87%，最終的航線規(guī)劃結(jié)果如圖6。

圖6 蟻群算法最終結(jié)果Fig. 6 Final result for ant colony algorithm

圖6中虛線就是由蟻群算法求得的航線，它不僅在總距離上和初始解相比減少了7.87%，而且航線在轉(zhuǎn)向次數(shù)和轉(zhuǎn)向角上也得到了優(yōu)化。為了滿足數(shù)學(xué)模型約束條件公式(6)的約束，觀察圖6可以發(fā)現(xiàn)蟻群算法優(yōu)化得到的船舶航線與施工水域的危險(xiǎn)水域都保持了一定安全距離，特別是圖中的A、B兩個(gè)地方，船舶航線與施工水域的危險(xiǎn)區(qū)域的距離已經(jīng)是算法設(shè)置的極限安全距離Dmin，本算法實(shí)驗(yàn)中Dmin取值為0.2 km。

為了觀察蟻群算法在搜索最優(yōu)航線時(shí)的收斂過程，圖7是蟻群算法優(yōu)化初始解的收斂過程。

圖7 蟻群算法收斂過程Fig. 7 Ant colony algorithm convergence process

觀察圖7可以發(fā)現(xiàn)蟻群算法在迭代的初期航線總長度迅速收斂，當(dāng)?shù)?00次以后算法收斂趨于平穩(wěn)，說明此時(shí)已經(jīng)得到了較為優(yōu)秀的滿意解，但是為了減小陷入局部最優(yōu)解的概率，提升算法搜索全局最優(yōu)解的能力,可以發(fā)現(xiàn)蟻群算法在迭代的中后期出現(xiàn)了幾次明顯的波動(dòng)。這是因?yàn)樵诶孟伻核惴ㄇ蠼饽Ｐ蜁r(shí)，螞蟻種群始終會(huì)以一個(gè)較小的概率去探索未知的領(lǐng)域，以便于發(fā)現(xiàn)全局最優(yōu)解。

5 結(jié) 論

1)為了解決船舶在施工水域中的船舶航線規(guī)劃問題，以圖論的相關(guān)理論為基礎(chǔ)，在海上施工水域中構(gòu)建起了可自由航行網(wǎng)絡(luò)圖，利用Dijkstra算法求得了航線的初始可行解。

2)在得到初始可行解的基礎(chǔ)上，利用蟻群算法對(duì)初始解進(jìn)行了智能優(yōu)化，使最終得到的海上施工水域的船舶航線得到了較大幅度的優(yōu)化。

3)通過兩階段的優(yōu)化求解求得了既滿足航線盡可能短又能夠避開海上施工水域危險(xiǎn)區(qū)域的航線。