摘要：和成熟的FE（有限元）方法相比，SPH方法可以较为真实地描述超高速碰撞过程的物理现象。本文利用SPH方法对Whipple防护结构在空间碎片超高速碰撞下的物理过程进行了数值模拟。在数值模拟中，为了充分发挥SPH方法和FE方法的优点，利用有限元单元和SPH节点混合建模，计算结果和试验结果相吻合，且计算效率高，可用于空间飞行器防护结构的设计。

1 SPH方法

SPH（Smoothed Particle Hydrodynamics）方法的全称为“光滑质点流体动力学”方法，是Lucy于1977年提出的用于天体物理计算的一种无网格化Lagrange计算方法[1,2]，其特点是在模拟物体大变形时，既可以克服Euler方法难于跟踪物质变形和不能识别材料界面位形的缺点，同时也解决了传统Lagrange方法在大变形下的网格扭曲（Distortion）问题，因此在工程上有很大的应用潜力。目前，SPH方法的应用领域主要有：天体物理，燃气动力学，结构动力学(包括超高速碰撞)等等[1,2]。

SPH方法的核心是一种插值技术。在SPH方法中，任一宏观变量（如密度、压力、温度、内能等）A（r）能方便地借助于一组无序点（disordered points）上的值表示成积分插值计算得到，其形式为：

（1）

其中，Ω为整个求解区域，W为核函数，h为“光滑长度”（smoothing length）,用于定义一个质点对周围区域的影响程度。

SPH方法的主要特点有：利用最接近相邻搜索算法(Nearest Neighbor Search)，可选择核函数W，提供可变的“光滑长度”选择以及可以和有限元方法混合建模进行求解（如图1所示）等等。

图1 通过定义接触条件，利用SPH节点和有限元单元实现混合建模

。由于算法差异，SPH方法相对有限元方法其求解时间更长，因此对于一些特殊问题的数值模拟，建模时可以将有限元单元和SPH节点（SPH Nodes）通过定义接触条件相结合，在大变形和飞溅区域采用SPH节点建模，而小变形区域则采用有限元单元建模，从而大大节省求解时间，提高计算效率。有关SPH方法的详细理论介绍，可参考文献[1,2]。

2 超高速碰撞的数值模拟

2.1 防护结构

随着人类空间发射次数的日益增加，超高速飞行的空间碎片（相对平均速度在10km/s左右[3]）对航天器的碰撞问题已成为航天器设计中必须考虑的因素。典型的Whipple防护结构是天体物理学家Whipple在20世纪40年代提出的一种双层板防护系统[3]。该结构由前板、后板和两板间一定的间隙组成。由于前板能起到将高速碎片充分粉碎甚至尽量融化或汽化的作用，从而大大消耗了碎片的动能和破坏了碎片的完整性，因而减轻了对于后板的破坏，起到防护作用。

为了和文献[3]的试验数据相对比，此次数值模拟的Whipple防护结构如图2所示，其中，空间碎片用圆柱形弹丸来模拟，试验中，前后板在距离中心RB处进行固定，其他结构特点及具体数值参见表1。

图2 Whipple结构尺寸示意 图3 碎片及防护板的离散模型

表1 Whipple防护结构参数

对上述的Whipple防护结构进行离散化建模，考虑到弹丸以及前后板中间部位的变形较大，且空间碎片碰撞后飞溅的特点，在建模时采用SPH节点，而前后板的其他部位则采用有限元单元（六面体8节点）进行建模，整个系统离散化后的计算模型如图3所示，放大后如图4所示。采用混合建模，既可以充分利用SPH方法的优点，同时又可以大大节省计算时间。

图4 碎片及防护板的离散模型放大图(a 防护板 b弹丸/碎片)

2.2 材料模型

超高速碰撞过程是一个极其复杂的过程，对其完善的描述应包括连续介质力学的各个方面，而且还要考虑到溅射、相变等问题，因此材料参数和模型的选取是超高速碰撞下非常重要的环节。在材料参数和模型选取时，主要考虑以下三方面的内容：①材料强度本构关系模型，即大应变率下材料的应力应变关系描述；②物态方程，即高压下材料的流体压力和密度、内能等的关系描述；③材料的失效准则，对于有限元计算主要是考虑单元的失效标准。由于SPH方法和有限元求解原理的差异，利用SPH方法求解时不需要定义材料的失效准则，材料是否失效，完全根据SPH节点的求解结果来确定。因此，和有限元方法相比，SPH方法减少了对试验数据的依赖性。

弹丸和前后板的材料参数可参见文献[3]。计算时，弹丸的材料为聚乙烯（Polyethylene），其材料参数均和文献[3]相同；前后板的材料均为Al-2024，其材料模型选取见表2，其中，强度模型选用Johnson-Cook模型，物态方程选用多项式形式来描述。

表2 前后板材料Al_2024在计算中使用的材料模型

在SPH方法求解中，选用的核函数W为三次样条函数，节点光滑长度h和半径的比值取为2.1，通过定义接触对来实现有限元单元和SPH节点的结合。

3 计算结果

选用PAM-SHOCK(ESI Group)提供的SPH求解器来进行此次超高速碰撞的求解。求解后，空间碎片的整个碰撞过程如图5(a)所示，从图中可以看到，弹丸在穿透前板后，形成二次碎片，碎片云经膨胀和拉长，对后板造成轻微的损伤。计算结果形象地描述了碎片云的形成和飞溅过程，符合有关超高速碰撞试验的图象记录[3,5]。

图5 空间碎片碰撞过程的各个阶段（时间单位 μs a SPH方法 b 有限元方法）

针对同一个算例，将SPH方法的计算结果和有限元方法相比较[4]，图5(b)是利用有限元方法建模和求解相同算例的碰撞过程。对比两种方法的计算结果可以看到：在描述主碎片云（碎片云的主体部分）的形成和运动轨迹上，两种方法的计算结果是基本一致的；SPH方法的优点是更为形象地描述了二次碎片云的膨胀和拉长过程（图6，注：为了清楚易见，且考虑到有限元单元建模部分位移很小，图中仅显示了SPH节点建模部分）。

图6 t=95us，碎片云撞击后板（仅显示SPH节点）

将SPH方法计算得到的前后板破坏情况和试验结果进行对比：后板有轻微破坏，与试验结果一致（图7）；前板破坏孔的直径为42.4mm，和试验结果（约42mm）几乎一致（图8）。此外，从SPH的计算结果的云图中，还可以看到各节点碎片云的密度分布情况，从而判断碎片云所处的状态（固、液、气）（参见图6至图8）。

图7 后板破坏情况（仅显示SPH节点）

图8 前板破坏情况（仅显示SPH节点 a 后视图 b 前视图）

4 结论

综上所述，利用SPH方法对空间碎片的超高速碰撞过程进行数值模拟，不仅很好地预测了Whipple防护结构的破坏情况，而且对整个碰撞过程，包括碎片云的形成、膨胀和拉长过程都有形象的描述，符合超高速碰撞的试验测试结果；利用SPH方法无需确定材料的失效准则，减少了对试验数据的依赖性。此外，鉴于SPH方法在求解时计算时间较长，在数学建模时，可以利用SPH节点和有限元单元混合建模，便于在工程设计中使用和推广。

非常感谢ESI Group中国公司对本文工作的支持！

参考文献
1. E. Haug, P. Groenenboom, A. Kamoulakos, etc. , Application of SPH Techniques in the PAM-SCL Code Family, PUCA'97, pp273-282, ESI Group.
2. Moreno Faraud, Roberto Destefanis, David Palmieri and Mario Marchetti, SPH Simulation of Debris Impacts Using Two Different Computer Codes, Int. J. of Impact Engineering, Vol.23, pp.249-260,1999.
3. M. Katayama, S. Toda and Kibe, Numerical Simulation of Space Debris Impacts on the Whipple Shield, Acta Astronautica, Vol. 40, No12, pp.859-869, 1997.
4. 闫晓军，张玉珠，聂景旭，超高速碰撞下Whipple防护结构的数值模拟，宇航学报，Vol23，No.5, pp81－84。
5. 闫晓军，超高速碰撞、复合冲蚀、盘/叶耦合振动的数值模拟与分析，北京航空航天大学博士后研究报告，2002。

闫晓军 聂景旭：北京航空航天大学能源与动力工程学院
张玉珠：北京航空航天大学宇航学院