微型固体姿控发动机微喷管内气粒两相流动规律的CFD-DSMC研究①

夏广庆，张 斌，孙得川，陈茂林，3

(1.大连理工大学工业装备结构分析国家重点实验室，大连 116024;2.中国航天科工集团公司第九总体设计部，武汉 430040;3.西北工业大学燃烧、热结构与内流场重点实验室，西安 710072)

微型固体姿控发动机微喷管内气粒两相流动规律的CFD-DSMC研究①

夏广庆1，张 斌2，孙得川1，陈茂林1，3

(1.大连理工大学工业装备结构分析国家重点实验室，大连 116024;2.中国航天科工集团公司第九总体设计部，武汉 430040;3.西北工业大学燃烧、热结构与内流场重点实验室，西安 710072)

微型固体姿控发动机在航天领域具有广泛的应用前景。以基于MEMS技术的微喷管为研究对象，首先通过计算微喷管中的克努森数，得到了微喷管中的气相流动状态;然后，采用CFD-DSMC方法，模拟了微喷管中的气粒两相流动，并研究了颗粒相质量分数和粒径对气相流动的影响。结果表明，在所研究的来流条件下，微喷管中的连续介质假设是成立的;气相与颗粒相间的动量和能量交换，导致气相马赫数降低、温度升高，同时也导致颗粒相速度增加、温度降低;颗粒相质量分数和粒径均能显著影响气相的马赫数和温度。

微型固体姿控发动机;微喷管;CFD;DSMC

0 引言

固体姿控发动机结构简单、可靠性高，在有效提高飞行器机动性能的同时，还能缩短控制系统的响应时间。因此，在国内外已得到了广泛的研究和使用［1-2］，其典型代表如使用180个固体姿控发动机作为末端控制装置的美国增程拦截弹ERINT-1。微型固体姿控发动机除具备以上优势外，以其在小型航天器上所具有的广泛应用前景，更成为当前的一个研究热点，其代表如基于MEMS技术的固体微推力器及其阵列［3-6］。

微型固体姿控发动机微喷管中的气粒两相流动一方面引发了比冲损失，同时又加剧了微喷管壁面的烧蚀。另一方面，对微喷管内气粒两相流动规律的研究是开展微型固体姿控发动机气粒两相羽流污染、气动力和气动热效应研究的首要前提［7-9］。由此可见，开展微型固体姿控发动机微喷管内气粒两相流动规律的研究工作，不仅对其优化设计具有十分重要的理论价值，同时对微型航天器的结构优化也具有十分重要的现实意义。

本文中，以基于MEMS技术加工的微喷管为研究对象。首先，开展了微喷管中纯气相流动规律研究，在得到气相主要特性参数分布规律的同时，通过计算流场中的克努森数，验证了气相计算方法的有效性;然后，采用CFD-DSMC方法，开展了微喷管中的气粒两相流动规律研究，得到了气相和颗粒相主要特征参数的分布规律，并研究了颗粒相浓度和粒径对气相相关参数的影响。

1 数学模型

气相控制方程组为不计体积力和无内热源的三维非定常N-S方程组。在笛卡尔坐标系下，其积分形式由式(1)表示为

颗粒相的控制方程组由式(4)～式(6)表示。

式中Vcell为网格单元体积;下标k表示第k束颗粒;Nk为第k束颗粒的真实颗粒数;Dk为第k束颗粒的直径;为第k束颗粒的滑移速度;CDk为第k束颗粒的滞止系数;Cp和C分别为颗粒和气体的比热容;Nuk为第k束颗粒的努谢尔特数;Pr为普朗特数。

2 数值方法

气相控制方程组采用中心有限体积法求解。其中，对流项采用添加各向异性人工粘性的中心差分格式求解，扩散项采用中心差分格式求解，时间项采用简化五步Runge-Kutta法求解，并采用隐式残值光顺法加速收敛。

对颗粒相方程组的求解，也采用时间推进法。在每一个时刻步开始时，类似经典的DSMC方法，在入口边界处加入新的颗粒，其速度和温度与气相一致，位置随机给定。在该时刻步内，计算颗粒的位置、速度和温度，并计算颗粒相对气相作用的动量和能量源项，同时对超出计算区域的颗粒施加相应的边界条件。在一个时刻步结束后，对颗粒的位置进行重新编号。在流场稳定后，统计流场信息。

对气相和颗粒相之间的耦合计算，采用PSIC算法。计算对象为方形截面的硅质微喷管，采用湿法化学腐蚀技术加工(腐蚀夹角约为54.7°)，其结构示意图如图1所示(单位μm)。计算中，采用结构网格离散计算区域，其示意图如图2所示。

图1 微喷管几何结构示意图Fig.1 Structure diagram of micro nozzle

图2 计算网格示意图Fig.2 Grid diagram for computation

气相入口给定来流的总温、总压、马赫数和方向角，气相出口处各流动参数均由内场按二价外插获得，气相壁面采用绝热、无滑移固壁边界条件。

颗粒相入口给定颗粒的质量分数、温度、速度和相关物性参数(密度为 4 004.6 kg/m3，比热容为1 380 J/(kg·K))，颗粒相出口处施加吸收边界，颗粒相壁面边界按Tabakoff［11］的经验公式，确定颗粒与壁面碰撞后的速度大小和方向。

3 计算结果与分析

3.1 微喷管中纯气相流动规律研究

在来流总温2 982 K、总压2 MPa条件下，微喷管中气相主要特性参数的分布规律如图3所示。由图3可知，由于微喷管加工工艺的限制，微喷管扩张段中壁面附近气流的膨胀并不充分，尤其是边角处气相的能量转换效率最低，导致气相的马赫数较低、温度较高。在微喷管的出口截面上，气相亚声速区域面积与该截面总面积的比值在扩张段内所有截面中位居最高，约为微喷管出口截面总面积的1/6。同时，气相在微喷管出口截面边角处的温度也在该截面中最高，其值约为2 600 K。图4为来流总温2 982 K、总压分别为1、2、3 MPa 3种情况下微喷管中的克努森数分布情况。由图4可知，来流总压越大，喷管中的克努森数越小。克努森数较大的区域位于喷管出口处的壁面附近。当来流总压1 MPa时，其最大值约0.008 5，微喷管中的流动为符合连续介质假设的无滑移流。可见，以上3种工况下，在全场采用无滑移边界的N-S方程进行流场求解是完全合理的。

图3 微喷管中气相主要特性参数分布规律Fig.3 Distribution contours of gas-phase characteristic parameters in the micro nozzle

图4 不同来流总压下微喷管中流动状态Fig.4 Flow state in the micro nozzle under different total pressures

3.2 微喷管中气粒两相流动规律研究

3.2.1 气相和颗粒相主要特性参数的分布规律

当颗粒质量分数为9.4%、粒径为1 μm时，气相主要特性参数的分布规律如图5所示(气相来流条件同3.1节)。由图5可知，颗粒相的存在对气相压强场和密度场的影响不大，但对气相温度场和马赫数场的影响较大，尤其是在扩张段中，在颗粒浓度较大的位置，气相场的温度升高，马赫数降低。与3.1节中纯气相流动计算结果对比可知，在微喷管出口截面的中心线处，颗粒相的存在导致气相在该处的马赫数降低约6%，温度升高约8%。

图5 气相主要特性参数分布规律Fig.5 Distribution contours of gas-phase characteristic parameters

颗粒相主要特性参数的分布规律如图6所示。由图6可知，沿流动方向颗粒相的速度逐渐增大，温度逐渐降低。在喷管的扩张段中，大部分区域为无粒子区(图6中参数无梯度变化的区域为无粒子区)。

图6 颗粒相主要特性参数分布规律Fig.6 Distribution contours of particle-phase characteristic parameters

3.2.2 颗粒质量分数对气相马赫数和温度的影响

在不同颗粒质量分数下，气相流场中心线处的马赫数和温度变化分别如图7和图8所示。由图7和图8可知，随颗粒质量分数的增加，颗粒对气相运动的阻碍效应增强，导致气相的马赫数呈下降趋势。同时，由于颗粒质量分数的增加，也直接导致了颗粒相与气相间换热效应的增强，从而导致气相温度呈上升趋势。

3.2.3 颗粒粒径对气相马赫数和温度的影响

在不同颗粒粒径下，气相流场中心线处的马赫数和温度变化分别如图9和图10所示。由图9和图10可知，随颗粒粒径的增大，颗粒相与气相接触的表面积减小，对气相的阻力下降，从而导致中心线上气相的马赫数和速度均呈上升趋势。同时，相间接触面积的减小，导致了颗粒相与气相间换热效应的减弱，从而导致气相温度呈下降趋势。

图7 不同颗粒质量分数下微喷管中心线上气相马赫数变化Fig.7 Mach number distribution of gas-phase along the micro nozzle centerline under different particle mass fractions

图8 不同颗粒质量分数下微喷管中心线上气相温度变化Fig.8 Temperature distribution of gas-phase along the micro nozzle centerline under different particle mass fractions

图9 不同颗粒粒径下微喷管中心线上气相马赫数变化Fig.9 Mach number distribution of gas-phase along the micro nozzle centerline under different particle diameters

图10 不同颗粒粒径下微喷管中心线上气相温度变化Fig.10 Temperature distribution of gas-phase along the micro nozzle centerline under different particle diameters

4 结论

(1)在气相来流总温2 982 K、总压1 MPa的条件下，微喷管中气相克努森数的最大值约为0.008 5，气相流动为符合连续介质假设的无滑移流。

(2)颗粒相的存在对气相温度场和马赫数场影响较大，尤其是在气相特性参数大幅变化的扩张段中，颗粒相的存在导致气相马赫数降低、温度升高。在气相来流总温2 982 K、总压2 MPa、颗粒质量分数9.4%、粒径1 μm的条件下，颗粒相的存在导致微喷管出口截面中心线处的气相马赫数降低约6%，温度升高约8%。

(3)随来流颗粒质量分数的增加，气相流场中心线处的马赫数和速度呈下降趋势，温度呈上升趋势。随着来流颗粒粒径的增大，气相流场中心线处的马赫数和速度均呈上升趋势，而温度呈下降趋势。

［1］邢继发，刘国球，黄坚定，等.世界导弹与航天发动机大全［M］.北京:军事科学出版社，1999.

［2］张钢锤，叶定友.一种新概念固体姿控发动机简介［J］.固体火箭技术，2001，24(3).

［3］Muller J，Marrese C，Polk J，et al.An overview of MEMS-based micropropulsion developments at JPL［J］.Acta Astronautica，2003，52(9).

［4］Rossi C，Orieux S，Larangot B，et al.Design fabrication and modeling of solid propellant microrocket-application to micropropulsion［J］.Sensors and Actuators，2002，99(1).

［5］周海清，张高飞，尤政.固体微推力器设计与数值分析［J］.推进技术，2007，28(3).

［6］刘旭辉，方蜀州.微型固体推力器阵列寻址点火控制系统研究［J］.固体火箭技术，2010，33(6).

［7］张建华，贺碧蛟，蔡国飙，等.卫星姿控发动机喷管羽流撞击效应试验［J］.空气动力学报，2007，25(2).

［8］Park J H，Baek S W，Kim J S.DSMC analysis of the interaction between thruster plume and satellite components［R］.AIAA 2002-0794.

［9］程晓丽，王强，阎喜勤.卫星姿控发动机高空羽流场工程分析［J］.空间科学学报，2006，26(1).

［10］陈伟芳，常雨.三维管道超声速气固两相流动的CFD/DSMC 仿真［J］.推进技术，2005，26(3).

［11］Tabakoff W，Sugiyama Y.Experimental method of determining particle restitution coefficients［C］//New York:Proc.Symp.on Polyphase Flow and Transport，ASME，1980，8:2032210.

Research on the gas-particle two-phase flow in the micro nozzle of attitude control micro solid rocket motor

XIA Guang-qing1，ZHANG Bin2，SUN De-chuan1，CHEN Mao-lin1，3
(1.State Key Laboratory of Structural Analysis for Industrial Equipment，Dalian University of Technology，Dalian 116024，China;2.The Ninth System Design Department of China Aerospace Science and Industry Corporation，Wuhan 430040，China;3.Science and Technology on Combustion，Internal Flow and Thermal-Structure Laboratory，Northwestern Polytechnical University，Xi＇an 710072，China)

Attitude control micro solid rocket motor has wide application potential in the aerospace field.The gas-particle twophase flow in the micro nozzle based on the MEMS technology was investigated.Firstly，through calculating the Knudsen number of the micro nozzle，the gas phase flow state in the micro nozzle was obtained.Then the gas-particle two-phase flow in the micro nozzle was simulated by using the method of CFD-DSMC.The influence of particle mass fraction and particle diameter on the gas phase flow was studied.The result shows that the continuum assumption in the micro nozzle is established under the conditions of the defined flow in the study.The exchange of momentum and energy between the gas phase and the particle can not noly reduce the gas phase Mach number and raise the temperature，but also increase the particle phase velocity and decrease the temperature.The particle phase mass fraction and particle diameter can significantly influence the Mach number and temperature of gas phase.

attitude control micro solid rocket motor;micro nozzle;CFD;DSMC

V435

A

1006-2793(2012)03-0356-05

2012-01-09;

2012-04-08。

国家自然科学基金(11105023);新世纪优秀人才支持计划(NCET-11-0054);工业装备结构分析国家重点实验室开放基金(GZ1101);中央高校基本科研业务费专项资金(DUT11ZD(G)02)。

夏广庆(1979—)，男，讲师，研究方向为发动机总体设计。E-mail:gq.xia@dlut.edu.cn

(编辑:崔贤彬)