考虑靶体自由表面和开裂区影响的可变形弹体斜侵彻脆性材料的终点弹道分析

孔祥振，方秦，吴昊

(解放军理工大学 国防工程学院，江苏 南京210007)

0 引言

动能弹斜侵彻问题是防护工程领域和武器研发领域的研究重点之一。弹体斜侵彻过程中，由于初始不对称力和自由表面的影响，弹体产生弯曲，且弹道发生偏转。对该问题的主要研究方法可分为经验公式法、理论分析方法和数值模拟法。经验公式和理论分析法可以较好地预测弹体的终点弹道参数，但若要研究斜侵彻过程中弹体的结构响应，则须采用数值模拟方法。传统的数值模拟方法基于守恒定律和控制方程，在离散的空间网格和时间点上求出所有控制方程的解，根据坐标系选取的不同可以分为Lagrangian 法、Euler 法和ALE 法，而这些方法都必须基于对靶体材料的准确描述，且需要考虑复杂的接触算法，计算成本较高，此外，Lagrangian 方法可能由于大变形产生网格畸变而使计算终止［1］。

基于球形空腔膨胀理论，Warren 等［2-3］将靶体对弹体的作用用阻力函数代替，采用弹靶分离的方法，使用Sandia 国家实验室编制的PRONTO 3D［4］，分别对可变形钢弹正侵彻铝合金靶体和石灰岩靶体进行了数值模拟，模拟中同时考虑了微小偏航角的影响，该方法不需要对靶体进行网格划分，避免了复杂的接触算法，从而节约了计算时间和成本，计算结果表明，该方法可以较好地预测弹体的运动及变形情况。何涛等［5-6］基于同样的思想，采用ABAQUS/Explicit，结合自主开发的子程序，分别对球形弹和卵形弹正侵彻铝合金靶体进行了数值模拟，模拟得到的侵彻深度和弹体变形情况与实验结果吻合较好。

上述方法能否准确模拟弹体的运动变形情况取决于阻力函数的准确性。而空腔膨胀理论得到的阻力函数基于半无限大靶体假设，当弹体斜侵彻时，由于自由表面的影响，需要对偏向自由表面一侧的阻力函数进行修正。Jung 等［7］认为弹体表面到靶体自由表面的距离存在一个阀值，若距离小于阀值，则弹体表面点不受力，并经验地确定了冻土的该阀值。然而，这忽略了自由表面效应随着空腔膨胀速度的增大而增大［8-9］。将靶体视为不可压缩材料，Macek等［10］研究了有限空腔膨胀理论，确定了考虑自由表面效应的阻力函数，但将靶体视为不可压缩材料，会高估阻力函数［11-12］。Warren 等［8-9］将靶体视为不可压缩材料，假定空腔膨胀产生塑性-弹性响应分区，通过比较不可压缩条件下得到的有限空腔表面应力与无限空腔表面应力，得到了自由表面效应衰减函数，然后将其乘以可压缩模型下得到的阻力函数，将可压缩性引入到考虑自由表面效应的阻力函数中，并分别对4340 高强钢弹体斜侵彻铝合金靶体和石灰岩靶体进行了数值模拟。对于铝合金靶体，模拟结果与实验结果吻合较好;而对于石灰岩靶体，入射速度较高时模拟结果与实验结果吻合较好，而速度较低时，模拟得到的弹体转动与实验结果差别较大。何涛等［13］将上述衰减函数推广应用于弹体贯穿金属靶体的数值模拟中(同时考虑两个自由表面的影响)，并对4340 高强钢弹体斜贯穿铝合金靶体进行了数值模拟，模拟得到的弹体最终形态与实验结果吻合较好。

本文基于Warren 等［9］的工作，将靶体视为不压缩的Mohr-Coulomb 材料，得到了基于塑性-开裂-弹性响应分区的自由表面效应衰减函数，并将脆性材料半经验阻力函数乘以该衰减函数进行修正，然后采用弹靶分离的方法将修正后的阻力函数作为边界条件施加于弹体表面，对4340(RC44.5)高强钢弹体斜侵彻石灰岩靶体进行了数值模拟，并分别将数值模拟结果与实验测得结果及Warren 等［9］模拟结果进行了对比，验证了本文提出方法的正确性和优越性。

1 考虑自由表面效应的脆性材料半经验阻力函数

本节将靶体视为不可压缩的Mohr-Coulomb 材料，假定空腔膨胀产生塑性-开裂-弹性响应分区，构造了自由表面效应衰减函数，并将其乘以脆性材料半经验阻力函数，得到了考虑自由表面效应的脆性材料半经验阻力函数。此外，由于弹靶之间摩擦机理尚不明确，而且弹靶之间摩擦系数是速度和压力相关的，难以测量和确定，并且弹靶间的摩擦效应已经拟合到了动态抗压强度R 中，因此，本文计算中忽略了弹靶间摩擦的影响。

1.1 半经验阻力函数

诸多混凝土和岩石等脆性材料靶体的正侵彻实验表明［14-15］，靶体的破坏形态为一个近似倒锥形的冲击坑和一个直径近似等于弹体直径的钻孔区组成，如图1 所示，其中冲击坑深度为Hc，弹体侵彻深度为Hp. 在钻孔区，Frew 等［14］、Forrestal 等［15］基于动态球形空腔膨胀理论和对实验的回归分析，提出了半经验的脆性材料靶体阻力函数公式，并被以后较多学者所采用，即

式中:σn为靶体阻力函数;ρ0为靶体初始密度;R 为靶体的动态抗压强度;vn为瞬时弹体表面法向速度。对于混凝土靶体，R=Sfc，其中S 为无量纲靶体经验参数，包含了高温、高压、高应变率、高静水压力等影响，Frew 等［16］通过对正侵彻实验数据的回归得到S=82.6 ×fc-0.544，式中fc单位为MPa. 对于石灰岩靶体，Frew 等［14］观察到R 随弹径的增大而减小，因此建议R =Φ +φ(D0/D)，其中D 为弹体直径，Φ、φ 和D0为经验系数，通过拟合实验数据，Frew等［14］建议Φ=607 MPa，φ=86 MPa，D0=25.4 mm.

图1 侵彻两阶段模型Fig.1 Two-stage penetration model

在开坑区，由于应力波在靶体自由表面的反射引起靶体破坏，动态球形空腔膨胀理论不再适用，Frew 等［14］、Forrestal 等［15］认为弹体从自由表面到开坑区和钻孔区交界处阻力呈线性增长，并认为开坑区深度为2 倍弹体直径。在数值方法中，将开坑区等分为10 层，则每一层的靶体阻力函数［3］为

式中:σn是按(1)式求得的靶体阻力函数。在开坑区，由于靶体破坏产生的材料飞散，弹身位置并不受力，靶体阻力函数只作用于弹头上。

1.2 考虑自由表面效应的衰减函数

如图2 所示，假设空腔在脆性材料中膨胀产生塑性-开裂-弹性响应分区，各响应分区的分界面分别为rb、rc，空腔半径为ra，空腔膨胀速度为，rd处为自由表面。将靶体视为不可压缩的Mohr-Coulomb 材料，以下推导构造塑性-开裂-弹性响应分区下的自由表面效应衰减函数。

图2 空腔膨胀响应区Fig.2 Response regions for cavity expansion

球对称欧拉坐标系下动量和质量守恒方程分别为

式中:σr和σθ分别为径向和环向应力，受压为正，受拉为负;ρ0和ρ 分别为靶体初始状态和最终状态的密度，当靶体不可压缩时有ρ=ρ0;粒子径向位移s 和粒子径向速度v 有以下关系:

对(4)式积分可以得粒子径向位移为

将(6)式代入(5)式可得粒子径向速度为

1.2.1 弹性区的求解

由于弹性区粒子位移很小，粒子径向位移［17］可近似为

则弹性应变为

对于不可压缩材料，由胡克定律:

将(10)式代入动量守恒方程(3)式，得

式中:Λ1为积分常数，可通过σr在自由表面处(r =rd)为0 求得。

1.2.2 开裂区的求解

当弹性区环向应力σθ到达材料最大抗拉强度时，靶体由弹性区进入开裂区，则开裂-弹性交界面上满足:

式中:ft为靶体单轴抗拉强度。由Hugoniot 跳跃条件，开裂-弹性交界面左右径向应力连续，则

利用牛顿迭代法解(14)式即可得开裂-弹性交界面位置rc.

在开裂区，σθ=0，将其代入动量守恒式(3)式，则

式中:Λ2为积分常数，可由(16)式与(12)式在r =rc处相等求得。

1.2.3 塑性区的求解

脆性材料塑性区由Mohr-Coulomb 破坏准则描述:

式中:λ、τ 为Mohr-Coulomb 常数。

塑性-开裂交界面左右径向应力连续，则有

式中:Y 为靶体单轴抗压强度。利用牛顿迭代法解(18)式即可得塑性-开裂界面位置rb.

将(17)式代入(3)式得

式中:Λ3为积分常数，可由(20)式与(16)式在r=rb处相等得到。

当空腔膨胀速度大于一定值时，开裂区消失(rb＞rc)，响应分区变为弹塑性分区，此时，塑性-弹性交界面位置可由下式［9］确定:

1.2.4 空腔表面径向应力的求解

由上述讨论，空腔表面径向应力可通过以下步骤进行求解:

1)由(14)式和(18)式分别求得rc、rb. 若rc＞rb，则空腔膨胀产生塑性-开裂-弹性响应分区;若rc＜rb，则空腔膨胀产生塑性-弹性响应分区。

2)若rc＞rb，依次利用(12)式、(16)式、(20)式即可求得空腔表面应力。但若rc＞rd时，即不存在弹性区时，脆性材料受拉破坏，此时，令空腔表面径向应力σr(ra)=0.

3)若rc≤rb，依次利用(12)式、(20)式即可求得空腔表面应力。但若rb＞rd时，即不存在弹性区时，脆性材料受压破坏，此时令空腔表面径向应力σr(ra)=0.

1.2.5 衰减函数的构造

若空腔在半无限靶体中膨胀，则产生塑性-开裂-弹性响应分区，此时空腔表面径向应力(ra)可由(12)式，(16)式和(20)式求得(令rd→∞)。则衰减函数构造为

图3 分别给出了不同空腔膨胀速度下，rd/ra与基于塑性-开裂-弹性响应分区得到的衰减函数和基于弹塑性响应分区得到的衰减函数关系曲线。其中，靶体为Salem 石灰岩，经验参数如表1 所示。由图3 可见:1)当空腔膨胀速度较低时，基于塑性-开裂-弹性响应分区得到的衰减函数比基于弹塑性响应分区得到的衰减函数小，即采用塑性-开裂-弹性响应分区时，自由表面效应影响较大;而当空腔膨胀速度较高时，由于开裂区消失，采用两种响应分区得到的衰减函数相同。2)随着空腔膨胀速度的增大，自由表面效应影响距离也越来越大。

图3 衰减函数与自由表面距离的关系曲线Fig.3 Decay function vs. distance to the free surface

由上述分析，考虑自由表面效应的半经验阻力函数可确定为

表1 石灰岩靶体经验参数［9］Tab.1 Empirical parameters of limestone target［9］

2 有限元模型

由第1 节分析，侵彻过程中靶体对弹体的作用可以用阻力函数σn来代替，这样，在进行有限元模拟时，就可以将靶体的响应作为边界条件施加于弹体表面，从而避免了靶体网格的划分和复杂的接触算法。图4 给出了弹体表面单元所受压力边界条件示意图，各单元所受压力可以表示为

式中:v 为结点速度矢量;n 为弹体表面在结点处的外法线方向。由上式可见，随着侵彻过程中v 和n不断变化，每个单元的瞬时压力也在不断变化，并且只有当弹体表面单元位于靶体内部时，压力边界条件才起作用。当弹体斜侵彻靶体时，还需要考虑弹靶之间的分离效应，本文采用的方法是当(v·n)＜0 时，令压力为0. 通过ABAQUS 用户子程序［18］将上述方法植入，在每一个时间步子程序都会被调用一次，用来计算施加于弹体表面的压力。

图4 弹体表面压力定义Fig.4 Definition of pressure boundary condition

3 数值模拟与实验结果对比

为验证本文提出方法的正确性和优越性，本节首先对4340 高强钢(RC44.5)弹体正侵彻石灰岩靶体进行了数值模拟，并对比了模拟得到的侵彻深度和实验测得值;然后对4340 高强钢弹体斜侵彻石灰岩靶体进行了数值模拟，并分别对比分析了模拟得到的弹体最终形态与弹尖最终位置和实验测得值及Warren 等［9］的模拟结果。

本文对4340 高强钢弹体采用幂次硬化本构模型［1］，即

3.1 正侵彻条件下侵彻深度

Frew 等［14］分 别 进 行 了7.1 mm、12.7 mm、25.4 mm不同直径4340 高强钢(RC44.5)弹体正侵彻石灰岩靶体的实验，图5 分别给出了数值模拟得到的侵彻深度和实验测得值。由图5 可见，数值模拟结果与实验吻合非常好，验证了本文的计算方法及其子程序的正确性。

图5 正侵彻侵彻深度模拟值与实验值［14］对比Fig.5 Comparison of the depths of numerical simulation and experimental［14］penetrations

3.2 斜侵彻条件下弹体变形

Warren 等［9］进行了4340 钢(RC44.5)弹体以不同倾角和初速斜侵彻石灰岩靶体的实验，表2 分别给出了不同工况下实验得到的弹体最终形态和依本文方法模拟得到的弹体最终形态及Warren 等［9］模拟得到的弹体最终形态图。由表2 可见，当倾角较小且初速较大时，本文模拟得到的结果和Warren等［9］模拟得到的结果与实验吻合均较好;而当倾角较大或初速较小时，本文得到的模拟结果较Warren等［9］模拟得到的结果与实验更为接近，验证了本文提出基于脆性材料塑性-开裂-弹性响应分区得到的自由表面效应模型的优越性。

表2 不同工况下数值模拟得到的弹体变形情况与实验结果的对比Tab.2 Comparison of the projectile deformations in numerical simulation，experiment and Ref.［9］

?

图6 分别给出了不同工况下实验得到的弹尖最终位置和依本文方法模拟得到的弹尖最终位置及Warren 等［9］模拟得到的弹尖最终位置。由图6 可见，依本文方法得到的弹尖最终位置较Warren 等［9］模拟得到的弹尖最终位置与实验更为接近，进一步验证了本文提出自由表面效应模型的优越性。

图6 不同工况下数值模拟得到的弹尖最终位置与实验结果的对比Fig.6 Comparison of the positions of the tips in numerical simulation，experiment and Ref. ［9］

4 结论

本文将脆性靶体视为不可压缩Mohr-Coulomb材料，假定空腔膨胀产生塑性-开裂-弹性响应分区，构造了自由表面效应的衰减函数，并将衰减函数乘以基于半无限大靶体假设的半经验阻力函数，得到了用于斜侵彻条件下半经验阻力函数。

基于弹靶分离的方法，将靶体对弹体的作用用阻力函数代替，避免了靶体网格的划分和复杂的接触算法，对可变形弹斜侵彻石灰岩靶体进行了数值模拟。数值模拟结果表明，当倾角较小且初速较大时，本文模拟得到的弹体最终形态和弹尖最终位置和Warren 等［9］模拟得到的结果与实验结果吻合均较好;而当倾角较大或初速较小时，本文得到的弹体最终形态和弹尖最终位置较Warren 等［9］模拟得到的结果与实验结果更为接近，验证了本文提出的基于脆性材料塑性-开裂-弹性响应分区的自由表面效应模型的优越性。

本文方法基于弹体无明显磨蚀假设，而当弹体高超声速侵彻混凝土或岩石靶体时(着靶速度达到1 500 m/s)，高温高压以及骨料的摩擦和切削等作用会产生弹体头部磨蚀效应。因此进一步的研究应在于将弹体磨蚀嵌入本文方法中，从而实现对弹体磨蚀和运动的耦合分析。

References)

［1］Camacho G T，Ortiz M. Adaptive Lagrangian modeling of ballistic penetration of metallic targets［J］. Computer Methods in Applied Mechanics and Engineering，1997，142:269 -301.

［2］Warren T L，Tabbara M R. Simulations of the penetration of 6061-T65ll aluminum targets by spherical-nosed VAR4340 steel projectiles［J］. International Journal of Solids and Structures，2000，37(44):4419 -4435.

［3］Warren T L. Simulations of the penetration of limestone targets by ogive-nose 4340 steel projectiles［J］. International Journal of Impact Engineering，2002，27:475 -496.

［4］Taylor L M，Flanagan D P. PRONTO 3D:a three-dimensional transient solid dynamics program［R］. Albuquerque，NM:Sandia National Labs，1989.

［5］何涛，文鹤鸣. 球形弹对金属靶板侵彻问题的数值模拟［J］.爆炸与冲击，2006，26(5):456 -461.HE Tao，Wen He-ming. Computer simulations of the penetration of metal targets by spherical nosed projectiles［J］. Explosion and Shock Waves，2006，26(5):456 -461. (in Chinese)

［6］何涛，文鹤鸣. 卵形钢弹对铝合金靶板侵彻问题的数值模拟［J］. 高压物理学报，2006，20(4):408 -414.HE Tao，WEN He-ming. Numerical simulations of the penetration of aluminum targets by ogive-nosed steel projectiles［J］. Chinese Journal of High Pressure Physics，2006，20(4):408 -414. (in Chinese)

［7］Jung J，Longcope D B，Tabbara M R. Modeling of oblique penetration into geologic targets using cavity expansion penetrator loading with target free-surface effects［R］. Albuquerque，NM:Sandia National Labs，1999.

［8］Warren T L，Poormon K L. Penetration of 6061-T6511 aluminum targets by ogive-nosed VAR 4340 steel projectiles at oblique angles:experiments and simulations ［J］. International Journal of Impact Engineering，2001;25:993 -1022.

［9］Warren T L，Hanchak S J，Poormon Kevin L. Penetration of limestone targets by ogive-nosed VAR 4340 steel projectiles at oblique angles:experiments and simulations［J］. International Journal of Impact Engineering，2004;30:1307 -1331.

［10］Macek R W，Duffey T A. Finite cavity expansion method for near-surface effects and layering during earth penetration ［J］.International Journal of Impact Engineering，2000，24 (3):239 -258.

［11］Forrestal M J，Luk V K. Dynamic spherical cavity-expansion in a compressible elastic-plastic solid ［J］. Journal of Applied Mechanics，1988，55(2):275 -279.

［12］Forrestal M J，Tzou D Y. A spherical cavity-expansion penetration model for concrete targets［J］. International Journal of Solids and Structures，1997，34(31):4127 -4146.

［13］何涛，文鹤鸣. 可变形弹丸贯穿铝合金靶的数值模拟［J］.高压物理学报，2008，22(2):153 -159.HE Tao，WEN He-ming. Numerical simulations of the perforation of aluminum targets struck by a deformable projectile［J］. Chinese Journal of High Pressure Physics，2008，22(2):153 -159. (in Chinese)

［14］Frew D J，Forrestal M J，Hanchak S J. Penetration experiments with limestone targets and ogive-nose steel projectiles［J］. Journal of Applied Mechanics，2000，67(4):841 -845.

［15］Forrestal M J，Altman B S，Cargile J D，et al. An empirical equation for penetration depth of ogive-nose projectiles into concrete targets［J］. International Journal of Impact Engineering，1994，15(4):395 -405.

［16］Frew D J，Hanchak S J，Green M L. Penetration of concrete targets with ogive-nose steel rods［J］. International Journal of Impact Engineering，1998，21(6):489 -497.

［17］Hopkins H G. Dynamic expansion of spherical cavities in metals［J］. Progress in Solid Mechanics，1960，1(3):5 -16.

［18］SIMULIA. ABAQUS，analysis user’s manual，Version 6. 5［M］. France:SIMULIA，2004.