地面力學在火星壤力學參數估計研究中的進展與展望

薛 龍，黨兆龍，陳百超，李建橋，鄒 猛

(1. 江西農業大學工學院，南昌 330045；2. 江西省現代農業裝備重點實驗室，南昌 330045；3. 中國空間技術研究院北京空間飛行總體設計部，北京 100094；4. 吉林大學工程仿生教育部重點實驗室，長春 130022)

0 引 言

為了尋找地外生命和可宜居星球，從登陸月球到對火星及太陽系中其他的天體的探測，人類的深空探測腳步從未停歇。但是由于地外行星惡劣的自然環境，使得人類必須使用可移動的星球車來完成對地外行星的研究和探測。即使人類可以在地外行星行走，也需要星球車來輔助人類完成各項科學探測及擴大人類探索的范圍[1-2]。這就要求星球車具有良好的機動性、平順性和通過性，使之能夠順利通過松軟、顛簸的地外行星表面[3-5]。星球車的越野機動性不僅與地面障礙物和地形特性有關，還與車輛外形、車輪形狀[6]和地外行星壤的特性有關[7]。而全面準確的了解行星表面星壤的物理力學參數，對于輪式星球車的路徑規劃、通過性能評估、自主導航、風險預測、仿真計算和脫困具有十分重要的意義。

由于受到火箭發射載荷的限制，目前的火星車并沒有攜帶專用的測量火星壤力學參數的設備，因而無法實時準確的獲取火星壤的力學參數。但根據車輛地面力學理論，基于輪壤作用數學模型可對輪下及火星車周邊火星壤的力學參數進行辨識。

為了更好的支持火星探測任務，本文探討了復雜火星地面對火星車移動性能的影響，分析了基于地面力學的火星壤力學參數辨識方法，總結了基于車輪的在軌火星壤力學參數辨識方法，給出火星壤的力學參數測試結果，并對未來地面力學在火星巡視探測領域的發展方向進行展望。

1 反演的目的和意義

目前，美國是火星探測最成功的國家，共有八個著陸器和火星車已經成功登陸火星表面并進行科學探測研究，其登陸時間和地點以及運行狀態如表1所示[8]。從美國火星著陸器和火星車發回的照片和數據可知，火星表面遍布著隕石坑、火山、沙丘與峽谷，且大部分被一層松散的火星壤覆蓋，一旦火星車車輪打滑并陷入其中，將導致巡視探測任務失敗[9]。例如，2009年5月7日，勇氣號的輪子陷入火星特洛伊沙地，車輪刨開了地表的一層硬殼，硬殼大約有2.45 cm(1英寸)厚，硬殼下面露出了松軟的含沙物質。當勇氣號開始嘗試擺脫困境時，它的車輪又攪起了更多的含沙物質。這些物質含有高濃度的硫酸鹽，比在火星上其他地區所看到的物質中硫酸鹽的濃度要高得多。根據車輪的直徑(約26 cm)，車輪所陷入的深度大約有30 cm。為了幫助勇氣號脫離困境，地面控制人員在土槽中對其狀態進行了模擬，以期找出使其脫困的方法。但隨著火星嚴冬的臨近，勇氣號的太陽能帆板已經無法吸收到足量的陽光供發電之用，它最終沒能挨過殘酷的冬天[10]，勇氣號在火星表面共運行了7.7 km。

勇氣號并不是唯一一個陷入松軟沙地的火星車，在梅里迪亞尼平原機遇號的六個車輪也曾陷入沙地達30 cm。為使機遇號脫困，科研人員同樣利用模擬火星壤進行地面土槽試驗，以圖尋找最佳的方法使機遇號脫離束縛。所幸的是經過23個火星日的努力，在地面指令的控制下機遇號終于擺脫了束縛[11]，而此處沙地也被正式的命名為“Purgatory Ripple”如圖1所示。勇氣號和機遇號車輪的沉陷主要是因為松軟的火星壤加劇了車輪的滑轉和車輪的沉陷，即車輪與火星壤所構成的地面車輪系統所能提供的驅動力不能克服地面的阻力。

表1 美國著陸器和火星車登陸地點和運行狀態

圖1 機遇號車輪沉陷于Purgatory Ripple處的軟沙中[11]

通過機遇號陷入沙地的深刻教訓，當地面人員控制好奇號在火星表面前進時，會定期的通過拍攝圖像的方式來觀察其是否確實在前進，并根據在線傳感器計算車輪實際的移動距離來確定火星車的總體滑轉率，并確保車輪的滑轉率低于限定值[12]。但是，由于地球和火星之間通訊時差的限制，這種方法也有其缺陷。例如2015年5月7日，好奇號就不得不停下其探索的“腳步”，其原因是由于好奇號車輪打滑，使得實際的行駛距離遠遠低于預計的19 m。實際上，從2015年5月7日至2015年5月13期間，好奇號試圖到達觀測目標進行研究的行進過程中，共行駛四次并且其中三次都因為車輪打滑而不得不中止行駛，以避免好奇號陷入沙地中。通過好奇號相機拍攝的圖像顯示，火星表面的沙丘橫亙在好奇號和觀測目標之間，因此在路徑規劃時選擇了沙丘周圍的一處看起來更堅固的地面，但是當好奇號達到這個區域時，意外發現這里竟然也是由松散物質構成的，說明火星表面具有非常大的欺騙性。圖2為好奇號在沙丘表面留下的車轍痕跡。

圖2 好奇號在火星表面留下的車轍[13]

可見，采用火星車進行探索和科學研究時，了解和掌握其周邊火星壤的力學狀態尤為重要，特別是在路徑規劃、風險評估、牽引控制等方面有著不可替代的作用。

地面力學是研究車輛與松軟地面之間的相互作用關系的一門學科[14-20]，主要用于設計車輛移動系統、評估車輛通過性能和分析車輛在非常規路面的機動性能。地面力學研究首先需要了解土壤的力學參數，例如：承壓模量kφ、內聚模量kc、沉陷指數n、內聚力c、內摩擦角φ和剪切變形模量K，其中前三個土壤參數稱為承壓力學參數(kc，kφ和n)，后三個土壤參數稱為剪切力學參數(c，φ和K)，而這些參數的測定需要專用的試驗設備對土壤進行測試才能獲得。但是在地外行星探測過程中，由于發射重量和其它研究項目等方面的各種限制，目前在軌的火星車都沒有配備專門的設備測量火星壤的物理力學參數，而是需要應用火星車或登陸器的觸土部件(車輪或者采樣鏟)，基于地面力學理論對火星壤的力學參數進行估算。車輪作為火星車重要組成部分，承擔載重、牽引和移動的功能，并且始終與地表接觸，因此根據輪壤相互接觸模型對火星壤的力學參數進行估計成為可能。

2 輪壤接觸模型簡介

當剛性車輪以一定速度向前運動時，在輪上載荷W和車輪扭矩T的共同作用下，車輪與松軟地面間產生的正應力σ(θ)并不是線性分布，而是以最大正應力為分界線劃分成兩個區域，以最大應力角θm和入土角θ1之間的區域為第一個應力區域，第二個應力區域包含離去角θ2和最大應力角θm，其正應力大小分別是σ1和σ2。車輪與土壤之間的相互作用關系如圖3所示，進而得到車輪的受力平衡方程。

(1)

(2)

(3)

式中：W為輪上載荷(N)，FDP為掛鉤牽引力(N)，T為車輪扭矩(N·m)，θ1為入土角(rad)，θ2為離去角(rad)，θm為最大應力角(rad)，σ(θ)為正應力(Pa)，τ(θ)為剪應力(Pa)，r為車輪半徑(m)，b為車輪輪寬(m)。

從車輪的受力平衡方程可以看出，土壤的力學特性對車輪牽引通過性有著重要的影響，其中承壓特性和剪切特性最為重要。承壓特性表示的是壓力和沉陷的關系，不僅可以對車輪沉陷和阻力進行評估，還可以評估其牽引特性，常用的承壓模型主要有Bekker模型[21]、Reece模型[22]、Ishigami模型[23]和滑轉模型[24-26]等。剪切特性表示驅動輪輪下土壤剪切應力和剪切位移之間的關系，常用的剪切模型主要有Janosi模型[18]、指數模型[27]和Wong-Reece應力分布模型[18-20]等?？梢娏私饣鹦侨懒W特性，即可計算車輪的牽引力、輪上載荷和車輪扭矩，也可以為車輪設計、風險評估和通過性估計提供依據。但是在火星探測過程中，由于缺乏火星壤力學參數，因此可以通過車輪的受力平衡方程，同時應用已知條件(如輪上載荷、車輪扭矩和滑轉率以及車輪沉陷等)，在不增加火星車自身重量和專用檢測設備的前提下，實現地面力學參數的辨識。

圖3 剛性車輪與土壤相互作用關系

3 輪壤接觸模型的火星壤力學參數反演方法

3.1 模型簡化

由于輪壤表達式(1)～式(3)中包含復雜的積分，很難求得精確的解析解，因此需要對求解表達式進行簡化。通過對四種不同類型土壤(干砂、砂壤土、粘土和雪)對比分析，輪下的正應力和剪切應力的分布是對稱于最大應力角，并且其分布曲線可以近似的看作線性分布，即以最大應力角為分界的兩個簡化的線性方程表示[28-31]，然后對式(1)～式(3)進行求解，得到包含正應力、剪切應力和入土角的簡化解。同時以Mohr-Coulomb失效準則為基礎，應用式(4)建立了正應力和剪切應力及土壤力學參數之間的關系，該模型也稱為Wong-Reece應力分布模型：

(4)

式中：s為滑轉率。通過式(4)、承壓模型以及簡化解，即可建立輪壤之間的動力學模型，對火星壤的力學參數進行辨識。

3.2 牛頓迭代法

應用Bekker承壓模型，根據圖3的輪壤關系，建立FDP的解析表達式[32-35]。該方法可應用復合辛普森公式、兩點高斯數值積分法對解析表達式進行簡化求解。由于內聚力變化相對于掛鉤牽引力不敏感，即當內聚力變化范圍很大時，其計算掛鉤牽引力的變化很小。相反，當掛鉤牽引力的值變化很小時，對反演的內聚力值的變化影響卻很大。另外，當把內聚力和其它土壤參數同時進行反演時，將導致其它土壤參數值與其真實值不符。因此為了避免這種現象的發生，在求解過程中需要設定內聚力為一個合理的固定值。

該方法對φ，ks和K進行求解，其中ks為組合內聚模量，ks=kc/b+kφ。由于該模型共包含有三個未知變量，因此至少需要三組試驗數據作為輸入數據才能求得收斂解，這三組數據中分別包括掛鉤牽引力、車輪滑轉率和車輪沉陷，即[FDP,s,z]。對于土壤剪切參數中的內聚力，在求解過程中設置為一個合理的平均值(c=3 kPa)，參與土壤力學參數反演。

3.3 最小二乘法

最小二乘法是基于經典的地面力學數學模型，通過求解[c,φ]=f(W,T,z,s)來估計內聚力和內摩擦角，其中z為車輪沉陷。該求解方法令θm=θ1/2，并應用簡化的線性模型，得到最大正應力σm和最大剪切應力τm，進而得到包含c和φ的表達式[28-31]，如式(5)所示：

(5)

在土壤力學參數反演的過程中，可以通過火星車車載傳感器測量記錄其行駛時的數據(例如z，T，W及s等)。z可以通過導航相機和避障相機拍攝的車轍圖像進行估計[36]，T可以通過車輪驅動電機電流計算獲得，W可以通過火星車的位姿進行估計，s可以通過車轍圖像輪刺的間距或者機器學習回歸算法進行估計得到[37]。由于火星車行駛速度較慢，因此在一段較短的距離內就可以收集多對行駛數據，通過收集的數據即可利用最小二乘法求解火星車車輪周邊行星壤的力學特性參數。

3.4 耦合法求解

耦合法是一種循環迭代的地面力學參數識別方法，分別利用掛鉤牽引力、車輪沉陷和前進阻力矩對接觸角系數、承壓特性參數和剪切特性參數進行辨識求解[38]。正應力σ(包括前部正壓力σ1和后部正壓力σ2兩部分)的計算模型中，用輪地作用沉陷指數N代替了傳統模型中的土壤沉陷指數n，用于反映車輪的滑轉沉陷，N由沉陷指數系數n0和n1計算得到。c1和c2為最大應力角系數，通常離去角很小，因此在計算過程中設為零，即c3=0。可計算得到車輪前部和后部的正應力表達式[24]。應用式(4)、簡化的線性模型[28-31]和正應力表達式[24]可得輪上載荷、掛鉤牽引力和車輪扭矩的解析表達式：

(6)

W，FDP和T是關于所有土壤力學參數的函數，是高度耦合的方程組，因此對眾多土壤參數同時求解是不合實際的。進行地面力學參數求解時需要先設定初始變量，給定輪上載荷，可以消除剪切特性參數PIII，同時參數PI對于沉陷量的影響很小，因此可以令c1=0.5，c2=0作為初始值，然后根據測量得到的沉陷量對承壓特性參數PII進行反演。根據反演得到的PII和給定的角接觸參數，及車輪扭矩對剪切特性參數PIII進行反演。最后，利用反演的PII和PIII進行角接觸參數PI的反演。通過上述步驟可以近似實現對土壤參數的估計，但是模型的簡化帶來了誤差，并且辨識得到的參數PI與進行參數PII和PIII辨識時所采用的參數PI不一致。因此在獲得初步結果的基礎上采用循環迭代的方法進一步進行參數反演。然后按照PII=fW(W,T,s,z)，PIII=fT(W,s,z)和PI=fDP(W,T,FDP,s,z)的順序循環求解[38-43]。

3.5 神經網絡法

當輪壤相互作用時，根據式(1)和式(3)建立輪上載荷和車輪扭矩的數學模型，利用簡化的線性模型對正應力和剪切應力簡化。而最大應力角為θm=(θ1+θ2)/2，車輪與地面的接觸角θc=|θ1-θ2|，由于離去角很小因此令θ2=0。對式(1)和式(3)求解，得到W和T的解析表達式[44]，如式(7)和式(8)所示：

(7)

(8)

可見，該公式建立了火星壤力學參數c,φ與s,T,W的理論模型，如式(9)：

[c,tanφ]T=f(s,T,W)

(9)

當給定的四種地面數據類型的剪切力學參數和車輪半徑和輪寬，根據特定車型的試驗數據分別選定z(0.005,0.0075,0.01,0.0125,0.015,0.0175,0.02) m，s(0.01,0.02,0.03)和I(50,60,70,80,90,100) mA電機電流可以轉換為車輪扭矩T，特征值的選擇是應用重量為30 kg的微型探測車Kapvik在松軟土壤行駛時的試驗數據為依據選定的[45]。對應給定的z，求解滑轉沉陷zss與接觸角θc[46]。根據s求解j=r(θc+ssin(θc)-sin(θc))，I得到車輪扭矩T。然后計算τm和σm。通過式(7)計算對應的W，同時建立對應的輸入數據對[W,T,s]和輸出數據對[c,tanφ]。遍歷所有給定的四種地面數據類型，形成一系列的輸入數據對和輸出數據對共504組。然后，應用504組數據訓練神經網絡，建立預測土壤內聚力和內摩擦角的評估模型，輸入變量和輸出變量分別為x=[W,T,s]和y=[c,tanφ]。

該方法可以在不需要獲得車輪沉陷的前提下，對火星壤的剪切力學參數進行在線估計，并且可對每個車輪周圍的火星壤進行同時估計。不僅有利于為充分了解巡視探測器的行駛狀態，也有利于巡視探測器的自主導航。

3.6 支持向量機法

最小二乘支持向量機(Least-squares SVM， LS-SVM)法通過訓練數據集建立土壤剪切力學參數辨識模型，對火星壤的剪切力學參數進行就位估計。

通常車輪行駛地面特別是地外行星表面，其土壤為塑性地面，車輪行駛過后，土壤的回彈量非常小，基本為零。當應力角θ分別等于0和θm時，可得最大剪切應力和剪切應力偏移量，通過簡化的線性模型[28-31]得到輪上載荷、車輪扭矩和掛鉤牽引力的表達式[47]，如式(10)～(12)所示。

(10)

(11)

(12)

建立火星壤剪切力學參數與s,T,W的理論模型，[c,φ,K]T=f(s,T,W)。

針對定型車輛的車輪，選定s,z和T的合理變化范圍[47]。確保所選的值包含盡可能多的數據，并且其值的范圍包含了輕載-土槽試驗臺所采集的s,z和T的最大值和最小值。給定車輪的幾何參數，通過解析模型計算入土角，最大剪切應力和最大正應力，計算輪上載荷，建立輸入變量x=[s,T,W]和輸出變量y=[c,φ,K]組成的訓練集數據對共計4237對。建立LS_SVM模型時，應用參數尋優函數Tunelssvm進行網格搜索，優化gam和sig2這兩個重要的參數[48]。

當火星車在火星土壤表面行駛時，通過車載傳感器獲取滑轉率、車輪扭矩和輪上載荷后，就可以根據建立的土壤參數辨識數學模型，得出車輪所處位置火星壤的剪切力學參數。

3.7 遺傳算法

遺傳算法(Genetic algorithm, GA)用于反演行星壤的承壓力學參數：組合變形模量和承壓指數。根據簡化的線性模型[28-31]，建立輪壤關系模型的矩陣表達式[49]，式(13)：

(13)

可以計算得到輪上載荷和掛鉤牽引力的表達式如式(14)和式(15)：

(14)

(15)

(16)

因此由式(14)、式(15)和式(16)可得，W和FDP的表達式，如式(17)和式(18)所示：

(17)

(18)

對于火星車而言，W,T,s和θ1可由火星車上攜帶的傳感器測量或者估計獲得[37, 51]。W與T,s,z和ks,n1,n2成函數關系，其表達式為：

Fz=fz(W,s,ks,n1,n2)

(19)

FW=fW(T,s,z,ks,n1,n2)

(20)

3.8 其他算法

自適應卡爾曼濾波算法(Adaptive robust extended Kalman filter，AREKF)，該算法可對預測沉陷指數N和內摩擦角φ進行反演，而其他力學參數需要給定一個經驗值[52-53]，該方法計算精度高，可適用于在線檢測。粒子濾波法(Particle filtering method，PFM)假設土壤參數為具有一定概率密度的隨機變量，對土壤的承壓力學參數(n，ks)和剪切力學參數進行估計[54]，該方法耗時較多。

表2總結了基于輪壤模型的火星壤力學參數辨識的方法，并對各種方法的輸入參數、輸出參數和優勢進行了總結，其中耦合法具有較高的計算精度和土壤力學參數的全面辨識，可實現離線測量。牛頓法具有計算速度快和辨識精度高的優點，適合于在線分析，可實現剪切力學參數的辨識和部分承壓力學參數的辨識。神經網絡法和支持向量機法，具有計算速度快和不需了解車輪沉陷的條件下進行土壤剪切力學參數在線辨識的優點。遺傳算法計算精度高，可實現對承壓力學參數的辨識，可適用于在線辨識。

表2 基于輪壤模型的反演方法

4 火星壤參數辨識在火星車上的應用

4.1 索杰納車輪反演

美國國家航空航天局(NASA)和噴氣推進實驗室(JPL)的科研人員利用索杰納火星車單輪旋轉剪切土壤，根據庫侖定律，對火星土壤的內聚力c和內摩擦角φ進行估測，進而判定該區域的可通過性。結果顯示，在軌火星車周邊火星表面的火星壤內摩擦角均值為35.0°，內聚力均值為0.01 kPa[56]。索杰納火星車在火星表面行走了105 m，進行了14次關于火星壤力學的試驗、6次空載電流標定和11次車輪磨損試驗[57]。當車輪在無負載條件下，車輪扭矩M與車輪驅動電機的電流和電機溫度呈線性關系[56-57]。

(21)

式中：M為扭矩，y為變量，I電機電流，x為電機無負載時的電機電流，t車輪驅動電機的溫度。y和x均隨著溫度t的變化而改變。在火星表面，通過抬起火星車兩個前輪和一個右邊中輪，然后驅動抬起的車輪旋轉，同時記錄電機電流，該電流為無負載電機電流。

土壤的強度采用Mohr-Coulomb失效準則，建立正應力和剪切應力之間的關系：

S=c+Ntanφ

(22)

式中：S為剪切應力，N為正應力。

在分析過程中，建立后輪(正轉)和前輪(逆轉)的輪上載荷與剪切應力、火星車重量、車輪沉陷深度、懸架力臂和彈簧剛度數據之間的關系方程[57]。正應力通過估計的輪上載荷、火星重力加速度和車輪在火星表面的投影面積估算，剪切應力通過車輪扭矩、車輪半徑和車輪投影面積估算。車輪投影面積通過車輪的寬度和車輪沉陷深度位置對應的弦長計算。應用最小二乘法和Mohr-Coulomb失效準則，對剪切應力和正應力數據對進行擬合，計算內聚力和內摩擦角。

在探路者號登陸火星后的第13火星日，在火星表面“Yogi”附近，索杰納進行右后輪試驗，車輪沉陷深度為1.3 cm，得到內聚力和內摩擦角分別為0.06 kPa和40.9°[57]。根據試驗結果，“cloddy deposits”類型火星壤的內聚力和內摩擦角的平均值分別為0.17±0.18 kPa和37.0°±2.6°。

4.2 基于MER車輪的反演

勇氣號和機遇號(Mars Exploration Rovers, MER)也通過車輪與地面的相互作用，進行了7次挖掘試驗和20次刮擦試驗，對古謝夫坑(Gusev Crater)和梅里迪亞尼平原(Meridiani Planum)進行了火星壤的內聚力和內摩擦角的估算，其中對內聚力的估算采用電機功的比例模型和輪緣剪切應力模型[58-60]。進行車輪挖掘和刮擦試驗時，應用電機電流、電壓和車輪持續于土壤相互作用的運行時間計算電機所做的功[58]。計算內摩擦角時，假定內聚力為零，則內摩擦角僅與車輪的正應力和剪切應力有關，因此可得

φ=arctan((A·V·t)/(RθFN))

(23)

式中：θ為車輪挖掘旋轉角度(rad)，實驗時前車輪順時針旋轉。電機做功需要減去電機內耗。車輪每旋轉60°(30個數據點)，計算一次內摩擦角，然后計算平均值，R為車輪半徑R=0.125 m，FN為法向力。

內聚力的計算有兩種方法，一種是電機功的比例模型(Ratios of electromechanical work)，另一種是輪緣剪切應力模型(Modeling shear stress along the wheel rim)，這兩種方法都需要已知內摩擦角和接觸面積。根據挖掘試驗得出火星壤的內摩擦角30°～37°，內聚力的值為0～2 kPa。而根據車輪與地面的刮擦試驗，得出火星壤的內聚力為0～11 kPa，不確定性的范圍為±0.7～±3.9 kPa。

車輪挖掘的方法被用于在軌估計土壤的力學特性，但是其缺點是針對每個不同的方法都需要建立校正模型，并且結果誤差大，而且應用車輪作為土壤力學測試的部件一旦出現沉陷等意外狀況就會使火星車無法移動，嚴重威脅火星車的行駛通過性。另外用車輪旋轉剪切火星壤的方法僅僅是在剪切應力最大的條件下獲得，不能獲取其它的地面力學參數，特別是承壓力學參數。

5 結論與展望

應用輪壤作用模型，可以對火星車移動系統進行設計和優化提供理論支撐，還可以利用輪壤相互作用模型對火星壤的力學參數進行辨識。

通過不同的分析方法得到的土壤力學參數具有差異性，其準確性和適用性還需要進一步的試驗驗證。因此，針對我國火星探測計劃的巡視探測器，開展整車與不同火星地表極端環境下的試驗，找出適合探測任務并可評估整車在極端環境下的通過性的方法是一個重要的研究方向，還需要進行大量的基礎性的研究工作。

另外，由于好奇號在火星地表運行過程中出現車輪磨損的情況，因此設計新型的具有一定剛度的彈性車輪，以適應苛刻的火星地表環境發揮最大的行駛性能，同時建立相適應的輪壤相互作用模型用于預測整車的通過性，也是地面力學領域中熱門研究方向之一。