气动加热的数值仿真及其地面试验模拟技术

摘 要： 气动加热仿真与相关试验研究是高超声速飞行器发展中的重要内容。 对超声速气动加热现象中的边界层流动和传热过程进行分析， 从物理角度探讨了气动加热仿真方法中工程模型、 气体模型、 湍流模型和网格尺度， 以及计算格式和边界条件对仿真的影响。 介绍了气动加热研究常用的氧-乙炔火焰烧蚀、 辐射加热器、 高焓风洞等试验设备及其优缺点， 讨论了气动加热环境地面模拟试验中的相似参数， 并分析了冷壁热流过高的原因。

关键词： 气动加热; 数值仿真; 地面试验; 试验设备; 相似参数

中图分类号：       TJ763

文献标识码：     A

文章编号：      1673-5048（2023）03-0011-09

DOI：  10.12132/ISSN.1673-5048.2022.0253

0 引 言

飞行器在大气层内做超声速或高超声速飞行时， 由于飞行器与环境气体之间存在巨大动能差以及气体黏性的作用， 环境气体会对飞行器表面产生剧烈的加热作用， 这种现象称为气动加热。

近几十年来， 对高超声速气动加热的研究非常多， 包括流场的数值仿真、 流场与结构的耦合仿真、 热防护材料的烧蚀、 以及试验模拟等［1-4］。 为了使同行能够从浩如烟海的文献中快速找到关键问题并理清研究思路， 本文将高超声速气动加热的数值仿真和地面模拟试验研究进行归纳总结， 希望能够引起同行更加深入的思考。

1 超声速飞行时的气动加热现象

从本质上看， 超声速飞行器的气动加热问题其实是飞行器和大气之间的能量交换问题。 图1以高超声速飞行的弹头为例进行说明。

从固定在弹体上的坐标系来观察， 气流以超声速冲击弹体， 首先形成脱体激波。 因为气流经过激波时有能量耗散（非等熵流动）， 且激波强度越高耗散越剧烈， 所以气流的动能经过激波后会衰减。 气流经过激波后状态发生变化， 如温度、 密度、 压强升高， 速度降低。 改变状态后的气流在弹体表面形成边界层流动， 且边界层厚度不断变化。 在边界层内， 气体黏性形成的摩擦作用使气体温度进一步升高， 最终在黏性底层以导热的形式将能量传递给壁面， 形成壁面热流。

当来流速度足够高、 使得激波后的气体温度升高到一定程度时， 气体分子本身的振动能被激发， 气体就会发生电离反应， 由分子状态变成离子状态， 进而引起物理性质的变化和激波位置的改变， 从而影响最终的气动加热效果。

而对于复杂几何外形的飞行器， 飞行器局部部件引起的激波可能会与飞行器其他部位发生交汇， 激波之间或者激波对边界层等局部流动产生干扰， 进而影响这些部位的换热。

而且， 高超声速飞行器在其所处的飞行环境下，  边界层内也会出现边界转捩现象。 边界层转捩通常是指边界层内流动由层流状态发展为湍流状态的过程。 在转捩起始点处飞行器表面热流密度会有一个陡升， 造成气动加热加剧， 致使对转捩区域内气动加热更难以预测［5］。

根据气动加热的物理过程可知， 要准确得到最终施加给壁面的热流， 就需要准确描述气流经历的各个环节。 Wuster等详细描述了NASA兰利研究中心对再入飞行器HL-20进行气动热评估和热防护系统选型的过程［6］， 指出确定气动热环境的典型方法是首先确定飞行弹道， 然后针对飞行弹道中的关键状态进行气动热环境评估； 而气动热环境的评估则是从简单的热分析方法（工程方法）到逐渐复杂的CFD（计算流体力学）方法。 随着设计过程的深入， CFD的方法也更加细致， 同时还需要进行关键热环境状态的风洞试验验证。 尽管这篇文献较早， 但是其详细阐述了气动热防护系统研究的各个环节， 研究路径基本相同。

无论超声速远场的流动结构如何， 气动加热终究是发生在壁面。 图2示意了壁面边界层内的气流速度分布和温度幅值的变化情况， 其中， Tw=T0是假定壁面温度恒为T0（常取值300 K）的状态， 对应的热流称为冷壁热流， 而Tw>>T0是指壁面温度在持续加热条件下的状态， 对应的热流称为热壁热流。 从局部传热的角度看， 气体传向壁面的热流qw可用傅里叶导热定律来表示：

式中： Tnw是指壁面处的法向温度梯度； k（T）是壁面处气体的导热系数（即温度的函数）。 可见， 壁面受热后Tw会一直升高， 对应的热壁热流不断减小， 若流场保持稳定， 则Tw最终会达到稳定值。 但是， 即使Tw达到稳定值， 也与表面材料相关。 受结构材料导热系数的影响， 即不同的材料在相同的来流条件下会达到不同的平衡温度， 因此， 为了便于评估不同飞行条件下的气动热， 通常采用冷壁热流作为衡量参数， 其排除了材料的影响， 并且在计算时不必考虑流场和结构之间的耦合。

2 气动加热的数值仿真

从仿真的角度看， 气动热问题的求解就是如何准确地给出壁面处的温度梯度（见图2）， 因为其中既涉及到流动问题， 也涉及到高温气体参数和壁面传热的问题， 所以物理模型方面就包括流动控制方程、 气体模型（含化学反应模型）、 湍流模型及壁面网格等问题， 而数学方法则需要关注高精度、 低数值耗散计算格式， 以及边界精度问题和化学反应求解的效率问题。

2.1 工程模型

因为气动热的起因就是气体的黏性， 所以在气动热的数值模拟中应采用黏性流体的N-S方程组。 但即使在计算技术得到了很好发展的今天， 求解整个飞行器的黏性绕流仍然是很艰巨的任务。 因此， 在早期的气动热问题求解， 以及目前复杂外形飞行器的气动热问题求解过程中， 采用纯工程算法、 或者求解无黏流动（Euler方程组）并结合附面层理论或工程方法求解壁面热流， 仍然是一种非常实用的方法。

纯工程算法实际上是对边界层内参数求解的一种最简单的等效方法。 因为不同外形的结构所对应的边界层不同， 所以纯工程算法一般是把复杂的飞行器外形分成球体、 锥体、 后掠圆柱体和楔形体等部分， 分别利用这几种形体各自经过试验验证的表面热流公式进行计算［7-11］。 例如， 驻点热流的工程计算方法通常采用的Fay-Riddell公式， 是对驻点区域的高温气体边界层方程进行简化得到的， 适用于计算来流总焓在1 549～24 158 kJ/kg、 壁温在300～3 000 K之间的驻点热流。 公式如下：

式中：  Pr为普朗特数； Le为路易斯数； ρw， μw和ρs， μs分别为壁面以及驻点的密度和粘性系数； hD， hs， hw分别为离解焓、 驻点焓以及壁面焓。 这些工程算法经过试验的检验， 具有较高的精度。

但是对于较复杂的飞行器外形、 或者某些较难近似的飞行姿态， 纯工程算法针对某些部位的气动热计算偏差就会变大， 因为毕竟这些算法是针对特定部位和飞行姿态下的边界层而提出的。 因而， 当气动热不涉及到激波/边界层干扰或缝隙流动等黏性起主要作用的问题时， 采用求解无黏流动+边界层分析（或工程算法）的方式是较好的选择， 既可以满足精度的要求， 又可以在一定程度上满足快速性的要求。

这方面的较早研究可参考Riley的论文［12］， 其思路是先利用无黏流动计算获得边界层外缘的参数， 再利用这些外缘参数来分析边界层流动， 并计算出边界层的位移厚度， 然后考虑位移厚度对飞行器外形的影响， 再根据修正的外形重新进行边界层外缘参数计算， 最终用这些参数进行边界层内的传热分析， 获得较为准确的当地壁面热流。 图3给出了一个比较的算例［12］， 其中AEROHEAT是工程算法的软件， 可见， 通过无黏流动计算边界层外的参数再结合边界层分析壁面热流， 可以得到更符合试验数据的结果。

在缺乏可靠的计算资源（包括硬件和软件）时， 这种无黏流计算+边界层分析的方法非常具有实用价值， 因此得到了充分的发展和应用。 例如文献［13］所描述的HEAT2D程序将该方法与二维结构传热模型相耦合（见图4）， 可以分析飞行器表面材料在整个飞行过程中的受热以及温度变化， 并可获得足够的精度。

2.2 气体模型

由图2和式（1）可知， 准确计算气动热的前提就是准确计算壁面法向的温度分布曲线， 当然也包括边界层外的温度。 而在气体流动中， 影响温度的参数是比热容， 描述比热容的模型就是气体模型， 故气体模型在气动热计算中特别重要。

CFD中常用的气体模型有量热完全气体模型（比热容为常数）、 单一组分或多组分热完全气体模型（其中某一组分的定压比热是温度T的函数， 且cp-cv=R， R为气体常数）。 其中热完全气体模型的函数系数是由试验和理论确定的， 比较符合实际情况； 而量热完全气体模型的定比热容假设只适用于温度变化不大的情况， 用该模型计算气动热的偏差可高达44%［14］； 另外， 定压比热的变化对边界层内速度的影响相对较小， 而对温度（平均温度、 温度波动）的影响非常明显［15］， 因此在气动热计算中不应再使用量热完全气体模型。 文献［16］通过数值计算验证了热完全气体模型的有效性， 指出采用热完全气体模型可以提高高温、 高速流动状态下气动加热的计算精度。 当气流总温不太高时（如低于1 500 K）， 可以采用单一组分的热完全气体模型。 但是， 当飞行器的速度足够高， 以致于激波后气体温度升高导致激发了气体分子的振动能、 并使气体发生电离后， 就应该采用考虑电离反应或化学非平衡过程的多组分热完全气体模型［17］。 与（单一组分的）热完全气体模型相比， 考虑高温空气的化学非平衡特性以后， 一个最明显的特征就是仿真得到的脱体激波更靠近钝体， 而且波后气体的计算温度会更低一些。 这是因为气体的电离反应是吸能过程， 电离反应后气体组分的增加， 使得混合气体的平均分子量降低、 比热容增大、 比热比减小， 这必然会导致混合气体的温度更低一些。 从气动热计算的角度考虑， 尽管激波位置更接近壁面， 但是电离使混合气体温度降低， 这必然会使壁面热流减小。 图5给出了文献［17］的计算结果， 从数值模拟的角度给出了例证。

2.3 湍流模型与网格尺度

因为湍流模型对流动的影响主要体现在壁面边界层内以及旋涡运动过程中， 而气动加热的理化过程是发生在边界层内， 所以湍流模型的选择对于气动加热仿真来说尤其重要。 只有当仿真所采用的湍流模型能够准确反映边界层内的流体参数变化时， 才能保证壁面热流计算的准确度。 文献［18］关于不同边界层壁面（第一层）网格厚度条件下， 剪应力传输（SST）湍流模型对气动加热计算结果的影响， 也恰恰说明这一点。 另外， 不仅湍流模型本身对边界层壁面网格厚度有要求， 而且壁面传热过程也对网格厚度有要求， 所以需将湍流模型和壁面网格厚度一起讨论。

从图2可以看到， 在气动加热的情况下， 边界层内的温度曲线比速度轮廓线要复杂。 速度轮廓通常是单调变化的， 无论是层流状态还是湍流状态均是如此； 但是温度曲线则不同， 这是因为黏性作用（无论是分子黏性亦或是湍流脉动）使得气体的动能转化为内能， 故边界层内的气体温度在靠近壁面时， 必然逐渐升高到大于壁面温度， 再通过导热将热量传递给壁面。 显然， 当假设壁面温度Tw=T0（如300 K）， 即不考虑与结构的热耦合而只计算冷壁热流时， 壁面处的温度梯度更大。 对于数值计算而言， 参数梯度越大， 达到精度所需的网格尺度越小， 这就是很多关于气动热计算对壁面网格的尺度有要求的根本原因。 从这个角度来说， 如果不考虑耦合计算的复杂度， 那么计算冷壁热流比计算热壁热流还要困难（需要更小的壁面网格和更小的时间步长， 或更长的收敛时间）。

关于气动热计算的网格尺度研究表明， 壁面第一层网格应落在线性（层流）层内； 另外， 与压强（以声速）传播的过程相比， 传热过程本身就是慢过程， 所以气动热的收敛比气动力的收敛慢得多。 在气动热计算中， 应直接观察热流数据的收敛， 以确保得到真正的收敛解［19］。