秦 操
(上海交通大學 高新船舶與深海開發裝備協同創新中心,上海 200240)
船舶運動數學模型是船舶運動仿真與控制問題的核心[1]。運動數學模型可以提供船舶操縱性限制條件,輔助人或系統進行路徑規劃,可以嵌入控制系統中以提高控制的精度,更是船舶運動模擬器的基礎。
傳統的建立船舶運動數學模型的方法是根據MMG模型結構[2]或者Abkowitz模型結構[3],求解其中的水動力導數。求解水動力導數的方法中,經驗公式估算十分快速,但是準確性依賴于所積累的數據和船型;約束模試驗數據準確,但是成本太高,也存在尺度效應問題。系統辨識可以從包含噪聲的數據中建立數學模型,并且可以不依賴于經典船舶運動數學模型結構,因此成為目前主流的建模方法。傳統的辨識方法如最小二乘法、極大似然法等由于精度過度依賴模型和初值,逐漸被以神經網絡為代表的新型辨識方法取代。近年來,擴展卡爾曼濾波器(EKF)[4]、撐向量回歸(SVR)[5]、螢火蟲算法[6]、人工蜂群算法[7]、多目標群智能算法[8]等都被應用在了船舶系統辨識領域。以神經網絡為代表的不依賴經典船舶數學模型結構的辨識方法要想得到準確的辨識結果,需要輸入輸出數據中包含船舶完整的操縱性信息。以EKF為代表的對經典船舶數學模型中的水動力導數進行辨識的方法由于已經引入了數學模型結構這一先驗信息,往往需要更少的數據進行辨識,但辨識的精度十分依賴經典數學模型結構,并且存在參數漂移和動力相消[9]的問題。
擴展卡爾曼濾波器利用一階泰勒級數將非線性項線性化,因此適用于非線性系統,也可以進行參數辨識,但是因為忽略了高階項,在非線性較強時誤差大。無跡卡爾曼濾波器[10](Unscented Kalman Filter,UKF)是無跡變換[11](Unscented Transformation,UT)與標準卡爾曼濾波器的結合,目的是為了提高非線性系統濾波的精度和效率。UKF利用UT來進行均值和協方差的非線性傳遞,其思想是近似非線性函數的概率分布而不是近似非線性函數。UKF相比EKF有著更高的精度,并且不用計算雅可比矩陣,因此被廣泛應用于導航、目標跟蹤等領域。
本文采用精度更高的UKF,結合成本相對較低的自航試驗對目標三體船進行辨識。針對三體船自身特點對MMG模型結構進行化簡,去掉了作用不大的水動力導數,并設計了不同的控制方式進行分步辨識,以減輕動力相消的問題。通過參數的試驗測量值與辨識值的對比,驗證了本文辨識出參數的準確性;通過模型預報數據與試驗觀測數據的對比,驗證了本文建立的三體船數學模型的準確性,具有較高的工程意義和參考價值。
隨機向量X的均值和協方差分別為和Pxx,對于非線性變換Y=f(X),無跡變換是利用采樣法求取Y的均值和協方差Pyy。首先構造Sigma點集 {χi},如下式:

式中:n為系統的維數; λ =α2(n+κ)-n為尺度參數; α 決定了Sigma點圍繞的分布,通常取一個小正值,1 >α>e-4; κ為第2個尺度參數,通常取為0;表示矩陣開方運算后選擇其中第i行并轉置后得到的n維向量。
對 {χi}中的每一個點進行非線性變換Yi=f(χi),變換后的點集 {Yi}可以近似表示Y=f(X)的分布,對{Yi}進行加權處理,得到Yˉ 和Pyy,如下式:

為敘述方便,將上述無跡變換過程用函數UT表示。對于隨機向量X和一個非線性變換Y=f(X),以及另一個非線性變換Z=h(X),互協方差Pzy的計算用函數UTc表示,如下式:

對n維離散的非線性系統,狀態空間模型如下式:

式中:Xk和Xk-1分別為k時刻和k-1時刻的系統狀態向量;Zk為系統觀測向量;f為狀態轉移函數;h為系統觀測函數;Uk為輸入控制;wk為服從高斯分布的過程噪聲,其均值為0,協方差矩陣為Qk;vk為服從高斯分布的觀測噪聲,其均值為0,協方差矩陣為Rk。
UKF的基本過程與標準卡爾曼濾波器一致,在誤差協方差矩陣預報和更新時都用無跡變換來計算,遞推公式如式(5)~式(8)所示。
狀態估計預報、誤差協方差矩陣預報:

濾波增益:

狀態估計更新:

誤差協方差矩陣更新:

安裝完成的目標三體船如圖1所示,3個船體都采用了S標準型,其主要參數如表1所示。三體船的推進裝置為船體兩側的1對螺旋槳,采用差分推力來控制前進和轉向。

圖1 安裝完成的目標三體船Fig.1 Installed trimaran

表1 船體主尺度信息Tab.1 Main parameters of trimaran
基于MMG模型建立目標三體船的運動數學模型結構,三自由度操縱性方程如下式:

式中:m,Izz分別為船舶的質量和繞z軸的轉動慣性矩;mx,my,Jzz分別為船舶縱向附加質量,橫向附加質量,附加轉動慣性矩;u,v,r分別為縱向速度,橫向速度,轉首角速度;XH,YH,NH分別為船體受到的縱向力,橫向力和轉首力矩;XP,YP,NP分別為槳受到的縱向力,橫向力和轉艏力矩。
XH,YH,NH是船體受到的粘性阻力和力矩,由于三體船沒有橫向主動力,v的產生靠的是坐標變換的科里奧利力,即r和u是產生v的必要條件,所以與v和r相關的水動力導數耦合在一起,產生了動力相消的問題,無法分開辨識。另外,由于v很小,橫向信噪比高,XH,YH,NH中與v耦合的高階項作用很小。因此可以在井上模型[1]的基礎上對與v相關的高階項和耦合項進行化簡,如下式:

式 中 :Xuu,Xrr,Yv,Yr,Nv,Nr,Nrr為 水 動 力 導 數 。
假設推力系數為進速系數的一次函數關系,如下式:

式中:Jp為 進速系數;wp為伴流分數;u為船速;n為螺旋槳轉速;D為螺旋槳直徑;a0,a1常系數。
單個螺旋槳的推力T(n,u)如下式:

式中:n表示螺旋槳轉速;t為推力減額; ρ 為液體密度;Tnn,Tun為系數;Tnn也是螺旋槳系柱推力系數。
對三體船,2個螺旋槳X方向位置與船重心X位置一致,距離船中的距離為l=0.45m,則螺旋槳的推力和力矩XP,YP,NP的計算如下式:

式中:nl,nr分別表示左右螺旋槳的轉速。
上述模型結構中,需要確定的參數為mx,my,。
構建系統狀態向量X是縱向速度u,橫向速度v,轉艏角速度r,大地坐標系下正北坐標x,大地坐標系下正東坐標y,艏向角 ψ 與待辨識的參數組合而成的向量;輸入控制U是左右螺旋槳轉速組成的向量;系統觀測向量Z是X前六維變量的測量值,如下式:

狀態轉移函數和系統觀測函數如下式:

其 中dt=0.2s,u˙ ,v˙ ,r˙由 式 (9) 計 算 得 到 。 無 跡 變換參數取值為 α =0.3, β =2, κ =-12。Qk為對角矩陣,u,v,r分量分別取 0.012, 0.012, 0.052,其余分量取0;Rk為對角矩陣,各分量按順序依次取 0.05542, 0.05542,0.00452,0.22,0.22,0.12。
系統辨識過程中,可能會出現參數漂移或動力相消[9]的情況,例如在式(9)中如果m+mx,m+my,Xuu,Xrr,Tnn,Tun同時增大同樣的倍數,縱向的運動方程并不會發生變化。為避免參數漂移,mx,my,Jzz不參與辨識,而是由由周昭明回歸公式[12]計算得到,如表2所示。分步辨識法也可以有效地緩解動力相消問題,其思想是分批辨識參數,將前步辨識出的參數作為已知值再進行后面其他參數的辨識。對于目標船,由于對稱性,在進行直航運動時v,r為0,Xrr,Yv,Yr,Nv, ,Nr,Nrr不產生作用,因此可以優先辨識出Xuu,Tnn,Tun這3個參數,然后將其設為已知值,再進行剩下參數的辨識。

表2 附加質量、轉動慣量值Tab.2 Added mass and moment of inertia
為獲取目標船參數辨識所需的觀測數據序列,進行目標船靜水中的自航試驗。目標船上搭載了AHRS,GPS傳感器,可以實時地獲取位置、速度、姿態信息。目標船由岸上的遙控器進行操控,控制信號通過WIFI發送給船上的電機驅動器,以控制2個螺旋槳的轉速。首先進行多組目標船直航試驗,再進行多組目標船S型回轉試驗,將傳感器觀測數據記錄下進行離線辨識。
直航試驗和S型回轉試驗螺旋槳轉速信號如圖2所示。將多組直航試驗觀測數據辨識的得到的參數Xuu,Tnn,Tun取均值,然后將其作為已知值,再利用S型回轉試驗觀測數據對Xrr,Yv,Yr,Nv,Nr,Nrr進行辨識,結果如表3所示。
將辨識出的參數代入運動數學模型,并將模型預報的結果與實際觀測數據進行比較,直航試驗對比結果如圖3所示,S型回轉試驗對比結果如圖4所示。由對比結果可以得出,直航試驗觀測的u和模型預報值比較吻合;S型回轉試驗觀測的u和r與模型預報值十分吻合,v存在一定的誤差,原因是三體船沒有橫向主動力,橫向速度比較小,因此信噪比較大。

圖2 直航試驗和S型回轉試驗螺旋槳轉速信號Fig.2 Propeller speed on straight sailing test and S rotation test

表3 參數辨識結果Tab.3 Results of parameter identification
辨識的參數中,Tnn是可以由螺旋槳系柱推力試驗獲得精確值的,為進一步驗證參數辨識的準確性,在上海交通大學循環水槽進行了目標船螺旋槳的系柱推力標定試驗。將左右螺旋槳轉速設為一樣,轉速依次為 5,10,15,20,25,30 r/min,測量對應的系柱推力,將系柱推力對轉速平方進行線性擬合,擬合直線的斜率即為2Tnn,如圖5所示。Tnn由系統辨識得到的結果與系柱推力測量的結果對比如表4所示,結果驗證了辨識參數的準確性。
本文針對目標三體船,基于MMG模型建立了其運動數學模型結構,利用無跡卡爾曼濾波器結合自航試驗數據對模型中的參數進行辨識。為減輕動力相消帶來的影響,結合三體船特點對模型結構進行了化簡,并設計了不同的航行方式和分步辨識方法。模型預報的速度、轉首角速度與實際觀測值十分吻合,對其中可以由試驗測量的系柱推力系數,試驗測量值與辨識值誤差很小,驗證了本文辨識方法的有效性和優越性。

圖3 直航試驗模型預報的結果與實際觀測數據對比Fig.3 Comparison between model prediction result and observed data on straight sailing test

圖4 S型回轉試驗模型預報的結果與實際觀測數據對比Fig.4 Comparison between model prediction result and observed data on S rotation test

圖5 系柱推力測量結果線性擬合Fig.5 Linear regression of propeller bollard thrust

表4 試驗結果與系統辨識結果比較Tab.4 Comparison between test result and identification result