吳利紅， 張愛鋒， 李一平， 封錫盛， 王詩文

(1. 大連海事大學(xué) 船舶與海洋工程學(xué)院，遼寧 大連 116026；2. 中國科學(xué)院沈陽自動化研究所 機器人學(xué)國家重點實驗室，遼寧 沈陽 110016)

預(yù)報水下機器人(autonomous underwater vehicle，AUV)定常航行的巡航速度和非定常運動的操縱性能不僅是評估其推進系統(tǒng)性能的重要指標，還有助于提高AUV航程估計，能量估算和減少定位誤差和AUV丟失的風險。采用類物理的數(shù)值方法進行海洋載體操縱預(yù)報是數(shù)值模擬的目標。此方法應(yīng)用真實模型模擬螺旋槳高速旋轉(zhuǎn)，能反映船槳舵內(nèi)在作用機理，再現(xiàn)流場特征。但它具有模型復(fù)雜，網(wǎng)格尺度相差大，載體各部件相對運動復(fù)雜，高頻與低頻運動混合的模擬難點[1]。自20世紀90年代美國海軍提出這個概念以來，直至2010年左右陸續(xù)有研究成果。如Chase等[1]、Pankajakshan等[2]、Poremba等[3]對潛艇自航運動進行了數(shù)值模擬；Mofidi等[4]對集裝箱船KCS帶全附體槳舵的ZigZag操縱運動的模擬；Carrica等[5]對水面戰(zhàn)艦的隨浪橫甩運動進行自航模擬；沈志榮[6]對水面船舶操縱運動進行數(shù)值模擬。以上研究均采用動態(tài)重疊網(wǎng)格法，在邊界移動時無需網(wǎng)格更新，但須實時確定移動邊界與背景網(wǎng)格的重疊區(qū)域，這是數(shù)值誤差的源頭[7]。而移動網(wǎng)格法在邊界移動時網(wǎng)格也隨著運動，運動更符合拉格朗日描述法，但邊界移動常導(dǎo)致復(fù)雜的網(wǎng)格運動，從而引起網(wǎng)格畸變和網(wǎng)格質(zhì)量降低，影響數(shù)值計算精度。受此限制，移動網(wǎng)格法主要應(yīng)用在二維物體的運動：圓球的振動，魚類的擺鰭運動[8]和三維物體的有限區(qū)間或簡單的運動，如翼的偏轉(zhuǎn)[9]，預(yù)定義AUV直航對接運動[10]。而本文則采用多塊動態(tài)混合網(wǎng)格方法[11]解決移動網(wǎng)格的網(wǎng)格畸變和數(shù)值精度差的問題[12]，實現(xiàn)AUV自航試驗的數(shù)值模擬，獲得了AUV巡航速度預(yù)報，再現(xiàn)了AUV自航過程非定常運動的航行特性。

1 水下機器人自航模型建立

AUV為改進型REMUS AUV，L=1.7 m，D=0.191m[13]。具有橢球型艏、圓柱形平行中體和圓錐型艉。艉部有“十字”型舵翼以及一個MAU4-40，直徑d=0.156 m的螺旋槳。AUV艏、舯和艉3段以及舵翼采用結(jié)構(gòu)網(wǎng)格繪制，螺旋槳采用非結(jié)構(gòu)網(wǎng)格離散。圖1為AUV帶舵翼和螺旋槳的網(wǎng)格圖。大地坐標系為oxyz,ox為AUV對稱軸，向艏為正(向前)；oy垂直向上為正；oz指向右舷。

圖1 帶槳和舵的AUV網(wǎng)格模型Fig.1 Grid model of AUV appended propeller and rudders

圖2為AUV帶槳舵流域在對稱面上網(wǎng)格分區(qū)，分區(qū)劃分是根據(jù)AUV的運動趨勢確定的。流域包含 6個分區(qū)。分區(qū)Ⅰ包含AUV和舵翼表面的四邊形網(wǎng)格以及近場邊界層區(qū)域的六面體結(jié)構(gòu)網(wǎng)格；分區(qū)Ⅱ包括螺旋槳表面的三角形網(wǎng)格以及螺旋槳擾動子區(qū)域的四面體非結(jié)構(gòu)網(wǎng)格；分區(qū)Ⅲ包含AUV、螺旋槳和舵翼擾動區(qū)域的六面體結(jié)構(gòu)網(wǎng)格；分區(qū)Ⅳ為AUV從艏部延伸到遠場的航行區(qū)域網(wǎng)格，靠近艏部加密，遠場較為稀疏，由六面體結(jié)構(gòu)網(wǎng)格構(gòu)成；分區(qū)Ⅴ為從AUV艉部延伸到遠場的結(jié)構(gòu)網(wǎng)格區(qū)域，靠近AUV艉部為加密區(qū)，遠場為稀疏區(qū)，也由六面體結(jié)構(gòu)網(wǎng)格構(gòu)成；外圍非結(jié)構(gòu)網(wǎng)格區(qū)域為分區(qū)Ⅵ。

圖2 動網(wǎng)格分區(qū)Fig.2 Dynamic meshing zones

在AUV自航推進數(shù)值模擬過程中，網(wǎng)格運動區(qū)域限制在圖2的AUV直航窄條帶區(qū)域，這有利于局部擾動區(qū)域梯度的捕捉。其中分區(qū)Ⅱ既有隨螺旋槳的旋轉(zhuǎn)運動，又有隨AUV的直航運動；分區(qū)Ⅰ和Ⅲ僅隨AUV作直航運動；分區(qū)Ⅳ和Ⅴ分別為網(wǎng)格壓縮區(qū)域和拉伸區(qū)域，在分區(qū)Ⅳ的最前端為靜止壁面，采用動態(tài)層方法對區(qū)域Ⅳ的最前端網(wǎng)格進行壓縮合并，因此最前端為網(wǎng)格消亡區(qū)域；分區(qū)Ⅴ的最后端為動態(tài)層網(wǎng)格的拉伸區(qū)域，為網(wǎng)格增長區(qū)域。外圍區(qū)域Ⅵ為網(wǎng)格變形區(qū)。各個分區(qū)之間的連接采用非一致連接。這種多塊結(jié)構(gòu)網(wǎng)格和非結(jié)構(gòu)網(wǎng)格混合的網(wǎng)格構(gòu)成方法，以及采用動態(tài)層網(wǎng)格更新方法，是借鑒了動態(tài)重疊網(wǎng)格的子區(qū)域網(wǎng)格運動的思想，從網(wǎng)格構(gòu)建質(zhì)量、網(wǎng)格更新質(zhì)量和網(wǎng)格更新速度等方面著手，提高數(shù)值求解的速度和計算精度。

AUV自航推進過程數(shù)值模擬中，AUV和螺旋槳之間力和速度的傳遞可以采用用戶自定義函數(shù)UDF(user defined function)進行編寫，嵌入到程序中[14]：首先，AUV靜止，螺旋槳以恒定轉(zhuǎn)速旋轉(zhuǎn)，產(chǎn)生推力，將此推力保存同時將推力傳遞給AUV；AUV讀取螺旋槳推力，同時實時求解URANS(unsteady reynolds averaged navier-stokes)方程，獲得AUV的阻力。AUV在推力和阻力的聯(lián)合作用下，通過在6-DOF方程加載推力和實時阻力，求解載體的加速度，并采用龍格庫塔積分獲得AUV實時速度V；此速度傳遞給螺旋槳，使AUV拖曳螺旋槳以V航速運動；隨著AUV航行，AUV的尾流場更新，影響螺旋槳的進流場，螺旋槳在新的進速和轉(zhuǎn)速下產(chǎn)生新的推力，重新記錄此推力并傳遞給AUV；AUV在阻力作用下，獲得新的航速，如此重復(fù)，直到AUV的推力T與其航行阻力R接近平衡，速度曲線趨于穩(wěn)定，加速度接近為零?；?qū)Ρ群剿倥c多參考系坐標[15](multiple frames of reference, MFR)的定常自航點航速，兩者誤差在5%內(nèi)，也認為達到自航點。記錄此時的推力、阻力、速度歷史數(shù)據(jù)，保存和分析。

2 水下機器人自航數(shù)值結(jié)果與分析

2.1 受力和速度

AUV 以不同轉(zhuǎn)速從靜止自航推進，最終會趨于勻速狀態(tài)。分別假定AUV以300、450、600 r/min 3種轉(zhuǎn)速航行，每種轉(zhuǎn)速自航一次，迭代步長為螺旋槳旋轉(zhuǎn)1°對應(yīng)的時間，每個外循環(huán)包含內(nèi)循環(huán)20次，收斂精度為10-4。如圖3分別為轉(zhuǎn)速為300、450和600 r/min的自航螺旋槳推力和AUV阻力的歷史變化曲線。推力和阻力曲線在數(shù)值求解中有點振蕩，這是動網(wǎng)格數(shù)值求解引起的不穩(wěn)定現(xiàn)象。阻力部分振蕩主要是由于載體低速巡航采用了和螺旋槳高速旋轉(zhuǎn)的小時間步引起的。隨著速度增加，推力減小，阻力增大，最終推力和阻力平衡，速度會趨于一個恒定的速度，這個速度即為該轉(zhuǎn)速下的自航點。3種轉(zhuǎn)速可得3次航速歷史變化，如圖4所示，可得3個不同的自航點。轉(zhuǎn)速越大，達到勻速的時間越短。600 r/min，需要約9 s, 可達到勻速1.55 m/s；450 r/min，需要約10.53 s, 可達到勻速1.15 m/s； 300 r/min，需要約14.39 s, 可達到勻速0.74 m/s。同時，轉(zhuǎn)速越大，起始加速度越大。

圖3 不同轉(zhuǎn)速下AUV自航阻力和推力變化Fig.3 The resistance (R) and Thrust (T) of AUV self-propulsion at different rotation speeds

圖4 不同轉(zhuǎn)速下AUV自航過程的航速變化Fig.4 Velocity history (V) in AUV self-propulsion at different rotation speeds (n)

2.2 數(shù)值驗證

為了驗證數(shù)值模擬結(jié)果的準確性，將AUV自航的數(shù)值結(jié)果與水池試驗(試驗在中國科學(xué)院沈陽自動化所的深水水池：L×B×D=25 m×15 m×9 m)結(jié)果進行對比。圖5記錄了螺旋槳轉(zhuǎn)速為300 r/min的AUV速度歷史數(shù)據(jù)。包括5次試驗結(jié)果，以離散數(shù)據(jù)點表示：Exp1、Exp2、Exp3、Exp4 和Exp5；數(shù)值結(jié)果以實線表示。數(shù)值模擬不受試驗水池邊界的限制，因此數(shù)值模擬的時間較長。數(shù)值結(jié)果和試驗結(jié)果吻合良好。此外，為了驗證不同轉(zhuǎn)速下AUV自航點的可靠性，將本文的結(jié)果與多參考坐標系MFR(multi-frame of reference)的結(jié)果進行對比驗證，如圖6所示，獲得的自航點速度略低于MFR，誤差在5%以內(nèi)。在一段轉(zhuǎn)速范圍內(nèi)，航速和轉(zhuǎn)速成正比。這是因為相近的一段轉(zhuǎn)速范圍內(nèi)，推力系數(shù)接近常數(shù)，阻力系數(shù)接近常數(shù)。螺旋槳推力與轉(zhuǎn)速成平方，勻速后推力與阻力相等，阻力與速度成平方，則自航點速度和轉(zhuǎn)速成正比。這個結(jié)論有助于更快找到不同轉(zhuǎn)速下的自航點。由轉(zhuǎn)速和航速的關(guān)系，根據(jù)設(shè)計航速1.5 m/s的要求，對圖6進行插值，可得對應(yīng)自航點的轉(zhuǎn)速為570 r/min。

圖5 AUV自航速度變化的試驗和數(shù)值結(jié)果 (n=300 r/min)Fig.5 The velocity history from experimental and numerical results in AUV self-propulsion at n=300 r/min

圖6 自航點對應(yīng)的轉(zhuǎn)速Fig.6 Approaching velocities vs. rotating speeds

2.3 流場分析

圖7為轉(zhuǎn)速n為300、 450和600 r/min時達到勻速直航的AUV對稱面上的速度云圖(注：遠場速度為零)?？梢姡D(zhuǎn)速越大，螺旋槳的艉跡越明顯，梢渦強度越大，向后推水的水速度越大，反作用力越大，則推力越大。且轉(zhuǎn)速越大，AUV艏部的速度擾動流場范圍越大，表明AUV達到自航點對應(yīng)的速度越大。3種勻速狀態(tài)時，AUV艏部最大截面處都產(chǎn)生了艏渦，且隨著速度增加，艏渦范圍加大，表明AUV艏部粘壓阻力增加。隨著AUV速度增加，AUV艉部截面變化區(qū)域開始產(chǎn)生艉渦(n=600 r/min)，表明艉部速降增加，艉部粘壓阻力增加。同時，3種航速時，AUV表面都有較明顯的邊界層，從艏至艉逐步增厚。

圖7 AUV不同轉(zhuǎn)速對應(yīng)的自航點速度云圖Fig.7 Velocity contours of AUV at different self-propulsion points with varying rotation speeds

圖8為設(shè)計航速1.5 m/s對應(yīng)轉(zhuǎn)速570 r/min 的AUV自航過程4個典型時刻(t為0.1，3,4.9和6.9 s)對稱面上AUV的速度分布(右圖)和螺旋槳艉跡局部放大圖(左圖)。相應(yīng)的壓力云圖如圖9所示?？梢姡斵D(zhuǎn)速一定時，AUV速度較低時，螺旋槳艉跡強，表明螺旋槳此時的推力較大(這與螺旋槳的敞水性能曲線一致，和圖3的推力曲線一致)。螺旋槳葉梢曳出的梢渦出現(xiàn)重疊；隨著AUV航速增加，梢渦強度降低，表明螺旋槳推力減少。但是梢渦分離清晰可見，這是由AUV航行時的艉跡反向作用引起的。

圖8 AUV自航對稱面上不同時刻的速度云圖(n=570 r/min)Fig.8 Velocity contours of AUV self-propulsion on the symmetry plane at different times with n=570 r/min

圖9 AUV自航對稱面上不同時刻的壓力云圖(n=570 r/min)Fig.9 Pressure contours of AUV self-propulsion at the symmetry plane at different times with n=570 r/min

隨著AUV航速的增加， 螺旋槳產(chǎn)生的推力在變小，這點在圖9的左側(cè)壓力云圖中清晰可見。低速時，螺旋槳葉面的高壓值很大，葉背的低壓值也較大，導(dǎo)致兩者的壓差較大，因此推力較大；隨著AUV速度增加，葉面高壓值減小，葉背低壓值也減小，兩者的壓差變小，則推力也在減少。但是，隨著螺旋槳持續(xù)旋轉(zhuǎn)，只要AUV推力大于阻力，AUV就獲得加速，隨著AUV速度增加，AUV艏部速度繞流場范圍加大，艏部駐點的壓強增加，導(dǎo)致AUV艏艉壓差阻力增加；同時艏肩部速度增加，導(dǎo)致摩擦阻力增加，AUV總阻力也在增加，直至與逐漸減小的推力平衡，達到勻速直航，即達到自航點。隨著螺旋槳推進AUV自航，AUV艉跡延長。螺旋槳艉跡梢渦強度隨AUV速度增加弱化。同時，在螺旋槳轂帽附近有一股根渦，和AUV航行方向一致，流速大小隨著AUV速度增加而增加。

3 結(jié)論

1) 水下機器人不同轉(zhuǎn)速對應(yīng)的最終巡航速度不同。300、450和600 r/min分別對應(yīng)的巡航速度為0.74、1.15和1.55 m/s。AUV的巡航速度1.5 m/s對應(yīng)的轉(zhuǎn)速是570 r/min。

2) AUV自航穩(wěn)定所需時間與AUV轉(zhuǎn)速大小有關(guān)。轉(zhuǎn)速越大，加速度越大，穩(wěn)定需要時間越短。300、450和600 r/min穩(wěn)定分別需要14.39、10.53和 9 s。

3) 在螺旋槳轉(zhuǎn)速變化較小的范圍內(nèi)，AUV的自航點速度和轉(zhuǎn)速成正比關(guān)系。

4) AUV自航時，螺旋槳旋轉(zhuǎn)運動曳出梢渦和根渦。梢渦強度隨著AUV航速增加而降低，運動方向與AUV航行方向相反。根渦大小隨AUV航速增加而增加，方向與AUV航行方向一致。

5) 螺旋槳推力的變化起源于其葉面和葉背壓差變化。AUV以恒定轉(zhuǎn)速自航，推力隨速度增加而減小。