朱昱雯，肖健梅，王錫淮

上海海事大學(xué)物流工程學(xué)院，上海201306

0 引 言

隨著航運業(yè)的不斷發(fā)展，人們對海洋資源的依賴程度越來越高，傳統(tǒng)的錨泊定位系統(tǒng)已經(jīng)無法滿足船舶的深海作業(yè)需求。目前，很多船舶都配置了動力定位系統(tǒng)，可以精準(zhǔn)定位船舶，從而控制船舶在期望方位上航行。作為動力定位系統(tǒng)的重要執(zhí)行模塊，推力分配模塊可以根據(jù)動力定位控制器給出的期望推力和力矩指令，對船舶推進器的推力和力矩進行合理分配，同時抵消環(huán)境擾動，從而使船舶在期望的艏向和位置上航行。此外，通過合理分配各個推進器的推力和方向角，可以解決推進器過度旋轉(zhuǎn)導(dǎo)致的磨損問題和推力損失引起的船舶電力系統(tǒng)功率波動問題。

在保持期望艏向和位置的前提下，推力分配模塊還需保證推進器產(chǎn)生的功率是最小推力組合，因此，推力分配問題可以視為具有約束條件的單目標(biāo)優(yōu)化問題。目前，大多采用序列二次規(guī)劃方法來解決推力分配問題，但該方法較為繁瑣且不利于編程實現(xiàn)，同時受約束條件的影響較大，故其實際可應(yīng)用的范圍非常有限［1］。徐云成等［1］利用了序列二次規(guī)劃來解決推力分配問題，結(jié)果表明：推進器推力和方向角的變化較為平緩，且誤差較??；但該研究未結(jié)合推力分配模塊與動力定位控制器進行實時仿真，即仿真系統(tǒng)沒有根據(jù)上一時刻船舶的位置和航速來改變推力力矩，而是僅對已知的推力力矩指令進行推力分配測試，故其應(yīng)用范圍較窄。而智能算法可以很好地解決推力分配問題，其過程也比序列二次規(guī)劃方法更簡明。目前，適用于推力分配問題的智能算法主要有遺傳算法、人工魚群算法和粒子群算法等［2］。Zhao 等［3］將遺傳算法應(yīng)用于推力分配問題，結(jié)果表明，遺傳算法可以響應(yīng)動力定位控制器根據(jù)上一時刻船舶位置而產(chǎn)生的期望指令，但其推進器推力和方向角的變化存在震蕩，所以增加了推進器的磨損。劉明等［4］結(jié)合偽逆法與混沌粒子群算法來解決推力分配問題，結(jié)果表明，推進器推力和方向角的前期變化存在波動，但后期變化較為平緩。雖然智能算法在處理速度、簡便程度、可實施性方面均優(yōu)于序列二次規(guī)劃方法，但在算法后期易陷入局部最優(yōu)問題，從而導(dǎo)致期望推力與實際推力的力矩指令之間存在一定誤差。

因此，本文擬采用基于成功歷史記錄的自適應(yīng)參數(shù)差分進化（success-history based parameter adaptation for differential evolution，SHADE）改進算法來解決船舶的推力分配問題。差分進化算法在前期收斂速度較快，后期因種群多樣性減少而過早停滯，傳統(tǒng)SHADE 算法已針對這一問題進行了改進。然而，對于存在環(huán)境擾動的推力分配問題，SHADE 算法的誤差較大且響應(yīng)速度較慢，故本文擬基于傳統(tǒng)SHADE 算法進行改進：根據(jù)當(dāng)前迭代次數(shù)，在動態(tài)范圍內(nèi)選擇較優(yōu)的個體進行變異操作，以提高后期的收斂速度；根據(jù)算法的迭代程度，對用以產(chǎn)生交叉率CR 的正態(tài)分布方差和用以動態(tài)調(diào)整位置參數(shù)的權(quán)重c 進行自適應(yīng)優(yōu)化，從而提高算法的穩(wěn)定性。最后，通過與傳統(tǒng)SHADE算法進行對比，證明改進SHADE 算法的可行性和有效性。

1 動力定位數(shù)學(xué)模型

典型的動力定位控制系統(tǒng)結(jié)構(gòu)如圖1 所示，本文將采用PID 控制器作為動力定位系統(tǒng)的控制器，用以控制船舶縱蕩、橫蕩和艏搖3 個自由度的運動。首先，PID 控制器通過對比實際位置η、實際航速v 和參考位置ηd，從而產(chǎn)生期望推力力矩指令τd；然后，基于改進SHADE 算法，即可得到推進器總功率最小時各個推進器的推力f 和方位角α 的最優(yōu)組合；最后，將該組合分配給各個推進器以產(chǎn)生實際合力τ ，根據(jù)動力學(xué)模型估算船舶的當(dāng)前位置和速度信息并反饋至動力定位控制器，從而實時控制船舶定位。

圖1 控制系統(tǒng)結(jié)構(gòu)圖Fig.1 Structure diagram of control system

PID 控制器的控制律為［5］

式中：RT(ψ)為旋轉(zhuǎn)矩陣的轉(zhuǎn)置，其中ψ 為船舶的艏搖；Kp，Ki，Kd均為PID 控制器的參數(shù)；s 為拉氏變換中的復(fù)變量。

1.1 船舶數(shù)學(xué)模型

船舶的動力定位系統(tǒng)通常僅考慮縱蕩、橫蕩和艏搖3 個自由度，這類船舶在北東坐標(biāo)系NED下的位置和船體坐標(biāo)系xyz 下的速度關(guān)系為［6］

1.2 推力分配數(shù)學(xué)模型

1.2.1 推進器數(shù)學(xué)模型

其中

1.2.2 推力分配目標(biāo)函數(shù)

推力分配問題不僅需考慮每個推進器產(chǎn)生的推力和方向角是否滿足動力定位控制器的期望推力和力矩指令，還需考慮推進器推力和方向角的變化率。在每次更新期望推力和力矩指令之后，推力分配模塊均需快速響應(yīng)每個推進器需產(chǎn)生的推力和方向角，但推進器在短時間內(nèi)的推力增加量有限。同時，全回轉(zhuǎn)推進器的過度旋轉(zhuǎn)將導(dǎo)致推進器磨損和功耗增加。此外，當(dāng)多個推進器同時工作時，相鄰?fù)七M器將互相影響，一個推進器螺旋槳產(chǎn)生的尾流會進入相鄰?fù)七M器的螺旋槳工作區(qū)域，對相鄰?fù)七M器產(chǎn)生水動力干擾，從而導(dǎo)致推進器的推力損失，因此需限制全回轉(zhuǎn)推進器的方向角，從而避免推力損失并減緩?fù)七M器方向角的變化速率［8］。

推力分配的目標(biāo)函數(shù)為

約束條件為

2 改進SHADE 算法

SHADE 算法是一種利用前期成功試驗向量的控制參數(shù)來調(diào)整當(dāng)前迭代控制參數(shù)的改進差分進化算法，其計算思路為：首先，利用存儲空間記錄多組可生成較優(yōu)控制參數(shù)的均值參數(shù)；然后，利用這些均值參數(shù)來生成每次迭代過程的控制參數(shù)，并在每次迭代結(jié)束之后利用生成成功試驗向量的控制參數(shù)來調(diào)整存儲空間的均值參數(shù)。這種自適應(yīng)調(diào)整控制參數(shù)的方法不僅加快了差分進化算法的收斂速度，還提高了魯棒性。

SHADE 算法在解決無約束問題時的表現(xiàn)較好，但處理復(fù)雜問題時將過于依賴隨機生成的初始種群。例如，對于推力分配問題，如果動力定位控制器產(chǎn)生的期望指令過小，則SHADE 算法的魯棒性將有所下降，且易出現(xiàn)過早停滯的現(xiàn)象，從而導(dǎo)致偏差較大。針對這些問題，本文將對傳統(tǒng)SHADE 算法進行改進。

2.1 算法原理

2.1.1 變異操作

式中：g 為當(dāng)前迭代次數(shù)，對應(yīng)第g 代種群；G 為最大迭代次數(shù)。

基于式（9），SHADE 算法的后期迭代過程可以減少在非精英區(qū)域搜索的范圍，從而加快后期的收斂速度并避免陷入局部最優(yōu)。本文將參數(shù)p的最小值設(shè)定為round(0.025×NP) ，當(dāng)種群規(guī)模NP=200 時，參數(shù)p 隨迭代次數(shù)的變化曲線如圖2 所示。

圖2 當(dāng)G=500，NP=200 時，參數(shù)p 隨迭代次數(shù)的變化曲線Fig.2 The curve of parameter p with the number of iterations at G=500，NP=200

2.1.2 交叉操作

由于推力分配優(yōu)化問題含有2 個變量：推進器推力f 和方向角α，因此需對2 個變量均進行交叉操作。本文將采用二項式交叉方式對目標(biāo)向量和變異向量進行交叉操作，具體如下：

通過對每個變量向量中的至少一維分量進行交叉操作，即可增加推進器推力和方向角的組合方式，用以保持種群的多樣性，從而提高算法的前期收斂速度。

2.1.3 選擇操作

分別計算試驗向量Ui，g和目標(biāo)向量Xi，g的適應(yīng)度值，從而選取可以產(chǎn)生更優(yōu)適應(yīng)度值的個體進入下一代，其選擇操作如下：

2.1.4 參數(shù)自適應(yīng)調(diào)整

由于差分進化算法依賴于控制參數(shù)（例如，縮放因子F 和交叉率CR）的設(shè)置，并且動力定位控制器每一次更新期望推力和力矩指令之后的工況均有所不同，因此，如果人為設(shè)置固定的控制參數(shù)，則將導(dǎo)致船舶航行過程中某一次期望推力及力矩指令與實際推力及力矩指令之間的誤差過大。雖然在SHADE 算法中已采用自適應(yīng)的方式對控制參數(shù)進行了調(diào)整，但用于生成控制參數(shù)分布的某些參數(shù)仍依賴于人為設(shè)定，因此，為了避免該誤差導(dǎo)致船舶偏離航線，本文將對這部分人為設(shè)定參數(shù)進行自適應(yīng)調(diào)整。

本文將分別采用柯西分布和正態(tài)分布隨機生成每個個體的縮放因子Fi，g和交叉率，即

式中：Fi，g∈［0.2，0.8］；CRi，g∈［0，1］；randc 為柯西分布函數(shù)；randn 為正態(tài)分布函數(shù)；μF為柯西分布的位置參數(shù)；μCR為正態(tài)分布的均值；σ2CR為正態(tài)分布的方差。

如果Fi，g和CRi，g超出閾值范圍，則將采用截斷的方式對其進行限制［11］。如果f（Ui，g）＜f（Xi，g），即表示生成了成功試驗向量Ui，g，其對應(yīng)的Fi，g和CRi，g將被分別存儲至空間SF和SCR中，用以在每一代結(jié)束后動態(tài)調(diào)整μF和μCR。同時，調(diào)整之后的μF和μCR將存儲至空間W 中，用以在每一個目標(biāo)向量變異之前從W 中隨機選取μF和μCR，其中空間W 的大小為

將空間W 中μF和μCR的初始值分別設(shè)為0.1和0.5，并在算法初始化時設(shè)定一個常數(shù)k，用以控制μF和μCR在空間W 中的存儲位置；同時，利用μF和μCR更新空間W 中第k 個位置的數(shù)值［12］，當(dāng)k ＞W(wǎng)size時，令k=1，重新開始計數(shù)。

針對柯西分布中的位置參數(shù)μF，本文基于SHADE 算法對其更新方式增加了一個權(quán)重參數(shù)c，用以控制μF的更新速率。雖然空間W 的存在可以預(yù)防μF更新過快所導(dǎo)致的算法陷入局部最優(yōu)的問題，但空間W 的大小并不容易控制，不同的期望推力指令需要不同大小的空間W，才能達到期望范圍內(nèi)的偏差。當(dāng)船舶需要穩(wěn)定在某一點時，過大的空間W 將導(dǎo)致算法迭代緩慢；但在船舶航行過程中，過小的空間W 將容易導(dǎo)致算法陷入局部最優(yōu)。因此，本文將通過新增的權(quán)重參數(shù)c來進行平衡，從而減小空間W 對船舶航行過程中各種工況下推力分配的影響。位置參數(shù)μF的動態(tài)調(diào)整方法具體如下：

對于控制μF更新的權(quán)重參數(shù)c，本文將根據(jù)迭代次數(shù)進行動態(tài)調(diào)整，即在算法前期更新速度較慢以保證算法的穩(wěn)定性，而在后期則將提高更新速度以滿足精度需求。權(quán)重c 的更新方式為

正態(tài)分布雖與柯西分布類似，但其存在均值和方差，故可較好地反映交叉率CR 的變化情況，因此，生成成功試驗向量的交叉率CR 可以引導(dǎo)均值和方差變化，從而自適應(yīng)調(diào)整交叉率CR。但是，如果在均值μCR的動態(tài)調(diào)整過程中增加了權(quán)重c，反而將降低最優(yōu)均值μCR的尋找速度。因此，本文將在正態(tài)分布的參數(shù)調(diào)整過程取消權(quán)重c，則均值μCR和方差σ2CR的動態(tài)調(diào)整方式如下［13］：

2.2 改進SHADE 算法的流程及步驟

在基本差分進化算法最優(yōu)值的基礎(chǔ)上，改進SHADE 算法增加了存儲歷史參數(shù)的空間和利用歷史參數(shù)調(diào)整控制參數(shù)的步驟，并在變異操作中增加了隨機選擇的精英個體。利用改進SHADE算法計算最優(yōu)推進器推力和方向角的流程如圖3所示。

圖3 改進SHADE 算法程序框圖Fig.3 Block diagram of improved SHADE algorithm program

采用改進SHADE 算法對推力分配問題進行優(yōu)化的具體步驟如下：

1）根據(jù)動力定位控制器給出的期望推力和力矩指令τd以及推進器在前一時刻的推力和方向角，在相應(yīng)的可行解空間中產(chǎn)生種群規(guī)模為NP的初始種群（包含每個推進器的推力和方向角）。根據(jù)式（6）計算松弛變量s0，然后計算每個個體的初始適應(yīng)度值并進行排序。

2）進行參數(shù)初始化，將空間SF和SCR清零。

3）根據(jù)式（9）計算精英個體數(shù)p。從空間W中隨機選取μF和μCR，根據(jù)式（12）和式（13）計算Fi，g和CRi，g。

4）根據(jù)式（8）的變異規(guī)則，選擇3 個不同個體，對當(dāng)前目標(biāo)向量按式（8）進行變異操作，生成變異個體Vi，g。

5）根據(jù)式（10）進行交叉操作，生成試驗向量Ui，g，并采用截斷的方 式對Ui，g的數(shù)值 進行邊界處理。

6）根據(jù)式（11）進行選擇操作，生成下一代的目標(biāo)向量Xi，g+1，并將失敗的目標(biāo)向量Xi，g存入集合A 中，將生成成功試驗向量的Fi和CRi分別存儲至空間SF和SCR。

7）判斷第g 代種群的個體是否全部完成變異、交叉和選擇操作，如果已完成，執(zhí)行步驟8；否則，返回步驟3。

8）判斷是否滿足終止條件，即迭代次數(shù)為最大迭代次數(shù)或s0達到最小要求誤差。如果滿足，輸出實際推力和力矩指令；否則，按式（15）～（21）更新空間W 中的μF，μCR，σ2CR，并返回步驟2。

3 仿真分析

為了驗證改進的SHADE 算法解決推力分配問題的可行性，本文將在相同的仿真環(huán)境和參數(shù)設(shè)置下，分別應(yīng)用改進SHADE 算法和傳統(tǒng)SHADE 算法對推進器的推力進行優(yōu)化分配。本文將采用文獻［14］中帶有動力定位系統(tǒng)的船舶進行仿真驗證，其中整個系統(tǒng)采用Simulink 進行仿真，而數(shù)值優(yōu)化部分則采用Matlab Function 模塊對算法進行編程實現(xiàn)。該船長76.2 m，寬18.8 m，質(zhì)量為4.591×106kg，配置了4 個推進器（2 個全回轉(zhuǎn)推進器和2 個隧道推進器），本文將每個推進器產(chǎn)生的最大推力設(shè)為船舶干重的1/60。推進器的具體布置如圖4 所示，其相關(guān)參數(shù)如表1 所示。

圖4 推進器的詳細(xì)布置圖Fig.4 Specific layout position of thruster

表1 推進器參數(shù)Table 1 Parameters of thruster

設(shè)定船舶從初始位置η0=［0，0，0］T航行至期望位置ηd=［60，40，30］T，本文將采用PID 控制器實現(xiàn)動力定位模式仿真，該控制器所產(chǎn)生的期望推力和力矩指令如圖5 所示。

在仿真過程中，本文為船舶施加了一個環(huán)境擾動b，該環(huán)境擾動由白噪聲經(jīng)過二階近似海浪模型而得，如圖6 所示。

圖5 動力定位控制器產(chǎn)生的期望推力和力矩指令Fig.5 Expected thrust and torque instruction generated by dynamic positioning controller

圖6 環(huán)境擾動圖Fig.6 Environmental disturbance

由圖5 和圖6 可見，經(jīng)仿真驗證，改進SHADE算法可以對推力進行合理分配，保證船舶由初始位置航行至期望位置，并在期望位置保持穩(wěn)定的環(huán)境擾動抑制能力。其航行路線和動力定位響應(yīng)如圖7 和圖8 所示。由圖8 可以看出，當(dāng)船舶到達期望位置后，在環(huán)境擾動下基本可以穩(wěn)定在原位，并且在3 個自由度上幾乎無波動。

圖9 所示為改進SHADE 算法和傳統(tǒng)SHADE算法的作用效果對比圖，在相同的仿真環(huán)境下，雖然2 種算法求解推力分配問題時的船舶動力定位響應(yīng)無明顯區(qū)別，但傳統(tǒng)SHADE 算法所計算的各個推進器推力均較大，這將導(dǎo)致功率消耗和方向角變化也較大。

圖7 船舶的航行路線Fig.7 Navigation route of ship

圖8 船舶動力定位響應(yīng)Fig.8 Dynamic positioning response of ship

由圖9（a）可以看出，基于改進SHADE 算法的推進器推力的變化曲線較為平緩，沒有大幅度振蕩的情況，并且未出現(xiàn)2 個相鄰?fù)七M器所產(chǎn)生的推力相互影響而導(dǎo)致的推力損失問題。由圖9（b）可以看出，當(dāng)船舶到達期望位置后，較小的推力即可使船舶在環(huán)境擾動下基本穩(wěn)定在原位，其中隧道推進器的推力值在-0.1～0.1 kN 范圍內(nèi)，全回轉(zhuǎn)推進器則為-0.5～0.5 kN。由圖9（c）中可以看出，基于傳統(tǒng)SHADE 算法的推進器推力的變化波動較大，并且在船舶航行過程中出現(xiàn)了相鄰?fù)七M器推力方向相反所導(dǎo)致的推力損失問題。由圖9（d）可以看出，當(dāng)船舶到達期望位置后，推力變化曲線的振蕩較大，且每個推進器均需產(chǎn)生較大的推力才能使船舶在環(huán)境擾動下基本穩(wěn)定在原位。

圖9 基于2 種算法的推進器推力變化對比圖Fig.9 Thrust change comparison diagram of thruster based on two algorithms

圖10 所示為基于2 種算法的全回轉(zhuǎn)推進器方向角的變化對比圖。由圖10 可以看出，基于改進SHADE 算法的推進器方向角的變化曲線較為平緩，沒有出現(xiàn)突變情況，并且當(dāng)船舶到達期望位置后，全回轉(zhuǎn)推進器的方向角趨于穩(wěn)定。而基于傳統(tǒng)SHADE 算法的推進器方向角在時刻變化，當(dāng)船舶到達期望位置后也無法抑制環(huán)境擾動的影響，故其角度無法穩(wěn)定。由此可見，改進SHADE 算法可以滿足船舶的基本航行要求和推力、方向角變化的約束要求，從而降低全回轉(zhuǎn)推進器的磨損。

圖11 所示為基于2 種算法的期望推力與實際推力力矩指令的誤差對比圖。由圖11 可以看出，由動力定位控制器產(chǎn)生的期望推力力矩指令與改進SHADE 算法得到的實際推力力矩指令之間的誤差非常小，其中：仿真前期的縱蕩和橫蕩誤差在±5×10-5kN 范圍內(nèi)，艏搖力矩誤差在±3×10-5kN·m 范圍內(nèi)；仿真后期 的縱蕩和橫蕩誤 差 在±1×10-7kN 范 圍 內(nèi)，艏 搖 力 矩 誤 差在±1×10-8kN·m 范圍內(nèi)。然而，傳統(tǒng)SHADE 算法的誤差則較大，且其在船舶到達期望位置后的誤差高于航行階段。因此，利用改進SHADE 算法來解決推力分配問題時，可以滿足較高的精度要求，并可將整個過程的誤差控制在較小范圍內(nèi)。

圖10 基于2 種算法的全回轉(zhuǎn)推進器方向角變化對比圖Fig.10 Direction change comparison diagram of full-rotating propulsion based on two algorithms

圖11 基于2 種算法的期望推力與實際推力力矩指令的誤差對比圖Fig.11 Error comparison diagram of torque instruction of expected thrust and actual thrust based on two algorithms

圖12 所示為基于2 種算法的推進器總功率變化對比圖。由圖12 可以看出，基于改進SHADE算法的推進器總功率變化較為平緩，未出現(xiàn)大幅振蕩，其在環(huán)境擾動下可以使船舶穩(wěn)定在期望位置，同時保證較小的功率消耗。而基于傳統(tǒng)SHADE 算法的推進器方向角和推力變化均存在振蕩現(xiàn)象，無法在較短時間內(nèi)達到最優(yōu)，故其推進器總功率的振幅和波動較大。因此，基于改進SHADE 算法對船舶進行推力分配，可以較好地抑制環(huán)境擾動對船舶航行的影響。

圖12 基于2 種算法的推進器總功率變化對比圖Fig.2 Total power change comparison diagram of thruster based on two algorithms

4 結(jié) 語

本文提出了一種運用改進SHADE 算法解決船舶推力分配問題的方法，仿真驗證對比結(jié)果表明，在整個航行過程中，推進器的推力和方向角變化較為平緩，可以降低推進器磨損和相互作用影響；推進系統(tǒng)的功率消耗變化幾乎無波動，期望推力和實際推力的力矩指令誤差可以控制在較小的范圍內(nèi)，其中縱蕩和橫蕩誤差在±5×10-5kN 范圍內(nèi)，艏搖力矩誤差在±3×10-5kN·m 范圍內(nèi)。

由此可見，該算法可以較好地解決船舶推力分配問題，且在環(huán)境擾動下可以保持航向穩(wěn)定；當(dāng)船舶到達期望位置后，較小的推力即可使船舶基本穩(wěn)定在原位。此外，該成果可為如何利用推力分配維持船舶電力系統(tǒng)的穩(wěn)定性研究奠定基礎(chǔ)。