多层防护结构广泛应用于军事与民用领域，如军舰、航天器等。其中，航天器大多采用多层间隔式防护结构^[1]。已有研究^[2-3]表明，以前用于拦截导弹等航天器的防空反导战斗部主要采用破片杀伤式，目前多采用动能块、动能拦截器等直接碰撞杀伤式。与破片类似，动能块也具有质量与形状可控的优势，并且加工相对简单方便，但动能块的质量可达两百多克甚至更大，而破片的质量一般不会很大。此外，动能块的质量和速度等存在差异，对目标的毁伤效果不尽相同，因此研究不同动能块质量和速度对目标毁伤的影响具有十分重要的意义。

目前，在动能毁伤元对多层防护结构的超高速撞击特性研究中，国内外研究主要集中在不同材料、破片形状以及层间距对多层防护结构抗侵彻性能的影响。Gupta等^[4]、董永香等^[5]对不同金属材料的多层间隔靶的侵彻响应进行了分析，得出了弹丸与间隔靶作用过程的物理图像和演变规律。吕珮毅等^[6]通过数值模拟研究了球形、圆柱形、立方体破片高速侵彻4层铝合金板，指出立方体破片的速度和动能的衰减最大，侵彻能力最小。关于层间距对多层防护结构抗侵彻性能的影响，不同学者针对不同形状弹体侵彻间隔式和叠合式多层防护结构进行了研究，得出了不同的结论。屈科佛等^[7]对球形、圆柱形、立方体破片分别高速正侵彻多层间隔式和叠合式靶板进行了数值模拟，指出间隔式多层靶板的抗侵彻能力强于叠合式多层靶板。Deng等^[8]对尖卵形弹体侵彻多层钢板进行了实验研究，指出叠合式多层结构的抗侵彻能力优于间隔式。针对动能块质量对多层防护结构抗侵彻性能的影响，大多数研究集中在小质量破片侵彻多层结构上。赵小峰^[9]对2.9、3.4和8.0 g 3种不同质量的立方体钨合金破片侵彻靶板进行了数值模拟、理论计算和实验研究，结果表明破片的侵彻能力随着质量的增加而增强，实验结果均略低于理论计算结果，且质量越低，偏差越大。在动能块速度对多层防护结构抗侵彻性能影响的研究中，侵彻速度大多集中在3 km/s以下，很少达到超高速。

相对于实验方法，数值模拟方法更加方便快捷，外界限制因素少，成本低，可重复性高，不仅能够弥补实验研究的不足，处理实验无法解决的工况，还能为实验提供指导。采用数值模拟方法研究多层防护结构的侵彻问题时，要求能够有效地再现碰撞过程中产生的碎片。有限元方法（finite element method，FEM）具有计算效率高、适用范围广的特点，但在模拟侵彻问题时几何大变形容易导致网格畸变，从而终止计算。为了避免网格畸变，商用软件引入了单元侵蚀算法，然而，该算法将畸变单元直接删除，因而无法模拟侵彻过程中产生的碎片。光滑粒子流体动力学（smoothed particle hydrodynamics，SPH）方法能有效地避免网格畸变，并且能自然地模拟材料的大变形、飞溅等现象，但是其计算效率相对较低。FEM和SPH方法在模拟强冲击问题时各有优缺点，为了集两者优势于一体，近年来，诸多学者^[10-18]对FEM-SPH耦合算法及其在高速冲击问题中的应用进行了系统的研究，但较少涉及多层防护结构的超高速碰撞问题。

本研究在He等^[19]研究的基础上，基于FEM-SPH自适应耦合算法，开展不同质量和撞击速度的动能块超高速碰撞多层防护结构的毁伤特性数值模拟，以期再现碰撞过程中产生的碎片，并应用量纲分析方法讨论动能块的质量和撞击速度对靶板穿孔直径的影响规律。

1.   数值模拟算法及模型验证

1.1   FEM-SPH自适应耦合算法

利用LS-DYNA软件对动能块超高速碰撞多层防护结构的穿孔特性进行数值模拟。LS-DYNA软件中，可以直接使用关键字*DEFINE_ADAPTIVE_SOLID_TO_SPH实现FEM-SPH自适应耦合方法。采用FEM建模，用六面体单元划分计算模型。对单元定义材料属性、失效判据、边界条件及接触条件等计算参数。此外，提前对计算中可能发生的失效转变为粒子的单元定义相应的关键字。当某单元失效后，删除该单元，并用SPH粒子代替，SPH粒子继承对应单元的质量、速度、材料、力等参数，继续参与后续计算。在每个计算时间步内，均进行一系列相同的操作，具体过程如图1所示。

图  1  FEM-SPH自适应算法的计算过程^[19]

Figure  1.  Calculation process of FEM-SPH adaptive algorithm^[19]

下载: 全尺寸图片 幻灯片

1.2   材料本构模型及参数

钨合金弹丸和Al2024-T351铝合金靶板均选用Johnson-Cook本构模型和Grüneisen状态方程描述^[19]。Johnson-Cook本构模型适用于高应变率、大变形以及高温情况，具有形式简单、精度高、各项物理意义明确等优点，其屈服应力的影响因素主要包括应变、应变率和温度，表达式为

式中：$\sigma$为屈服应力；A、B、n、C、m分别为材料的静态屈服极限、应变硬化模量、应变硬化指数、应变率强化参数、温度软化指数，由实验确定；${\overline{\varepsilon }}_{\mathrm{p}}$为有效塑性应变；${\dot{\varepsilon }}^{\mathrm{*}}$为归一化的有效应变率；${T}^{\mathrm{*}}$为相对温度。

Johnson-Cook本构模型采用累积损伤来考虑材料的破坏，单元的损伤度D的定义为

式中：${\Delta \varepsilon }_{\mathrm{p}}$为一个时间步长的等效塑性应变增量，${\varepsilon }_{\mathrm{f}}$为当前时刻的破坏应变。D的取值介于0～1之间，初始未损伤时，D = 0；当D = 1时，材料发生失效，单元发生破坏。

式中：D₁～D₅为材料失效参数，ε_f,min为破坏应变的最小值。

由于撞击速度较高，压缩波在动能块和靶板的自由表面反射之后产生的拉伸波强度远大于材料所能承受的极限强度，动能块和靶板内部会出现层裂现象，仅考虑Johnson-Cook失效模型尚不能模拟出层裂现象，可同时采用最大拉应力失效判据来实现材料层裂现象的模拟。此外，发生层裂时，动能块和靶板内部剩余单元之间会产生接触，在LS-DYNA软件中使用关键字*CONTACT_INTERIOR实现动能块和靶板内部单元的自接触^[19]。

钨合金和Al2024-T351铝合金材料的Johnson-Cook本构模型参数、Grüneisen状态方程参数和Johnson-Cook失效模型参数如表1^[19–22]所示，其中：ρ、μ、E、G分别为材料的密度、泊松比、弹性模量和剪切模量，${T}_{\mathrm{m}}$、${T}_{\mathrm{r}}$分别为材料熔点和环境温度，$c$为材料的体积声速，S₁为无量纲系数，${\gamma }_{0}$为Grüneisen系数，$a$为${\gamma }_{0}$的无量纲修正系数，${\sigma }_{\mathrm{p}}$为材料的最大失效拉应力。

表  1  钨合金和铝合金的材料模型参数^[19–22]

Table  1.  Material parameters of aluminum and tungsten alloy^[19–22]

Material	$\, \rho$/(g·cm–3)	$\,\mu$	E/GPa	G/GPa	A/GPa	B/GPa	n
Tungsten alloy	17.000	0.28	409.6	160.0	1.506	0.177	0.12
Al2024-T351	2.785	0.33	73.4	27.6	0.265	0.426	0.34
Material	C	m	Tm/K	Tr/K	D1	D2	D3
Tungsten alloy	0.016	1.0	1723	300	1.5	0	0
Al2024-T351	0.015	1.0	775	300	1.0	0	0
Material	D4	D5	c/(km·s−1)	S1	${\gamma }{_{0}}$	a	${\sigma }{_{\mathrm{p} }}$/GPa
Tungsten alloy	0	0	4.029	1.237	1.54	0.134	3.5
Al2024-T351	0	0	5.328	1.338	2.00	0.875	2.6

下载: 导出CSV 

| 显示表格

1.3   模型验证

数值模拟所采用的算法、材料模型及参数在很大程度上影响数值模拟结果的准确性和可信度，因此，对所采用的算法和材料模型进行了有效性验证。考虑到关于钨合金材料的动态本构参数研究较为成熟，可以直接参考，为此本研究仅对Al2024-T351铝合金材料的超高速碰撞数值模拟结果进行验证^[23]。

图2显示了直径为10 mm的Al2024-T351铝合金球弹以5.941 km/s的速度正撞击2 mm厚的Al2024-T351铝合金板的数值模拟结果与实验结果^[24]的比较。从图2可以看出，弹丸撞击靶板后15.9 µs时的数值模拟结果与实验测得的碎片云形态吻合较好。表2给出了15.9 µs时穿孔直径d_h、碎片云轴向速度v_a和径向速度v_r、碎片云长度L_d和直径d_d的数值模拟结果和实验数据，相对误差均小于6%。由此验证了本研究所采用的数值模拟算法、材料模型和参数的准确性和可信度。

图  2  15.9 µs时刻碎片云数值模拟与实验结果^[24]的比较

Figure  2.  Comparison between numerical simulation and experimental result^[24] of debris cloud at 15.9 µs

下载: 全尺寸图片 幻灯片

表  2  15.9 μs时的数值模拟结果与实验数据^[24]的比较

Table  2.  Comparison between numerical simulation and experimental result^[24] at 15.9 μs

Method	dh/cm	va/(m·s−1)	vr/(m·s−1)	Ld/cm	dd/cm
Simulation	2.00	5320	1860	8.16	6.57
Experiment	1.89	5296	1913	8.11	6.56
Error/%	5.82	0.45	2.77	0.62	0.15

下载: 导出CSV 

| 显示表格

2.   数值模拟结果及分析

为了研究动能块的质量和撞击速度对多层防护结构抗侵彻性能的影响，保持其他参数不变，动能块形状选择球形，材料为钨合金，动能块质量范围为5～250 g，撞击速度范围为3～7 km/s，靶板设置为10 mm间隔式多层防护结构。多层防护结构均由17层280 mm×280 mm的Al2024-T351铝合金靶板构成，第1层厚度为20 mm，其余厚度均为2 mm，靶板中心大变形区域网格加密，建立1/4对称模型，靶板外围采用全自由度固定约束条件，对称面采用对称约束。接触方式采用自动面面和自动点面接触，SPH的人工黏性选择Monaghan形式，采用正撞击方式。有限元模型如图3所示，数值模拟工况设置见表3，其中：m_p、d_p、v_p分别为动能块的质量、直径和撞击速度。

图  3  数值模拟的有限元模型

Figure  3.  Finite element model of numerical simulation

下载: 全尺寸图片 幻灯片

表  3  数值模拟工况

Table  3.  Conditions of numerical simulation

Group	No.	mp/g	dp/mm	vp/(km·s−1)		Group	No.	mp/g	dp/mm	vp/(km·s−1)
A	1	5	8.16	3		B	1	5	8.16	3
	2	45	16.96	3			2	5	8.16	4
	3	85	20.98	3			3	5	8.16	5
	4	125	23.84	3			4	5	8.16	6
	5	165	26.16	3			5	5	8.16	7
	6	205	28.12	3
	7	250	30.04	3

下载: 导出CSV 

| 显示表格

2.1   动能块质量的影响

通过A组工况研究不同动能块质量对多层防护结构超高速碰撞特性的影响。图4为工况A1～A7中动能块刚好穿透多层防护结构时刻的数值模拟结果。

图  4  工况A1～A7的数值模拟结果

Figure  4.  Numerical simulation results of different conditions (Case A1–A7)

下载: 全尺寸图片 幻灯片

从图4可以看出：在所有工况中，动能块均能成功贯穿所有17层铝合金板，并在靶后形成碎片云；随着动能块质量的增加，靶板的穿孔形貌近似为纺锤形；当动能块的质量较小时，贯穿所有铝合金板之后，动能块几乎完全破碎成小碎片，并且在撞击初始阶段，动能块内部均出现不同程度的层裂现象。以工况A7为例，在动能块与第1层铝合金板撞击阶段，动能块和第1层铝合金板内部的冲击波传播过程如图5所示。

图  5  动能块和铝合金板内部的冲击波传播（剖视图）

Figure  5.  Propagation of shock wave in kinetic energy block and aluminum alloy plate (sectional view)

下载: 全尺寸图片 幻灯片

从图5可以看出，当动能块撞到第1层铝合金板时会产生冲击波，冲击波同时向动能块和铝合金板中传播。对于球形动能块来说，冲击波首先在动能块侧向反射形成稀疏波，将侧向的冲击波卸载；而轴向冲击波继续向动能块背面传播，到达动能块背部自由面后反射形成稀疏波；在第1层铝合金板内部，当轴向冲击波到达铝合金板背部自由表面时，反射形成的稀疏波会追赶动能块内冲击波并将其卸载；之后铝合金板内冲击波沿着径向传播，并在铝合金板四周固定端不断反射形成稀疏波，从而引起铝合金板的破碎^[25]。由图5可知，在此工况下，第1层铝合金板和球形动能块内部出现层裂的时间大约在弹靶撞击后4和7 µs。

为了更直观地观测铝合金板的破坏情况，对铝合金板的穿孔情况进行分析，工况A1～A7中各层铝合金板的穿孔直径d_h如图6所示。从图6可以看出，每层靶板的穿孔直径大致都随着动能块质量的增加而增大。选取工况A1～A7中第1层铝合金板的穿孔直径数据，对其进行无量纲处理，如表4所示，其中：d_h,max为最大穿孔直径，d_h,min为最小穿孔直径。

表  4  第1层铝合金板的穿孔直径

Table  4.  Perforation diameter of the first layer of aluminum alloy plate

Case	${d}{_{\mathrm{p} } }$/mm	${d}{_{\mathrm{h},\mathrm{m}\mathrm{a}\mathrm{x} }}$/mm	${d}{_{\mathrm{h},\mathrm{m}\mathrm{i}\mathrm{n} } }$/mm	${d}{_{\mathrm{h} } }$/mm	${d}{_{\mathrm{h} } }/{d}{_{\mathrm{p} } }$
A1	8.16	22.72	22.68	22.70	2.7819
A2	16.96	39.08	39.00	39.04	2.3019
A3	20.98	47.12	45.98	46.55	2.2188
A4	23.84	51.34	50.96	51.15	2.1456
A5	26.16	56.80	56.72	56.76	2.1697
A6	28.12	58.50	57.74	58.12	2.0669
A7	30.04	62.08	59.62	60.85	2.0256

下载: 导出CSV 

| 显示表格

图  6  各层铝合金板的穿孔直径统计（工况A1～A7）

Figure  6.  Statistics of perforation diameter of each layer of aluminum alloy plate (Case A1–A7)

下载: 全尺寸图片 幻灯片

对于第1层铝合金板来说，影响其穿孔直径的物理量主要有两类：

(1) 动能块参数，包括直径${d}_{\mathrm{p}}$、撞击速度${v}_{\mathrm{p}}$、密度${\rho }_{\mathrm{p}}$、屈服强度${\sigma }_{\mathrm{p}}$、材料声速${c}_{\mathrm{p}}$；

(2) 铝合金板参数，包括厚度${T}_{\mathrm{t}}$、密度${\,\rho }_{\mathrm{t}}$、屈服强度${\sigma }_{\mathrm{t}}$、材料声速${c}_{\mathrm{t}}$。

基于以上各物理量参数，第1层铝合金板的穿孔直径${d}_{\mathrm{h}}$与各物理量之间的函数关系可写为

根据量纲分析的 $\text{π }$ 定理^[26]，取动能块直径${d}_{\mathrm{p}}$、铝合金板密度${\rho }_{\mathrm{t}}$、铝合金板屈服强度${\sigma }_{\mathrm{t}}$作为基本量，可得到以下无量纲函数关系

穿孔直径${d}_{\mathrm{h}}$与动能块直径${d}_{\mathrm{p}}$之比称为无量纲穿孔直径。在A组工况中，变化的参数为动能块质量，即动能块直径，动能块的撞击速度和材料以及铝合金板的材料和厚度均不发生变化，不考虑材料可压缩性的影响，即

则式(5)可简化为幂次关系，即

式中：$\alpha$、$\beta$均为待定常数。

利用表4中的数值模拟数据对式(7)进行回归分析，可得到球形动能块撞击第1层铝合金板的穿孔直径经验公式

式(8)的适用范围：动能块质量5～250 g，撞击速度3 km/s。

根据式(8)计算得到的第1层铝合金板穿孔直径与数值模拟结果见表5，可见，相对误差均小于5%，且正负比相对均衡，经验公式与数值模拟结果吻合较好。

表  5  第1层铝合金穿孔直径的计算数据与数值模拟结果的对比

Table  5.  Comparison between calculation and simulation of perforation diameter of the first layer of aluminum alloy plate

mp/g	${d} {_{\mathrm{h} } }/{d}{_{\mathrm{p} } }$		Error/%
	Calc.	Sim.
25	2.4437	2.4835	−1.60
65	2.2683	2.2320	1.63
105	2.1854	2.0982	4.16
145	2.1312	2.0571	3.60
185	2.0912	2.1155	−1.15
225	2.0594	2.0255	1.67

下载: 导出CSV 

| 显示表格

2.2   动能块撞击速度的影响

通过B组工况研究不同撞击速度对多层防护结构超高速碰撞特性的影响，图7为工况B1～B5中刚好穿透多层防护结构时刻的数值模拟结果。从图7可以看出，在所有工况下，动能块均能成功贯穿所有17层铝合金板，并在靶后形成碎片云。当撞击速度v_p为3、4、5、6、7 km/s时，碎片云头部速度v_d分别为2.06、2.53、2.57、2.88、3.05 km/s。在其他条件不变的情况下，随着动能块撞击速度的提升，形成的碎片云头部速度近似呈线性增大，拟合曲线如图8所示，拟合公式为${v}_{\mathrm{d}}=1.453+0.233{v}_{\mathrm{p}}$。

图  7  工况B1～B5的数值模拟结果

Figure  7.  Numerical simulation results of different conditions (Case B1−B5)

下载: 全尺寸图片 幻灯片

图  8  碎片云头部速度随动能块速度的变化曲线

Figure  8.  Variation of debris cloud head velocity with impact velocity

下载: 全尺寸图片 幻灯片

工况B1～B5中，各层铝合金板的穿孔直径的数值模拟结果如图9所示。从图9可以看出：在保持其他参数不变的情况下，铝合金板的最大穿孔直径均出现在第2层铝合金板上，不随撞击速度的改变而改变；对于每种工况来说，铝合金板的穿孔直径整体上呈现先增后减最后基本保持平稳的趋势。

图  9  各层铝合金板的穿孔直径统计（工况B1～B5）

Figure  9.  Statistics of perforation diameter of each layer of aluminum alloy plate (Case B1−B5)

下载: 全尺寸图片 幻灯片

选取每种工况中出现最大穿孔直径位置的第2层铝合金板的穿孔直径数据，对其进行无量纲处理，结果如表6所示。

表  6  第2层铝合金板的穿孔直径（${d}_{\mathrm{p}}$=8.16 mm）

Table  6.  Perforation diameter of the second layer of aluminum alloy plate (${d}_{\mathrm{p}}$=8.16 mm)

Case	${v}{_{\mathrm{p} } }$/(km·s−1)	${d}{_{\mathrm{h},\mathrm{m}\mathrm{a}\mathrm{x} } }$/mm	${d}{_{\mathrm{h},\mathrm{m}\mathrm{i}\mathrm{n} } }$/mm	${d}{_{\mathrm{h} } }$/mm	${d}{_{\mathrm{h} }}/{d}{_{\mathrm{p} } }$
B1	3	36.50	36.48	36.49	4.4718
B2	4	46.08	45.94	46.01	5.6385
B3	5	54.86	52.84	53.85	6.5993
B4	6	61.16	57.84	59.50	7.2917
B5	7	65.94	61.54	63.74	7.8113

下载: 导出CSV 

| 显示表格

影响第2层铝合金板穿孔直径的物理量主要有两类：

(1) 动能块参数，包括直径${d}_{\mathrm{p}}$、撞击速度${v}_{\mathrm{p}}$、密度${\rho }_{\mathrm{p}}$、屈服强度${\sigma }_{\mathrm{p}}$、材料声速${c}_{\mathrm{p}}$；

(2) 铝合金板参数，包括第1层铝合金板厚度${T}_{\mathrm{t}}$、第2层铝合金板厚度${T}_{\mathrm{t}1}$、密度${\rho }_{\mathrm{t}}$、屈服强度${\sigma }_{\mathrm{t}}$、材料声速${c}_{\mathrm{t}}$、靶板层间距L。

可以看出，分析第2层铝合金板所考虑的参数比分析第1层铝合金板时多2个物理量，分别是第2层铝合金板的厚度以及两层铝合金板的层间距。与2.1节的量纲分析过程类似，根据 $\text{π}$ 定理可得

对于B组工况，变化的参数为撞击速度，动能块直径和材料、各层铝合金板的材料和厚度以及靶板层间距均不发生变化，不考虑材料可压缩性的影响，即

则式(9)可简化为幂次关系，即

式中：$\gamma$、$\delta$均为待定常数。

利用表6中的数值模拟数据对式(11)进行回归分析，可得到球形动能块撞击第2层铝合金板的穿孔直径经验公式

式(12)的适用范围：动能块质量5 g，撞击速度3～7 km/s。

根据式(12)计算得到的第2层铝合金板穿孔直径与数值模拟结果见表7，可见，相对误差均小于10%，且正负比相对均衡，经验公式的计算结果与数值模拟结果吻合较好。

表  7  第2层铝合金板穿孔直径的计算结果与数值模拟结果的对比

Table  7.  Comparison between calculation and simulation of perforation diameter of the second layer of aluminum alloy plate

${v}{_{\mathrm{p} } }$/(km·s−1)	${d}{_{\mathrm{h} }}/{d}{_{\mathrm{p} }}$		Error/%
	Calc.	Sim.
3.5	5.0670	4.7426	6.84
4.5	5.9888	6.0515	−1.04
5.5	6.8405	6.4877	5.44
6.5	7.6414	7.3627	3.79

下载: 导出CSV 

| 显示表格

3.   结　论

利用LS-DYNA有限元软件，采用FEM-SPH自适应算法，对多层防护结构在动能块超高速碰撞下的特性进行了数值模拟，结合量纲分析方法，讨论了动能块质量和撞击速度对铝合金板穿孔直径的影响，在本研究的动能块质量和撞击速度范围内，保持其他参数不变，可得出以下结论。

(1) 不同质量的钨合金动能块均能够贯穿所有17层铝合金板，在靶后形成碎片云，并且在撞击作用过程中动能块和铝合金板内部均出现层裂现象。

(2) 对于不同质量的动能块，铝合金板的最大穿孔直径出现的位置无明显规律，第1层铝合金板的穿孔直径随着动能块质量的增大近似以幂函数${d}_{\mathrm{h}}/{d}_{\mathrm{p}}=2.246\,3\;{\left({T}_{\mathrm{t}}/{d}_{\mathrm{p}}\right)}^{0.233\,3}$的形式增大。

(3) 不同撞击速度的钨合金动能块在贯穿所有17层铝合金板的过程中，动能块和铝合金板内部均出现层裂现象，所形成的靶后碎片云的头部速度随着撞击速度的提升近似呈线性增大，拟合的函数关系为${v}_{\mathrm{d}}=1.453+0.233{v}_{\mathrm{p}}$。

(4) 在不同的撞击速度下，铝合金板的最大穿孔直径位置均出现在第2层，并且第2层铝合金板的穿孔直径随着撞击速度的提升近似以幂函数${d}_{\mathrm{h}}/{d}_{\mathrm{p}}=1.013\,8\;{\left({v}_{\mathrm{p}}/\sqrt{{\sigma }_{\mathrm{t}}/{\rho }_{\mathrm{t}}}\right)}^{0.662\,7}$的形式增大。

(5) 考虑到靶后的碎片云主要由动能块和破碎的铝合金板组成，铝合金板的穿孔直径越大，表明消耗的动能越大，参与生成碎片云的铝合金板质量越大。因此，研究铝合金板的穿孔特性可为后期分析靶后碎片云的质量和速度分布、建立靶后的冲击载荷模型、研究靶后带板装药的冲击起爆规律奠定基础。