陳 杰，賀碧蛟，蔡國(guó)飆

（北京航空航天大學(xué)宇航學(xué)院，北京 100191）

火星環(huán)繞器羽流效應(yīng)仿真研究

陳 杰，賀碧蛟?，蔡國(guó)飆

（北京航空航天大學(xué)宇航學(xué)院，北京 100191）

火星探測(cè)器發(fā)動(dòng)機(jī)羽流效應(yīng)影響飛行任務(wù)的安全，采用差分求解NS方程與直接模擬蒙特卡洛耦合的方法對(duì)火星環(huán)繞器姿軌控發(fā)動(dòng)機(jī)羽流力、熱以及污染效應(yīng)進(jìn)行了仿真研究。結(jié)果表明：發(fā)動(dòng)機(jī)羽流對(duì)火星環(huán)繞器氣動(dòng)力、熱及污染效應(yīng)影響主要在發(fā)動(dòng)機(jī)儲(chǔ)箱上，最大熱流值可以達(dá)到1789 W／m2，需要注意；發(fā)動(dòng)機(jī)羽流對(duì)展開(kāi)太陽(yáng)能電池翼的氣動(dòng)力和氣動(dòng)力矩相對(duì)于發(fā)動(dòng)機(jī)推力和推力力矩來(lái)說(shuō)，小3～4個(gè)量級(jí)，影響較小；已分離著陸器在距發(fā)動(dòng)機(jī)出口300 m以內(nèi)時(shí)，羽流對(duì)已分離著陸器的影響較大，氣動(dòng)熱流在100 W／m2以上；在100 m以內(nèi)時(shí)，氣動(dòng)熱流在1000 W／m2以上，需要對(duì)敏感設(shè)備采取一定的防護(hù)措施。不同熱適應(yīng)系數(shù)氣動(dòng)熱流仿真對(duì)比計(jì)算結(jié)果表明，隨著熱適應(yīng)系數(shù)的增加，熱流密度值呈線性增加。數(shù)值模擬結(jié)果可以為火星環(huán)繞器工程設(shè)計(jì)提供參考。

火星探測(cè)；火星環(huán)繞器；羽流效應(yīng)；DSMC

1 引言

2020年我國(guó)將進(jìn)行首次火星探測(cè)，火星環(huán)繞器作為火星探測(cè)系統(tǒng)的一部分，將要完成地火轉(zhuǎn)移飛行和軌道修正、完成火星捕獲、進(jìn)入火星停泊軌道、與巡視器可靠分離、攜帶載荷完成環(huán)火遙感等科學(xué)探測(cè)任務(wù)。由于發(fā)動(dòng)機(jī)羽流會(huì)對(duì)火星環(huán)繞器及其上敏感設(shè)備、太陽(yáng)能帆板等產(chǎn)生氣動(dòng)力和力矩、氣動(dòng)熱以及污染效應(yīng)，對(duì)火星環(huán)繞器的正常工作造成影響。為了使火星環(huán)繞器能夠順利完成任務(wù)，需要對(duì)火星環(huán)繞器發(fā)動(dòng)機(jī)羽流效應(yīng)進(jìn)行仿真研究。

直接模擬蒙特卡洛（DSMC）方法［1］是目前模擬發(fā)動(dòng)機(jī)羽流效應(yīng)較為成熟的方法?；贒SMC方法，北京航空航天大學(xué)的蔡國(guó)飆等人開(kāi)發(fā)了PWS（Plume Work Station）軟件［2-3］，驗(yàn)證算例表明軟件具有較高的置信度［4］。PWS軟件已成功應(yīng)用在神舟飛船、嫦娥3號(hào)著陸器等航天器的羽流仿真研究中［5-6］。

本文采用差分求解NS（Navier-Stokes）方程與DSMC相耦合的方法［7］，使用PWS對(duì)不同工況下發(fā)動(dòng)機(jī)羽流對(duì)環(huán)繞器表面及敏感元件的氣動(dòng)力、熱和污染效應(yīng)進(jìn)行仿真分析，對(duì)比不同熱適應(yīng)系數(shù)對(duì)環(huán)繞器氣動(dòng)熱流仿真結(jié)果的影響。對(duì)于火星環(huán)繞器與火星著陸巡視器分離的過(guò)程，由于已釋放的著陸巡視器與環(huán)繞器的距離不斷改變，著陸巡視器與羽流作用是一個(gè)復(fù)雜的非穩(wěn)態(tài)過(guò)程。本文采用將羽流流場(chǎng)與已釋放著陸巡視器解耦的辦法對(duì)羽流對(duì)已釋放著陸巡視器的影響進(jìn)行估算。

2 仿真模型

PWS軟件采用二級(jí)構(gòu)造法來(lái)構(gòu)建復(fù)雜仿真模型［2］。首先在不影響羽流效應(yīng)仿真分析結(jié)論的前提下，對(duì)實(shí)際模型進(jìn)行適當(dāng)簡(jiǎn)化。然后對(duì)PWS軟件定義的幾種常用的平面（直線）和曲面（曲線）等基本元素，使用“與”和“或”運(yùn)算最終組合成擁有復(fù)雜邊界的仿真模型。圖1為使用PWS軟件建立的火星環(huán)繞器羽流仿真簡(jiǎn)化模型。

3 數(shù)值方法

考慮到發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)以及羽流核心區(qū)流場(chǎng)具有較高的壓強(qiáng)和密度，使用DSMC方法求解，計(jì)算量過(guò)大。本文采用差分求解NS方程與DSMC相耦合的方法，對(duì)火星環(huán)繞器發(fā)動(dòng)機(jī)羽流效應(yīng)進(jìn)行仿真研究。

對(duì)于火星環(huán)繞器與火星著陸巡視器分離的過(guò)程，本文采用解耦的辦法對(duì)羽流對(duì)已釋放著陸巡視器的影響進(jìn)行估算。即首先計(jì)算發(fā)動(dòng)機(jī)羽流場(chǎng)，得到羽流場(chǎng)壓強(qiáng)、密度、溫度、速度分布，然后根據(jù)著陸巡視器與環(huán)繞器距離和與發(fā)動(dòng)機(jī)推力夾角，確定火星探測(cè)器分離軌跡，并得到軌跡上羽流壓強(qiáng)、密度、溫度和速度參數(shù)，最后根據(jù)局地平衡假設(shè)及分子動(dòng)力學(xué)理論，計(jì)算分離軌跡上的羽流質(zhì)量流率、壓強(qiáng)和熱流密度，估算羽流效應(yīng)。

1）差分求解NS方程

目前計(jì)算連續(xù)流域流動(dòng)的最常用方法是數(shù)值求解NS方程的方法。矢量形式的NS方程如式（1）所示［8］，其中，Q為守恒變量，F(xiàn)、G為無(wú)粘通量，F(xiàn)v、Gv為粘性通量，S為源項(xiàng)。

一般真實(shí)的流場(chǎng)都是以湍流狀態(tài)存在的，因此，控制方程中的粘性項(xiàng)和熱傳導(dǎo)項(xiàng)中的系數(shù)需由層流和湍流共同確定［9］，如式（2）、式（3）：

式中，μ、μl和μt分別為粘性系數(shù)、層流和湍流粘性系數(shù)，Pr、Prl和Prt分別為普朗特?cái)?shù)、層流和湍流普朗特?cái)?shù)。層流粘性系數(shù)一般隨溫度而變化，可由Sutherland公式較為精確地給出［9］。在求解雷諾平均控制方程時(shí)，為了使其封閉，必須引入計(jì)算湍流粘性系數(shù)的湍流模型［8］，計(jì)算中采用SST k-ω兩方程模型。

2）DSMC方法

由于DSMC計(jì)算域的流場(chǎng)密度和氣體溫度較低，計(jì)算時(shí)做以下假設(shè)［10］：（1）假設(shè)流場(chǎng)中分子間的碰撞為二體碰撞；（2）僅考慮分子的轉(zhuǎn)動(dòng)內(nèi)能；（3）氣體流動(dòng)為定常流動(dòng)；（4）假設(shè)分子為可變徑硬球模型（VHS）。

對(duì)于分子與壁面的作用模型，目前應(yīng)用較為廣泛的是Maxwell反射模型。Maxwell反射模型采用熱適應(yīng)系數(shù)α來(lái)表示分子在與壁面作用的過(guò)程中，鏡面反射與漫反射各占的比例［10］，計(jì)算公式如式（4）：

其中，qi和qr是來(lái)流和反射流的能流量，qw是反射溫度為壁溫時(shí)完全Maxwell漫反射下的能流量。在完全Maxwell漫反射時(shí)，α＝1，完全鏡面反射時(shí)，α＝0。熱適應(yīng)系數(shù)表征了反射分子的溫度在多大程度上已經(jīng)“適應(yīng)”了壁面的溫度。

4 羽流效應(yīng)計(jì)算

火星環(huán)繞器在不同階段將執(zhí)行不同的飛行任務(wù)，需要不同數(shù)量、不同推力的發(fā)動(dòng)機(jī)組合工作。對(duì)每個(gè)發(fā)動(dòng)機(jī)工況，首先使用差分求解NS方程的方法得到發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)及部分外流場(chǎng)，再使用得到的發(fā)動(dòng)機(jī)出口截面參數(shù)作為發(fā)動(dòng)機(jī)二維自由羽流場(chǎng)DSMC計(jì)算的入口條件，然后從得到的二維流場(chǎng)計(jì)算結(jié)果中截取流場(chǎng)參數(shù)作為三維羽流場(chǎng)入口條件，對(duì)三維羽流場(chǎng)進(jìn)行仿真，得到不同工況下發(fā)動(dòng)機(jī)羽流對(duì)火星環(huán)繞器及敏感元件的氣動(dòng)力、熱及污染效應(yīng)。采用不同熱適應(yīng)系數(shù)對(duì)發(fā)動(dòng)機(jī)儲(chǔ)箱的氣動(dòng)熱流密度進(jìn)行了對(duì)比計(jì)算。

1）入口條件

對(duì)火星環(huán)繞器進(jìn)行全流場(chǎng)計(jì)算時(shí)，計(jì)算域較大，發(fā)動(dòng)機(jī)推力較大以及發(fā)動(dòng)機(jī)數(shù)量較多。在不影響計(jì)算結(jié)果的情況下，為減小計(jì)算量，將發(fā)動(dòng)機(jī)噴管向外延伸一個(gè)母線傾角為15°的圓錐，以此圓錐面上的流場(chǎng)參數(shù)作為三維DSMC計(jì)算入口條件。

圖2為分別以發(fā)動(dòng)機(jī)出口和發(fā)動(dòng)機(jī)噴管向外延伸一個(gè)母線傾角為15°的圓錐為DSMC入口條件得到的羽流場(chǎng)壓力對(duì)比云圖?？梢钥闯?，兩種入口條件得到的羽流場(chǎng)基本吻合，只有在離發(fā)動(dòng)機(jī)較遠(yuǎn)處、壓力較小的流場(chǎng)區(qū)域有些許偏差，且偏差在工程允許范圍內(nèi)。結(jié)果說(shuō)明，將發(fā)動(dòng)機(jī)噴管向外延伸一個(gè)母線傾角為15°的圓錐作為DSMC入口條件，對(duì)羽流流場(chǎng)計(jì)算結(jié)果影響不大，此做法較為合理。對(duì)120 N和25 N發(fā)動(dòng)機(jī)羽流場(chǎng)分析可以得到類似的結(jié)論。

為了得到此圓錐面上的流場(chǎng)參數(shù)，首先差分求解NS方程，得到發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)及部分外流場(chǎng)結(jié)果，然后以發(fā)動(dòng)機(jī)出口截面參數(shù)作為入口條件用DSMC方法對(duì)二維發(fā)動(dòng)機(jī)羽流場(chǎng)進(jìn)行計(jì)算，最后在得到的二維羽流場(chǎng)中截取圓錐面母線上的流場(chǎng)參數(shù)，作為三維羽流場(chǎng)DSMC計(jì)算入口條件。

圖3為從3000 N發(fā)動(dòng)機(jī)出口開(kāi)始的5 m×5 m的外流場(chǎng)壓力分布云圖，圖中紅線表示作為三維DSMC計(jì)算入口的圓錐面位置。圓錐面所在位置的粒子壓強(qiáng)和密度與發(fā)動(dòng)機(jī)出口處相比，已大大降低，達(dá)到了減小計(jì)算量的目的。

2）邊界條件

計(jì)算域邊界：當(dāng)粒子經(jīng)計(jì)算域邊界運(yùn)動(dòng)出計(jì)算域時(shí)，按逃逸處理，即將粒子刪除。

模型壁面：當(dāng)粒子撞擊到模型邊界時(shí)，將粒子運(yùn)動(dòng)按Maxwell反射處理。模型壁面的溫度設(shè)置為常溫300 K（模型壁面溫度的變化，對(duì)羽流仿真結(jié)果影響不大）［2］，發(fā)動(dòng)機(jī)壁面溫度設(shè)置為1000 K，熱適應(yīng)系數(shù)設(shè)定為1。在對(duì)比算例中，分別計(jì)算了熱適應(yīng)系數(shù)為0、0.25、0.5、0.75和1的對(duì)比算例。

入口條件邊界：將入口條件圓錐面設(shè)置為吸收邊界條件，即在入口條件邊界內(nèi)部不存在粒子，粒子運(yùn)動(dòng)穿過(guò)吸收邊界時(shí)即被刪除。

5 結(jié)果分析

5.1 差分求解NS方程仿真結(jié)果

首先對(duì)發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)及部分外流場(chǎng)使用差分求解NS方程的方法進(jìn)行仿真計(jì)算。

圖4為3000 N發(fā)動(dòng)機(jī)內(nèi)流場(chǎng)及部分外流場(chǎng)流場(chǎng)溫度分布云圖。可以看出，溫度在燃燒室內(nèi)最高，噴管擴(kuò)張段為特形面，氣體充分?jǐn)U張。高溫高速噴流噴入真空環(huán)境中，羽流邊界形成了拋物線形的膨脹波，噴管型面作用使噴流形成兩個(gè)連在一起的錐形激波。3000 N發(fā)動(dòng)機(jī)流場(chǎng)壓強(qiáng)、密度、Ma數(shù)流場(chǎng)分布云圖也呈現(xiàn)出類似的規(guī)律。

5.2 羽流氣動(dòng)力、熱以及污染效應(yīng)的影響

在差分求解NS結(jié)果的基礎(chǔ)上，使用DSMC方法對(duì)發(fā)動(dòng)機(jī)二維羽流場(chǎng)和三維羽流場(chǎng)進(jìn)行仿真計(jì)算。對(duì)所有工況仿真結(jié)果分析總結(jié)，發(fā)現(xiàn)在3000 N+4×120 N+2×25 N發(fā)動(dòng)機(jī)組合工作時(shí)，火星環(huán)繞器及其上設(shè)備所受到的羽流效應(yīng)影響最大。

圖5為該工況下的三維流場(chǎng)切面溫度分布圖?？梢钥闯觯诎l(fā)動(dòng)機(jī)儲(chǔ)箱上方流場(chǎng)溫度最高，這是由于3000 N與120 N發(fā)動(dòng)機(jī)噴流相互阻滯的緣故，這與發(fā)動(dòng)機(jī)儲(chǔ)箱上的熱流密度最大的仿真結(jié)果相一致。流場(chǎng)密度、壓強(qiáng)分布圖等也呈現(xiàn)出相似的規(guī)律。

圖6為該工況下環(huán)繞器底板及其上設(shè)備羽流影響分布云圖?？梢钥闯觯谠摴r下，發(fā)動(dòng)機(jī)羽流對(duì)發(fā)動(dòng)機(jī)儲(chǔ)箱的羽流效應(yīng)影響最大，儲(chǔ)箱上的最大熱流值可以達(dá)到1789 W／m2，CO2質(zhì)量到達(dá)率可以達(dá)到9.12×10-5kg／m2／s。對(duì)其它結(jié)果云圖分析發(fā)現(xiàn)，羽流氣動(dòng)力壓強(qiáng)、H2O質(zhì)量到達(dá)率等分布呈現(xiàn)出類似的規(guī)律，各最大值均在發(fā)動(dòng)機(jī)儲(chǔ)箱上；而對(duì)于離3000 N主發(fā)動(dòng)機(jī)距離相對(duì)較遠(yuǎn)的火星環(huán)繞器上的其它敏感元器件來(lái)說(shuō)，其附近的羽流流場(chǎng)密度壓強(qiáng)、溫度、密度等變得相對(duì)較小，對(duì)敏感元器件的羽流氣動(dòng)力、熱以及污染效應(yīng)也因此變得較小。因此，在該火星環(huán)繞器設(shè)計(jì)方案下，主要考慮發(fā)動(dòng)機(jī)羽流對(duì)發(fā)動(dòng)機(jī)儲(chǔ)箱所造成的氣動(dòng)力、熱及污染效應(yīng)。

5.3 太陽(yáng)能電池翼上羽流氣動(dòng)力和力矩的影響

由于太陽(yáng)能電池翼展開(kāi)半徑較大，需要考慮各工況下展開(kāi)太陽(yáng)能電池翼上的羽流氣動(dòng)力以及其對(duì)探測(cè)器中心點(diǎn)坐標(biāo)（0，0，0.705）的力矩值。圖7為0°、45°和90°的太陽(yáng)翼所受到的發(fā)動(dòng)機(jī)羽流氣動(dòng)力壓強(qiáng)分布云圖。

對(duì)不同角度的太陽(yáng)能電池翼上的力和力矩的計(jì)算結(jié)果總結(jié)分析，得到不同角度（°）下的太陽(yáng)能電池翼所受到的氣動(dòng)力在坐標(biāo)軸上的分力值如表1所示。

表1 不同角度太陽(yáng)能電池翼上羽流氣動(dòng)力Table 1 Aerodynamic force on solar wings with differrent angles

表2為工作的兩臺(tái)發(fā)動(dòng)機(jī)A1和A2的各自推力在坐標(biāo)軸上的分力值，以及兩臺(tái)發(fā)動(dòng)機(jī)推力合力在坐標(biāo)軸上的推力值。表3為不同角度下的太陽(yáng)能電池翼所受到的羽流氣動(dòng)力相對(duì)于參考點(diǎn)坐標(biāo)的力矩值。表4為A1、A2發(fā)動(dòng)機(jī)推力相對(duì)于參考點(diǎn)坐標(biāo)的力矩值以及它們的和推力相對(duì)于參考坐標(biāo)的力矩值。

表2 發(fā)動(dòng)機(jī)推力Table 2 Thrust of motors

表3 太陽(yáng)能電池翼上羽流氣動(dòng)力相對(duì)于參考點(diǎn)的力矩Table 3 Torque of aerodynamic force on solar wings relative to the reference point

表4 發(fā)動(dòng)機(jī)推力相對(duì)于參考點(diǎn)的力矩Table 4 Torque of motor thrust relative to the reference point

由表1和表2可以看出，發(fā)動(dòng)機(jī)羽流對(duì)太陽(yáng)能電池翼的羽流氣動(dòng)力在0.001 N量級(jí)，而發(fā)動(dòng)機(jī)的推力則在10 N量級(jí)，即與發(fā)動(dòng)機(jī)推力相比，太陽(yáng)能電池翼上的羽流氣動(dòng)力影響可以忽略。

由表3和表4可以看出，發(fā)動(dòng)機(jī)羽流對(duì)太陽(yáng)能電池翼的羽流氣動(dòng)力相對(duì)于參考點(diǎn)（0，0，0.705）的力矩在0.01 N·m量級(jí)，而發(fā)動(dòng)機(jī)推力相對(duì)于參考點(diǎn)的力矩值在10 N·m量級(jí)，與發(fā)動(dòng)機(jī)推力矩相比，羽流氣動(dòng)力所造成的力矩影響也較小。

圖8為太陽(yáng)能電池翼上的氣動(dòng)力隨著角度變化的變化趨勢(shì)圖。圖9為太陽(yáng)能電池翼上的氣動(dòng)力相對(duì)于參考點(diǎn)（0，0，0.705）的氣動(dòng)力矩隨著角度變化的變化趨勢(shì)圖。

由圖8可以看出，45°太陽(yáng)能電池翼所受到的x軸向的氣動(dòng)力最大，即45°太陽(yáng)翼受到的繞火星環(huán)繞器中心軸的扭轉(zhuǎn)力最大。隨著太陽(yáng)能電池翼偏轉(zhuǎn)角度的增加，太陽(yáng)能電池翼上所受到的沿z軸的氣動(dòng)力逐漸增大，即90°時(shí)最大，0°時(shí)最??；沿y軸的氣動(dòng)力先減小后增大，在45°時(shí)最小，90°時(shí)最大。由圖9可以看出，沿x軸方向的羽流氣動(dòng)力矩隨著太陽(yáng)能電池翼轉(zhuǎn)角的增大而增大，這是因?yàn)殡S著轉(zhuǎn)角的增大，太陽(yáng)能電池翼受到的沿x軸方向的力越來(lái)越大，而力臂基本相同；沿y軸方向的羽流氣動(dòng)力矩在各角度下都基本為0；沿z軸方向的力矩隨著太陽(yáng)能電池翼轉(zhuǎn)角的增大先減小后增大，在轉(zhuǎn)角為45°時(shí)絕對(duì)值最大，方向與z軸正向相反；合力矩隨著太陽(yáng)能電池翼轉(zhuǎn)角的增大而增大。

5.4 發(fā)動(dòng)機(jī)羽流對(duì)已分離著陸器的影響

采用解耦的方法估算不同推力發(fā)動(dòng)機(jī)羽流對(duì)距離發(fā)動(dòng)機(jī)出口不同距離的已分離火星著陸器的影響。圖10為距發(fā)動(dòng)機(jī)出口不同距離下，3000 N發(fā)動(dòng)機(jī)羽流場(chǎng)氣動(dòng)熱流分布對(duì)比圖?？梢钥闯觯谂c發(fā)動(dòng)機(jī)軸線夾角相同時(shí)，熱流密度隨著與發(fā)動(dòng)機(jī)出口距離的增加而減小；而在距離發(fā)動(dòng)機(jī)出口相同距離時(shí)，隨著與發(fā)動(dòng)機(jī)軸線夾角度數(shù)的增加，熱流密度也呈下降趨勢(shì)。仿真結(jié)果符合羽流場(chǎng)結(jié)構(gòu)特點(diǎn)，即隨著流場(chǎng)位置與發(fā)動(dòng)機(jī)軸線夾角的不斷增大，流場(chǎng)逐漸由核心流域過(guò)渡到稀薄流域，密度、壓強(qiáng)、溫度等流場(chǎng)參數(shù)值不斷降低；而在與發(fā)動(dòng)機(jī)軸線同一夾角下的位置上，噴流逐漸向真空膨脹，密度、壓強(qiáng)、溫度等流場(chǎng)參數(shù)值隨著噴流的膨脹而降低。對(duì)其它結(jié)果云圖分析發(fā)現(xiàn)，羽流場(chǎng)質(zhì)量通量分布以及正應(yīng)力分布也呈現(xiàn)出類似的規(guī)律。表5為3000 N發(fā)動(dòng)機(jī)工作時(shí)，在已分離著陸器距發(fā)動(dòng)機(jī)出口不同距離下的羽流對(duì)著陸器表面的最大質(zhì)量通量、最大正應(yīng)力和最大熱流密度值?？梢钥闯?，當(dāng)已分離著陸器距發(fā)動(dòng)機(jī)出口300 m時(shí)，質(zhì)量通量達(dá)到10-5kg／m2量級(jí)，正應(yīng)力達(dá)到0.01 Pa量級(jí)，熱流達(dá)到100 W／m2量級(jí)；當(dāng)已分離著陸器距發(fā)動(dòng)機(jī)出口100 m時(shí)，質(zhì)量通量達(dá)到10-4kg／m2量級(jí)，正應(yīng)力達(dá)到0.1 Pa量級(jí)，熱流達(dá)到1000 W／m2量級(jí)，即對(duì)于分離著陸器來(lái)說(shuō)，發(fā)動(dòng)機(jī)羽流對(duì)已分離著陸器的影響在300 m內(nèi)應(yīng)該考慮，在100 m內(nèi)應(yīng)該對(duì)相關(guān)敏感元器件采取一定的防護(hù)措施。

表5 3000 N發(fā)動(dòng)機(jī)各計(jì)算結(jié)果的最大值Table 5 Maximum value of plume effects of 3000 N motor

5.5 熱適應(yīng)系數(shù)對(duì)羽流氣動(dòng)熱流仿真結(jié)果的影響

圖11為不同熱適應(yīng)系數(shù)下發(fā)動(dòng)機(jī)儲(chǔ)箱上最大熱流密度值的變化趨勢(shì)圖?？梢钥闯?，隨著熱適應(yīng)系數(shù)的增加，發(fā)動(dòng)機(jī)儲(chǔ)箱上最大熱流密度值基本呈線性增加。

當(dāng)熱適應(yīng)系數(shù)為0時(shí)，發(fā)動(dòng)機(jī)儲(chǔ)箱上的熱流密度為0，即此時(shí)粒子與發(fā)動(dòng)機(jī)儲(chǔ)箱表面沒(méi)有發(fā)生能量交換，仿真結(jié)果符合理論實(shí)際。當(dāng)熱適應(yīng)系數(shù)由0.25增大到1時(shí)，發(fā)動(dòng)機(jī)儲(chǔ)箱上熱流密度分布基本相同，只是熱流密度值隨著熱適應(yīng)系數(shù)的增大而增大。發(fā)動(dòng)機(jī)儲(chǔ)箱上的熱流密度最大值由熱適應(yīng)系數(shù)α＝0.25時(shí)的464 W／m2，增大到熱適應(yīng)系數(shù)α＝1時(shí)的1789 W／m2。

6 結(jié)論

對(duì)火星環(huán)繞器在入軌及器器分離過(guò)程中的羽流效應(yīng)進(jìn)行了仿真分析；對(duì)同一發(fā)動(dòng)機(jī)工況下不同角度太陽(yáng)能電池翼羽流氣動(dòng)熱流進(jìn)行了對(duì)比計(jì)算；對(duì)熱適應(yīng)系數(shù)對(duì)羽流氣動(dòng)熱流仿真結(jié)果的影響進(jìn)行了研究，結(jié)果表明：

1）受羽流氣動(dòng)力、熱及污染影響較大的為發(fā)動(dòng)機(jī)儲(chǔ)箱，氣動(dòng)熱流值可以達(dá)到1789 W／m2；羽流對(duì)環(huán)繞器表面及其他敏感元件的氣動(dòng)力、熱及污染影響相對(duì)較小，氣動(dòng)熱流值基本都在100 W／m2量級(jí)和10 W／m2量級(jí)。

2）發(fā)動(dòng)機(jī)羽流對(duì)太陽(yáng)能電池翼的氣動(dòng)力以及氣動(dòng)力相對(duì)于參考點(diǎn)的力矩?cái)?shù)值較小，氣動(dòng)力在0.001 N量級(jí)；相對(duì)于參考點(diǎn)的氣動(dòng)力矩值在0.01 N·m量級(jí)，比發(fā)動(dòng)機(jī)推力和發(fā)動(dòng)機(jī)推力相對(duì)于參考點(diǎn)的力矩小3～4個(gè)量級(jí)；太陽(yáng)能電池翼上羽流氣動(dòng)力和力矩相對(duì)于發(fā)動(dòng)機(jī)合力和合力矩可以忽略。

3）對(duì)于分離著陸器來(lái)說(shuō)，發(fā)動(dòng)機(jī)羽流對(duì)距發(fā)動(dòng)機(jī)出口100 m的已分離著陸器的質(zhì)量通量達(dá)到10-4kg／m2量級(jí)，正應(yīng)力達(dá)到0.1 Pa量級(jí)，熱通量達(dá)到1000 W／m2量級(jí)；對(duì)距發(fā)動(dòng)機(jī)出口300 m的已分離著陸器質(zhì)量通量達(dá)到10-5kg／m2量級(jí)，正應(yīng)力達(dá)到0.01 Pa量級(jí)，熱通量達(dá)到100 W／m2量級(jí)，即對(duì)于分離著陸器來(lái)說(shuō)，發(fā)動(dòng)機(jī)羽流對(duì)已分離著陸器的影響在300 m內(nèi)應(yīng)該考慮，在100 m內(nèi)應(yīng)該對(duì)相關(guān)敏感元器件采取一定的防護(hù)措施。

4）當(dāng)太陽(yáng)能電池翼處于不同角度時(shí)，45°太陽(yáng)能電池翼所受到的x軸向的氣動(dòng)力最大。隨著太陽(yáng)能電池翼偏轉(zhuǎn)角度的增加，太陽(yáng)能電池翼上所受到的沿z軸的氣動(dòng)力逐漸增大；沿y軸的氣動(dòng)力先減小后增大，在45°時(shí)最小。對(duì)于相對(duì)于參考點(diǎn)的氣動(dòng)力矩來(lái)說(shuō)，沿x軸方向的羽流氣動(dòng)力矩隨著太陽(yáng)能電池翼轉(zhuǎn)角的增大而增大；沿y軸方向的羽流氣動(dòng)力矩在各角度下都基本為0；沿z軸方向的力矩隨著太陽(yáng)能電池翼轉(zhuǎn)角的增大先減小后增大，在轉(zhuǎn)角為45°時(shí)絕對(duì)值最大，方向與z軸正向相反；合力矩隨著太陽(yáng)能電池翼轉(zhuǎn)角的增大而增大。

5）熱適應(yīng)系數(shù)對(duì)羽流氣動(dòng)熱流值的仿真結(jié)果的影響是線性的，即隨著熱適應(yīng)系數(shù)的增大，羽流氣動(dòng)熱流值仿真結(jié)果呈線性增加。

（References）

［1］ Bird G.Molecular Gas Dynamics［M］.Oxford：Clarendon Press，1976.

［2］ 蔡國(guó)飆，賀碧蛟.PWS軟件應(yīng)用于探月著陸器羽流效應(yīng)數(shù)值模擬研究［J］.航天器環(huán)境工程，2010，27（01）：18-23.Cai G B，He B J.PWS software applied to numerical simulation on plume effect of lunar exploration lander［J］.Spacecraft Environment Engineering，2010，27（01）：18-23.（in Chinese）

［3］ 侯鳳龍，蔡國(guó)飆.應(yīng)用于DSMC方法的直角網(wǎng)格技術(shù)［J］.北京航空航天大學(xué)學(xué)報(bào)，2010，36（06）：694-699.Hou F L，Cai G B.Cartesian grid techmology applied in DSMC method［J］.Journal of Beijing University of Aeronautics and Astronautics，2010，36（06）：694-699.（in Chinese）

［4］ He B，He X，Zhang M，et al.Plume aerodynamic effects of cushion engine in lunar landing［J］.Chinese Journal of Aeronautics，2013，26（2）：269-278.

［5］ 蔡國(guó)飆，賀碧蛟.神舟飛船平移發(fā)動(dòng)機(jī)對(duì)太陽(yáng)能帆板羽流效應(yīng)數(shù)值模擬研究［J］.載人航天，2011，17（04）：48-53.Cai G B，He B J.DSMC simulation for thee plume effect of the shenzhou spacecraft attitude control thrusters on the solar board［J］.Manned Spaceflight，2011，17（04）：48-53.（in Chinese）

［6］ 張熇，蔡國(guó)飆，許映喬，等.嫦娥三號(hào)著陸器軟著陸過(guò)程中羽流仿真分析及試驗(yàn)研究［J］.中國(guó)科學(xué)：技術(shù)科學(xué)，2014，44（04）：344-352.Zhang H，Cai G B，Xu Y Q，et al.Simulation and experimental study of the plume during the Chang’E-3 lunar landing［J］.Scientia Sinica Technologica，2014，44（04）：344-352.（in Chinese）

［7］ 唐振宇，賀碧蛟，蔡國(guó)飆.解耦N-S／DSMC方法計(jì)算推力器真空羽流的邊界條件研究［J］.推進(jìn)技術(shù)，2014，35（07）：897-904.Tang Z Y，He B J，Cai G B.Investigation on boundary conditions in decoupled N-S／DSMC method for vacuum plume simulation of thrusters［J］.Journal of Propulsion Technology，2014，35（07）：897-904.（in Chinese）

［8］ 約翰D.安德森.計(jì)算流體力學(xué)基礎(chǔ)及應(yīng)用［M］.吳頌平，劉趙淼，譯.北京：機(jī)械工業(yè)出版社，2015.JohnD.Anderson.Foundation and Application of Computational Fluid Dynamics［M］.Wu S P，Liu Z M，translate.Beijing：China Machine Press，2015.（in Chinese）

［9］ 李斌.超燃沖壓發(fā)動(dòng)機(jī)進(jìn)氣道的數(shù)值模擬研究［D］.湘潭：湘潭大學(xué)，2013.Li B.Numerical Simulation of Scramjet Inlet／Isolator［D］.Xiangtan：Xiangtan University，2013.（in Chinese）

［10］ 沈青.稀薄氣體動(dòng)力學(xué)［M］.北京：國(guó)防工業(yè)出版社，2003.Shen Q.Rarefied Gas Dynamics［M］.Beijing：National Defense Industry Press，2003.（in Chinese）

Simulation Study of Plume Effect on Mars Probe

CHEN Jie，HE Bijiao?，CAI Guobiao
（School of Astronautics，Beihang University，Beijing 100191，China）

To mitigate the threat of the motor plume on the flight mission of the Mars probe，the hybrid NS equations／DSMC method was used to simulate the plume effect on the Mars probe.The simulation results showed that the plume of the motor mainly affected the storage box and the maximum heat flow reached 1789 W／m2.The effects of the aerodynamic force and the torque relative to the reference point on the expended solar wing were 3～4 orders’smaller than that of the motor.The plume effect on the separated Mars lander within 300 m was large，and the value of the aerodynamic heat flux was more than 100 W／m2.Within 100 m，the value of the aerodynamic heat flux was more than 1000 W／m2.So some protective measures for the sensitive equipment was necessary.The simulation with different thermal accommodation coefficients showed that with the increase of the thermal accommodation coefficient，the heat flow value increased linearly.The numerical simulation results can provide references for the engineering design of the Mars probs.

Mars exploration；Mars orbiters；plume effect；DSMC

V476.4

A

1674-5825（2017）06-0743-08

2017-04-17；

2017-10-02

陳杰，男，碩士研究生，研究方向?yàn)樽塑壙匕l(fā)動(dòng)機(jī)羽流效應(yīng)仿真。E-mail：fyaycj＠163.com

?通訊作者：賀碧蛟，男，博士，講師，研究方向?yàn)楹教炱髡婵沼鹆餍?yīng)。E-mail：hbj＠buaa.edu.cn

（責(zé)任編輯：龐迎春）