枪管内膛损伤对弹头外弹道过程的影响研究

2020-07-15沈超周克栋赫雷陆野李峻松

北京理工大学学报 2020年6期

沈超，周克栋，赫雷，陆野，李峻松

(1.南京理工大学机械工程学院，江苏,南京 210094；2.中国兵器工业第208研究所，北京 102202)

枪弹发射过程中，膛线与弹头壳之间紧密接触以减少火药气体的泄露，并且阳线会在弹头上形成若干道刻槽，具有一定缠度的阳线与刻槽之间的接触力，可以分解为沿枪管轴向的挤进阻力和沿周向的导转力，后者使得弹头在沿枪管轴向运动的同时产生高速自转，以保证弹头出膛后的飞行稳定性[1-2].线膛的结构参数(如阴、阳线宽度，阳线高度，缠度等)直接决定了弹头的表面形貌.随着射弹量的增加，枪管内膛磨损导致的阳线高度、宽度等的降低会使得弹头壳上的刻槽高度、宽度等随之降低.此外，内膛表面还分布着裂纹，烧蚀坑和铬层剥落等损伤形式，在弹头高速运动过程中会划伤弹头壳，破坏弹头壳的表面完整性[3].因此，不同寿命阶段枪管所发射弹头的表面形貌是有差异的，弹头的气动特性必然会随之发生改变.研究枪管内膛损伤对弹头气动特性的影响规律，是获得射击精度、横弹孔率等(枪管寿终与否的直接判据)随射弹量增加的变化规律的基础，对于加深对枪管寿终机理的认识有重要意义.

国内外学者们对旋转弹头或弹丸气动特性和外弹道过程的实验测试与数值计算方面进行了较为深入的研究，并取得了一定的成果.De Spirito等[4-5]使用雷诺平均和大涡模拟(RANS/LES)混合模型获得了稳定旋转的M910弹丸在亚声速、跨声速及超声速来流下的气动特性参数，并通过与实验结果的对比，得到了不同的数值模拟方法对弹丸各气动特性参数预测的准确度.Silton等[6]采用实验测试与数值模拟相结合的方法研究了膛线在小口径弹丸上留下的刻槽对弹丸气动特性的影响，发现刻槽对不同气动参数的影响程度也存在较大差异.Jacobson 等[7]采用有限差分法计算了旋成体在超声速流中的马格努斯效应，研究了各种边界层状态下马赫数、攻角和转速等变化对马格努斯效应的影响.孟鹏等[8]基于SSTk-ω湍流模型对3维Navier-Stokes方程进行求解，通过数值模拟研究了弹带对高速旋转弹丸气动特性的影响规律；陈亮等[9]结合数值仿真及近似理论方法分析了弹箭高速滚转运动对其气动参数的影响，并获得了不同马赫数下的平衡转速随尾翼结构参数的变化规律.赵清杰等[10]提出了一种新的大迎角非线性非定常气动力及气动力矩建模方法，降低了气动参数的预测结果与真实数据的误差.

对于旋转弹头或弹丸的气动特性，现有的试验分析及数值模拟大多研究的是无膛线刻槽、无损伤的弹头或弹丸的气动特性，也有部分学者考虑了膛线在弹头表面形成的刻槽对弹头气动特性的影响，但均未对不同寿命阶段的含内膛损伤枪管所发射的弹头的气动特性及外弹道过程受内膛损伤的影响规律进行深入分析，而是用初始状态下弹头的气动力参数去计算损伤枪管所发射弹头的外弹道过程参数，这种处理方法的准确度及可信度均较低.

为了更深入地研究枪管内膛损伤导致枪管寿终的机理，本文基于某12.7 mm机枪枪管寿命试验，研究对象为12.7 mm重机枪枪管在1个寿命周期中的4个不同的寿命阶段所发射的具有不同表面形貌的弹头，通过对4个阶段枪管的完整内弹道过程进行仿真计算，获得了各阶段枪管所发射弹头出膛时的表面形貌状态.基于SSTk-ω湍流模型对3维Navier-Stokes方程进行雷诺时间平均，采用滑移网格技术处理弹头的高速旋转所引起的动边界问题，对4种不同表面形貌的弹头的气动特性进行了数值模拟；采用均匀设计方法安排随机因素影响下的弹头内弹道计算过程，获得了各寿命阶段枪管发射弹头的膛口扰动随机响应状态；基于弹头的6自由度刚体外弹道方程对其外弹道过程进行了求解，结合枪管寿命试验结果对内膛损伤导致弹头椭圆弹孔率超标进而导致枪管寿终的过程进行了深入研究.

1 数值计算方法

本文以3维Navier-Stokes方程(简称N-S方程)为基础，采用剪应力输运湍流模型(shear stress transportk-ω,简称SSTk-ω)，针对枪管寿命试验中弹头初速变化较小，而运动姿态变化较大的情况[11]，对4个寿命阶段的枪管所发射的具有不同表面形貌的弹头的气动特性及外部流场分布情况进行数值模拟.由于膛线对弹头的导转力作用，弹头出膛时转速可以达到甚至超过10 000 rad/s，为准确模拟弹头旋转导致的非定常流对其气动特性的影响，本文采用滑移网格方法处理弹头高速旋转所引起的动边界问题.

采用滑移网格技术需将计算域分为两个区域，即包围弹头并随弹头以相同角速度旋转的内部旋转区和该内部旋转区之外的外部固定区(内部旋转区相对于外部固定区产生滑移运动)，两区域之间有一对交界面，交界面上的网格节点不需要相互重合，在两区域的滑移交界面上采用数值插值以保证区域间的通量守恒.这样处理既便于在流动核心区域建立更密的网格，又可以降低远离核心区域的网格数量，且内部旋转区网格在运动过程中不发生变形，求解精度高，计算速度快.

1.1 控制方程

采用滑移网格技术时，内部旋转区相对外部固定区做定常运动以模拟弹头的旋转，内部旋转区网格转速与弹头转速相等，外部固定区网格保持不动，且由于旋转弹头飞行速度可达数倍声速，空气的压缩性不可忽略，因此N-S方程应用于求解滑移网格流场域时的基本控制方程组可表示为

(1)

式中：V为控制体体积；ρ为空气密度；∂V为控制体的边界面；u,ug分别为流速矢量和内部旋转区旋转速度矢量；A为控制体的表面区域矢量；S为比焓源项，由黏性耗散、可压缩性和辐射等引起.

1.2 湍流模型

本文采用剪应力输运SSTk-ω湍流模型，该模型是标准k-ω模型修正后的两方程模型，综合了标准k-ω模型能够很好地模拟边界层内的低雷诺数流动,以及标准k-ε模型能够很好地模拟边界层外的充分发展的湍流流动的优势，能适应压力梯度变化的各种复杂物理现象.该模型对湍流黏性系数进行了修正，且考虑了湍流剪切应力从而不会对涡流黏度进行过度预测，因而能更好地适用于跨声速及超声速来流下的弹头扰流场的模拟.

SSTk-ω湍流模型的数学表达式为

(2)

(3)

1.3 计算模型及边界条件

1.3.1计算模型

计算模型为带膛线刻槽的某12.7 mm穿甲燃烧弹，该弹头主要包括弧形部，圆柱部和尾锥部3个部分，如图1所示，圆柱部直径D=13 mm.如图2所示，将计算域分为3部分，从内到外分别为边界层、内部旋转区和外部固定区.

由于弹头表面刻槽增加了其几何形状的复杂性，为增强解的收敛性与稳定性，全弹面网格及全计算域空间网格分别采用正交性较好的四边形网格和六面体网格，并对有刻槽的弹头圆柱部网格进行了局部加密.综合多次计算并收敛后的结果，对边界层划分了20层网格，第一层网格厚度满足y+≤0.5，相邻网格厚度沿径向的增加率小于1.15，边界层、内部旋转区和外部固定区网格数分别为100万、140万和90万.

1.3.2边界条件

根据枪管寿命试验及文献[11]，一个枪管寿命周期内，弹头的初速变化很小，试验时采用的100 m靶道内弹头的速度下降值也很小，而弹头出膛时的运动姿态改变更为明显，后者是导致100 m靶处椭圆弹孔率超标(长轴与短轴比大于1.25)进而判断枪管寿终的主要原因.弹孔的长短轴之比由入射角、攻角和弹头的外形尺寸决定，自动武器多采用低伸弹道，100 m范围内几乎相当于是平射，在入射角很小时，根据几何关系可得，当弹头攻角大于3°时，就已满足弹孔长轴与短轴之比大于1.25.因此，本文采用的计算条件为：来流马赫数Ma为2.2，攻角α为0°～3°，量纲一化转速ω*=0.215(ω*=ωD/v∞)，其中：ω为弹头旋转角速度；D为弹头圆柱部直径；v∞为来流速度.

弹头附近网格划分如图3所示，弹头表面采用无滑移壁面边界条件，内部旋转区与弹头壁面相关联，随壁面以相同角速度旋转；内部旋转区与外部固定区通过交界面传递数据，交界面采用滑移边界条件；外部固定区的外边界采用压力远场边界条件.

1.4 6自由度刚体外弹道方程

由于12.7 mm弹头的长径比较小(弹长与弹径之比约为5)，其外弹道过程可以近似为6自由度的刚体外弹道过程，弹头在外弹道飞行过程中受到的外力主要有重力、空气阻力、升力、马格努斯力等，外力矩主要有俯仰力矩、马格努斯力矩、极阻尼力矩、赤道阻尼力矩等.弹头在诸外力、外力矩的作用下运动过程满足的方程为

(4)

(5)

(6)

式中：m为弹头质量；v为弹头质心速度矢量；x为弹头质心空间位置矢量；F为合外力矢量；G为弹头角动量矢量；M为合外力矩矢量.建立合适的辅助坐标系并将各矢量分解即可建立弹头的6自由度刚体外弹道运动方程，将弹头气动特性数值模型计算得到的各气动参数代入式(4)～(6)，结合弹头的初始运动状态即可基于4阶龙格-库塔法编程求解弹头的外弹道运动过程.

2 损伤枪管所发射弹头表面形貌的获得方法

要研究枪管内膛损伤对弹头气动特性的影响，首先需要获得各寿命阶段枪管所发射弹头出膛时的表面形貌，具体流程如下：

① 根据系统的枪管寿命试验中获得的多根相同材料、工艺和结构的枪管在相同射击规范下完成不同射弹量的射击试验过程后的内膛尺寸数据、内窥和解剖过程的视频及照片等，总结枪管在各寿命阶段其内膛损伤的主要形式及损伤沿枪管轴向的分布规律；

② 根据枪管的内膛损伤形式及其分布规律因射弹量增加和轴向位置不同而导致的差异性，在枪管内膛预置诸如裂纹、烧蚀坑及镀铬层剥落等损伤，完成不同寿命阶段损伤枪管有限元模型的建立.

③ 考虑内膛损伤导致的弹后空间的改变，建立弹头在损伤枪管内完成挤进和沿线膛运动直至出膛的热力耦合模型，并将计算结果与试验结果进行对比以验证模型的准确性.

④ 基于该热力耦合模型对弹头沿射弹量分别为0,1 400,3 000和6 000发(分别对应了枪管的无损伤状态、寿命中前期状态、寿命中期状态和寿终状态)的枪管内膛运动直至出膛的过程进行仿真分析，获得4个寿命阶段枪管所发射弹头出膛时的表面形貌状态和运动姿态等.

内膛损伤数据的分析及损伤形式和分布规律的总结过程、有限元模型的建立过程及弹头沿4个寿命阶段枪管运动的内弹道计算过程详见参考文献[11].最终获得了上述4个寿命阶段的枪管所发射弹头的表面形貌状态如图4所示(右侧为弹头头部方向)以及弹头表面刻槽尺寸如表1所示.

表1 弹头出膛时表面刻槽尺寸

将有限元分析结果中的弹头网格导入Hypermesh或ICEM-CFD等软件进行流场计算域网格划分后，即可利用Fluent对4种不同表面形貌的弹头的气动特性进行仿真分析.

3 数值计算模型的实验验证

为验证本文所建立的高速旋转弹头气动特性数值计算方法的有效性和准确性，选取DeSpirito在文献[4]及Plostins在文献[12]风洞实验中采用的弹头外形，对其进行气动特性的数值仿真，并将仿真结果与风洞实验结果进行对比.DeSpirito及Plostins实验中采用的弹头外形结构与本文中的12.7 mm弹头相似，如图5所示.

数值计算条件为：来流马赫数Ma在0.6～4.5之间，攻角分别为0°和3°，弹头量纲一化转速ω*=0.17，力矩参考点取在弹头质心，距弹头头部49.9 mm.选取数值计算结果中攻角为0°时的阻力系数以及攻角为3°时的俯仰力矩系数导数及法向力压力中心位置等气动特性参数，其数值模拟结果与实验结果的比较情况如图6～图8所示.从3个气动特性参数的对比曲线图上可以看出，数值计算得到的气动参数的变化规律与风洞试验得到的结果相一致，并且数值计算结果与风洞试验所获得的相应气动参数值的最大误差均在10%以内，这说明本文所建立的高速旋转弹头的气动特性数值计算模型是准确可信的.

4 气动特性计算结果及外弹道过程分析

为了更深入地研究枪管内膛损伤对弹头外弹道过程的影响，本节基于上述数值计算方法获得了4个寿命阶段的枪管所发射的具有不同表面形貌的弹头的气动特性参数之间的差异及变化规律，并结合6自由度刚体外弹道方程，分析了随机扰动作用下，12.7 mm机枪枪管处于4个寿命阶段时，其发射的弹头在100 m靶处的散布圆半径和椭圆弹孔率(弹孔的长轴与短轴之比大于1.25即为椭圆弹孔).

4.1 内膛损伤对弹头的气动特性的影响

将4个不同表面形貌的弹头依次编号为1，2，3和4号，分别对应于射弹量0，1 400，3 000和6 000发的枪管所发射的弹头.

图9给出了4个弹头在Ma=2.2，ω*=0.215的飞行状态下，法向力压力中心距离弹头质心的距离(以质心为原点，指向弹头头部方向为正)随攻角α的变化规律.从图9中可以看出，4个不同表面形貌的弹头的压心均在质心前方，且随着攻角的增加，压心总体均呈现后移的趋势.此外还可以看出，损伤程度越严重的枪管所发射的弹头，其压心位置就越靠近弹头头部，这一现象在身管寿命的中后期(对应于3，4号弹头)尤为明显，根据文献[13]中弹丸急螺稳定系数计算公式(见式(7))可知，由内膛损伤的加剧导致的弹头压心的前移会加大压心质心的距离，降低弹头的飞行稳定行，使散布直径和椭圆弹孔率增加.

(7)

式中：σ为急螺稳定系数；h为压心距质心距离；D为弹头直径；H(y)为气重函数，取决于弹道高度；v为弹头速度；kmz(Ma)为翻转力矩的速度函数；α为攻角；A、B分别为弹头的极转动惯量和赤道转动惯量；ω为弹头转速；g为重力加速度.

图10～图12分别对应为来流马赫数为2.2、量纲一转速为0.215下的4个不同表面形貌弹头的阻力系数、升力系数和马格努斯力矩系数(规定马格努斯力矩方向竖直向上为正)随攻角变化规律.由图10可以看出，随着弹头表面形貌较初始弹头的改变量的增大，阻力系数也逐渐增大，且攻角越大，阻力系数的增大量也越大，这主要是因为随着枪管内膛损伤的不断发展，其发射的弹头表面刻槽的深度在不断降低[11]，导致了弹头阻力面积的增大，进而引起了阻力系数的上升.由图11、图12可以看出，从绝对数值上来看，各寿命阶段枪管所发射弹头的升力系数和马格努斯力矩系数在各攻角下均随射弹量增加而增加，且攻角越大不同弹头间气动参数差值越大，其中升力系数的增加原因与阻力系数相同，即刻槽变浅导致的受力面积的增加导致了气动参数的增加；马格努斯力矩系数增加的主要原因为内膛损伤枪管所发射的弹头，其表面除规则的阳线刻槽外，还会出现由内膛损伤导致的不规则分布的局部刻痕，内膛损伤还会划破弹头被甲材料，破坏其表面完整性，这会使得弹头高速旋转时边界层厚度的非对称畸变、径向压力梯度的非对称畸变等增大，从而引起马格努斯力矩系数的增大.

结合图10～图12及寿终枪管和无损伤枪管所发射弹头气动参数在α=3°时对比表(表2)可以看出，4号弹头的各气动参数相对1号弹头的改变量在α=3°达到最大，相对变化量均在10%以上，对弹头外弹道过程的影响不可忽略.

表2 攻角为3°时1号、4号弹头气动参数对比

4.2 外弹道过程仿真分析

弹头发射过程中，造成弹头出膛时运动状态产生随机波动的原因有很多，对于某一支特定的自动武器来说，将其架在固定枪架上进行射击时，枪械自身的随机因素对弹头的影响已降至较低，使弹头出膛状态产生随机波动进而产生射弹散布的主要原因在于弹药自身参数在一定范围内的波动，主要有弹头外形尺寸偏差、弹头质量及质量偏心偏差、弹头壳材料力学性能偏差、发射药主控参数偏差等.

本文针对上述5个会对弹头出膛扰动随机响应产生较大影响的因素，依据设计要求和工程实践经验，给出了各因素随机波动的范围如表3所示.5个因素分别取11个水平值，若按照正交试验方法进行因素不同水平值间的组合，则每个寿命阶段的枪管均需要进行112次理论及仿真计算.为避免参数组合情况太多带来的时间成本的增加，同时又要充分考虑各因素的各水平值对弹头初始扰动的影响，采用均匀设计方法进行各因素间的组合[14]，均匀设计能在较好地反映试验体系主要特征的前提下大大地减少试验次数，且相对于正交设计不会产生大的偏差.基于均匀设计方法获得了11组弹头主要随机因素的组合状态，并将这11个不同初始状态的弹头分别代入第2节建立的4个寿命阶段枪管内弹道过程计算模型中，获得了4个寿命阶段的枪管所发射弹头出膛扰动的随机响应结果.

表3 弹药随机因素取值范围表

将1～4号弹头的各11组膛口随机扰动状态作为初值代入式(4)～(6)，结合气动特性参数对各阶段弹头的外弹道过程进行求解，获得了4个寿命阶段枪管所发射的各11发弹头在100 m处的散布密集度R50和椭圆弹孔率.判断弹孔是否为椭圆弹孔的方法如图13所示，图13截面为弹头轴线和质心速度所在平面(椭圆弹孔长轴所在平面)，由6自由度外弹道方程可以得到100 m处弹头的入射角值(质心速度和靶面法线夹角)及攻角值(质心速度和弹头轴线夹角)，根据几何关系即可求得弹孔的长轴长度，短轴长与弹头圆柱部直径相等，二者之比大于1.25即为椭圆弹孔.

各阶段枪管R50、椭圆弹孔率与枪管寿命试验结果的对比情况如表4所示(寿命试验的R50和椭圆弹孔率取多次试验的平均值，每次发射20发弹头)，以寿终枪管为例，其计算得到100 m靶处弹着点分布情况如图14所示，可以看出椭圆弹孔往往也偏离散布中心较远.枪管寿终的判断标准为R50≥30 cm或者椭圆弹孔率≥50%.

表4 各阶段枪管散布密集度和椭圆弹孔率

Tab.4 Dispersion intensity and ratio of elliptical bullet hole at different life period

枪管射弹量/发R50/cm椭圆弹孔率计算值试验平均值计算值试验平均值012.814.600140013.013.400300015.115.600600021.722.854.5%60%

从表4可以看出，散布密集度R50和椭圆弹孔率的计算值均与枪管寿命试验实测值较好地符合，这进一步验证了损伤弹头的获得方法及气动特性参数的计算模型的准确性.从表4中还可以看到，试验过程中无损伤枪管的R50值反而大于射弹量1 400发的枪管的相应值，而数值仿真结果则恰恰相反，这也导致了无损伤枪管R50的数值计算结果与试验结果误差较大，这是因为试验时的无损伤枪管内膛表面存在制造和镀铬时留下的细小毛刺等缺陷，导致了弹头出膛扰动增大，经过一定次数的射击后，细小毛刺被磨去，而且内膛的损伤还未充分形成，因此在身管寿终的中前期，R50往往会出现不升反降的现象.计算值与试验情况均表明，该12.7 mm重机枪枪管寿终表现为椭圆弹孔率的超标，枪管在寿命中前期R50的增大量较小，且基本不会出现椭圆弹孔，在枪管寿命的后期，R50显著增加以及椭圆弹孔率的急剧增加导致了枪管的寿终，这是因为枪管寿命后期弹头出膛时的初始攻角等膛口扰动量会显著增加[11]，各气动特性参数在大攻角时的明显改变(如压心前移，马格努斯力矩系数的增大等)也会显著降低弹头外弹道过程的飞行稳定性，从而进一步增大了弹头的攻角，形成了恶性循环，导致了枪管的迅速寿终.