張 峰,楊旭鋒,劉永壽,南 華,2
(1.西北工業大學力學與土木建筑學院飛行器可靠性工程研究所,陜西 西安710129;2.蘭州萬里航空機電有限責任公司,甘肅 蘭州730700))
緩沖系統是飛機起落架的重要組成部分。飛機在著陸過程中受到很大的沖擊載荷,這要求緩沖系統能消耗掉飛機著陸沖擊帶來的能量。緩沖系統性能的好壞直接影響到飛機著陸過程的安全性,其故障輕者導致飛機受損,重者釀成機毀人亡的慘劇。因此,緩沖系統的設計技術在起落架設計中占有重要的地位。緩沖器是緩沖系統中重要的組成部件,其參數配置的恰當與否,對緩沖系統性能有著決定性的影響[1,2]。
目前,起落架緩沖系統的設計仍基于安全系數法的設計理念[3~6]。文獻[4]建立了前起落架單腔定油孔和單腔變油孔的虛擬樣機、主起落架單腔定油孔和雙腔定油孔的虛擬樣機,分別開展在使用功和最大功情況下的落震仿真分析。文獻[5]建立某型飛機主起落架緩沖器的虛擬樣機模型,采用基于正交試驗設計的響應面方法獲得優化參數和優化目標間的數學關系,并通過優化來提高多支柱起落架的著陸緩沖性能。文獻[6]借助 MSC.ADAMS軟件建立了起落架緩沖系統多體運動學仿真模型,研究緩沖系統充填參數對著陸性能的影響。
需要指出的是,這些研究忽略了緩沖器在設計、制造過程參數的隨機不確定性[7]。研究表明,緩沖器一些關鍵設計參數的微小波動會導致緩沖系統性能有很大的變化。可靠性靈敏度采用失效概率對基本變量分布參數的偏導數來度量,得到的靈敏度序列能快速識別出影響系統可靠性的關鍵參數,對系統改設計提供指導[8,9]。因此,本文開展起落架緩沖器參數的靈敏度分析。首先采用MSC.ADAMS軟件建立起落架系統落震仿真模型,通過多體動力學分析對緩沖系統的效率、緩沖器的行程以及起落架最大過載等指標進行考核;考慮設計、制造、操作以及環境因素帶來的影響,將緩沖器的初始充填壓力、初始充氣體積、壓氣面積、活塞桿面積、氣體多變指數、壓油面積、主油孔面積、主油液流量系數和油液密度等參數處理為隨機變量,并以緩沖系統的效率、緩沖器的行程以及起落架最大過載等指標建立多模式下緩沖器參數的靈敏度分析模型,采用基于Kriging方法的代理模型建立緩沖器輸出響應與設計參數的函數關系,再采用重要抽樣法求解緩沖器參數的靈敏度。
某型飛機搖臂式起落架采用的是油氣式緩沖器,其結構模型如圖1所示。
在MSC.ADAMS軟件中建立起落架系統落震仿真模型如圖2所示,通過多體動力學理論來分析緩沖系統著陸過程中的性能。

圖1 起落架緩沖器結構模型Fig.1 Shock absorber model of a landing gear system

圖2 某型飛機起落架系統落震仿真模型Fig.2 The drop simulation model of an aircraft landing gear system
通過起落架系統落震仿真分析得到緩沖器行程為294.7mm,緩沖系統的效率為79.97%。起落架輪胎所受地面沖擊載荷變化曲線如圖3所示,所受地面最大載荷為340.7kN,著陸過載為1.66。緩沖系統的效率η、緩沖器的行程Stroke以及起落架最大過載max(Force)等指標均滿足設計要求。

圖3 起落架所受沖擊載荷的變化曲線Fig.3 The change curve of the landing gear impact load
由于設計、制造、操作以及環境因素帶來的影響,導致緩沖器的初始充填壓力、初始充氣體積、壓氣面積、活塞桿面積、氣體多變指數、壓油面積、主油孔面積、主油液流量系數和油液密度等參數為隨機變量。當隨機變量之間具有相關性時,需要預先通過Nataf變換將相關變量轉化為相互獨立的隨機變量[10]。
緩沖器參數的分布類型及其均值、變異系數、標準差如表1所示,參數的變異系數通過查閱相關手冊或憑工程經驗獲得[2]。

表1 隨機變量的分布類型及參數Tab.1 Distribution type and parameters of random variables
假定隨機變量xk(k=1,2,…,9)的概率密度函數為fk(xk),向量x=(x1,x2,…,x9),則系統不確定性變量的聯合概率密度
考核緩沖器結構設計是否合理,重在于考核緩沖系統的效率、緩沖器的行程以及起落架的最大過載等三個指標。由于(x1,x2,…,x9)均為隨機變量,導致η,Stroke和 max(Force)均為隨機變量,是向量x的隱函數,需要借助MSC.ADAMS軟件進行起落架系統落震仿真分析來獲得。
在評估緩沖系統功能可靠性時,采用系統的效率η(x)、緩沖器的行程Stroke(x)以及起落架的最大過載max(Force(x))等三個指標分別建立相應的功能函數。
(1)效率失效:該緩沖器作為起落架系統主要吸收能量的裝置,其效率η(x)不能低于預定值η*,由此定義功能函數g1(x)如下式

(2)行程失效:該緩沖器的使用行程Stroke(x)不得大于0.9倍的最大行程S*,由此定義功能函數g2(x)如下式

(3)過載失效:該起落架著陸載荷實際值的最大值max(Force(x))不能大于最大著陸載荷容許值F*,由此定義功能函數g3(x)如下式

三個指標中的任一一個不滿足設計要求,都會影響緩沖系統的功能,所以這三個指標為串聯關系。定義緩沖系統功能可靠性的功能函數為G(x),G(x)可表示為

緩沖系統的失效概率Pf用下式來計算

式中D為整個變量空間,I[·]為指示函數,滿足。
根據文獻[8,9],定義緩沖器參數的可靠性靈敏度為失效概率對基本變量分布參數的偏導數。式(5)對第k個變量xk的分布參數θk(θk為μk或σk)求導數,得到可靠性靈敏度如下式所示

式(6)中,當θk=μk時,?fx(x)/?θk可寫為

當θk=σk時,?fx(x)/?θk可寫為

在緩沖器參數可靠性靈敏度分析中,存在兩個難題。
(1)由于緩沖器的復雜性,緩沖系統的效率η(x)、緩沖器的行程Stroke(x)、著陸載荷實際值的最大值max(Force(x))等指標均是x的隱函數,需要借助MSC.ADAMS軟件進行起落架系統落震仿真分析才能得到響應值。可靠性靈敏度分析多次調用MSC.ADAMS仿真軟件進行功能函數的計算,這需要很大的計算量。因此,緩沖器參數的可靠性靈敏度分析需要借助高精度的代理模型對功能函數進行擬合[11]。Kriging方法從變量相關性和變異性出發,在有限區域內對區域化變量的取值并進行無偏和最優估計,具有很高的擬合精度[12,13]。因此,本文采用Kriging代理模型來擬合緩沖系統的三個功能函數gi(x)(i=1,2,3)。
(2)緩沖系統的功能涉及到η(x),Stroke(x),max(Force(x))三個指標,其參數的可靠性靈敏度是一個多模式耦合問題。在多模式的參數可靠性靈敏度分析中,一次二階矩方法適合于解決線性極限狀態方程這類問題,對緩沖器這類復雜系統進行可靠性靈敏度分析估算誤差較大。Monte Carlo法是對變量的個數、變量的概率分布及極限狀態的形式等均無限制,但該方法收斂緩慢、計算量巨大。重要抽樣法采用重要抽樣函數進行抽樣,能使更多的樣本點落入失效域,從而提高抽樣效率,減少計算量[14,15]。傳統重要抽樣法需要預先計算極限狀態方程的設計點,而多極限狀態方程的設計點難以求解。構建重要抽樣函數的另一種思路是保持抽樣中心在均值點不變,將標準差擴大λ倍來構建重要抽樣函數,該策略不需要求解極限狀態方程的設計點,也不受極限狀態方程的個數限制[16]。因此本文通過該方式構建重要抽樣函數來求解緩沖器參數的可靠性靈敏度。
在Kriging代理模型中,首先需要產生插值試驗點。拉丁超立方抽樣設計(LHS)方法能在預定的抽樣空間抽出較為均勻的實驗點[17],因此本文采用LHS方法產生m個試驗樣本,記為Xt=(Xt1,Xt2,…,Xt9)(t=1,2,…,m),對應的功能函數值為gi(Xt)(t=1,2,…,m;i=1,2,3)。利用功能函數值和試驗點建立Kriging的近似模型,如下式所示

式中C= [c1(X),c2(X),…,c9(X)],α= [α1,α2,…,α9]T,ck(X)(k=1,2,…,9)是預先給定的多項式函數,αk為回歸模型的待定系數,Z(X)是一均值為0、方差為σ2的高斯過程,其協方差如下式所示

式中R為試驗點Xj和Xp的相關函數,它與試驗點之間在空間的位置密切相關。R選取Gaussian函數形式,如下式所示

式中d表征了試驗點Xj、Xp之間的距離,如下式

當j=p時,

從式(13)可以看出,R(Xj,Xp)構成的對角元為1。R是m×m的相關系數陣,形式如下式

假設采用重要抽樣法獲得的樣本為x,通過式(9)得到基于Kriging模型樣本x的預測響應值~gi(x),如下式

式中Gi=[gi(X1),gi(X2),…,gi(Xm)],r(x)為預測點x與m個插值點之間的相關系數構成的列向量,如下式

x的回歸模型系數α采用最小二乘法來求解,其估計式為

對式(9)中Z(X)的方差進行估計,估計值如下式

在式(18)中,R滿足下式

采用極大似然算法求解γi即可建立式(9)所示的Kriging代理模型。
在本文中,通過LHS方法產生500個試驗點,隨機抽取400個試驗點建立緩沖系統三個功能函數g1(x),g2(x),g3(x)基于 Kriging方法的代理模型,將余下的100個試驗點進行模型的預測精度測試,各響應的均方根誤差(RMSE)和最大相對誤差(MRE)均小于3%,滿足工程要求。
對于式(6)所示的可靠性靈敏度,引入重要抽樣函數h(λ,x)(λ>1),則參數的可靠性靈敏度改寫為下式

推導式(20)基于重要抽樣法的無偏估計?^Pf/?θk及其方差Var(?^Pf/?θk)、變異 系數CV(?^Pf/?θk)如下式所示


式中xj為重要抽樣函數h(λ,x)抽取的第j個子樣,N為樣本總數目。
在本文中,變量xk(k=1,2,…,9)是服從N(μk,σ2k)正態分布,重要抽樣函數h(λ,x)可表示為下式

在緩沖器參數可靠性靈敏度分析中,需要對靈敏度結果進行無量綱的歸一化處理,如下式所示

基于Kriging模型和重要抽樣法的緩沖器參數可靠性靈敏度分析結果如圖4所示。

圖4 緩沖器參數的可靠性靈敏度Fig.4 Reliability sensitivities of shock absorber parameters
從圖4可以看出,壓氣面積、壓油面積、主油孔面積、主油孔流量系數的均值和標準差對起落架緩沖系統的可靠性影響都非常大,在緩沖器改設計中要重點關注這些參數。
(1)在多體動力學仿真的基礎上,建立了起落架緩沖器參數的可靠性靈敏度分析模型;
(2)針對緩沖器參數可靠性靈敏度分析中的多模式及隱式功能函數問題,給出了該問題基于Kriging代理模型和重要抽樣法的求解過程;
(3)通過參數可靠性靈敏度分析,確定壓氣面積、壓油面積、主油孔面積、主油孔流量等參數為關鍵參數,為緩沖器改設計提供指導。
[1] 航空航天工業部科學技術委員會.飛機起落架強度設計指南[M].成都:四川科學技術出版社,1989.Science and technology committee of aeronautics and astronautics.Introduction to design for airplane landing gear[M].Chengdu:Sichuan Science and Technology Press,1989.
[2] 《飛機設計手冊》總編委員會.飛機設計手冊(第14冊)[M].北京:航空工業出版社,2003.Aircraft design handbook general editorial board.Aircraft Design Handbook (Volume 14)[M].Beijing:Aviation Industry Press,2003.
[3] Zdravko T,Milan V,Hinko W.Numerical simulation of landing aircraft dynamics[J].Strojarstvo,2009,51(6):657—665.
[4] 張璞,劉向堯,聶宏,等.民用飛機起落架緩沖器選型與對比分析研究[J].航空工程進展,2011,2(3):287—291.Zhang P,Liu XY,Nie H,et al.Research on types and comparisons of landing gears shock absorber for civil airplane[J].Advances in Aeronautical Science and Engineering,2011,2(3):287—291.
[5] 劉小川,馬曉利,孫俠生,等.基于響應面方法的多支柱起落架著陸緩沖性能優化[J].振動工程學報,2010,23(3):305—309.Liu X C,Ma X L,Sun X S,et al.Performance optimization of shock absorber for multi-strut landing gear based on RSM[J].Journal of Vibration Engineering,2010,23(3):305—309.
[6] 侯赤,萬小朋,趙美英.虛擬樣機技術在飛機地面載荷分析中的應用[J].哈爾濱工業大學學報,2009,41(11):134—137.Hou C,Wan X P,Zhao M Y.Application of virtual prototyping technology in the analysis of aircraft ground load[J].Journal of Harbin Institute of Technology,2009,41(11):134—137.
[7] 趙世春,喻天翔,王慧,等.仿真技術在飛機起落架可靠性分析中的應用研究[J].航空工程進展,2011,2(4):444—448.Zhao S C,Yu T X,Wang H.Application of simulation technology on reliability analysis for aircraft landing gear system[J].Advances in Aeronautical Science and Engineering,2011,2(4):444—448.
[8] Wu Y T.Computational methods for efficient structural reliability and reliability sensitivity analysis[J].AIAA J,1994:32(8):1 717—1 723.
[9] Wu Y T,Sitakanta M.Variable screening and ranking using sampling-based sensitivity measures[J].Reliability Engineering and System Safety,2006,91(6):634—647.
[10]Li H S,Lu Z Z,Yuan X K.Nataf transformation based point estimate method[J].Chinese Science Bulletin,2008,53(17):2 586—2 592.
[11]張崎,李興斯.基于Kriging模型的結構可靠性分析[J].計算力學學報,2006,23(2):176—179.Zhang Q,Li X S.Analysis of structural reliability based on Kriging model[J].Chinese Journal of Computational Mechanics,2006,23(2):176—179.
[12]Echard B,Gayton N,Lemaire M.AK-MCS:An active learning reliability method combining Kriging and Monte Carlo Simulation[J].Structural Safety,2011,33:145—154.
[13]Vincent D,Bruno S,Bourinet J M.Reliability-based design optimization using Kriging surrogates and subset simulation[J].Struct.Multidisc Optim,2011,44:673—690.
[14]Melchers R E.Search-based importance sampling[J].Structure Safety,1990(9):117—128.
[15]呂震宙,岳珠峰.多個模式聯合失效的設計點及概率的計算[J].計算力學學報,1998,15(4):442—447.Lu Z Z,Yue Z F,Joint design point and failure probability calculations of multi-mode[J].Chinese Journal of Computational Mechanics,1998,15(4):442—447.
[16]吳斌,歐進萍,張紀剛.結構動力可靠度的重要抽樣法[J].計算力學學報,2001,18(4):478—482.Wu B,Ou J P,Zhang J G.Importance sampling techniques in dynamical structural reliability[J].Chinese Journal of Computational Mechanics,2001,18(4):478—482.
[17]Helton J C,Davis F J.Latin hypercube sampling and the propagation of uncertainty in analysis of complex systems[J].Reliability Engineering and System Safety,2003,81(1):23—69.