介绍
在航空航天工业中，鸟击通常是指空中动物（鸟或蝙蝠）与飞机之间的碰撞。每年都会报告数以万计的鸟击事故，但其中只有不到 5% 的事故可被视为 破坏性撞击。相关货币成本估计超过 12 亿美元。
鸟击物理测试项目包括将一只真正的鸟射入正在运行的喷气发动机并确定发动机的损坏程度。提高测试程序效率的有效方法是进行数值模拟。高保真数值模型意味着可以在物理测试的实际最后阶段之前优化发动机设计，从而节省大量资源。
在本文中，受到感恩节和我对喷气发动机的热爱的启发，我使用全尺寸高旁通喷气发动机和实际的火鸡模型进行了鸟击模拟。鸟撞击的动态性质以及鸟在撞击时的高变形能力使得经典有限元方法受到限制，因此选择平滑粒子流体动力学方法。
背景：平滑粒子流体动力学
鸟击测试是一个高度动态的问题。研究表明（Heimbs）经典拉格朗日 框架无法产生可靠的结果。鸟在撞击过程中的高度变形表明，可以更有效地使用其他建模技术来模拟鸟与喷气发动机风扇叶片的撞击。平滑粒子流体动力学可以轻松满足这些建模需求。
拉格朗日 方法
欧拉方法
任意拉格朗日-欧拉 (ALE) 方法
平滑粒子流体动力学 (SPH) 方法
图1。不同的建模方法
在深入探讨如何实施 SPH 技术并进行鸟击模拟之前，我想快速概述和比较可用于解决此类问题的不同有限元方法。
在图 1中，我试图呈现一个直观的示意图，说明上述技术之间的差异 。在拉格朗日 描述（通常是结构分析的首选方法）中，模型的节点与质点精确移动。在此框架内，可以更直观地跟踪自由表面并应用不同类型的载荷和边界条件。然而，这种方法的缺点是，面对高应变梯度时，网格更容易变形。
与拉格朗日描述相反，在欧拉框架中，空间网格的节点在材料流过时保持固定。在欧拉描述中跟踪自由表面很困难，但由于网格在空间中是固定的，所以不存在网格变形的问题。欧拉方法可在 Abaqus/CAE 中使用耦合欧拉-拉格朗日 (CEL) 功能来实现。
在任意拉格朗日-欧拉（ALE）方法中，网格运动仅在必要时才被约束到材料上，主要是在自由边界处，否则，材料运动和网格是独立的。这种方法可以最大限度地减少网格变形，同时保持相同的网格拓扑。
平滑粒子流体动力学 (SPH) 方法可解决其他经典方法失败或效率低下的建模需求。传统 CFD 无法应对的剧烈流体流动（例如浅水流）或 CEL 难以/不可能实施的极高变形（例如弹道和积雪压实）都是 SPH 可以高效使用的例子。
图2 . 内核函数
平滑粒子流体动力学方法的新颖之处在于在移动宏观粒子的不规则网格内进行平滑插值和微分的特定方法。由于节点连通性不固定，因此避免了严重的单元变形，因此该方法允许非常高的应变梯度。
SPH 方法只是连续偏微分方程的离散化方法，而不是模拟不同粒子之间相互作用（接触）的方法（如离散元法或 DEM）。该方法使用演化插值方案来近似任意空间点的场变量。
详细介绍数学公式和方法实现的基础知识超出了本文的范围。我可以向读者推荐 JJ Monaghan 的“平滑粒子流体动力学代码基础知识”以及 Abaqus 文档。
Abaqus 中的平滑粒子流体动力学
SPH 方法是 仅适用于Abaqus/Explicit 的一项功能，只能针对 3-D 模型实现。SPH 粒子是使用单节点 PC3D 元素定义的，这些粒子元素使用 Abaqus 的现有功能来引用与单元相关的特征，例如材料定义、载荷和边界条件。
激活 SPH 粒子可以使用两种不同的转换技术来实现；可以为每个父元素生成粒子，也可以基于统一的网格生成粒子。
当预期有高应变梯度时，通常使用从有限元网格到每个父单元的粒子的转换技术，使得严重扭曲单元的经典单元删除技术受到限制。每个父元素生成的粒子数量以及转换标准都可以控制。根据模拟的性质，可以指定四种不同的转换标准；基于时间、基于应变、基于应力和基于用户子例程。此功能可通过 Abaqus/CAE 的“单元类型”控制菜单下的网格模块来使用。
图3 . 每个父元素内部生成的粒子，针对每个等参方向三个粒子进行说明
基于均匀背景网格方法生成粒子的目的是在分析开始时生成均匀分布的粒子。该方法对于需要粒子均匀分布的应用非常有用，例如在模拟流体时。背景网格的原点放置在要转换的元素集的质量中心在分析开始时，集 合中的有限元将转换为 SPH 粒子。粒子仅在父元素体积内的网格线交叉处生成。所有生成的粒子都具有相同的体积和质量。任何不包含网格交集的父元素都不会生成与其关联的粒子。
图4 . 基于均匀背景的内部生成粒子。
在本文中，我将使用第一种方法将有限元网格转换为粒子。我将使用基于时间的转换标准。
让火鸡做好准备
在深入研究模型设置之前，我想先谈谈鸟击模拟的撞击物几何形状。知道没有标准化的人造鸟形状（Heimbs），可以在文献中找到各种几何形状。通常选择不同的形状来反映鸟类躯干的质量和几何形状。文献中发现的主要几何形状是圆柱体、两端为半球形的圆柱体、椭球体、球体。
图5 . 鸟击模拟过程中鸟类的不同替代几何形状
在本文中，我在模拟过程中使用了鸟的真实几何形状，在这个场景中是一只火鸡。土耳其模型是通过真正的火鸡扫描得到的，然后得到.stl文件。CATIA 用于生成土耳其模型的 CAD 文件。然后将火鸡模型导入 Abaqus/CAE。
图6 . 火鸡 3D 模型
我使用 C3D4 元素（4 节点线性四面体）对土耳其模型进行网格划分。选择网格元素类型时，我从“转换为粒子”子菜单中切换为“是” ，如图8所示。我选择了基于时间的转换方法，并选择 0 作为Threshold，这意味着土耳其模型中的所有元素都将在第一步增量开始时转换为粒子。
图7 . 火鸡网格模型
图8 . 单元类型选择和转换为粒子定义
我还选择了PPD 为 1，这是每个等参方向生成的粒子数。SPH模型如图9所示。
图9 . 将模型转换为 SPH 粒子
一般的假设是鸟主要由水组成。弹塑性材料模型甚至超弹性材料模型的总体框架似乎是有限的。因此，通常使用状态方程（EOS）来构建鸟类的本构模型。EOS 表示鸟类撞击器的压力和体积与室温下水的参数之间的关系。常见的状态方程包括线性方程、多项式方程、表格方程和 Mie-Grüneisen 方程。我在 Abaqus/CAE 中使用了激光状态方程定义，它还指定了冲击速度与粒子速度之间的线性差来定义 Hugoniot 压力。火鸡材料参数包括质量密度以及 EOS 参数，其中包括c0和 Γ0。
喷气发动机模型准备
为鸟击模拟准备了虚拟涡轮风扇喷气发动机模型。图 10显示了发动机的等视图和侧视图。本文的主要目的是展示 SPH 方法的正确实现技术，而不是提供任何工程见解。然而，根据我自己的文献研究，材料参数和模拟边界条件被选择为最接近的现实场景。
新型发动机风扇叶片一般采用钛合金Ti-6Al-4V，这是一种α-β钛合金，具有高强度重量比和优异的耐腐蚀性能。冲击的高动态效应要求考虑应变率对风扇材料的影响。风机所用材料参数汇总如表1和表2所示。
图10 . 虚拟喷气发动机的等视图和侧视图
为了使模拟更轻松，我使用了 Display Body 约束，并将模拟域仅限制为风扇（图 10中以绿色突出显示）。
表1：钛合金的材料特性（Lesuer，2000）。
表2：钛合金的 Johnson-Cook 损伤参数（Lesuer，2000）。
鸟击模型设置
如上所述，在模拟过程中，只有风扇被视为连续体，喷气发动机的其余部分是使用无运动的显示体约束进行建模的。添加耦合约束以将风扇旋转限制到位于主轴中心的参考点的旋转。然后使用预定义场将旋转速度边界条件应用于该参考点。选择的角速度为 41.88 弧度/秒。
图11 . 火鸡与喷气式飞机发动机组件的 ISO 视图
图12 . 带有风扇叶片的火鸡的放大视图
火鸡初始撞击速度约为 400 英里/小时。我使用无摩擦公式的一般接触方法来模拟撞击时火鸡和风扇之间的摩擦。使用 Abaqus/Explicit 运行仿真 0.01 秒或 10 毫秒。火鸡和风扇之间的撞击角度为 90 度。
结果
图13 . 从左上到右下的冲击进展（0 秒、2.5 毫秒、5 毫秒和 7.5 毫秒）
如上所述，土耳其有限元模型在第一次增量开始时转换为 SPH 粒子。火鸡和风扇叶片之间的撞击发生在风扇旋转时。
图14 . Von Mises 和风扇上的压应力分布
大部分冲击发生在风扇叶片上，产生约 155 MPa 的局部应力，以及 31 MPa 的局部压力。然而，这种特定的模型设置并未在任何受影响的叶片处产生任何局部塑性变形。
图15 . 压力应力随时间的变化
图 15显示了冲击期间和冲击之后压力应力的演变。这表明压力在撞击时立即增加，然后随着时间的推移而减少。该模型的动能显示出类似的趋势。
图16 . 模型的动能
由于模型主要受鸟的速度控制，因此模型的动能随着时间显着减小。我们可以看到，撞击时的动能显着减少，因为它以应变能的形式传递到风扇叶片。
图17 . 模型的应变能
总结
本文进行了鸟击分析。
选择平滑粒子流体动力学方法来执行这种模拟
使用每父元素内部生成粒子方法将有限元网格转换为 SPH 粒子
使用基于时间的转换标准将 FEA 网格转换为SPH 粒子
