基于五阶WENO格式的燃气在药床中流动过程二维两相流研究

程申申; 王浩; 薛绍; 陶如意

doi:10.11883/bzycj-2020-0200

基于五阶WENO格式的燃气在药床中流动过程二维两相流研究

doi: 10.11883/bzycj-2020-0200

南京理工大学能源与动力工程学院，江苏南京 210094

详细信息

作者简介:
程申申（1994－　），男，博士研究生，1447202650@qq.com

通讯作者:
陶如意（1978－　），女，博士，副研究员，硕士生导师，tao801801@163.com

中图分类号: O359
计量
- 文章访问数: 699
- HTML全文浏览量: 361
- PDF下载量: 50
- 被引次数: 5
出版历程
- 收稿日期: 2020-06-16
- 修回日期: 2020-10-09
- 网络出版日期: 2021-05-31
- 刊出日期: 2021-06-05

Two-dimensional numerical simulation on gas-solid two-phase flow induced by combustion gas flow in a chamber based on a fifth-order WENO scheme

School of Energy and Power Engineering, Nanjing University of Science and Technology, Nanjing 210094, Jiangsu, China

摘要

摘要: 为研究内弹道初始阶段中心点火管燃气在膛内药床中的流动特性和传播规律，设计了可视化点传火实验平台，并进行了膛内假药床的点传火实验。基于加权本质无震荡(weighted essentially non-oscillatory, WENO)格式，构造了膛内轴对称二维内弹道两相流模型，对膛内燃气在假药床中的流动过程进行数值模拟。计算结果与可视化实验结果符合较好，全局压力平均误差为5.35%。表明数值计算准确地描述了燃气流动特性，完整地呈现了点火管燃气在假药床中的发展过程。在点火初始阶段，膛内压力径向效应明显，气相沿径向传播较快，药床药粒基本不会发生运动；随着燃气逐渐在膛内传播，膛内压力呈现径向一致、轴向梯度分布的特征，在压力梯度作用下，气相轴向速度开始占据主导，径向速度在膛底和中部区域减小为零，而固相速度随气相速度变化而变化；气相在到达弹底前，由于固相颗粒的壅塞，会提前出现速度反向波动现象。
- 内弹道 /
- 两相流 /
- WENO格式 /
- 点传火实验
Abstract: In order to explore the flow characteristics and propagation law of combustion gas from the central ignition tube in the initial stage of internal ballistic, a visual experimental platform was designed to carry out ignition experiments with the substituted particle bed in the chamber. A high-speed camera system was used to capture the gas flow and flame propagation in the chamber, and a dynamic data acquisition and analysis system with pressure sensors was applied to record the pressure data at characteristic positions in the chamber. A two-dimensional, axisymmetrical, two-phase flow model of internal ballistics was constructed to simulate the flow process of the gas in the substituted particle bed based on a weighted essentially non-oscillatory (WENO) scheme, and the time term was determined by the third-order TVD Runge-Kuta method. The calculated results are in good agreement with the visual experimental results, and the global pressure average error is 5.35%. It indicates that the numerical simulation can accurately describe the gas flow characteristics and present the development process of the gas from the ignition tube in the substituted particle bed. The radial effect of the chamber pressure is obvious, and the gas moves rapidly along the radial direction, and the substituted particle basically does not move in the initial stage of ignition. Moreover, with the gradual propagation of the gas in the chamber, the chamber pressure is characterized by a radial uniformity and an axial gradient distribution. Under the action of the pressure gradient, the axial velocity of the gas phase begins to dominate, and furthermore, the radial velocity decreases to zero in the bottom and the middle region of the chamber, while the solid phase velocity varies with the gas phase velocity. In addition, before the gas reaches the bottom of the right end, the inverse velocity fluctuation appears in advance due to the solid particle congestion.
- internal ballistic /
- two-phase flow /
- WENO scheme /
- ignition experiment

HTML全文

膛内点传火过程是一个复杂的高温高压传热传质过程，包括点火管火药燃烧、点火管燃气流出以及膛内药床的传热燃烧过程，对火炮的药床燃烧以及后续推进过程有重要的影响。在实际实验中，现有测试手段很难得到实际状态下这一阶段的实验现象。因此，开展药床燃气流动的内弹道模拟对整个火炮发射过程的研究具有重要意义^[1-3]。

近年来针对点火管和膛内火药床燃烧的研究成果有很多。王浩等^[4]通过点传火特性实验，对比了不同装填条件下的传火管多点压力曲线及火焰传播速率，分析了底火焰在传火管中的传输特性，是早期对点传火进行的实验研究。王珊珊等^[5-6]对点火管和中心管建立了一维两相流模型进行内弹道两相流研究，分析了不同结构点火管的点传火性能。郭梦婷等^[7]对中心炸管定容阶段进行了二维两相流研究，探究了中心管式抛撒机构的内弹道流场特性。但是以上研究都着眼于点火管内火药燃烧和传播，对于燃气在膛内流动及火药燃烧的研究较少。Miura等^[8]使用带破裂膜片的燃烧装置进行实验，并采用两相流体动力学代码进行模拟计算，用于探索膛内负压差产生的机理,并得出结论,压差是由固体推进剂运动导致的，压差波动取决于底火喷口的流量和喷口的位置，但实验和模拟过程未对火焰发展及流场规律进行详细描述和观测。Jaramaz等^[9]对整个火炮膛内流动过程进行实验和一维两相流研究，但是未考虑径向效应。程诚等对某制导炮弹二维两相流内弹道性能进行研究，并采用高阶黎曼近似模型求解两相流方程，详细研究了带推进系统的整个内弹道过程中膛内流场特性及其发展^[10-11]，但主要研究火药床燃烧及弹丸推进过程，并未仔细对初始阶段点传火过程进行流场分析。因此，开展中心传火管装药点火过程的理论与试验研究，对进一步分析火炮发射起始阶段膛内的点火过程燃气流动规律，降低起始压力波的产生，保证火炮射击的安全性至关重要。

本文中，拟设计低压可视化实验装置，观察点火药燃气在药床中的流动过程。考虑到可视化装置的承压能力，采用假药粒装填结构，并结合二维两相流数值模拟来预测这一阶段的实际过程。使用五阶精度的WENO格式，相比于经典的Mac Cormack格式与TVD格式，其具有精度高、间断分辨清晰，而且效率高的优点，是十分有效的算法^[12-14]。以期研究结果对进一步开展火炮起始阶段燃气流动及压力变化规律研究具有应用价值。

1. 数理模型

1.1 物理模型

对膛内初始阶段点传火过程建立的计算模型如图1所示。膛内填充由尼龙材料仿真实药粒大小制成的假药粒，中心为点传火管，内装小粒黑火药，左端用端盖密封，装有底火以点燃点火管内火药，右端用铝合金膜片密封。

图 1 装药结构示意图

Figure 1. Schematic of the charge structure

下载: 全尺寸图片幻灯片

根据上述模型建立真实药粒运动燃烧和假药粒运动的两相流模型，为了简化内弹道数值模拟的复杂性，需对实际物理过程进行适当假设，基本假设如下：

（1）采用双流体模型假设，将发射药颗粒群作为具有连续介质特性的拟流体来处理;

（2）忽略湍流影响和气固两相间黏性；

（3）火药颗粒服从几何燃烧和指数燃烧定律^[2]；

（4）气相状态方程采用Nobel-Abel状态方程^[2]；

（5）假设底火内点火药完全以燃气形式加入点火管内，点火管采用一维两相流计算，破膜后源项在膛内当地释放；

（6）假药粒不产生燃气，尺寸与真实药粒相同。

1.2 膛内轴对称二维两相流方程组

对于内弹道膛内气-固两相流过程，由于膛内沿身管轴线方向呈对称形式，因此在内弹道膛内气-固两相流二维计算时，往往采用轴对称模型，守恒形式可以表示如下：

$\frac{\partial {{U}}}{\partial t}+\frac{\partial {{F}}}{\partial r}+\frac{\partial {{G}}}{\partial z}+\frac{{{H}}}{r}={{S}}$

(1)

式中：U为时间守恒矢量，F、G分别为径向和轴向的矢通量，H为轴对称变换矢量，S为源项。具体可表示如下：

$\begin{split} &{{U}}=\left[ \begin{array}{c}\phi {\rho }_{\rm{g}}\\ \left(1-\phi \right){\rho }_{\rm{p}}\\ \phi {\rho }_{\rm{g}}{v}_{\rm{g}}\\ \phi {\rho }_{\rm{g}}{u}_{\rm{g}}\\ \left(1-\phi \right){\rho }_{\rm{p}}{v}_{\rm{p}}\\ \left(1-\phi \right){\rho }_{\rm{p}}{u}_{\rm{p}}\\ \phi {\rho }_{\rm{g}}{E}_{\rm{g}}\end{array}\right],\qquad {{F}}=\left[ \begin{array}{c}\phi {\rho }_{\rm{g}}{v}_{\rm{g}}\\ \left(1-\phi \right){\rho }_{\rm{p}}{v}_{\rm{p}}\\ \phi \left({\rho }_{\rm{g}}{v}_{\rm{g}}^{2}+p\right)\\ \phi {\rho }_{\rm{g}}{u}_{\rm{g}}{v}_{\rm{g}}\\ \left(1-\phi \right)\left({\rho}_{\rm{p}}{v}_{\rm{p}}^{2}+p+{R}_{\rm{p}}\right)\\ \left(1-\phi \right){\rho }_{\rm{p}}{u}_{\rm{p}}{v}_{\rm{p}}\\ \phi {\rho }_{\rm{g}}{v}_{\rm{g}}\left({E}_{\rm{g}}+\dfrac{p}{{\rho }_{\rm{g}}}\right)\end{array}\right],\qquad {{G}}=\left[ \begin{array}{c}\phi {\rho }_{\rm{g}}{u}_{\rm{g}}\\ \left(1-\phi \right){\rho }_{\rm{p}}{u}_{\rm{p}}\\ \phi {\rho }_{\rm{g}}{v}_{\rm{g}}{u}_{\rm{g}}\\ \phi \left({\rho }_{\rm{g}}{u}_{\rm{g}}^{2}+p\right)\\ \left(1-\phi \right){\rho }_{\rm{p}}{v}_{\rm{p}}{u}_{\rm{p}}\\ \left(1-\phi \right)\left({\rho}_{\rm{p}}{u}_{\rm{p}}^{2}+p+{R}_{\rm{p}}\right)\\ \phi {\rho }_{\rm{g}}{u}_{\rm{g}}\left({E}_{\rm{g}}+\dfrac{p}{{\rho }_{\rm{g}}}\right)\end{array}\right] \\ & {{H}}=\left[ \begin{array}{c}\phi {\rho }_{\rm{g}}{v}_{\rm{g}}\\ \left(1-\phi \right){\rho }_{\rm{p}}{v}_{\rm{p}}\\ \phi {\rho }_{\rm{g}}{v}_{\rm{g}}^{2}\\ \phi {\rho }_{\rm{g}}{u}_{\rm{g}}{v}_{\rm{g}}\\ \left(1-\phi \right){\rho }_{\rm{p}}{v}_{\rm{p}}^{2}\\ \left(1-\phi \right){\rho }_{\rm{p}}{u}_{\rm{p}}{v}_{\rm{p}}\\ \phi {\rho }_{\rm{g}}{v}_{\rm{g}}\left({E}_{\rm{g}}+p/{\rho }_{\rm{g}}\right)\end{array}\right],\qquad {{S}}=\left[ \begin{array}{c}{\dot{m}}_{{\rm{i}}{\rm{g}}{\rm{n}}}\\ 0\\ -{f}_{r}+{\dot{m}}_{{\rm{i}}{\rm{g}}{\rm{n}}}{v}_{{\rm{i}}{\rm{g}}{\rm{n}}}+p\dfrac{\partial \phi }{\partial r}\\ -{f}_{z}+{\dot{m}}_{{\rm{i}}{\rm{g}}{\rm{n}}}{u}_{{\rm{i}}{\rm{g}}{\rm{n}}}+p\dfrac{\partial \phi }{\partial z}\\ {f}_{r}-p\dfrac{\partial \phi }{\partial r}\\ {f}_{z}-p\dfrac{\partial \phi }{\partial z}\\ -{f}_{r}{v}_{\rm{p}}-{f}_{z}{u}_{\rm{p}}-{A}_{{\rm{t}}}{Q}_{\rm{p}}+{m}_{{\rm{i}}{\rm{g}}{\rm{n}}}{H}_{{\rm{i}}{\rm{g}}{\rm{n}}}-p\dfrac{\partial \phi }{\partial t}\end{array}\right] \end{split}$

式中： $\phi$ 为孔隙率，ρ_g为气相密度, ρ_p为固相密度, u_g、v_g分别为气相轴向和径向速度, u_p、v_p分别为固相轴向和径向速度，p为气相压力；气相总能E_g=(v_g²+u_g²)/2+e_g, e_g为气相内能；m_ign为点火管流入气相单位体积质量流量，u_ign、v_ign分别为点火燃气轴向和径向速度，f_z、f_r分别为轴向和径向气-固相间阻力，R_p、Q_p分别为颗粒间应力和相间传热，H_ign为点火燃气的滞止焓, A_t为单位体积内固相颗粒等效比表面积。

提供辅助方程来封闭基本方程。膛内点传火过程源项由点火管提供，点火管为一维线源，从破膜后点火燃气流入膛内过程视为源项当地释放。源项与膛内燃烧过程进行耦合计算，具体一维点火管控制方程及流出流量公式可参考文献[5]，其他如状态方程、燃烧规律、相间阻力、相间传热、颗粒间应力、颗粒表面温度等辅助方程等参见文献[2]。

2. 数值计算方法

基于有限差分法，对上述模型在时间项上采用三阶TVD型龙格库塔法进行时间推进，空间项采用五阶WENO格式进行空间离散。

2.1 五阶WENO格式

WENO格式基本思路是将基架点分为多个组，每个组独立计算某一点的导数逼近，得到多个差分，由于流场中间断不可能处处存在，所以根据每个模板的光滑程度，设定权重，如果光滑度越高，则权重越大，存在间断则权重趋近零。

对于向量形式的守恒方程：

$\frac{\partial {{U}}}{\partial t}+\frac{\partial {{F}}}{\partial r}+\frac{\partial {{G}}}{\partial z}+\frac{{{H}}}{r}={{S}}$

(2)

采用全局Lax-Friedrichs分裂对空间项F、G进行流通矢量分裂：

${{{F}}}_{i,j}^{(\pm )}=\left[{{{F}}}_{i.j}\pm {\lambda }_{F}^{*}{{{U}}}_{i,j}\right]/2$

(3)

${{{G}}}_{i,j}^{(\pm )}=\left[{{{G}}}_{i.j}\pm {\lambda }_{G}^{*}{{{U}}}_{i,j}\right]/2$

(4)

式中：λ^*为全局最大特征值，为保险起见通常取 ${\lambda }_{F}^{*}=1.1\; {\rm{max}}\left[\left|{v}_{\rm{g}}\right|+c\right]$ ， ${\lambda }_{G}^{*}=1.1\; {\rm{max}}\left[\left|{\rm{u}}_{\rm{g}}\right|+{\rm{c}}\right]$ 。

利用分裂后的矢通量构造矢通量空间导数项的差分逼近式，以F项为例：

$\frac{\partial {{F}}}{\partial r}=\frac{\partial {{{F}}}^{(+)}}{\partial r}+\frac{\partial {{{F}}}^{(-)}}{\partial r}\approx \frac{{\rm{d}}{{{F}}}^{(+)}}{{\rm{d}}r}+\frac{{\rm{d}}{{{F}}}^{(-)}}{{\rm{d}}r}$

(5)

采用Jiang等^[14]的五阶精度WENO格式构造矢通量导数逼近式的通量形式为（以正特征值为例）：

$\frac{{\rm{d}}{{{F}}}_{j}^{(+)}}{{\rm{d}}r}=\frac{\left({h}_{j+1/2}-{h}_{j-1/2}\right)}{{\rm{d}}r}$

(6)

h_j+1/2为正通量，由3个通量加权平均而得：

${h}_{j+1/2}={w}_{1}{h}_{j+1/2}^{\left(1\right)}+{w}_{2}{h}_{j+1/2}^{\left(2\right)}+{w}_{3}{h}_{j+1/2}^{\left(3\right)}$

(7)

${h}_{j+1/2}^{\left(1\right)}$ 、 ${h}_{j+1/2}^{\left(2\right)}$ 、 ${h}_{j+1/2}^{\left(3\right)}$ 的表达式分别为：

${h}_{j+1/2}^{\left(1\right)}=\frac{1}{3}{{{F}}}_{j-2}^{(+)}-\frac{7}{6}{{{F}}}_{j-1}^{(+)}+\frac{11}{6}{{{F}}}_{j}^{(+)}$

(8)

${h}_{j+1/2}^{\left(2\right)}=-\frac{1}{6}{{{F}}}_{j-1}^{(+)}-\frac{5}{6}{{{F}}}_{j}^{(+)}+\frac{1}{3}{{{F}}}_{j+1}^{(+)}$

(9)

${h}_{j+1/2}^{\left(3\right)}=\frac{1}{3}{{{F}}}_{j}^{(+)}+\frac{5}{6}{{{F}}}_{j+1}^{(+)}-\frac{1}{6}{{{F}}}_{j+2}^{(+)}$

(10)

加权系数为（k=1,2,3）：

${\omega }_{k}=\frac{{\alpha }_{k}}{{\alpha }_{1}+{\alpha }_{2}+{\alpha }_{3}}, {\alpha }_{k}=\frac{{C}_{k}}{{\left({\varepsilon }_{0}+{I}_{k}\right)}^{\beta }}$

(11)

式中：C₁=1/10，C₂=6/10，C₃=3/10，ε₀=10⁻⁶，β=2，I_k为光滑度量因子。I_k的表达式为：

${I}_{1}=\frac{1}{4}{\left({h}_{j-2}-4{h}_{j-1}+3{h}_{j}\right)}^{2}+\frac{13}{12}{\left({h}_{j-2}-2{h}_{j-1}+{h}_{j}\right)}^{2}$

(12)

${I}_{2}=\frac{1}{4}{\left({h}_{j-1}-{h}_{j+1}\right)}^{2}+\frac{13}{12}{\left({h}_{j-1}-2{h}_{j}+{h}_{j+1}\right)}^{2}$

(13)

${I}_{3}=\frac{1}{4}{\left(3{h}_{j}-4{h}_{j+1}+{h}_{j+2}\right)}^{2}+\frac{13}{12}{\left({h}_{j}-2{h}_{j+1}+{h}_{j+2}\right)}^{2}$

(14)

当特征值为负时，利用对称性，将上式中所有下标“j+m”换成“j−m”即可。对G项的离散同理。

2.2 时间方向离散

空间导数离散化后的半离散方程为：

$\frac{\partial {{U}}}{\partial t}={{Q}}_{n}+{{R}}_{n}$

(15)

式中：Q_n为离散后空间导数项，R_n为源项。使用三阶TVD龙格库塔方法进行离散：

${{{U}}}^{\left(1\right)}={{{U}}}^{(n)}+{\Delta }t{{Q}}\left({{{U}}}^{(n)}\right)+{{R}}_{n}$

(16)

${{{U}}}^{\left(2\right)}=3/4{{{U}}}^{\left(n\right)}+3/4{{{U}}}^{\left(1\right)}+1/4{\Delta }t{{Q}}\left({{{U}}}^{\left(1\right)}\right)+{{R}}_{n}$

(17)

${{{U}}}^{(n+1)} =1/3{{{U}}}^{\left(n\right)}+2/3{{{U}}}^{\left(2\right)}+2/3{\Delta }t{{Q}}\left({{{U}}}^{\left(2\right)}\right)+{{R}}_{n}$

(18)

∆t为时间步长，可以根据下式确定：

${\Delta }t=\frac{C{\rm{d}}x}{{\left|{\lambda }^{*}\right|}_{{\rm{m}}{\rm{a}}{\rm{x}}}}$

(19)

式中：C为柯朗数，dx为计算网格尺度；λ^*为全局最大特征值，由上述算法可计算得出。

2.3 初始条件和边界条件

膛内装填条件采用如图1所示的装填结构，并结合实验装填数据和环境参数，确定模拟过程中的初始条件。膛底和出口处边界条件采用镜面反射方法，计算至膜片破裂。由于WENO格式采用多个基架点进行计算，在边界外设置3个虚拟点可将边界点转化为内点参与计算，还可以保持边界点的计算精度。中心轴线边界满足轴对称关系式。

3. 实验验证与数值结果分析

3.1 假药粒实验

为了验证所建立模型和数值方法的可靠性，对假药粒情况进行实验验证。采用了针对105 mm口径火炮药室缩比的可视化实验装置，如图2所示。装置一端端盖内侧连接点火管，外侧可以装底火，用以激发点火点燃点火管内的点火药；另一端端盖为中通结构，用以锁紧闭压膜片。装置中间装有耐高压玻璃管，管内装有定量的假药粒，通过装置中部预留的矩形窗口可以有效观测装置内部的实验现象。装置上方设置3个传感器接头，分别对应药室的膛底、中部和弹底（膜片）处，可以通过压力传感器测得这几处的压力。

图 2 可视化实验装置

Figure 2. Visual experimental system

下载: 全尺寸图片幻灯片

整体实验装置如图2所示，可视化窗口为400 mm×120 mm，高速摄像系统可通过此窗口拍摄实验过程中火焰发展传播过程。传感器距膛底分别为40、247.5和455 mm，可通过数据采集系统获取各位置压力。整体实验装填点火药为2#小粒黑火药，质量为30 g；假药粒为PA尼龙，质量为3 kg；膜片为铝合金，厚度为0.5 mm。

3.2 数值计算与实验数据比对

以上述实验装置进行可视化模拟实验的内弹道过程二维两相流动为研究对象，其具体结构可参见图1中的简化物理模型，计算参量设置为：点火药质量，30.0 g；点火药密度，1650 kg/m³；点火药燃速指数，0.98；点火药燃速系数，3.68 mm·s⁻¹·MPa^−0.98；火药力，300 J/g；点火药平均直径，1 mm；假药粒密度，1100 kg/m³；假药粒质量，3.0 kg；药室长度，495 mm；药室直径，110 mm；点火管长度，25 mm；点火管直径，19 mm；气体绝热指数，1.3；网格长度，5 mm；气体余容，1 cm³/g。

为了简化计算，忽略点火管的二维效应，即在传火孔破膜后将点火管作为线性点火源加入膛内进行两相耦合计算，弹底与膛底为垂直于轴向的固壁平面，采用反射壁面条件。通过压力传感器对膛内各位置采集的压力与数值模拟获取压力进行比对，由于模拟计算只到破膜后截止，因此只取破膜前的压力进行对比分析。

图3所示为中间测点压力随时间变化的实验结果与模拟结果的对比曲线。可以看出，模拟结果与实验结果符合较好，全局平均误差为5.35%，最高峰值处（即2.5 ms后区域）误差仅为1.42%，说明所建立的模型和采用的数值方法合理有效。

图 3 计算压力曲线与实验对比图

Figure 3. Comparison of pressure-time curves between numerical simulation and experiment

下载: 全尺寸图片幻灯片

整个膛内点传火过程从点火管火药点燃开始，高温高压燃气经由点火管小孔流入膛内，点燃主装药。燃气在点火管燃烧破膜以及在膛内流动需要一定的时间，而压力传感器测点在管壁处，所以初始阶段压力基本没有变化，从点火开始到1 ms时间段内，均保持在一个大气压左右；在随后的1～2 ms时间段内，压力快速上升，这是由于假药粒装填几乎占满整个药室，药室透气性很差，药室自由空间很小，膛内压力上升也较快；压力会在到达峰值时，会产生波动现象，这种波动现象主要由于假药粒装填密度过高，装药床透气性差，且假药粒不会燃烧导致燃气流动性变差，局部压力会产生振荡。

3.3 膛内两相流动过程分析

通过高速摄像拍摄可视化实验装置内火药燃气流动过程，与计算得到的压力场进行对比，可以清晰反映火药燃气在膛内的发展状况。由于可视化装置只能观察到膛内部分区域火焰传播过程，而火焰反映了高温燃气的流动过程，进而反映膛内的压力场的变化过程，因此对计算得到的压力云图进行相对应位置截取对比，如图4所示。

图 4 数值模拟压力云图与高速摄影图像的对比

Figure 4. Comparison between pressure clouds by numerical simulation and images by high-speed photography

下载: 全尺寸图片幻灯片

整个点传火可分为3个阶段；初始破孔阶段，在压力云图中可以看出随着点火管依次破孔，小孔破孔位置出现不同程度压力梯度，而在高速摄影中，由于点火管被假药粒包围，很难直接观测到初始阶段的破孔现象；随着点火管高温燃气不断向膛内喷射，高温燃气在膛内向膛壁和药室两端传播，在压力云图中可以观察到其压力梯度沿径向和轴向都有明显的变化，相对应的高速摄影图像中，火焰传播也呈现与压力云图同样的特征；到达2.5 ms时，点火管对应的膛内区域燃气流动在径向上基本一致，即沿径向的二维效应基本消失，而在远离点火管的区域仍存在一定的径向效应，从点火管开孔处向膛底和弹底形成明显的轴向压力梯度。

为有效呈现膛内压力随时间的发展变化，在径向上取r=20 mm处观察压力沿轴向的变化，如图5所示。图5(a)展示点火破膜初期阶段（0.5～1.2 ms）膛内压力沿轴向的变化，这一阶段膛内径向效应明显，压力从破孔区域不断上升，且随着点火管依次破孔，压力出现了相对应的峰值，并逐渐向周围发展。随着点火管燃气不断向膛内喷射，点火管对应的膛内区域压力基本趋于稳定，如图5(b)所示：2.5 ms后中心位置压力约18 MPa，径向压力梯度基本消失，此时在轴向压力梯度作用下，燃气向药室两端传播，膛底和弹底的压力不断升高，沿轴向压力变化区域也不断增大。

图 5 r=20 mm处沿轴向方向的压力分布

Figure 5. Pressure distributions along the axis at r=20 mm

下载: 全尺寸图片幻灯片

图6和图7描述了r=20 mm处膛内气相的径向速度和轴向速度沿药室轴向的变化规律。从图6可以看出，点火管破孔初始阶段，燃气不断从点火管内流向膛内，此时点火管附近区域燃气的径向速度要远大于其轴向速度，沿径向传播也比较快。初始破孔阶段，由于膛内压力远低于管内压力，气体喷射出来速度极大（由于压力差导致由点火管流入膛内气相速度很大），随着燃气不断流入，径向压力差变小，径向的燃气速度也逐渐降低。

图 6 0.5～1.0 ms不同时刻r=20 mm处沿轴向方向的气相速度分布

Figure 6. Gas phase velocity distributions along the axis at r=20 mm and different times from 0.5 ms to 1.0 ms

下载: 全尺寸图片幻灯片

图 7 1.5～3.2 ms不同时刻r=20 mm处气相速度沿轴向的分布

Figure 7. Gas phase velocity distributions along the axis at r=20 mm and different times from 1.5 ms to 3.2 ms

下载: 全尺寸图片幻灯片

在点火1.5 ms后，大量燃气由点火管破孔处向两端流动，并逐渐堆积形成高压区，在1.5～2.0 ms内，燃气逐渐向膛底和弹底运动，气相轴向速度渐渐升高，如图7(a)所示；待燃气运动到膛底或弹底时遇到阻力后又气相轴向速度慢慢降低至反向运动，如图7(b)所示，由于高密度装填，在弹底形成了药粒壅塞，所以气相速度在接近弹底处已开始有反向的波动；随着燃气轴向运动加剧，燃气的径向速度在点火管对应区域基本为零，在2.0 ms前燃气向膛底和弹底方向传播，由于燃气压力较低，初始传播时受颗粒扰动较大，径向速度波动剧烈，如图7(c)所示；2.2 ms后，由于燃气充盈膛底和中部区域，燃气径向速度除在膛底有小范围波动外，其他区域燃气径向速度降为零，燃气的径向扰动逐渐向弹底传播，如图7(d)所示。

火药颗粒运动对初始阶段膛内点传火研究很重要，为进一步研究药床中药粒运动规律，图8～9展示了r=20 mm膛内固相速度沿药室轴向发展变化。从图8可以看出，破孔初始阶段火药颗粒轴向和径向速度非常低，最高速度不超过0.5 m/s，说明膛内的火药颗粒基本保持不动，呈现一种多孔介质状态。从实验实际拍摄的高速影像来看，固相颗粒确实未发生大范围运动。

图 8 0.5～1.2 ms不同时刻r=20 mm处固相速度沿轴向的分布

Figure 8. Solid phase velocity distributions along the axis at r=20 mm and different times from 0.5 ms to 1.2 ms

下载: 全尺寸图片幻灯片

图 9 2.5～3.2 ms不同时刻r=20 mm处固相速度沿轴向的分布

Figure 9. Solid phase velocity distributions along the axis at r=20 mm and different times from 2.5 ms to 3.2 ms

下载: 全尺寸图片幻灯片

随着火药燃气在药床内不断流动，势必会带动火药颗粒运动。在高密度装填的药床中，火药颗粒可运动的范围十分有限，且模拟和实验采用的假药粒不会燃烧且变形较小，所以颗粒间相互作用力十分大，导致药床内的固相速度非常低。如图9(a)所示，火药颗粒径向速度随时间发展逐渐为零，从破孔位置沿轴向至弹底方向有一定的径向速度，至膛底处径向速度也逐渐趋向于零，与径向气相速度相对应，说明气相速度也是影响颗粒运动的重要原因。图9(b)中火药颗粒轴向速度也验证了这一点，由于燃气轴向压力梯度大，气相轴向速度高，所以导致火药颗粒向膛底和弹底两侧运动速度较高。在真实装药情况下，会造成膛底火药粒挤压破碎，形成小范围高压区域，继而推动药床向弹底运动。

4. 结　论

以火炮膛内初期点传火性能研究为背景，建立了轴对称二维两相流数理模型，对假药粒进行相应合理化假设，在空间上采用五阶WENO格式进行空间离散，时间上采用三阶TVD型龙格库塔法进行时间推进，并进行实验与数值模拟对比研究。数值分析了点火初期膛内压力的发展及气固两相速度，具体结论如下。

（1）可视化实验装置与高精度数值方法，有效模拟了初始阶段点火燃气在药床中的传播规律，通过与可视化实验平台测试数据和现象进行对比，验证了数理模型和数值方法的有效性。

（2）点火初期膛内径向效应明显，随着传火不断进行，点火管对应膛内区域径向效应逐渐消失，燃气的轴向传播开始占据主导地位，膛内压力呈现径向一致、轴向梯度分布的特征。

（3）通过数值模拟与实验观测可以得出，在高密度装填情况下，初始阶段药床药粒基本不会发生运动，最高速度为0.5 m/s，后续固相速度的发展与气相速度有紧密的关系；气相速度在初始阶段沿径向传播较快，后续在轴向压力梯度作用下，轴向速度开始占据主导，在膛底和中部区域气相速度降低为零，径向速度扰动也逐渐向弹底传播，在到达弹底之前由于药粒壅塞现象导致阻力变大，气相速度会出现振荡。

图 1 装药结构示意图

Figure 1. Schematic of the charge structure

下载: 全尺寸图片幻灯片

图 2 可视化实验装置

Figure 2. Visual experimental system

下载: 全尺寸图片幻灯片

图 3 计算压力曲线与实验对比图

Figure 3. Comparison of pressure-time curves between numerical simulation and experiment

下载: 全尺寸图片幻灯片

图 4 数值模拟压力云图与高速摄影图像的对比

Figure 4. Comparison between pressure clouds by numerical simulation and images by high-speed photography

下载: 全尺寸图片幻灯片

图 5 r=20 mm处沿轴向方向的压力分布

Figure 5. Pressure distributions along the axis at r=20 mm

下载: 全尺寸图片幻灯片

图 6 0.5～1.0 ms不同时刻r=20 mm处沿轴向方向的气相速度分布

Figure 6. Gas phase velocity distributions along the axis at r=20 mm and different times from 0.5 ms to 1.0 ms

下载: 全尺寸图片幻灯片

图 7 1.5～3.2 ms不同时刻r=20 mm处气相速度沿轴向的分布

Figure 7. Gas phase velocity distributions along the axis at r=20 mm and different times from 1.5 ms to 3.2 ms

下载: 全尺寸图片幻灯片

图 8 0.5～1.2 ms不同时刻r=20 mm处固相速度沿轴向的分布

Figure 8. Solid phase velocity distributions along the axis at r=20 mm and different times from 0.5 ms to 1.2 ms

下载: 全尺寸图片幻灯片

图 9 2.5～3.2 ms不同时刻r=20 mm处固相速度沿轴向的分布

Figure 9. Solid phase velocity distributions along the axis at r=20 mm and different times from 2.5 ms to 3.2 ms

下载: 全尺寸图片幻灯片

参考文献(14)

[1]	翁春生, 王浩. 计算内弹道学[M]. 北京: 国防工业出版社, 2006: 113−123.
[2]	袁亚雄, 张小兵. 高温高压多相流体动力学基础[M]. 哈尔滨: 哈尔滨工业大学出版社, 2005: 140−203.
[3]	王升晨, 周彦煌, 刘千里, 等. 膛内多相燃烧理论及应用[M]. 北京: 兵器工业出版社, 1994: 160−197.
[4]	王浩, 梁世超, 张莺, 等. 火焰在传火管装药床中的传输特性研究 [J]. 爆炸与冲击, 1999, 19(1): 66–71. WANG H, LIANG S C, ZHANG Y, et al. Study of flame propagation characteristics in fire transfering tube [J]. Explosion and Shock Waves, 1999, 19(1): 66–71.
[5]	王珊珊, 张玉成, 王浩, 等. 大长径比点火管高密实火药床点传火过程两相流的数值模拟 [J]. 爆炸与冲击, 2013, 33(4): 444–448. DOI: 10.11883/1001-1455(2013)04-0444-05. WANG S S, ZHANG Y C, WANG H, et al. Two-phase flow in ignition process of consolidated charge bed within a large length-to-diameter ratio igniter tube [J]. Explosion and Shock Waves, 2013, 33(4): 444–448. DOI: 10.11883/1001-1455(2013)04-0444-05.
[6]	王珊珊, 王浩, 黄明, 等. 药床透气性对大长径比点传火管点传火性能的影响研究 [J]. 弹道学报, 2013, 25(4): 95–99. DOI: 10.3969/j.issn.1004-499X.2013.04.019. WANG S S, WANG H, HUANG M, et al. Research on influence of charge permeability on ignition characteristics on large length-diameter ratio igniter tube [J]. Journal of Ballistics, 2013, 25(4): 95–99. DOI: 10.3969/j.issn.1004-499X.2013.04.019.
[7]	郭梦婷, 陶如意, 李子杰. 大长径比中心炸管式抛撒定容阶段两相流模拟 [J]. 弹道学报, 2017, 29(1): 68–72; 92. DOI: 10.3969/j.issn.1004-499X.2017.01.012. GUO M T, TAO R Y, LI Z J. Simulation of two-phase flow in the constant volume of dispersal system using centralizing blast-tube with large length-diameter ratio [J]. Journal of Ballistics, 2017, 29(1): 68–72; 92. DOI: 10.3969/j.issn.1004-499X.2017.01.012.
[8]	MIURA H, MATSUO A, NAKAMURA Y. Three-dimensional simulation of pressure fluctuation in a granular solid propellant chamber within an ignition stage [J]. Propellants, Explosives, Pyrotechnics, 2011, 36(3): 259–267. DOI: 10.1002/prep.201000058.
[9]	JARAMAZ S, MICKOVIĆ D, ELEK P. Two-phase flows in gun barrel: theoretical and experimental studies [J]. International Journal of Multiphase Flow, 2011, 37(5): 475–487. DOI: 10.1016/j.ijmultiphaseflow.2011.01.003.
[10]	程诚, 张小兵. 某制导炮弹二维两相流内弹道性能分析与数值模拟研究 [J]. 兵工学报, 2015, 36(1): 58–63. DOI: 10.3969/j.issn.1000-1093.2015.01.009. CHENG C, ZHANG X B. Two-dimensional numerical simulation on two-phase flow interior ballistic performance of a guided projectile [J]. Acta Armamentarii, 2015, 36(1): 58–63. DOI: 10.3969/j.issn.1000-1093.2015.01.009.
[11]	程诚, 张小兵. 高阶近似黎曼解模型在火炮内弹道两相流中的应用研究 [J]. 兵工学报, 2010, 32(10): 1200–1205. CHENG C, ZHANG X B. Research and application of higher-order approximate Riemann solver to two-phase flow in gun [J]. Acta Armamentarii, 2010, 32(10): 1200–1205.
[12]	张德良. 计算流体力学教程[M]. 北京: 高等教育出版社, 2010: 431 439.
[13]	李新亮. 高超声速湍流直接数值模拟技术 [J]. 航空学报, 2015, 36(1): 147–158. DOI: 10.7527/S1000-6893.2014.0233. LI X L. Direct numerical simulation techniques for hypersonic turbulent flows [J]. Acta Aeronautica et Astronautica Sinica, 2015, 36(1): 147–158. DOI: 10.7527/S1000-6893.2014.0233.
[14]	JIANG G S, SHU C W. Efficient implementation of weighted ENO schemes [J]. Journal of Computational Physics, 1996, 126(1): 202–228. DOI: 10.1006/jcph.1996.0130.

施引文献

期刊类型引用(2)

1.	马天一，陶如意，周思渊，王金龙，王浩. 模块装药短管炮实验及数值模拟研究. 弹道学报. 2024(02): 44-52 . 百度学术
2.	马天一，陶如意，森思义，王浩. 单模块装药膛内流场特性试验及数值仿真. 含能材料. 2023(01): 70-82 . 百度学术