超高速碎片云撞击下薄板变形与撕裂建模研究

刘泽荣; 龙仁荣; 张庆明; 陈利

doi:10.11858/gywlxb.20210811

超高速碎片云撞击下薄板变形与撕裂建模研究

doi: 10.11858/gywlxb.20210811

北京理工大学爆炸科学与技术国家重点实验室，北京　100081

基金项目: 国防科工局“十三五”碎片专项（KJSP2016030301）；国防科工局民用航天技术预研项目（D020304）

详细信息

作者简介:
刘泽荣（1995－），男，硕士研究生，主要从事冲击动力学研究. E-mail：rongzixu@126.com

通讯作者:
龙仁荣（1982－），男，博士，副教授，主要从事冲击动力学研究. E-mail：longrenrong@bit.edu.cn

中图分类号: O347; V423.42
计量
- 文章访问数: 878
- HTML全文浏览量: 604
- PDF下载量: 57
出版历程
- 收稿日期: 2021-06-11
- 修回日期: 2021-06-29

Theoretical Study of Deflecting and Petalling of Thin Plate under Debris Cloud Loading Induced by Hypervelocity Impact

State Key Laboratory of Explosion Science and Technology, Beijing Institute of Technology, Beijing 100081, China

摘要

摘要: 为研究多层板结构中薄板在碎片云作用下的变形与破坏问题，开展了超高速撞击多层板实验。实验结果表明，薄板在高速碎片云冲击下的典型破坏特征为中央穿孔及环孔凹陷变形与花瓣型撕裂。在此基础上，考虑弯矩和膜力作用，建立了描述薄板在轴对称分布强冲击载荷作用下大变形的理想刚塑性环板模型，据此可以计算环板变形过程的横向与径向速度场，结合Grady破碎理论，可以计算花瓣型撕裂的花瓣数，理论计算值与实验比较吻合。研究结果可以为多层板结构在超高速弹丸撞击下的毁伤评估提供理论基础。
- 超高速撞击 /
- 碎片云 /
- 环板动态变形 /
- 花瓣型撕裂
Abstract: In order to study the deformation and failure of the thin plate under impact of debris cloud in the multi-layer structure, a series of hypervelocity impact experiments on multi-layer structure are carried out. The experimental results showed that typical failure characteristics of the thin plate under impact of the high-speed debris cloud are central perforation, deflection and petalling. Considering the effect of bending moment and membrane force, an ideal rigid-plastic annular plate deflecting model is established, which can be used to describe large deflection of a thin plate under the strong axisymmetric impact load. The transverse and radial velocity fields of the thin plate under the strong axisymmetric impact load can also be calculated by the established ideal rigid-plastic annular plate deflecting model. With the Grady fragmentation theory, the number of petals can be calculated. The theoretical calculation result is consistent with the experiment. The research results can provide a basic theory for the damage assessment of the multilayer structure under the impact of hypervelocity projectiles.
- hypervelocity impact /
- debris cloud /
- dynamic deflection of annular plate /
- petalling

HTML全文

多层板结构广泛存在于航天器、舰船等结构中^[1]，然而航天器会遭受在轨超高速空间碎片的撞击，舰船等目标会遭受超高速动能武器的撞击^[2]。在超高速弹丸撞击下，多层板结构中的第1层板形成穿孔，第2层及其后各层板则承受碎片云的冲击，即分布式载荷。碎片云的破坏形式与超高速弹丸明显不同，其破坏特征与尺寸是毁伤评估的重要参数，研究碎片云撞击下薄板的变形与破坏十分必要。

碎片云碰撞下薄板的变形与破坏可以看作薄板在具有空间分布的瞬态载荷下发生塑性动态大变形问题。薄壁结构是一种常见的结构，薄板变形理论由来已久，其中薄板小变形理论研究相对完善^[3]，而薄板动态塑性大变形理论仍在发展中。在研究薄板的大变形过程时，一方面需要同时考虑弯矩和膜力的作用，所建模型相对复杂，另一方面其几何方程具有非线性项，使得位移求解非常困难。Craggs^[4]提出用薄膜理论分析薄板大变形，在模型中不考虑弯矩作用而只考虑膜力，这虽然可以简化理论模型，但是却只适用于很薄且强度较低的薄板大变形挠曲问题，难以描述常用金属薄板大变形问题^[5]。Reissner^[6]首先提出板壳大变形基本理论，随后Jones^[7-8]应用此理论求解了理想刚塑性板在动载荷下的挠曲变形，通过忽略面内变形，从而简化方程得到理论解。Yu等^[9]引入膜力因子，将膜力对板变形的影响通过膜力因子乘以平衡方程中的弯矩体现出来，由于方程中不包含膜力项，求解薄板挠曲变形方程的难度大大降低。然而，这些模型忽略或简化了膜力作用，难以实现大变形高精度的计算。

本研究利用二级轻气炮开展球形弹丸超高速正撞击多层板实验，分析总结薄板在碎片云撞击下的变形与破坏特征，在此基础上建立考虑弯矩和膜力作用的环板变形模型，将大变形的薄板材料运动分解为横向的外凸运动和环向的膨胀运动，并结合Grady破碎理论，建立环向膨胀花瓣撕裂模型。

1. 超高速撞击多层板实验及结果分析

采用二级轻气炮将直径d_p=9.13 mm的LY12铝合金球形弹丸加载至速度v₀=5 km/s后撞击多层板防护结构。多层板防护结构由4块尺寸为100 mm×100 mm的LY12铝合金薄板组成，相邻两板间距为：S₁=17 mm，S₂=50 mm，S₃=38 mm；第1层至第4层板的厚度分别为：h_b1=1.0 mm，h_b2=2.0 mm，h_b3=2.0 mm，h_b4=1.5 mm，如图1所示。

图 1 弹丸超高速碰撞多层板结构

Figure 1. Multi-shock shield structure subjected to hypervelocity impact

下载: 全尺寸图片幻灯片

实验后各层板的照片如图2所示，第1层板只发生穿孔破坏，靶板中心有一个规整圆孔。第1层与第2层板的间距较小，弹丸与第1层板碰撞形成的一次碎片云没能充分膨胀，碎片云载荷依然较为集中，因此第2层板中央形成了一个大的穿孔，孔四周出现一个由小碎片撞击形成的环状小坑区域。第3层板与第4层板的破坏特征相似，板中间有穿孔破坏，孔周边材料发生拉伸完全变形，并且发生撕裂破坏，如图2(c)、图2(d)、图2(e)所示。第3层板弯曲撕裂至靶板边界，边沿夹持的卡具也发生了弯曲大变形。由于第3层板的阻挡，作用到第4层板的碎片云能量大幅减小，第4层板的弯曲变形主要集中在中部区域，边沿和卡具均没有明显的弯曲变形。由于裂纹扩展所需要的能量较小，第4层板的裂纹依然扩展至靶板边沿。第4层板表面覆盖一层白色物体，弹坑也相对较少，经分析，白色物是熔融铝冷却凝固物，熔融态的碎片云对第4层板的作用更趋于整体作用，而不是大尺寸固态碎片的局部作用。整个多层板结构在超高速弹丸撞击下的变形特征如图1所示。

图 2 高速冲击后各层板的破坏照片

Figure 2. Damage results of each layer after hypervelocity impact

下载: 全尺寸图片幻灯片

实验结果表明，弹丸与集中分布的碎片云撞击薄板会造成穿孔破坏，随着撞击次数的增加，弹靶材料破碎得更细，且熔融态的比例增加，碎片云载荷分布越分散，薄板的变形与破坏模式转变为中心穿孔及环孔区域大凹陷变形并伴有撕裂。第3层板阻挡了大部分碎片云能量，第4层板的弯曲变形比第3层板小很多，其中心穿孔直径的测量平均值为46 mm，凹陷变形区域直径为穿孔直径的2～3倍，凹陷板与靶板原平面的夹角较小，约为12°，说明其凹陷变形的速度较小。薄板穿孔四周材料的凹陷变形可简单类比成圆环的膨胀变形，当膨胀速度大于某值后，圆环在膨胀过程会断裂成多段。材料起裂后裂纹在板中传播形成长裂纹，与实验中薄板撕裂成多块的结果一致。根据前述分析，第4层板环孔凹陷区域的变形速度较低，其环向膨胀更趋近于薄板起裂的临界速度。

2. 碎片云作用下环形薄板变形模型

由实验结果总结得到在分散碎片云载荷下薄板的响应特征如下：薄板中间区域材料在碎片云中部高幅值冲击载荷下发生穿孔，环孔区域的环向板承受的碎片云载荷稍微小一点，环板材料发生凹陷变形，孔边材料在凹陷膨胀过程中破碎成几段，裂纹扩展直至薄板裂成几块。碎片云速度分布与环板变形如图3所示。碎片云撞击薄板穿孔尺寸的计算已有研究^[10]，在此不做考虑，本研究侧重碎片云载荷下环形薄板的变形建模。

图 3 环板变形侧视图

Figure 3. Lateral schematic of annular plate deformation

下载: 全尺寸图片幻灯片

2.1 薄板变形假设与一般方程

根据实验结果可知，碎片云分布载荷作用下的变形属于典型的薄板大变形问题。薄板在变形过程中发生弯曲，需要考虑弯矩作用，另外薄板面内还有膜应力。要建立薄板变形的一般方程，需作如下基本假设：

(1) 考虑到薄板挠度是板厚的几倍甚至几十倍，为此忽略板的弹性小变形，采用理想刚塑性体模型描述薄板变形；

(2) 薄板变形服从Love-Kirchhoff假设^[11]，即中面法线变形后仍保持为中面法线，厚度变化忽略不计，同时略去中面法向应力；

(3) 不区分变形前后参考位形的差异，采用变形前的位形作为参考基准建立平衡方程；

(4) 薄板内力除径向弯矩M_r和环向弯矩M_θ外，还需要考虑径向膜力N_r和环向膜力N_θ。

基于以上假设，可以建立轴对称薄板变形过程的横向、环向和径向动力学方程^[12]

$\frac{\partial }{{\partial r}}(r{Q_r}) + \frac{\partial }{{\partial r}}\left(r{N_r}\frac{{\partial w}}{{\partial r}}\right) + pr - \mu r\frac{{{\partial ^2}w}}{{\partial {t^2}}} = 0$

(1)

$\frac{\partial }{{\partial r}}(r{N_r}) + {N_\theta } - \mu r\frac{{{\partial ^2}u}}{{\partial {t^2}}} = 0$

(2)

$\frac{\partial }{{\partial r}}(r{M_r}) + {M_\theta } - r{Q_r} = 0$

(3)

式中：r为径向坐标，θ为环向坐标，u为径向位移，w为挠度，Q_r为剪力，μ为板的面密度，p为轴对称的横向分布载荷。

考虑面内变形与出面变形的薄板轴对称几何方程，其内力（或称广义应力^[3]）径向弯矩M_r、环向弯矩M_θ、径向膜力N_r和环向膜力N_θ对应的广义应变分别为径向膜应变ε_r、环向膜应变ε_θ、径向曲率κ_r和环向曲率κ_θ。几何方程如下

${\varepsilon _r} = \frac{{\partial u}}{{\partial r}} + \frac{1}{2}{\left( {\frac{{\partial w}}{{\partial r}}} \right)^2}$

(4)

${\varepsilon _\theta } = \frac{u}{r}$

(5)

${\kappa _r} = - \frac{{{\partial ^2}w}}{{\partial {r^2}}}$

(6)

${\kappa _\theta } = - \frac{1}{r}\frac{{\partial w}}{{\partial r}}$

(7)

采用理想刚塑性体模型，材料发生刚性位移或塑性变形，则薄板发生塑性变形时膜力N_r、N_θ达到N₀，弯矩M_r、M_θ达到M₀

${N}_{0}={\sigma }_{\rm{y}}h ,\;\;\;\;{M}_{0}={\sigma }_{\rm{y}}{h}^{2}/4$

(8)

式中：h为板厚度，σ_y为屈服强度。

薄板材料是理想刚塑性体模型，薄板的应力状态位于屈服面上，为此采用Hodge^[13]提出的包含弯矩、膜力作用的近似屈服准则，如图4所示。可以看出，屈服面上的膜力和弯矩独立作用，各自位于Tresca屈服面上。这样的近似处理使求解过程变得非常简便。另外，横向剪力Q_r出现在平衡方程中，但是由于不考虑其产生的剪切变形，因此不出现在屈服条件中，并且后面计算中可以看到，化简的平衡方程不包含剪力。

图 4 Tresca屈服条件与正交法则

Figure 4. Tresca condition and flow rule

下载: 全尺寸图片幻灯片

图4中m_r、m_θ、n_r、n_θ为无量纲内力，m_r = M_r/M₀，m_θ = M_θ/M₀，n_r = N_r/N₀，n_θ = N_θ/N₀；A、B、C、D、E、F是屈服条件上的角点，各线段代表材料的某一内力状态。正交法则规定一点广义塑性应变增量的方向与屈服面外法线方向相同，如图4所示。对于具体的问题，在约束允许的机动场中，根据广义塑性应变增量（或广义塑性应变率）的方向，可以确定塑性变形由哪一部分屈服面控制，从而为具体问题的求解确定内力大小。

2.2 环板变形方程简化

根据薄板变形方程，结合屈服条件和正交法则，就可以求解薄板变形的初边值问题，从而得到环板的变形。然而多个偏微分方程的相互耦合以及非线性项的引入，使方程很难得到解析解，需要根据本问题的特征对方程进行简化。理想刚塑性体的使用、碎片云撞击载荷的特性以及环板构型的特殊性会使方程变量减少，复杂程度降低。

如图5所示，薄板中部发生穿孔（0 < r < a，其中r为径向坐标，a为穿孔半径），形成圆环薄板，其中剪力做功转化为环板动能；内环之外的环状区域发生塑性凹陷（a ≤ r ≤ l，其中l为塑性凹陷半径），当环板的环向变形达到一定程度时，内环边发生破碎，成为薄板花瓣型撕裂破坏的起始点；其他区域（r > l）在本问题中则可认为是未变形区。可以认为碎片云冲击载荷为瞬态载荷作用于环板面上。根据已有的碎片云模型^[10]，得到作用在环板上的冲量面密度分布i = i(r)，碰撞过程瞬时完成，将冲量载荷转化为环板横向速度初始条件v(r) = i/μ(t = 0)，如图3所示。这样将碎片云碰撞薄板问题简化为圆环板在初始横向速度条件下的变形问题。

图 5 环板变形俯视图

Figure 5. Top schematic of annular plate deformation

下载: 全尺寸图片幻灯片

刚塑性环板挠曲变形采用锥形机动场^[14]，即a ≤ r < l处挠度w满足

$w = \frac{{{w_a}}}{{1 - \alpha }}\left(1 - \frac{r}{l}\right)$

(9)

$\dot w = \frac{{{{\dot w}_a}}}{{1 - \alpha }}\left(1 - \frac{r}{l}\right)$

(10)

式中：w_a和 $\dot{{w}_{a}}$ 分别为内环边挠度和横向速度； $\alpha$ 为无量纲长度，α = a/l；为方便表示，以下统一规定 $\dfrac{\partial f}{\partial t}=\dot{f},\;\dfrac{\partial f}{\partial r}={f}{'}$ 。如图5所示，理想刚塑性体变形仅发生在塑性铰上，分别是位于r = l的塑性铰圆<1>和a ≤ r < l的径向塑性铰线<2>。

在a ≤ r < l区域，由式(6)和式(10)计算得到径向曲率变化率 $\dot{{\kappa }_{r}}$

${\dot \kappa _r} = - {(\dot w)^{\prime \prime }} = 0$

(11)

面内变形只在塑性铰线<2>上产生，对于环板，令<2>上的径向应变率 $\dot{{\varepsilon }_{r}}=0$ ^[15]，由式(4)得

${\varepsilon _r} = {u^\prime } + \frac{1}{2}{({w^\prime })^2} = 0$

(12)

${\dot \varepsilon _r} = {\dot u^\prime } + {w^\prime }{\dot w^\prime } = 0$

(13)

对径向坐标在[0, r]积分式(13)，得

$\dot u = {\dot u_a} + \frac{{l - r}}{{{{(l - a)}^2}}}{w_a}{\dot w_a}$

(14)

式中： ${\dot{u}}_{a}$ 为环板内环边径向速度。未知的速度变量 $\dot{w}(r,t)$ 和 $\dot{u}(r,t)$ 变为内环边速度 $\dot{{w}_{a}}\left(t\right)$ 和 ${\dot{u}}_{a}\left(t\right)$ 。下面进一步化简以减少未知数。

图6显示了刚塑性环板位移场，可得径向位移u_a

图 6 挠度与径向关系示意图

Figure 6. Schematic of relationship between deflection and radial displacement

下载: 全尺寸图片幻灯片

${u_a} = (l - a)(1 - \cos \theta )$

(15)

式中：θ为环板塑性铰处转角。对时间求导，得到速度 ${\dot{u}}_{a}$

${\dot u_a} = (l - a)\sin \theta \,\dot \theta$

(16)

代入 ${w}_{a}=(l-a)\sin \theta$ 以及线性的 ${\dot{w}}_{a}=(l-a)\dot{\theta }$ ，得

${\dot u_a} = \frac{{{w_a}{{\dot w}_a}}}{{l - a}}$

(17)

代入式(14)，得

$\dot u = \frac{{{{\dot u}_a}}}{{1 - \alpha }}\left(1 - \frac{r}{l}\right)$

(18)

则速度未知变量数降为1，从而方便之后求解。

下面分析环板变形内力。塑性铰线<2>上有 $\dot{{\varepsilon }_{r}}=0$ 和 ${\dot{\kappa }}_{r}=0$ ，如图4所示：广义塑性应变率 $(\dot{{\kappa }_{\theta }},\dot{{\kappa }_{r}})$ 和 $(\dot{{\varepsilon }_{\theta }},\dot{{\varepsilon }_{r}})$ 的方向分别为A₁B₁和A₂B₂线段上的箭头方向。根据正交法则，<2>上的应力状态由屈服曲线的一部分A₁B₁和A₂B₂决定，即M_θ = M₀，N_θ = N₀，M_r = 0，N_r = 0。

由驻定塑性铰上的边界条件，塑性铰圆<1>上有 $\dot{{\varepsilon }_{\theta }}=0$ 和 ${\dot{\kappa }}_{\theta }=0$ ，如图4所示：广义塑性应变率 $(\dot{{\kappa }_{\theta }},\dot{{\kappa }_{r}})$ 和 $(\dot{{\varepsilon }_{\theta }},\dot{{\varepsilon }_{r}})$ 的方向分别为E₁F₁和B₂C₂线段上的箭头方向。根据正交法则，<1>上的应力状态由屈服曲线的一部分E₁F₁和B₂C₂决定，即M_θ = 0，N_θ = 0，M_r = M₀。刚塑性环板应力状态已知，唯一未知变量为速度变量，则通过动力学方程即可求解。

2.3 动力学方程求解

由于只有单一未知量，使用能量方程微分形式求解更方便，即塑性铰上塑性功率与动能变化率平衡

${\dot W\rm_p} + {\dot E\rm_d} = 0$

(19)

${\dot W\rm_p} = {\dot W_{{M_r}}} + {\dot W_{{M_\theta }}} + {\dot W_{{N_\theta }}}$

(20)

式中： ${\dot{W}}_{\rm p}$ 为塑性功率， ${\dot{E}}_{\rm d}$ 为动能变化率， ${\dot{W}}_{{M}_{r}}$ 为<1>径向弯矩功率， ${\dot{W}}_{{M}_{\theta }}$ 为<2>环向弯矩功率， ${\dot{W}}_{{N}_{\theta }}$ 为<2>环向膜力功率。在r = l的圆周上，径向弯矩功率 ${\dot{W}}_{{M}_{r}}$ 为

$\begin{split} {{\dot W}_{{M_r}}} &= \displaystyle\oint_l {{M_r}{{\dot \kappa }_r}{\rm d}r} = \displaystyle\int_0^{2{\text{π}} } {{M_0}\frac{{{{\dot w}_a}}}{{l - a}}l{\rm d}\theta } = \frac{{2{\text{π}} }}{{1 - \alpha }}{M_0}{{\dot w}_a} \end{split}$

(21)

对于a < r < l上的环向弯矩功率，由式(7)得

$\begin{split} {{\dot W}_{{M_\theta }}} &= \iint_A {{M_\theta }{{\dot \kappa }_\theta }{\rm d}A} = \int_a^l { - {M_0}\frac{{{{\dot w}^\prime }}}{r}2{\text{π}} r{\rm d}r} = 2{\text{π}} {M_0}{{\dot w}_a} \end{split}$

(22)

对于a < r < l上的环向膜力功率，由式(5)得

$\begin{split} {{\dot W}_{{N_\theta }}} &= \iint_A {{N_\theta }{{\dot \varepsilon }_\theta }{\rm d}A} = \int_a^l {{N_0}\frac{{\dot u}}{r}2{\text{π}} r{\rm d}r} = 2{\text{π}} {N_0}\frac{{\rm{d}}}{{{\rm{d}}t}}\frac{{w_a^2}}{4} \end{split}$

(23)

环板整体的动能

$\begin{split} {E\rm_d} &= \iint_A {\left( {\frac{1}{2}{{\dot u}^2} + \frac{1}{2}{{\dot w}^2}} \right)\mu {\rm d}A} = \int_a^l {\left[{\frac{1}{2}{{\dot u}_a}^2{{\left( {\frac{{1 - r/l}}{{1 - \alpha }}} \right)}^2} + \frac{1}{2}{{\dot w}_a}^2{{\left( {\frac{{1 - r/l}}{{1 - \alpha }}} \right)}^2}} \right]\mu 2{\text{π}} r{\rm d}r} \\ &= \frac{{2{\text{π}} \mu h{l^2}}}{6} \, \frac{{1 - 3{\alpha ^2} + 2{\alpha ^3}}}{{1 - \alpha }}\left( {\frac{1}{2}\dot u_a^2 + \frac{1}{2}\dot w_a^2} \right) \end{split}$

(24)

对时间求导，可得动能变化率

${\dot E\rm_d} = \frac{{2{\text{π}} \mu h{l^2}}}{6} \, \frac{{1 - 3{\alpha ^2} + 2{\alpha ^3}}}{{1 - \alpha }}\frac{{\rm{d}}}{{{\rm{d}}t}}\left( {\frac{1}{2}\dot u_a^2 + \frac{1}{2}\dot w_a^2} \right)$

(25)

将式(21)～式(23)代入式(20)，并将式(20)、式(25)代入式(19)，得到关于w_a的微分方程

$\frac{{{N_0}}}{2}\frac{{\rm{d}}}{{{\rm{d}}t}}\left( {\frac{1}{2}{w_a^2}} \right) + A{M_0}{\dot w_a} + B\frac{{\rm{d}}}{{{\rm{d}}t}}\left[ {\frac{1}{2}{{\left( {\frac{{\rm{d}}}{{{\rm{d}}t}}\frac{{{w_a^2}}}{{2l(1 - \alpha )}}} \right)}^2} + \frac{1}{2}{{\dot w}_a^2}} \right] = 0$

(26)

式中： $A = \dfrac{{2 - \alpha }}{{1 - \alpha }}$ ， $B = \dfrac{{\mu h{l^2}}}{6} \, \dfrac{{1 - 3{\alpha ^2} + 2{\alpha ^3}}}{{1 - \alpha }}$ 。联立本问题的初始条件，即可解出w_a(t)，代入式(10)可以得到横向位移场 $w(r,t)$ 。

可以看出，式(26)为非线性微分方程，不易获得解析形式解。然而对于环形板在碎片云载荷作用下的变形问题，更多地关注内环环向应变率，因此可不直接求解式(26)，由式(17)可知，需要求出内环横向速度 ${\dot{w}}_{a}$ 。求式(26)对时间t在[0, t]上的积分，并代入初始条件：W_p = 0，E_d = E_d0，t = 0，则

${W\rm_p} + {E\rm_d} - {E\rm_{d0}} = 0$

(27)

从而推导出

$\frac{{{N_0}}}{2} {\frac{1}{2}{w_a^2}} + A{M_0}{w_a} + B\frac{1}{2}{{\dot w}_a}^2\left[ {1 + {{\left( {\frac{{{w_a}}}{{l(1 - \alpha )}}} \right)}^2}} \right] - \frac{{E\rm_d}_0}{2{\text{π}} } = 0$

(28)

由此得到横向速度 ${\dot{w}}_{a}$ 关于挠度w_a的关系

${\dot w_a}^2 = - \frac{2}{B}\frac{{({N_0}/4)w_a^2 + A{M_0}{w_a} - {E_{\rm d0}}/(2{\text{π}} )}}{{1 + {{\left[ {\dfrac{{{w_a}}}{{l(1 - \alpha )}}} \right]}^2}}}$

(29)

式(29)说明，当环板内环边变形至挠度为w_a时，内环边横向速度为 ${\dot{w}}_{a}$ 。式(29)中塑性铰圆半径l由碎片云模型作用于靶板的范围确定^[10]，初始动能E_d0由碎片云冲量分布得出，通过能量守恒可求出初始动能E_d0。碎片云冲击薄板时，环板获得的初始动能E_d0由内环边剪切力做功W_s和环板受冲击获得的动能 ${E}_{\rm d0}^{i}$ 组成，则

${E_{\rm d0}} = {W\rm_s} + E_{\rm d0}^i$

(30)

${W\rm_s} = \mathop \int \nolimits_0^h \frac{{{\sigma _y}}}{{\sqrt 3 }}\left( {h - \textit{z}} \right)2{\text{π}} a{\rm d}\textit{z}$

(31)

$E_{\rm d0}^i = \int_0^l {\frac{1}{2}{{\left( {\frac{{i\left( r \right)}}{\mu }} \right)}^2}2{\rm{{\text{π}} }}r\mu {\rm{d}}r}$

(32)

其中剪切力^[16]可认为是 $\dfrac{\sigma_{y}}{\sqrt{3}}$ 。令式(28)中 ${\dot{w}}_{a}=0$ 作为环板变形终态，将求得的塑性铰圆半径l和初始动能E_d0代入式(28)，求出环板终态内环边挠度 ${w}_{a}^{\rm t}$

$w_a^{\rm t} = \frac{{ - A{M_0} + \sqrt {{A^2}M_0^2 - {N_0}{E_{\rm d0}}/\left( {2{\rm{\text{π} }}} \right)} }}{{{N_0}/2}}$

(33)

代入式(34)得到环板凹陷变形终态角度 $\theta ^{\rm{ t}}$

$\sin {\theta \rm^t} = w_a^{\rm t}/\left( {l - a} \right)$

(34)

由式(33)和式(34)即可求出环板变形最终挠度 ${w}_{a}^{\rm t}$ 和终态角度θ^t，由式(29)可得到环板变形与速度的关系式。下面将计算环板在碎片云撞击下的具体变形问题，得出环板径向速度，进而计算环板花瓣型撕裂。

3. 计算验证

环板在碎片云冲击作用下既发生弯曲变形又有拉伸变形，环板内环边沿（r = a）以一定横向速度 ${\dot{w}}_{a}$ 和径向速度 ${\dot{u}}_{a}$ 扩张，如图7所示，环向应变率 ${\dot{\varepsilon }}_{\theta a}$ 达到某临界值 ${\dot{\varepsilon }}_{\theta a}^{\rm c}$ 时，内环发生破碎，进而裂纹扩展形成花瓣型撕裂。采用本研究提出的环板变形模型，计算与实验工况相同的靶板穿孔形成的环板变形，得到内环边变形速度，将环板内环边起裂视为膨胀环破碎问题，采用Grady破碎模型^[17]计算撕裂花瓣数，并与实验得到的第4层板结果比对。

图 7 冲量载荷下环板变形示意图

Figure 7. Deflection of annular plate subjected to impulse loading

下载: 全尺寸图片幻灯片

为计算环板撕裂起始时的变形，首先计算出环板在碎片云撞击下达到的最大径向速度。为得到环板变形时材料受力状态，将板厚h = 1.5 mm、材料屈服强度 $\sigma$ _y = 265 MPa代入式(8)，求得刚塑性弯矩M₀和刚塑性膜力N₀。根据多次碎片云模型^[10]，可以计算本研究实验工况下3次碎片云对第4层环板施加的冲量面密度i（见图8）以及碎片云作用范围的半径l = 47.2 mm，再通过式(30)～式(32)求得碎片云撞击环板时环板的初始动能E_d0 = 97.5 J。由式(29)和式(17)得到内环边径向速度 ${\dot{u}}_{a}$ 随内环边挠度w_a变化曲线，如图9所示，得到最大内环边径向速度 ${\dot{u}}_{a}^{\rm m}$ 为25.99 m/s。

图 8 冲量面密度载荷分布曲线

Figure 8. Distribution of the impulse per area

下载: 全尺寸图片幻灯片

图 9 环板内环边径向速度与挠度的关系

Figure 9. Relationship between inner edge radial velocity and deflection

下载: 全尺寸图片幻灯片

将计算得到的环板变形所达到的最大径向速度代入圆环的Grady破碎模型，计算环板的撕裂情况。将速度代入式(5)，得到内环边最大环向应变率 ${\dot{\varepsilon }}_{\theta a}^{\rm m}$ 为1.13×10³ s⁻¹。当内环边径向速度最大时，即膨胀环达到最大膨胀速度时，材料更易破碎，将此时的最大环向应变率 ${\dot{\varepsilon }}_{\theta a}^{\rm m}$ 代入Grady膨胀环破碎模型，求得撕裂花瓣数n

$n={{\left( \frac{4{{\text{π} }^{3}}{{\rho }^{2}}c_{0}^{2}{{a}^{3}}\dot{\varepsilon }_{\theta a}^{\rm m}}{3K_{\rm c}^{2}} \right)}^{1/3}}$

(35)

式中：K_c为材料断裂韧性，c₀为材料声速， $\,\rho$ 为薄板材料密度，a为环板内环半径。由材料与实验工况参数：K_c = 3 × 10⁷N·m^−3/2，c₀ = 5328 m/s，a = 23 mm， $\,\rho$ = 2.78 g/cm³，计算得到n = 5.95，即实验中第4层板花瓣数为6，实验结果与理论计算结果比较符合。将破碎时的环板变形状态作为变形终态，由式(34)可以求出环板凹陷变形终态角度 $\theta ^{\rm{t}}$ = 14.36°，与实验测得的12°的相对误差不超过20%，说明本研究的理论模型是有效的。

上述模型可以用于描述碎片云载荷分布较分散的环形薄板变形，其环板凹陷变形区域集中在碎片云载荷分布区域内，由于载荷不太集中，环板的变形不太快，变形区与未变形区分界的塑性铰采用驻定塑性铰，其计算结果与第4层板的实验结果相符合。实验中二次碎片云撞击第3层板时的冲量面密度相对较高，第3层板凹陷变形的速度较快，拉着未变形区域发生凹陷，其变形区与未变形区的塑性铰是移动的，因此第3层板的变形需采用移行的塑性铰圆来描述。对于移行铰，l = l(t)， $\dot{l}\ne 0$ ，由图6可知，这时内环边挠度w_a和径向位移u_a为

$\begin{cases} {w_a} = (l - a)\sin \theta \\ {u_a} = (l - a)(1 - \cos \theta ) \end{cases}$

(36)

求导得到内环速度，并取二阶近似

$\begin{cases} {{\dot w}_a} = \dot l\theta + (l - a)\left( {1 - \dfrac{{{\theta ^2}}}{2}} \right)\dot \theta \\ {{\dot u}_a} = \dfrac{{\dot l{\theta ^2}}}{2} + (l - a)\theta \dot \theta \end{cases}$

(37)

由式(16)得到驻定铰工况内环边径向速度 ${\dot{u}}_{a}$ 近似为 $\left(l-a\right)\theta \dot{\theta }$ ，比式(37)计算得到的移行铰工况计算结果小，二者之差为弯曲转角 $\theta$ 的二阶量： ${\text{Δ}}{\dot{u}}_{a} \approx O\left({\theta }^{2}\right)$ 。薄板移行铰变形过程中会发生靶板撕裂，撕裂后各块板的速度依然较大，其塑性铰圆仍在移动中，整个薄板的变形建模和求解更为复杂。

4. 结　论

通过超高速撞击实验总结了薄板在轴对称碎片云载荷作用下的变形与破坏模式，并对此完成了建模研究，结论如下：

(1) 实验中薄板在分散的碎片云载荷下经历中部穿孔，形成环板，环板在剩余碎片云载荷下发生凹陷，孔边材料继而发生拉伸破碎，形成花瓣型撕裂；

(2) 基于实验现象与经典理论，建立了环形薄板在碎片云冲击下的刚塑性变形模型，该模型同时考虑弯矩和膜力，通过简化得到了环板变形微分方程，可以求解得到板横向速度场和径向速度场；

(3)将孔周边的薄板材料膨胀视为圆环的膨胀，根据模型计算得到环板内环最大环向应变率，代入Grady破碎模型，计算实验中第4层板撕裂的花瓣数，另外根据模型还计算得到环板凹陷变形终态角度，计算结果均与实验结果吻合；

(4) 对于高幅值碎片载荷的环形薄板变形，宜采用移行塑性铰描述变形区域的扩展过程，其变形的理论建模和求解还需要在本研究的基础上进一步开展。

图 1 弹丸超高速碰撞多层板结构

Figure 1. Multi-shock shield structure subjected to hypervelocity impact

下载: 全尺寸图片幻灯片

图 2 高速冲击后各层板的破坏照片

Figure 2. Damage results of each layer after hypervelocity impact

下载: 全尺寸图片幻灯片

图 3 环板变形侧视图

Figure 3. Lateral schematic of annular plate deformation

下载: 全尺寸图片幻灯片

图 4 Tresca屈服条件与正交法则

Figure 4. Tresca condition and flow rule

下载: 全尺寸图片幻灯片

图 5 环板变形俯视图

Figure 5. Top schematic of annular plate deformation

下载: 全尺寸图片幻灯片

图 6 挠度与径向关系示意图

Figure 6. Schematic of relationship between deflection and radial displacement

下载: 全尺寸图片幻灯片

图 7 冲量载荷下环板变形示意图

Figure 7. Deflection of annular plate subjected to impulse loading

下载: 全尺寸图片幻灯片

图 8 冲量面密度载荷分布曲线

Figure 8. Distribution of the impulse per area

下载: 全尺寸图片幻灯片

图 9 环板内环边径向速度与挠度的关系

Figure 9. Relationship between inner edge radial velocity and deflection

下载: 全尺寸图片幻灯片

参考文献(17)

[1]	BALOS S, HOWARD D, BREZULIANU A, et al. Perforated plate for ballistic protection: a review [J]. Metals, 2021, 11(4): 526. doi: 10.3390/met11040526
[2]	LU Y Y, ZHANG Q M, XUE Y J, et al. High-velocity impact performance of aluminum and B₄C/UHMW-PE composite plate for multi-wall shielding [J]. Applied Sciences, 2020, 10(2): 721. doi: 10.3390/app10020721
[3]	JONES N. Structural impact [M]. 2nd ed. Cambridge: Cambridge University Press, 2011: 1−10.
[4]	CRAGGS J W. The normal penetration of a thin elastic-plastic plate by a right circular cone [J]. Proceedings of the Royal Society of Edinburgh, 1952, 63(4): 359–370.
[5]	钱伟长. 穿甲力学[M]. 北京: 国防工业出版社, 1984: 177−215. QIAN W C. Perforation mechanics [M]. Beijing: National Defense Industry Press, 1984: 177−215.
[6]	REISSNER E. On finite deflections of circular plates [C]//REISSNER E, PRAGER W, STOKER J J.Proceedings of a symposia in applied mathematics (Volume Ⅰ): nonlinear problems in mechanics of continua. New York: American Mathematical Society, 1949: 213−219.
[7]	JONES N. Impulsive loading of a simply supported circular rigid plastic plate [J]. Journal of Applied Mechanics, 1968, 35(1): 59–65. doi: 10.1115/1.3601174
[8]	JONES N. A theoretical study of the dynamic plastic behavior of beams and plates with finite-deflections [J]. International Journal of Solids and Structures, 1971, 7(8): 1007–1029. doi: 10.1016/0020-7683(71)90078-3
[9]	YU T X, CHEN F L. Analysis of the large deflection dynamic plastic response of simply-supported circular plates by the “membrane factor method” [J]. Acta Mechanica Sinica, 1990, 6(4): 333–342. doi: 10.1007/BF02486892
[10]	龙仁荣. 超高速碰撞多层板结构的碎片云运动破坏模型 [D]. 北京: 北京理工大学, 2008. LONG R R. Description and damage effects of debris cloud formed by impact of hypervelocity projectile on multi-plate structure [D]. Beijing: Beijing Institute of Technology, 2008.
[11]	陈至达. 板、壳有限变形分析 [J]. 力学进展, 1983, 13(2): 1–2. CHEN Z D. Nonlinear analysis of large deformation of plates and shells [J]. Advances Mechanics, 1983, 13(2): 1–2.
[12]	TIMOSHENKO S, WOINOWSKY-KRIEGER S. Theory of plates and shells [M]. New York: McGraw-Hill, 1959: 100−112.
[13]	HODGE P G JR. Yield conditions for rotationally symmetric shells under axisymmetric loading [J]. Journal of Applied Mechanics, 1960, 27(2): 323–331. doi: 10.1115/1.3643960
[14]	AGGARWAL H R, ABLOW C M. Plastic bending of an annular plate by uniform impulse [J]. International Journal of Non-Linear Mechanics, 1971, 6(1): 69–80. doi: 10.1016/0020-7462(71)90035-7
[15]	欧阳成生. 薄膜力作用下钢筋混凝土环板的塑性分析 [J]. 土木工程学报, 1983, 3: 59–70. OUYANG C S. A plastic analysis of R. C. ring slabs under membrane action [J]. China Civil Engineering Journal, 1983, 3: 59–70.
[16]	ZUKAS J A. High velocity impact dynamics [M]. New York: Wiley-Interscience, 1990: 65−126.
[17]	GRADY D E. Local inertial effects in dynamic fragmentation [J]. Journal of Applied Physics, 1982, 53(1): 322–325. doi: 10.1063/1.329934