再入飛行器端頭燒蝕的耦合計算方法

龍麗平 韓俊 萬田 仲峰泉, 田保未

再入飛行器端頭燒蝕的耦合計算方法

龍麗平1韓俊2萬田1仲峰泉1,2田保未2

（1 中國科學(xué)院力學(xué)研究所高溫氣體動力學(xué)國家重點實驗室，北京 100190；2 中國科學(xué)院大學(xué)工程科學(xué)學(xué)院，北京 100049）

端頭燒蝕是超高速飛行器再入過程中非常關(guān)心的問題。端頭材料在高溫高壓環(huán)境中，因燒蝕其原有氣動外形和結(jié)構(gòu)傳熱邊界不斷變化，而氣動外形和傳熱邊界的變化又反過來影響端頭熱流、溫度分布和燒蝕量。它們之間表現(xiàn)出復(fù)雜的強耦合、非線性特征。本文以碳基材料端頭帽燒蝕過程為例，發(fā)展了端頭帽繞流、燒蝕和結(jié)構(gòu)傳熱耦合計算方法。通過氣動、燒蝕和結(jié)構(gòu)熱響應(yīng)計算程序的耦合和迭代，實現(xiàn)了對端頭帽再入燒蝕過程的實時動邊界模擬，并在飛行試驗條件下，得到了與測量數(shù)據(jù)基本吻合的結(jié)果。

再入過程；端頭帽；燒蝕；移動邊界；多學(xué)科耦合計算方法

0 引言

熱防護(hù)一直是高超聲速飛行器及各類飛行器研制的核心技術(shù)之一。對于高馬赫數(shù)再入大氣層或在大氣層內(nèi)長時間巡航的飛行器來說，結(jié)構(gòu)的熱防護(hù)及其對飛行器氣動性能的影響始終是一個關(guān)鍵的問題[1-3]。端頭帽作為鈍頭體飛行器的頭部，其氣動加熱環(huán)境最為惡劣，也是飛行器熱防護(hù)的重點區(qū)域之一。端頭帽的選材經(jīng)歷了難熔金屬、陶瓷、石墨、C/C以及在C/C材料中添加金屬材料的低燒蝕C/C等一系列材料應(yīng)用過程[4-6]。國內(nèi)外許多學(xué)者針對各種材料進(jìn)行了燒蝕性能的研究[7-14]和燒蝕模型的發(fā)展[15-18]，對端頭結(jié)構(gòu)設(shè)計及優(yōu)化起到了重要的指導(dǎo)作用。但以往報道的再入燒蝕分析大多以理論分析和解析關(guān)系為主，并且大多研究針對某一來流條件（飛行速度、高度）研究燒蝕特性，缺乏對再入燒蝕全程以及對燒蝕參數(shù)多維分布的系統(tǒng)性研究。因此本文旨在針對端頭帽結(jié)構(gòu)和環(huán)境條件，考慮燒蝕與結(jié)構(gòu)傳熱、燒蝕外形與氣動力/熱的耦合作用關(guān)系，編制了氣動—燒蝕—熱傳導(dǎo)耦合計算程序，實現(xiàn)了對再入過程中端頭帽燒蝕與結(jié)構(gòu)演化的全程模擬。該方法為高超聲速飛行器端頭帽的防熱設(shè)計與優(yōu)化以及評估燒蝕對飛行器氣動性能的影響提供了有力工具。

1 計算流程

本文發(fā)展的端頭帽外流、燒蝕以及結(jié)構(gòu)傳熱耦合計算方法的思路如下：首先將再入過程按照飛行軌跡離散成若干準(zhǔn)定常的飛行狀態(tài)，對于每一個飛行狀態(tài)及其時間，將外流與燒蝕/結(jié)構(gòu)計算解耦。首先進(jìn)行外流計算，將外流計算獲得的邊界層傳熱系數(shù)、摩阻系數(shù)、壓力等參數(shù)作為輸入條件，求解碳燒蝕速率、組分等。燒蝕與端頭結(jié)構(gòu)傳熱實時耦合計算，兩者自動完成溫度、熱流和燒蝕量等參數(shù)的傳遞并實時更新燒蝕外形進(jìn)行結(jié)構(gòu)熱傳導(dǎo)計算。燒蝕/結(jié)構(gòu)計算獲得的燒蝕結(jié)構(gòu)外形（由燒蝕后退量決定）、溫度、碳蒸汽引射速率等參數(shù)再作為外流計算的邊界條件，重新計算外流。如此，形成了外流與燒蝕/結(jié)構(gòu)的耦合計算過程。

圖1 繞流—燒蝕—結(jié)構(gòu)耦合計算流程圖

2 計算模塊

2.1 氣動計算

連續(xù)介質(zhì)計算方法基于N-S方程求解，包含熱化學(xué)非平衡效應(yīng)。質(zhì)量守恒方程為

振動能守恒方程為

總能守恒方程為

化學(xué)反應(yīng)機理采用Park 1990模型（REF），逆向反應(yīng)速率常數(shù)由公式（5）給出。

其中k為反應(yīng)平衡常數(shù)，由NASA CEA程序求出?；旌蠚怏w的輸運系數(shù)由Blottner公式和Wilke混合模型給出，其中擴(kuò)散系數(shù)采用定Lewis數(shù)假設(shè)，Lewis數(shù)采用標(biāo)準(zhǔn)值0.7，振動松弛模型采用Landau-Teller公式

控制方程采用有限體積方法離散，其中無粘項采用三階迎風(fēng)MUSCL TVD格式，限制器為minmod。粘性項采用二階中心差分。整個方程組采用DPLR方法[19]隱式求解。湍流模型采用k-omega SST雙方程模型。圖2為典型再入鈍頭體壁面的湍流熱流分布圖，通常層流的熱流峰值在頭部駐點，而層流轉(zhuǎn)捩為湍流后，熱流峰值則向下游移動到肩部，這會導(dǎo)致鈍頭體肩部出現(xiàn)局部燒蝕量增大的后果。

2.2 燒蝕計算

本文建立了碳材料的熱化學(xué)燒蝕模型，給出不同溫度、壓力下燒蝕速率與產(chǎn)物組分，發(fā)展了燒蝕與端頭結(jié)構(gòu)非穩(wěn)態(tài)熱傳導(dǎo)耦合的求解方法。

圖2 典型再入鈍頭體壁面的湍流熱流分布圖

圖3 燒蝕表面熱流、質(zhì)量交換模型

進(jìn)一步的，我們根據(jù)燒蝕表面的質(zhì)量守恒條件，建立壁面的元素濃度與質(zhì)量損失因子的相容關(guān)系

有關(guān)燒蝕花紋與機械剝蝕的預(yù)測模型仍在發(fā)展中，隨著對燒蝕機理的深入認(rèn)識與測試手段的改進(jìn)，燒蝕花紋與機械剝蝕模型的準(zhǔn)確性及適用范圍將會有很大改善。

2.3 結(jié)構(gòu)傳熱計算

結(jié)構(gòu)熱傳導(dǎo)控制方程為瞬態(tài)熱傳導(dǎo)方程

邊界條件

基于ABAQUS二次開發(fā)平臺，同時還實現(xiàn)了燒蝕—熱傳導(dǎo)模塊與氣動計算模塊的非實時耦合計算，根據(jù)燒蝕量的情況，每隔一定時間間隔進(jìn)行數(shù)據(jù)交換。以下為利用本文發(fā)展的燒蝕多學(xué)科耦合計算程序測試的一個算例。

3 算例

以某飛行試驗為例，端頭為石墨材料，其外形如圖4所示，從高度70km，以5.8km/s的速度、零度攻角再入飛行，計算總時間約28.25s。圖5為軸對稱結(jié)構(gòu)計算模型。=0.032m，錐角9.8°。經(jīng)本文的流/固耦合燒蝕程序計算，給出圖6的端頭帽溫度隨時間變化云圖，圖7為駐點附近各個時間段的燒蝕量云圖。圖8為該算例端頭帽燒蝕最終外形與初始外形對比。如圖8所示，燒蝕外形近似為典型的雙錐結(jié)構(gòu)[24]。其中燒蝕錐角為51.72度。這與該端頭帽燒蝕的測量數(shù)據(jù)（錐角約48度）吻合的較好（誤差為7.75%）。計算給出的駐點燒蝕后退量為0.306R（約9.79mm），測量數(shù)據(jù)為0.328R（約10.51 mm）。針對端頭帽再入燒蝕的復(fù)雜問題，計算數(shù)據(jù)與試驗數(shù)據(jù)取得了令人滿意的相符程度。

圖4 端頭帽幾何模型

圖5 計算模型

圖6 隨時間變化的溫度云圖（單位: K）

圖8 端頭帽燒蝕外形與初始外形對比

從細(xì)節(jié)看，計算和實驗結(jié)果的偏差主要出現(xiàn)在壁面附近速度大的地方，這可能是因為在肩點之后，邊界層與壁面相互作用，出現(xiàn)了湍流的緣故，端頭燒蝕花紋對于氣動熱流產(chǎn)生了放大作用，導(dǎo)致這部分的計算后退量相較于實驗值偏小。經(jīng)過分析，后退量的差異主要是表面粗糙度模型以及碳基材料高溫?zé)崃W(xué)參數(shù)（如導(dǎo)熱系數(shù)、比熱容）的誤差所導(dǎo)致的。這說明燒蝕花紋與氣動計算的相互作用影響燒蝕后退量，也是本文計算方法今后需要重點改進(jìn)之處。

4 結(jié)論

本文從端頭燒蝕的多學(xué)科耦合物理過程出發(fā)，將自編氣動計算程序、燒蝕計算程序和商業(yè)軟件ABAQUS相結(jié)合，實現(xiàn)了基于移動網(wǎng)格技術(shù)的多學(xué)科耦合燒蝕計算程序，并把該程序應(yīng)用來模擬某飛行試驗再入過程中端頭帽的燒蝕，順利實現(xiàn)了氣動、燒蝕、傳熱三個模塊的數(shù)據(jù)傳輸和聯(lián)合求解，給出了整個再入過程中端頭帽燒蝕與結(jié)構(gòu)傳熱的時間、空間演化過程，并且具有較好的計算效率和計算精度。同時，結(jié)合流動穩(wěn)定性理論、材料屬性等，該方法初步考慮了湍流轉(zhuǎn)捩、碳—碳表面花紋演化以及機械剝落等因素。計算結(jié)果總體上與試驗數(shù)據(jù)相符，基本模擬了再入過程中燒蝕外形影響氣動特性，氣動特性又返回來影響端頭熱流和燒蝕量的耦合效應(yīng)。本文建立的高超聲速飛行器端頭燒蝕多學(xué)科耦合計算方法，較之前的單學(xué)科計算燒蝕方法有較大的進(jìn)步，為高超聲速飛行器端頭防熱設(shè)計與優(yōu)化提供了有力工具。

[1] 賈洲俠, 吳振強, 吳建國, 等.飛行器氣動熱與結(jié)構(gòu)傳熱雙向耦合研究[J]. 強度與環(huán)境, 2019, 46(6): 16-23.[JIA Zhouxia, WU Zhenqiang, WU Jianguo, et al. Study on two-way coupled fluid-structure-thermal analysis for hypersonic vehicles[J].Structure & Environment Engineering, 2019, 46(6): 16-23.]

[2] 何欽華, 劉佳, 田玉坤, 等. 高超聲速飛行器低氣壓環(huán)境結(jié)構(gòu)熱試驗控制技術(shù)研究[J]. 強度與環(huán)境, 2019, 46(6): 9-15.[HE Qinhua, LIU Jia, TIAN Yukun, et al. Research on thermal test control technology of hypersonic vehicle in low atmosphere pressure environment[J]. Structure & Environment Engineering, 2019, 46(6): 9-15.]

[3] 張東, 任東鴻, 趙春章, 等. 高速飛行器流固界面熱力載荷轉(zhuǎn)換方法及其應(yīng)用[J]. 強度與環(huán)境, 2019, 46(4): 54-58.[Zhang Dong, Ren Donghong, Zhao Chunzhang, et al. The high speed aircraft fluid-solid interface thermal load conversion and its application [J]. Structure & Environment Engineering, 2019, 46(4): 54-58.]

[4] 楊雯, 霍浩亮, 李海波, 等. 航天多功能熱控材料及結(jié)構(gòu)研究進(jìn)展[J]. 強度與環(huán)境, 2020, 47(2): 1-12.[Yang Wen, Huo Haoliang, Li Haibo, et al. Research progress of multifunctional thermal control materials and structures of aerospace vehicles[J]. Structure & Environment Engineering, 2020, 47(2): 1-12.]

[5] 巴曉蕾, 郭文婧, 羅嘉, 等. 航天復(fù)合材料機翼疲勞試驗加載技術(shù)研究[J]. 強度與環(huán)境, 2019 (3): 14-19.[Ba Xiaolei, Guo Wenjing, Luo Jia, et al. Research on loading technology of the fatigue life testing for composite wings[J]. Structure & Environment Engineering, 2019 (3): 14-19.]

[6] 高魁垠, 李海波, 吳建國, 等. 2D-C/SiC復(fù)合材料螺栓連接結(jié)構(gòu)可靠性分析[J].強度與環(huán)境, 2020, 47(1): 33-40.[Gao Kuiyin, LI Haibo, Wu Jianguo, et al. Reliability analysis of 2D-C/SiC composite bolted joints[J]. Structure & Environment Engineering, 2020, 47(1): 33-40.]

[7] 黃海明, 杜善義, 施惠基. C/TaC/C復(fù)合材料的燒蝕性能研究[J]. 航空兵器, 2002, 5:5-7. [Huang H M, Du S Y, Shi H J. The ablation property of C/TaC/C [J]. Hang Kong Bin Qi, 2002, 5: 5-7.]

[8] 黃海明, 杜善義, 吳林志, 等. C/C復(fù)合材料燒蝕性能分析[J]. 復(fù)合材料學(xué)報, 2001, 18(3):76-80. [Huang H M, Du S Y, Wu L Z, et al. Analysis of the ablation of C/C composites [J]. Acta Materiae Compositae Sinica, 2001, 18(3): 76-80.]

[9] 王德文, 楊月誠, 查柏林, 等. 炭/炭復(fù)合材料耐高溫和高速粒子燒蝕/侵蝕研究[J]. 彈箭與制導(dǎo)學(xué)報, 2013, 33(1): 132-134. [Wang D W, Yang Y C, Cha B L., et al. The study on anti-ablative property of C/C Composites [J]. Journal of Projectiles, Rockets, Missiles and Guidance, 2013, 33(1): 132-134.]

[10] R R Williams. ATJ-S Graphite Ablation at Angle of Attack[C]. AIAA 8th Thermophysics Conference, Palm Spring, California, July 16-18, 1973.

[11] W R Grabowsky, H P Eldridge, R C Wheaton. Reentry vehicle instrumentation survey [C]. AIAA 1st Flight Testing Conference, Las Vegas, Nevada, November 11-13, 1981.

[12] J W Metzger, M J Engel, N S Diaconis. Oxidation and sublimation of graphite in simulated re-entry environments [J]. AIAA Journal, 1976, 5(8): 451-460.

[13] 尹健, 張紅波, 熊翔, 等. 燒蝕條件對混合基體C/C復(fù)合材料燒蝕性能的影響[J]. 材料科學(xué)與工程學(xué)報, 2007, 25(1): 26-29. [Yin J, Zhang H B, Xiong X, et al. Influences of ablation conditions on the ablation performances of Mix-matrix C/C composites [J]. Journal of Materials Science and Engineering. 2007, 25(1): 26-29.]

[14] L W Hunter, L L Perini, D W Conn, et al. Calculation of carbon ablation on a re-entry body during supersonic/subsonic flight [J]. Journal of Spacecraft, 1986, 23(5): 487-492.

[15] 王臣, 梁軍, 吳世平, 等. 高溫?zé)g條件下C/C材料熱力耦合場模擬[J]. 復(fù)合材料學(xué)報, 2006, 23(5): 143-148. [Wang C, Liang J, Wu S P, et al. Numerical simulation of C/C composites coupled thermos-mechanical field under the condition of high temperatures and ablation[J]. Acta Materiae Compositae Sinica, 2006, 23(5): 143-148.]

[16] 王中原, 史金光. 超高速飛行彈箭氣動燒蝕數(shù)值模擬研究[J]. 南京理工大學(xué)學(xué)報, 2003, 27(5): 595-602. [Wang Z Y, Shi J G. The simulation of aerodynamic heating and ablation for a hypervelocity projectile [J]. Journal of Nanjing University of Science and Technology, 2003, 27(5): 595-602.].

[17] 王松柏. 耐熱燒蝕材料燒蝕和侵蝕速度的計算方法[J]. 兵工學(xué)報, 1996, 17(4): 374-377. [Wang Songbai. The calculated method of the property and velocity of erosion with hot preservative materials [J]. Acta Armamentarii, 1996, 17(4): 374-377]

[18] 楊瓊梁, 史曉鳴, 許斌, 等. 燒蝕防熱層與結(jié)構(gòu)耦合熱傳導(dǎo)分析的交替計算法[J]. 宇航學(xué)報, 2011, 32(8): 1854-1858. [Yang Q L, Shi X M, Xu B, et al. Alternative algorithm for heat transfer analysis of ablative protection layer coupled with structure [J]. Journal of Astronautics, 2011, 32(8): 1854-1858.]

[19] Michael J W, Graham V C, Deepak B. Data-parallel line relaxation method for the navier-stokes equations[J]. AIAA Journal, 1998, 36(9): 1603-1609.

[20] 姜貴慶, 劉連元. 高速氣流傳熱與燒蝕熱防護(hù)[M]. 北京：國防工業(yè)出版社, 2003.

[21] 劉志剛, 張巨生, 韓杰才. 炭基復(fù)合材料超高溫?zé)峄瘜W(xué)燒蝕產(chǎn)物分析[J]. 炭素, 2005, 124(4): 21-25.[Liu Z G, Zhang J S, Han J C. Results analyzed of super-high temperature thermal chemical ablation for carbon-based composite material[J]. CARBON, 2006, 124(4): 21-25.]

[22] R M Grabow, C O White. Surface roughness effects on nose tip ablation characteristics [J]. AIAA Journal, 1975, 13(5): 605-610.

[23] M B Ziering, V Di Cristina, Thermo-mechanical erosion of ablative plastic composites [C]. AIAA 7th Thermo-physics Conference, San Antonio, April 10-12, 1972.

[24] G J Hochrein, G F Wright. Analysis of the tater nose-tip boundary layer transition and ablation experiment [C]. AIAA 14th Aerospace Sciences Meeting, Washington D.C., 1976.

Research of Coupling Computational Method about the Nose-Tip Ablation of Reentry Vehicle

LONG Li-ping1HAN Jun2WAN Tian1ZHONG Feng-quan1,2TIAN Bao-wei2

（1 Institute of Mechanics, Chinese Academy of Sciences, Be?ing 100190, China；2 School of Engineering Science, University of Chinese Academy of Sciences, Be?ing 100049, China)

Nose-tip ablation is a very important problem in the reentry process of hypersonic vehicle. In the high temperature and high pressure environment, the original aerodynamic shape and structure heat transfer boundary are constantly changing because of the ablation, and the change of aerodynamic shape and heat transfer boundary affects nose-tip heat flux, temperature distribution and ablation quantity in turn. They exhibit complex strong coupling and non-linear characteristics. Taking the ablation process of nose-tip of carbon based materials as an example, the coupling calculation method of flow, ablation and structure heat transfer of nose-tip are developed. Through the coupling and iteration of the calculation program of the aerodynamic, ablation and structure thermal response, the real time moving boundary simulation of the nose-tip reentry ablation process is realized, and the results are basically in agreement with the measured data under the flight test conditions.

Reentry process; nose tip; ablation; moving boundary; coupling simulation

V416.5

A

1006-3919(2021)02-0008-07

10.19447/j.cnki.11-1773/v.2021.02.002

2020-12-03；

2021-01-05

中國科學(xué)院戰(zhàn)略性先導(dǎo)科技專項（XDA17030100）

龍麗平（1977—），女，博士，副研究員，研究方向：高超聲速飛行器結(jié)構(gòu)優(yōu)化設(shè)計，可變形高超飛行器總體優(yōu)化設(shè)計；（100190）北四環(huán)西路15號中國科學(xué)院力學(xué)研究所.