水陸兩棲機器人水下動力學建模及操縱性研究

陳 龍，童保成，祖斯羽，張 豪，胡 橋，2，*

（1．西安交通大學 機械工程學院，陜西 西安 710049；2．西安交通大學 陜西省智能機器人重點實驗室，陜西 西安 710049）

0 引言

海洋一直以來都是國家高質(zhì)量發(fā)展的戰(zhàn)略要地，建設海洋強國是實現(xiàn)中華民族偉大復興的重大戰(zhàn)略任務。研究一種既能實現(xiàn)海洋與陸地環(huán)境探測，又具有高效推進能力、環(huán)境適應性強的兩棲機器人，對海洋資源探測和近海區(qū)域警戒具有重要的現(xiàn)實意義[1-2]。針對現(xiàn)有水陸兩棲機器人環(huán)境適應性差和穩(wěn)定性弱等問題，本文將仿生推進機構與傳統(tǒng)推進機構結合設計出一種高性能組合推進型水陸兩棲機器人。機器人的水下操縱性影響其工作狀態(tài)，因此研究機器人的水動力系數(shù)以及水下操縱性有重要意義。

通過CFD方法建立機器人水下動力學模型的方法已通過標模試驗數(shù)據(jù)的驗證，得到了國內(nèi)外學者的廣泛應用。2021年，F(xiàn)RANCESCHI等人利用開源CFD代碼OpenFoam求解了船體本地和附加質(zhì)量的水動力系數(shù)，完成了對雙軸艦船的操縱性質(zhì)量評估[3]，研究表明，應用CFD流體分析方法可以得到較為準確的結果。2019年，GO等人提出了一種利用計算流體動力學技術求解拖魚水下航行器水動力導數(shù)的新方法[4]，并進行了不同速度下拖魚在3種運動場景的六自由度仿真，最后通過虛擬仿真結果與理論對比分析，驗證了方法的有效性[3]。2017年，意大利DUBBIOSO等人利用非定常CFD求解器詳細分析了潛艇的機動特性，針對十字形和X形舵進行水平面的三自由度運動分析，并通過實驗進行驗證，研究結果表明X形舵具備更加優(yōu)異的轉彎能力[5]。上海交通大學萬教授團隊基于OpenFOAM仿真平臺，開發(fā)了一款船舶自航行求解器naoe-FOAM-SJTU[6]，實現(xiàn)了船舶航行的直接數(shù)值模擬。結合上述研究現(xiàn)況分析可知，基于CFD的間接求解方法計算效率和精度高，且易于參數(shù)化，便于進行動力學模型的修正和機器人優(yōu)化設計。

本文以一種利用波動鰭推進的輪–鰭–槳機器人為研究對象，通過CFD技術針對機器人水下多種運動模態(tài)進行數(shù)值求解，分析得到機器人不同運動模態(tài)下的水動力特性，求解水動力參數(shù)，建立機器人的動力學模型。構建水下運動仿真平臺，通過機器人直航、回轉等運動模態(tài)研究其水下操縱性，并驗證動力學模型的準確性和有效性[6]，為進一步研究多模態(tài)水陸兩棲機器人操縱特性奠定了理論與技術基礎。

1 機器人坐標系建立及數(shù)值求解方法

1.1 坐標系建立

如下圖1所示，本文采用了國際水池會議（ITTC）推薦的體系，建立水陸兩棲機器人的坐標系：

圖1 水陸兩棲機器人坐標系示意Fig.1 Schematic diagram of the coordinate system for amphibious robots

大地坐標系Oexeyeze–固定坐標系。該坐標系Oe為空間中任一點，且Oeze正向指向地心，Oexe和Oeye位于水平面內(nèi)且互相垂直。

隨體坐標系Obxbybzb–運動坐標系。原點與機器人浮心重合，Obxb軸沿機器人前進方向為正，Obyb軸沿機器人右側方向為正，Obzb軸垂直于Obxbyb平面指向下為正。

鰭面坐標系Ofxfyfzf–鰭面坐標系。各坐標軸正向與波動鰭動力學建模分析中坐標系相同。

1.2 計算模型

本文以機器人三維模型為基礎構建仿真計算模型，計算模型將其內(nèi)部填充為實體，保留整體的外形輪廓；同時針對局部細節(jié)進行簡化處理，提高計算速率。最終構建的機器人水下仿真模型如下圖2所示，其尺寸為1 150 mm×607 mm×280 mm（長×寬×高）。

1.3 網(wǎng)格劃分

根據(jù)機器人尺寸劃分外流場和密度盒內(nèi)外2部分計算域。其中：外流場尺寸為8L×6L×6L，密度盒尺寸為1.3L×0.8L×0.5L，L為機器人體長。計算模型放置于密度盒中，距離流體域速度入口的距離為3L，尾部距離壓力出口4L，如圖3所示。

圖3 流場范圍及邊界條件Fig.3 Flow field range and boundary conditions

本文采用ICEM軟件進行網(wǎng)格劃分。外流場計算域尺寸較大，且無計算模型，故采用結構六面體網(wǎng)格進行劃分；由于計算模型在部分工況中非靜止，故密度盒內(nèi)網(wǎng)格劃分采用適應性強的四面體非結構網(wǎng)格，該網(wǎng)格支持smooth和remesh共2種動網(wǎng)格節(jié)點更新方法，在機器人計算模型表面進行加密，保證計算精度；2個計算域間采用interface邊界條件連接，并在接觸面處進行加密，如圖4所示。

圖4 計算域網(wǎng)格劃分Fig.4 Computational domain grid partitioning

1.4 仿真參數(shù)設定

機器人運動的外流場通常呈現(xiàn)大雷諾數(shù)湍流，因此采用可準確描述流場演化規(guī)律的N-S方程進行求解，笛卡爾坐標系下RANS方程如下：

經(jīng)無關性驗證，確定密度盒內(nèi)流場網(wǎng)格數(shù)量為2640 000和時間步長為1/40 T時滿足精度要求，同時計算效率較高。

機器人仿真計算參數(shù)設定如表1所示。

表1 機器人運動仿真參數(shù)Table 1 Robot motion simulation parameters

2 水動力系數(shù)求解與動力學建模

根據(jù)受力分析可知，若要建立機器人水下動力學模型，首先需要求解機器人在水下運動過程中受到的水動力[6]。求解機器人水動力系數(shù)成為建立動力學模型的關鍵。出于研究周期和成本考慮，常采用CFD技術求解機器人不同運動狀態(tài)下的水動力系數(shù)。

2.1 純縱蕩運動

進行純縱蕩運動仿真可求解機器人沿x軸方向速度u的慣性水動力系數(shù)Xu′˙，其運動方程如式（2）所示，根據(jù)運動方程編寫流體仿真UDF文件，實現(xiàn)計算模型在流場內(nèi)的運動。

式中：U為流場流速，m/s；A為機器人正弦運動的幅值，設置A=40 mm；f為機器人純縱蕩運動頻率，取f=0.3 Hz、f=0.4 Hz、f=0.5 Hz共3組工況進行求解。

純縱蕩運動中機器人所受力的表達式為

計算收斂之后，機器人純縱蕩運動中受到的x軸方向的阻力X如圖5所示，由于數(shù)值仿真前期數(shù)據(jù)波動較大，可靠性低，因此本文選擇10～20 s仿真時間內(nèi)的數(shù)據(jù)進行擬合，此時計算已穩(wěn)定。

圖5 純縱蕩運動中x軸受力曲線圖Fig.5 X-axis force curves in pure surge motion

通過曲線擬合，得到公式系數(shù)如表2所示。

由于機器人在水下運動過程中很少進行反向運動，若采用純縱蕩運動求解粘性水動力系數(shù)會造成較大誤差，因此該工況僅用于求解慣性水動力系數(shù)，而粘性水動力系數(shù)由穩(wěn)態(tài)直航運動進行求解。公式系數(shù)與水動力系數(shù)的關系式為

通過建立公式系數(shù)Xa和-1/2ρL3Aω2之間的線性關系求解縱蕩運動的水動力系數(shù)，曲線的斜率即為水動力系數(shù)，而從小幅度振蕩的線性理論可知，由多工況繪制的曲線應該具備過原點特征[9]，因此建立的過原點曲線如圖6所示。

圖6 過原點擬合X軸加速度系數(shù)曲線Fig.6 A curve of fitting X-axis acceleration coefficient through the origin

由圖6可以看出，公式系數(shù)Xa和-1/2ρL3Aω2之間整體呈線性關系，且數(shù)據(jù)擬合效果較好，證明了動態(tài)運動仿真求解機器人水動力系數(shù)的有效性，求得機器人純縱蕩運動無因次化水動力系數(shù)：

2.2 水動力系數(shù)求解匯總

以純縱蕩運動求解過程為例，可以求解與無因次化處理其它水動力系數(shù)。具體可參考文獻[10]。機器人水動力系數(shù)計算結果匯總如表3所示。

表3 水動力系數(shù)計算結果匯總Table 3 Summary of hydrodynamic coefficient calculation results

2.3 五自由度動力學及運動學建模

結合機器人受力分析及水動力系數(shù)求解結果，可建立其動力學模型，忽略機器人橫搖小擾動的影響[7]，將機器人運動方程簡化為五自由度動力學方程：

同理，可將機器人運動學方程進一步簡化，得到其平動方程：

轉動方程為

基于數(shù)值仿真結果建立了機器人五自由度運動學和動力學模型，為機器人水下操縱性分析研究奠定了基礎。

3 機器人水下操縱性運動特性研究

機器人水下操縱性運動特性分析是基于機器人運動仿真方程，通過設定初始運動參數(shù)及機器人推力，以數(shù)值積分的方式求解機器人下一時刻的運動狀態(tài)，可以直觀地看到水動力對機器人運動的影響?？梢灶A報機器人運動，驗證運動模型的準確性；同時可以總結機器人運動規(guī)律，為機器人的優(yōu)化和研制提供指導。

3.1 運動仿真數(shù)值計算方法

為了保證仿真結果的準確性，機器人的仿真方程計算過程中，本文采用精度高、易收斂的四階經(jīng)典龍格庫塔法進行積分。計算方法如下：

其中：

3.2 機器人運動仿真結果及分析

基于操縱性仿真運動模型，本文針對機器人直航運動及水平面回轉運動開展仿真研究，分析機器人運動狀態(tài)，并對水動力系數(shù)和動力學模型的準確性進行驗證。

3.2.1 機器人直航運動仿真

初始狀態(tài)下，機器人的初始運動速度和姿態(tài)角均為0，首先設置推進器推力為10 N，接著每隔30 N進行一次運動仿真，求解不同推力下機器人航速和最大推力下的極限速度，結果如圖7所示。

圖7 不同推力下機器人航速曲線Fig.7 Robot speed curves under different thrusts

由圖7中機器人運動航速曲線可知：

1）當推進器推力P＝10N時，機器人航速u＝0.36m/s；當P＝40N時，u＝0.72m/s；P＝70N ，u＝0.95m/s；而推進器最大合力106 N左右，機器人最大航速約為1.17 m/s。

2）機器人在螺旋槳推進力的作用下，由靜止開始加速，且加速時間隨著推進力的增大而縮短，最終速度趨于穩(wěn)定。

3）當水動力與推力平衡時，機器人進入勻速行駛階段。隨著推力增大，機器人航速越大，阻力逐漸增大，而航速增長速率逐漸減小，符合機器人水下直航的一般規(guī)律。

3.2.2 機器人水平面回轉運動

設置機器人左右推進器推力P＝10N，搖艏力矩N＝-3N?m ，機器人由靜止開始運動，仿真時間設置為70 s。其運動軌跡、運動過程中機器人在隨體坐標系中的縱向速度u、橫向速度v、垂向速度w的變化曲線如圖8所示。

圖8 機器人水平面回轉運動Fig.8 Robot horizontal rotation motion

由圖8可知：

1）機器人在推力P和搖艏力矩N的作用下，由靜止開始加速，在10 s之后開始進入勻速回轉階段。

2）當機器人開始勻速回轉時，此時縱向速度u＝0.32m/s ，橫向速度v＝0.035m/s，回轉半徑R約為3.095 m。

3）機器人在推進器推力和搖艏力矩的作用下可實現(xiàn)水下回轉運動，符合水下機器人水平回轉的一般規(guī)律。

結果表明，機器人水下運動基本符合一般水下機器人的運動規(guī)律，驗證了水動力系數(shù)求解的準確性和五自由度運動方程的有效性，該方程可通過一系列仿真預報機器人運動狀態(tài)，總結運動規(guī)律，為機器人閉環(huán)運動控制研究奠定了基礎。

4 機器人水下操縱性運動實驗

為了進一步研究機器人水下運動操縱性，本文開展水下操縱性實驗，與數(shù)值仿真結果對比分析，從而驗證機器人水動力系數(shù)和動力學模型的準確性。

4.1 直航運動實驗

為了驗證機器人直航運動仿真結果的準確性，本實驗通過改變推進器推力進行機器人不同航速測試。機器人在推進器作用下從靜止開始加速，加速完成后通過10 m標準航段，監(jiān)測所需時長，從而計算出機器人運動速度，實驗結果如下圖9所示。

圖9 不同推力下機器人水下航行速度Fig.9 Underwater navigation speed of robots under different thrusts

將實驗結果與仿真數(shù)據(jù)對比分析，其結果如下表4所示。

表4 不同推力下機器人航行速度Table 4 Robot navigation speed under different thrusts

由實驗結果可知，機器人水下航行速度隨著推力增大而不斷增大，而阻力也顯著提高，因此機器人速度增長速率逐漸減小，符合水下機器人的一般運動規(guī)律；且速度誤差小于6%，可見實驗結果與仿真數(shù)據(jù)吻合較好，初步驗證了機器人水動力系數(shù)和動力學模型的準確性，說明了利用動力學模型研究機器人水下運動特性的可行性。

4.2 水平面回轉運動

為了進一步驗證水動力系數(shù)求解和動力學模型的準確性，進行機器人水平面回轉運動實驗。實驗通過UWB（Ultra Wide Band，超寬帶）高精度導航定位系統(tǒng)采集機器人的位置信息，如圖10所示。UWB定位系統(tǒng)由基站、標簽和控制臺3部分組成，基站的4個Link Track P-B模塊分別位于矩形實驗區(qū)域四角；標簽同樣采用Link Track P-B模塊，安裝在機器人上，標簽可測量到基站的距離并解算出坐標信息；控制臺可與標簽、基站通信，在上位機軟件中實時顯示機器人位置信息。UWB定位系統(tǒng)的2維定位精度為10 cm；當定位頻率為200 Hz時，延遲為5 ms；最遠通信距離為1 500 m。

圖10 UWB定位系統(tǒng)原理示意圖Fig.10 Schematic diagram of UWB positioning system principle

實驗過程中，初始條件設置機器人推進力合力為10 N，接著基于初始條件分別設置兩側的推進器推力，在滿足推力要求的同時使機器人搖艏力矩為3 N·m。實驗過程如圖11所示，機器人從靜止加速進行水平面回轉運動。

實驗數(shù)據(jù)由于噪聲干擾出現(xiàn)小范圍波動，因此將實驗數(shù)據(jù)通過局部二次回歸進行平滑處理，機器人水平面回轉運動的仿真數(shù)據(jù)、實驗數(shù)據(jù)及擬合數(shù)據(jù)如圖12所示。

圖12 機器人水平面回轉運動軌跡Fig.12 Robot horizontal rotation motion trajectories

由于游泳池水流循環(huán)系統(tǒng)的影響，實驗場地并非靜水環(huán)境，存在隨機水流的干擾。由運動軌跡可知，在第2圈回轉運動過程中，機器人于1/4航程處偏航角發(fā)生了明顯變化，導致后續(xù)運動軌跡發(fā)生了偏轉，因?qū)嶒灜h(huán)境與仿真環(huán)境一致性不高，誤差較大。故以第1圈機器人軌跡進行數(shù)據(jù)分析，在推進力和搖艏力矩作用下機器人進行水平面回轉運動，運動軌跡近似為圓，與仿真數(shù)據(jù)重合度較高，規(guī)律一致，最大偏移量為0.15 m，偏差小于5%，滿足精度要求。實驗結果表明機器人水平面回轉運動特性與仿真結果相符，進一步驗證了動力學模型的準確性。

5 結束語

仿生波動鰭水陸兩棲機器人憑借隱蔽性強、機動性好等優(yōu)點成為海洋設備的研究熱點，然而由于現(xiàn)階段仿生機器人仿生程度低、推進性能差、作業(yè)環(huán)境復雜等客觀因素，使其在執(zhí)行任務中面臨著推進速度慢、效率低等困境。針對上述問題，本文圍繞新型水陸兩棲機器人水下動力學建模、水下操縱性數(shù)值仿真及實驗驗證等方面開展研究，仿真結果與實驗結果誤差均不超過10%，驗證了動力學建模的準確性。相關成果為進一步研究多模態(tài)水陸兩棲機器人運動性能奠定了理論與技術基礎。