趙知辛，王 琨，汪 杰，張昌明

（1.陜西理工大學(xué)機(jī)械工程學(xué)院，陜西 漢中723000；2.陜西省工業(yè)自動化重點(diǎn)實(shí)驗(yàn)室，陜西 漢中723000）

1 引言

飛機(jī)起落架是實(shí)現(xiàn)飛機(jī)著陸、滑跑以及降落的重要裝置，并在起降、滑跑過程中用于支撐飛機(jī)的重量［1，2］。因此對起落架進(jìn)行落震仿真分析與結(jié)構(gòu)優(yōu)化設(shè)計(jì)尤為重要。在飛機(jī)飛行過程中，起落架是收在飛機(jī)機(jī)身內(nèi)部的，如果起落架的體積很大，意味著收放空間就要越大，這就會造成起落架的收放困難［3，4］。并且隨著起落架的重量增加，對燃油量的消耗也會增多，實(shí)驗(yàn)數(shù)據(jù)表明，每減少10kg的飛機(jī)重量，就能夠節(jié)省3925kg的燃油量并少排放4噸CO2。因此對起落架進(jìn)行減重是當(dāng)前的研究熱點(diǎn)。國內(nèi)外已經(jīng)有學(xué)者對飛機(jī)起落架進(jìn)行研究，在國外，文獻(xiàn)［5］提出一種輕型緊湊型直升機(jī)著陸支柱的優(yōu)化設(shè)計(jì)方法。文獻(xiàn)［6］提出了一種考慮強(qiáng)度、穩(wěn)定性和撓度限制最小重量的飛機(jī)壁板層復(fù)合材料結(jié)構(gòu)設(shè)計(jì)方法。文獻(xiàn)［7］提出了一種解決航天設(shè)計(jì)問題的拓?fù)鋬?yōu)化算法。在國內(nèi)，文獻(xiàn)［8］建立了飛機(jī)起落架帶鎖撐桿支桿的三維模型，利用Matlab對支桿的截面尺寸進(jìn)行了優(yōu)化。文獻(xiàn)［9］利用Adams軟件對飛機(jī)起落架的支柱進(jìn)行了參數(shù)化建模，并以油孔面積為設(shè)計(jì)變量，分析了油孔面積的變化對緩沖性能的影響。文獻(xiàn)［10］借助Optistruct結(jié)構(gòu)優(yōu)化平臺，利用變密度法，在給出的多種工況下對飛機(jī)起落架的外筒進(jìn)行了拓?fù)鋬?yōu)化設(shè)計(jì)。文獻(xiàn)［11］提出了一種施加攔阻載荷的方法，對艦載無人機(jī)攔阻著艦進(jìn)行剛?cè)狁詈隙囿w動力學(xué)仿真分析。然而，上述學(xué)者僅僅是在給定工況條件下，對其進(jìn)行了優(yōu)化設(shè)計(jì)，工況條件一般是通過試驗(yàn)、與理論計(jì)算得來的，用這兩種方法獲得的工況條件不僅耗時費(fèi)力，而且大都是針對特定機(jī)型而獲得的工況條件，用這種工況條件進(jìn)行優(yōu)化設(shè)計(jì)往往不準(zhǔn)確。

通過Adams軟件建立了飛機(jī)前起落架的動力學(xué)仿真模型，然后對其進(jìn)行落震仿真分析，并利用動力學(xué)分析所獲得的工況條件，利用變密度法與優(yōu)化準(zhǔn)則法以扭力臂的單元相對密度作為設(shè)計(jì)變量、以柔度最小作為目標(biāo)函數(shù)對飛機(jī)前起落架的扭力臂結(jié)構(gòu)進(jìn)行了拓?fù)鋬?yōu)化設(shè)計(jì)。

2 計(jì)算方法

在以往對飛機(jī)前起落架的研究基礎(chǔ)上，利用Adams軟件建立了飛機(jī)前起落架的落震仿真模型，并對其進(jìn)行了動力學(xué)分析，最后基于變密度法與優(yōu)化準(zhǔn)則法對飛機(jī)前起落架扭力臂結(jié)構(gòu)進(jìn)行了拓?fù)鋬?yōu)化設(shè)計(jì)。

2.1 基于變密度法的拓?fù)鋬?yōu)化數(shù)學(xué)模型

連續(xù)體結(jié)構(gòu)拓?fù)鋬?yōu)化模型的本質(zhì)分為兩類：一類是尋求0-1離散變量的組合優(yōu)化，然后用離散變量0和1來描述優(yōu)化模型，采用組合優(yōu)化方法求解拓?fù)鋬?yōu)化問題。對于規(guī)模較小的問題，這種組合優(yōu)化方法擁有較強(qiáng)的全局尋優(yōu)能力，但是當(dāng)連續(xù)體結(jié)構(gòu)拓?fù)鋬?yōu)化的規(guī)模達(dá)到一定程度時，結(jié)果就截然相反了，不僅優(yōu)化問題的求解效率大大降低，而且由于設(shè)計(jì)變量數(shù)量的增多，就會出現(xiàn)“組合爆炸”問題；

另一類方法，變密度法是將離散變量的優(yōu)化問題松弛為連續(xù)變量的優(yōu)化問題，從而使原本離散變量的設(shè)計(jì)模型變?yōu)檫B續(xù)設(shè)計(jì)變量的優(yōu)化模型，這樣就避免了“組合爆炸問題”，因此在拓?fù)鋬?yōu)化設(shè)計(jì)中選用了變密度法。變密度法首先選定設(shè)計(jì)域與非設(shè)計(jì)域，然后將設(shè)計(jì)域進(jìn)行有限元離散，以每個單元介于0、1之間的單元相對密度作為設(shè)計(jì)變量，將一個離散的優(yōu)化問題轉(zhuǎn)變成更易求解的連續(xù)性優(yōu)化問題［12］。在運(yùn)用變密度法對研究對象進(jìn)行拓?fù)鋬?yōu)化設(shè)計(jì)時，需要引入懲罰因子對出現(xiàn)的中間密度值進(jìn)行懲罰，使中間密度值向0/1兩端聚集，在此次優(yōu)化設(shè)計(jì)中采用的材料插值模型為SIMP法，SIMP法用公式表示為：

式中：E（bi）、E1、E0、bi、p-插值之后的彈性模量、孔洞部分材料的彈性模量、實(shí)體部分材料的彈性模量、單元相對密度，取值為0的代表無材料，取值為1的代表有材料和懲罰因子。以單元相對密度作為設(shè)計(jì)變量，結(jié)構(gòu)的最小柔度作為目標(biāo)函數(shù)，基于變密度法SIMP法的拓?fù)鋬?yōu)化數(shù)學(xué)模型可以表示為：

式中：B、C、A、S、K、f、V、V0、bmin、bmax-設(shè)計(jì)變量、結(jié)構(gòu)柔度、載荷矢量、位移矢量、結(jié)構(gòu)剛度矩陣、保留的體積分?jǐn)?shù)、優(yōu)化后的體積、結(jié)構(gòu)的初始體積、設(shè)計(jì)變量的下限值和設(shè)計(jì)變量的上限值。

2.2 優(yōu)化算法

在拓?fù)鋬?yōu)化問題中，用來求解的常用兩種算法是優(yōu)化準(zhǔn)則法與移動漸近線法。優(yōu)化準(zhǔn)則法是以目標(biāo)函數(shù)和約束函數(shù)的一階導(dǎo)數(shù)之商來更新單元密度變量，此法最大的特點(diǎn)是適用于設(shè)計(jì)變量比較多與約束條件較少的情況，且迭代次數(shù)少，計(jì)算效率高，通用性較強(qiáng)［13，14］。移動漸近線法是利用設(shè)計(jì)點(diǎn)目標(biāo)函數(shù)與其一階導(dǎo)數(shù)構(gòu)建一個近似凸函數(shù)去逼近實(shí)際的隱函數(shù)，此法適用于求解較為平滑的非線性優(yōu)化問題，缺點(diǎn)是尋找近似函數(shù)比較困難。此次拓?fù)鋬?yōu)化問題中，選擇更易求解的優(yōu)化準(zhǔn)則法。

式中：w、n、η、λ-正的移動步長、迭代次數(shù)、阻尼系數(shù)和拉格朗日乘子。

3 落震仿真模型的建立

根據(jù)某軍用A型飛機(jī)前起落架的實(shí)際尺寸，建立了前起落架的三維模型，并對此模型進(jìn)行了合理的簡化，然后將模型保存為X_T格式并將其導(dǎo)入到Adams中，接著對其添加必要的運(yùn)動副與力的加載。用一個小球來模擬機(jī)身的重量，并將小球與機(jī)身用固定副相連接，因?yàn)榛钊麠U與外筒之間既可以相互轉(zhuǎn)動，又可以相互軸向移動，故在二者之間添加圓柱副；考慮到防扭臂可相對于外筒及活塞桿轉(zhuǎn)動，通過旋轉(zhuǎn)副分別與之連接；起落架在著陸過程中，斜撐桿處于鎖死狀態(tài)，故使用固定副與外筒連接；在落震仿真過程中，考慮到輪胎有一定的轉(zhuǎn)速，故在輪胎與輪軸之間通過旋轉(zhuǎn)副相連接；在跑道與機(jī)輪之間定義了一個接觸，并用impact函數(shù)來計(jì)算接觸力；在起落架緩沖系統(tǒng)中主要是對空氣彈簧力、油液阻尼力與結(jié)構(gòu)限制力的加載，由于在Adams中無法模擬油液與空氣，故在外筒軸線上一點(diǎn)MARKER_7與活塞桿最頂端MARKER_8之間通過添加單作用力來分別模擬空氣彈簧力、油液阻尼力與結(jié)構(gòu)限制力。根據(jù)空氣彈簧力的公式，在空氣彈簧力的函數(shù)定義窗口輸入其表達(dá)式：initial_air_pressure*pneumatic_area*（initial_air_volume/（initial_air_volume-pneumatic_area*（uppercylinder_length-DM（MARKER_7，MARKER_8））））**polytropic_exponent，由于油液阻尼力與活塞桿的運(yùn)動方向、運(yùn)動速度都有關(guān)系，因此利用VR函數(shù)與IF函數(shù)模擬了油液阻尼力，部分函數(shù)表達(dá)式為：IF（VR（MARKER_7，MARKER_8）：表達(dá)式1，0，表達(dá)式3），其中表達(dá)式1與表達(dá)式3為油液阻尼力的公式；最后對結(jié)構(gòu)限制力進(jìn)行模擬，從結(jié)構(gòu)限制力的公式可以看出，結(jié)構(gòu)限制力是一個分段函數(shù)，因此利用IF函數(shù)模擬了結(jié)構(gòu)限制力，其函數(shù)表達(dá)式為：

stopper_stiffness*（uppercylinder_length-DM（MARKER_7，MARKER_8））*IF（uppercylinder_length-DM（MARKER_7，MARKER_8）：1，0，0），如圖1所示飛機(jī)前起落架的落震仿真模型。

圖1 前起落架的落震仿真模型Fig.1 Simulation Model of Nose Landing Gear

4 仿真分析

在飛機(jī)前起落架著陸之前需做一個靜平衡，讓飛機(jī)前起落架在前兩秒時間里是靜力學(xué)分析，兩秒之后飛機(jī)前起落架才開始下落，初始仿真參數(shù)由某軍用A型飛機(jī)前起落架設(shè)計(jì)參數(shù)通過計(jì)算得到［15］，部分?jǐn)?shù)據(jù)如表1所示。然后設(shè)定仿真時間為4s，迭代子步為400步，最大壓縮行程為600mm，如圖2所示為飛機(jī)前起落架空氣彈簧力隨時間的變化曲線。

表1 初始仿真參數(shù)及取值Tab.1 Initial Simulation Parameters and Values

圖2 空氣彈簧力隨時間變化規(guī)律Fig.2 Variation Law of Air Spring Force with Time

圖2 中可以看出，由于前2.74s屬于飛機(jī)空中降落階段，此時空氣彈簧力是一個定值約為17792 N；在2.74s時，此時飛機(jī)機(jī)輪剛好接觸地面，緩沖器開始壓縮，經(jīng)過0.32s后空氣彈簧力達(dá)到最大值約為84293 N，之后緩沖器開始回彈，空氣彈簧力逐漸下降，在3.58s時達(dá)到最小值約為37600 N，完成了一次循環(huán)往返，之后空氣彈簧力震蕩幅度越來越小，最終趨于穩(wěn)定，如圖3所示為緩沖器行程隨時間變化曲線。

從圖3中可以看出，在前2.74s屬于飛機(jī)空中降落階段，此時緩沖器并沒有開始壓縮，緩沖器的行程為0，在2.74s后緩沖器開始壓縮且在3.06s時緩沖器的行程達(dá)到最大值，最大值約為439 mm；隨后緩沖器開始釋放能量，行程開始回彈，由圖可知當(dāng)行程回彈到約239 mm時，完成了第一次往返行程，又開始第二次的循環(huán)往返行程時間約為0.84s。如圖4所示，為油液阻尼力隨時間變化曲線。

圖3 緩沖器行程隨時間變化規(guī)律Fig.3 Variation of Buffer Stroke with Time

圖4 油液阻尼力隨時間變化規(guī)律Fig.4 Oil Damping Force Changing with Time

從圖4可以看出，曲線在2.74s后油液阻尼力迅速增大，與緩沖器行程、空氣彈簧力增大的時刻一致，在t=2.78s時油液阻尼力達(dá)到最大值為220300 N，之后迅速減小，成正負(fù)交替震蕩，最終趨于零，說明外筒與活塞桿之間相對運(yùn)動速度迅速趨于減小。如圖5所示，為緩沖器軸向力的變化規(guī)律。

圖5 緩沖器軸向力的變化規(guī)律Fig.5 Variation Law of Buffer Axial Force

在圖5中，可以看出t=2.18s時，最大軸向力為Py=123770 N，且軸向力最終趨于穩(wěn)定，在拓?fù)鋬?yōu)化設(shè)計(jì)中選取最大軸向力作為危險(xiǎn)工況對扭力臂進(jìn)行拓?fù)鋬?yōu)化。最大起轉(zhuǎn)載荷為［16］：

式中：Ay、Ax-最大起轉(zhuǎn)垂直載荷和最大起轉(zhuǎn)水平載荷。

由此計(jì)算出緩沖器所受扭矩為Mt=Ax*r=6280N*m，r-扭力臂與活塞桿連接處軸線到活塞桿軸線的距離。

5 拓?fù)鋬?yōu)化設(shè)計(jì)

在進(jìn)行拓?fù)鋬?yōu)化之前，先對前起落架進(jìn)行靜力學(xué)分析，前起落架材料為40CrMnSiMoVA，強(qiáng)度極限為1760 MPa，彈性模量為200 GPa，泊松比為0.3。由于上扭力臂與外筒之間是相互鉸接在一起，而外筒是固定不動的，故在上扭力臂兩端施加固定約束，然后在活塞桿上施加大小為6280N*m的扭矩，通過有限元方法得到上下扭力臂的變形及應(yīng)力分布情況，如圖6所示。

圖6 優(yōu)化前扭力臂變形、應(yīng)力圖Fig.6 Deformation and Stress Diagram of Torsion Arm Before Optimization

從（a）圖中可以看到扭力臂的最大變形量約為0.84mm，變形量很小，且從（b）圖中可以看到，扭力臂的最大應(yīng)力約為769MPa，最大應(yīng)力值出現(xiàn)在扭力臂的邊角部分，中部區(qū)域應(yīng)力值較小，遠(yuǎn)小于材料的強(qiáng)度極限，這說明中間區(qū)域扭力臂的剛度冗余很大，具有很大的優(yōu)化設(shè)計(jì)空間。考慮到扭力臂的優(yōu)化設(shè)計(jì)空間較大，進(jìn)一步對扭力臂進(jìn)行拓?fù)鋬?yōu)化設(shè)計(jì)，在拓?fù)鋬?yōu)化設(shè)計(jì)過程中需要選擇設(shè)計(jì)域與非設(shè)計(jì)域，此次優(yōu)化設(shè)計(jì)過程中設(shè)計(jì)域選擇為除去扭力臂吊耳部分的中間區(qū)域，體積保留分?jǐn)?shù)選擇為0.7，優(yōu)化后的拓?fù)鋬?yōu)化結(jié)果，如圖7所示。

從圖7中可以看到，拓?fù)鋬?yōu)化結(jié)果比較清晰，扭力臂中間部分形成了明顯的開槽區(qū)域，按照拓?fù)鋬?yōu)化生成的結(jié)果，在三維軟件中，重新對起落架進(jìn)行精確建模。

圖7 拓?fù)鋬?yōu)化結(jié)果云圖Fig.7 Cloud Chart of Topology Optimization Results

接下來對重新建立后的模型，進(jìn)行強(qiáng)度設(shè)計(jì)校核，利用有限元方法得到了扭力臂優(yōu)化后的變形與應(yīng)力分布結(jié)果，如圖8所示。

圖8 優(yōu)化后扭力臂變形、應(yīng)力圖Fig.8 Deformation and Stress Diagram of Torsion Arm after Optimization

從圖8中可以看到優(yōu)化后的扭力臂最大變形量約為0.88mm，比優(yōu)化前的變形量僅增加了0.04mm；優(yōu)化后的最大應(yīng)力值約為842MPa，比優(yōu)化前的應(yīng)力值增加了73MPa，增加的百分比約為9.4%，且優(yōu)化后的應(yīng)力值仍在材料的強(qiáng)度極限范圍內(nèi)，故優(yōu)化后的結(jié)構(gòu)滿足強(qiáng)度要求；優(yōu)化后扭力臂的質(zhì)量為7.96kg，優(yōu)化前扭力臂的質(zhì)量為10.68kg，優(yōu)化后扭力臂的質(zhì)量比優(yōu)化前的質(zhì)量減少了2.72kg，減少的百分比為25.4%，輕量化效果明顯。

6 結(jié)論

利用Adams軟件建立了飛機(jī)前起落架的落震仿真模型，并進(jìn)行仿真分析，得到了飛機(jī)前起落架緩沖系統(tǒng)軸向載荷及行程的變化規(guī)律，最后在動力學(xué)分析的基礎(chǔ)上基于變密度法與優(yōu)化準(zhǔn)則法對飛機(jī)起落架進(jìn)行了拓?fù)鋬?yōu)化設(shè)計(jì)，得到如下結(jié)果：

（1）在降落階段，此時空氣彈簧力是一個定值，當(dāng)飛機(jī)著陸后緩沖器開始壓縮，空氣彈簧力迅速增大，增大到峰值后，緩沖器開始回彈，空氣彈簧力迅速下降，之后緩沖器重復(fù)此過程，空氣彈簧力來回震蕩，且震蕩幅度越來越小，最終趨于穩(wěn)定；緩沖器行程在空中降落階段為0，此后變化規(guī)律與空氣彈簧力變化規(guī)律一致；在空中降落階段，油液阻尼力也為0，在之后的0.3s內(nèi)，油液阻尼力震蕩幅度很大，此后震蕩幅度很小，直至趨于穩(wěn)定；

（2）在動力學(xué)分析的基礎(chǔ)上，基于變密度法與優(yōu)化準(zhǔn)則法以單元相對密度作為設(shè)計(jì)變量，以結(jié)構(gòu)柔度最小作為目標(biāo)函數(shù)，在滿足強(qiáng)度與剛度的條件下，對扭力臂結(jié)構(gòu)進(jìn)行了拓?fù)鋬?yōu)化設(shè)計(jì)，得到的拓?fù)鋬?yōu)化結(jié)果清晰明了且易于加工，優(yōu)化后扭力臂質(zhì)量減小了2.72kg，減少的百分比為25.4%，輕量化設(shè)計(jì)效果明顯。