Kriging气动力/热代理模型与热防护壁板热气动弹性分析

易子钧，冀春秀，谢丹，唐硕; YI Zijun; JI Chunxiu; XIE Dan; TANG Shuo

doi: 10.11887/j.issn.1001-2486.25010025

易子钧，冀春秀，谢丹，唐硕

西北工业大学航天学院, 陕西西安 710072

基金项目: 国家自然科学基金面上基金资助项目（52572457）；西北工业大学博士论文创新基金资助项目（CX2025037）

详细信息

作者简介

易子钧（2000—），男，河南郑州人，博士研究生，E-mail: yizijun@mail.nwpu.edu.cn

通讯作者

谢丹（1986—），女，陕西西安人，教授，博士，博士生导师，E-mail: dxie@nwpu.edu.cn

中图分类号: V211.47

文献标识码: A

文章编号: 1001-2486(2026)01-139-11

Kriging-based aerodynamic/aerothermal surrogate models and aerothermoelastic analysis of TPS panel

YI Zijun ， JI Chunxiu ， XIE Dan ， TANG Shuo

School of Astronautics, Northwestern Polytechnical University, Xi′an 710072 , China

摘要

针对高速气动力/热不同计算方法存在的效率和精度之间的矛盾，本文以典型高速飞行器的热防护壁板为研究对象，基于Kriging方法建立了气动力/热代理模型，将计算效率提升4个数量级。基于气动力/热代理模型，并采用有限元方法和自编热传导程序搭建了热防护壁板的热气动弹性计算框架，开展了热防护壁板的热气动弹性分析。该研究将为高速飞行器的气动力/热载荷快速准确计算、热防护系统设计和飞行安全评估提供重要理论基础。

关键词

高速飞行器 / 热防护壁板 / 热气动弹性 / 气动力/热代理模型 / Kriging代理模型

Abstract

To address the contradiction between efficiency and accuracy that exists among different computational methods for high-speed aerodynamic/thermal loads, a typical TPS (thermal protection system) panel of a high-speed vehicle was focused, and aerodynamic and aerothermal surrogate models based on the Kriging method were developed, which achieved a four-orders-of-magnitude improvement in computational efficiency. Based on these surrogate models, a computational framework for the aerothermoelastic analysis of the TPS panel was established using the finite element method and a self-developed heat conduction program. The aerothermoelastic behavior of the TPS panel was then analyzed within this framework. This research will provide an important theoretical foundation for the rapid and accurate prediction of aerodynamic and thermal loads, the design of thermal protection systems, and the flight safety assessment of high-speed vehicles.

Keywords

high-speed vehicles / TPS panels / aerothermoelasticity / aerodynamic/aerothermal surrogate model / Kriging surrogate model

1 研究对象与研究方法 1.1 样本数据获取方法 1.2 Kriging代理模型 1.3 热气动弹性耦合计算方法 2 基于Kriging的壁板气动力/热代理模型 2.1 气动力/热样本数据计算方法验证 2.2 气动力/热代理模型的建立与评估方法 2.2.1 状态空间参数化 2.2.2 代理模型的生成输入与预测输出 2.2.3 代理模型的精度评估方法 2.3 基于Kriging的气动力代理模型 2.3.1 气动力代理模型精度评估 2.3.2 非定常气动力代理模型 2.3.3 非定常气动力模型验证 2.4 基于Kriging的气动热代理模型 2.4.1 基于遗传算法的气动热代理模型参数优化 2.4.2 气动热代理模型精度评估 3 基于代理模型的热防护壁板热气动弹性分析 3.1 热屈曲 3.2 极限环运动 3.3 倍周期运动 3.4 准周期运动 3.5 混沌运动 3.6 结果分析 4 结论

高速飞行器作为新的战略制高点，在科研和国防方面有广阔应用前景，成为各国研究热点^[1-2]。其研制和装备涉及多领域、多学科的交叉前沿技术。高速飞行器与大气的剧烈摩擦，同时强激波下的高压力梯度导致气体温度急剧上升，产生气动加热效应，空气流动状态和物理性质剧烈改变^[3]，气动力/热为研究中的重点。另外，高速流中飞行器的气动热会引发材料温变与温结构热应力两方面影响，因此在流固耦合问题的基础上，由流-固-热三个物理场耦合引发的热气动弹性问题越发不可忽视。

典型的气动力/热计算方法包括实验测量、工程算法和数值计算。然而工程算法效率高但精度不足，数值计算具有较高精度但计算资源消耗大，为解决这种矛盾，发展出模型降阶方法^[4-6]，以建立兼有精度和计算效率的气动力/热模型。

典型的基于系统辨识的气动力降阶模型方法有Volterra级数法、Kriging模型、人工智能模型、模块式模型等^[7]。Kriging模型作为一种典型非线性系统辨识和插值方法，同时描述系统全局特性和局部细节。研究表明，使用Kriging模型可以建立高效、准确的气动力降阶模型，应用于气动弹性问题研究^[8]。该方法参数均由数学方法回归得到，无须参数调整和长时间数据训练。全局气动热降阶方法模型^[9]直接以气动热分布为输出。王洋等^[10]使用本征正交分解（proper orthogonal decomposition，POD）结合Kriging代理模型建立快速求解气动热参数的降阶模型。

本文基于Kriging代理模型，建立高速来流中弹性壁板的气动力/热代理模型。结合有限元方法构建气动弹性子系统，编写结构热传导程序构建气动热子系统，搭建热防护（thermal protection system，TPS）壁板热气动弹性计算框架，使用气动力/热代理模型开展热防护壁板的热气动弹性计算与分析。

1 研究对象与研究方法

典型高速飞行器多采用超燃冲压发动机及前体-进气道-推进等诸多一体化设计，存在诸多薄壁结构，其热气动弹性效应十分显著。本文针对典型高速飞行器前体-进气道腹板处的热防护壁板作为研究对象，飞行器外形简化为如图1所示的高速流场中的半顶角θ=5°的尖楔形结构，飞行高度为30 km。在高速来流（1区）下，尖楔前缘产生斜激波，波后气流（2区）平行于尖楔表面流动，弹性壁板前缘（3区）的气体参数与2区一致。2区参数通过斜激波理论得到。

尖楔结构前体表面的弹性板如图1所示，两端简支，壁板长度为a，壁板厚度为h。此外，假定气流在壁板前缘上游方向1 m处发生层流到湍流的转捩^[11]。

图1热防护壁板

Fig.1TPS panel

本文所研究的弹性壁板为飞行器表面具有热防护功能的热防护壁板。由上中下三层构成，最上层为PM-2000蜂窝夹层的辐射层，中间层为内部多屏隔热（internal multiscreen insulation，IMI）层，最下层为钛合金（Ti-6Al-2Sn-4Zr-2Mo）金属板，如图1所示。

热防护壁板三层结构之间为理想热接触，没有热量损失。另外，由于辐射层为很薄的金属镀层，隔热层为稀松的复合材料，忽略其刚度影响。热防护壁板相关参数见表1。

表1热防护壁板相关参数

Tab.1Relevant parameters of the TPS panel

材料参数随温度变化显著，热防护壁板各层材料比热c、材料导热系数k以及金属板弹性模量E和热膨胀系数α等材料属性随温度的变化见文献^[12-14]，随温度升高，弹性模量降低，表现为材料温变。

1.1 样本数据获取方法

本文使用计算流体力学（computational fluid dynamics，CFD）方法作为生成气动力/热样本数据的方法。

本文采用三维可压缩非定常Navier-Stockes（N-S）方程作为流动控制方程。直角坐标系下，方程组的守恒积分形式为：

(1)

式中：Q=（ρ_F，ρ_Fu，ρ_Fv，ρ_Fw，ρ_Fe）^T为守恒能量，ρ_F为流体密度，（u，v，w）为直角坐标系下沿三个坐标轴方向的速度分量，e为单位质量气体的总能量；

\partial Ω

为某一固定流体域Ω的边界；V为该控制体的体积，S为其边界；n是边界的外法向量；F为矢通量。

空气属性动力黏性系数μ通过Sutherland公式计算，根据空气属性表数据进行插值确定导热系数k的值。

本文使用国家数值风洞NNW-HyFLOW软件进行CFD计算，采用k-ω SST湍流模型。数值计算相关边界条件、网格设置等参考文献^[11]进行设置，以满足转捩位置、计算收敛性等方面的要求。

1.2 Kriging代理模型

Kriging代理模型的数学表达式可写成全局估计和当地差异的线性叠加：

f (x) = R (x) + C (x)

(2)

式中：x为输入参数向量； R（x）为代表全局近似估计的回归函数；C（x）为表征当地差异的相关函数。C（x）的协方差矩阵可表示为：

C o v [C (x_{i}), C (x_{j})] = σ^{2} R [R (x_{i}, x_{j})]

(3)

Kriging的预测输出可写为：

\hat{y} (x) = \hat{β} + r^{T} R^{- 1} (y - 1 \hat{β})

(4)

式中：

\hat{β}

为β的估计值；R为m×m维相关矩阵；r为m维相关向量；y为采样样本点响应值组成的m维列项量。未知参数β和σ²通过最小二乘意义下的估计得到。

1.3 热气动弹性耦合计算方法

热气动弹性问题是典型的流-固-热耦合问题，将三场耦合系统分解为两个耦合子系统模块（热模块和气动弹性模块）进行求解。

热模块是气动热和结构瞬态热传导耦合系统，将气动热计算得到的热流传输给热传导模块，再将热传导计算获得的结构瞬态温度分布传输给气动热计算，更新下一时间步的热流；气动弹性模块将气动力模型计算的非定常气动力传给结构动力学模块，再将由结构动力学模型计算的弹性变形传输给气动力计算，更新下一时间步的非定常气动力。本文采用气动热-气动弹性的单向耦合策略，将热模块得到的结构温度场传输给气动弹性模块。

气动弹性模块中使用有限元方法进行结构动力学求解，瞬态动力学分析中的基本运动方程为：

M \ddot{u} + C \dot{u} + K u = F (t)

(5)

式中：M为质量阵；C为阻尼阵；K为刚度阵；

\ddot{u}

为节点加速度向量；

\dot{u}

为节点速度向量；u为节点位移向量；F（t）为随时间t变化的载荷矩阵。

热模块中，由于辐射层的存在，上表面为对流边界，同时有气动热流Q_aero输入和辐射热流Q_rad输出，因此，热防护壁板上表面的净热流通量为Q=Q_aero-Q_rad。辐射热流与壁面温度T_wall和环境温度T_env之间的温差有关，计算式为：

Q_{r a d} = σ_{S - B} ε (T_{w a l l}^{4} - T_{e n v}^{4})

(6)

式中：σ_S-B为Stefan-Boltzmann常数，取σ_S-B=5.669×10^-8 W·m^-2·T^-4；辐射发生率ε=0.7。

结构内部考虑热传导，二维壁板的瞬态热传导偏微分方程为：

ρ_{S} c \frac{\partial T}{\partial t} = k (\frac{\partial^{2} T}{\partial x^{2}} + \frac{\partial^{2} T}{\partial z^{2}})

(7)

式中，x为壁板弦向坐标，z为壁板厚度方向坐标， ρ_S为材料密度，c为材料比热，T为温度，k为材料导热系数，t为时间。采用中心差分进行空间离散，向前差分进行时间推进，前后缘端面与下表面为绝热边界。基于Kriging代理模型的热防护壁板热气动弹性分析流程如图2所示。

图2热防护壁板热气动弹性分析流程

Fig.2Aerothermoelastic analysis flow of TPS panels

2 基于Kriging的壁板气动力/热代理模型

图3展示了建立代理模型的一般流程，包括：①确定系统的状态空间，将系统的输入输出变量参数化，并确定参数边界；②参数状态空间取样，使用一定的采样方法，从状态空间中选取样本集和验证集；③样本点响应，高阶方法对样本集工况进行计算，生成训练数据；④生成代理模型，基于样本状态点与响应结果的输入输出，建立代理模型；⑤精度评估，根据交叉验证思想在评估集对代理模型精度进行评估；⑥应用代理模型进行预测。

图3代理模型建立流程

Fig.3Flow for building the surrogate model

2.1 气动力/热样本数据计算方法验证

参考高速圆管绕流实验^[15]以及对给定形变和表面温度的壁板进行计算来验证CFD方法的准确性，以保证样本的精度。

图4（a）展示了高速流绕过圆管后的流场马赫数云图以及1/4截面圆管表面位置在极坐标轴下的角坐标表达。图4（b）、图4（c）为CFD计算圆管表面的无量纲压强、热流与实验数据的对比^[15]，CFD计算结果与实验结果吻合。CFD计算驻点热流为506.26 kW/m²，文献实验值为483 kW/m²，两者相差4.82%。

图4高速圆管绕流马赫数云图和表面压强及热流分布

Fig.4Mach number cloud chart of high-speed flow around a tube & surface pressure and heat flux distributions

使用壁板的模态表示壁板的结构位移，式（8）验证工况壁板位移；式（9）验证工况壁板表面温度分布。

w (x) = \overset{3}{\underset{i = 1}{Σ}} η_{i} s i n (i π \frac{x}{a})

(8)

式中，η₁=η₂=0，η₃=3.5。

\begin{matrix} T (x) = 701.6 x^{4} - 2087 x^{3} + 1475 x^{2} - \\ 193 . 9 x + 1246 \end{matrix}

(9)

表2展示了CFD方法的网格无关性验证参数，2号网格计算壁板表面最大压强与最大热流与密网格结果相对误差都在0.5%以下，为兼顾计算精度与效率，本文使用2号网格设置。图5展示了使用CFD方法计算表面的气动压强和表面热流分布，与文献^[11]结果的对比。图中同样展示了气动力、热工程算法，即活塞理论、参考焓法的计算结果。

表2网格无关性验证参数

Tab.2Grid independence study parameter

图5CFD方法进行气动力/热验证

Fig.5Validation of aerodynamic forces & aerodynamic heating by CFD method

CFD计算结果与文献^[11]结果吻合，气动压强平均相对误差为3.4%，而气动热流平均相对误差为15.6%，考虑到采用的湍流模型与文献不同，且不同数值计算模型的误差大于计算本身的误差^[11]，该偏差是可接受的。

综上，认为CFD方法气动力/热计算精度满足要求。此外，CFD方法计算的气动压强、气动热流与工程算法的计算结果在壁板上的分布不同，凸显了建立保留高保真特性的代理模型带来的预测精度提升的意义。

2.2 气动力/热代理模型的建立与评估方法

2.2.1 状态空间参数化

气动力/热代理模型的输入为来流马赫数Ma₁、壁板形变以及壁板表面温度。马赫数作为标量可直接作为模型输入参数，壁板形变和表面温度通过空间坐标进行参数化处理。

壁板形变w（x）用其前六阶模态函数叠加表示，如式（10）所示，模态函数Φ（x）=sin（ixπ/a），h为壁板厚度；表面温度分布T（x）用壁板弦向的坐标x的n阶多项式表达，如式（11）所示，已有研究表明，三次多项式可满足精度要求^[11]。

w (x) = η_{1} Φ_{1} (x) + η_{2} Φ_{2} (x) + \dots + η_{6} Φ_{6} (x)

(10)

T (x) = T_{1} + T_{2} x + T_{3} x^{2} + T_{4} x^{3}

(11)

参考已有的研究工作^[11]，参数空间如表3所示。此外，为优化样本点对应壁板变形和壁板表面温度在参数空间内的分布，设置位移振幅和表面温度极值约束：壁板的无量纲形变幅值≤10，200 K≤壁板表面温度极值≤2 000 K。本文使用有约束的拉丁超立方取样方法进行样本点和评估点的选取。

表3代理模型参数空间

Tab.3Parameter space of the surrogate model

2.2.2 代理模型的生成输入与预测输出

参考式（2），基于Kriging的气动力/热代理模型可以表示为局部偏差C（d，X）和全局估计R（d，X）的和：

y (d) = R (d, X) + C (d, X)

(12)

(13)

其中： y（d）为代理模型对参数空间内某状态点的预测输出；d为某状态点所对应的输入参数组成的向量；X为训练样本点的输入参数矩阵；n为建立模型所用样本数量。

为减少数据量级对模型的干扰，使用无量纲参数气动压强系数C_P和斯坦顿数St作为模型的输出。通过式（14）和式（15）得到气动压强和气动热流。

P_{steady} (x) = \frac{1}{2} y_{C_{P}} (d) ρ_{3} U_{3}^{2}

(14)

式中，

y_{C_{P}} （ d ）

是壁板表面压强系数的预测值，ρ₃和U₃为3区密度和速度。

Q_{aero} = c_{p} U_{3} ρ_{3} (T_{s} (x) - T_{0,3}) y_{S t} (d)

(15)

式中，y_St（d）是壁板表面斯坦顿数的预测值，c_p为恒压比热，T_s（x）为壁板表面温度分布，T_0，3为3区总温。

2.2.3 代理模型的精度评估方法

本文基于交叉验证的思想，选取与构建代理模型所用样本点不重合的300个评估点，通过将代理模型预测结果与CFD方法计算结果进行对比，对模型精度进行量化评价。

归一化均方根误差（normalized root-mean-square error，NRMSE）和最大值误差为：

e_{N R M S E}^{j} = \frac{\sqrt{\frac{1}{p} \overset{p}{\underset{i = 1}{Σ}} {(S_{i}^{j} - C_{i}^{j})}^{2}}}{C_{m a x}^{j} - C_{m i n}^{j}}

(16)

L_{\infty}^{j} = \frac{{|S^{j} - C^{j}|}_{m a x}}{C_{m a x}^{j} - C_{m i n}^{j}}

(17)

其中，S为代表模型预测结果，C为CFD方法计算结果，i为壁板弦向节点标号，j为评估点标号，p为板上节点总数。在以上指标的基础上，计算四个值作为量化评价整个代理模型的评价指标：平均均方根误差

{\bar{e}}_{NRMSE}

、平均最大值误差

{\bar{L}}_{\infty}

、最大最大值误差L_∞max以及L_∞小于10%的评估点数量占比

p_{L_{\infty}}

≤10%。同时，检验上述四个评价指标随样本数量变化情况，以验证模型的收敛性。

此外，由式（18）计算预测值与CFD方法计算结果的平均相对误差在壁板弦向位置的分布：

e_{a r e (i)} = \frac{1}{n} \overset{n}{\underset{j = 1}{Σ}} \frac{|S_{i}^{j} - C_{i}^{j}|}{C_{i}^{j}}

(18)

2.3 基于Kriging的气动力代理模型

2.3.1 气动力代理模型精度评估

表4给出了使用不同样本数量所建立的气动力代理模型的精度评价指标和前人研究评价指标^[11]作为参考。由表4可知，使用1 200个和1 500个样本点得到的气动力代理模型，评价指标全部优于文献^[11]。

表4不同样本数量的气动力代理模型评价指标

Tab.4Evaluation metrics of aerodynamic surrogate models with different sample sizes

图6展示了气动力代理模型评价指标随样本数量的变化。样本数量从200增加到1 500，代理模型的精度评价指标不再提升，验证了生成代理模型所用样本数量的收敛性。

图6气动力代理模型评价指标随样本数量的变化

Fig.6Variation of aerodynamic surrogate models evaluation metrics with sample size

图7展示了1 500个样本点生成的气动力代理模型的气动压强系数计算结果与CFD方法计算结果的相对误差分布。差异较大的地方发生在壁板前后缘处，可见Kriging代理模型作为广义插值算法存在插值模型在边缘区域预测精度下降的通病。有待后续进一步探索以提高代理模型在边缘的预测精度。

2.3.2 非定常气动力代理模型

对于非定常气动力的预测，在代理模型预测稳态气动力的基础上，结合三阶活塞理论非定常效应项，建立热防护壁板的非定常气动力代理模型。代理模型非定常压强系数计算公式为：

图7气动力代理模型平均相对误差位置分布

Fig.7Position distribution of average relative error for aerodynamic surrogate models

C_{P} (x, t) = C_{P, vel} (x) + {\bar{C}}_{P} (x, t) + y_{C_{P}} (d)

(19)

式中，

y_{C_{P}} （ d ）

为代理模型输出的壁板表面压强系数的预测值，C_P_，vel与

{\bar{C}}_{P}

为基于三阶活塞理论气动阻尼项。

最终，壁板表面非定常气动压强的计算式为：

q_{a} (x, t) = \frac{1}{2} C_{P} (x, t) ρ_{3} U_{3}^{2}

(20)

2.3.3 非定常气动力模型验证

使用气动力代理模型对指定运动形式的壁板的非定常广义气动力进行预测，验证基于代理模型的非定常气动力模型的准确性。

广义气动力计算式为：

F_{G} (t) = \int (Φ (x) \cdot q_{a} (x, t)) d x

(21)

式中，Φ（x）为模态函数。壁板振动方程为：

y (x, t) = 0.005 s i n (2 π f t) \overset{n}{\underset{i = 1}{Σ}} η_{i} s i n (i π \frac{x}{1.5})

(22)

图8展示了运动壁板非定常广义气动力的结果。选取一阶模态作为壁板的运动形式，

{\bar{η}}_{1}

=12，频率f=140 Hz。计算工况：来流马赫数为10，壁板表面温度均为1 200 K。结果表明，代理模型预测的非定常广义气动力结果与CFD非定常计算结果一致，验证了基于代理模型计算非定常气动力的准确性。图8同样展示了文献^[11]和本文使用非定常三阶活塞理论的计算结果，与CFD非定常计算的结果存在明显差异，再一次展示了建立具有高精度预测的代理模型的重要意义。

注意到壁板振动最大无量纲幅值为12，超过了代理模型的状态空间，而代理模型仍可以进行准确非定常广义气动力预测。这说明Kriging模型作为典型插值模型，具有一定的外插泛化能力。图9显示了三个不同壁板振动幅值对应的代理模型预测结果和CFD计算非定常广义气动力的结果，

{\bar{η}}_{1}

=15时误差最大为13.4%，而

{\bar{η}}_{1}

=20时代理模型预测误差最大为17.2%，超过参数边界越多，泛化能力越差，因此实际应用时可在参数空间基础上适当拓展应用范围，但是不建议太多。

图8运动壁板非定常F_G

Fig.8Unsteady F_G of moving panel

图9不同振动幅值代理模型和CFD非定常F_G结果

Fig.9Surrogate model and CFD unsteady F_G results under different vibration amplitudes

计算效率方面，使用三阶活塞理论以及气动力代理模型，计算验证工况所用时间均在0.5 s以内。而使用CFD方法进行非定常气动力计算的耗时在8 h以上，是代理模型的6万倍（CPU：Gold 6152 2.10 GHz）。因此代理模型方法相较于CFD方法，将计算效率提高了4个数量级。

2.4 基于Kriging的气动热代理模型

2.4.1 基于遗传算法的气动热代理模型参数优化

本文在建立气动热代理模型时，使用遗传算法（genetic algorithm，GA），以模型精度评价指标为优化目标，对Kriging模型中的相关参数θ进行寻优。

表5展示了使用遗传算法前后，由1 200个样本数据生成的气动热代理模型的评价指标。由表5可知，使用遗传算法后，

{\bar{e}}_{NRMSE}

、

{\bar{L}}_{\infty}

、L_∞max减小，

p_{L_{\infty}}

≤10%增大，模型的精度得到了提升。

表5遗传算法优化前后的模型评价指标

Tab.5Model evaluation metrics before and after optimization by genetic algorithm

2.4.2 气动热代理模型精度评估

表6给出了使用不同样本数量所建立的气动热代理模型的精度评价指标，并给出文献^[11]中的评价指标作为参考。由表6可知，使用2 000个和2 400个样本点得到的气动力代理模型的评价指标全部优于文献^[11]。

表6不同样本数量的气动热代理模型评价指标

Tab.6Evaluation metrics of aerothermal surrogate models with different sample sizes

图10展示了气动热代理模型评价指标随样本数量的变化。样本数量从200增加到2 400，代理模型的精度评价指标不再提升，验证了生成代理模型所用样本数量的收敛性。图11展示了2 200个样本点生成的气动热代理模型的斯坦顿数计算结果与CFD方法计算结果的相对误差分布。差异较大的地方发生在壁板后缘处，同样可解释为是由插值模型对边缘区域预测能力不足引起，同样有待进一步探索。

图10气动热代理模型评价指标随样本数量的变化

Fig.10Variation of aerothermal surrogate models evaluation metrics with sample size

图11气动热代理模型平均相对误差位置分布

Fig.11Position distribution of average relative error for aerothermal surrogate models

由表5和表6可知，生成相似精度的气动热代理模型相比气动力代理模型需要更多的样本数量。此外，生成气动热代理模型的过程还用到了遗传算法进行参数优化。结果表现出生成气动热代理模型的复杂性。主要原因在于气动热的计算不仅与马赫数和壁板形变有关，还受到壁板表面温度分布的影响，影响参数更多，生成代理模型的各项成本也就越高。

3 基于代理模型的热防护壁板热气动弹性分析

考虑典型气动响应条件下的壁板形变快照和壁板刚体外形（即平板）两种外形作为气动热模块的输入。使用气动热代理模型，在气动热模块进行热防护壁板结构瞬态热传导计算，将某时刻的温度场快照作为热载荷输入气动弹性求解模块，使用气动力代理模型，进行热气动弹性计算。典型形变快照^[16]为：

{\bar{η}}_{1}

=6，

{\bar{η}}_{2}

=-3，

{\bar{η}}_{3}

=1。

本文进行了如表7所示的四组工况的热气动弹性计算，其中来流条件均为1区参数Ma₁。由于气动力和气动热非定常特征时间差异较大，模拟高速飞行器进行诸如加速突防机动动作的情景时，飞行速度大幅增大而飞行器结构温度变化不大，工况3相比工况2模拟了类似的情形，热传导计算来流马赫数相同，而热气动弹性计算来流动压更大。四组热气动弹性工况得到以下典型非线性响应。

表7热防护壁板热气动弹性计算工况

Tab.7Aerothermoelastic calculation conditions of TPS panel

3.1 热屈曲

工况2中，在传热60 s、100 s和180 s时刻温度场下，壁板热气动弹性响应形式均为热屈曲。图12（a）展示了60 s时刻温度场下，壁板响应时程图。图12（b）展示了0 s、60 s、100 s和180 s时刻温度场下，壁板热屈曲响应收敛达到稳定后，壁板位移形变挠曲线，可见在加入温度场之后，壁板产生了失稳。图13展示了该四个温度场下壁板的静力学分析得到的热形变和金属板中性层无量纲热应力分布。0 s时刻结构没有升温，热应力为0。随着热防护壁板热传导的进行，温度场板内温度越高，热膨胀所引起的热应力越大，也导致热屈曲位移极值越大。

图12屈曲时程图和稳定时刻壁板位移

Fig.12Buckling time history & panel displacement at response stabilization

3.2 极限环运动

工况3中，在传热300 s时刻温度场下，壁板热气动弹性响应形式为极限环运动。图14展示了壁板响应时程图和相平面图。极限环运动的时程图中，壁板振动形式呈现一定的周期性规律，图14（a）中将部分响应曲线进行了放大以直观展示。相应相平面图中只有单条封闭的相平面曲线。

图13热形变和中性层无量纲应力

Fig.13Thermal deformation and dimensionless stress of the neutral layer

图14极限环运动响应图

Fig.14Limit cycle oscillation response diagram

3.3 倍周期运动

工况3中，在传热100 s时刻温度场下，壁板热气动弹性响应形式为倍周期运动。图15展示了壁板响应时程图、频谱图和相平面图。倍周期运动的时程图中，壁板振动形式呈现一定的周期性规律，图15（a）中同样进行了放大，但比极限环运动的振动形式更加复杂。频谱图中各波峰对应的频率为10.22 Hz、20 Hz、30.22 Hz和40 Hz，为1、2、3、4倍整数关系。相平面图呈多条封闭曲线绕成的带状。

图15倍周期运动响应图

Fig.15Period-doubling motion response diagram

3.4 准周期运动

工况4中，在传热300 s时刻温度场下，壁板热气动弹性响应形式为准周期运动。图16展示了壁板响应相关结果图。准周期运动与倍周期运动类似，振动呈现一定的周期性（时程图中进行了放大），但频谱图中各波峰对应的频率之间不存在整倍数关系，相平面图呈多条封闭曲线绕成的带状。

图16准周期运动响应图

Fig.16Quasi-periodic motion response diagram

3.5 混沌运动

工况1中，在传热100 s时刻温度场下，壁板热气动弹性响应形式为混沌运动。图17展示了壁板响应相关结果图。混沌运动响应中，时程图表现的壁板振动没有明显规律，频谱图中多个峰值在一定频率范围内连续分布，相平面曲线分布杂乱。

3.6 结果分析

表8展示了本文计算工况下的热气动弹性响应形式。除表8中的算例外，所有工况加载热传导0 s时刻的温度场的响应结果均为静稳定收敛。这说明壁板内温度升高而引起的热应力和材料温变是热防护壁板热气动弹性失稳的主要因素。此外，随着所加载的温度场时刻越来越靠后，壁板的热气动弹性响应形式一般呈现由简单到复杂再到简单的发展规律。在响应幅值方面，随着壁板内温度越来越高，热应力增大，响应幅值增大。

图17混沌运动响应图

Fig.17Chaotic motion response diagram

如表8中工况4相对工况3，热混沌发展为极限环运动出现所对应的温度场时刻更靠前。这说明来流马赫数越高，热防护壁板温度升高越快，周期运动、热混沌、极限环运动等失稳时刻提前。如工况1和其余工况的响应结果所示，得到的壁板热气弹响应形式随温度场时刻增大的变化趋势相同。因此对于单向耦合策略，壁板响应形式随壁板温度的变化规律并不受热防护壁板热传导计算时的壁板形变形式的影响。

表8各工况热气动弹性响应形式

Tab.8Response forms of aerothermoelastic under various cases

4 结论

本文针对高速来流中的热防护壁板，基于Kriging代理模型方法，对建立高速流中的热防护壁板非定常气动力/热代理模型开展了研究。使用基于有限元法的结构动力学建立气动弹性模块，编写结构瞬态热传导程序建立气动热模块，搭建热防护壁板热气动弹性分析框架。最终，使用基于Kriging的气动力/热代理模型，开展热防护壁板热气动弹性计算与分析。得到以下结论：

1）对于气动力/热代理模型：使用Kriging代理模型方法建立的气动力/热代理模型，除了能够准确、高效预测气动力/热，还具有一定的泛化能力；对于气动力代理模型，建立气动热代理模型需要更多的样本和方法成本；代理模型方法相比传统CFD方法，计算效率提高了4个数量级。

2）对于热气动弹性分析：气动热效应引起的面内热应力和材料温变是热防护壁板在高速流中失稳的主要因素；随着壁板内温度的升高，热防护壁板的热气动弹性响应幅值逐渐增大，响应形式一般呈现由简单到复杂再到简单的发展规律。

综合来看，本文所建立的基于Kriging的气动力/热代理模型，能够弥补工程算法（如活塞理论、参考焓法等）在气动力/热预测精度方面的不足，且保留了工程算法效率高的特点，极大提升了效率。该方法为非定常气动力/热的快速准确计算提供了有效手段，能够为高速飞行器的气动力/热载荷准确预估、热防护结构设计、飞行安全评估等提供重要理论依据。

图1热防护壁板

Fig.1TPS panel

下载: 全尺寸图片

图2热防护壁板热气动弹性分析流程

Fig.2Aerothermoelastic analysis flow of TPS panels

下载: 全尺寸图片

图3代理模型建立流程

Fig.3Flow for building the surrogate model

下载: 全尺寸图片

图4高速圆管绕流马赫数云图和表面压强及热流分布

Fig.4Mach number cloud chart of high-speed flow around a tube & surface pressure and heat flux distributions

下载: 全尺寸图片

图5CFD方法进行气动力/热验证

Fig.5Validation of aerodynamic forces & aerodynamic heating by CFD method

下载: 全尺寸图片

图6气动力代理模型评价指标随样本数量的变化

Fig.6Variation of aerodynamic surrogate models evaluation metrics with sample size

下载: 全尺寸图片

图7气动力代理模型平均相对误差位置分布

Fig.7Position distribution of average relative error for aerodynamic surrogate models

下载: 全尺寸图片

图8运动壁板非定常F_G

Fig.8Unsteady F_G of moving panel

下载: 全尺寸图片

图9不同振动幅值代理模型和CFD非定常F_G结果

Fig.9Surrogate model and CFD unsteady F_G results under different vibration amplitudes

下载: 全尺寸图片

图10气动热代理模型评价指标随样本数量的变化

Fig.10Variation of aerothermal surrogate models evaluation metrics with sample size

下载: 全尺寸图片

图11气动热代理模型平均相对误差位置分布

Fig.11Position distribution of average relative error for aerothermal surrogate models

下载: 全尺寸图片

图12屈曲时程图和稳定时刻壁板位移

Fig.12Buckling time history & panel displacement at response stabilization

下载: 全尺寸图片

图13热形变和中性层无量纲应力

Fig.13Thermal deformation and dimensionless stress of the neutral layer

下载: 全尺寸图片

图14极限环运动响应图

Fig.14Limit cycle oscillation response diagram

下载: 全尺寸图片

图15倍周期运动响应图

Fig.15Period-doubling motion response diagram

下载: 全尺寸图片

图16准周期运动响应图

Fig.16Quasi-periodic motion response diagram

下载: 全尺寸图片

图17混沌运动响应图

Fig.17Chaotic motion response diagram

下载: 全尺寸图片

表1热防护壁板相关参数

Tab.1Relevant parameters of the TPS panel

下载: 全尺寸图片

表2网格无关性验证参数

Tab.2Grid independence study parameter

下载: 全尺寸图片

表3代理模型参数空间

Tab.3Parameter space of the surrogate model

下载: 全尺寸图片

表4不同样本数量的气动力代理模型评价指标

Tab.4Evaluation metrics of aerodynamic surrogate models with different sample sizes

下载: 全尺寸图片

表5遗传算法优化前后的模型评价指标

Tab.5Model evaluation metrics before and after optimization by genetic algorithm

下载: 全尺寸图片

表6不同样本数量的气动热代理模型评价指标

Tab.6Evaluation metrics of aerothermal surrogate models with different sample sizes

下载: 全尺寸图片

表7热防护壁板热气动弹性计算工况

Tab.7Aerothermoelastic calculation conditions of TPS panel

下载: 全尺寸图片

表8各工况热气动弹性响应形式

Tab.8Response forms of aerothermoelastic under various cases

下载: 全尺寸图片

图1热防护壁板

Fig.1TPS panel

图2热防护壁板热气动弹性分析流程

Fig.2Aerothermoelastic analysis flow of TPS panels

图3代理模型建立流程

Fig.3Flow for building the surrogate model

图4高速圆管绕流马赫数云图和表面压强及热流分布

Fig.4Mach number cloud chart of high-speed flow around a tube & surface pressure and heat flux distributions

图5CFD方法进行气动力/热验证

Fig.5Validation of aerodynamic forces & aerodynamic heating by CFD method

图6气动力代理模型评价指标随样本数量的变化

Fig.6Variation of aerodynamic surrogate models evaluation metrics with sample size

图7气动力代理模型平均相对误差位置分布

Fig.7Position distribution of average relative error for aerodynamic surrogate models

图8运动壁板非定常F_G

Fig.8Unsteady F_G of moving panel

图9不同振动幅值代理模型和CFD非定常F_G结果

Fig.9Surrogate model and CFD unsteady F_G results under different vibration amplitudes

图10气动热代理模型评价指标随样本数量的变化

Fig.10Variation of aerothermal surrogate models evaluation metrics with sample size

图11气动热代理模型平均相对误差位置分布

Fig.11Position distribution of average relative error for aerothermal surrogate models

图12屈曲时程图和稳定时刻壁板位移

Fig.12Buckling time history & panel displacement at response stabilization

图13热形变和中性层无量纲应力

Fig.13Thermal deformation and dimensionless stress of the neutral layer

图14极限环运动响应图

Fig.14Limit cycle oscillation response diagram

图15倍周期运动响应图

Fig.15Period-doubling motion response diagram

图16准周期运动响应图

Fig.16Quasi-periodic motion response diagram

图17混沌运动响应图

Fig.17Chaotic motion response diagram

表1热防护壁板相关参数

Tab.1Relevant parameters of the TPS panel

表2网格无关性验证参数

Tab.2Grid independence study parameter

表3代理模型参数空间

Tab.3Parameter space of the surrogate model

表4不同样本数量的气动力代理模型评价指标

Tab.4Evaluation metrics of aerodynamic surrogate models with different sample sizes

表5遗传算法优化前后的模型评价指标

Tab.5Model evaluation metrics before and after optimization by genetic algorithm

表6不同样本数量的气动热代理模型评价指标

Tab.6Evaluation metrics of aerothermal surrogate models with different sample sizes

表7热防护壁板热气动弹性计算工况

Tab.7Aerothermoelastic calculation conditions of TPS panel

表8各工况热气动弹性响应形式

Tab.8Response forms of aerothermoelastic under various cases

引用提醒

图(17) / 表(8)

引用本文

易子钧, 冀春秀, 谢丹, 等. Kriging气动力/热代理模型与热防护壁板热气动弹性分析[J]. 国防科技大学学报, 2026, 48(1): 139-149.

复制

YI Z J, JI C X, XIE D, et al. Kriging-based aerodynamic/aerothermal surrogate models and aerothermoelastic analysis of TPS panel[J]. Journal of National University of Defense Technology, 2026, 48(1): 139-149.

Copy

计量

图1热防护壁板

Fig.1TPS panel

图2热防护壁板热气动弹性分析流程

Fig.2Aerothermoelastic analysis flow of TPS panels

图3代理模型建立流程

Fig.3Flow for building the surrogate model

图4高速圆管绕流马赫数云图和表面压强及热流分布

Fig.4Mach number cloud chart of high-speed flow around a tube & surface pressure and heat flux distributions

图5CFD方法进行气动力/热验证

Fig.5Validation of aerodynamic forces & aerodynamic heating by CFD method

图6气动力代理模型评价指标随样本数量的变化

Fig.6Variation of aerodynamic surrogate models evaluation metrics with sample size

图7气动力代理模型平均相对误差位置分布

Fig.7Position distribution of average relative error for aerodynamic surrogate models

图8运动壁板非定常F_G

Fig.8Unsteady F_G of moving panel

图9不同振动幅值代理模型和CFD非定常F_G结果

Fig.9Surrogate model and CFD unsteady F_G results under different vibration amplitudes

图10气动热代理模型评价指标随样本数量的变化

Fig.10Variation of aerothermal surrogate models evaluation metrics with sample size

图11气动热代理模型平均相对误差位置分布

Fig.11Position distribution of average relative error for aerothermal surrogate models

图12屈曲时程图和稳定时刻壁板位移

Fig.12Buckling time history & panel displacement at response stabilization

图13热形变和中性层无量纲应力

Fig.13Thermal deformation and dimensionless stress of the neutral layer

图14极限环运动响应图

Fig.14Limit cycle oscillation response diagram

图15倍周期运动响应图

Fig.15Period-doubling motion response diagram

图16准周期运动响应图

Fig.16Quasi-periodic motion response diagram

图17混沌运动响应图

Fig.17Chaotic motion response diagram

表1热防护壁板相关参数

Tab.1Relevant parameters of the TPS panel

表2网格无关性验证参数

Tab.2Grid independence study parameter

表3代理模型参数空间

Tab.3Parameter space of the surrogate model

表4不同样本数量的气动力代理模型评价指标

Tab.4Evaluation metrics of aerodynamic surrogate models with different sample sizes

表5遗传算法优化前后的模型评价指标

Tab.5Model evaluation metrics before and after optimization by genetic algorithm

表6不同样本数量的气动热代理模型评价指标

Tab.6Evaluation metrics of aerothermal surrogate models with different sample sizes

表7热防护壁板热气动弹性计算工况

Tab.7Aerothermoelastic calculation conditions of TPS panel

表8各工况热气动弹性响应形式

Tab.8Response forms of aerothermoelastic under various cases

陈浩宇, 王彬文, 宋巧治, 等. 高超声速飞行器热颤振研究现状与展望[J]. 航空工程进展,2022,13(1):19-27.CHEN H Y, WANG B W, SONG Q Z,et al. Research progress and prospect of thermal flutter of hypersonic vehicles[J]. Advances in Aeronautical Science and Engineering,2022,13(1):19-27.(in Chinese)

樊会涛, 段鹏飞, 袁成. 航空颠覆性技术初探[J]. 航空学报,2024,45(5):529893.FAN H T, DUAN P F, YUAN C. Disruptive technologies in aviation:preliminary study[J]. Acta Aeronautica et Astronautica Sinica,2024,45(5):529893.(in Chinese)

王泽, 宋述芳, 王旭, 等. 数据驱动的气动热建模预测方法总结与展望[J]. 气体物理,2024,9(4):39-55.WANG Z, SONG S F, WANG X,et al. Summary and prospect of data-driven aerothermal modeling prediction methods[J]. Physics of Gases,2024,9(4):39-55.(in Chinese)

JI C X, YI Z, XIE D,et al. Exploration on flutter mechanism of a damaged transonic rotor blade using high-fidelity fluid-solid coupling method[J]. Journal of Fluid Mechanics,2024,998: A15.

李凯, 杨静媛, 高传强, 等. 基于POD和代理模型的静气动弹性分析方法[J]. 力学学报,2023,55(2):299-308.LI K, YANG J Y, GAO C Q,et al. Static aeroelastic analysis based on proper orthogonal decomposition and surrogate model[J]. Chinese Journal of Theoretical and Applied Mechanics,2023,55(2):299-308.(in Chinese)

王泽, 王梓伊, 王旭, 等. 一种数据驱动的气动热预示模型[J]. 空气动力学学报,2023,41(5):12-19.WANG Z, WANG Z Y, WANG X,et al. A data-driven aeroheating prediction model[J]. Acta Aerodynamica Sinica,2023,41(5):12-19.(in Chinese)

张伟伟, 王旭, 寇家庆. 面向流体力学的多范式融合研究展望[J]. 力学进展,2023,53(2):433-467.ZHANG W W, WANG X, KOU J Q. Prospects of multi-paradigm fusion methods for fluid mechanics research[J]. Advances in Mechanics,2023,53(2):433-467.(in Chinese)

LIU H J, HU H Y, ZHAO Y H,et al. Efficient reduced-order modeling of unsteady aerodynamics robust to flight parameter variations[J]. Journal of Fluids and Structures,2014,49:728-741.

晏筱璇, 韩景龙, 马瑞群. 高超声速气动热弹性分析降阶研究[J]. 振动工程学报,2022,35(2):475-486.YAN X X, HAN J L, MA R Q. Reduced-order modeling research for hypersonic aerothermoelastic analysis[J]. Journal of Vibration Engineering,2022,35(2):475-486.(in Chinese)

王洋, 袁军娅, 王洪兴. 基于代理模型和线性近似的快速气动热边界求解方法[J]. 导弹与航天运载技术,2018(4):11-17.WANG Y, YUAN J Y, WANG H X. Fast method to determine thermal boundary based on surrogate model and linear approximation[J]. Missiles and Space Vehicles,2018(4):11-17.(in Chinese)

CROWELL A R, MCNAMARA J J, MILLER B A. Hypersonic aerothermoelastic response prediction of skin panels using computational fluid dynamic surrogates[J]. Journal of Aeroelasticity and Structural Dynamics,2011,2(2):3-30.

谢丹, 冀春秀, 景兴建. 高超声速典型弹道下的壁板热气动弹性动力学分析[J]. 航空学报,2021,42(11):524843.XIE D, JI C X, JING X J. Dynamics analysis of panel aerothermoelasticity in typical hypersonic trajectories[J]. Acta Aeronautica et Astronautica Sinica,2021,42(11):524843.(in Chinese)

XIE D, DONG B, JING X J. Effect of thermal protection system size on aerothermoelastic stability of the hypersonic panel[J]. Aerospace Science and Technology,2020,106:106170.

JI C X, XIE D, ZHANG S H,et al. Spectral proper orthogonal decomposition reduced-order model for analysis of aerothermoelasticity[J]. AIAA Journal,2023,61(2):793-807.

WIETING A R. Experimental study of shock wave interference heating on a cylindrical leading edge: NASA-TM-100484[R]. Hampton, Virginia: NASA Langley Research Center,1987.

CULLER A J, MCNAMARA J J. Studies on fluid-thermal-structural coupling for aerothermoelasticity in hypersonic flow[J]. AIAA Journal,2010,48(8):1721-1738.

Home

About Journal

Guide for Authors

Editorial Board

Publication Statement

Open Access Statement

Contact

Journal Subscription

Rss

AI assistant

Chinese

1 研究对象与研究方法

2 基于Kriging的壁板气动力/热代理模型

3 基于代理模型的热防护壁板热气动弹性分析

4 结论