巡航导弹，是在打击目标过程中，大部分时间都以巡航状态飞行的有翼制导弹药，是一种成本低、通用性强、作战使用灵活的武器系统[1]。只有将敌方巡航导弹拦截于其到达目标并对目标造成摧毁之前，才能奠定战争或战斗胜利的基础[2]。
在对导弹目标的毁伤及易损性研究方面，窦丽华等[3]提出了一种单发高炮弹药对巡航导弹目标毁伤概率的数学计算方式。常二莉等[4]提出了巡航导弹目标在高炮弹药打击下的毁伤模型；陈健等[5]分析了导弹制导系统及其部分构件的易损性；姜颖资[6]等利用数值模拟方法，研究了钨穿甲弹对反舰导弹发动机舱的毁伤规律；卢明章等[7]通过有限元仿真得到弹丸对导弹燃料舱毁伤的过程。在实验研究方面，曹兵[8]研究了破片在不同条件下撞击模拟油箱毁伤情况；李世纪等[9]分析了杆式射流对屏蔽PBX炸药的起爆能力。上述研究表明，一般而言对导弹毁伤和易损性问题的研究比较复杂，尤其是使用有限元模拟或试验方法时，时间成本或试验成本较高。为此，本文基于MATLAB软件，通过将破片对导弹毁伤问题转化为数学模型，完成对目标模型、破片射击线、弹目交会等问题的数值仿真，通过建立破片对目标的穿透毁伤判据，得到一种简化的能够计算破片对导弹目标舱段毁伤概率的方法，能够节省时间成本，可以较快速地对导弹舱段易损性作出分析，可为破片战斗部对导弹目标的毁伤评估提供一定的参考。
1 毁伤概率计算方法
研究单枚破片对目标舱段的毁伤概率主要使用蒙特卡洛方法，即同一破片在一定范围内的随机位置以随机角度打击目标，以对目标的毁伤频率去代替毁伤概率，进一步地完成对目标易损性的分析。
1.1 目标模型仿真
使用数值仿真方法研究破片打击巡航导弹毁伤概率问题需要借助MATLAB软件。在利用数值仿真方法进行毁伤概率研究之前，需要对目标进行处理，即让目标转化为由“面元”组成。这里需要将“面元”取为三角形面元，因为三角形自身在一个平面上，在计算过程中更为容易。通过有限元处理软件hypermesh将目标划分为由三角形网格组成的有限元模型，可以得到目标节点数据以及三角形网格节点信息。继而，可以使用MATLAB仿真程序对目标节点数据进行处理，得到由三角形面元组成的目标仿真模型。得到的巡航导弹制导舱节点模型和三角形面元仿真模型如图1所示，巡航导弹其他舱段面元模型同样可以依照此方法获得。
图1 巡航导弹制导舱节点模型及三角形面元仿真模型示意图
1.2 破片毁伤元及其射击线
影响毁伤概率的破片参数主要有破片材料、质量、初始速度、形状等等。在破片材料和形状选定的前提下，可研究破片初始速度和破片质量对毁伤概率的影响和变化规律。
对破片打击过程的仿真计算需要借助破片射击线，破片射击线是包括破片质量、破片初始速度、破片初始位置和破片打击方向的用来模拟破片弹道轨迹的模型。若破片初始位置为(x0,y0,z0)，破片初始速度v0为(vx,vy,vz)，三个方向分速度vx、vy、vz的合速度方向即为v0方向，则该破片空间射击线方程可以表示为：
(1)
本文研究的目标易损性主要是破片毁伤元命中后目标的毁伤概率，因此可以简化射击线模型，将破片简化为一个点，且假定在破片射击线上，破片打击速度保持为初速度为v0。为了得到某质量破片对目标的毁伤概率，应该使得破片对目标的打击方向随机，应用蒙特卡洛方法，通过大量随机试验，进而得到毁伤概率，因此，需要数值仿真破片随机射击线。
由于使用蒙特卡洛方法求单个破片对目标的毁伤概率的基本思想是用频率代替概率，即命中并造成毁伤的破片数目与命中总数目的比值来代替毁伤概率，因此，为了使得到的所有数据都有用，不妨使得每次破片都可以命中目标。同时，为了最大可能保证打击角度的随机性，可以使得每次运算破片初始位置也在一定范围内随机。因此，通过这种方法可以得到运算10次的随机试验射击线结果示意图如图2。
图2 随机产生的10条射击线示意图
从图2中可以看出，这种方法产生的随机破片射击线，每次运算破片都可以命中目标，能够保证计算过程的高效性和有效性。
1.3 弹目交会计算
假定破片初始位置为O，方向为D(D的模为1)，则该破片射击线可以表示为方程：
R(t)=O+tD
(2)
t表示射击线上一点到破片初始位置O的距离，t≥0。根据t值的不同，D方向上所有点构成过破片初始位置O点沿着方向D的破片射击线。
构成目标三角形面元的三个顶点由P1、P2、P3定义。在弹目交会问题当中，需要确定破片是否命中目标，即需要确定射击线是否与目标某个三角形面元相交。求解问题示意图如图3。
图3 破片命中目标三角面元问题示意图
任意一个三角形面元上的点可以由以下方程(3)表示。
T(u,v)=(1-u-v)P1+uP2+vP3
(3)
其中，u、v、1-u-v分别表示三个顶点P2、P3和P1的“权重”，并且同时满足：
计算破片射击线与三角形面元交点，即求解：
R(t)=T(u,v)
(4)
其中，参数t、u、v为未知参量。通过求解可得，
(5)
其中：E1=P2-P1,E2=P3-P1,T=O-P1,P=D×E2,Q=T×E1，这样，如果方程有实数解，即说明破片射击线与三角形面元相交，说明破片命中目标，可求解出破片射击线与目标三角形面元交点。但是，当求解破片射击线与多个三角形面元构成目标的交点时，很有可能求解出两组解，因为，在破片射击线方向上，如果射击线与某个面元相交，则在几何模型上该射击线可能与目标另一侧某个三角形面元也相交，即该射击线“贯穿”了目标。而在实际当中破片贯穿巡航导弹目标的可能性往往比较小，因此，还要找出该问题的唯一解，由于解得的t是破片初始位置到命中点的距离，因此，该解可取t取值最小的一组解为唯一解。图4为破片命中目标示意图，其中与射击线相交的三角形面元为命中单元，单元上标记的“*”为命中点。
图4 破片命中目标示意图
1.4 穿透毁伤判据
根据修正后的THOR侵彻方程[10-11]，球形破片侵彻靶板极限速度公式为：
(6)
式中: Vj为破片穿透靶板临界速度(m/s)；h为靶板厚度(m)；d为球形破片直径(m)；ρt为靶板材料密度(kg/m3)；ρp为破片材料密度(kg/m3)；θ为破片入射角；σt为靶板材料强度极限(MPa)；a、b为与弹靶材料有关的常数，通常由弹靶侵彻试验确定。
由于巡航导弹目标材料组成复杂，因此在对其进行毁伤研究时，常对目标进行简化处理，通过对各个舱段建立等效模型，将各个舱段等效为具有一定厚度的LY-12硬铝材料模型。当破片为球形钨破片，靶板为铝靶板时，a=4.24，b=0.75[10]。又由于钨密度为ρp=17 600 kg/m3，LY-12硬铝密度为ρt=2 730 kg/m3，强度极限为460 MPa，故球形钨破片打击等效为LY-12硬铝材料目标的极限速度公式为：
(7)
当破片打击具有一定厚度的某一舱段时，如果打击速度超过了弹道极限，即可以将目标穿透，可认为对目标舱段造成了毁伤。
2 破片作用下巡航导弹舱段的易损性分析
由于使用蒙特卡洛方法对目标易损性的计算是用频率替代概率完成的，因而要考虑数值仿真运算次数对结果的影响；同时，目标的仿真计算模型由一个个三角形面元所构成，因此，构成目标模型的三角形面元数量同样会对计算结果精度产生影响。通过进一步试验性计算得到，运算次数取 1 000 次，三角形面元数量无需过多，即可以在保证运算精确度的基础上较快速地完成巡航导弹舱段在单枚破片作用下的易损性计算。
巡航导弹目标的毁伤建立在巡航导弹舱段毁伤的基础之上，在研究破片对巡航导弹毁伤概率之前，有必要对单枚破片对巡航导弹各个舱段的毁伤概率进行研究。
上一节中提到，在对巡航导弹各舱段毁伤进行研究时，可以将各舱段等效为一定厚度的LY-12硬铝材料模型，某巡航导弹整体上可以划分为5个舱段或结构，将某巡航导弹各个舱段进行等效之后，可以得到各个舱段等效厚度值如表1所示[12]。
表1 巡航导弹不同舱段等效LY-12硬铝材料厚度
巡航导弹舱段等效LY-12硬铝厚度/mm制导舱15.00战斗部20.00燃油舱7.00发动机舱14.60弹翼8.00
巡航导弹各个舱段数值仿真模型如图5所示。
图5 巡航导弹各舱段数值仿真模型示意图
根据导弹近炸引信制导规律[13-15]，一般对空导弹在距离目标5 m左右时，引信开始引爆战斗部，即对空导弹脱靶量一般为5 m左右，近似可取5 m，因此，根据蒙特卡洛方法，可取随机破片产生位置距离目标小于5 m，以随机方向打击目标舱段。
通过改变破片质量和速度对单枚破片打击各个舱段进行数值仿真计算，可以获得不同舱段毁伤概率随着破片打击速度的变化曲线，即目标舱段在破片毁伤作用下的毁伤概率随破片速度变化曲线(即易损性曲线)如图6所示。
图6 巡航导弹不同舱段在破片毁伤作用下的毁伤概率随破片速度变化曲线
从图6可以看出，对于某一质量的破片，打击巡航导弹某个舱段时，随着破片速度增大，对该舱段毁伤概率也增大；在某一打击速度下，质量大的破片对舱段毁伤概率要高于质量小的破片；对于不同舱段而言，一般情况下即使是同一质量破片在同一速度下进行打击，对舱段毁伤概率也不同。至于破片对各个舱段毁伤的难易程度，整体上看，破片的质量和速度较低时，就有较大概率能够对燃油舱造成毁伤，说明对燃油舱毁伤相对而言比较容易，相反，能对战斗部舱段造成较大毁伤概率时，破片质量和速度普遍较高，则说明破片对战斗部舱段的毁伤难度比较大。
通过对计算方法的加工处理，通过整合将算法集成为一个具有前台界面的能够对单个破片打击巡航导弹舱段毁伤概率进行计算的独立完整程序，程序计算界面如图7所示。
图7 毁伤概率计算程序计算界面
通过程序面板左上方下拉菜单可以选择不同舱段作为计算目标，通过将破片质量等初始参数输入后，可以对破片打击某一舱段的毁伤概率进行计算。并且，为了对程序计算过程和结果有更直观的理解，可以对一系列过程及计算结果图像显示和输出，程序也可以直接计算输出某个舱段在某一质量破片作用下的易损性曲线。
3 结论
1) 利用数学方法完成了巡航导弹目标模型仿真、破片毁伤元及其射击线仿真、弹目交会过程计算以及破片对目标造成穿透毁伤的判断方法等，得到一种巡航导弹舱段在破片毁伤作用下的易损性快速分析方法，能够节省时间成本，用于导弹舱段性研究。
2) 根据计算结果，作出了1 g、2 g、3 g、5 g、10 g破片打击下巡航导弹各个舱段目标的易损性曲线。
3) 对算法进行了整合处理，将算法开发为具有清楚计算界面的独立计算程序，使用该程序可计算不同条件下单枚破片对巡航导弹各个舱段的毁伤概率以及对目标舱段易损性快速分析。
[1] 陈胜政,高晓波.小型飞航导弹发展综述[J].兵器装备工程学报,2020,41(01):94-99.
[2] 刘桐林,何立萍,高雁翎.巡航导弹防御技术发展[J].飞航导弹,2000(06):36-43.
[3] 窦丽华,王高鹏,张娟.高炮对巡航导弹毁伤概率仿真[J].火力与指挥控制,2007,32(12):49-51,58.
[4] 常二莉,韩亮.高炮对巡航导弹拦阻射击毁伤计算模型[J].系统仿真学报,2006(S2):83-84.
[5] 陈健,米双山,易胜,等.巡航导弹制导系统毁伤模型研究[J].弹道学报,2004(03):44-48.
[6] 姜颖资,王伟力,黄雪峰,等.钨穿甲弹对超音速反舰导弹发动机舱的毁伤效应[J].科技导报,2014,32(13):19-22.
[7] 卢明章,杜茂华,王伟力,等.基于AYSYS/LS-DYNA的破片碰撞作用下导弹燃料舱毁伤机理仿真研究[J].海军航空工程学院学报,2011,26(02):185-188.
[8] 曹兵.破片对巡航导弹燃油舱冲击引燃的实验研究[J].火炸药学报,2008,31(06):45-49.
[9] 李世纪，魏锦.杆式射流冲击屏蔽PBX炸药数值模拟研究[J].兵器装备工程学报，2020,41(09):135-138,143.
[10] 黄长强,朱鹤松.球形破片对靶板极限穿透速度公式的建立[J].弹箭与制导学报,1993(02):58-61.
[11] THOR P.The Resistance of Various Metallic Materials to Perforation by Steel Fragments Empirical Relationships for Fragment Residual Velocity and Residual Weight.[R].Maryland:Ballistic Analysis Laboratory,Johns Hopkins Univ,1961.
[12] 任丹萍.破片和冲击波复合作用下对导弹的毁伤[D].南京：南京理工大学,2006.
[13] 王少峰.某型空空导弹精度分析的数字仿真技术研究[D].西安： 西北工业大学,2006.
[14] 张志鸿.防空导弹引信技术的发展[J].现代防御技术,2001,29(4):26-31.
[15] 周凤岐,李言俊.用蒙特卡罗法进行导弹脱靶量计算时的精度估计[J].西北工业大学学报,1988(03):271-280.