何成陳國平何歡姜金輝

（1機械結(jié)構(gòu)力學(xué)及控制國家重點實驗室，南京 210016）（2南京航空航天大學(xué)振動工程研究所，南京 210016）

1 引言

探測器的著陸過程主要分為著陸下降減速和著陸接地緩沖兩個階段。在前一個階段，大多是利用降落傘的氣動阻力對登陸系統(tǒng)進行減速。然而，在實際工程中，降落傘降速能力的提高是以總質(zhì)量和總體積為代價，因此受降落傘系統(tǒng)總質(zhì)量和體積的限制，必須采取進一步的接地緩沖措施以確保登陸艙內(nèi)攜帶的精密、貴重儀器設(shè)備等有效負載的安全。

緩沖氣囊緩沖原理是利用系統(tǒng)著陸時，氣囊內(nèi)氣體壓縮來吸收著陸沖擊能量，從而達到緩沖的目的。由于其良好的緩沖性能，并且具有質(zhì)量輕，可折疊等特點，因此在汽車、航空、航天領(lǐng)域都得到了廣泛應(yīng)用[1,2]。近些年隨著深空探測發(fā)展，大量精密探測器的著陸安全對緩沖系統(tǒng)的緩沖特性提出了更高的要求，因此如何合理設(shè)計緩沖系統(tǒng)，進一步提高其緩沖性能就顯得尤為重要。

目前對氣囊緩沖系統(tǒng)的沖擊動力學(xué)研究方法主要分為試驗分析方法和有限元仿真分析方法兩類。由于試驗方法最能直觀反映出緩沖系統(tǒng)的真實特性，在國內(nèi)外開展了大量研究工作[3-5]。但是受實際試驗場地、試驗技術(shù)水平、試驗周期以及研究經(jīng)費等條件限制，無法開展大量、大規(guī)模驗證試驗和研究。因此采用有限元仿真分析方法成為近些年進行氣囊系統(tǒng)緩沖特性研究的熱點，文獻［6］以“火星探路者”著陸系統(tǒng)的全向式氣囊緩沖裝置為對象，利用有限元分析技術(shù)，實現(xiàn)氣囊系統(tǒng)的建模與著陸全過程的仿真研究；文獻［7］采用歐拉數(shù)值模型，研究了氣囊初始內(nèi)壓、材料參數(shù)以及地面參數(shù)對著陸系統(tǒng)過載的影響。

雖然在氣囊式著陸緩沖系統(tǒng)方面的研究論文較多，但在氣囊系統(tǒng)著陸緩沖性能優(yōu)化方面的研究開展相對較少。本文基于D最優(yōu)試驗設(shè)計建立的移動最小二乘（Moving Least Square，MLS）代理模型對類似于“獵兔犬”式氣囊緩沖系統(tǒng)進行沖擊動力學(xué)優(yōu)化研究，將首次沖擊最大過載、囊體織物最大應(yīng)力作為優(yōu)化目標，以囊內(nèi)初始氣壓、繩剛度作為設(shè)計變量，引入理想點法完成了緩沖系統(tǒng)的多目標優(yōu)化研究，并給出了多目標函數(shù)之間的Pareto前沿，為設(shè)計者使用封閉氣囊作為緩沖系統(tǒng)提供參考依據(jù)。

2 代理模型的構(gòu)造

2.1 D最優(yōu)試驗設(shè)計

由于實際試驗需要考慮到研究周期以及試驗成本等問題，因此在試驗前進行合理的試驗設(shè)計，可以在取盡可能少的試驗樣本點的情況下獲得較好的試驗效果，試驗設(shè)計的任務(wù)就是為了解決這一問題。在構(gòu)造全局近似函數(shù)時通常也借助試驗設(shè)計來減少總的響應(yīng)計算量。

D最優(yōu)設(shè)計[8]作為一種試驗設(shè)計方法在工程界得到了廣泛應(yīng)用，它是基于回歸模型的一類試驗設(shè)計方法，該方法對設(shè)計點的選擇自由，在樣本空間中無明顯的分布規(guī)律。假設(shè)研究的回歸模型為 ′，即

式中 x為設(shè)計空間中一點，設(shè)計空間為m維歐式空間；ε為服從方差為σ2的標準正態(tài)分布的測量噪聲；fi（x）為形式已知的函數(shù)；a為回歸模型系數(shù)。

將式（1）表示為矩陣形式

式中 F為模型的結(jié)構(gòu)矩陣，通過最小二乘法可獲得回歸系數(shù)向量a的估計值。估計參數(shù)協(xié)方差矩陣為

通常定義信息矩陣M=FTF，D最優(yōu)設(shè)計即為使得協(xié)方差矩陣取極小時所對應(yīng)的設(shè)計點集，或等價于使得信息矩陣的行列式取極大。

2.2 移動最小二乘法基本格式

MLS法最早是無網(wǎng)格法（SPH）逼近函數(shù)的方法之一，是無網(wǎng)格伽遼金法的核心。隨著學(xué)科間的交叉與滲透，MLS法在圖形處理和數(shù)據(jù)擬合等領(lǐng)域都獲得了較廣泛的研究與應(yīng)用[9]。近些年，在工程界逐漸將MLS法應(yīng)用于代理模型構(gòu)建中[10]，并很快形成了一種新的代理模型構(gòu)造方法。另外由于MLS法采用回歸近似技術(shù)并且系數(shù)具有可變性，從而使得它不僅對工程問題中的數(shù)值噪聲有很好的過濾作用而且在擬合高度非線性模型時具有良好的精度。

假定真實響應(yīng)函數(shù)為 ，由MSL法回歸得到的近似值為 ′，MLS法將待求函數(shù)表示為一系列基函數(shù)的線性疊加形式，即

式中 a=[a1（x），a2（x），…，am（x）]T為待定系數(shù)，它是 l維設(shè)計參數(shù)向量 x=｛x1，x2，…，xl｝的函數(shù)。p=[p1（x），p2（x），…，pm（x）]T為基函數(shù)，它是一個完備多項式，m為基函數(shù)項數(shù)。二維空間中單項式基函數(shù)為

為度量代理模型對目標函數(shù)的逼近程度，定義逼近函數(shù)殘差的離散加權(quán)L2范數(shù)J為：

生態(tài)環(huán)境遭受人為的嚴重破壞，生態(tài)環(huán)境日益惡化，自然森林植被、灌木、樹少林稀，極大地削弱了植被涵養(yǎng)水源、保持水源、調(diào)節(jié)氣候、減少蒸發(fā)的作用，致使本區(qū)地表、地下水資源短缺，地表河流及山間溪溝水流量逐年減少。該區(qū)內(nèi)許多巖溶泉水正因為周邊生態(tài)環(huán)境遭受破壞，導(dǎo)致流量逐年減少，甚至枯竭斷流。此外，由于生態(tài)環(huán)境的破壞，每逢豐水雨季徑流量大還造成水土流失，直接淤積江河及巖溶地下空間，從而引發(fā)嚴重的洪澇災(zāi)害。

式中i=‖x-xi‖／r，（r是影響半徑）；n是影響區(qū)域內(nèi)設(shè)計樣本點的個數(shù)；i是 xi處樣本點處響應(yīng)值，w（i）是樣本點xi的權(quán)函數(shù)，常用的有3次樣條權(quán)函數(shù)與高斯權(quán)函數(shù)，在本文中采用工程中常用的高斯權(quán)函數(shù)：

將式（2）寫成矩陣形式：

式中

由最小二乘法求得系數(shù)

式中

移動最小二乘法可以通過取不同階的基函數(shù)獲得不同的擬合精度，取不同的權(quán)函數(shù)以改變擬合曲面光滑度，這是其他擬合方法無法比擬的。需要注意的是，A矩陣可逆影響半徑必須滿足一定的條件：p（x）為二維線性基時，影響半徑內(nèi)至少包含3個不共線的樣本點；p（x）為二維二次基時，影響半徑內(nèi)至少6個樣本點，且這些點不在任何1條二次曲線上。

2.3 模型誤差評價準則

ANOVA統(tǒng)計分析技術(shù)不僅可以作為評價響應(yīng)面模型適合性的一種方法，并且能夠?qū)υO(shè)計變量的影響水平進行評價。本文選取常用與評價模型適合性的統(tǒng)計學(xué)判定系數(shù)R2與均方誤差根（RMSE）作為代理模型精度的評價參數(shù)。其中

式中 N為總的樣本點個數(shù)；i與i′是樣本點空間中各點的響應(yīng)真值和代理模型計算值；是樣本點空間中各點響應(yīng)真值的均值。RMSE的大小表示代理模型誤差值所占真實平均值的百分比；R2判定系數(shù)代理模型與真值之間的總體差異程度，取值范圍為0～1之間，一般來說，R2值越大，RMSE值越小，模型的合適性就越高。

3 封閉式氣囊系統(tǒng)的著陸緩沖性能多目標優(yōu)化

3.1 有限元模型及其優(yōu)化問題描述

本文以類似于“獵兔犬”緩沖登陸系統(tǒng)為對象，以真實火星環(huán)境為輸入條件，對該系統(tǒng)進行緩沖性能優(yōu)化。有限元模型如圖1所示，該登陸系統(tǒng)由氣囊織物、登陸艙以及收縮繩3部分組成。收縮繩采用MSCDytran所支持的“安全帶”一維單元進行模擬，并假設(shè)繩索的力與應(yīng)變呈線性關(guān)系。初始模型織物參數(shù)為，密度875kg/m3，彈性模量6.43GPa，泊松比0.3，阻尼系數(shù)0.1，收縮繩剛度為1.5×104N。

著陸緩沖系統(tǒng)以60°的入射角，著陸速度為18m/s進行著陸，接觸分析中采用的靜摩擦系數(shù)和動摩擦系數(shù)均取0.3。由于目前尚缺乏一些真實參數(shù)的數(shù)據(jù)支持，本文所取的分析參數(shù)有可能與真實情況不符，但并不妨礙所研究方法的普遍性。由于在整個著陸緩沖過程中，首次著陸沖擊過載遠遠大于后續(xù)彈跳過程中產(chǎn)生的沖擊過載，因此只以首次著陸的一個完整沖擊過載的峰值作為優(yōu)化目標。另外，在著陸沖擊過程中，氣囊織物承受較大的面內(nèi)張力，對系統(tǒng)進行合理的設(shè)計可以防止氣囊因沖擊載荷作用而產(chǎn)生爆破，因此對氣囊織物承受的最大應(yīng)力也需要進行優(yōu)化控制。

本文選取初始氣壓Pr與收縮繩剛度k（每單位應(yīng)變所需要的力）作為設(shè)計變量，以初始模型參數(shù)作為優(yōu)化空間中心，將設(shè)計變量Pr、k設(shè)計范圍分別設(shè)定為[3 000,6 000]Pa，[0,3]×104N，則優(yōu)化問題可描述為

式中f1為最大過載；f2為氣囊織物最大應(yīng)力；Φ為多目標優(yōu)化的目標函數(shù)。

3.2 樣本點的確定及其對應(yīng)的響應(yīng)

采用D最優(yōu)試驗設(shè)計方法，選取樣本點數(shù)目為22個。通過樣本點對應(yīng)的樣本值，采用非線性瞬態(tài)有限元分析程序MSC.Dytran作為核心求解器，計算得到在樣本點處的最大過載與最大應(yīng)力。由于計算得到的加速度信號中含有大量的高頻成分，掩蓋了對登陸艙設(shè)備過載起決定作用的低頻脈沖，故采用文獻[11]的方法對過載信號進行濾波，本文選取的濾波截止頻率為60Hz。表1給出了樣本點的空間位置及其對應(yīng)的目標函數(shù)值。

表1 樣本點及對應(yīng)的目標函數(shù)值Tab.1 Design matrix and values of objective functions

3.3 優(yōu)化分析

利用MSL構(gòu)造設(shè)計變量與各個目標函數(shù)之間的映射關(guān)系，其對應(yīng)的響應(yīng)面如圖2、3所示。對每個目標函數(shù)的代理模型適合性進行評價，計算發(fā)現(xiàn)：最大過載的響應(yīng)面RMSE值僅為0.001 8，判定系數(shù)R2最小值為0.980 7，最大應(yīng)力的響應(yīng)面RMSE值僅為0.006 3，判定系數(shù)R2最小值為0.967 2，可見每個模型完全滿足精度要求，可以進行后續(xù)優(yōu)化工作。

從圖2、3中可以看到，在設(shè)計空間內(nèi)最大過載極大值為78.2gn，而最大過載極大值與極小值之間相差18.5gn，接近極大值的24%。同樣的最大應(yīng)力極大值為310MPa,極大值與極小值之間相差154MPa，占到極大值的49%，可見各個目標函數(shù)都具有較大的可優(yōu)化空間；設(shè)計變量對目標函數(shù)的影響程度有差別，目標函數(shù)對收縮繩剛度的變化更為敏感，而隨初始氣壓的變化相對平緩；系統(tǒng)的最大過載隨囊內(nèi)部初始氣壓和收縮繩剛度增大而增大，在低壓力區(qū)收縮繩剛度影響較為復(fù)雜；總體而言最大應(yīng)力呈現(xiàn)出隨收縮繩剛度增加而增大的趨勢，但是在設(shè)計空間的頂點處最大應(yīng)力有不同程度的升高。

引入理想點法，通過構(gòu)造統(tǒng)一目標函數(shù)使多目標問題轉(zhuǎn)化為單目標優(yōu)化問題進行優(yōu)化。令各分目標的理想最優(yōu)值為（i=1，2），定義

式中 q為目標函數(shù)個數(shù)；wi是第i個目標函數(shù)所對應(yīng)的加權(quán)因子，文中采用歸一權(quán)處理。由于響應(yīng)面表現(xiàn)出的非線性程度不同，特別是針對響應(yīng)面是多峰函數(shù)的情況，因此采用全局尋優(yōu)能力較強的遺傳算法進行尋優(yōu)。

通過變化不同的權(quán)系數(shù)可以得到側(cè)重不同目標分量的Pareto解集。解集是一個全局的相對最優(yōu)解，設(shè)計者可以在解集中根據(jù)不同的優(yōu)化要求與側(cè)重點進行選取。本文采用Monte-Carlo模擬結(jié)合MSL代理模型，計算獲得了目標函數(shù)的Pareto前沿，如圖4和圖5所示。為便于區(qū)分，在圖中進行了編號處理。

圖4 緩沖系統(tǒng)設(shè)計參數(shù)Pareto前沿Fig.4 The Pareto front of design parameters

圖5 緩沖系統(tǒng)優(yōu)化目標的Pareto前沿Fig.5 The Pareto front of objective functions

從圖中可以看出：1）Pareto解集多集中于特定域內(nèi)，且多分布于設(shè)計空間的邊界，呈現(xiàn)出聚類分布特征。2）在設(shè)計空間頂點“1”區(qū)域，可以獲得設(shè)計域內(nèi)最大過載的極小值，若決策者更加注重對登陸艙過載的控制，可以選擇該區(qū)域內(nèi)的設(shè)計點。3）當初始氣壓處于3 000Pa附近的低壓力區(qū)時（對應(yīng)圖4與圖5中“1”-“3”區(qū)域），收縮繩的剛度是對目標函數(shù)起主要作用的因素，最大過載與最大應(yīng)力表現(xiàn)出相反的變化趨勢，即織物最大應(yīng)力隨著首次最大過載的增大而減小。4）當系統(tǒng)初始氣壓在高壓力區(qū)時（對應(yīng)圖4與圖5中“4”區(qū)域），可以得到設(shè)計域內(nèi)織物最大應(yīng)力的全局最優(yōu)值。若決策者較為關(guān)注氣囊織物強度，可以選擇該區(qū)域內(nèi)的設(shè)計點。

分別從“1”和“4”區(qū)域內(nèi)選出優(yōu)化解分別定義為“優(yōu)化系統(tǒng)1”與“優(yōu)化系統(tǒng)2”，圖6和圖7分別給出了優(yōu)化后系統(tǒng)的變形和首次沖擊過載的時域曲線。圖7結(jié)合表2可以看出優(yōu)化系統(tǒng)不僅最大程度的優(yōu)化了重點關(guān)注的單個目標函數(shù)，而且其他優(yōu)化目標性能也獲得了不同程度的提高，采用此種優(yōu)化方法達到了預(yù)期效果。

圖6 優(yōu)化系統(tǒng)2的應(yīng)力與變形圖Fig.6 The stress of second landing system after optimization

圖7 優(yōu)化前后系統(tǒng)沖擊過載時域曲線Fig.7 Time history of overloading before and after optimization

表2 優(yōu)化前后系統(tǒng)緩沖性能比較Tab.2 Cushioning performance of landing system before and after optimization

4 結(jié)束語

本文研究了基于D最優(yōu)試驗設(shè)計結(jié)合MLS法的氣囊著陸緩沖動力學(xué)多目標優(yōu)化問題。優(yōu)化結(jié)果表明，利用MSL構(gòu)造的代理模型能夠非常好的逼近強沖擊載荷下的沖擊響應(yīng)響應(yīng)面，在此基礎(chǔ)上進行優(yōu)化分析能夠給出合理的優(yōu)化分析結(jié)果。

文章針對氣囊著陸緩沖系統(tǒng)的囊內(nèi)部初始氣壓和收縮繩剛度進行了優(yōu)化分析，有效的降低了氣囊著陸沖擊最大過載和織物上承受的最大應(yīng)力，證明了所提出的優(yōu)化設(shè)計方法的有效性。研究結(jié)果還表明，收縮繩的引入使得初始氣壓對目標函數(shù)影響規(guī)律更加復(fù)雜，最大過載和最大應(yīng)力呈現(xiàn)出較高的非線性特征?？傮w而言最大過載表現(xiàn)為隨囊內(nèi)部初始氣壓和收縮繩剛度增大而增大的趨勢，最大應(yīng)力隨收縮繩剛度增大而增大，且收縮繩剛度對最大過載和織物最大應(yīng)力的影響更加顯著。

（References）

[1] 戈嗣誠.火星全向氣囊的著陸緩沖特性研究[J].航天返回與遙感,2012,33(2)∶16-22.GE Sicheng.Study on Cushioning Characteristic of Mars Omni-direction Airbag[J].Spacecraft Recovery and Remote Sensing,2012,33(2)∶16-22.(in Chinese)

[2] 萬志敏,謝志民.氣囊緩沖飛行器模型著陸特性的試驗研究[J].試驗技術(shù)與試驗機,2003,43(1)∶9-12.WANG Zhimin,XIE Zhiming.Study on Cushioning Characteristic of Airbag Landing System by Experiment[J].Test Technology and Testing Machine,2003,43(1)∶9-12.(in Chinese)

[3] Gardinier D J,Taylor A P,.Design and Testing of the K-1 Reusable Launch Vehicle Landing System Airbag[R].AIAA,1999,1757∶1-10.

[4] Cadogan D,Sandy C,Grahne M.Development and Evaluation of the Mars Pathfinder Inflatable Airbag Landing System[J].Acta Astronautica,2002，50(10)∶633-640.

[5]DaveN,Chris M.Improved Inflatable Landing Systems for Low Cost Planetary Landers[J].Acta Astronautica,2006，59∶726-733.

[6] 鄧春燕,裴錦華.全向式氣囊著陸裝置緩沖過程的仿真研究[J].中國空間科學(xué)技術(shù),2010,30(1)∶78-83.DENG Chunyan,PEI Jinhua.Simulation about Buffer Process of Omni-directional-type Airbag Landing Device[J].Chinese Space Science and Technology,2010,30(1)∶78-83.(in Chinese)

[7] 戴華杰,胡振東,咸奎成，等.火星探測器氣囊緩沖系統(tǒng)著陸過程仿真[J].力學(xué)季刊,2010,31(4)∶555-561.DAI Huajie,HU Zhendong,XIAN Kuicheng,et al.Simulation Analysis for Landing Process of A Mars Detector with Airbag Buffer System[J].Chinese Quarterly of Mechanics， 2010,31(4)∶555-561.(in Chinese)

[8] Bu S L,Wang S Q,Ye H Q.An Algorithm Based on Variable Feedback to Synchronize Chaotic Systems[J].Physica D,2002,164∶45-52.

[9] 邵衛(wèi)云,張雄.水泵水輪機全特性曲線的擬合-移動最小二乘近似[J].水力發(fā)電學(xué)報.2004,23(5)∶102-106.SHAO Weiyun,ZHANG Xiong.A New Simulation Method of Complete Characteristic Curves of Reversible Pump Turbine-moving Least Square Approximation[J].Journal of Hydroelectric Engineering.2004,23(5)∶102-106.(in Chinese)

[10] 張勇.基于近似模型的汽車輕量化優(yōu)化設(shè)計方法[D].湖南大學(xué),2008.ZHANG Yong.Optimization Design Method of Vehicle Lightweight Based on Approximate Model[D].Hunan University,2008.(in Chinese)

[11] Burrow L T.Verification Testing of a UH-1 Wire Strike Protection System(WSPS)[M].S Army Applied Technology Laboratory(AVRDCOM),1982.