0.5 N模型发动机羽流撞击效应的D S MC模拟

王 倩 何小英 蔡国飙

（北京航空航天大学宇航学院）

1 引言

航天器的姿态控制往往利用姿态控制发动机来进行调节，发动机喷出的羽流会对流场中的航天器带来羽流污染、气动力、气动热效应等影响。羽流污染可能会导致航天器上敏感原件无法正常工作，甚至影响航天器的寿命，所以研究羽流的撞击效应就显得非常重要。

国外对羽流撞击效应的研究开展较早，1986年德国的Dettleff试验研究了模型发动机羽流场[1]；1984年法国的Allegre，在SR3风洞中测量了羽流撞击平板的气动力效应[2]；1990年Legge测量了羽流撞击倾斜挡板的气动力[3]；Boyd等人在1992年对低密度喷管流动和羽流场进行了试验研究和数值模拟[4]。在实验研究的基础上，国外还进行了相关数值模拟研究，1985年Legge等人对卫星模型发动机的羽流场和羽流撞击效应进行了数值模拟[5]；2002年Park数值模拟了羽流和卫星的相互作用[6]。国内由于试验条件的限制，研究工作主要集中在数值计算方面。

本文采用自开发的PWS软件，对推力0.5N的模型发动机进行了羽流撞击效应DSMC模拟。工质为CO2。首先，对流场的压强、温度、马赫数、流线图等进行研究；然后，对挡板压强及羽流污染情况进行了分析。

2 DSMC方法及PWS软件介绍

DSMC方法的基本思想是：利用少量的模拟分子代替真实流场内数目众多的气体分子，用计算模拟实际物理过程，即对模拟分子进行跟踪，记录它们的位置和速度的变化，计算分子之间，分子和物面间的碰撞，再经统计平均，获得所需的各种流动参数。DSMC的基本假设是：模拟每个网格中的分子间的碰撞时，忽略分子的具体位置。

基于DSMC的PWS软件的主程序包括网络模块、粒子参数模块、粒子进入模块、粒子运动及边界模块、粒子碰撞模块、并行模块、统计输出模块共七大模块。在这七大基本模块的基础上，针对各种实际问题，PWS羽流计算软件都可以方便的搭建数值模拟应用平台[7]。

PWS软件是根据模块化和面向对象的原则对DSMC计算程序进行改进，具体如下：

（1）使用松散的软件架构保证各个计算模块的相对独立性，即改进其中一个模块时，基本上不需要改动其他模块；

（2）对各个模块的功能实现采用通用性设计，即在使用单/多组分分子、不同分子模型、不同固体边界条件时不需要改动源代码。

使用图1的步骤生成所需边界条件，这样可以将粒子-固体边界处理归结到几种基本类型，可以方便的实现通用性设计。

（1）定义几种常用的平面和曲面（直线）构成基本元素；

（2）使用多个平面或曲面相交构成凸多面体；

（3）使用多个凸多面体来构成最终的复杂边界。

图1 边界生成过程说明

3 计算模型

本文采用推力0.5N模型发动机，边长300mm方形挡板，工质为CO2，总压Pc=8000Pa，总温Tc=900K，环境温度为300K。

3.1 物理模型

发动机与挡板空间位置如图2所示，设定X轴为对称轴。发动机轴线经过挡板中心，a=0.1m，角度ω=0°，n=m=300mm，h=0.1mm，Dt=0.6mm，De=4.7mm。

图2 发动机与挡板空间位置

3.2 边界条件

DSMC粒子入口条件：以0.5N发动机出口截面为粒子入口截面，模拟粒子选取CO2。真空边界：认为粒子逃逸，即粒子经过边界后注销。粒子入口截面参数分布如图3至图6所示。图3为粒子入口密度分布图，图4为粒子入口温度分布图，图5为粒子入口轴向速度分布图，图6为粒子入口径向速度分布图。

图3 粒子入口密度分布图

图4 粒子入口温度分布图

图5 粒子入口轴向速度分布图

图6 粒子入口径向速度分布图

3.3 网格划分

X，Y，Z方向网格数目分别为60,120,120，且网格宽度为等比数列，X方向比例因子为1.035，Y,Z方向网格宽度从中心向两侧等比增加，比例因子为1.046，Y-Z平面网格视图如图7所示。

图7 网格划分

4 计算结果及分析

使用PWS软件对羽流撞击效应进行计算。采用64台节点机并行计算约2.9h，总共计算模拟分子数达到950万个以上，迭代到10,000步时开始采样，继续统计20,000步得出羽流流场。流场参数分布如图8-图 11。

图8 流场压强分布图

图9 流场温度分布图

图10 流场马赫数分布图

由图8至图10可以看出，发动机羽流场的压强在粒子入口处最高，最大压强约为950Pa，在轴线和挡板附近一定区域内比较高，随着气流向真空中扩散而降低。在轴线附近一椭球形区域内，温度较低，马赫数较高；而在这个椭球形边缘的附近，由于从发动机喷出的气流和反射气流相遇，造成局部温度的升高，马赫数减小；在椭球形以外的流场中，则以反射气流为主导，压强减小，温度下降，马赫数增加。最高温度约为450K，最大马赫数约为12。根据图11可以看出，气流沿着发动机出口呈放射状流向挡板，喷流撞击到挡板上发生反射，挡板表面附近的气流向四周反射流走。

挡板压强分布如图12所示，由于核心气流对挡板的撞击，在发动机轴线与挡板的交点处压强较高，以该交点为圆心，板上的压强分布具有对称性，且沿径向逐渐减小。平板中心最高压强可达110Pa左右。由此可见，挡板中心处羽流效应影响最严重，影响程度由中心沿径向逐渐减小。

图12 挡板压强分布图

5 结论

本文采用自开发的PWS软件，对推力0.5N的模型发动机进行了羽流撞击效应进行DSMS模拟，工质为CO2，得到了流场的压强、温度、马赫数、流线图；然后，对挡板压强和羽流效应情况进行了分析，得到以下结论：

（1）发动机羽流场压强在粒子入口处最高，在轴线和挡板附近一定区域内比较高，随着气流向真空中扩散而降低。

（2）在轴线附近一个椭球形区域内，温度较低，马赫数较高；而在这个椭球形边缘的附近，由于从发动机喷出的气流和反射气流相遇，造成局部温度的升高，马赫数减小；在椭球形以外的流场中，则以反射气流为主导，压强减小，温度下降，马赫数增加。

（3）气流沿着发动机出口呈放射状流向挡板，喷流撞击到挡板上发生反射，挡板表面附近的气流向四周反射流走。

（4）由于核心气流对挡板的撞击，在发动机轴线与挡板的交点处压强较高，以该交点为圆心，板上的压强分布具有对称性，且沿径向逐渐减小。挡板中心处羽流效应影响最严重，影响程度由中心沿径向逐渐减小。 ◇

［1］ DETTLEFF G,BOETTCHER R D,DANKERT C,KPPPENWALLNER G,LEGGE H.Attitude Control Thruster Plume Flow Modeling and Experiments［J］.Journal of Spacecraft and Rockets,1986,23（5）:476-481.

［2］ALLEGRE J,RAFFIN M,LENGRAND J C.Forces Induced by a Simulated Rocket Exhaust Plume Impinging upon a Flat Plate［A］.Proceedings of the 14th International Symposium on Rarefied Gas Dynamics［C］,1984:287-294.

［3］LEGGE H.Plume impingement forces on Inclined Flat Plates［A］.Proceedings of the 17th International Symposium on Rarefied Gas Dynamics［C］,1990:955-962.

［4］BOYD I D,PENKO P F,MEISSNER D L,DEWITT K J.Experimental and Numerical Investigations of Low-density Nozzle and Plume Flows of Nitrogen ［J］.AIAA Journal,1992,30（10）:2453-2461.

［5］LEGGE H,BOETTCHER R D.Modelling Control Thruster Plume Flow and Impingement［A］.Proceedings of the 13thInternational Symposium on Rarefied Gas Dynamics［C］,1985:983-992.

［6］PARK J H,BAEK S W,KIM J S.DSMC Analysis of the Interaction between Thruster Plume and Satellite Components［R］.AIAA 2002-0794.

［7］蔡国飙，贺碧蛟，PWS软件应用于探月着陆器羽流效应数值模拟研究［A］.航天器环境工程［J］.2010:18-23.