小发动机羽流的试验和数值模拟研究张建华+ 田辉蔡国飙陈爱国杨彦广( 北京航空航天大学宇航学院，北京1 0 0 0 8 3 )( 中国空气动力研究与发展中心超高速空气动力研究所四川绵阳6 2 2 6 6 3 )摘要：试验研究和数值模拟了小发动机的羽流场，对数值计算和试验测量结果进行了对比，得出了一致的结果．试验在低密度风洞中进行测量．试验气体选为氮气，对两个不同几何形状的小发动机，共五种试验状态进行了测量，测出了小发动机的流场压力等参数．比较测量结果得出了符合实际工作情况的规律．同时对每种试验状态进行了数值模拟．模拟时采取N s 方程和D S M C 方法耦合求解．对喷管内部流动采用差分求解N s 方程的方法进行数值计算，将喷管出口截面的流动参数作为羽流场计算的进口边界条件．对羽流场采用直接模拟的蒙特卡罗方法( D S M C ) 进行数值模拟．考虑了环境背压对流场参数的影响．最后将计算结果和试验值进行了对比，得出了一致性的结论．本文进一步验证了D S M C 方法能够模拟低密度稀薄流动并保证了数值上的准确性，为工程应用尤其是8 0 k m 高度左右的发动机的流场特性提供了有力证据．主题词；小发动机；羽流；数值模拟；试验；N S 方程；D S M C 方法1 引言人造卫星、飞船、空间站以及深空探测器上大量使用的姿态控制发动机等推进系统，在高空稀薄环境下，姿态控制发动机工作时，其喷流向外部环境自由膨胀形成真空羽流，会对卫星、飞船、空间站以及深空探测器产生羽流污染、羽流干扰力和羽流热效应等羽流影响，这是进行卫星、飞船、空间站、深空探测器设计时必须解决的重要问题。
对丁：高空羽流问题的研究，从6 0 年代以来就各受人们的关注各国均开展了相关研究早期主要进行了大量的高空羽流试验研究⋯”2 ，随着计算条件的日益改善，数值模拟的方法成了研究这一问题的有力工具““”同时对于羽流问题的研究，也向这羽流碰撞、气动力和气动热研究方面发展”1 针对羽流问题，本文所开展的工作主要是：利用现有条件对小发动机进行了低密度条件下的试验研究，并针对试验状态进行了数值模拟验证一方面得到了试验结果规律，另一方面验证了数值模拟方法的可行性和准确性文中设计了两个几何形状的小发动机，对试验设备、试验过程和试验仪器进行了介绍对试验结果进行了分析试验在高超声速低密度风洞设备上进行工质选为氮气数值模拟时，采取计算喷管内部流动和外部流动耦合计算的方法来进行采用差分求解N s 方程的方法对喷管内部流动进行数值计算，得到喷管内部流动参数的变化规律及喷管出口截面的流动参数，为进行羽流场计算提供进口边界条件——2 4 9 ——在喷管外部，由于羽流为稀薄流，采用数值模拟稀薄气体的D S M C 方法进行羽流场的计算数值模拟时候考虑了环境压力对流场参数的影响针对试验状态，对同样条件的羽流场进行数值计算。
最后的数值模拟结果和试验测量值的对比表明非常好的一致性，验证了程序的可靠性本文的研究工作对深入开展羽流污染、气动力和气动热效应的研究打下了良好基础进一步的工作将开展对卫星上应用的实际发动机进行流场测量和气动力、气动热的试验研究和数值模拟研究2 试验测量2 ．1试验设备和测试仪器2 ．1 _ 1 试验设备试验是在高超声速低密度风洞上进彳丁= 的该风洞有气源系统、调压系统、前室、喷管、试验段、扩压段、冷却段、真空系统和测试系统构成试验时利用了该风洞的试验段，羽流的气源由自由射流装置提供小发动机的参数如下，如表1 所示：表1 喷管参数N o z z l e lN o z z l e 2N o m i n a lt h r u s tf N2 ．0 ．5I n l e td i a m e t e rD pm m1 0 ．1 0 ．T h r o a td i a m e t e rD f ，m Ⅲ1 ．20 ．6E x i td i a m e t e rD m m8 ．54 ．7 5I n l e th a l f a n g l e 巩d e g2 3 ．7 54 21 lE x i th a J 二a n g l e 如d e g1 51 52 ．1 ．2 测试仪器a ．试验段环境压力用Z R - 3 J 电阻真空计测量；b ．羽流前室压力用C Y G 一0 2 型( 量程为O 一1 M K P a ) 差压式固态压阻传感器测量：C ．皮托压力用C Y G 一1 9 型( 量程为0 , , , 4 0 K P a ) 差压式固态压阻传感器测量。
2 ．2 试验内容、方法和试验状态2 ．2 ．1 试验内容( 1 ) 羽流流场参数测量，测出皮托管的压力，然后确定流场的马赫数M 分布和当地静压；( 2 ) 辉光放电显示流场[ 1 1 02 ．2 ．2 试验方法用L Z - 1 5 真空泵将试验段的环境压力抽至1 0 P a 左右时．打开快速阔，调节调压阀，使发动机喷嘴的前室压力为O ．5 M P a 或1 M P a 时，气流通过喷管加速，在真空环境里膨胀形成羽流用辉光放电显示羽流流场结构用平头总压探头测量羽流流场压力张建华，博士生，北京航空航天大学2 - 6 5 信箱，1 0 0 0 8 3 ，0 1 0 - 8 2 3 1 7 2 0 6 ，z j h @ b u a a ．e d uc T I——2 5 0 ——2 ．2 ．3 试验状态试验状态如表2 所示表2 试验状态状志喷管总压，M P a总温．KC a 5 e l1I ．O3 0 0C a s e 21053 0 0C a s e 32l03 0 0C a s c 420 ．53 0 0C a s c 520 ．1 l3 0 02 ．3 试验结果与分析2 ．3 ．1 羽流流场结构图1 为羽流流场辉光显示结果。
其中( a ) 为F = 2 N 喷管在总压为1 ．0 M P a 时候的显示结果，( b ) 为F = 0 ．5 N 喷管在总压为1 ．0 M P a 时的显示结果从图1 看出：羽流流场结构较为简单，状似羽毛N 2 从喷嘴喷出，在低压力环境里膨胀，周边形成抛物面形的膨胀波构成羽流边界 a ) 试验状态1( b ) 试验状态3图1 羽流试验辉光显示图2 ．3 ．2 羽流流场压力测量结果在流场中测量羽流轴线不同位置处的压力皮托压力的测量结果见图2 图2 试验测量皮托压力一2 5 I 一2 ．4 数据处理流场中马赫数由前室总压R 和皮托压力只：求得，具体公式为： 弘《笔肛筹，一击r 蒜，击㈨此公式求解方法较多超声速流场中皮托管测压公式为：等划{ 篆2 ．／C ，| V ／a 毪r∞ n{1 十卜‘8I其中：K ‘= 生笋，下标l 为波前；2 为波后；k 为比热比为试验时候皮托管测定的总压通过公式( 2 ) ，由皮托压力可以求山当地静压3 数值计算进行计算时，对喷管内部流动采用差分求解N S 方程的方法进行数值计算，将喷管山口截面的流动参数作为羽流场计算的进口边界条件对羽流场采用直接模拟的蒙特卡罗方法( D S M C ) 进行数值模拟。
3 ．1 赜管内流场计算3 ．1 ．1 控制方程及坐标变换在忽略了彻体力、热辐射、且认为喷管为一绝热体系的前提条件下，轴对称、可压、冻结流、湍流流动雷诺平均的N a v i e r - - S t o k e s 方程形式为： 旦+ 堡+ 堡+ H ：0( 3 )a岔丹进一步把物理平面的控制方程变换到计算平面上， 署㈢+ 磬D 删小知D 玑删+ ㈢= ④{ ，r ／为计算坐标系下的坐标坐标变换式为：善：= 也，r ／，= 一％孝，= 一止q ，r ／，= d z #( 5 )其中J 为雅可比矩阵：，= 善：叩，一善，印：= ————= _ 一z { ～一z 口-3 ．1 ．2 湍流模型( 6 )由于湍流的复杂性，人们提出了许多种湍流模型，例如：0 方程模型、1 方程模型、2 方程模型等等，但其应用范围都是有限的在所有的湍流模型中，最简单的就是代数湍流模型在本研究的计算中，就选用这种模型米确定湍流粘性系数“的值——2 5 2 —较为著名的代数湍流模型有：C e b e c i 和s m i t h 的C —s 模型，B a l d w i n 和L o m a x [ f ib —L 模型。
1 等本研究在计算中选用B - - L 两层代数湍流模型，以给出壁面区的涡粘性系数3 ．1 I3 网格划分本文采用代数方法生成计算网格坐标变换公式如下：』善= 石( 7 ) I 叩= y ／f ( x )，' ( J ) 是喷管型面函数 为了反映壁面和喉部参数的剧烈变化程度，在壁面和喉部分别进行了网格加密，壁面加密体现在7 7方向，喉部加密体现在善方向a壁面加密公式：喉部加密公式％小鲁J l - J‰吐卜刽1 ≤／≤^( 8 )l 蔓f ≤f ，( 9 )岛出旷斗一箸{㈡⋯1 0 )L，式中：f ．是} 方向的网格总数，i ：是{ 方向喉部之前的网格数，尾．、卢f 2 分别是舌方向喉部前后的加密囡子，J l 是叩方向的网格总数，岛是，7 方向加密因子，工是喷管长度，L t 是喉部之前的长度·本文采用M a c C o r m a c k 时间推进预报校正二步格式在计算平面求解方程预报步和校正步分别采用向前差分和向后差分，计算平面网格间距取△善= A t ／= l ，所有差分格式中上标n 表示时问步数，下标‘' ，表示网格点差分格式如下：预测步：． ( 盼( 盼觚．水望笋] 0 ( 型笋] ：，J也( ( 哨0 ( 华H 圳Ⅲ，——2 5 3 —校正步：( ∥州，咄“[ ( 等≈一[ 型笋) ：，] 咄。
( 竽) H 华] “剐∞，( 等) j ：1 = 三f ( 等) i + ( 等) ：：『圹]c 切3 ．1 ．5 边界条件本文所进行模拟的流场是一个轴对称的喷管内流场，由于对称性，只需求解半个计算域的流场即可在本研究中，需给出如下边界条什：喷管内场的进I ：1 边界条件，出口边界条什，对称轴边界条件，I 古| 体壁面边界条件3 ．2 羽流场计算由于喷管外部是低密区域，从喷管喷射出来的气流很快会膨胀，因此采用求解稀薄气体的直接模拟的蒙特卡罗方法来进行数值模拟3 ．2 ．1D S M C 方法基本思想D S M C 方法的基本思想是：利用少量的模拟分子代替真实流场内数目众多的气体分子，用计算模拟实际物理过程，即对模拟分子进行跟踪，记录它们的位置和速度的交化，计算分子之间，分子和物面问的碰撞，再经统计平均，获得所需的各种流动参数在真实气体流动中，气体分子的运动与碰撞总是在同一时间进行的，它们相互耦合，相互影响在D S M C 方法中采用了将气体分子运动与碰撞解耦的方法，认为气体分子的碰撞是在瞬问完成的，它不改变分子的运动轨迹，而在相邻的两次碰撞之间，气体分子作匀速直线运动在计算模拟分子的碰撞时，首先随机选取两个可能的碰撞分子，然后计算该分子对的碰撞几率函数，这里我们采用分子对的相对运动速度作为几率函数。
如果几率函数大于某个随机函数，则认为该分子对将发生碰撞，于是计算碰撞后分子的速度如果儿率函数小于随机数，则认为该分子对将不会碰撞，需要重新选取分子对3 ．2 ．2D S M C 方法不同于传统C F D 方法的特点D S M C 方法不同于基于连续介质的传统C F D 方法的最主要特点在于：C F D 问题基于连续介质，其物理过程的基本特征是物理作用的邻接性，任何物理作用直接来自相邻的介质，而D S M C 基于离散介质，其物理过程的基本特征是分子的运动而且两者处理的对象也不同，C F D 计算的是网格点参数．而D S M C计算的网格内分子的运动参数D S M C 方法不同于基于连续介质的传统C F D 方法的特点还有：①D S M C 方法具有随机性②D S M C 方法不受边界形状的影响，可以求解具有任何复杂形状边界的物理问题③D S M C 方法不存在收敛与发散的问题。