某型液体火箭发动机喷注参数与尾焰红外图像特征研究

申帅帅，张 楠，王大锐

0 引 言

发动机尾焰红外图像的研究主要在以下两个方面：a）利用尾焰红外辐射能量大、红外图像突出的特点，来识别目标飞行器；b）对尾焰红外图像本身特征的研究。胡炳梁等[1]成功研制了中国第1个超声速机载红外测量吊舱系统，为测量目标机尾焰红外特性提供了有效的测量手段；史丽芳等[2]完成了飞机尾焰红外图像识别软件系统，该系统利用低空卫星携带的红外热像仪来探测目标，并根据各个阶段尾焰呈现的不同特性（形体特征、温度）来检测与识别来袭飞机，及时预警并主动攻击，具有较高的辨识率和较好的发展前景；张硕等[3]开发了红外图像识别系统，将小波分析、中值滤波以及数学形态学的方法、密度分割法和频率阈法应用于固体火箭发动机排气羽焰图像的识别处理中，得到了较好的火焰轮廓及其基本参数，如火焰量值、长宽比、均值对比度、亮度最大像素点、复杂度等参数等，具有较好的识别效果。

尾焰红外图像中包含丰富的燃烧信息，较好地反映发动机工作状态。但是，目前通过建立尾焰红外图像特征与发动机工作状态之间的关系对发动机工作状态进行监测诊断的研究还比较少。由于室压、壁温等常规监测参数不对发动机燃烧过程中所有状态变化都敏感，而尾焰流场的亮度、特定区域的形状等特征参数对于不同的燃烧工况反应明显[4]，且尾焰红外图像的测量方法具有非接触、无干扰、精度高等特点，因此对发动机尾焰光学测量监测的应用研究越来越广泛。

本文采用多学科仿真技术，耦合蒸发、燃烧、传热过程，建立了三维非稳态推力室蒸发-燃烧-尾焰场一体化仿真模型；采用尾焰图像识别方法，开展液体火箭发动机喷注参数与尾焰红外图像特征之间关系的研究。

1 物理模型和计算方法

1.1 物理模型及边界条件

本文仿真对象为某型液体火箭发动机，燃料组合为甲基肼（MMH）和绿色四氧化二氮（NTO）。该仿真模型采用三维非稳态仿真模型，耦合了蒸发、燃烧、传热等多个过程。

三维仿真计算模型网格划分如图 1所示，整个计算域采用多层O形网格划分，壁面附近划分边界层网格，网格精度较高。计算域由流体区域Fluid块和固体区域Solid块组成，内壁面采用Coupled壁面条件，外壁面采用固定热流壁面条件，尾焰场区域边界采用Pressure-Out边界条件。

图1 三维仿真计算模型示意Fig.1 Schematic Diagram of Three-dimensional Simulation Model

1.2 计算模型

该推力室-尾焰场仿真模型的计算采用Euler-Lagrange方法，控制方程包括 Euler坐标系下的N-S方程和 Lagrange方程坐标系下的离散颗粒轨道方程。气液两相之间的耦合通过源项实现。将液滴简化为质点，用 Lagrange法进行追踪，液滴状态由其中心位置XP、直径DP、速度UP和密度ρp4个物理量确定；对于气相流体采用Euler坐标系下的N-S方程进行描述。

为提高计算精度，在计算允许的条件下，应尽量增加喷注器数量，本次仿真中布置了 352个喷注器，共计1958束液滴流束。喷注器布置面垂直于推力室轴线方向，其示意如图2所示。

图2 喷注器布置示意Fig.2 Schematic Diagram of Injector Layout

各个喷注器的位置、液滴喷射速度、粒径、流量均由相位多普勒粒子动态分析仪（Phase Doppler Particle Analyzer，PDPA）雾化试验测量获得。由于仿真过程中采用的工质为MMH和NTO，而雾化试验采用的工质是水，所以布置喷注器时，液滴的喷注参数需要进行换算。

蒸发模型中采用的颗粒类型为液滴，遵循加热、冷却、蒸发、沸腾规律，当液滴剩余质量大于挥发极限质量，且温度分别达到汽化和沸腾温度时，液滴便会以相应的速率发生蒸发和沸腾，液滴表面质量输运方程为

式中 mp(t)为 t时刻液滴质量；mp(t+Δt)为 t+Δt时刻液滴质量；Ni为蒸发速率系数，与饱和蒸汽压、蒸汽扩散系数、运动粘度系数有关；Ap(t)为t时刻液滴表面积。

化学反应采用MMH和NTO一步总包反应，化学反应模型为Finite-Rate/Eddy-Dissipation燃烧模型[5]，假定化学反应速率由化学动力学模型和涡耗散模型共同控制，取两模型各自计算反应速率的较小值作为化学反应速率。其中，化学动力学控制的反应速率采用阿累尼乌斯（Arrhenius）公式计算，涡耗散模型假定化学反应速率由涡团破碎的速率计算。

MMH和NTO液体推进剂的物性参数（密度、比热容、潜热、饱和蒸汽压等）采用《液体推进剂》[6]提供的数据；壁面材料的铌钨合金物性参数（热导率、比热容等）采用《航天用铌钨合金棒材规范》[7]提供的数据；各气相物质的物性参数，如热导率、粘度等为常数，定压比热容为关于温度的多项式。

1.3 尾焰红外图像识别方法

本文采用液体火箭发动机尾焰红外图像特征识别方法[4]作为尾焰场特征识别的工具，该方法具体流程如图3所示。

基本步骤如下：

a）对图像进行灰度化处理，并设定合适阈值获取需要的灰度图像，过滤掉部分干扰点后得到边缘点图像；

b）将边缘点连接成为线段，过滤掉干扰线段，并为线段排序编号；

c）将各边缘线段和构造的区域边界连接起来，构成指定区域（尾焰出口区、尾焰射流区、尾焰核心区）；

d）计算各区域周长、面积、灰度熵，第1马赫盘位置等特征参数值。

图3 尾焰红外图像特征识别方法流程Fig.3 Schematic Diagram of the Process of Infrared Flame Image Identification

通过红外图像特征识别方法识别出2个区域，如图4所示，区域1为尾焰出口区，区域2为尾焰射流区。

图4 仿真尾焰红外图像识别分区示意Fig.4 Schematic Diagram of Simulation of Tail Flame Infrared Image Recognition

将2个区域的周长、面积、灰度熵、第1马赫盘位置作为尾焰场红外图像的特征指标。

假定喷管出口直径为D ，图像识别出其像素个数为n，则单个像素尺寸为D/n。

a）区域周长Li：识别区域的边缘像素点个数与单个像素尺寸的乘积：

式中 n1i为横或者竖方向线段数量；n2i为对角线方向的线段数；下标i为识别区域代号，i=1，2。

b）区域面积Si：识别区域的像素点个数与单个像素尺寸的平方的乘积：

式中 Ni为区域内的像素点数，下标i为识别区域代号，i=1，2。

c）灰度熵：灰度熵是图像分割中常用的特征参数，它表示图形中像素灰度的不均匀程度或者复杂程度。灰度熵越小，表明图像均匀性越好[8]。

对于一幅256级灰度图像而言，每个像素的灰度值可能是0～255级灰度中的某一级，因此像素灰度值的符号集为{0，1，2，…，255}。假设图像中灰度值为j的像素个数为kj，则图像中灰度值为j的像素出现的概率为

整幅图像的灰度熵值为

式中 I为灰度。

d）第1马赫盘位置：如图5所示，喷管出口射流经强烈压缩形成的第 1道正激波所在的位置即为第 1马赫盘位置。用第1马赫盘所在位置与喷口之间的像素点数与单个像素尺寸的乘积表示，计算表达式如下：

式中 y1为第1马赫盘所在位置的像素坐标；y0为喷管出口所在位置的像素坐标；y1-y0为第 1马赫盘距离喷管出口的像素个数。

图5 喷管出口附近波系示意Fig.5 Schematic Diagram of the Waves Construction out of the Nozzle

2 仿真和热试车结果对比分析

2.1 仿真和热试车室温、室压、壁温、尾焰场对照

多对互击式发动机仿真和热试车外壁面、尾焰场红外图像特性如图6所示。

图6 仿真和热试车外壁面、尾焰场红外特性云图Fig.6 Simulation and Hot Firing Test of Outer Wall Surface and Tail Flame Infrared Image

续图6

由图6可知：图6a和图6b的仿真结果较为一致，仿真尾焰场出现明显的马赫盘结构，且第1马赫盘位置清晰；图6c和图6d尾焰轮廓相差较大，分析可知：由于图6d中的尾焰图像是发生复燃反应的发动机的试车数据，所以仿真结果尾焰边缘未出现外扩现象是正常的，并发现仿真和热试车尾焰场温度波动幅度不一样，这是因为不同温度下流场的辐射系数不同，而图6d是根据红外图像对全尾焰场取同一辐射系数求解所得温度场，导致结果有偏差。

多对互击式发动机推力室内温度、压力、马赫数特性如图7所示。

图7 推力室中轴面温度、压力、马赫数仿真云图Fig.7 Simulation of Axial Plane in the Thrust Chamber Temperature, Pressure and Mach Number Cloud

由于离散化喷嘴布置在喷注面板处，液滴汽化吸热，使该处温度较低；推力室直线段压力变化不大，再经过喉部和喷管扩张段，燃气压力逐渐减小，马赫数逐渐增大，且在喉部附近达到声速，与试车结果和气体动力学理论符合较好。仿真与热试车结果对照如表1所示。

表1 仿真、热试车结果对照Tab.1 Comparison Table of Simulation and Hot Firing Test

如表1可知，仿真和热试车的室压和外壁最高温度偏差均小于5%，该仿真模型结果较为准确。

仿真和热试车外壁面温度沿推力室轴向无量纲位置变化规律如图8所示。

图8 仿真和热试车外壁面温度沿轴向变化规律Fig.8 Figure of Simulation and Hot Firing Test of Outer Wall Surface Temperature Changes Along the Axis

由图 8可知，在距离喷注面板一定距离内（区域1），由于边区液膜、气膜冷却作用，壁面温度较低；距离喷注面板一定距离（区域2）时，气膜消失，高温燃气直接冲刷内壁，壁面温度迅速升高；随着距离的增加（区域3），流动截面积变小，热流密度增加，壁面温度继续上升（一直持续到喉部前端）；随着轴向位置的进一步增加（区域4），由于流动截面积变小，燃气内能转换为动能加剧，气流温度降低，壁面温度也随之降低；过了喉部（区域5），流动截面积变大，热流密度变小，而且燃气内能进一步转化为动能，燃气温度迅速降低，壁面温度随之急剧降低。

对比分析仿真结果可知：仿真和热试车外壁面温度变化曲线相吻合，验证了推力室蒸发-燃烧模型的计算精度可用该型发动机的仿真分析。

2.2 仿真和热试车尾焰识别参数对照

为了验证该仿真模型可以应用上述尾焰红外图像特征识别方法进行后续研究，将仿真和热试车尾焰红外图像识别结果进行对比分析，如表2所示。

表2 额定工况下仿真和热试车尾焰特征参数对照Tab.2 Comparison Table of Simulation and Hot Firing Test the Tail Flame Characteristic Parameters

由表2可知，除了灰度熵，尾焰图像的其他特征，第1马赫盘位置、区域周长和面积的偏差均小于10%，由于热试车采用的是红外热像仪拍摄的尾焰红外图像，而仿真结果采用的是尾焰场温度云图近似代替尾焰红外图像，其颜色梯度上有差别，导致灰度熵有所偏差。综上所述，该尾焰红外图像特征识别方法对仿真尾焰红外图像的识别结果准确可靠。

3 喷注参数与尾焰特征关系

采用上述仿真模型和尾焰红外图像特征识别方法研究了不同的喷注参数（液滴喷射速度、液滴平均粒径、边区冷却流量百分比）条件下尾焰红外图像特征的变化规律。

3.1 液滴平均喷射速度与尾焰特征关系

通过对不同液滴喷射速度下尾焰场红外图像进行特征识别，得到不同喷注速度下第1马赫盘位置和区域1、2灰度熵随液滴平均喷注速度变化图（见图9）、区域1、2周长随液滴平均喷注速度变化图（见图10）、区域1、2面积随液滴平均喷注速度变化图（见图11）。

由图9～11可知，随着液滴喷射速度由20 m/s变化到40 m/s，尾焰场红外图像中第1马赫盘位置、尾焰出口区和尾焰射流区的灰度熵和面积均出现明显变化，且出现极值，而尾焰出口区和尾焰射流区的周长没有明显变化，说明这两个特征参数对于液滴喷射速度变化不敏感。

图9 第1马赫盘位置和区域灰度熵随液滴平均喷注速度变化Fig.9 Figure of the Area of the First Mach Locus and Gray Value Entropy Variation with Droplet Ejection Average Speed

图10 区域1、区域2周长随液滴平均喷注速度变化Fig.10 Figure of the Zone 1, 2 Perimeter Variation with the Average Jetting Speed

3.2 液滴平均粒径与尾焰特征关系的研究

通过对不同的液滴平均粒径下尾焰场红外图像进行特征识别，得到不同液滴平均粒径下第1马赫盘位置和区域1、区域2灰度熵随液滴平均粒径变化图（见图12）、区域1、区域2周长随液滴平均粒径变化图（见图13）、区域1、区域2面积随液滴平均粒径变化图（见图 14）。

由图12～14可知，随着液滴平均粒径由15 μm变化到35 μm，尾焰场红外图像中第1马赫盘位置、尾焰出口区和尾焰射流区的灰度熵和面积均出现明显变化。在平均粒径为25 μm附近，尾焰出口区和尾焰射流区灰度熵均出现极小值，尾焰出口区、尾焰射流区面积和第1马赫盘位置均出现极大值，而尾焰出口区和尾焰射流区的周长没有明显变化，说明这两个特征参数对于液滴平均粒径变化不敏感。

图12 第1马赫盘位置和区域灰度熵随液滴平均粒径变化Fig.12 Figure of the Area of the First Mach Locus and Gray Value Entropy Variation with Droplet Average Particle Size

图13 区域1、区域2周长随液滴平均粒径变化Fig.13 Figure of the Zone 1, 2 Perimeter Variation with roplet Average Particle Size

图14 区域1、区域2面积随液滴平均粒径变化Fig.14 Figure of the Zone 1, 2 Area Variation with Droplet Average Particle Size

3.3 边区冷却流量百分比与尾焰特征关系

通过对不同的边区冷却流量百分比下尾焰场红外图像进行特征识别，得到不同边区冷却流量百分比下第1马赫盘位置和区域1、区域2灰度熵随边区冷却流量百分比变化图（见图15）、区域1、区域2周长随边区冷却流量百分比变化图（见图16）、区域1、区域2面积随边区冷却流量百分比变化图（见图17）。

图15 第1马赫盘位置和区域灰度熵随边区冷却流量百分比变化Fig.15 Figure of the Area of the First Mach Locus and Gray Value Entropy Variation with Side to Cool the Flow Percentage

图16 区域1、区域2周长随边区冷却流量百分比变化Fig.16 Figure of the Zone 1, 2 Perimeter Variation with Side to Cool the Flow Percentage

图17 区域1、区域2面积随边区冷却流量百分比变化Fig.17 Figure of the Zone 1, 2 Area Variation with Side to Cool the Flow Percentage

由图15～17可知，随着边区冷却流量百分比由22%变化到30%，尾焰场红外图像中第1马赫盘位置、尾焰出口区和尾焰射流区的灰度熵和面积均出现明显变化，且第1马赫盘位置、尾焰出口区和尾焰射流区的灰度熵出现双峰变化规律，尾焰出口区和尾焰射流区的面积变化明显，但没有发现特定的变化规律；尾焰出口区和尾焰射流区的周长没有明显变化，说明这两个特征参数对于边区冷却流量百分比变化不敏感。

观察还发现，随着这3个喷注参数的变化，区域灰度熵和第1马赫盘距离喷管出口的位置均呈现出很强的负相关性。分析知，随着燃烧效率增大，喷管出口马赫数增大，同时，燃烧更充分，尾焰红外图像灰度熵减小。说明喷管出口马赫数与第1马赫盘距离喷管出口的距离成正相关，与文献[9]中结论一致。

4 结 论

采用蒸发、燃烧、传热多学科算法和尾焰场红外图像识别方法，利用雾化试验液滴分布规律，建立的三维非稳态推力室蒸发-燃烧-尾焰场一体化仿真模型，对某型发动机推力室喷嘴雾化特性和边区冷却流量与尾焰场红外图像特征的关系仿真分析结论如下：

a）本文建立的模型与热试车结果吻合较好，具有较高的计算准确性；

b）尾焰射流场的周长对喷注参数不敏感，而灰度熵、马赫盘位置和图像面积的敏感性较强，能够用来监测喷注特性；

c）采用马赫盘位置、尾焰出口区和尾焰射流区的灰度熵和面积峰值可以评估推力室喷嘴雾化参数和边区冷却流量的选取。

[1] 胡炳梁, 刘学斌，等. 超音速机载吊舱红外测量系统研究[J]. 红外与激光工程, 2003, 32(5)： 441-444.

Hu Bingliang, Liu Xuebin, et al. Research on an ultrasonic pod infrared measurement system[J]. Infrared and Laser Engineering, 2003, 32(5)：441-444.

[2] 史丽芳. 基于图像处理的飞机尾焰检测与识别技术的研究[D]. 西安：西安电子科技大学, 2015.

Shi Lifang. A research of detection and recognition technology of the plane tail flame based on image processing[D]. Xi’an： Xidian University,2015.

[3] 张硕, 王宁飞, 张平. SRM羽焰图像的检测与参数预估[J]. 导弹与航天运载技术, 2008 (4)： 47-50, 55.

Zhang Shuo, Wang Ningfei, Zhang Ping. Image detection and parameter estimation on flame images of solid rocket motor[J]. Missiles and Space Vehicles, 2008(4)： 47-50, 55.

[4] 王大锐, 张楠, 葛明和. 液体火箭发动机喷管燃气外流场红外图像研究[J]. 导弹与航天运载技术, 2016(2)： 26-30.

Wang Darui, Zhang Nan, Ge Minghe. Infrared image study on gas external flow field of liquid rocket engine nozzle[J]. Missiles and Space Vehicles,2016(2)： 26-30.

[5] 刘涛. 燃烧理论[M]. 长沙： 国防科学技术大学出版社, 2010.

Liu Tao. Combustion Theory[M]. Changsha： University of Defense Science and Technology University Press, 2010.

[6] 李亚裕. 液体推进剂[M]. 北京： 中国宇航出版社, 2011.

Li Yayu. Liquid propellant[M]. Beijing： China Aerospace Press, 2011.

[7] 周小军. GJB8507-2015. 航天用铌钨合金棒材规范[S]. 北京： 总装备部军标出版发行部, 2015.

Zhou Xiaojun. GJB8507-2015. Specification for niobium-tungsten alloy bars for aerospace[S]. Beijing： General Armament Department Military Standards Publication and Distribution Department, 2015.

[8] 田相军, 罗琳. 综合图像灰度熵和灰度值的人脸识别方法[J]. 现代电子技术, 2005, 29(24)： 46-48.

Tian Xiangjun, Luo Lin. Face recognition of grayhound of entropy and grayhound of value[J]. Modern Electronic Technique, 2005, 29(24)： 46-48.

[9] 胡俊, 邱明勇, 郭绍刚. 激光切割气体出口马赫数对马赫盘的影响[J].中国激光, 2008, 35(11)： 1789-1794.

Hu Jun, Qiu Mingyong, Guo Shaogang. Effect of gas outlet mach number on mach shock disk in laser cutting[J]. Chinese Journal of Lasers, 2008,35(11)： 1789-1794.