祝世興，李佩窈，祝恒佳，魏戩

(1.中國民航大學(xué)航空工程學(xué)院，天津 300300；2.昆明理工大學(xué)機(jī)電工程學(xué)院，云南昆明 650000)

0 前言

緩沖器是所有現(xiàn)代飛機(jī)不可或缺的部分，目前油氣緩沖器在現(xiàn)代飛機(jī)上應(yīng)用最為廣泛[1-2]。雙腔式油氣緩沖器包含低壓腔和高壓腔2個(gè)空氣腔，擁有單腔緩沖器優(yōu)點(diǎn)的同時(shí)在應(yīng)對粗糙路況時(shí)的吸能特性也得到明顯改善。雙腔式油氣緩沖器剛度特性對改善飛機(jī)在沖擊振動(dòng)、高速滑跑時(shí)的運(yùn)行品質(zhì)，提升乘客的舒適度，保護(hù)飛機(jī)零部件以及儀表設(shè)備等方面影響較大。目前，緩沖器的剛度特性一般采用激振試驗(yàn)的方法獲得，但是由于試驗(yàn)自身的滯后性、周期長以及成本高等缺點(diǎn)，采用仿真模擬的方法精確地模擬其動(dòng)態(tài)特性成為研究趨勢[3-5]。雙腔式油氣緩沖器內(nèi)部結(jié)構(gòu)復(fù)雜，一般的仿真模擬并沒有考慮到流體在工作過程中的流動(dòng)損失對緩沖器整體剛度特性的影響。因此，為更接近緩沖器實(shí)際工作情況，采用流固耦合仿真方法對雙腔式油氣緩沖器整體的剛度特性進(jìn)行研究是非常重要的。

早期的雙腔式緩沖器多應(yīng)用于直升機(jī)上，并且大部分是對起落架整體進(jìn)行研究，而針對雙腔式緩沖器的單獨(dú)研究較少。文獻(xiàn)[6]中研究了符合直升機(jī)地面共振穩(wěn)定性要求的起落架剛度及阻尼的優(yōu)化設(shè)計(jì)方法。文獻(xiàn)[7]中利用動(dòng)力學(xué)分析軟件ADAMS對起落架進(jìn)行了落震動(dòng)力響應(yīng)的模擬與分析。文獻(xiàn)[8]中對大柔性飛機(jī)起落架緩沖器參數(shù)進(jìn)行了優(yōu)化設(shè)計(jì)。文獻(xiàn)[9]中利用ADAMS建立了柔性活塞桿的緩沖器動(dòng)力學(xué)模型并進(jìn)行了落震仿真模擬，探討活塞桿剛度對緩沖器吸能特性的影響。大多數(shù)研究中均將起落架簡化為二自由度彈簧質(zhì)量塊模型，將緩沖器簡化為所需力的形式加載到質(zhì)量塊上，但這并未考慮到起落架結(jié)構(gòu)的彈性特性和緩沖器的空氣彈簧和阻尼的非線性特性。因此，文獻(xiàn)[10]中利用哈密頓變分原理，建立了彈性結(jié)構(gòu)的起落架動(dòng)力學(xué)模型，以彌補(bǔ)上述缺點(diǎn)。對于直升機(jī)起落架的防“地面共振”設(shè)計(jì)和抗墜毀能力，文獻(xiàn)[11]中設(shè)計(jì)了一種單級雙腔緩沖器，其中高壓腔與低壓腔分離，并且保留了單級雙腔緩沖器的優(yōu)點(diǎn)，串聯(lián)雙腔緩沖器目的是將飛機(jī)的起飛著陸和耐墜毀緩沖行程進(jìn)行分離。文獻(xiàn)[12]中利用等效線性化方法，分析了串聯(lián)雙腔緩沖器在多種頻率下的緩沖器等效剛度和當(dāng)量阻尼分析值。文獻(xiàn)[13]中利用流固耦合及動(dòng)網(wǎng)格方法，對緩沖器的阻尼特性進(jìn)行了仿真。但是到目前為止，仍缺少考慮緩沖器復(fù)雜的流固耦合問題尤其是緩沖器工作過程中流體的流動(dòng)損失有限元模擬問題對緩沖器剛度特性影響的相關(guān)算例與研究。

本文作者以某大型民航客機(jī)主起落架為例，建立適用于ABAQUS的雙腔油氣緩沖器仿真模型。利用工程估算方法模擬緩沖器沖擊過程并與ABAQUS中的仿真結(jié)果進(jìn)行比對，驗(yàn)證ABAQUS中模型的可用性；采用正弦激勵(lì)，分析緩沖器隨幅值和頻率變化時(shí)的剛度特性；采用隨機(jī)激勵(lì)，模擬實(shí)際工況，分析緩沖器的動(dòng)態(tài)特性，從而為雙腔式緩沖器結(jié)構(gòu)性能的優(yōu)化提供依據(jù)，并為進(jìn)一步利用ABAQUS軟件開發(fā)緩沖器樣機(jī)設(shè)計(jì)奠定基礎(chǔ)。

1 雙腔緩沖器模型描述

如圖1所示，雙腔式緩沖器低壓腔與高壓腔通常充入一定初始?xì)鈮旱母稍锏獨(dú)饣蚩諝猓ぷ鲿r(shí)活塞運(yùn)動(dòng)使油腔中的油液流向低壓腔壓縮氣體，同時(shí)吸收能量。油液在外力的作用下，流經(jīng)油孔耗散掉起落架一部分能量。當(dāng)?shù)蛪呵坏目諝鈮毫εc油液阻尼力增加到等于高壓腔初始壓力時(shí)，位于高低壓腔之間的活塞開始移動(dòng)，高壓腔開始被壓縮[14]。壓縮緩沖器使得活塞上移稱為“正行程”。經(jīng)過初始沖擊后，空氣腔氣體作為彈性體開始膨脹釋放能量，迫使油液又反流回油腔，活塞向下移動(dòng)，此過程稱為“反行程”。經(jīng)過一正一反2個(gè)行程完成一個(gè)循環(huán)，經(jīng)多個(gè)循環(huán)，可耗散掉絕大部分飛機(jī)降落時(shí)產(chǎn)生的能量，進(jìn)一步達(dá)到減震的效果。飛機(jī)的著陸撞擊時(shí)間在0.6～0.9 s之間[15]。正反行程總時(shí)間不超過0.8 s。

圖1 緩沖器結(jié)構(gòu)

2 緩沖器數(shù)學(xué)建模

緩沖器軸向力Fs是由空氣彈簧力Fa、油液阻尼力Fh、緩沖器結(jié)構(gòu)限制力Fl以及內(nèi)部摩擦力Ff構(gòu)成[16]。

Fs=Fa+Fh+Fl+Ff

(1)

結(jié)合氣體多變方程式可以得到空氣彈簧力的力學(xué)方程為

(2)

利用MATLAB擬合緩沖行程和活塞行程為

s=8.582×104t5-1.064×105t4+4.769×

104t3-1.025×104t2+1 180t

(3)

s′=1.276×104t4-1.203×104t3+3 180t2-124t

(4)

油液阻尼力是油液流過油孔時(shí)，油孔兩側(cè)存在壓力差而產(chǎn)生的[17]。假設(shè)油液為不可壓縮流體，由Bernoulli方程和連續(xù)性方程[18]可以得到油液阻尼力為

(5)

其中：Ah是緩沖器有效壓油面積；Aoil是油孔截面積；Cd是油孔縮流系數(shù)；ρ是油液密度。

緩沖器內(nèi)部庫侖摩擦力假定為一個(gè)由緩沖器彎曲位移引起的內(nèi)部摩擦力和內(nèi)部壓力引起的內(nèi)部摩擦力的合力，表達(dá)式為

Ff=Ff1+Ff2

(6)

緩沖器彎曲位移引起的內(nèi)部摩擦力為

(7)

其中：μb是庫侖摩擦因數(shù)；Nu和NL分別是活塞對缸體的作用力與柱塞對缸體的作用力。

緩沖器內(nèi)部壓力引起的內(nèi)部摩擦力為

(8)

其中：Db是軸套直徑；Hb是軸套高度。

緩沖器的設(shè)計(jì)要求指出:緩沖器系統(tǒng)應(yīng)適當(dāng)保留某些行程余量，吸收給定過載下的使用功量[1]。只有在極端條件下，才允許緩沖系統(tǒng)吸收最大的功量。正常使用情況下,緩沖器結(jié)構(gòu)限制力為0，只有在過度伸長與極限壓縮的情況下才考慮結(jié)構(gòu)限制力。文中仿真緩沖器正常的使用情況，所以結(jié)構(gòu)限制力為0。

緩沖器的模型參數(shù)如表1所示。

表1 緩沖器模型參數(shù)

3 有限元模型建立

緩沖器在工作過程中，大致分為兩類接觸問題：第一類是缸體與活塞和柱塞的摩擦接觸，ABAQUS中對于此類接觸問題常采用拉格朗日法進(jìn)行計(jì)算；第二類是氮?dú)夂陀鸵号c缸體、活塞和柱塞的接觸，應(yīng)采用流固耦合的方法進(jìn)行計(jì)算。對于流固耦合問題，ABAQUS中通常采用歐拉-拉格朗日(CEL)接觸算法來模擬流體與固體邊界處高速的相互作用。CEL方法中固體采用拉格朗日單元、流體采用歐拉單元，歐拉單元可以克服流體大變形下網(wǎng)格嚴(yán)重扭曲的現(xiàn)象，并且與拉格朗日單元不同之處在于:拉格朗日單元隨材料移動(dòng)而移動(dòng)，歐拉單元是材料在單元內(nèi)流動(dòng)。如圖2所示。

圖2 連續(xù)體變形

在默認(rèn)情況下，歐拉體內(nèi)部是沒有任何材料的，因此采用CEL方法的關(guān)鍵在于材料的填充，應(yīng)在Load模塊下對相應(yīng)的歐拉體區(qū)域進(jìn)行操作。該模型材料的屬性定義如表2所示，其中Eos: US-UP為Mie-Grüneisen狀態(tài)方程中定義流體的材料參數(shù)。

表2 油氣材料屬性定義

在通常情況下，建立油氣緩沖器有限元模型時(shí)需要考慮氮?dú)馇?、油腔、柱塞、活塞、缸體等，為提升計(jì)算速度，應(yīng)將模型進(jìn)行必要的簡化。假定油液密封性良好，將密封圈簡化為柱塞的一部分，不予以考慮,去除不必要的轉(zhuǎn)角與溝槽,忽略緩沖器自身重力的影響等。

該模型中，緩沖器柱塞與缸體采用C3D8R六面體單元，氣體與液體采用EC3D8R單元，進(jìn)行材料的填充。利用ABAQUS/CAE建立緩沖器三維實(shí)體模型，而后建立離散單元模型，離散剛體單元有41 339個(gè)，其中C3D8R有26 438個(gè)、EC3D8R有14 901個(gè)，如圖3所示。

圖3 網(wǎng)格劃分

應(yīng)用ABAQUS/Explicit模塊進(jìn)行動(dòng)態(tài)分析，緩沖器有限元模型計(jì)算過程分為2個(gè)部分：(1)在特定的頻率與幅值下給予上端一個(gè)沖擊位移，模擬緩沖器的沖擊過程，輸出隨時(shí)間變化的緩沖器軸向力曲線；(2)給予緩沖器不同的動(dòng)態(tài)激勵(lì)，模擬它的運(yùn)動(dòng)狀態(tài)與工作特性。

4 模型仿真與結(jié)果分析

有限元模型建立后，利用軟件按照設(shè)定的要求進(jìn)行仿真分析，直到設(shè)定的仿真時(shí)間為止。緩沖器的正反行程總時(shí)間不超過0.8 s，模型設(shè)定的仿真時(shí)間到0.5 s為止。

4.1 沖擊過程

給予緩沖器上端240 mm的垂向位移，模擬緩沖器沖擊過程。圖4和圖5分別為壓縮過程中緩沖器的內(nèi)部流速與內(nèi)部液體壓力。

圖4 緩沖器的內(nèi)部流速

圖5 緩沖器內(nèi)部液體壓力

基于理論基礎(chǔ)建立數(shù)學(xué)模型，采用MATLAB軟件編寫相應(yīng)的計(jì)算程序，得到空氣彈簧力、油液阻尼力、緩沖器軸向力隨時(shí)間變化的曲線，并將ABAQUS仿真結(jié)果與MATLAB結(jié)果進(jìn)行比對，分析二者之間的誤差。

緩沖器正行程時(shí)間約為0.23 s，圖6所示是空氣彈簧力的理論值與仿真值隨時(shí)間變化曲線?？芍涸趖=0.13 s之前只有低壓腔被壓縮，t=0.13 s后高壓腔處于壓縮階段;在t=0.19 s后，由于空氣腔的氣體在高壓狀態(tài)下液化導(dǎo)致氣體溫度升高，內(nèi)能增加并伴隨了一定的動(dòng)能損失，所以在t=0.19 s之后理論值與仿真值產(chǎn)生了較大誤差。

圖6 空氣彈簧力隨時(shí)間變化曲線 圖7 油液阻尼力隨時(shí)間變化曲線

圖7所示為油液阻尼力理論值與仿真值隨時(shí)間變化曲線?？芍河鸵鹤枘崃﹄S時(shí)間的增加先增加后減小，這是因?yàn)榫彌_器內(nèi)部回油孔前后壓差處于先增加后減小的趨勢，使得油液流速處于先增加后減慢的情況。

數(shù)值仿真和有限元仿真的誤差可以通過均方根值來定義：

(9)

式中：m表示離散力信號的數(shù)據(jù)個(gè)數(shù)；FNUM和FSIM分別表示數(shù)值仿真和有限元仿真的輸出力。

表3所示為空氣彈簧力和油液阻尼力的數(shù)值仿真和有限元仿真之間的誤差大小，可以看出：空氣彈簧力的誤差較大，但由于空氣彈簧力的平均值為17 609 N，因此誤差是可以接受的，驗(yàn)證了文中建立的緩沖器有限元模型的可用性。

表3 不同仿真模型之間的誤差 單位：N

圖8所示為緩沖器軸向力理論值與仿真值隨時(shí)間變化的曲線，可得：理論值與仿真值相對誤差為0.96%，從而驗(yàn)證了仿真模型的準(zhǔn)確性。

圖8 緩沖器軸向力隨時(shí)間變化曲線

經(jīng)比對，空氣彈簧力、油液阻尼力和軸向力仿真值與理論值誤差均小于1%，證明了該模型的可用性，為后續(xù)研究緩沖器特性提供了參考。

4.2 不同激振對緩沖器等效動(dòng)剛度影響

在緩沖系統(tǒng)運(yùn)行過程中，外界環(huán)境溫度與充氣壓力的改變會直接影響緩沖系統(tǒng)的性能。初始擾動(dòng)的大小會影響起落架的振動(dòng)幅值，飛機(jī)裝載不同會引起機(jī)體轉(zhuǎn)動(dòng)慣量的變化，從而導(dǎo)致起落架固有頻率的改變[19]。因此，本文作者給予緩沖器不同位移激勵(lì)方式，改變頻率與幅值，觀察其變化情況。

圖9所示為沖擊激勵(lì)下緩沖器內(nèi)部流體壓力，顯示了緩沖器在工作過程中油液流動(dòng)時(shí)的噴濺現(xiàn)象。由于空氣壓縮的速度很快，以至于在壓縮過程中空氣所產(chǎn)生的熱量來不及向外界傳遞，因此空氣壓縮的過程可近似看成絕熱過程。但油液高速流動(dòng)時(shí)會發(fā)生噴濺現(xiàn)象，與低壓腔的氣體進(jìn)行強(qiáng)烈的熱交換，因此空氣壓縮的過程為多變過程。圖9也從側(cè)面說明了在壓縮過程中低壓腔的多變過程。

圖9 沖擊激勵(lì)下緩沖器內(nèi)部流體壓力

設(shè)對緩沖器施加的位移激勵(lì)為簡諧激勵(lì)，表達(dá)式為

x(t)=H0+Hsinωt

(10)

其中：ω為激振力的圓頻率，ω=2πf；f為激勵(lì)頻率；H為激勵(lì)幅值；H0為初始壓縮量。

根據(jù)雙腔式磁流變減振器的軸向力與位移功量圖，結(jié)合式(11)求得緩沖器的等效動(dòng)剛度K為

(11)

假定緩沖器的振動(dòng)幅值保持12 mm不變，緩沖器等效動(dòng)剛度隨振動(dòng)頻率的變化如圖10所示?？梢钥吹剑弘S頻率的增加等效動(dòng)剛度降低，當(dāng)頻率超過8 Hz時(shí)，等效動(dòng)剛度降低速度減慢。

圖10 頻率對等效動(dòng)剛度的影響

假定緩沖器的振動(dòng)頻率保持2 Hz不變，以不同振幅的正弦位移信號作為輸入，圖11所示為緩沖器在激振幅值為6 mm時(shí)的力-位移遲滯曲線，其中位移為正時(shí)表示壓縮，為正行程；位移為負(fù)時(shí)表示拉伸，為負(fù)行程。圖12所示為幅值對等效動(dòng)剛度的影響,可以看出：幅值與等效動(dòng)剛度呈現(xiàn)嚴(yán)重的非線性特性，隨振動(dòng)幅值增加等效動(dòng)剛度減?。慌c頻率對等效動(dòng)剛度的影響相比，幅值的變化對等效動(dòng)剛度的影響較大。

圖11 激振幅值為6 mm時(shí)的力-位移曲線 圖12 幅值對等效動(dòng)剛度的影響

然而，正弦激勵(lì)在實(shí)際工作過程中很少遇到，作為隔振元件的油氣緩沖器被廣泛應(yīng)用于隨機(jī)振動(dòng)的環(huán)境中。隨機(jī)信號可以通過離散傅里葉變換分解為一系列具有不同頻率和幅值的簡諧波。諧波疊加法基本原理是將路面的不平度表示為大量具有隨機(jī)相位的正弦或余弦諧波之和[20]。

在ABAQUS中，可以考慮采用周期形幅值曲線來近似代替隨機(jī)激勵(lì)，以產(chǎn)生偽隨機(jī)激勵(lì)曲線。周期形幅值曲線用傅里葉級數(shù)表示，表達(dá)式為

(12)

式中：N為傅里葉級數(shù)項(xiàng)的個(gè)數(shù)；t0為起始時(shí)刻；A0為初始幅值；An為cos項(xiàng)的系數(shù)；Bn為sin項(xiàng)的系數(shù)。

取A0=0.01、t0=0，頻率與幅值按照表4取值。

表4 偽隨機(jī)激勵(lì)曲線幅值與頻率部分取值

圖13所示為偽隨機(jī)位移激勵(lì)時(shí)域曲線，圖14所示為0.2～1 s緩沖器輸出的軸向力。選取此時(shí)間段是整個(gè)仿真過程中的隨機(jī)選擇，以方便觀察軸向力變化趨勢。由圖14可以看出：有限元方法輸出緩沖器軸向力略小于數(shù)值方法輸出的緩沖器軸向力，原因是油液與筒壁之間產(chǎn)生摩擦從而引起沿程損失，但總體趨勢一致。

圖13 偽隨機(jī)位移激勵(lì)時(shí)域曲線 圖14 0.2～1 s內(nèi)的軸向力

5 結(jié)論

通過建立起落架雙腔緩沖器的數(shù)學(xué)模型與有限元模型，模擬不同工況下的緩沖器工作特性，得到以下結(jié)論：

(1)與數(shù)學(xué)模型相比，基于流固耦合有限元模型可以較好地處理油液在流動(dòng)過程中產(chǎn)生的壓力損失，與緩沖器的實(shí)際工況更為接近；有限元仿真分析結(jié)果與數(shù)學(xué)模型計(jì)算結(jié)果在大位移沖擊作用下相對誤差均小于1%，進(jìn)一步說明了有限元模型的可行性；

(2)在正弦激勵(lì)工況下，隨激振頻率與幅值的改變，緩沖器等效動(dòng)剛度表現(xiàn)出分段非線性動(dòng)力學(xué)特性，隨頻率與幅值的增加等效動(dòng)剛度減小，并且幅值的改變對等效動(dòng)剛度產(chǎn)生的影響較大；

(3)在隨機(jī)激勵(lì)工況下，由于油液的流動(dòng)損失導(dǎo)致了仿真結(jié)果小于理論結(jié)果，但仿真結(jié)果的總體趨勢與理論結(jié)果基本保持一致，吻合度較高。