表面热流辨识技术在边界层转捩位置测量中的应用初步研究

钱炜祺，周 宇，何开锋，石友安

（中国空气动力研究与发展中心空气动力学国家重点实验室，四川绵阳 621000）

0 引 言

飞行器飞行过程中，边界层转捩位置对飞行器的摩阻、压阻及飞行性能有很大的影响，转捩位置的确定是飞行器设计的关键技术之一。除了航空航天领域，边界层转捩位置的确定在船舶、机械制造、降噪、化工等领域也具有实用价值。然而，由于转捩过程相当复杂，包含着层流失稳、T－S波的出现和放大、三维扰动的发展和流向涡的形成、流动二次失稳、湍流斑生成等一系列复杂的流动现象，使得目前对转捩机理的认识还不完善，采用流体力学理论方法来计算转捩位置较为困难。因此，迄今为止，实验仍是做出转捩判断的最可靠的手段。

1 实验中转捩测量的主要方法

目前，实验中转捩的测量方法主要有升华法、脉动压力测试技术、红外成像热图技术、磷光粉热图显示技术、液晶涂层法等［1－2］。其中，升华法是利用某些有升华现象的固态物质在边界层层流区和湍流区的升华速度不同的原理来显示转捩位置，但由于需要较长的吹风时间，限制了该方法在高速风洞中的使用。脉动压力测试技术是利用脉动压力传感器测量压力脉动，根据脉动压力的变化来确定转捩位置，但风洞背景噪声的干扰对测量结果有较大影响。由于层流和湍流区域的壁面温度不一样，理论上所有测壁面温度的方法都可以用来测量转捩，红外成像热图技术、磷光粉热图显示技术和液晶涂层法正是基于这一思想测量表面温度的技术，目前在风洞中使用较广，而在飞行试验中由于需要获取飞行器表面的图像信息，限制了这些技术在飞行试验中的使用。

在飞行试验中测量转捩的主要方法有：近壁热电偶测温法、光电二极管指示法、边界层声监测法等［3］。由于转捩前后边界层内脉动幅值增加，后两种方法分别对边界层内的光学信息和压力脉动进行监测来判断转捩的发生，近壁热电偶测温法则是利用边界层转捩前后壁面温度的变化来判断转捩。由于温度测量技术相对较为成熟，因而目前在工程实际飞行试验中主要采用的是近壁热电偶测温法。美国20世纪60年代进行的Reentry－F飞行试验项目［4－5］和Sandia实验室近年开展的SHARP－B2飞行试验项目［3］，都采用的是近壁热电偶测温法来判断转捩。尤其是Reentry－F飞行试验项目，该飞行试验的主要目的就是进行飞行器表面的转捩试验。试验飞行器是长3.962m、半锥角5°的尖锥体，最大飞行速度为20马赫。尖锥体的弹体由铍合金构成，在弹身上布置了21个温度测点（如图1所示），在每个测点沿厚度方向布置了4个热电偶传感器，其中距离表面最近的传感器距离表面0.254mm。飞行试验结束后，希望利用测点温度测量数据来得出表面热流沿飞行器表面的分布，进而确定转捩位置，文献［5］中采用了两种方法来进行温度测量数据处理，一种是用类似顺序函数法的反演方法由测点温度确定表面热流［6］；另一种是认为测点温度近似等于表面温度，由温度历程对时间求导，导数最大的点对应的就是转捩位置。结果表明，利用两种方法确定出的转捩位置基本一致（如图2所示，图中“↑”箭头表示采用温度历程求导方法得出的转捩位置；“↓”箭头表示采用分析表面热流辨识结果得出的转捩位置），同时，这一结果和地面试验结果也吻合较好。

图1 Reentry－F飞行器及温度测点布置Fig.1 Reentry－F spacecraft and its measurement positions

图2 对飞行试验数据用两种方法确定出转捩位置的比较Fig.2 Comparison of transition location

目前，虽然近壁热电偶测温法已在工程上实际应用，但该方法还存在两个不足：

（1）该方法要求热电偶安装位置靠近表面，但靠近表面位置测点的温升较高，对热电偶传感器的量程和精度要求较高。

（2）只有通过在防热层内钻孔的方式才能将温度传感器安装在靠近表面的位置，这势必对防热层的结构强度产生较大影响。

因此，有必要探索新的处理方法和思路对近壁热电偶测温法进行改进与完善，可对热电偶测点位置进行适当选取，使其既远离表面，又能够从测点温度有效地辨识反演出表面热流。

2 基于表面热流辨识技术确定转捩位置的处理方法

气动热参数辨识技术是热传导逆问题在气动热环境分析中的实际应用，包含防热材料热传导系数辨识和表面热流辨识两方面研究内容［7］。转捩位置的确定属于表面热流辨识问题，其基本原理如图3所示，在转捩测试平板内壁或材料内部沿流向布置多个温度测点Pi，通过这些测点的温度历程来辨识表面热流随流向位置和时间的变化历程。该方法与近壁热电偶测温法的主要差异在于测点离表面相对较远。但是，温度传感器远离表面后会带来表面热流辨识的不适定性问题。体现在两方面：一是在厚度方向，由于热传导的延迟性和耗散性，测点温度信息可能无法反映表面热流在时间方向上的某些变化；二是在流动方向，同样由于热传导在x方向的耗散，表面热流在x方向上的变化可能在测点信息中体现不出来。对于这两个问题，需要通过仿真计算辨识的方法来确定传感器安装位置，即首先根据理论计算或风洞实验事先确定转捩位置附近表面热流的大致形式，然后利用热传导正问题计算出测点上的温度历程，并在该温度历程上叠加随机误差来模拟测点温度的测量值，用来辨识出表面热流，最后将热流辨识结果与给定值比较，如果二者符合，则当前的测点位置是合理的，否则需调整测点位置。

图3 基于表面热流辨识技术的流动转捩测试示意图Fig.3 Sketch of transition determination method based on surface heat flux estimation

上述仿真辨识技术的核心是表面热流的辨识算法，对于一维传热问题的表面热流辨识算法，文献［8－9］中有详细的介绍，下面简要给出二维传热的表面热流辨识算法。仍以图3为例，设长度方向（流向）的范围为x∈［0，H］，厚度方向的范围为y∈［0，L］，时间方向为t∈［0，tf］。二维传热的控制方程为：

初始条件：t＝0：T＝T0。

式中ρ、Cp、k分别为材料的密度、比热和热传导系数，Q为随时间和空间变化的表面热流函数，就辨识问题而言，该函数是未知函数，通过如下的观测信息来确定。测点Pi（i＝1，M）的坐标为（xi，0），对应的温度测量历程：

式中，T表示测点温度的真值，v（t）表示测量误差噪声。因此，辨识问题等价于求合适的Q（x，t）使如下目标函数达极小的优化问题：

其中，T（xi，t，Q）表示表面热流为Q时的测点温度历程，通过有限元计算方法来求取；˜T为相应的实测值。由于Q（x，t）是空间和时间的函数，因此式（3）实际上是一个泛函表达式，考虑到有限元计算方法对计算空间域和时间域的离散，此时的待优化变量取为各时刻受热边界节点上的热流值，即

其中Nj为受热边界离散单元节点数，Np为时间方向离散层数。则式（3）的优化变成了参数优化问题，可以采用梯度类优化算法来进行处理。采用伴随方程推导的方法来求梯度，首先将目标函数写为如下的扩展形式：

式中λ称为伴随变量。对（5）式进行分部积分再做变分后［8，10］得伴随变量应满足：

时域边界条件：t＝tf，λ（tf）＝0。

同时，可导出目标函数对Q的导数为：

由式（6）、（7）求出梯度值后，即可用最速下降法或共轭梯度法进行优化，注意为克服不适定性的影响，当计算结果满足如下收敛准则时终止优化迭代计算：

式中σ为温度测量值的标准差。

3 仿真算例

3.1 X－43飞行器的流动转捩分析

X－43飞行器马赫7飞行试验中，T021点位置表面热流随时间变化历程如图4所示［11］。从图中可以看出T021点在哪些时间段处于层流状态，哪些时间段处于湍流状态，由此可将多个测点的信息组合起来进行转捩位置判断。现考虑采用表面热流辨识方法，在飞行器壳体内壁测量温度变化历程，将壳体内的传热简化为一维传热模型，由温度历程来辨识出表面热流，进而对表面流态做出判断。设防热层是厚度10mm的不锈钢，利用图4中的热流仿真计算出内壁测点温度，考虑2%的测量噪声，得到测点温度历程如图5所示；再由此温度历程来辨识表面热流，得到的辨识结果与图4中给定热流的比较如图6所示，从中可以看到，辨识结果较好地再现了表面热流的变化情况，能够用于各时刻表面流态判读。

3.2 二维平板流动转捩位置测量仿真分析

二维平板模型及其表面上有转捩时的表面热流分布示意图如图3所示，设平板模型为长1m厚 10mm的合金，给定时空变化热流如图7所示。可以看到，t＝30s后表面热流逐渐达到稳定状态。针对计算域生成如图8所示的有限元网格，利用给定热流可以计算出内壁单元节点的温度变化历程。

图7 二维平板转捩分析的给定热流时空分布Fig.7 Temporal and spatial function of twodimensional plate surface heat flux

图8 有限元计算网格Fig.8 Computational grid of Finite Element Method

分别考虑内壁x方向均匀安装10个和20个测点的情况，利用这些测点的温度测量历程，用上一节中的方法来辨识表面热流，在t＝30s时刻得出的x方向热流辨识结果如图9中的“Estimated Q”所示，并与给定值“Exact Q”进行了比较。可以看到，在安装20个测点的情况下，辨识结果能较好地反映出表面热流的变化情况，对转捩位置进行较好的预测。

图9 给定热流与辨识结果比较（不考虑测量误差）Fig.9 Comparison of surface heat flux（no measurement noise）

进一步在测量结果中考虑1%的随机测量误差，图10给出了此时的辨识结果。可以看出：在此情况下，安装20个测点的辨识结果仍能较好地再现出表面热流的变化情况，对转捩位置进行较好预测。

此外，针对这一表面热流辨识问题，也可以考虑在20个测点位置都假设为一维传热问题，用一维辨识方法来进行计算，图11中的虚线给出了这样处理得出的t＝30s时刻x方向表面热流分布情况，与给得出的t＝30s时刻x方向表面热流分布情况，与给定值、二维传热辨识结果有较明显差异，得出的转捩区域明显大于真实情况。因此，此时采用基于二维传热模型的表面热流辨识方法是有必要的。

图10 给定热流与辨识结果比较（考虑1%测量误差）Fig.10 Comparison of surface heat flux（with measurement noise of 1%relative error）

图11 一维和二维传热模型辨识结果比较Fig.11 Comparison of estimated heat flux for one－and two－dimensional heat conduction model

4 结 论

介绍了风洞实验和飞行试验中测量转捩的主要方法及其优缺点。提出了基于表面热流辨识技术确定转捩位置的基本思想和处理方法，与目前工程上主要采用的近壁热电偶测温法相比，表面热流辨识方法的优势在于测点可以距离受热表面相对较远，从而放宽对温度传感器量程和结构强度的要求。但测点的位置也不是没有限制，因为测点越远离受热面，辨识问题的不适定性越强，因此需要采用仿真辨识方法对传感器安装位置进行合理选取。在给出二维传热模型表面热流辨识算法的基础上，对两个算例进行了仿真辨识分析。结果表明：基于表面热流辨识技术确定转捩位置是可行的，能给出较为准确的转捩区域判断。

［1］ 卞於中，徐铁军，范洁川，等.边界层综合诊断技术研究［J］.流体力学实验与测量，2000，14（4）：7－13.

［2］ 胡成行，黄叙辉，李红梅，等.应用脉动压力测试技术探测边界层转捩［J］.流体力学实验与测量，2002，16（2）：67－71.

［3］ KUNTZ D W，POTTER D L.Boundary layer transition and hypersonic flight testing［R］.AIAA－2007－308，2007.

［4］ ILIFF K W，SHAFER M F.A comparison of hypersonic vehicle flight and prediction results［R］.NASA TM－104313，1995.

［5］ WRIGHT R L，ZOBY E V.Flight boundary layer transition measurements on a slender cone at Mach 20［R］.AIAA－77－719，1977.

［6］ HOWARD F G.Single－thermocouple method for determining heat flux to a thermally thick wall［R］.NASATND－4737，1968.

［7］ 钱炜祺，何开锋，桂业伟，等.飞行器气动热参数辨识研究进展［A］.第六届全国流体力学青年研讨会论文集［C］，浙江杭州，2009：79－84.

［8］ ALIFANOV O M.Inverse heat transfer problems［M］.Springer－Verlag，Berlin，1994.

［9］ 钱炜祺，蔡金狮.再入航天飞机表面热流密度辨识［J］.宇航学报，2000，21（4）：1－6.

［10］钱炜祺，何开锋，桂业伟，等.非稳态表面热流反演算法研究［J］.空气动力学学报，2010，28（2）：155－161.

［11］BERRY S，DARYABEIGI K，WURSTER K，et al.Boundary layer transition on X－43A［R］.AIAA－2008－3736，2008.