月球着陆器着陆腿非线性有限元建模与仿真

梁东平 柴洪友 曾福明

(北京空间飞行器总体设计部，北京100094)

月球着陆器着陆缓冲分系统的主要功能是缓冲着陆器月面软着陆时的冲击载荷，保证结构和有效载荷不因着陆冲击而破坏，并吸收着陆过程中的动能，使着陆器安全着陆到月面上.由于着陆冲击试验成本高、周期长、难度大，只能进行有限次数的特定工况地面模拟验证试验，因此必须开展大量着陆冲击动力学仿真分析，验证着陆缓冲分系统的性能，获得相关性能参数，从而为着陆缓冲分系统的优化设计和试验验证提供指导，并为着陆器结构和仪器设备的设计提供输入条件［1-4］.

着陆缓冲分系统着陆冲击仿真大多数是基于多刚体动力学模型，多刚体模型比较简单，计算效率高，但没有考虑结构柔性［5-7］.为研究结构柔性对着陆腿性能的影响，通常采用部件模态综合的方法.通过有限元分析软件生成着陆腿的模态中性文件，并导入多体动力学软件中，建立柔性体模型［8］.由于着陆冲击过程中涉及材料、几何和接触等非线性问题，采用多体动力学建模方法很难将这些非线性因素考虑在内.例如考虑结构柔性时内外套筒接触面法向力很复杂，而多体动力学模型内外套筒之间往往采用滑动副的约束形式，因此很难准确地计算着陆冲击过程中内外套筒之间的摩擦力.着陆腿多柔性体动力学模型的柔性信息通常是以冲击前着陆腿的构型和边界条件给出的，在着陆冲击过程中，着陆腿支柱内外套筒会产生压缩行程，着陆腿构型和边界条件是不断变化的，用着陆冲击前的构型和边界条件生成的模态中性文件来分析着陆冲击过程也会产生一定的模型误差.此外，多体动力学模型中着陆腿足垫与月壤的接触力只能根据具体工况采用经验公式计算，通用性较差.

本文以非线性有限元理论为基础，充分考虑着陆器着陆腿在着陆冲击过程中的几何、材料和接触等非线性因素，建立了详细的着陆腿有限元模型，并采用ABAQUS/Explicit作为求解器进行着陆冲击非线性动力学分析，并与试验数据进行了对比.仿真结果表明，非线性有限元动力学模型能较准确地预测着陆腿的着陆缓冲性能.

1 着陆腿的构型及缓冲原理

月球着陆器多采用腿式着陆缓冲机构，着陆腿围绕着陆器中心轴均匀分布，如图1所示.单个着陆腿缓冲机构由主支柱、辅助支柱和足垫组成.主支柱外筒上端通过万向节与着陆器主体连接，主支柱内筒下端通过球铰与足垫连接.辅助支柱通过万向节着陆器主体连接，通过球铰与主支柱连接.主、辅支柱内都安装有铝蜂窝芯子缓冲器，利用蜂窝芯子塑性压缩变形吸收着陆冲击过程中的动能，起到缓冲着陆冲击载荷的作用.主支柱内的缓冲器具有单向缓冲功能，辅助支柱缓冲器具有双向缓冲功能.

图1 软着陆机构示意图

2 显式非线性有限元分析理论

显式非线性有限元动力学模型的数值解法，是将结构在空间域进行离散，把连续的微分方程转换成有限阶的代数方程组:

对求解的时间历程在时间域内离散，并采用显式积分方法进行求解.中心差分法是最常用的显式积分方法，积分步骤如下:

其中，u·为速度列阵;u为位移列阵;上标(i)表示增量步的次数，(i-1/2)和(i+1/2)是指增量步的中间值.这样，在增量步开始时提供了满足动力学平衡条件的加速度.得到了加速度后，在时间上“显式地”前推速度和位移.显式积分算法中，不需要形成整体的剪切刚度矩阵，由于是显式的前推模型的状态，所以不需要迭代和收敛准则.显式动力学求解方法是条件稳定的，最大稳定时间积分步长Δts与单元的特征尺寸Le和波速cd相关:

对于线弹性材料，cd的计算公式如下:

其中，E为弹性模量;ν为泊松比;ρ为密度.

3 着陆冲击有限元建模

3.1 着陆腿有限元模型

在ABAQUS/CAE平台上建立着陆腿有限元模型，模型主要包括铝蜂窝缓冲器、主支柱、辅助支柱和足垫，如图2所示.着陆腿支柱内外套筒采用S4R减缩积分壳单元及S3壳单元建模，并对不同区域赋予不同的单元截面属性(section property).足垫采用解析刚体建模.辅助支柱和主支柱之间、主支柱与足垫之间采用球铰连接单元(connector element)连接.

着陆腿中对仿真结果影响不大的细节已经做了简化处理.由于着陆腿一些局部区域的单元尺寸较小，这些单元控制最大稳定时间增量步长，对这些单元较小的区域采用质量缩放技术，将这些局部区域单元的质量增加5倍，人为地增大材料的密度，由式(4)、式(5)可知，稳定时间步长可增加约2倍.这些局部区域质量的增大不会对着陆腿系统的总体动力学响应产生太大影响，但却可以明显地提高计算效率，减少计算量.

图2 着陆腿有限元模型

在着陆冲击过程中，着陆腿各支柱内外套筒接触面之间以及足垫与着陆面之间的相互作用采用接触算法来实现.模型采用ABAQUS/Explicit中的通用接触算法定义着陆腿构件间的接触，同时对不同的接触区域定义不同的接触属性.采用库仑摩擦模型来描述接触面间摩擦的相互作用.摩擦系数定义如下:

其中，μ为摩擦系数;μk为动摩擦系数;μs为静摩擦系数;dc为指数衰减系数;γ·eq为接触面间的相对滑移速度.模型中着陆腿支柱内外套筒之间的摩擦参数取值均由试验标定.

3.2 铝蜂窝缓冲器模型

采用Crushable Foam体积硬化模型作为铝蜂窝缓冲器的材料模型.根据铝蜂窝缓冲器的额定压缩载荷行程特性，定义材料模型的弹性及硬化参数.初始屈服应力由铝蜂窝的压缩载荷和铝蜂窝芯子的横截面积确定.将弱蜂窝极限压缩塑性应变对应的屈服应力设为一个大值，当弱蜂窝达到极限压缩行程后将会硬化，从而促使强蜂窝开始压缩.

为了验证模型的正确性，本文以主支柱二级铝蜂窝缓冲器的载荷行程特性进行了仿真分析，分析模型如图3所示.底部刚性板固支，蜂窝缓冲器固定在底板上，顶部刚性板以4 m/s的速度向下运动，达到蜂窝的极限压缩行程时停止.底部刚性板的支反力(80 Hz滤波)以及铝蜂窝的压缩行程曲线如图4所示.曲线中较平滑的两段对应的压缩载荷与设计值基本一致.整个压缩过程中铝蜂窝的总塑性耗散能为7 548 J，与额定值7 520 J相当.仿真结果证明，该模型能很好地模拟蜂窝缓冲器的缓冲及吸能特性.

图3 蜂窝缓冲器有限元模型

图4 铝蜂窝缓冲器载荷行程曲线

4 单腿着陆冲击仿真及试验验证

单腿着陆冲击试验装置如图5所示.吊篮用于模拟着陆腿的真实安装情况，吊篮安装在导轨上，只能沿导轨在竖直方向运动，吊篮上安装有配重用于模拟不同的着陆工况.试验时，由起吊装置将吊篮及着陆腿吊至指定高度后由释放装置释放，吊篮及着陆腿在重力作用下冲击到着陆平台上.由安装在吊篮上的加速度仪测量吊篮中心处的加速度.本文对表1中的着陆工况进行仿真分析，以验证着陆腿模型的正确性.

图5 单腿着陆冲击试验示意图［7］

表1 着陆冲击试验工况

单腿着陆冲击有限元仿真模型如图6所示.模型由吊篮、着陆腿和着陆面组成.由于吊篮为钢结构，其刚度比着陆腿大，且不是仿真研究所关注的重点，因此采用离散刚体模型建模，其边界条件施加在离散刚体的质心参考点上，约束参考点上的3个旋转自由度和y向、z向的平移自由度.着陆腿支柱和吊篮之间采用万向节连接单元连接.着陆平台分为刚性月面和模拟月壤两种情况.采用Cap Drucke-Prager(CDP)模型作为模拟月壤的材料模型，该模型是在子午线为线性的Drucker-Prager模型上加上一个帽盖屈服面形成的，这个帽盖屈服面设置了静水压力的上限，因此提供了一个非弹性硬化机制用于描述材料的塑性压缩.模型参数由JLU2模拟月壤［9］的三轴压缩试验数据及其相关研究成果标定，具体的模型参数如表2所示.表中，d为材料粘性;β为材料摩擦角;R为帽盖偏心率;p0为帽盖初始屈服面位置;r为过度面半径参数;K为流动应力比.模拟月壤的压缩系数Cc和膨胀系数Cs分别取0.02和0.005，月壤的压缩硬化特性由压缩系数、膨胀系数和初始孔隙率定义.

图6 单腿着陆冲击仿真模型

表2 模拟月壤材料属性

着陆冲击仿真结果与试验结果的对比如表3所示，着陆腿主支柱和辅助支柱的缓冲器的缓冲行程仿真结果与试验结果基本一致，两个工况的误差都在工程上可接受的范围内.这说明非线性有限元模型可以较准确地预测着陆腿的缓冲和吸能特性.工况2的主支柱行程仿真值要比工况1小，这是因为模拟月壤的塑性变形会吸收一部分动能;而工况2主支柱行程试验值却比工况1大，原因主要有以下两方面:蜂窝芯子由于存在制造误差，工况2所用蜂窝芯子实际压缩载荷可能要比工况1所用的蜂窝芯子小;试验中工况1所用铝蜂窝芯子存在压偏的现象，这会增大着陆腿套筒间摩擦力，使得工况1的压缩行程偏小.

为了研究着陆腿的柔性对着陆冲击性能的影响，在工况1条件下［8］，对着陆腿施加刚性约束进行着陆冲击仿真分析，吊篮中心加速度(80 Hz滤波)时间历程如图7所示.吊篮中心加速度试验峰值约为9.18 g，柔性着陆腿模型和刚性着陆腿模型仿真结果分别为8.0 g和7.67 g.由于试验过程中加速度测量点的测试结果受安装刚度影响较大，而仿真模型的吊篮采用刚体建模，导致加速度试验值和仿真值匹配性不是很好，但总体变化趋势基本一致.同时由图7可知，柔性腿模型比刚性腿模型的仿真值更接近试验测试值.

柔性腿和刚性腿主支柱内外套筒之间的摩擦力(80 Hz滤波)如图8所示.柔性腿套筒摩擦力峰值约为5.4 kN，刚性腿峰值约为4.36 kN.刚性腿主支柱套筒间的摩擦力变化比较平缓，柔性腿主支柱摩擦力的波动性较大，这主要是主支柱本身，以及辅助支柱的弹性变形引起的.在着陆冲击前0.15 s过程中，柔性腿摩擦力曲线大部分都位于刚性腿摩擦力曲线上方，可见着陆腿的柔性变形会增大套筒间的摩擦力.

工况2着陆冲击仿真能量时间历程如图9所示.在着陆冲击后0.15 s时刻，系统动能变化量、重力做功及人工应变能总和为4 135 J，塑性变形和摩擦力耗散能总和约为4 077 J，相差的58 J由数值求解过程中引入的粘性耗散掉，由此可见，冲击过程中能量的变化符合能量守恒定律.人工应变能很小说明仿真模型没有发生明显的沙漏现象，仿真结果是可靠的.系统着陆时刻理论初始动能为3367 J，仿真分析初始时刻动能为3384 J，这是由于采用了局部质量放大技术引起的，这个微小的差异对着陆腿的整体冲击动力学响应的影响不大，但使稳定时间步长由原来7.08×10-7s增大到1.46×10-6s，显著地提高了计算效率.

表3 仿真结果与试验结果对比

图7 单腿着陆冲击仿真模型

图8 主支柱内外套筒间的摩擦力时间历程

图9 工况2能量时间历程

5 结论

1)非线性有限元动力学建模方法可以充分考虑着陆腿着陆冲击过程中的材料、几何、接触等非线性因素，能较准确地预测着陆腿的载荷缓冲和吸能特性，可用于指导着陆腿设计及优化.

2)着陆腿柔性对其内外套筒之间的摩擦力有较大的影响，分析时应予以考虑;增加着陆腿刚度可以减小摩擦力.

3)采用局部质量缩放技术可以显著地增大显式动力学程序的稳定时间增量步长，提高计算效率，且不会对总体动态响应产生太大影响.

References)

［1］ William F R.Apollo experience report-lunar module landing gear subsystem［R］.NASA TN D-6850，1972

［2］ Merchant D H，Sawdy D T.Monte Carlo dynamic analysis for lunar module landing loads［J］.Journal of Spacecraft and Rockets，1970，8(1):48-55

［3］ Doiron H H，George A Z.Apollo lunar module landing dynamics［C］//The 41st Structures，Structural Dynamics，and Materials Conferenceand Exhibit.Washington:AIAA/ASME/ASCE/AHS/ASC，2000

［4］ Andrew L B，James R C，Ralph D L，et al.Planetary landers and entry probes［M］.New York:Cambridge University Press，2007

［5］ Masahiro N，Akira M.Modeling for lunar lander by mechanical dynamics software［C］//Modeling and Simulation Technologies Conference and Exhibit.Washington:AIAA，2005

［6］ Hilchenbach M，Kuchemann O，Rosenbauer H.Impact on a comet:Rosetta lander simulations［J］.Planetary and Space Science，2000，48(5):361-369

［7］朱汪，杨建中.月球探测器软着陆机构着陆腿模型与仿真分析［J］.宇航学报，2008，29(6):1723-1728 Zhu Wang，Yang Jianzhong.Modeling and simulation of landing leg for lunar landing gear system ［J］.Journal of Astronautics，2008，29(6):1723-1728(in Chinese)

［8］ 曾福明，杨建中，朱汪，等.月球着陆器着陆缓冲性能研究［J］.航天器工程，2010，19(5):43-49 Zeng Fuming，Yang Jianzhong，Zhu Wang，et al.Research on landing impact attenuation performance of lunar lander［J］.Spacecraft Engineering，2010，19(5):43-49(in Chinese)

［9］杨艳静.模拟月壤和月球车车轮相互作用的数值仿真和试验研究［D］.北京:中国空间技术研究院，2010 Yang Yanjing.Numerical and experimental research of lunar soil stimulant and lunar rover wheel interaction［D］.Beijing:China Academy of Space Technology，2010(in Chinese)