葉青,饒煒,劉鋒,孫澤洲,劉國(guó)強(qiáng),王闖,董捷,韓泉東,繆遠(yuǎn)明,譚志云,陸慧林

1.北京空間飛行器總體設(shè)計(jì)部，北京 100094 2.上?？臻g推進(jìn)研究所，上海 201112 3.哈爾濱工業(yè)大學(xué) 能源科學(xué)與工程學(xué)院，哈爾濱 150006

航天器著陸過(guò)程發(fā)動(dòng)機(jī)羽流與地外天體土壤的相互作用過(guò)程的研究，對(duì)于保障和提高探測(cè)器的可靠性具有重要意義。國(guó)內(nèi)外針對(duì)該問(wèn)題開(kāi)展了一系列的數(shù)值模擬與試驗(yàn)測(cè)試研究。

Hutton等對(duì)海盜號(hào)的著陸發(fā)動(dòng)機(jī)羽流氣體與顆粒作用進(jìn)行了研究，并針對(duì)發(fā)動(dòng)機(jī)對(duì)地外天體土壤表面的侵蝕作用進(jìn)行了評(píng)估。Scott和Immer等針對(duì)Apollo登月艙著陸過(guò)程的羽流與月壤開(kāi)展了研究工作，給出了Apollo登月艙最終著陸段的揚(yáng)塵情況和受侵蝕遷移的月壤厚度。Liever等采用整體聯(lián)合流動(dòng)求解器(Unified FlowSolver，UFS)對(duì)Apollo月球探測(cè)器羽流與月壤顆粒的作用過(guò)程和羽流流動(dòng)過(guò)程進(jìn)行預(yù)測(cè)，研究發(fā)現(xiàn)羽流出口速度會(huì)直接影響羽流流動(dòng)過(guò)程，包括羽流的壓力和速度分布。

Kuang等將著陸發(fā)動(dòng)機(jī)與陸面的作用比擬為高速射流與顆粒層的作用，采用離散顆粒方法(DEM)模擬高速氣體射流沖擊顆粒時(shí)的氣體-顆粒流動(dòng)過(guò)程，研究發(fā)現(xiàn)坑的速度和寬度與射流速度有關(guān)，但是當(dāng)坑形成后，再繼續(xù)增加射流速度，坑的寬度和深度不會(huì)繼續(xù)增加，基本維持不變。

中國(guó)在發(fā)動(dòng)機(jī)羽流和顆粒方面也做了大量的研究。王淑彥建立了二維拱形模擬體模型，開(kāi)展了其作用下的月壤顆粒運(yùn)動(dòng)過(guò)程研究。陸鑫等提出了一種離散單元接觸力模型，利用該模型模擬計(jì)算了不同速度和不同剛度范圍內(nèi)，月壤顆粒漂浮顆粒的最大高度隨時(shí)間的變化等。耿動(dòng)梁和崔玉紅等通過(guò)將CFD(Computational Fluid Dynamics)與DSMC(Direct Simulation Monte Carlo)方法進(jìn)行分區(qū)耦合，詳細(xì)地模擬出發(fā)動(dòng)機(jī)羽流場(chǎng)與月壤相互作用過(guò)程。

除了數(shù)值模擬研究，人們還建立了地球模擬試驗(yàn)，將發(fā)動(dòng)機(jī)置于地面的真空艙內(nèi)，在真空艙底部鋪設(shè)對(duì)應(yīng)的模擬顆粒。采用此類(lèi)模擬試驗(yàn)來(lái)測(cè)量空間的顆粒濃度分布，獲取顆粒激揚(yáng)和羽流作用后表面坑特性的數(shù)據(jù)。

然而，當(dāng)前的數(shù)值仿真缺乏與在軌飛行數(shù)據(jù)的對(duì)比和分析，而地面試驗(yàn)難以模擬地外天體表面重力，且成本較高。因此本文主要采用數(shù)值仿真結(jié)合火星在軌飛行數(shù)據(jù)來(lái)開(kāi)展火星著陸發(fā)動(dòng)機(jī)羽流與火壤相互作用的研究。

1 模型與方法

1.1 物理模型

針對(duì)火星探測(cè)器的狀態(tài)和特征，本文建立如圖1所示的物理模型，其中航天器和發(fā)動(dòng)機(jī)按真實(shí)外形尺寸建模，并按1.5 m/s的速度垂直勻速下降，直到落火。計(jì)算區(qū)域長(zhǎng)和寬均為50 m、高度為35 m，其中底部為厚5 m的火星土壤顆粒，火星土壤顆粒底部和兩側(cè)在高度5 m范圍內(nèi)假設(shè)為固體壁面，其余為出口。

發(fā)動(dòng)機(jī)采用火星典型著陸變推力發(fā)動(dòng)機(jī)，安裝位置與航天器底部平面平齊，發(fā)動(dòng)機(jī)高度即為羽流出口與火壤的距離，發(fā)動(dòng)機(jī)推力即為變推力發(fā)動(dòng)機(jī)當(dāng)前的羽流推力。不同發(fā)動(dòng)機(jī)推力對(duì)應(yīng)的出口燃?xì)鈪?shù)如表1所示。

火星土壤顆粒由多種不同直徑的顆粒組成，各自有不同的質(zhì)量分?jǐn)?shù)?？紤]到國(guó)外調(diào)研火壤結(jié)果中細(xì)顆粒比例最大，且顆粒粒徑越小越容易被激揚(yáng)，后續(xù)仿真統(tǒng)一選取細(xì)顆粒進(jìn)行計(jì)算，其中顆粒填充體積濃度為0.65，顆粒密度為1 100 kg/m，粒徑為15 μm。

圖1 物理模型與計(jì)算區(qū)域Fig.1 Physical model and computation area

表1 發(fā)動(dòng)機(jī)燃?xì)鈪?shù)Table 1 Engine gas parameters

1.2 數(shù)值仿真方法

1) 數(shù)值方法介紹

本文主要開(kāi)展火星探測(cè)器著陸期間羽流氣體與火壤顆粒的相互作用，以及火壤顆粒濃度的大范圍變化特征研究。采用雙流體模型進(jìn)行研究，即將發(fā)動(dòng)機(jī)的羽流假設(shè)為連續(xù)的流體，同時(shí)采用擬流體處理方法將大量離散火壤顆粒的統(tǒng)計(jì)行為用廣義連續(xù)介質(zhì)來(lái)描述。此時(shí)流動(dòng)模型的基本方程由質(zhì)量守恒方程(連續(xù)性方程)和動(dòng)量守恒方程組成。

① 質(zhì)量守恒方程

連續(xù)性方程采取的形式為

② 動(dòng)量守恒方程

從微觀連續(xù)假設(shè)出發(fā)，動(dòng)量方程為

(-)+

式中：為重力加速度；為羽流氣體與火星土壤顆粒的作用系數(shù)或者曳力系數(shù)；為火星土壤顆粒體積分率；為羽流氣體體積分率；為火星土壤顆粒密度；為羽流氣體密度；和分別為羽流方程和顆粒方程中的火星土壤顆粒速度；和分別為羽流方程和顆粒方程中的羽流氣體速度；為火星土壤顆粒相壓力；為羽流氣體相壓力；為羽流氣體相應(yīng)力；為火星土壤顆粒相應(yīng)力。

火星土壤顆粒流動(dòng)產(chǎn)生的剪切應(yīng)力由火星土壤顆粒瞬間碰撞產(chǎn)生的動(dòng)力剪切應(yīng)力和火星土壤顆?；瑒?dòng)和滾動(dòng)形成的半持續(xù)接觸產(chǎn)生的摩擦剪切應(yīng)力組成。

采用Navier-Stokes方程進(jìn)行求解，湍流模型選擇Shear Stress Transport(SST)-湍流模型，其參數(shù)為FLUENT默認(rèn)。Courant數(shù)設(shè)置為0.4，通量格式使用Roe格式。

2) 網(wǎng)格劃分及網(wǎng)格無(wú)關(guān)性檢驗(yàn)

將計(jì)算區(qū)域劃分為兩個(gè)子區(qū)域：在計(jì)算高度范圍內(nèi)以噴管軸線為中心，直徑為1.5區(qū)域?yàn)榧用苤骶W(wǎng)格區(qū)(其中，為噴管出口直徑)，其他區(qū)域?yàn)榇尉W(wǎng)格區(qū)。次網(wǎng)格平均尺寸為加密主網(wǎng)格尺寸的2倍。對(duì)4組不同的主網(wǎng)格尺寸(細(xì)網(wǎng)格尺寸為37.5 mm，中網(wǎng)格I尺寸為52 mm，中網(wǎng)格II尺寸為75 mm，粗網(wǎng)格尺寸為175 mm)進(jìn)行網(wǎng)格敏感性分析。圖2給出了不同網(wǎng)格密度下羽流氣體軸向速度的變化。不同計(jì)算網(wǎng)格下羽流氣體軸向速度的變化趨勢(shì)相同，然而，粗計(jì)算網(wǎng)格下羽流氣體速度最低，細(xì)網(wǎng)格和中網(wǎng)格I尺寸預(yù)測(cè)的羽流氣體軸向速度相差不大，中網(wǎng)格II尺寸預(yù)測(cè)的結(jié)果介于中間，表明粗網(wǎng)格和中網(wǎng)格II低估了羽流氣體速度分布。因此，綜合考慮對(duì)計(jì)算精度的要求以及計(jì)算資源的消耗，數(shù)值模擬選取中網(wǎng)格I尺寸進(jìn)行計(jì)算。

3) 計(jì)算工況參數(shù)

航天器在進(jìn)行地外天體著陸時(shí)，發(fā)動(dòng)機(jī)羽流與顆粒相互作用，在這一過(guò)程中，發(fā)動(dòng)機(jī)推力和高度以及顆粒狀態(tài)將直接影響坑的大小和顆粒激揚(yáng)的高度。

結(jié)合火星探測(cè)器著陸過(guò)程和發(fā)動(dòng)機(jī)的推力變化，開(kāi)展了仿真分析，工況如表2所示。其中將Case 1、Case 2和Case 3對(duì)比可以得出發(fā)動(dòng)機(jī)推力對(duì)顆粒的影響，將Case 2、Case 4和Case 5對(duì)比可以得出發(fā)動(dòng)機(jī)高度對(duì)顆粒的影響，Case 6用于與在軌數(shù)據(jù)對(duì)比。

圖2 不同網(wǎng)格密度下羽流氣體軸向速度的變化Fig.2 Change of plume gas axial velocity under different grid densities

表2 計(jì)算工況參數(shù)Table 2 Calculated operating parameters

2 結(jié)果分析與討論

根據(jù)以上因素，基于雙流體計(jì)算方法，以期得到發(fā)動(dòng)機(jī)推力和發(fā)動(dòng)機(jī)高度對(duì)顆粒影響的變化規(guī)律，火坑的寬度和深度通過(guò)火壤的顆粒濃度范圍來(lái)判讀，顆粒濃度為0即代表顆粒已經(jīng)被吹走，如圖3所示。

圖3 火坑的寬度和深度示意圖Fig.3 Sketch of width and depth of Mars crater

2.1 仿真結(jié)果分析

2.1.1 發(fā)動(dòng)機(jī)推力對(duì)顆粒的影響

發(fā)動(dòng)機(jī)高度=1 m時(shí)不同發(fā)動(dòng)機(jī)推力對(duì)應(yīng)的顆粒濃度云圖和三維圖分別如圖4和圖5所示。分析數(shù)據(jù)可知：

1) 發(fā)動(dòng)機(jī)高度=1 m時(shí)，隨著發(fā)動(dòng)機(jī)推力的增大，發(fā)動(dòng)機(jī)羽流吹掃火壤顆粒的范圍越來(lái)越大，火坑的寬度和深度越來(lái)越大。

2) 發(fā)動(dòng)機(jī)高度=1 m時(shí)，隨著發(fā)動(dòng)機(jī)推力的增大，顆粒與發(fā)動(dòng)機(jī)羽流的作用越來(lái)越激烈，火壤顆粒激揚(yáng)的高度越來(lái)越高。

圖4 不同發(fā)動(dòng)機(jī)推力對(duì)應(yīng)的顆粒濃度云圖(H=1 m)Fig.4 Contours of particle concentration corresponding to different engine thrusts (H=1 m)

3) 綜合分析火壤顆粒的激揚(yáng)高度、火坑深度和火坑寬度，如表3所示，數(shù)據(jù)表明為了確保航天器著陸火面的安全和穩(wěn)定，最佳方案是盡量用較小的推力狀態(tài)來(lái)著陸。

圖5 不同發(fā)動(dòng)機(jī)推力對(duì)應(yīng)的顆粒濃度三維圖(H=1 m)Fig.5 Three-dimensional diagram of particle concentration corresponding to different engine thrusts (H=1 m)

表3 不同發(fā)動(dòng)機(jī)推力下火壤激揚(yáng)高度、火坑深度和火坑寬度統(tǒng)計(jì)

2.1.2 發(fā)動(dòng)機(jī)高度對(duì)顆粒的影響

發(fā)動(dòng)機(jī)推力=5 000 N時(shí)不同發(fā)動(dòng)機(jī)高度對(duì)應(yīng)的顆粒濃度云圖和三維圖分別如圖6和圖7所示。分析數(shù)據(jù)可知：

1) 發(fā)動(dòng)機(jī)推力=5 000 N時(shí)，=3 m時(shí)發(fā)動(dòng)機(jī)羽流與火壤顆粒的相互作用可忽略，隨著發(fā)動(dòng)機(jī)高度的下降，發(fā)動(dòng)機(jī)羽流吹掃火壤顆粒的范圍越來(lái)越大，火坑的寬度和深度越來(lái)越大。

2) 發(fā)動(dòng)機(jī)推力=5 000 N時(shí)，隨著發(fā)動(dòng)機(jī)高度的下降，顆粒與發(fā)動(dòng)機(jī)羽流的作用越來(lái)越激烈，火壤顆粒激揚(yáng)的高度越來(lái)越高。

3) 綜合分析火壤顆粒的激揚(yáng)高度、火坑深度和火坑寬度，如表4所示，表明為了確保航天器著陸火面的安全和穩(wěn)定，盡量在發(fā)動(dòng)機(jī)高度較高的位置處關(guān)機(jī)(約1～3 m之間)。

圖6 不同發(fā)動(dòng)機(jī)高度對(duì)應(yīng)的顆粒濃度云圖(F=5 000 N)Fig.6 Contours of particle concentration corresponding to different engine heights (F=5 000 N)

圖7 不同發(fā)動(dòng)機(jī)高度對(duì)應(yīng)的顆粒濃度三維圖(F=5 000 N)Fig.7 Three-dimensional diagram of particle concentration corresponding to different engine heights (F=5 000 N)

表4 不同發(fā)動(dòng)機(jī)高度下火壤激揚(yáng)高度、火坑深度和火坑寬度統(tǒng)計(jì)

2.2 飛行數(shù)據(jù)驗(yàn)證

中國(guó)火星探測(cè)器在火星車(chē)上安裝相機(jī)對(duì)著陸器進(jìn)行了拍攝，通過(guò)相片的像素、灰度以及相機(jī)相對(duì)于航天器的相對(duì)位置關(guān)系可以得到火坑的寬度和深度，如圖8所示。

圖8 火星探測(cè)器著陸后照片F(xiàn)ig.8 Photos of Mars Rover after landing

火星探測(cè)器在軌采用觸火關(guān)機(jī)，關(guān)機(jī)時(shí)發(fā)動(dòng)機(jī)推力約為3 000 N，發(fā)動(dòng)機(jī)距離火面的高度約為0.35 m，通過(guò)相機(jī)對(duì)火星探測(cè)器著陸后火坑情況進(jìn)行了拍攝。首先根據(jù)相片中灰度識(shí)別出火坑的邊緣，然后通過(guò)相片的像素和相機(jī)相對(duì)于航天器的位置可以得出火坑的深度和寬度，與Case 6的仿真結(jié)果對(duì)比可知：

1) 相機(jī)拍攝的火坑形態(tài)與顆粒濃度仿真云圖基本一致，如圖9所示，其中黑色箭頭代表氣流方向，表明采用雙流體模型進(jìn)行仿真的機(jī)理基本準(zhǔn)確，且火面實(shí)際狀態(tài)接近細(xì)顆粒狀態(tài)。

2) 從火坑深度和寬度的量級(jí)來(lái)看，如表5所示，仿真值大于在軌實(shí)測(cè)值，且標(biāo)稱(chēng)誤差在30%以內(nèi)。表明仿真算法可信，偏差的來(lái)源包括火星表面的背壓不均勻以及火壤特性選取等。

圖9 相機(jī)照片與顆粒濃度仿真云圖比對(duì)Fig.9 Comparison of camera photos with simulated contours of particle concentration

表5 火坑深度和火坑寬度統(tǒng)計(jì)Table 5 Statistics of Mars crater depth and width

3 結(jié) 論

本文針對(duì)火星探測(cè)器地外天體著陸出現(xiàn)的發(fā)動(dòng)機(jī)羽流與火壤表面相互作用的問(wèn)題，針對(duì)不同的發(fā)動(dòng)機(jī)高度和推力，基于CFD方法給出了顆粒激揚(yáng)、火坑深度和寬度的仿真分析，得到了如下結(jié)論：

1) 隨著發(fā)動(dòng)機(jī)推力的增大和發(fā)動(dòng)機(jī)高度的降低，發(fā)動(dòng)機(jī)羽流吹掃火壤顆粒的范圍越來(lái)越大，火坑的寬度和深度越來(lái)越大。

2) 綜合分析火壤顆粒的激揚(yáng)高度、火坑深度和火坑寬度，為了確保航天器著陸火面的安全和穩(wěn)定，盡量在發(fā)動(dòng)機(jī)高度較高的位置處關(guān)機(jī)(約1～3 m之間)，同時(shí)盡量用較小的推力狀態(tài)來(lái)著陸。

3) 相機(jī)拍攝的火坑形態(tài)與仿真基本一致，表明采用雙流體模型進(jìn)行研究的機(jī)理基本準(zhǔn)確，且火面實(shí)際狀態(tài)接近細(xì)顆粒狀態(tài)。

4) 從火坑深度和寬度的量級(jí)來(lái)看，仿真值大于在軌實(shí)測(cè)值，且標(biāo)稱(chēng)誤差在30%以內(nèi)。表明仿真算法可信，偏差的來(lái)源包括火星表面的背壓不均勻、火壤特性難以準(zhǔn)確選取參數(shù)以及圖片識(shí)別和濃度表征之間的差別等因素。