吳宏宇，王春潔，丁宗茂，丁建中，董 洋，滿劍鋒

(1. 北京航空航天大學機械工程及自動化學院，北京 100083； 2. 北京航空航天大學虛擬現(xiàn)實技術與系統(tǒng)國家重點實驗室，北京 100083;3. 清華大學機械工程系, 北京 100084)

0 引 言

隨著科學技術的進步，人類的研究空間已不再局限于地球，探索月球、火星、乃至距地球更遠星球上的奧秘已成為當今科學界的研究熱點之一。而探測任務中，將機器人或宇航員安全送往外星球表面是探測任務的關鍵環(huán)節(jié)[1]。縱觀世界各國的深空探測任務，常用腿式著陸器搭載探測設備軟著陸外星球表面，如美國的阿波羅號載人登月艙、海盜2號火星探測器、鳳凰號火星探測器及中國的嫦娥三號月球探測器等。腿式著陸器作為一種簡單、可靠的著陸緩沖裝置，需要持續(xù)開展其研究工作，而研究工作的重點內(nèi)容包括著陸器軟著陸過程的動力學仿真與著陸器的性能優(yōu)化。

關于腿式著陸器軟著陸過程仿真計算與性能優(yōu)化的研究，已有較多學者開展了相關的工作。針對軟著陸過程仿真分析的研究，文獻[2-6]利用多體動力學分析軟件ADAMS建立了腿式著陸器軟著陸過程的動力學仿真模型，分析了著陸器在典型工況下的軟著陸性能，其中文獻[5-6]通過對比仿真結果與試驗數(shù)據(jù)證明了建模方法的合理性；文獻[7]基于基本力學原理，研究了著陸器沖擊動力學相似性技術，提出了等效低重力工況的新方法，并利用物理樣機試驗驗證了方法的有效性；文獻[8-10]基于Monte Carlo法研究了不確定著陸工況下的月球和火星著陸器的軟著陸性能；文獻[11-12]則考慮多種工況因子，重點針對腿式著陸器的抗翻倒能力與載荷緩沖性能繪制了穩(wěn)定性邊界；文獻[13-14]建立了腿式著陸器軟著陸過程的動力學方程，并對方程進行求解，研究了影響著陸器軟著陸性能的關鍵因素。針對腿式著陸器優(yōu)化設計的研究，文獻[15]基于響應面法與單一惡劣著陸工況，優(yōu)化了月球著陸器的緩沖機構；文獻[16]同時考慮多組惡劣著陸工況，直接利用動力學仿真模型參與優(yōu)化迭代計算，優(yōu)化了某型著陸器的緩沖機構；文獻[17]除考慮多組惡劣著陸工況外，同時考慮某型著陸器的兩種著陸模式，結合響應面法對某型著陸器的緩沖機構進行了優(yōu)化設計。

上述研究的集中點在于腿式著陸器軟著陸過程的仿真方法與特殊工況下的著陸器性能優(yōu)化方法，但鮮有文章報道著陸姿態(tài)不確定條件下的著陸器優(yōu)化方法。由于控制延遲、誤差等原因，著陸器與星球表面接觸時姿態(tài)具有隨機性，故而將著陸姿態(tài)的不確定性引入著陸器的優(yōu)化設計過程中非常必要。本文結合響應面法與分級優(yōu)化設計方法，以某型著陸器為例，重點研究了不確定著陸姿態(tài)下的著陸器緩沖機構優(yōu)化方法。

1 著陸器仿真模型建立

本文以圖1所示的某型著陸器為研究對象，著陸器由主體和四套完全相同的緩沖機構組成。緩沖機構均勻分布在主體四周，由主支柱、輔助支柱、足墊和緩沖桿組成，輔助支柱內(nèi)外筒間裝有吸收沖擊能量的緩沖元件[16-17]。

采用ADAMS軟件建立非0高度關機[17]軟著陸模式下的著陸器軟著陸過程仿真模型，模型從著陸器與星球表面接觸瞬間開始進行仿真分析，且各構件均簡化為剛體或利用剛體等效建模。

考慮大變形構件緩沖桿主要通過彎曲變形的方式緩沖載荷，其等效建模方式如圖2所示。首先，建立兩段軸線共線的剛體，其質量特性參考真實緩沖桿賦予。剛體1利用萬向節(jié)副與主體在O點聯(lián)接，并限制剛體2只能沿軸線移動。緩沖桿最初狀態(tài)時，軸線重合于X軸，當其末端受力導致其彎曲時，則等效于O點產(chǎn)生反作用彎矩M，將M投影在Y與Z軸上，得到兩個彎矩分量MY,MZ。反之，在O點處對剛體1施加兩個彎矩分量MY,MZ可以合成彎矩M。同時，為了模擬末端位置的變化，在剛體2上施加運動S。結合圖2可知，M,S均為O點與末端連線同X軸夾角α的函數(shù)。利用測量函數(shù)，求得緩沖桿末端相對于O點的坐標(x,y,z)，進而可得α的表達式如式(1)所示，并可得到M,S的表達式如式(2)和式(3)所示。

(1)

M=f(α)

(2)

S=g(α)

(3)

建立緩沖桿的有限元模型，末端施加彎矩，可仿真計算得到α-M，α-S函數(shù)關系f，g如圖3、圖4所示。

僅考慮受彎，緩沖桿的俯視圖為一條線段，如圖5所示。結合圖5，根據(jù)幾何關系，可計算出M在Y軸、Z軸的分量MY，MZ的表達式如式(4)、式(5)所示。

(4)

(5)

進而，利用式(3)～(5)結合圖2所示的等效模型可模擬緩沖桿的受力彎曲情況。

輔助支柱內(nèi)緩沖元件的吸能過程通過在內(nèi)外筒之間施加與壓縮行程D相關的緩沖力FD來模擬，F(xiàn)D與D的函數(shù)關系如圖6所示。

足墊和星球表面之間的垂向力和法向力分別采用非線性阻尼彈簧模型和庫倫摩擦模型來模擬[4,16]。

2 著陸工況參數(shù)及軟著陸性能判據(jù)

2.1 著陸工況參數(shù)

設著陸器與星球表面接觸瞬間的整機質量為800 kg，豎直速度vx為3.7 m/s、水平速度與三軸角速度為0，星球表面的坡度為10°、動摩擦系數(shù)為0.4[9,17]，不考慮著陸器足墊撞擊巖石或落入凹坑，此時著陸工況參數(shù)主要為著陸器的姿態(tài)角，而通過調整俯仰角θf與偏航角θp即可獲得評判第2.1節(jié)提出的軟著陸性能所需的著陸工況?？紤]著陸器整機的對稱性，設定θf與θp的取值范圍如表1所示。

表1 著陸姿態(tài)角取值范圍Table 1 Range of parameters of landing attitude angle

θf與θp的幾何含義如圖7所示，圖中xM，yM，zM表示著陸器的機械坐標系。

2.2 軟著陸性能判據(jù)

本文主要針對著陸器的軟著陸穩(wěn)定性，提出以下三點軟著陸性能判據(jù)：

1)著陸器抗翻倒能力判據(jù)。定義包含著陸器任意相鄰兩足墊中心點的豎直平面為翻倒平面，著陸器的抗翻倒能力與著陸器在著陸過程質心距翻倒平面的最小距離LD成正相關，且LD應始終大于0，否則認為著陸器翻倒[18]。

2)輔助支柱緩沖行程約束判據(jù)?？紤]輔助支柱內(nèi)部緩沖元件的承載能力有限，需要保證輔助支柱在正常工作狀態(tài)下的最大緩沖行程DM不大于80 mm。

3)著陸器底面抗損壞能力判據(jù)。著陸器的底面抗損壞能力與著陸器尾噴管底面中心點在著陸過程中距星球表面的最小距離HM成正相關，且HM過小不利于探測車的釋放，在此約束HM不小于350 mm。

3 著陸器軟著陸性能分析與緩沖機構優(yōu)化

3.1 不確定著陸姿態(tài)下的軟著陸性能分析

本文基于ADAMS模型，采用Monte Carlo法分析著陸器在不確定著陸姿態(tài)下的軟著陸性能，且認為各種著陸姿態(tài)對軟著陸性能的影響同等重要。首先利用描述性采樣方法在表1所示的姿態(tài)角取值范圍內(nèi)等概率抽取樣本點，之后將樣本點代入仿真模型進行仿真計算，最后對仿真模型輸出結果進行統(tǒng)計，以定量分析著陸器的軟著陸性能。本文抽取2000個樣本點進行仿真分析，根據(jù)2000組仿真結果繪制性能指標值的散點圖如圖8所示，并得到表2所示的軟著陸性能分析結果。

圖8中，方形散點表示HM小于350 mm的樣本點。

表2 Monte Carlo模擬分析結果Table 2 Analysis results of Monte Carlo simulation

分析結果顯示，著陸器底面抗損壞能力不理想，有必要對著陸器進行優(yōu)化設計以提升其軟著陸性能。

3.2 緩沖機構優(yōu)化設計

3.2.1優(yōu)化設計變量

為保證著陸器在運載火箭中的安裝要求，輔助支柱與主體的連接點位置、主支柱與緩沖桿的連接點位置、輔助支柱與主支柱的連接點相對于足墊的位置均難以調整[17]。由緩沖機構構型可知，著陸器足墊的中心點至主體中軸線的水平距離dz、主支柱長度LM與輔助支柱緩沖元件的緩沖力大小可作為優(yōu)化設計變量。dz，LM的幾何表示如圖9所示，緩沖元件緩沖力的大小利用其增減倍數(shù)kz來表達。

結合圖9，可直接分析出dz，LM和kz對著陸器軟著陸性能有以下三方面的影響：

1)在其他兩個參數(shù)不變的前提下，dz的增加會直接加大質心距翻倒平面的初值，減小主體的高度，故dz與著陸器的抗翻倒能力成正相關、與底面抗損壞能力成負相關。

2)在其他兩個參數(shù)不變的前提下，LM的增加會直接提升著陸器主體的高度，故dz與著陸器底面抗損壞能力成正相關，而主體提高的同時著陸器質心勢必增高，故LM的增加會間接導致著陸器抗翻倒能力下降。

在其他兩個參數(shù)不變的前提下，kz的增加會直接導致輔助支柱緩沖行程與主體下降距離減小，不難得知kz與著陸器的抗翻倒能力和輔助支柱最大緩沖行程成負相關、與底面抗損壞能力成正相關。

dz，LM和kz的初值與范圍如表3所示。

表3 設計變量范圍及初值Table 3 Initial values and range of design variable

3.2.2響應面代理模型

考慮結合Monte Carlo法，使用ADAMS模型參與緩沖機構的優(yōu)化迭代計算效率極低。對此，本文基于響應面方法建立LD，DM，HM與dz，LM，kz，θf，θp的映射關系，響應面模型的引入將大幅度提高優(yōu)化計算效率。本文采用二階多項式建立響應面模型，包含dz，LM，kz，θf，θp五個參數(shù)的二階多項式表達式：

φ=β0+β1dz+β2LM+β3kz+β4θf+β5θp+

β12dzkz+β13dzθf+β14dzθp+β15LMkz+β16LMθf+

β17LMθp+β18kzθf+β19kzθp+β20θfθp

(6)

式中：φ表示LD，DM，HM中任意一項的近似函數(shù)、βi表示多項式的各項系數(shù)。根據(jù)文獻[19]，利用最小二乘法確定βi的取值，進而可得響應面模型的具體表達式。

本文根據(jù)均方根相對誤差值RMSE和決定系數(shù)R2來判斷響應面模型的擬合精度，對應表達式如式(7)、式(8)所示[19]。

(7)

(8)

利用優(yōu)化拉丁超立方實驗設計在表1和表3中的著陸姿態(tài)與緩沖機構構型參數(shù)范圍內(nèi)抽取30個樣本點，利用ADAMS模型計算出30組LD,DM,HM值如表4所示，進而利用以上30組數(shù)據(jù)擬合響應面模型。

為避免多項式項數(shù)過多引起的多項式擺動等誤差源，本文引入了響應面模型關鍵項篩選策略，即對響應面模型的多項式進行逐項篩選，保留10項來擬合樣本點，選取RMSE最小的項組合方式來確定響應面模型的最終表達式如式(9)～式(11)所示。

2.732θpLM-121.739dzkz+1148.830dzLM

(9)

0.685θpLM-120.259dzkz+205.839dzLM

(10)

HM= -21936.723+23353.525dz+13572.952LM-

0.0345θfθp+2.612θfLM+402.763dzkz-

183.107kzLM

(11)

利用RMSE和R2分析響應面模型對30組樣本點的擬合精度，分析結果如表5所示。

Dynamic Analysis of Semi-type Floating Offshore Wind Turbine with Failure Conditions Under Metocean Conditions in South China Sea SHI Wei,ZHENG Kan,REN Nianxin (12)

檢驗結果顯示響應面模型可準確描述LD，DM，HM與dz，LM，kz，θf，θp的映射關系。

3.2.3優(yōu)化數(shù)學模型的建立及求解

上文的分析結果顯示，著陸器在初始構型參數(shù)下的底面抗損壞能力不理想，故本節(jié)將提出優(yōu)化方法計算dz，LM和kz的最佳取值，以提升著陸器的軟著陸性能。為使各種著陸姿態(tài)在優(yōu)化過程中起同等重要的作用以保證優(yōu)化結果的穩(wěn)健性，本文以最大化不確定著陸姿態(tài)下LD與HM的平均值LD-ave與HM-ave為目標，以不確定著陸姿態(tài)下DM最大值DM-max與HM的最小值HM-min不超過許用值作為約束條件建立優(yōu)化數(shù)學模型，為了防止優(yōu)化后DM-max與HM-min與許用邊界重合，優(yōu)化過程中將DM許用值降為75 mm,HM許用值提升為360 mm。綜上，得到優(yōu)化數(shù)學模型如式(12)所示。

表4 實驗設計結果Table 4 Result of design of experiment

表5 響應面精度檢驗結果Table 5 Accuracy analysis table of response surface models

max{LD-ave,HM-ave} s.t.DM-max≤75HM-min≥360x(L)O≤xO≤x(U)O(12)

本文結合式(9)與式(11)，利用多元積分法求解LD-ave與HM-ave，鑒于θf與θp相互獨立，多元積分求解如下所示：

(13)

122.948LM-121.739dzkz+1148.830dzLM

(14)

HM-ave=-21936.723+23353.525dz+13572.952LM-

402.763dzkz-183.107kzLM

(15)

本文結合式(10)與式(11)，在給定dz，LM，kz取值下，利用MIGA算法求解以θf,θp為設計變量的DM-max與HM-min，算法參數(shù)配置如表6所示。

表6 MIGA算法參數(shù)配置Table 6 Parameters configuration of MIGA

對于整個多目標優(yōu)化數(shù)學模型，本文利用NSGA-II算法求解，算法參數(shù)配置如表7所示。

表7 NSGA-II算法參數(shù)配置Table 7 Parameters configuration of NSGA-II

綜上，得到不確定著陸姿態(tài)下的著陸器緩沖機構分級優(yōu)化設計流程如圖10所示。

經(jīng)過優(yōu)化迭代計算，得到帕累托最優(yōu)解集，并繪制帕累托前沿曲線如圖11所示。

由于著陸器的翻倒更能直接導致探測任務的失敗，故在帕累托最優(yōu)解集中選取著陸器抗翻倒能力最強的解作為本文的最優(yōu)解，dz，LM和kz的取值如表8所示。

表8 最優(yōu)解取值Table 8 Value of optimum solution

將最優(yōu)解代入ADAMS模型，再次抽取2000個樣本點進行Monte Carlo模擬，根據(jù)2000組仿真結果繪制性能指標值的散點圖如圖12所示，并得到表9所示的軟著陸性能分析結果。

將表2、表9的分析結果做對比可知：相比于優(yōu)化前，優(yōu)化后著陸器的HM超出許用值的工況數(shù)由73組減為0組，LD與HM的平均值分別提升了3.546%和5.140%，優(yōu)化效果明顯。

4 結 論

本文在某型著陸器軟著陸動力學仿真模型的基礎上，結合Monte Carlo法研究了著陸姿態(tài)不確定條件下的著陸器軟著陸性能，并結合優(yōu)化拉丁超立方實驗設計得到樣本點，建立了以著陸器緩沖機構構型參數(shù)和著陸姿態(tài)參數(shù)為輸入、軟著陸性能指標值為輸出的響應面模型；針對著陸姿態(tài)不確定條件下的著陸器緩沖機構優(yōu)化問題，基于響應面模型，提出了結合NSGA-II算法與MIGA算法的分級優(yōu)化求解方法，并計算得到了使著陸器軟著陸性能最佳的緩沖機構構型參數(shù)；再次結合Monte Carlo法與動力學仿真模型校驗優(yōu)化結果，優(yōu)化后著陸器的抗翻倒能力和底面抗損壞能力分別提升了3.546%和5.140%，且尾噴管底面中心點距星球表面的最小距離與輔助支柱最大緩沖行程均不超出許用值。

指標超出許用值工況數(shù)平均值LD/mm0868.558DM/mm037.605HM/mm0374.169

本文所述的研究方法具有成本低、效率高等特點，也可為其他型號腿式著陸器的優(yōu)化設計提供參考。