蜂窝夹芯雷达罩抗鸟撞冲击性能仿真分析

李国举，张潜锐，侯晖东，张昕喆,4

（1.郑州航空工业管理学院 航空宇航学院，河南 郑州 450046；2.中国兵器工业标准化研究所，北京 100089；3.中国兵器科学研究院，北京 100089；4.郑州航空工业管理学院 河南省通用航空技术重点实验室，河南 郑州 450046）

0 引 言

现代飞机的机头外部均装备雷达天线罩以保护雷达天线安全。然而，雷达天线罩作为飞机主要的迎风面结构之一，在服役过程中不可避免地会受到鸟体、冰雹等的高速撞击，因此其抗冲击性能对飞机安全有着重要影响。

在所有高速撞击载荷中，鸟撞载荷是最为危险的极限载荷，因此，开展鸟体撞击飞机的研究与测试尤为重要。当前，对于飞机鸟撞研究主要分为鸟撞试验和鸟撞仿真模拟两类方法[1]。孔令勇等[2]采用鸟撞试验装置对碳纤维复合材料层合板进行撞击测试，重点对材料的铺层方式与其耗能之间的关系进行了系统研究。陈琨等[3]采用空气炮发射明胶鸟弹撞击复合材料蜂窝夹芯板，借助高速摄影研究了夹芯结构在软体高速冲击下的损伤情况。但由于高速摄影、空气炮等试验设备对场所有一定要求且价格高昂，导致研究成本过高，且操作过程复杂、耗时。随着计算机性能和计算力学的不断发展，有限元仿真方法逐渐成为分析鸟撞问题的主要手段[4-6]。郑涵天等[7]采用有限单元法（finite element method，简称FEM），研究了雷达罩抗鸟撞动力响应及结构受损情况。刘玄等[8]采用FEM 方法对无人机雷达天线罩在抗鸟撞冲击过程的耗能进行了研究。然而，FEM 方法在结构受到大变形时容易引起单元畸变，难以有效模拟高速鸟撞过程中的鸟体塑性流动和铺展过程，而光滑粒子流体动力学(smooth particle hydrodynamics，简称SPH)方法可以有效克服FEM 方法因结构大变形引起的计算不稳定问题，能够较好地模拟鸟体在高速冲击过程中的飞溅效果，因此已被广泛应用于鸟撞仿真分析中[9-12]。毋玲等[13]采用SPH粒子算法，对雷达罩抗高速鸟撞过程中的耗能随撞击位置的变化情况进行了系统研究。于永强等[14]同样采用SPH 粒子研究了复合材料层合板受鸟撞击的损失失效情况。研究发现，鸟体与天线罩之间的接触传力对鸟撞过程分析非常关键，而当前单一SPH方法对于接触模拟尚不够精确，限制了其在模拟鸟撞过程中的进一步应用。

为弥补上述研究中单一FEM 与SPH 方法的不足并充分发挥两种方法各自的优势，本文引入FEMSPH 耦合算法[15-18]，开展对无人机天线雷达罩不同位置的抗鸟撞冲击动态响应过程仿真模拟研究。该算法在鸟体材料失效前以拉格朗日有限元网格单元模拟鸟体与无人机雷达天线罩之间的接触传力过程；鸟体失效后有限元网格转换为SPH 粒子，以模拟其后续的飞溅和铺展等大变形过程。由此得到的仿真结果相较于单一FEM 与SPH 方法更加真实可信，可为鸟撞仿真模拟以及飞行器结构设计领域提供理论依据与参考数据。

1 雷达天线罩鸟撞有限元模型建立

1.1 简化鸟体模型选择及有限元模型建立

能够描述真实鸟体结构特征是开展鸟撞仿真分析的首要前提。在实际有限元仿真计算中，通常将真实鸟体形状简化为圆柱与其端部的两个半球形成的半球端圆柱体结构。其中，圆柱体的高度是半球体半径的两倍，半球端圆柱体的特征几何尺寸和鸟体的质量、密度的关系由式(1)给出。

式(1)中，L为圆柱体的高度，R为半球体半径，m为鸟体的质量，ρ为鸟体的密度。

国际上对鸟撞研究有着较多的规范，一般选取四类鸟体进行表征，这四类鸟体的质量和密度参数如表1 所示。因此，本文同样选择家鸡为研究对象，将其质量m=1.8Kg，密度ρ=977.71Kg/m3带入式(1)，计算可得家鸡简化模型中圆柱高度L=112mm，半球半径R=56mm，如图1(a)所示。最终通过CATIA 软件建模获得鸟体几何模型如图1( b)所示。采用LSPrepost 软件建立的鸟体有限元网格模型如图1(c)所示，并进一步基于有限单元中心点获得对应的SPH粒子模型如图1(d)所示。鸟体的FEM 单元为六面体实体单元，其网格尺寸控制在3mm—7mm。最终获得鸟体有限元网格单元数量为56 832个，节点数量为58 407个。

图1 鸟体形状尺寸及其计算模型

表1 四类鸟体质量和密度参数

本文将分别采用FEM、SPH 和FEM-SPH 耦合方法描述鸟体材料在冲击过程中的变形失效情况，并对比分析拉格朗日有限元网格模拟碰撞接触的精确性以及SPH 粒子允许大变形特性在鸟撞问题中的影响。

拉格朗日有限元网格单元采用式（2）给出的弹塑性硬化模型描述鸟体的动态力学行为。

式（2）中：σy为塑性应力；σs为屈服强度应力；EP和εP分别为切线模量和塑性应变。鸟体的失效准则采用最大塑性应变进行描述。

其中，鸟体的弹塑性硬化模型参数取值如下：E=1×104MPa，v=0.3，σs=1.0MPa，Ep=5.0MPa，εcr=1.22[19]。

考虑到拉格朗日有限元网格单元变形过大时会因产生负体积而导致计算中止，因此本文设置当塑性应变超过1.22 时单元被删除，从而保证鸟撞仿真过程中计算的稳定性。对于采用SPH 粒子算法的鸟体模型，则直接采用空材料（*MAT_NULL）描述其力学行为。而对于采用FEM-SPH 耦合算法的鸟体模型，单元未失效前采用拉格朗日有限元网格，此时用弹塑性硬化模型描述鸟体动态力学行为，而当鸟体有限元单元塑性应变达到1.22 时自动转换为SPH 粒子，并将失效单元节点位移、质量、速度、应力和声速等物理量赋给转换后的SPH 粒子，相应转换关系见式(3)—(7)。

式(3)—(7)中，角标p 和n 分别代表SPH 粒子和有FEM 的节点，而角标e 和g 则分别单元和高斯点；x和m为鸟体模型的位移和质量；Ne为与节点相关的FEM 单元总数，Nn为任意单元的节点数量；ρei、Vei、σei分别为单元材料的密度、体积和应力张量；wgj和σgj为单元序号i内部的高斯积分点g 处对应的加权系数和应力张量；Ng为任一FEM 单元所包含的高斯积分点数量。SPH 粒子采用空材料（*MAT_NULL）描述其随后的飞溅和铺展流动效果，鸟体失效后的偏应力由式(8)给出。

式(8)中，σ为鸟体失效后的偏应力;υd为是动态粘度;D为变形率张量;P为压力;II为张量不变量。

空材料模型一般配以多项式状态方程进行描述。鸟体体积变形与压力之间关系由式(9)给出，材料的压缩比由式(10)描述，对应的鸟体状态方程参数如表3所示。

式(9)—(10)中，U为鸟体单位体积的内能;μ为材料的压缩比;ρ和ρ0分别为材料的当前密度和初始密度;C0—C6均为材料参数。

其中，鸟体的状态方程参数取值如下：C0=0，C1=2.25×103MPa，C2—C6均为0。

1.2 无人机雷达天线罩有限元模型建立

本文的无人机雷达天线罩采用三明治夹芯结构，其结构示意图见图2。天线罩内外两层为0.8mm的2024T3Al 合金板，中间夹芯层为8.4mm 的Nomex蜂窝铝材料[19]。

图2 无人机雷达天线罩的三明治夹芯结构示意图

通过CATIA 软件建立三维实体模型并采用Hypermesh 对其进行网格剖分，获得无人机雷达天线罩有限元网格模型如图3 所示。为了体现天线罩的三维结构情况，图中分别给出了如图3(a)所示的正视图和如图3(b)所示的侧视图。其中夹芯层采用六面体实体网格单元(Solid 164)进行剖分，截面方向的网格尺寸设置为5mm,厚度方向尺寸为2.8mm，以确保厚度方向至少划分3 层实体单元。最终获得夹芯层结构的六面体单元数量为115 680 个，节点数量为155 848 个。对内、外两层铝合金面板采用壳单元(shell 163)进行面网格划分，单元尺寸同样设为5 mm，最终获得内、外层单元数量为38 880，节点数为39 284,并将内、外层壳单元的厚度均设置为0.8 mm。

图3 无人机雷达天线罩

雷达罩内、外层材料的动态力学行为由式(11)所给出的基于连续损伤理论的104 号材料模型(*MAT_DAMAGE_1)进行描述。

式(11)中，r为损伤累积的塑性应变判据为等效塑性应变率;σ0为材料初始屈服强度;Q1、Q2、C1、C2为材料各向同性硬化系数。

其中，天线罩外壳的材料模型参数取值如下：σ0=364.5MPa，Q1=334.7MPa，C1=6.16，Q2=0，C2=0，r=0.18。

天线罩夹芯层泡沫铝结构材料采用154 号材料（* MAT_DESH PANDE_FLECK_FOAM）进行描述。其中，材料屈服强度σy和材料参数的关系由式(12)给出，压实应变由式(13)定义给出。

式(12)—(13)中，α2、γ和β都为材料参数;σp为平台应力;εD为压实应变;ε̂为等效应变;ρf为泡沫材料当前密度;ρf0为泡沫材料初始密度。

其中，天线罩夹芯层泡沫铝结构的材料参数取值如下：ρf0=300Kg/m3，E=1.5×103MPa，α=2.1，γ=6.10MPa，εD=2.2，α2=38.1MPa，β=3.1，σp=4.41MPa，εcr=0.2。[20]

本文中，鸟体相对天线罩结构的速度设为150m/s，鸟体与天线罩正面相撞，撞击点沿高度方向（Z 轴正方向）排列，设定为5 组，每组撞击位置在Z 方向的位移间隔为100mm，如图4所示。天线罩模型边界采用固支约束，鸟体拉格朗日有限元网格与天线罩结构使用自动面面接触算法进行接触控制，鸟体网格失效后转换为SPH 粒子与天线罩结构采用自动点面接触算法进行接触控制。模型整体求解时间设置为6ms，每0.1ms 输出一个结果，时间步长系数设置为0.6。

图4 无人机雷达天线罩的鸟撞不同位置有限元模型示意图

2 天线罩鸟撞过程有限元仿真分析

2.1 天线罩正面中间鸟撞位置处鸟体不同算法仿真结果对比

天线罩正面中间位置处不同鸟体描述算法对应的天线罩撞击点整体位移—时间曲线及其在T=3ms时刻的鸟体变形铺展状态如图5所示。由图中可知，在时间为0ms—1.5ms 时，三种算法描述的鸟体撞击后获得的天线罩撞击点整体位移—时间曲线基本一致。在同一时刻对应的撞击点整体位移中，FEM 算法获得的位移最大，表明鸟体初始接触天线罩体时FEM 算法接触最为充分。随着时间的增加，同一时刻SPH算法获得的位移最大，FEM最小。观察T=3ms时鸟体变形铺展过程可知，随着FEM 算法中单元塑性变形增大并达到失效应变而被不断删除，最终由于质量损失过大而无法对天线罩继续作用。而SPH和FEM-SPH 算法中，不存在单元删除，因此两者获得位移均高于FEM 算法。对比SPH 算法和FEMSPH 算法在T=3ms时的变形铺展过程可知，SPH 算法鸟体已经完全铺展开并向四周飞散。而在FEM-SPH耦合算法中，由于未失效拉格朗日有限单元之间的相互连接，使得鸟体模型完全铺展后并未如SPH 算法一样快速飞散，这表明FEM-SPH 耦合算法更能精细描述鸟体撞击到天线罩上的变形、失效以及铺展过程。

图5 天线罩正面中间位置处不同鸟体算法对应的撞击点整体位移—时间曲线以及T=3 ms对应的鸟体变形铺展过程

2.2 天线罩正面中间鸟撞位置处鸟体铺展和天线罩变形过程

天线罩正面中间位置处随鸟撞不同时刻对应的鸟体铺展和天线罩整体位移云图如图6 所示。从图中可知，在T=2ms 时，鸟体与天线罩开始接触并产生局部变形，鸟体头部有限元网格单元已经转化为SPH粒子并开始在天线罩四周铺展，此时，天线罩结构整体的位移峰值为83mm，见图6(a)。在T=4ms时，天线罩位移已经达到最大值96mm，且鸟体大部分FEM 单元已经失效转换为SPH 粒子并向天线罩上部铺展，此时的鸟体铺展与飞溅效果已经比较明显，如图6(b)所示。在T=6ms时，鸟体网格单元基本上全部转换为SPH 粒子，鸟体飞溅效果更加显著，并开始整体滑离出天线罩表面，如图6(c)所示。

图6 天线罩正面中间位置鸟撞不同时刻鸟体铺展和天线罩位移变形云图

2.3 天线罩正面中间鸟撞位置处外层铝合金板的应力演变过程

天线罩正面中间位置处鸟撞不同时刻对应的天线罩外层铝合金板Von-Mises等效应力云图如图7所示。从图中可知，在T=0.1ms 时，鸟体与天线罩开始接触，此时撞击点位置等效应力最大值达到了383.0MPa，如图7(a)所示。随着时间的增加，当T=0.8ms时，鸟撞接触位置处产生应力集中，最大等效应力增加到427.3MPa，此时由于鸟体头部单元尚未失效，整个接触区域呈现出应力值集中现象，如图7(b)所示。当T=1.3ms 时，应力峰值增加到431.8MPa，此时与鸟体头部上侧位置转换为SPH 粒子，导致此处应力值降低，仅为316MPa，远小于其他接触区域应力值，如图7(c)所示。在T=3.3ms时，鸟体头部接触区域应力峰值达到整个鸟撞过程最大值452.0MPa，此时鸟体头部一部分单元已经转换为SPH 粒子并呈现出流体的铺展和飞溅等流动性特征，使得这些接触区域应力值稳定在300MPa 左右，如图7(d)所示。在T=5.4ms 时，鸟体几乎完全转化为SPH 粒子，即主要呈现其流动性特征，使得应力集中转移到鸟撞形成的坑体周围区域，最大值也降至408MPa，如图7(e)所示。在T=6.0ms时，由于大部分SPH 粒子已经从鸟撞过程中形成的坑体内滑移出去，因此整体接触区域应力值较小，仅为200MPa左右，如图7(f)所示。

图7 天线罩正面中间位置鸟撞不同时刻外层铝合金板Von-Mises等效应力云图

2.4 天线罩正面中间鸟撞位置处三明治结构的耗能情况

无人机雷达天线罩内外层以及夹芯层吸收的内能随时间变化对比曲线如图8所示。

图8 工况1条件下天线罩三明治结构内能—时间曲线对比

从曲线对比可知，在T=4.0ms 时，天线罩内外层及夹芯层内能均达到最大值，后续随时间变化保持不变。此时，泡沫铝夹芯层内能达到3.6×103J，而内、外铝合金面板内能仅为775J 和682J，这表明鸟撞过程中冲击主要由泡沫铝夹芯层吸收，同时也表明泡沫铝夹芯层材料在天线罩结构中起着吸收鸟撞冲击能力的重要作用。因此为了有效保护天线罩内部的雷达和精密设备，需要对夹芯层材料进行选择和设计。

2.5 不同鸟撞位置下天线罩位移形变情况

不同鸟撞位置下无人机雷达天线罩对应的最大整体位移云图如图9 所示。从图中数据可知，天线罩鸟撞位置从上部向下移动过程中，天线罩结构最大整体位移呈现先增加后减小的趋势。无人机天线罩最危险鸟撞位置为中间靠下的位置，此时对应的最大整体位移峰值为123.5mm，较之鸟撞最上部位置的67.0mm，位移增加了约84%。鸟撞下部位置最大位移为103.4mm，稍高于中间位置的99.7mm。这表明鸟撞位置变化对天线罩整体变形影响很大，因此需要重点加强危险位置的结构设置。

图9 不同鸟撞位置下对应的天线罩模型最大整体位移云图

3 结 论

（1）FEM-SPH 耦合方法能够充分发挥FEM 法和SPH法各自的优势，能更精细描述鸟体撞击到天线罩上的变形、失效以及铺展飞溅过程。在鸟体撞击天线罩正面中间位置的过程中，鸟体有限元网格单元不断失效并转变为SPH 粒子，在T=4ms 时出现明显的飞溅和铺展现象。进一步研究其外层铝合金板应力演化发现，在T=5.4ms 后由于SPH 飞溅出撞击凹陷，应力主要集中在凹陷周围，凹陷底部应力显著减小。

（2）在鸟体撞击天线罩正面中间位置的过程中，泡沫铝夹芯层材料耗能达到3666J，远大于内外层铝合金面板耗能，表明泡沫铝夹芯层材料在天线罩结构中起着吸收鸟撞冲击力的重要作用。

（3）鸟体以相同速度撞击天线罩不同位置处的有限元仿真结果表明，随着天线罩鸟撞位置从顶部向下移动，天线罩结构最大整体位移呈现先增加后减小的趋势。无人机天线罩最危险的鸟撞位置为中间靠下位置，此时最大位移为123.5 mm，为上部位置最大位移值的1.8倍。