刘宇^1,2, 赵淼^1,2, 张章^1,2, 贾贺^1,2, 黄伟^1,2
1. 北京空间机电研究所, 北京 100094;
2. 中国航天科技集团有限公司 航天进入减速与着陆技术实验室, 北京 100094
收稿日期：2021-12-18
基金项目：国家自然科学基金资助项目（11602018）
作者简介：刘宇(1986-), 男, 工程师。E-mail: 413908179@qq.com

摘要：该文对充气式再入与减速技术（inflatable reentry and descent technology，IRDT）中使用的锥形减速结构流场进行了气动热仿真计算。计算模型基于有限体积法对N-S（Navier-Stokes）方程进行求解，使用非平衡双温度模型计算流场热化学反应。为了验证算法准确性，对钝体标准模型ELECTRE进行仿真计算，计算结果与飞行试验和文献结果基本相符。锥形减速结构仿真工况高度为70 km，来流Ma为13，仿真结果表明：激波后振动温度被激活，并逐渐升高至平动温度，同时空气中离解组元浓度逐渐升高；结构表面热流与压强在驻点附近沿径向快速降低，随后热流呈线性下降，压强近似为常量。对4种不同半锥角的锥形减速结构仿真结果进行了对比，结果显示：50°、55°和60°半锥角激波位置及表面热流基本相同，65°半锥角激波距离前缘点更远，同时表面热流更低；4种半锥角驻点压强基本相同，外围压强随半锥角增加呈线性增加。仿真结果可为IRDT方案设计提供参考。
关键词：锥形减速结构充气式减速气动热非平衡反应双温度模型数值仿真
Simulation of thermochemical nonequilibrium flow around a conical deceleration structure
LIU Yu^1,2, ZHAO Miao^1,2, ZHANG Zhang^1,2, JIA He^1,2, HUANG Wei^1,2
1. Beijing Institute of Aerospace Mechanics and Electricity, Beijing 100094, China;
2. Laboratory of Aerospace Entry, Descent, and Landing Technology, China Aerospace Science and Technology Corporation, Beijing 100094, China

Abstract: Objective The conical deceleration structure is a typical shape in inflatable reentry and descent technology (IRDT). Compared with the traditional rigid deceleration structure, the inflatable deceleration structure represented by alumina fiber has lower heat resistance. Therefore, accurate thermal environment prediction is crucial for designing the IRDT system. Moreover, high pressure deforms the surface of the inflatable structure, so the surface pressure distribution is another issue that needs attention. The surface heat flux and pressure distribution of a conical deceleration structure under thermochemical reaction conditions are analyzed through numerical simulation. At the same time, the influence of different half-cone angles on the surface heat flow and pressure distribution is analyzed. Methods The numerical model is based on the integral Navier-Stokes (N-S) equation. The Park85 and the two-temperature nonequilibrium models are used to calculate the thermochemical reaction with a noncatalytic wall condition. The equations are solved using the finite volume method. The lower-upper symmetric Gauss-Seidel method is adopted for iteration. The blunt body standard model ELECTRE is used to validate the numerical model. The calculation case of a conical deceleration structure with a height of 70 km is investigated, and the inlet Mach number is 13. The variations in temperature and chemical component concentration along the stagnation line, as well as heat flow and pressure distributions on the structure surface are studied. In addition, the simulation of four conical deceleration structures with different half-cone angles is carried out to analyze the effect of the half-cone angle on the surface heat flow and pressure. Results The simulation results show that 1) the gas translational temperature after the shock wave is approximately 7000 K. Along the stagnation point line, the vibrational temperature gradually increases, and the two temperatures reach equilibrium near the stagnation point and decrease to the wall temperature. 2) The concentration of the N component in the shock layer is low and decreases to 0 at the stagnation point. The O and NO components gradually increase along the stagnation point line and reach the maximum near the stagnation point. 3) The surface heat flow and pressure are the highest at the stagnation point and decrease rapidly along the radial direction near the stagnation point. Then, the heat flow decreases linearly, and the pressure is approximately constant. 4) For different half-cone angle conical deceleration structures, the shock wave positions and surface heat flow distributions of the 50°, 55°, and 60° cases are basically identical. The shock wave position of the 65° case is farther from the leading edge, and the surface heat flux is lower. 5) Finally, the stagnation pressures of the four cases are basically identical, and the peripheral pressures increase linearly with increasing half-cone angle. Conclusions The surface heat flow and pressure distributions on the conical deceleration structure can be revealed by the numerical calculation. The change in the half-cone angle significantly impacts the surface heat flow and pressure distributions of the conical deceleration structure.
Key words: conical deceleration structureinflatable reentry and descent technologyaerothermodynamicsnonequilibrium reactiontwo-temperature modelnumerical simulation
充气式再入与减速技术(inflatable reentry and descent technology，IRDT)是一种新型的大气再入减速技术^[1]，该技术特点在于其减速结构是由耐高温柔性材料制成的充气式结构。与传统返回舱式减速技术相比，充气式减速技术具有重量轻、体积小、成本低等优点，拥有广泛的应用前景。充气减速结构作为IRDT系统的核心组成部分，有多种不同的外形设计方案^[2-4]，其中典型的充气减速结构为多个充气圆环组成的锥形结构^[5-6]。
在IRDT系统减速的过程中，会在弓形激波与锥形减速结构表面形成高温激波层，对外部热环境进行准确的预测是IRDT系统热防护设计的前提。以氧化铝纤维为代表的充气减速结构热防护材料的耐热性能要低于传统的返回舱热防护材料^[7]，这也对充气式减速结构热环境预测提出了更高的要求。对于航天器再入大气层的热环境问题，试验手段包括高空飞行试验^[8-9]和高焓风洞试验^[10-11]，高空飞行试验成本高且风险大，而高焓风洞试验数据的分析、外推和使用与常规风洞试验相比更加困难和复杂^[12]。目前高超声速流场气动热的计算方法主要有工程算法、工程数值耦合方法及数值仿真方法。工程算法^[13-15]借助经验公式求解驻点热流密度，计算效率高且可以保证一定的精度；工程数值耦合方法^[16]对边界层外无黏流场进行数值求解，对边界层内流场用工程算法求解，可应用于复杂外形；数值仿真方法主要包括求解连续流动的计算流体动力学(computational fluid dynamics, CFD)方法^[17-19]和求解稀薄流动的直接模拟Monte Carlo(direct simulation of Monte Carlo，DSMC)方法^[20]，虽然可以获得复杂流场的流动特征，但是计算量大，特别是在计算包含复杂的热化学反应时。
在高超音速流场中激波层温度非常高，分子内的振动能得到激发，同时会引起边界层内的气体出现离解甚至电离^[21]，因此想要得到正确的计算结果，数值仿真计算必须考虑激波层内的热化学反应。本文基于有限体积方法，采用双温非平衡热化学反应模型，对锥形减速结构的外部热环境进行仿真计算，分析了沿驻点线的温度变化与化学组元浓度变化，以及结构表面的热流及压强分布；同时，对4种不同半锥角的锥形减速结构进行仿真计算，对减速结构外部热流及压强受半锥角变化的影响进行研究。
1 计算方法1.1 非平衡热化学反应N-S方程气体的内能模式包含分子的平动能、转动能、振动能及电子能^[21]。在低温下，只有分子的平动能和转动能处于激发状态并达到平衡，气体温度才可以用单一平动温度T进行描述。在高温气体中，振动能及电子能被激发，其能量传递需要用振动温度T_v及电子温度T_e进行描述。在双温模型^[22]中，原子振动能及电子能处于平衡状态，其能量传递可以统一用T_v进行描述。基于双温模型，包含N种不同化学组元的气体在非平衡状态下的N-S(Navier-Stokes)方程积分形式可表示为

$\frac{\partial}{\partial t} \int_V \boldsymbol{Q} \mathrm{d} \sigma+\oint_{\partial V}\left(\boldsymbol{F}_{\mathrm{c}}-\boldsymbol{F}_{\mathrm{v}}\right) \mathrm{d} s=\int_V \boldsymbol{S} \mathrm{d} \sigma .$

(1)

其中：F_c和F_v分别为对流通量及黏性通量，S为源项，t为时间，V和?V分别为控制体体积及表面积，dσ和ds分别为体积微元和面积微元，Q为流场守恒变量。Q可表示如下：

$\boldsymbol{Q}=\left[\begin{array}{c}\rho \\\rho_1 \\\vdots \\\rho_{N-1} \\\rho u \\\rho v \\\rho w \\E \\E_{\mathrm{v}}\end{array}\right], $

(2)

$E=\sum\limits_{i=1}^N \rho_i h_i+\frac{1}{2} \rho|\boldsymbol{V}|^2, $

(3)

$h_i=\left[\frac{H_i^0(T)-R}{M_i}\right]+e_{\mathrm{v}, i}, $

(4)

$E_{\mathrm{v}}=\sum\limits_{i=1}^N \rho_i e_{\mathrm{v}, i}, $

(5)

$e_{\mathrm{v}, i}=\frac{R}{M_i} \frac{\theta_i}{e^{\theta_i / T_{\mathrm{v}}}-1} .$

(6)

其中：ρ为气体密度，ρ_i为第i种气体组元密度，i=1，2，…，N-1；速度矢量V=(u, v, w)，u、v和w为速度分量；E为气体的总内能，E_v为气体的振动内能；h_i为单位质量气体组元的焓；H_i⁰(T)为每摩尔气体组元在温度T下的焓值，可通过文[23-24]获得；R为理想气体常数；M_i为气体组元相对分子质量；e_{v, i}为单位质量气体组元的振动能；θ_i为气体组元振动特征温度。
对流通量F_c可表示如下：

$\boldsymbol{F}_{\mathrm{c}}=\left[\begin{array}{c}\rho V_n \\\rho_1 V_n \\\vdots \\\rho_{N-1} V_n \\\rho u V_n+p n_x \\\rho v V_n+p n_x \\\rho w V_n+p n_y \\(E+p) V_n \\E_{\mathrm{v}} V_n\end{array}\right] .$

(7)

其中：P为气体压强；n(n_x, n_y, n_z)为控制体表面?V法向量，x、y、z分别为3个坐标轴方向；V_n=V·n为流体速度在n方向投影，n为法向量n的方向。
扩散通量F_v可表示如下：

$\boldsymbol{F}_{\mathrm{v}}=\left[\begin{array}{c}0 \\D_1 \frac{\partial \rho_1}{\partial n} \\\vdots \\D_{N-1} \frac{\partial \rho_{N-1}}{\partial n} \\\tau_{n x} \\\tau_{n y} \\\tau_{n z} \\u \tau_{n x}+v \tau_{n y}+w \tau_{n z}+ \\q_n+q_{\mathrm{v} n}+\sum\limits_{i=1}^N D_i h_i \frac{\partial \rho_i}{\partial n} \\q_{\mathrm{v} n}+\sum\limits_{i=1}^M D_i e_{\mathrm{v}, i} \frac{\partial \rho_i}{\partial n}\end{array}\right], $

(8)

$\boldsymbol{\tau}_n=\boldsymbol{\tau} \cdot \boldsymbol{n}=2 \mu\left[{\mathit{\pmb{ε}}}-\frac{1}{3}(\nabla \cdot \boldsymbol{V}) \boldsymbol{I}\right] \cdot \boldsymbol{n}, $

(9)

${\mathit{\pmb{ε}}}=\frac{1}{2}\left[(\nabla \cdot \boldsymbol{V})+(\nabla \cdot \boldsymbol{V})^{\mathrm{T}}\right], $

(10)

$q_n=k \frac{\partial T}{\partial n}, $

(11)

$q_{\mathrm{v} n}=k_{\mathrm{v}} \frac{\partial T_{\mathrm{v}}}{\partial n}.$

(12)

其中：D_i为气体组元的质量扩散系数，假设各气体组元的二元扩散系数^[25]相同，则D_i=μ/ρSc，其中，Sc为Schimidt数，对于中性粒子组元可取0.5；τ为切应力张量, τ_n为τ在n方向投影，τ_n=(τ_nx, τ_ny, τ_nz)；μ为气体黏性系数；I为单位张量；ε为应变率张量；q_n和q_vn分别为气体平动和振动热流；k和k_v分别为气体平动导热系数和振动导热系数；对于混合气体，k、k_v和黏性系数μ可由Wilke^[26]半经验方法通过各组元的黏性系数和导热系数确定。
源项S表示如下：

$\boldsymbol{S}=\left[\begin{array}{c}0 \\\dot{\omega}_1 \\\vdots \\\dot{\omega}_i \\0 \\0 \\0 \\0 \\\sum\limits_{i=1}^N\left(\rho_i q_{\mathrm{T}-\mathrm{v}, i}+\dot{\omega}_i e_{\mathrm{v}, i}\right)\end{array}\right], $

(13)

$\dot{\omega}_i=M_i\left(\nu_i^{\prime \prime}-\nu_i^{\prime}\right)\left(k_{\mathrm{f}} \prod\limits_{i=1}^N\left(\frac{\rho_i}{M_i}\right)^{\nu_i^{\prime}}-k_{\mathrm{b}} \prod\limits_{i=1}^N\left(\frac{\rho_i}{M_i}\right)^{\nu_i^{\prime \prime}}\right) .$

(14)

其中：$\dot{\omega}$_i为组元生成速率；q_{T-V, i}为平动能-振动能转化速率；ν'_i和ν″_i分别为反应物和生成物的化学反应式计量系数；k_f和k_b分别为正向及逆向反应速率，可由Arrhenius公式^[27]表示，两者形式相同，以k_f为例，

$k_{\mathrm{f}}=A T^B \exp (-C / T) .$

(15)

其中：A、B和C为反应常数，可通过实验获得。k_f和k_b满足关系

$\frac{k_{\mathrm{f}}}{k_{\mathrm{b}}}=K_{\mathrm{eq}}=\exp \left(-\frac{\Delta_{\mathrm{f}} G^0}{R T}\right) .$

(16)

其中：K_eq为反应平衡常数，Δ_fG⁰为正向反应过程中Gibbs自由能的变化量^[21]。当k_f为已知时，可通过式(16)计算k_b。
q_{T－V, i}可通过Landau-Teller公式^[28]表示，即

$q_{\mathrm{T}-\mathrm{v}, i}=\frac{e_{\mathrm{v}, i}^*-e_{\mathrm{v}, i}}{\tau_{\mathrm{v}, i}} .$

(17)

其中：e_{v, i}^*为平衡振动能，τ_{v, i}为振动松弛时间。
本文采用有限体积法对式(1)进行离散求解，迭代采用LU-SGS方法，对流通量采用Roe通量分裂格式，通过MUSCL方法结合Van Leer限制器^[29]进行变量重构获得二阶精度，黏性通量采用中心格式。
1.2 热化学反应模型对于空气中不同组元热化学反应式的建立及式(15)中各化学反应常数的确定，有多种不同的模型供选择。常用的热化学反应模型包括Park85^[30]、Park93^[31]及Dunn-Kang^[32]模型等。其中，Park85和Park93模型只给出了正向反应的参数值，逆向反应需通过反应平衡常数求得，而Dunn-Kang同时给出了正向及逆向反应的参数值。在实际应用中，为提高计算效率，一般对上述化学反应模型组元进行简化。在工程计算满足精度需求的情况下，5组元Park85模型可以兼顾计算的准确度和效率^[33-34]。本文计算采用基于Park85模型的5组元模型，相应的反应式及反应参数如表 1所示。
表 1 5组元化学反应模型^[30]

化学反应	A/(m3·mol-1·s-1)	B	C/K
N2+M$\rightleftharpoons$2N+M	3.70×1015(M=N2, O2, NO) 1.11×1016(M=N, O)	-1.60	113 200
O2+M$\rightleftharpoons$2O+M	2.75×1013(M=N2, O2, NO) 8.25×1013(M=N, O)	-1.00	59 500
NO+M$\rightleftharpoons$N+O+M	2.30×1011(M=N2, O2, NO) 4.60×1011(M=N, O)	-0.5	75 500
N2+O$\rightleftharpoons$NO+N	2.16×102	1.29	19 220
NO+O$\rightleftharpoons$O2+N	3.18×107	0.10	37 700

表选项






热化学反应壁面边界条件可分为非催化壁及完全催化壁^[21]，本文采用非催化壁，即壁面化学反应速率为0。
2 算法验证为了验证计算方法的准确性，本文首先对钝体标准模型ELECTRE进行仿真计算，并将结果与文[11]飞行试验数据及计算结果进行对比。ELECTRE模型为钝锥形旋成体，总长度为2 m，球形头部半径R为0.175 m，半锥角为4.6°，如图 1所示。流场环境参数与文献中飞行试验一致，即高度为53.3 km，来流速度为4 230 m/s，温度为265 K，压强为53 Pa，Ma为13，密度为6.94×10^-4 kg/m³，壁面温度为343 K。计算网格如图 2所示，基于1/4外形划分结构网格，网格总数为36万，其中壁面法向网格节点数为80。为验证计算结果的网格相关性，壁面网格法向尺寸分别设置为1.7×10^-5 m和1.7×10^-4 m，对应网格Reynolds数Re_c分别为3和30。

图 1 ELECTRE模型外形

图选项

图 2 ELECTRE计算网格

图选项

部分变量的归一化残差迭代收敛曲线如图 3所示，各变量残差值均收敛于稳定结果。壁面热流计算结果与飞行试验数据及文[11]计算值对比如图 4所示，X即图 1中x方向坐标值, 其中文[11]使用模型为5组分Dunn-Kang模型。两套网格计算结果整体相差在3%以内，可认为计算结果满足网格无关性。在驻点附近本文结果与试验及文[11]结果一致，在流场下游3者变化趋势相同，本文计算结果低于文[11]及试验值。造成结果偏差的主要原因为两者选取的化学反应模型不同。在飞行试验中壁面催化效应及高温气体的热辐射效应都会导致壁面热流增加，但在计算中没有考虑这两方面因素。

图 3 残差收敛曲线(Rec=30)

图选项

图 4 壁面热流对比

图选项

沿驻点线平动温度和振动温度计算结果与文[11]对比，如图 5所示，X即图 1中x方向坐标值。定义组元浓度Y_i=ρ_i/ρ，沿驻点线组元O及NO浓度与文[11]对比如图 6所示，X为图 1中x方向坐标值。本文激波脱体距离与文[11]结果一致，激波后平动温度T略低于文[11]结果。T与T_v变化趋势及组元浓度变化趋势一致，因选取的化学反应模型不同，本文激波后组元离解速度更高，导致温度计算结果低于文[11]值。本文与文[11]计算结果对比及误差如表 2所示。

图 5 沿驻点线温度对比

图选项

图 6 沿驻点线组元浓度对比

图选项

表 2 本文与文[11]计算结果对比

	文[11]	本文	误差/%
激波脱体距离/m	0.015	0.015	0
激波后温度/K	7 997	7 572	5.3
NO组元壁面浓度	0.072	0.055	23.6
O组元壁面浓度	0.153	0.172	12.4

表选项






3 锥形减速结构流场仿真本文计算所用锥形减速结构如图 7所示，R为0.48 m，半锥角为60°，底边半径为2.3 m。流场环境为高度70 km，来流速度3 861 m/s，压强5.2 Pa，温度220 K, Ma为13，密度8.28×10^-5 kg/m³，壁面温度1 600 K。锥形减速结构计算网格如图 8所示，基于1/4外形划分结构网格，网格总数为21.6万，其中壁面法向网格节点数为80，壁面网格法向尺寸为6.0×10^-4 m，网格Re_c为3。

图 7 锥形减速结构

图选项

图 8 锥形减速结构计算网格

图选项

沿驻点线平动温度及振动温度变化曲线如图 9所示, X即图 7中x方向坐标值。激波距离前缘约0.05 m，激波后温度约7 000 K，激波后空气振动能被激活，振动温度逐渐升高，在前缘附近平动温度与振动温度达到平衡。沿驻点线组元O、NO及N浓度变化曲线如图 10所示, X即图 7中x方向坐标值。激波后空气中N₂及O₂发生离解反应生成O及N，进而通过置换反应生成NO。因N₂离解反应活化温度(即表 1中系数C)较高，N组元浓度较低且在驻点处几乎完全复合为N₂。O及N组元浓度持续升高，到驻点附近分别达到最大值。

图 9 沿驻点线温度

图选项

图 10 沿驻点线组元浓度

图选项

在IRDT方案设计中，柔性结构热响应及受力变形分析需要以表面热流及压强分布作为边界输入条件。锥体表面热流及压强沿展向变化曲线如图 11所示, Y即图 7中y方向坐标值。在驻点处热流及压强最大值分别为84 kW/m²和1 140 Pa，在驻点附近热流及压强沿展向快速降低，在距离驻点展向距离0.5 m以后，热流呈线性缓慢下降，压强近似为常量。在实际应用中，锥形结构中心部位设计为刚性防热结构，并与外部柔性结构相连接，根据上述计算结果可知，当柔性结构内缘尺寸满足条件时，其表面热流及压强分布可以认为是线性，柔性结构热响应及受力变形分析边界输入条件可以得到简化。

图 11 沿展向表面热流及压强

图选项

为了分析锥形减速结构半锥角变化对表面热流及压强分布的影响，本文比较了50°、55°、60°及65°共4种半锥角仿真结果，不同半锥角下沿驻点线平动温度变化对比如图 12所示。50°、55°及60°半锥角激波位置、激波后最高温度及温度变化基本相同；相比前3种半锥角，65°半锥角激波脱体距离增加了0.027 m，激波后最高温度略低，约为6 870 K。

图 12 不同半锥角沿驻点线平动温度对比

图选项

不同半锥角沿展向表面热流对比如图 13所示。50°、55°及60° 3种半锥角表面热流结果相近，当半锥角增加到65°以后，表面热流出现整体下降，幅度约为10 kW/m²。激波脱离距离对表面热流产生影响，一方面因为激波层厚度增加使温度梯度降低，另一方面激波层内反应时间增加导致空气离解更加充分。不同半锥角沿驻点线组元浓度对比如图 14所示，65°半锥角的2种离解气体组元浓度均高于60°半锥角。

图 13 不同半锥角沿展向表面热流对比

图选项

图 14 不同半锥角沿驻点线组元浓度对比

图选项

不同半锥角沿展向表面压强对比如图 15所示。除65°半锥角略高外，其他3种半锥角驻点压强基本相等。在锥体外围表面压强随半锥角增大而升高，分别为747、856、962和1 079 Pa，表面压强与半锥角变化近似为线性关系。

图 15 不同半锥角沿展向表面压强对比

图选项

4 结论本文对IRDT方案中使用的锥形减速结构流场进行了化学非平衡气动热仿真计算。计算模型基于有限体积法对N-S方程进行求解，利用5组元非平衡双温度模型计算流场热化学反应，壁面化学反应条件为无催化条件。对钝体标准模型ELECTRE进了仿真计算，仿真结果与飞行试验及文献计算结果基本相符。
对半锥角为60°的锥形减速结构在70 km高度、来流Ma为13工况下的流场进行了仿真计算，结果显示：激波后气体平动温度约为7 000 K，气体振动能被激活后沿驻点线平动能向振动能转化，振动温度逐渐升高，平动温度和振动温度在驻点附近达到平衡并降低至壁面温度。激波层内N组元浓度较低，并在驻点处重新降低为0；O及NO组元沿驻点线逐渐升高，在驻点处达到最大值；表面热流及压强在驻点处为最大值，在驻点附近沿径向快速降低，随后热流呈线性下降，压强近似为常量。
对50°、55°、60°及65°共4种半锥角锥形减速结构仿真结果进行对比，结果表明：50°、55°和60°半锥角激波位置及表面热流基本相同；相比前3种半锥角，65°半锥角激波距离前缘点更远，同时表面热流更低；4种半锥角驻点压强基本相同，外围压强随半锥角增加呈线性增加。

参考文献

[1]	黄伟, 曹旭, 张章. 充气式进入减速技术的发展[J]. 航天返回与遥感, 2019, 40(2): 14-24. HUANG W, CAO X, ZHANG Z. The development of inflatable entry decelerator technology[J]. Spacecraft Recovery & Remote Sensing, 2019, 40(2): 14-24. (in Chinese)
[2]	REZA S, HUND R, KUSTAS F, et al. Aerocapture inflatable decelerator (AID) for planetary entry[C]//Proceedings of the 19th AIAA Aerodynamic Decelerator Systems Technology Conference and Seminar. Williamsburg, USA: AIAA, 2007: 1-18.
[3]	卫剑征, 谭惠丰, 王伟志, 等. 充气式再入减速器研究最新进展[J]. 宇航学报, 2013, 34(7): 881-890. WEI J Z, TAN H F, WANG W Z, et al. New trends in inflatable re-entry aeroshell[J]. Journal of Astronautics, 2013, 34(7): 881-890. (in Chinese)
[4]	WALTHER S, THAETER J, REIMERS C, et al. New space application opportunities based on the inflatable reentry & descent technology (IRDT)[C]// Proceedings of AIAA/ICAS International Air and Space Symposium and Exposition: The Next 100 Years. Dayton, USA: AIAA, 2003: 1-7.
[5]	OLDS A, BECK R E, BOSE D M, et al. IRVE-3 post-flight reconstruction[C]// Proceedings of AIAA Aerodynamic Decelerator Systems (ADS) Conference. Daytona Beach, USA: AIAA, 2013: 1-24.
[6]	曹旭, 廖航, 许望晶, 等. 充气式减速技术试验器的设计和飞行试验[J]. 载人航天, 2018, 24(6): 802-808. CAO X, LIAO H, XU W J, et al. Design and flight test of demonstration aircraft with inflatable reentry and descent technology[J]. Manned Spaceflight, 2018, 24(6): 802-808. (in Chinese)
[7]	黄明星, 王伟志. 某型充气式再入减速热防护结构优化分析[J]. 航天返回与遥感, 2016, 37(1): 22-31. HAUNG M X, WANG W Z. Optimization on a flexible thermal protection structure of inflatable reentry system[J]. Spacecraft Recovery & Remote Sensing, 2016, 37(1): 22-31. (in Chinese)
[8]	LEE D B, BERTIN J J, GOODRICH W D. Heat-transfer rate and pressure measurements obtained during Apollo orbital entries[R]. Washington DC: National Aeronautics and Space Ministration, 1970.
[9]	KIMMEL R L, ADAMCZAK D W, BORG M P, et al. First and fifth hypersonic international flight research experimentation's flight and ground tests[J]. Journal of Spacecraft and Rockets, 2019, 56(2): 421-431. DOI:10.2514/1.A34287
[10]	MUYLAERT J, WALPOT L, SPEL M, et al. A review of European code validation studies in high enthalpy flow[C]// Proceedings of the 20th AIAA Advanced Measurement and Ground Testing Technology Conference. Albuquerque, USA: AIAA, 1998: 1-15.
[11]	MUYLAERT J, WALPOT L, H?USER J, et al. Standard model testing in the European High Enthalpy Facility F4 and extrapolation to flight[C]// Proceedings of the 28th Joint Propulsion Conference and Exhibit. Nashville, USA: AIAA, 1992: 1-16.
[12]	董维中, 乐嘉陵, 高铁锁. 钝体标模高焓风洞试验和飞行试验相关性的数值分析[J]. 流体力学实验与测量, 2002, 16(2): 1-8, 20. DONG W Z, LE J L, GAO T S. Numerical analysis for correlation of standard model testing in high enthalpy facility and flight test[J]. Experiments and Measurements in Fluid Mechanics, 2002, 16(2): 1-8, 20. (in Chinese)
[13]	HUGHES S, DILLMAN R A, STARR B R, et al. Inflatable re-entry vehicle experiment (IRVE) design overview[C]// Proceedings of the 18th AIAA Aerodynamic Decelerator Systems Technology Conference and Seminar. Munich, Germany: AIAA, 2005: 1-14.
[14]	高艺航, 贺卫亮. 充气式返回舱气动热特性研究[J]. 航天返回与遥感, 2014, 35(4): 17-25. GAO Y H, HE W L. Research on aerodynamic heating characteristics of inflatable reentry decelerator[J]. Spacecraft Recovery & Remote Sensing, 2014, 35(4): 17-25. (in Chinese)
[15]	周印佳, 张志贤, 付新卫, 等. 再入飞行器烧蚀热防护一体化计算方法[J]. 航空学报, 2021, 42(7): 124520. ZHOU Y J, ZHANG Z X, FU X W, et al. Integrated computing method for ablative thermal protection system of reentry vehicles[J]. Acta Aeronautica et Astronautica Sinica, 2021, 42(7): 124520. (in Chinese)
[16]	韩东, 方磊. 基于流线跟踪法的气动热工程计算研究[J]. 航空动力学报, 2009, 24(1): 65-69. HAN D, FANG L. Study of aerodynamic heating predictions for hypersonic aircrafts by tracking the surface stream trace[J]. Journal of Aerospace Power, 2009, 24(1): 65-69. (in Chinese)
[17]	粟斯尧, 石义雷, 柳森, 等. 有限催化对返回舱气动热环境影响[J]. 空气动力学学报, 2018, 36(5): 878-884. SU S Y, SHI Y L, LIU S, et al. Finite-rate surface catalysis effects on aero-heating environment of a reentry capsule[J]. Acta Aerodynamica Sinica, 2018, 36(5): 878-884. (in Chinese)
[18]	KOIKE T, TAKAHASHI Y, OSHIMA N, et al. Aerodynamic heating prediction of flare-type membrane inflatable reentry vehicle from low earth orbit[C]// Proceedings of 2018 AIAA Atmospheric Flight Mechanics Conference. Kissimmee, USA: AIAA, 2018: 1-21.
[19]	MAZAHERI A. High-energy atmospheric reentry test aerothermodynamic analysis[J]. Journal of Spacecraft and Rockets, 2013, 50(2): 270-281.
[20]	PALHARINI R C, AZEVEDO J L F, WHITE C. DSMC computations of SARA reentry capsule exposed to weakly ionized gas flow[C]// Proceedings of AIAA Scitech 2019 Forum. San Diego, USA: AIAA, 2019: 1-19.
[21]	卞荫贵, 徐立功. 气动热力学[M]. 2版. 合肥: 中国科学技术大学出版社, 2011. BIAN Y G, XU L G. Aerothermodynamics[M]. 2nd ed. Hefei: Press of University of Science and Technology of China, 2011. (in Chinese)
[22]	PARK C. Assessment of two-temperature kinetic model for ionizing air[C]// Proceedings of the 22nd Thermophysics Conference. Honolulu, USA: AIAA, 1987: 1-13.
[23]	MCBRIDE B J, ZEHE M J, GORDON S. NASA Glenn coefficients for calculating thermodynamic properties of individual species[R]. Cleveland: Glenn Research Center, 2002.
[24]	EHLERS J, GORDON S, HEIMEL S, et al. Thermodynamic properties to 6 000 K for 210 substances[R]. Cleveland: Lewis Research Center, 1963.
[25]	YOON B, RASMUSSEN M L. Diffusion effects in hypersonic flows with a ternary mixture[J]. KSME International Journal, 1999, 13(5): 432-442.
[26]	WILKE C R. A viscosity equation for gas mixtures[J]. The Journal of Chemical Physics, 1950, 18(4): 517-519.
[27]	FORD D I, JOHNSON R E. Dependence of rate constants on vibrational temperatures: An arrhenius description[C]// Proceedings of the 26th Aerospace Sciences Meeting. Reno, USA: AIAA, 1988: 1-5.
[28]	RIBEIRO V G, ZABADAL J G, MONTICELLI C O, et al. Estimating heat transfer coefficients for solid-gas interfaces using the Landau-Teller model[J]. Applied Mathematics and Computation, 2017, 301: 135-139.
[29]	ABU TALIP M S, AKAMINE T, OSANA Y, et al. Cost effective implementation of flux limiter functions using partial reconfiguration[C]// Proceedings of the 8th International Symposium on Applied Reconfigurable Computing. Hong Kong, China: Springer, 2012: 215-226.
[30]	PARK C. On convergence of computation of chemically reacting flows[C]// Proceedings of the Space Sciences Meeting. Reno, USA: AIAA, 1985: 1-17.
[31]	PARK C. Review of chemical-kinetic problems of future NASA missions, Ⅰ: Earth entries[J]. Journal of Thermophysics and Heat Transfer, 1993, 7(3): 385-398.
[32]	DUNN M G, KANG S. Theoretical and experimental studies of reentry plasmas[R]. Buffalo: Cornell Aeronautical Laboratory, 1973.
[33]	杨建龙, 刘猛, 阿嵘. 高超声速热化学非平衡对气动热环境影响[J]. 北京航空航天大学学报, 2017, 43(10): 2063-2072. YANG J L, LIU M, A R. Influence of hypersonic thermo-chemical non-equilibrium on aerodynamic thermal environments[J]. Journal of Beijing University of Aeronautics and Astronautics, 2017, 43(10): 2063-2072. (in Chinese)
[34]	张敏捷, 向树红. 高超声速三维热化学非平衡流场的数值计算对比研究[J]. 航天器环境工程, 2016, 33(1): 35-41. ZHANG M J, XIANG S H. A comparative study of the computation of 3D hypersonic flow in thermochemical nonequilibrium state[J]. Spacecraft Environment Engineering, 2016, 33(1): 35-41. (in Chinese)