0 引 言

子母弹武器系统通过在战区上方散布大量子弹药，有效地提高了杀伤范围，实现了打击纵深大面积多目标，被广泛应用于区域战争中^[1]。其抛撒分离过程是子母弹武器系统研制中的关键环节，特别是在超声速飞行中，母弹开舱后外形的改变、母弹与子弹间的激波干扰使得流场结构异常复杂。为此分析子弹药抛撒分离过程，探讨其分离过程中流场结构特性、气动特性及干扰流场生成机理，对保证子母弹安全有效分离具有重要意义。

国外对于多体分离相关问题研究工作开展的较早，并开展了一系列风洞试验研究工作^[2-5]。Cavallo Peter A和Lee Robert^[6]基于有限体积法采用非结构动网格方法，对子母弹分离进行了数值求解，计算结果与试验结果表现了良好的一致性。Magdi Rizk^[7]等结合嵌套网格技术编译了解决多体分离问题的专用代码，得到了较好的计算结果。Bhange^[8]等基于6DOF弹道方程对穿甲弹脱壳过程进行了数值仿真，捕捉了流场细节特征，并提出来优化设计准则。

国内，雷娟棉，苗瑞生^[9]等以三维N-S方程为基础，运用CFD数值模拟方法，对火箭子母弹分离初期不同位置上定常流场进行了数值计算，并进行了风洞试验研究^[10]，分析了子弹在穿越母弹激波过程中干扰流场特性。张玉东、纪楚群^[11-12]采用有限体积法，通过耦合求解三维非定常Euler方程和6DOF弹道方程，对子母弹分离非定常流场进行了系统研究，分析了子母弹干扰流场气动特性。随后，陶如意^[13-14]等对双层弹舱结构下子母弹时序抛撒定常流场进行了数值模拟，对其干扰流场结构及气动干扰特性进行了分析，并结合风洞试验进行了对比。

从国内外公开发行文献来看，对于子母弹分离气动干扰特性研究较多，但在数值模拟时主要针对受激波干扰严重的前舱子弹分离过程展开，而子母弹武器系统实际包含多层弹舱，且抛撒时通常采用时序控制，各舱段子弹按照时序间隔依次分离。各舱段子弹由于轴向排布位置不同，在不同舱段分离时弹体气动特性差异较大，因而各个舱段子弹药的运动特性及分离姿态并不相同。因此在前人研究基础上，为了获得子母弹在时序抛撒过程中子弹在不同舱段下分离流场特性及运动特性，为工程应用提供最优抛撒方案，本文基于有限体积法，通过耦合求解ALE法描述下的三维非定常可压缩N-S方程及6DOF刚体运动方程，首先对位于不同舱段下子弹分离非定常流场进行了数值求解，其中控制方程及湍流方程采用二阶迎风格式进行空间离散，时间离散采用隐式方案，揭示了不同舱段下子弹分离过程流场特性及气动干扰特性。

1 控制方程及耦合求解方法 1.1 控制方程

ALE有限体积法描述下的三维非定常可压缩N-S方程其积分形式表达式如下：

式中，Ω为控制体积，S为控制体表面边界，n为控制体边界外法向单位向量，守恒变量Q及对流项为：

控制方程右端为作用在控制体表面的粘性通量：

式中ρ、p、e分别为控制体内流体密度、压强及比内能，n_x，n_y，n_z分别为外法向单位向量n的三个分量，u、v、w分别为速度在三个方向上的分量，q_x、q_y、q_z为微元热流量。其中

式中，x_t、y_t、z_t分别为网格移动速度在三个方向的分量。为使方程组封闭，引入气体状态方程：

1.2 耦合求解方法

为实现子母弹分离过程流场区域的更新，将刚体六自由度动力学方程耦合流体控制方程求解，运用动网格技术以达到动边界区域网格的运动更新，六自由度刚体动力学方程组详见文献[15]。耦合求解过程如图 1所示，具体过程如下：

(1)为保证非定常计算的准确性，首先对子母弹初始流场进行定常求解，得到其流场气动载荷分布，作为非定常流场计算的初始流场，并沿弹体表面对压力进行积分求得子弹所受气动力。

(2)以子弹当前气动载荷以及子弹初始运动参数，包括抛撒速度、质心位置、角速度等作为流场计算的初始条件，调用六自由度动力学方程求解，通过对时间的积分获得下一时刻子弹速度、角速度及姿态角。

(3)运用弹簧光顺法和局部网格重构法，根据新的边界值对子弹动边界进行网格的更新与重构，并进行下一步迭代计算，从而实现整个流场区域的更新。

2 计算模型及条件 2.1 计算模型

图 2为子母弹双层弹舱抛撒计算模型，母弹包含前后两层弹舱，每层弹舱装配4枚子弹药，为节约计算时间，提高计算效率，取整个区域的1/4作为计算区域，如图所示。其中子弹直径D=170mm，长度L=595mm。母弹直径D=440mm，长度L=1200mm，中间隔板厚度δ=10mm。

外流场区域整体采用四面体非结构网格划分，整体网格数量约为200万。为保证计算过程中子弹动区域周围网格质量，在子弹表面生成棱柱型贴体网格，并对弹丸周围部分区域加密网格，以提高求解精度。图 3为计算模型局部网格放大图。

2.2 边界条件

边界条件包括固壁边界及流场外边界，其中母弹为静止固壁，子弹为移动固壁，物面上切向速度为零，法向为无穿透条件。来流边界取自由来流边界条件，出流边界采用场内外推处理。子弹初始抛速v=20m/s，角速度ω=200°/s。流场初始参数为：Ma=3、P=101325Pa、T=300K。

2.3 算例验证

为了验证上述非定常数值计算方法的可靠性，本文选择常用来检验动网格计算方法的NACA0012翼型俯仰振动绕流问题进行数值验证，翼型俯仰振动迎角随时间变化规律为:

式中，迎角平均值α₀=0.0016°，振幅α_m=2.51。无量纲角频率k=ωc/2U_∞，c为弦长，U_∞为自由来流速度。俯仰轴心x_m/c=0.25，来流马赫数M_∞=0.755。

计算结果如图 4所示，图 4(a)、(b)分别为升力系数及俯仰力矩系数随迎角变化的曲线，与文献[16]中实验值进行了对比。计算结果与实验结果吻合良好，验证了本文数值方法的可靠性。

3 计算结果与分析 3.1 流场结构

图 5、图 6分别为分离初期流场速度分布云图及流线图。在超声速子母弹分离过程中，其干扰流场结构主要有母弹激波、子弹激波及弹体间的反射波组成，由图 5、图 6可以看出，在分离初期前舱子弹头部逐渐进入母弹激波，并在子弹头部形成高压区。由于子弹弹体的阻塞作用，部分母弹激波绕过弹体上表面随主流激波向后发展。同时部分反射激波流入子弹与母弹分离形成的弹槽空腔中，此时由子弹反射激波、弹槽反射波在子弹弹体下表面形成一个三角形低压区域。并在母弹弹肩及隔板的影响下，分别在弹槽前端及隔板处形成低速度、低密度及低压力回流区，产生涡流，随着气流的进一步膨胀流出，在子弹尾部同时形成涡流。由于此时子弹与母弹间隙较小，导致壅塞现象出现，在子弹弹尾及隔板处形成喉状激波，使得母弹和子弹表面对应区域压力增加，其压力值随着子弹径向位移的增加而迅速减小。

图 7为前舱子弹分离过程压力分布云图。由图 7可看出，随着子弹径向位移的增加，前舱段子弹在分离时刻t=10ms后子弹头部摆脱母弹激波，在超声速来流下子弹头部形成激波干扰。此时母弹激波逐渐作用于子弹头部下表面，母弹激波与子弹激波产生相交波系，并在子弹下表面形成反射激波，反射波向下发展作用于母弹弹舱表面再次形成反射。随着子弹位移的进一步增加，母弹激波作用点扫过子弹弹身并向子弹尾部移动，子弹逐渐摆脱母弹激波影响。整个分离过程中，母弹激波的干扰作用不仅为子弹分离提供了轴向力与径向力，而且为子弹的俯仰运动提供了力矩。

由图 8后舱段子弹分离过程压力分布云图可以看出，与前舱子弹相比而言，由于后舱子弹轴向位置靠后，其干扰力相比前舱子弹较弱，母弹激波提供的轴向力与径向力较小，在分离时刻35ms时子弹头部才摆脱母弹激波影响，此时母弹激波与后舱子弹激波产生相交干扰，并作用于弹体下表面，为后舱子弹提供俯仰力矩。

图 9为不同舱段下子母分离轨迹图，由图中可以直观地看出，与后舱子弹相比，在整个分离过程中，前舱子弹姿态变化较为明显，能够迅速建立有效的径向分离距离及分离速度，而后舱子弹由于长期处于母弹激波影响下，姿态变化程度较低，整个分离过程较为缓慢。

3.2 弹体间干扰特性分析

子弹抛撒过程中，不仅受母弹气动力影响，子弹间也存在相互作用影响。当抛撒时序间隔较长时，后舱子弹抛出时前舱子弹已经分离，前后弹舱段子弹间干扰可以忽略，而弹舱间装配的四枚子弹干扰影响仍然存在，为此分别对前后舱段装配的四枚子弹相互间干扰特性进行分析，计算结果如图 10、图 11所示。

图 10、图 11分别为前舱子弹、后舱子弹截面处速度分布云图，由计算结果可以看出分离1ms时弹体间激波间产生相交干扰，10ms后前舱子弹激波锥面彼此相对独立，由于后舱子弹轴向位置靠后，5ms时后舱子弹间干扰已经相对很弱，由计算结果可以得出，子弹间干扰主要存在于分离初始阶段，此时周向排布子弹较为紧密，子弹间激波有相交干扰，随着弹体位移的增加，5~10ms后弹体间干扰减弱，此后可认为弹体间无相互干扰。

3.3 运动特性分析

图 12、图 13分别为子弹沿轴向X方向和径向Y方向相对速度变化曲线。由图中可以看出，由于分离过程中来流阻力作用，子弹轴向速度逐渐增加，整体呈线性增长趋势，其径向速度分布呈先减小后增加的趋势，这是由于子弹穿越母弹激波过程中，弹体上表面压力高于下表面以及克服自身重力影响导致。

结合图 7、图 8可以看出，前舱子弹摆脱母弹激波时间较快，其轴向和径向分离速度整体高于后舱子弹。

图 14为子弹分离过程中升力系数C_l变化曲线，由图中可以看出，初始阶段由于子弹弹身上表面区域压力高于子弹下表面与弹舱形成的低压区域以及子弹克服自身重力影响，分离初期其升力系数呈下降趋势。由前舱子弹升力系数变化曲线可以看出，随着子弹逐渐进入母弹激波，母弹激波作用于子弹下表面，为子弹提供了轴向力，使得子弹升力系数急剧上升。伴随着子弹位移的增加，母弹激波作用点扫过弹身，向子弹尾部移动为子弹提供了负向俯仰力矩，即使得子弹低头，子弹在t=15ms时前舱子弹升力系数达到最大值，随后其升力系数逐渐减小。在25ms时刻再次增加，这是由于子弹完全摆脱母弹激波影响，在自身角速度作用下，子弹迎角增加导致，其整个振荡过程表征了子弹在不同时刻所受气动干扰特性和姿态的变化的过程。与前舱子弹相比，在相同的分离时间内，后舱子弹长期处于母弹激波下，升力系数相比前舱较低，其整体变化过程较为缓慢。

图 15给出了各舱子弹阻力系数C_d变化曲线，可以看出，由于后舱子弹轴向位置靠后，母弹激波干扰力减弱，轴向及径向受力较小，因而在分离初始阶段其阻力系数相比前舱子弹较低。由前舱子弹变化曲线来看，随着子弹逐渐进入母弹激波，在5ms时刻其弹体头部完全进入母弹激波，其阻力系数达到最大值，之后头部摆脱母弹激波影响，其阻力值逐渐下降，并在25ms后再次增加，其变化过程同样是由于子弹摆脱母弹激波后自身翻转角速度导致姿态角增加所致。后舱子弹其阻力系数同样在弹体完全进入母弹激波时达到最大值，之后逐渐减小呈现先增加后减小的振荡过程。

图 16为子弹分离过程中角速度随时间变化曲线，通过前后舱角速度变化曲线可以看出，前舱子弹弹身在摆脱母弹激波过程中，母弹激波主要作用于子弹上表面，使得子弹产生沿Z轴正方向的俯仰力矩，即使子弹低头的趋势，其角速度值逐渐增加。结合图 7可以看出，在子弹头部摆脱母弹激波后，由于母弹激波干扰并作用于子弹下表面，为子弹提供了沿Z轴负方向俯仰力，使得姿态角增加，角速度值再次减小。后舱子弹由于位置靠后，受激波影响力较小，在自身角速度及气动干扰作用下，使得子弹沿Z轴负方向转动，子弹姿态角增加，角速度减小。在15ms后后舱子弹逐渐穿越母弹激波，其角速度同样呈现先增加再减小的趋势。通过角速度变化曲线可以看出，角速度振荡过程表征了子弹穿越母弹激波及子弹分离姿态的整个动态变化过程。

4 结 论

本文以双层弹舱轴向排布子母弹系统为模型，采用有限体积法耦合求解流体控制方程及6DOF刚体运动方程，运用动网格技术对时序抛撒方式下的不同舱段子母弹抛撒分离过程进行了数值模拟，获得了不同舱段下子母弹分离干扰流场特性，分析了子弹气动力系数随时间的变化关系。

研究结果表明，在分离初期，前舱子弹与母弹之间产生强耦合作用，在子弹与母弹弹舱间形成壅塞和激波的多次反射现象。随着子弹位移的增加，子弹与母弹弹舱间气动干扰作用逐渐降低。在分离的中后期，伴随着母弹激波逐步扫过子弹弹身表面，其气动干扰为子弹提供了俯仰运动力矩，造成子弹升、阻力系数曲线产生振荡。与前舱子弹相比，由于后舱子弹位置靠后，母弹激波干扰力减弱，造成分离初期子弹轴向力及径向力较小，气动参数相比前舱子弹较低，整体气动力参数变化过程较为缓慢。

文中通过对子弹在不同舱段下分离的三维非定常流场进行了数值仿真，得到了子弹在不同舱段下分离流场特性，但在计算过程中忽略了子母弹在短时序抛撒过程中前舱子弹对后舱子弹气动干扰特性的影响，描述不够细致。在下一步工作中，将针对子母弹在多舱段低时序抛撒过程中，前舱子弹对后舱子弹气动干扰特性及流场特性的影响进行进一步研究。