陶周亮，方建義，張 梅

(中國直升機設(shè)計研究所，江西 景德鎮(zhèn) 333001)

?

基于LS-DYNA的滑橇起落架落震分析及二次開發(fā)

陶周亮，方建義，張 梅

(中國直升機設(shè)計研究所，江西 景德鎮(zhèn) 333001)

滑橇起落架在著陸過程中通過結(jié)構(gòu)變形來吸收著陸功量而不產(chǎn)生結(jié)構(gòu)破壞，因此迫切需要開發(fā)一套高效可靠的直升機滑橇起落架落震分析方法?；陲@式動力學和接觸算法，在考慮旋翼升力的影響和機身與弓形梁連接點處的彎矩傳遞問題的基礎(chǔ)上，建立滑橇起落架落震分析模型。為提高仿真分析效率，搭建由二次開發(fā)軟件、ANSYS和LS-DYNA組成的落震分析系統(tǒng)。在兩種工況下對滑橇起落架進行了落震仿真分析，仿真結(jié)果與試驗結(jié)果吻合較好，表明了滑橇起落架落震仿真分析方法的有效性。

直升機；滑橇起落架；落震；顯式動力學；二次開發(fā)

0 引言

滑橇起落架結(jié)構(gòu)簡單，易于加工和維護。現(xiàn)代輕型直升機多采用滑橇式起落架[1]，在直升機著陸過程中通過前、后弓形梁的變形來吸收著陸能量，前后弓形梁會產(chǎn)生很大的塑性變形而要求結(jié)構(gòu)不產(chǎn)生破壞。

在滑橇起落架著陸沖擊的過程中，滑橇受到地面的垂直沖擊載荷和地面的摩擦力載荷。弓形梁產(chǎn)生垂直地面的垂向位移和沿地面方向的側(cè)向位移，隨著載荷增大，弓形梁的各個截面依次由彈性狀態(tài)進入塑性狀態(tài)?；疗鹇浼艿闹懶阅苡嬎闶且粋€求解滑橇在動載荷作用下的幾何非線性、材料非線性和狀態(tài)非線性的問題，難度較大且不可能用解析方法進行求解[2]。目前國內(nèi)外主要采用數(shù)值分析方法。貝爾直升機公司的Cheng-Ho Tho等[3]對帶圓角的矩形截面弓形梁進行簡化，完成了滑橇起落架的落震分析。黃生月、張紹儀等[4]采用“位移-增量迭代法”來處理滑橇起落架在墜毀過程中吸收著陸功量這一彈塑性、大變形問題。為防止地面共振，滑橇起落架前弓形梁上通常會安裝減振器，黃生月、彭宗梁等[5]在前面的工作基礎(chǔ)之上進行了帶減振器的滑橇式起落架耐墜毀研究。對于這種復雜的沖擊碰撞過程，以非線性有限元基本理論和接觸算法為基礎(chǔ)的非線性有限元軟件可以得到較好的結(jié)果[6]。本文在LS-DYNA中對滑橇起落架進行落震仿真分析，將計算結(jié)果與試驗對比以驗證準確性。為了克服在設(shè)計過程中模型反復修改帶來的重復仿真分析并提高效率，本文基于LS-DYNA和VC++進行滑橇起落架落震分析的二次開發(fā)，該研究對滑橇起落架的設(shè)計具有非常重要的工程實際應用價值。

1 滑橇起落架落震分析方案

滑橇起落架的落震分析方案由三部分組成：二次開發(fā)軟件、ANSYS和LS-DYNA。二次開發(fā)軟件基于VC++[7]開發(fā)，其功能為調(diào)用ANSYS、修改K文件以及調(diào)用LS-DYNA進行仿真分析。其中LS-DYNA為滑橇起落架落震分析方案的計算核心。

滑橇起落架落震分析流程如圖1所示，通過在二次開發(fā)軟件中輸入建模參數(shù)、實常數(shù)參數(shù)、材料參數(shù)和試驗參數(shù)，自動生成APDL代碼，代碼包含了建模、創(chuàng)建材料、網(wǎng)格劃分、加載、初始速度以及K文件生成等等。通過創(chuàng)建進程的方式調(diào)用ANSYS運行APDL代碼，生成供LS-DYNA執(zhí)行的K文件。軟件對成功生成的K文件進行修改，釋放機身與弓形梁連接點處繞航向的轉(zhuǎn)動自由度。軟件調(diào)用LS-DYNA運行成功修改后的K文件進行滑橇起落架的落震仿真分析。

由于在ANSYS中添加constrained_nodal_rigid_body約束時，無法釋放連接點處的自由度，通常是在LS-PrePost中完成自由度釋放，這也是二次開發(fā)軟件兩次分別調(diào)用ANSYS和LS-DYNA的原因。通過對滑橇起落架落震進行參數(shù)化，二次開發(fā)軟件根據(jù)輸入?yún)?shù)快速進行滑橇起落架落震的建模，極大地提高了落震分析的效率和質(zhì)量，尤其是對于設(shè)計的反復修改導致的迭代分析。

圖1 滑橇起落架落震分析流程圖

2 非線性有限元計算理論

非線性方程的求解方法分兩種：隱式方法和顯式方法。非線性有限元中一般采用Newton-Raphson的隱式求解方法，它是將非線性方程線性化進行求解，通過不斷迭代，在每次迭代步中修改剛度矩陣，最終達到平衡。而顯式方法主要指中心差分法。

2.1 顯式積分

對于非線性方程：

式中：U為位移矩陣，P(t)為外力向量矩陣，M為質(zhì)量矩陣，C為阻尼矩陣，K為剛度矩陣。

LS-DYNA采用顯式中心差分法來進行時間積分[8]，在已知0, …,tn時間步解的情況下，求解tn+1時間步的解，運動方程為：

式中：Fint(tn)為內(nèi)力矢量，H(tn)為沙漏阻力。

tn+1時刻的速度和位移由下面公式求得：

集中質(zhì)量矩陣M為對角矩陣，很容易進行求逆，運動方程的求解是非耦合的，不需要組成總體剛度矩陣和對剛度矩陣求逆，與隱式方法對比，大大節(jié)省了存儲空間和計算時間。

2.2 接觸算法

LS-DYNA處理接觸一般采用三種不同的算法：動力約束法、分配參數(shù)法和對稱罰函數(shù)法。對稱罰函數(shù)法為缺省使用算法[9]，其原理是每一時先檢查各從節(jié)點是否穿透主表面，沒有穿透則不對該從節(jié)點做任何處理。如果穿透，則在該節(jié)點與主表面間，主節(jié)點與從表面引入一個較大的界面接觸力，大小與穿透深度、接觸剛度成正比，稱為罰函數(shù)值。其物理意義相當于在其中放置一系列法向彈簧，限制穿透，接觸力的表達式為：

式中：K為接觸界面剛度，δ為穿透量。

LS-DYNA中接觸的摩擦系數(shù)由下式計算可得：

式中：fd為動摩擦系數(shù)，fs為靜摩擦系數(shù)，β為指數(shù)衰減系數(shù)，Vrel為接觸表面的相對速度。

3 有限元模型

考慮到仿真分析的計算效率，采用梁單元beam161對滑橇起落架的前、后弓形梁和滑管進行建模，而機身幾乎不產(chǎn)生吸能變形，故采用mass166來表示。通過約束constrained_nodal_rigid_body將機身質(zhì)量點和前、后弓形梁的四個連接點進行剛體連接。在沖擊過程中，對單個弓形梁而言，地面支撐力和摩擦力產(chǎn)生的彎矩在兩端傳遞到弓形梁的中部。而剛體連接會限制彎矩的傳遞，需釋放四個連接點繞航向的轉(zhuǎn)動自由度，在四個連接點處建立局部坐標系，把RRFLAG置為1。試驗臺采用殼單元shell163和剛體材料建模，限制其平動和轉(zhuǎn)動的六個自由度。有限元模型如圖2所示。

圖2 滑橇起落架落震有限元模型

在整個落震沖擊的過程中，滑橇起落架是吸收著陸功量的主要部件，在大的著陸功量下，前、后弓形梁會進入塑性變形，以Hughes-Liu非線性梁單元模型對前后弓形梁進行建模，前、后弓形梁總共采用了70個梁單元。為了模擬弓形梁進入塑性的狀態(tài)，弓形梁的材料采用雙線性等向強化(BISO)?？紤]到滑管有可能進入塑性變形狀態(tài)，故采用與弓形梁相同的材料模型，滑管的有限元模型由30個梁單元組成。試驗臺在落震沖擊中產(chǎn)生的變形很小，所以在建模時將試驗臺考慮為剛體，由100個殼單元組成，剛體單元幾乎不參與迭代計算，給剛體劃分網(wǎng)格時適當增加網(wǎng)格數(shù)量并不會增加計算時間，卻有助于提高主從面接觸計算的準確性。

滑管和試驗臺之間的接觸采用點面接觸contact_nodes_to_surface，定義滑管上所有的節(jié)點為從面，試驗臺上所有的節(jié)點為主面，輸出滑橇和試驗臺之間的接觸界面力(rcforc)。

而實際落震過程中，直升機旋翼處于旋轉(zhuǎn)狀態(tài)，其產(chǎn)生的旋翼升力是落震沖擊過程中不可忽視的關(guān)鍵影響因素，旋翼產(chǎn)生的升力如下式所示：

式中：L為旋翼產(chǎn)生的升力，KL為旋翼升力系數(shù)，一般取2/3，m為機身質(zhì)量，g為重力加速度。

在有限元模型施加重力場載荷時，考慮旋翼升力的影響下等效的重力加速度為：

式中：Kres為直升機升力卸載影響系數(shù)[10]，取1/3。

4 計算結(jié)果與試驗對比

在質(zhì)心處定義飛機的有效落震重量，為縮短仿真時間，不考慮自由落體階段，從滑管接觸試驗臺前一刻開始仿真，仿真時間為0.12s。為驗證仿真結(jié)果的準確性，在動摩擦系數(shù)為0.25時，對兩種工況(落震高度為64mm和落震高度為204mm)進行仿真，并與試驗數(shù)據(jù)進行對比。

4.1 工況一(落震高度為64mm)

圖3所示為落震高度為64mm時重心處垂直壓縮量曲線，曲線在同一時刻達到峰值，仿真結(jié)果的峰值略大于試驗結(jié)果。兩條曲線的上升以及下降保持一致，仿真與試驗結(jié)果整體吻合度較好。

圖3 重心處垂直壓縮量曲線(落震高度為64mm)

圖4所示為落震高度為64mm時滑管與試驗臺之間的接觸力曲線，經(jīng)輕微光滑處理的仿真結(jié)果曲線帶有略微的毛刺感。仿真結(jié)果的峰值略大于試驗結(jié)果的峰值，仿真結(jié)果的上升及下降均比試驗結(jié)果要快，這主要是試驗中存在的阻尼影響所致， 除此之外，兩種結(jié)果非常接近。

圖4 滑管與試驗臺之間的接觸力曲線(落震高度為64mm)

4.2 工況二(落震高度為204mm)

圖5所示為落震高度為204mm時重心處垂直壓縮量曲線，在整個落震過程中，仿真結(jié)果的峰值與試驗結(jié)果的峰值接近，兩種結(jié)果的吻合度非常好。

圖6所示為落震高度為204mm時滑管與試驗臺之間的接觸力曲線， 兩者的峰值接近。相對落震高度為64mm的結(jié)果來說，曲線在峰值附近“滯留”的時間較長，這是由于滑橇進入了塑性之后變得相對較“軟”。從整體上來看，兩種結(jié)果吻合度較好。

圖5 重心處垂直壓縮量曲線(落震高度為204mm)

圖6 滑管與試驗臺之間的接觸力曲線(落震高度為204mm)

4.3 結(jié)果對比分析

表1所示為仿真結(jié)果與試驗結(jié)果的對比，從表中數(shù)據(jù)可以看出，工況二下的仿真結(jié)果比工況一略微更貼近于試驗結(jié)果。相對于滑管與試驗臺之間的接觸力峰值，重心處的垂向壓縮量峰值的仿真結(jié)果與試驗結(jié)果的誤差要更小，都在2%以內(nèi)。而滑管與試驗臺之間的接觸力峰值的仿真結(jié)果與試驗結(jié)果的誤差最大值為5.83%，處于合理的范圍內(nèi)?？傮w上，兩種工況顯示仿真結(jié)果與試驗結(jié)果非常接近，表明落震分析方法準確有效。

表1 仿真結(jié)果與試驗結(jié)果對比

5 結(jié)論

1) 建立了滑橇起落架落震仿真模型，模型中考慮了旋翼升力的影響，考慮到弓形梁彎矩的傳遞問題，釋放了機身與弓形梁連接點處繞航向的轉(zhuǎn)動自由度，模型更加符合實際情況。

2) 采用梁單元對弓形梁和滑管進行建模，采用殼單元對試驗臺進行建模，把試驗臺考慮為剛體，在不影響精度的前提下大大提高了仿真分析效率。

3) 仿真結(jié)果與試驗結(jié)果對比顯示，仿真結(jié)果與試驗結(jié)果吻合度很好，表明本文提供的仿真分析方法能有效地對滑橇起落架落震進行計算。

4) 由二次開發(fā)軟件、ANSYS和LS-DYNA組成的滑橇起落架落震分析方案能快速提高落震計算效率及質(zhì)量，具有較大的工程實際應用價值。

[1] 楊俊嚴. 直升機載荷手冊[M]. 北京：航空工業(yè)出版社, 1991.

[2] 王新敏. ANSYS工程結(jié)構(gòu)數(shù)值分析[M]. 北京：人民交通出版社, 2007.

[3] Tho C H, Sparks C E, Sareen A K. efficient helicopter skid landing gear dynamic drop simulation using LS-DYNA[J]. Journal of the American Helicopter Society, 2004, 49(4): 483-492.

[4] 黃生月, 張紹儀, 王志榮. 耐墜毀滑橇式起落架的設(shè)計研究[J]. 中國航空科技文獻,1986.

[5] 黃生月, 彭宗梁, 王仁良, 等. 帶減震器的滑橇式起落架耐墜毀研究[J]. 中國航空科技文獻,1992.

[6] 王俊峰, 張志誼. 基于LS—DYNA的沖擊試驗機碰撞分析[J]. 噪聲與振動控制，2007(6): 10-12.

[7] Prata S. C++ Primer Plus(第五版)中文版[M]. 北京：人民郵電出版社, 2005.

[8] 趙海鷗. LS-DYNA動力學分析指南[M]. 北京：兵器工業(yè)出版社, 2003.

[9] Hallquist J 0. LS-DYNA Theory Manual[Z].1998．

[10] 高澤迥. 飛機設(shè)計手冊：第14分冊[M].北京：航空工業(yè)出版社, 2002.

Analysis and Secondary Development of Skid Landing Gear Shock Based on LS-DYNA

TAO Zhouliang, FANG Jianyi, ZHANG Mei

(China Helicopter Research and Development Institute, Jingdezhen 333001, China)

Skid landing gear absorb landing energy by the way of structural deformation in the landing process without structural failure, a efficient and reliable analysis method of skid landing gear shock is urgently needed. Considering the influence of rotor lift and the transfer problem of bending moment on the joints of the fuselage and the skid based on explicit dynamics and contact algorithm, the analysis model of skid landing gear shock was established. In order to raise the efficiency of simulation analysis, an analysis system of shock was established with secondary development software, ANSYS and LS-DYNA. The shock of skid landing gear was analyzed in two working conditions The simulation results were in good agreement with test results, which showed that the provided simulation analysis method of skid landing gear shock in this paper was effective.

helicopter; skid landing gear; shock; explicit dynamics; secondary development

2015-05-21

陶周亮(1989-)，男，江西九江人，碩士，助理工程師，主要研究方向：起落架設(shè)計。

1673-1220(2015)03-025-04

V214.1+3；V226

A