徐 顥，居鶴華，余 萌，顧 朋

（1.南京航空航天大學(xué)航天學(xué)院，南京210016；2.北京控制工程研究所，北京100190）

目前，火星已經(jīng)成為世界上各個(gè)國家主要的探索對象[1]?；鹦茄惨暺魃涎b載了多種科學(xué)儀器，是人類執(zhí)行火星表面探測任務(wù)的重要載體。由于火星表面環(huán)境具有復(fù)雜性與不確定性，巡視器應(yīng)具備穿越各種崎嶇地形的能力，從而保證火表探測工作的可靠性[2]。例如：美國的“機(jī)遇號”火星巡視器因車輪陷入沙土而無法動(dòng)彈[3]；“好奇號”火星巡視器自2012 年8 月登陸火星至2013 年10 月執(zhí)行科學(xué)探測任務(wù)期間，車輪本體結(jié)構(gòu)破損嚴(yán)重[4]。因此，針對防止巡視器發(fā)生結(jié)構(gòu)性破壞以及巡視器的運(yùn)行狀態(tài)等問題，國內(nèi)外學(xué)者對行星巡視器的運(yùn)動(dòng)學(xué)、動(dòng)力學(xué)等問題進(jìn)行了大量深入的研究。

Yoshida 等[5]從輪土力學(xué)和鉸接體動(dòng)力學(xué)兩方面出發(fā)研究了行星車的動(dòng)力學(xué)問題，建立了懸架機(jī)構(gòu)的動(dòng)力學(xué)仿真模型；以車輪滑移率為狀態(tài)變量，在行星車試驗(yàn)臺上進(jìn)行了實(shí)驗(yàn)研究。為了研究在崎嶇地形中的高速移動(dòng)機(jī)器人的行為控制方法，Iagnemma 等[6]開發(fā)了一套動(dòng)力學(xué)軟件，該軟件采用四階龍格庫塔法求解以“牛頓?歐拉”分析動(dòng)力學(xué)方法建立的非線性多階常微分方程，并以四輪移動(dòng)機(jī)器人為例，進(jìn)行了仿真和實(shí)驗(yàn)。Lindemann[7]通過動(dòng)態(tài)測試模型（Dynamic test model, DTM）評估了火星巡視器在模擬火星環(huán)境中的移動(dòng)特性，同時(shí)在ADAMS 軟件中建立了巡視器的動(dòng)力學(xué)模型并給出了仿真結(jié)果。居鶴華等[8]研究了六輪搖臂式月面巡視器在復(fù)雜月面環(huán)境下的動(dòng)力學(xué)建模問題，采用牛頓?歐拉法建立了巡視器動(dòng)力學(xué)模型，利用開發(fā)的三維仿真平臺對建立的動(dòng)力學(xué)模型進(jìn)行了驗(yàn)證。目前為止，國內(nèi)針對行星巡視器動(dòng)力學(xué)的相關(guān)研究較少，仍依賴于ADAMS、ODE、Vortex 等商業(yè)或開源動(dòng)力學(xué)軟件進(jìn)行二次開發(fā)，缺少自主理論和相關(guān)技術(shù)研發(fā)，難以滿足復(fù)雜的行星表面環(huán)境探測的工程需求。

針對復(fù)雜火星表面環(huán)境及高自由度的火星巡視器動(dòng)力學(xué)建模問題，本文引入了運(yùn)動(dòng)鏈符號系統(tǒng)，其具有偽代碼的功能，易于計(jì)算機(jī)編程，計(jì)算效率高；同時(shí)，本文方法可以動(dòng)態(tài)地建立巡視器的動(dòng)力學(xué)方程，其具有完全參數(shù)化、分層迭代式的優(yōu)點(diǎn)。

1 火星巡視器結(jié)構(gòu)

一類典型火星巡視器的三維模型如圖1 所示。巡視器采用主動(dòng)式搖臂結(jié)構(gòu)，是在被動(dòng)式搖臂的基礎(chǔ)上，在左、右主臂與差速器輸出端之間以及左、右副臂與各自的主臂之間各加裝一個(gè)控制關(guān)節(jié)，從而實(shí)現(xiàn)主臂與車體、副臂與主臂的角度控制功能。車體上的差速器保證其出軸與本體的角度大小相等，方向相反。巡視器6 個(gè)驅(qū)動(dòng)輪均具有舵機(jī)，則該主動(dòng)式搖臂結(jié)構(gòu)具有16 個(gè)獨(dú)立的自由度。6 個(gè)驅(qū)動(dòng)輪的舵機(jī)、主臂關(guān)節(jié)、副臂關(guān)節(jié)及差速器均具有絕對編碼器，驅(qū)動(dòng)輪轉(zhuǎn)動(dòng)速度由相對編碼器進(jìn)行測量，用作里程計(jì)。巡視器還具有4 個(gè)單自由度的太陽翼、1 個(gè)2 自由度的數(shù)傳天線及1 個(gè)3 自由度的桅桿等。

圖1 火星巡視器三維模型Fig.1 3D model of the Mars rover

在機(jī)器人工程中，首先定義坐標(biāo)系，然后通過工程測量確定坐標(biāo)系間的關(guān)系，最后以該坐標(biāo)系為參考進(jìn)行機(jī)器人的運(yùn)動(dòng)學(xué)及動(dòng)力學(xué)分析。巡視器的軸間運(yùn)動(dòng)關(guān)系同樣需要通過與它們固結(jié)的坐標(biāo)系進(jìn)行度量。

為了便于分析，本文采用自然坐標(biāo)系描述巡視器的軸間關(guān)系。巡視器移動(dòng)系統(tǒng)自然坐標(biāo)系如圖2 所示，當(dāng)系統(tǒng)處于初始位置時(shí)，只需建立慣性系以及各軸上的參考點(diǎn)，其他運(yùn)動(dòng)軸坐標(biāo)系也自然確定，而不必為每一運(yùn)動(dòng)軸確定各自的體坐標(biāo)系。

圖2 巡視器各部件的體坐標(biāo)系編號與部件英文縮略標(biāo)志符一一對應(yīng)，原點(diǎn)均位于部件質(zhì)量中心，x、y、z 表示坐標(biāo)軸，下標(biāo)i 和c 分別表示系統(tǒng)慣性系和車體系；下標(biāo)lr、rr 分別表示左、右主臂，lb、rb 分別表示左、右副臂；三位下標(biāo)編號的首位縮略標(biāo)志l 或r 分別表示左或右部件，第二位縮略標(biāo)志f、m、r 分別表示前、中、后部件，第三位縮略標(biāo)志w表示車輪、d 表示方向機(jī)，例如lfd 表示左前方向機(jī)，lfw 表示左前輪，以此類推；坐標(biāo)軸xc、yc軸及zc軸分別指向巡視器正前方、正左側(cè)及正上方。

圖2 巡視器移動(dòng)系統(tǒng)自然坐標(biāo)系Fig.2 Natural coordinate system of the rover mobile system

2 火星巡視器動(dòng)力學(xué)建模研究

2.1 基本符號定義

本文以運(yùn)動(dòng)鏈符號演算系統(tǒng)為基礎(chǔ)[9]，通過精確的符號定義以解決高自由度系統(tǒng)符號規(guī)范的問題。為方便閱讀，下面對本文出現(xiàn)的相關(guān)符號及公式做簡單說明。

基本定義：在運(yùn)動(dòng)鏈符號演算系統(tǒng)中，具有偏序的屬性變量或常量，在名稱上包含表示偏序的指標(biāo)，包含左上角及右下角指標(biāo)或者右上角及右下角指標(biāo)；而它們的方向總是由左上角指標(biāo)至右下角指標(biāo)。

符號定義如下：

圖3 軸l 的閉子樹Fig.3 Closed subtree of Axis l

（2）任何復(fù)合的運(yùn)動(dòng)副可由兩個(gè)基本運(yùn)動(dòng)副組成，即轉(zhuǎn)動(dòng)副與平動(dòng)副，分別用R 和P 來表示。

表1 符號定義Table 1 Symbol definition

2.2 火星巡視器動(dòng)力學(xué)建模

火星巡視器動(dòng)力學(xué)建模分為兩個(gè)步驟。首先，采用有向Span 樹[9]表征巡視器的拓?fù)潢P(guān)系，本質(zhì)上反映軸與軸的連接關(guān)系以及軸間運(yùn)動(dòng)量的參考關(guān)系與作用關(guān)系；隨后，遍歷火星巡視器的軸鏈拓?fù)浣Y(jié)構(gòu)有向Span 樹，完成巡視器的動(dòng)力學(xué)建模。

(1)建立火星巡視器拓?fù)潢P(guān)系

火星巡視器的運(yùn)動(dòng)學(xué)與動(dòng)力學(xué)依賴于其拓?fù)潢P(guān)系分析，本文采用有向Span 樹來表征巡視器的拓?fù)潢P(guān)系，即忽略巡視器各部分結(jié)構(gòu)的尺寸，僅考慮運(yùn)動(dòng)副和相互連接關(guān)系而組成的結(jié)構(gòu)，用Span樹的節(jié)點(diǎn)編號表達(dá)動(dòng)力學(xué)系統(tǒng)的運(yùn)動(dòng)副。火星巡視器的有向Span 樹如圖4 所示。

火星巡視器通過運(yùn)動(dòng)軸的連接確定了一個(gè)運(yùn)動(dòng)鏈，該鏈?zhǔn)怯行虻倪\(yùn)動(dòng)軸的集合，記22 軸巡視器移動(dòng)系統(tǒng)的運(yùn)動(dòng)鏈為A=(0,1,2,…,22]。

圖4 火星巡視器軸鏈有向Span 樹Fig.4 Directed Span-tree of the Mars rover axis-chain

(2)建立火星巡視器動(dòng)力學(xué)方程

一個(gè)多自由度的火星巡視器系統(tǒng)可視為一個(gè)多軸剛體系統(tǒng)即iL，由轉(zhuǎn)動(dòng)副與平動(dòng)副組成，下面將推導(dǎo)其拉格朗日遞歸方程。一般的拉格朗日方程表示為

式中

由式(1)得多軸系統(tǒng)拉格朗日方程

動(dòng)力學(xué)系統(tǒng)iL 的平動(dòng)動(dòng)能及轉(zhuǎn)動(dòng)動(dòng)能分別表示為

考慮式(2)和式(3)，即有

至此，已完成能量對關(guān)節(jié)速度及坐標(biāo)的偏速度計(jì)算。

至此，已完成對時(shí)間t 的求導(dǎo)。

由偏導(dǎo)數(shù)鏈規(guī)則可知[9]

根據(jù)式(18)偏導(dǎo)數(shù)關(guān)系，將式（12，15～17）代入式(9)，可以推導(dǎo)出軸l 的拉格朗日遞歸動(dòng)力學(xué)方程如下

由式(19)可知，對于多軸系統(tǒng)iL，遍歷火星巡視器的軸鏈拓?fù)浣Y(jié)構(gòu)有向Span 樹即可完成對巡視器的動(dòng)力學(xué)建模。由此可見，針對高自由度的火星巡視器動(dòng)力學(xué)建模問題，本文的動(dòng)力學(xué)方法可以動(dòng)態(tài)地建立基于運(yùn)動(dòng)軸的有向Span 樹，具有完全參數(shù)化、分層迭代式的優(yōu)點(diǎn)，較傳統(tǒng)的牛頓?歐拉法、拉格朗日法、凱恩法[11?13]及目前國際上的6D 空間算子代數(shù)法[14?15]具有更低的計(jì)算復(fù)雜度；同時(shí)，該方法引入了運(yùn)動(dòng)鏈符號系統(tǒng)，具備了偽代碼的功能，使得計(jì)算機(jī)編程更加簡單，計(jì)算效率高于傳統(tǒng)的拉格朗日遞歸法，適用于高自由度、可變拓?fù)涞臋C(jī)器人動(dòng)力學(xué)建模問題。

在建立系統(tǒng)動(dòng)力學(xué)方程后，緊接著就是方程求解的問題，此為本文后續(xù)的工作重點(diǎn)，根據(jù)文獻(xiàn)[9]，可將式(19)更進(jìn)一步表達(dá)為規(guī)范型方程，并進(jìn)行多軸系統(tǒng)動(dòng)力學(xué)求解，對此本文不再展開討論。

3 輪土力學(xué)與輪式系統(tǒng)移動(dòng)維度研究

輪土作用力是輪式多軸系統(tǒng)移動(dòng)的動(dòng)因，利用輪土力學(xué)知識預(yù)測火星巡視器的車輪?土壤相互作用性能，對巡視器的輪土作用力識別、動(dòng)力學(xué)控制具有重要意義。Bekker 輪土模型目前應(yīng)用廣泛，但該模型沒有定義驅(qū)動(dòng)輪與土壤相互作用過程中相關(guān)物理參數(shù)的參考系。因此，本文將鏈符號系統(tǒng)引入經(jīng)典Bekker 輪土模型，以提高其準(zhǔn)確性。

3.1 Bekker 輪土力學(xué)

Bekker[16]輪土力學(xué)模型是專門為彈塑性土和剛性驅(qū)動(dòng)輪設(shè)計(jì)的作用力模型，該模型建立基于水平地形及重力方向載荷。

(1)輪土作用幾何關(guān)系

圖5 輪土模型幾何關(guān)系Fig.5 Geometrical relationship of terramechanics model

根據(jù)圖5 所示幾何關(guān)系，可得

則輪速v 表達(dá)為

(2)改進(jìn)的輪土應(yīng)力模型

根據(jù)文獻(xiàn)[17]，可得正應(yīng)力σ(z)及剪應(yīng)力κ(z)的表達(dá)式，Yoshida 等[18]應(yīng)用法向最大應(yīng)力角θm、泥土進(jìn)入角θ1、輪半徑r、滑移率將σ(z)及κ(z)的表達(dá)式修改為

式中：S(θ)和C(θ)分別代表θ 的正弦與余弦；其他的輪土參數(shù)如表2 所示

表2 輪土參數(shù)Table 2 Terramechanics parameters

給定輪土參數(shù)及滑移率，由式(27)得到θm，由式(24)及 式(25)分 別 得σm及κm；然 后，由 式(26)得σ1(θ)，σ2(θ)，κ1(θ)及κ2(θ)。

由正應(yīng)力分量σ(θ) ?C(θ)、剪應(yīng)力分量κ(θ) ?S(θ)及式(24)至式(26)，得到輪地接觸點(diǎn)Ol對輪的正壓力表達(dá)式

可以證明正壓力Nllˉ是關(guān)于進(jìn)入角θ1的單調(diào)函數(shù)，即

(3)牽引力、牽引力矩及側(cè)向力

圖6 滑移角Fig.6 Slip angle

式中：A 為輪地 接觸面（m2）；為 滑 移 角；B 為 常 數(shù)（1/rad）；為正壓力（N）。

3.2 輪地接觸約束方程

理想的輪地接觸可視為接觸副，該接觸副僅約束了關(guān)節(jié)速度而不是位形，故接觸副存在單邊約束，即輪地正壓力總大于等于零。

記u″矢量空間為接觸副u″Ok″的約束矢量空間，它是u 的子空間，u″位形空間除維度不同外，與矢量空間u 基本一致。因此，可以定義約束副u″Cu。

接觸副u″Ok″的約束要求體u″不能刺入體k″中，故有

記接觸面法向支撐力為u″Sfk″S，接觸條件為

由式(33)至式(35)表明接觸副約束是不等式方程組。

下 面 建 立 接 觸 副u″SOk″S的 單 邊 約 束 方 程?？紤]接觸副u″SOk″S，接觸點(diǎn)記為u″S，接 觸力記為，接觸點(diǎn)位置記為uru″S。接觸副u″SOk″S是特殊的運(yùn)動(dòng)副u″Skk″S，其 中 一 個(gè) 平 動(dòng) 軸 受 單 邊 約 束，轉(zhuǎn) 動(dòng) 軸 受控。記u″對k″施 加 的i|u″Sτk″S或i|u″Sfk″S為 廣 義 約 束力，i|u″JuS及i|k″JkS為 廣 義 約 束 雅 可 比 矩 陣。考 慮 單位質(zhì)量及轉(zhuǎn)動(dòng)慣量，則有

式中

由式(33～39)得

至 此，式(40～42)建 立 了 接 觸 副u″SOk″S的 單 邊約束方程。顯然，輪地接觸的前向牽引力及側(cè)向力是正壓力的函數(shù)，輪地作用的合力與正壓力成正比。因此，該合力也是單邊作用力，通過一個(gè)單邊約束方程表征。

基于上述Bekker 輪土力學(xué)代數(shù)方程及式(19)，可得完整的巡視器動(dòng)力學(xué)模型，其中包含與25 個(gè)軸對應(yīng)的25 個(gè)動(dòng)力學(xué)方程、與一個(gè)差速器對應(yīng)的一個(gè)約束方程及與6 個(gè)輪地接觸點(diǎn)對應(yīng)的6 個(gè)輪地接觸單邊約束方程。

4 仿真平臺實(shí)現(xiàn)

基于VS2013 軟件中的MFC 框架，采用C++語言開發(fā)了火星巡視器三維可視化仿真平臺，在仿真平臺中，建立了模擬地形場景及一類典型火星巡視器三維模型。仿真平臺主要?jiǎng)恿W(xué)參數(shù)中，巡視器車體質(zhì)量127 kg，相對于本體系的三軸主轉(zhuǎn)動(dòng)慣量為[40,40,40] kg·m2，車輪質(zhì)量為3.674 kg，主轉(zhuǎn)動(dòng)慣量為[2,2,2] kg·m2，其余部件的質(zhì)量和轉(zhuǎn)動(dòng)慣量較小。以如圖7 和圖8 所示的火星表面模擬環(huán)境為場景，進(jìn)行可視化仿真分析，巡視器可根據(jù)期望的移動(dòng)路線、速度及加速度移動(dòng)，同時(shí)在期望路徑上完成原地轉(zhuǎn)向、直線移動(dòng)、圓弧運(yùn)動(dòng)等基本行為。

圖7 巡視器沿曲線行駛Fig.7 Rover trajectory:curvature steering

圖8 巡視器原地轉(zhuǎn)向Fig.8 Rover trajectory:circling around

圖9～14 給出了巡視器在圖7 模擬場景下的運(yùn)動(dòng)狀態(tài)及對應(yīng)的方向機(jī)力矩、角度，車輪角速度、力矩，巡視器姿態(tài)及速度曲線圖。其中，每幅圖的點(diǎn)劃線、虛線及實(shí)線分別與火星巡視器部分結(jié)構(gòu)的前、中及后部件對應(yīng)。

圖9 方向機(jī)力矩Fig.9 Steering gear torque

圖10 方向機(jī)角度Fig.10 Steering gear angle

圖11 車輪角速度Fig.11 Angular velocity of wheel

圖12 驅(qū)動(dòng)輪力矩Fig.12 Wheel torque

火星表面重力加速度設(shè)為3.72 m/s2，圖9～14表明巡視器由起始至5.0 s 時(shí)的動(dòng)力學(xué)過程，包含500 仿真步，每步10 ms。圖9 表明巡視器在行駛過程中方向機(jī)力矩變化較大，這是車輪接觸崎嶇地形導(dǎo)致的結(jié)果。圖10 表明兩側(cè)方向機(jī)角度變化幅度較小，因?yàn)榍? s 內(nèi)火星巡視器剛開始與地面接觸。圖11 表 明 左 側(cè) 車 輪 角 速 度 由-2.0 (o)/s 至2.1 (o)/s變 化，右 側(cè) 車 輪 角 速 度 由-1.2 (o)/s 至2.1 (o)/s 變化；且右側(cè)車輪角速度變化幅度大于左側(cè)，符合轉(zhuǎn)向特性。圖12 反映6 個(gè)車輪所受到輪土力矩隨著輪土作用而逐漸變化，且輪土支撐力之和小于巡視器在火表加速度下的總重力。圖13 和圖14 反映了巡視器的前5 s 姿態(tài)坐標(biāo)是(200,200,0)，且巡視器以約0.1 m/s 速度前進(jìn)時(shí)三軸的平動(dòng)速度及轉(zhuǎn)動(dòng)速度，其中巡視器的偏航角變化較大，這是由于巡視器在崎嶇的火表地形沿曲線行駛，前側(cè)驅(qū)動(dòng)輪較后側(cè)驅(qū)動(dòng)輪的沉陷量大，反映了地形模擬的真實(shí)性。

由于火星表面環(huán)境具有復(fù)雜性與不確定性，傳統(tǒng)的被動(dòng)六輪搖臂-轉(zhuǎn)向架火星巡視器只能被動(dòng)地適應(yīng)地形，在不確定強(qiáng)的地表不能進(jìn)行自我調(diào)節(jié)，例如：巡視器行駛在松軟地形時(shí)易發(fā)生車輪沉陷等問題。本文巡視器采用主動(dòng)式搖臂結(jié)構(gòu)，具有主臂與車體、副臂與主臂的角度控制功能，例如：在復(fù)雜的火表中，通過巡視器懸架的主動(dòng)變形，改變了巡視器移動(dòng)系統(tǒng)的質(zhì)心相對位置，可以實(shí)現(xiàn)車體高度調(diào)整、抬輪、“蟹形”移動(dòng)等，增加了對不確定性地表的適應(yīng)能力，更加適合巡視器探測崎嶇復(fù)雜地形。

圖13 巡視器姿態(tài)Fig.13 Attitude of the Mars rover

圖14 巡視器速度Fig.14 Speed of the Mars rover

5 結(jié) 論

針對中國未來行星表面探測任務(wù)的復(fù)雜性與不確定性，本文對火星巡視器動(dòng)力學(xué)原理進(jìn)行了研究。首先對火星巡視器的結(jié)構(gòu)進(jìn)行介紹并采用自然坐標(biāo)系統(tǒng)描述巡視器移動(dòng)系統(tǒng)。然后對本文的基本符號進(jìn)行定義，建立了巡視器的樹形拓?fù)浣Y(jié)構(gòu)，應(yīng)用鏈符號系統(tǒng)推導(dǎo)了多軸拉格朗日方程，計(jì)算能量對關(guān)節(jié)速度和坐標(biāo)的偏速度及對時(shí)間的導(dǎo)數(shù)，遍歷巡視器拓?fù)浣Y(jié)構(gòu)即可完成巡視器的動(dòng)力學(xué)建模，該動(dòng)力學(xué)方法具有完全參數(shù)化、分層迭代式的優(yōu)點(diǎn)。接著，通過改進(jìn)的輪土應(yīng)力模型，得到輪地作用力和力矩，同時(shí)建立了輪地作用接觸副的約束方程，用于后續(xù)巡視器動(dòng)力學(xué)建模。最后，在開發(fā)的可視化仿真平臺對巡視器模型進(jìn)行測試。結(jié)果表明，巡視器可以根據(jù)期望的速度及加速度移動(dòng)，同時(shí)可在平坦或崎嶇地形完成原地轉(zhuǎn)向、直線、圓弧移動(dòng)等基本行為，體現(xiàn)了本文動(dòng)力學(xué)方法的有效性與優(yōu)越性，對中國未來行星車動(dòng)力學(xué)研究及仿真軟件的開發(fā)具有借鑒意義。