衛(wèi)星回收艙再入過渡流區(qū)氣動(dòng)熱數(shù)值計(jì)算

賈居紅 胡麗杰

?

衛(wèi)星回收艙再入過渡流區(qū)氣動(dòng)熱數(shù)值計(jì)算

賈居紅1胡麗杰2

（1 海軍91267部隊(duì)，福州 350015）（2 中國人民大學(xué)信息學(xué)院，北京 100872）

返回式衛(wèi)星進(jìn)入地球大氣層95km高度左右時(shí)，速度達(dá)到>20，回收艙被高超聲速稀薄來流形成的弓形激波環(huán)繞，氣動(dòng)加熱問題非常明顯，準(zhǔn)確預(yù)測(cè)過渡流區(qū)的氣動(dòng)熱成為一個(gè)十分突出的問題。過渡區(qū)由于氣體空氣分子仍然比較密集，基于分子動(dòng)力學(xué)的直接蒙特卡洛模擬方法（DSMC）極其耗費(fèi)計(jì)算資源，而求解Navier-Stokes（N-S）方程的方法誤差較大。采用添加二階滑移條件的N-S方程求解過渡流區(qū)氣動(dòng)熱，并與開源DSMC2V程序計(jì)算結(jié)果對(duì)比，研究了艙體母線熱流、壓強(qiáng)系數(shù)變化；通過對(duì)比兩種方法下艙體前緣弓形激波及流場(chǎng)參數(shù)變化，分析了滑移條件影響壁面和流場(chǎng)參數(shù)的機(jī)理。結(jié)果表明：回收艙再入過渡區(qū)時(shí)，鈍頭駐點(diǎn)區(qū)為高壓、高熱流區(qū)，錐身區(qū)氣動(dòng)熱和壁面壓強(qiáng)保持在較低水平。帶滑移條件的計(jì)算機(jī)流體力學(xué)方法（Computational Fluid Dynamics，CFD）計(jì)算得到的熱流、壓強(qiáng)系數(shù)與DSMC結(jié)果吻合良好，具有計(jì)算效率高、精度較高的優(yōu)勢(shì)。對(duì)流場(chǎng)壓力、溫度、速度等參數(shù)分析顯示，滑移條件中壁面速度滑移和溫度跳躍的加入，改變了壁面流動(dòng)參數(shù)，進(jìn)而改善了壁面熱流和壓強(qiáng)的準(zhǔn)確模擬能力?；茥l件對(duì)外流場(chǎng)參數(shù)的影響極小，沒有改善稀薄流區(qū)流場(chǎng)的模擬能力和激波捕獲能力，模擬得到的溫度、速度等流場(chǎng)參數(shù)及激波位置、激波層厚度等與DSMC結(jié)果仍有差別?？梢哉J(rèn)為，文章采用方法能夠滿足工程上快速高效預(yù)測(cè)衛(wèi)星回收艙再入過渡流區(qū)氣動(dòng)熱需要。

高超聲速 氣動(dòng)熱 回收艙 計(jì)算流體力學(xué) 滑移條件 返回式衛(wèi)星

0 引言

返回式衛(wèi)星主要用于遙感和空間科學(xué)實(shí)驗(yàn)，具有顯著的經(jīng)濟(jì)和社會(huì)效益。其微重力條件可以達(dá)到10–5n甚至更高的水平，這是載人飛船與航天飛機(jī)無法比擬的[1-3]。近年來，返回式衛(wèi)星由于其優(yōu)越的微重力環(huán)境和低成本，而在空間科學(xué)試驗(yàn)領(lǐng)域受到各國青睞，最具代表的是俄羅斯和歐空局的“光子”號(hào)和“生物”號(hào)系列空間科學(xué)實(shí)驗(yàn)衛(wèi)星。

衛(wèi)星完成在軌任務(wù)后，經(jīng)過制動(dòng)減速返回地球，再入過程中速度達(dá)>20，高超聲速來流形成弓形激波包裹在衛(wèi)星周圍，流場(chǎng)溫度接近2 000K[4]。除了來流空氣與衛(wèi)星艙體摩擦生熱之外，高溫條件下空氣發(fā)生化學(xué)反應(yīng)、返回艙體表面材料催化特性也會(huì)給返回艙壁面熱載荷帶來影響[5-6]，準(zhǔn)確預(yù)測(cè)再入過程氣動(dòng)熱是回收艙再入需要解決的難題之一。地面預(yù)測(cè)高超聲速氣動(dòng)熱的方法主要有風(fēng)洞試驗(yàn)、數(shù)值模擬和工程預(yù)估，由于稀薄空氣高超聲速氣動(dòng)熱風(fēng)洞試驗(yàn)難度較大且費(fèi)用高昂，工程預(yù)估精度難于保證，采用數(shù)值方法預(yù)測(cè)再入過程氣動(dòng)加熱成為一種經(jīng)濟(jì)高效的選擇[7]。目前，數(shù)值模擬高超聲速氣動(dòng)熱的方法主要有基于分子動(dòng)力學(xué)的直接蒙特卡洛模擬法（DSMC）[8-9]和求解Navier-Stoke（N-S）方程的計(jì)算機(jī)流體力學(xué)方法（CFD）?？諝獾南”〕潭韧ǔＳ每伺瓟?shù)來表征[10]，定義為分子平均自由程與流動(dòng)特征尺度的比值。對(duì)于連續(xù)流區(qū)（<0.001），帶滑移條件的計(jì)算機(jī)流體力學(xué)方法（Computational Fluid Dynamics，CFD）方法能具有較高氣動(dòng)熱預(yù)測(cè)精度；對(duì)于自由分子流區(qū)（>10），DSMC方法是不二之選；對(duì)于過渡流區(qū)（0.1<<10）而言，由于氣體分子比較密集，DSMC方法極其耗費(fèi)計(jì)算資源，而基于連續(xù)性假設(shè)的N-S方程方法誤差過大。因此，添加滑移條件的N-S方程方法成為一種折中選擇[11-12]。本文以DSMC方法為對(duì)照，采用添加滑移條件的N-S方程求解過渡流區(qū)小型再入衛(wèi)星表面氣動(dòng)熱，為工程上多工況快速預(yù)測(cè)提供借鑒。

1 幾何模型與來流條件

衛(wèi)星返回艙外形如圖1（a），前部為回收艙，后部為制動(dòng)艙?；厥张擉w布局為球錐結(jié)構(gòu)形式[13]，球頭半徑為0.267 8m，后體半錐角11.4°，底部半徑0.503 5m，艙體總長(zhǎng)1.41m，詳細(xì)參數(shù)如圖1（b）所示。圖1中以返回艙頭部頂點(diǎn)為坐標(biāo)原點(diǎn)，橫坐標(biāo)為艙體軸向距離,縱坐標(biāo)為艙體切向距離。

圖1 衛(wèi)星回收艙體結(jié)構(gòu)參數(shù)

伴隨流條件為相應(yīng)高度的國際標(biāo)準(zhǔn)大氣[14]，衛(wèi)星返回至95km高空時(shí)速度7 860m/s，>28.45，以衛(wèi)星球頭半徑為特征尺度的克努森數(shù)（）為0.216，屬于過渡流區(qū)。詳細(xì)來流參數(shù)如表1所示。文獻(xiàn)[15]采用DSMC方法對(duì)三維鈍頭再入艙進(jìn)行了化學(xué)平衡和非平衡兩種工況計(jì)算，本文以此為參照，采用帶二階滑移條件的CFD方法進(jìn)行模擬，并與開源的DSMC2V程序[16]二維計(jì)算結(jié)果進(jìn)行對(duì)比分析。

表1 來流條件

Tab.1 Flow conditions

2 計(jì)算方法

2.1 數(shù)值方法

本文采用雷諾平均N-S方程求解，其形式可表示為[17]

2.2 滑移條件

采用Beskok等[21]提出的二階滑移條件，其壁面速度滑移和溫度跳躍公式為

2.3 網(wǎng)格劃分與邊界條件

圖2 對(duì)稱面網(wǎng)格及邊界條件

3 結(jié)果與討論

3.1 網(wǎng)格無關(guān)性驗(yàn)證

物體壁面法向網(wǎng)格尺度對(duì)氣動(dòng)熱模擬結(jié)果影響較大，為了消除網(wǎng)格因素的影響，采用由疏到密的5套網(wǎng)格，對(duì)二維軸對(duì)稱模型進(jìn)行了網(wǎng)格無關(guān)性測(cè)試。測(cè)試網(wǎng)格尺度、網(wǎng)格延伸率及計(jì)算獲得壁面+如表2所示。

表2 測(cè)試網(wǎng)格參數(shù)

Tab.2 Test grids parameters

表3 不同網(wǎng)格下熱流與壓強(qiáng)系數(shù)峰值（誤差）

Tab.3 Heat flux & pressure coefficient with different grid（Difference）

圖3 不同網(wǎng)格尺度下壁面熱流與壓強(qiáng)系數(shù)

3.2 壁面熱流與壓強(qiáng)分布

為了分析帶滑移條件的CFD方法模擬再入體氣動(dòng)熱的能力，以DSMC2V計(jì)算結(jié)果為對(duì)照，對(duì)壁面熱流、壓強(qiáng)和流場(chǎng)參數(shù)進(jìn)行對(duì)比分析。表4為駐點(diǎn)熱流與壓強(qiáng)系數(shù)峰值對(duì)比，可以看出：相對(duì)于DSMC非平衡流結(jié)果，CFD無滑移計(jì)算得到的量綱—熱流系數(shù)峰值誤差13.3%，添加滑移條件后，CFD滑移熱流峰值誤差降低至2.9%；CFD無滑移壓強(qiáng)系數(shù)峰值誤差2.8%，添加滑移條件后，壓強(qiáng)峰值誤差降低至0.5%。

表4 駐點(diǎn)熱流與壓強(qiáng)系數(shù)峰值（誤差）

Tab.4 Peak heat flux & pressure coefficient（Difference）

圖4（a）為三維壁面熱流分布，可以看到熱流峰值出現(xiàn)在球頭駐點(diǎn)區(qū)域，并沿球頭母線快速遞減，錐體部分熱流密度遠(yuǎn)低于駐點(diǎn)區(qū)，并且保持緩慢降低。圖4（b）為沿艙體對(duì)稱面母線的熱流系數(shù)，可以看出，無滑移CFD方法模擬壁面熱流明顯高于DSMC結(jié)果，誤差較大；添加滑移條件后模擬得到數(shù)據(jù)曲線與DSMC非平衡結(jié)果較為接近，誤差較小。圖4（c）給出了艙體表面三維壓強(qiáng)分布，可以看出壓強(qiáng)與熱流分布規(guī)律較為相似，駐點(diǎn)區(qū)壓強(qiáng)明顯大于錐身區(qū)域，錐體部分壓強(qiáng)保持在較低水平。圖4（d）給出沿艙體對(duì)稱面母線壓強(qiáng)系數(shù)分布，可以看出無滑移CFD與帶滑移條件CFD結(jié)果均與DSMC變化趨勢(shì)相同，但添加滑移條件后，模擬得到的壁面壓強(qiáng)更接近DSMC非平衡結(jié)果，誤差較小。綜合表4和圖4可以看出，添加滑移條件的CFD方法能較為精確模擬過渡流區(qū)衛(wèi)星回收艙表面熱流分布和壓強(qiáng)分布，滿足工程上氣動(dòng)熱預(yù)測(cè)需要。

圖4 壁面熱流與壓強(qiáng)分布

3.3 流場(chǎng)參數(shù)分布

綜合對(duì)比圖5（b）（d）（f）中是否添加滑移條件的兩種CFD方法模擬結(jié)果，可以看出：添加滑移條件后，對(duì)返回艙壁面附近的速度和溫度參數(shù)進(jìn)行了矯正，進(jìn)而造成壁面溫度梯度變化，進(jìn)一步改善壁面熱流和壁面壓強(qiáng)分布，使之更加接近真實(shí)狀態(tài)。綜合圖4、圖5可以看出，滑移條件對(duì)外流場(chǎng)壓強(qiáng)、速度帶來影響非常小，沒有改變?cè)辛鲌?chǎng)分布，也沒有改善CFD方法模擬稀薄流區(qū)流場(chǎng)能力。

最后，對(duì)比CFD方法與DSMC方法在計(jì)算資源和計(jì)算效率上的差別。計(jì)算機(jī)配置為英特爾Core i5-650，3.2G Hz，雙核4線程，內(nèi)存8G，本文計(jì)算均采用單線程。兩種方法計(jì)算至完全收斂的時(shí)長(zhǎng)如表4所示?？梢钥闯觯珼SMC由于依靠大量分子進(jìn)行模擬，其計(jì)算量巨大，即使二維結(jié)果也需要計(jì)算20小時(shí)以上才能收斂，而三維CFD方法僅需要2小時(shí)左右。此外，必須看到，當(dāng)進(jìn)入稠密大氣層時(shí)，DSMC方法效率大為降低，耗時(shí)將增加數(shù)十倍，但CFD方法計(jì)算效率變化不會(huì)太大。

表5 計(jì)算效率

Tab.5 Computational efficiency

4 結(jié)束語

本文以DSMC方法為對(duì)比，采用CFD方法模擬了克努森數(shù)不太大的過渡流區(qū)衛(wèi)星回收艙表面壓力和熱流分布，研究了艙體母線熱流、壓強(qiáng)系數(shù)變化，分析了艙體前緣弓形激波厚度及流場(chǎng)溫度、壓強(qiáng)分布。數(shù)值模擬結(jié)果顯示：

1）添加移條件后，CFD方法能夠有效求解過渡流區(qū)回收艙氣動(dòng)熱和壁面壓強(qiáng)分布，模擬的母線方向壁面熱流、壓強(qiáng)系數(shù)與DSMC結(jié)果吻合較好，具有計(jì)算效率高、精度較高的優(yōu)勢(shì)，滿足工程設(shè)計(jì)階段氣動(dòng)熱預(yù)測(cè)需求。

2）駐點(diǎn)前緣流線方向壓力、溫度、速度對(duì)比分析顯示，添加滑移條件改善了壁面熱流和壓強(qiáng)分布，但沒有改善CFD方法模擬稀薄流區(qū)流場(chǎng)的能力，得到的流場(chǎng)溫度、壓強(qiáng)、速度分布與DSMC結(jié)果有所差別。

[1] ANTONIO V, GIUSEPPE P. Nonequilibrium Aerothermodynamics for a Capsule Reentry Vehicle[J]. Engineering Applications of Computational Fluid Mechanics, 2009, 3(4): 543-561.

[2] 王希季. 中國返回式衛(wèi)星的發(fā)展—紀(jì)念中國第一顆返回式衛(wèi)星成功返回40周年[J]. 航天返回與遙感, 2015, 36(6): 1-5. WANG Xiji. Development of Chinese Returnable Satellite—To Commemorate the 40 Anniversary of the Successful Return of Chinese First Returnable Satellite[J]. Spacecraft Recovery & Remote Sensing, 2015, 36(6): 1-5. (in Chinese)

[3] 邢連群. 返回式衛(wèi)星燒蝕防熱結(jié)構(gòu)的工程計(jì)算[J]. 中國空間科學(xué)技術(shù), 1991, 4(2): 26-34. XING Lianqun. An Engineering Computation of Ablative Thermal Protection Structure of Returnable Satellite[J]. Chinese Space Science and Technology, 1991, 4(2): 26-34. (in Chinese)

[4] 唐伯昶. 中國返回式衛(wèi)星[J]. 中國航天, 2008(5): 8-17. TANG Bochang. Chinese Returnable Satellites[J]. Aerospace China, 2008(5): 8-17. (in Chinese)

[5] SCANLON T, WHITE J, BORG C, et, al. Open Source Direct Simulation Monte Carlo Chemistry Modeling for Hypersonic Flows[J]. AIAA Journal, 2015, 53(6): 1670-1680.

[6] MAZAHERI A R. Computational Aerothermodynamic Analysis of the Space Shuttle Orbiter Tile Overlay Repair with Different Geometries[R]. Reston, VA, US: AIAA, 2004.

[7] YANG G D, DUAN Y H, LIU C Z, et al. Approximate Prediction for Aerodynamic Heating and Design for Leading-edge Bluntness on Hypersonic Vehicles[R]. Reston, VA, US: AIAA, 2014.

[8] RODRIGO C, CRAIG W, MATTHEW K B, et al. Benchmark Numerical Simulations of Rarefied Non-reacting Gas Flows Using an Open-source DSMC Code[J]. Computers and Fluids, 2015(120): 140-157.

[9] MOSS J N, GLASS C E, GRANNE F A. DSMC Simulations of Apollo Capsule Aerodynamics for Hypersonic Rarefied Conditions[C]//Thermophysics and Heat Transfer Conference. San Francisco, CA, US, 2006.

[10] 賈居紅. 臨近空間高超聲速氣動(dòng)熱數(shù)值模擬研究[D]. 北京: 北京理工大學(xué), 2017. JIA Juhong. Numerical Study on Hypersonic Aerothermodynamic of Near Space Flight Vichle[D]. Beijing: Beijing Institute of Technology, 2017.

[11] SCANLON T J, WHITE C, BORG M K, et, al. Open-source Direct Simulation Monte Carlo Chemistry Modeling for Hypersonic Flows[J]. AIAA Journal, 2015, 53(6): 1670-1680.

[12] FREIDOONIMEHR N, RASHIDI M M, YANG Z G, et al. Velocity Slip and Temperature Jump Effects for an Unsteady Flow over a Stretching Permeable Surfaces[J]. Thermal Science, 2016(1): 1-8.

[13] SANTOS W F. Aero-thermodynamic Analysis of a Reentry Brazilian Satellite[J]. Brazilian Journal of Physics, 2012, 42(5): 373-390.

[14] Standard Atmosphere Computations[EB/OL].[2018-3-23]. http://www-mdp.eng.cam.ac.uk.

[15] PAHLARINI R C, AZEVEDO J L. Thermochemical Nonequilibrium Computations of a Brazilian Reentry Satellite[J]. Journal of Spacecraft & Rockets, 2017, 54(4): 1-10.

[16] BIRD G A. The DSMC Method[M]. Oxford: Claredon Press, 2013: 57-92.

[17] WILCOX D C. Turbulence Modeling for CFD[M]. 3rd ed. California: DCW lndustries Inc, 2006: 243-249.

[18] GOLDBERG U C, OTA D K. A k-e Near-wall Formulation for Separated Flows[C]//AIAA 21st Fluid Dynamics, Plasma Dynamicsand Lasers Conference. Seattle, WA, US, 1990.

[19] GOLDBERG U C. Hypersonic Turbulent Flow Predictions Using CFD++[R].Reston, VA, US: AIAA, 2005.

[20] FURFARO D, SAUREL R. A Simple HLLC-type Riemann Solver for Compressible Non-equilibrium Two-phase Flows[J]. Computer & Fluid, 2015(111): 159-178.

[21] BESKOK A, KARNIADAKIS G E, TRIMMER W. Rarefaction and Compressibility Effects in Gas Microflows[J]. Transactions of the ASME, 1996, 118(9): 448-456.

[22] 毛宏霞, 賈居紅, 傅德彬, 等. HIFIRE-1飛行器激波與邊界層干擾氣動(dòng)熱研究[J]. 兵工學(xué)報(bào), 2018(3): 528-536. MAO Hongxia, JIA Juhong, FU Debin, et al. Research on Aero Thermodynamics and Influencing Factors for HIFIRE-1[J]. Acta Armamentarii, 2018(3): 528-536. (in Chinese)

[23] NASA Viscous Grid Spacing Calculator[EB/OL].[2018-3-30]. http://geolab.larc.nasa.gov/APPS/YPlus.

Aerothermodynamics Computations of Reentry Satellite in Transition Region

JIA Juhong1HU Lijie2

（1 Unit 91267, PLA Navy, Fuzhou 350015, China）（2 School of Information, Renmin University of China, Beijing 100872, China）

When the reentry satellite comes into the altitude of 95 km, its reentry Mach number is more than 20, high temperature wave surrounded the vehicle, and aero-heating becomes a serious problem. the Navier-Stokes (N-S) equations based on second-order velocity slip and temperature jump boundary conditions suggested by Beskok and Karniadakis are used to calculate the heat flux and pressure of the capsule, the results are compared with the direct simulation Monte Carlo (DSMC) method result, the DSMC2V opened by Bird is used. The influence of slip effects on aerothermodynamics properties are analyzed, the pressure contours, Mach numbers and temperature contours of flow field are studied. The result shows that, N-S equations with velocity slip and temperature jump changed the near wall flow contours, then it can calculate the aerothermodynamics and aerodynamics correctly, and has the advantage of high efficiency and high accuracy. Havever, the analysis of flow contours depicts, slip boundary condition do not change the out flow conditions. As a result, the simulated flow contours, along with the position and thickness of shock wave is a little different from the DSMC results. It can be conclude that this method is satisfied with the prediction of heat flux for reentry satellite in transition regime.

hypersonic; aerothermodynamics; capsule; Computational fluid dynamics; slip boundary conditions; reentry satellite

V211

A

1009-8518(2018)06-0037-09

10.3969/j.issn.1009-8518.2018.06.005

2018-07-20

賈居紅，男，1989年生，2017年獲北京理工大學(xué)宇航學(xué)院兵器科學(xué)與技術(shù)專業(yè)碩士學(xué)位，工程師。研究方向?yàn)樵偃塍w氣動(dòng)熱問題。E-mail：jiajuhongbit@126.com。

胡麗杰，女，1995年生，現(xiàn)在中國人民大學(xué)信息學(xué)院攻讀碩士學(xué)位。主要研究方向?yàn)槲⒎址匠虜?shù)值解問題研究。E-mail：lily_950827@126.com。

（編輯：劉穎）