0   引 言

复杂外形飞行器的高超声速流场中存在复杂的流场结构,如激波、激波边界层相互作用等。对其进行深入研究,不仅有利于加深对其流动机理的认识,而且对工程应用有着重要的指导意义^[1-2]。

为了能在地面风洞试验中对这些现象进行深入研究,先进的流场显示和测量技术尤为重要。其中,随着激光技术、计算机技术和信息处理技术的发展,基于示踪物的光学非接触测量技术成为现代空气动力学研究的重要手段。PIV (Particle Image Velocimetry)速度场测试技术^[3-4]、平面激光诱导荧光(PlanarLaserInduced Fluorescence,PLIF)技术^[5-6]等由于各自的原因在超声速、高超声速流场测试中难以达到理想的效果。PIV速度场测试技术使用微米级尺寸的粒子作为示踪物，不能满足超声速、高超声速流场对Induced Fluorescence,PLIF)技术^[5-6]等由于各自的原因在超声速、高超声速流场测试中难以达到理想的效果。PIV 速度场测试技术使用微米级尺寸的粒子作为示踪物,不能满足超声速、高超声速流场对示踪粒子跟随性的要求。而PLIF技术使用气体分子散射则不存在粒子跟随性难题，但其信号非常弱，信噪比很低。这些都限制了上述测量技术在超声速流场测试中的应用。

为了在超声速流场中进行高分辨率的流场显示，易仕和等发展了1种基于纳米粒子的超声速流动成像NPLS(Nano-Tracer Planar Laser Scattering)流场显示和测量技术^[7-10]，NPLS以纳米粒子作为示踪粒子，以脉冲平面激光作为光源，通过CCD记录流场中的粒子图像实现超声速流动的高分辨率成像。NPLS 方法采用的纳米示踪粒子，体积远小于常规PIV 的示踪粒子,解决了示踪粒子在超声速流场中的跟随性难题。另一方面，纳米示踪粒子的散射信号比PLIF 分子散射要强几个量级，而且通过激光偏振调节技术，实现信噪比大幅增强。 所以，利用常规CCD 相机即可获取高质量的流场结构图像，实现了流动结构高分辨率可视化。通过对流场结构NPLS 图像的校准和计算，可以得到超声速流动高分辨率的密度场和速度场等流场参数。目前，NPLS 技术已成功应用于激波边界层干扰、超声速混合层、超声速绕流和超声速边界层等多种流场精细结构的测量研究。

李明等把NPLS技术推广应用到高超声速流场显示中^[11]，开展高超声速条件下激波/边界层干扰的试验研究。结果表明，在其试验条件下，粒子散射光的信号与模型表面反射光的信号处于相当的水平，在来流马赫数为10的条件下，部分区域表面反射光的信号更强，这给图像信号处理带来一定的困难。这种现象1个可能的原因是在高超声速条件下纳米粒子的跟随性不好，特别是在流场速度出现急剧变化的区域，纳米粒子并不能实时跟随气流变化。

对于纳米粒子的跟随性，一般采用粒子松弛时间或响应时间来表征^[12]:

式中：ρ_p为粒子密度；d_p为粒子直径；μ为流体的粘性系数；Kn_d=λ/d_p 为粒子克努森数，λ为分子平均自由程。

粒子的松弛时间越短，在超声速流动中的跟随性就越好。在多相流体动力学中，定义斯托克斯数St 来衡量粒子跟随性，St=τ_p/τ_f，τ_f为流动特征时间。当斯托克斯数很小时，粒子具有足够的时间去响应流场的变化，否则就可能无法准确反映实际的流场结构。 为了研究在不同稀薄度流场条件下纳米粒子的跟随性，本文采用数值仿真方法，把纳米粒子作为固体颗粒，仿真气-固两相流流场，研究纳米粒子在气体流场中的跟随性。

1   数值模拟算法

数值仿真采用基于DSMC^[13](Direct Simulation Monte Carlo)的稀薄两相流方法。该方法采用有限数目的仿真分子模拟实际流场中数目巨大的真实分子，通过跟踪流场中仿真分子的运动和分子间的碰撞达到流场模拟的目的。该方法广泛用于模拟稀薄气体的流动，具有很高的精确度。Gallis提出了1种改进的DSMC方法^[14]。他利用Green 函数发展的DSMC 方法适用于求解在任意分子速度分布的气相流场中颗粒所受力和热的特点，可以模拟包括稀薄和化学惰性固体颗粒相在内的稀薄流动。后来经Burt、Boyd等人的发展，建立了1种适用于DSMC方法的双向耦合算法(Two-way Coupled)，既考虑气相对固相的力和热的作用，又考虑固相颗粒对气相的作用，能够准确描述固相颗粒在稀薄过渡流中的输运过程^[14-16]。

1.1   DSMC方法

对气相流场，采用DSMC方法。在DSMC方法中，采用变刚球分子模型，能量交换采用L-B模型。网格采用二级直角网格，碰撞网格根据密度自适应，碰撞对的选取限制在碰撞网格内。为了提高计算效率，本文采用了基于MPI的并行计算技术，将计算区域分解为N个子区域，分别分配给N个节点。各节点进行独立且完整的DSMC仿真。仿真分子穿过边界区域时，计算区域之间进行数据交换^[13]。

1.2   两相流双向耦合算法

两相流DSMC模拟运算法则是基于相间动量和能量从微粒特性暂时变化传递的解耦，即把气相对固相的作用和固相对气相的作用分别处理。

第一步，考虑气相对固体颗粒的作用。假设固体颗粒处于当地自由分子流的状态，从颗粒表面反射的气体分子与颗粒周围的气体分子不发生碰撞(碰撞在气相DSMC方法中处理(见1.1节))，在颗粒周围不会形成流场结构。同时不考虑多原子气体的振动激发，在同一个网格里每个DSMC气体仿真分子作用到一个固体颗粒上的动量和能量分别为^[4]：

式中：R_p为等效颗粒半径；N_g为每个计算分子所模拟的真实气体分子数目；τ为颗粒表面导热调节系数；V_c为网格体积；m为单个气体分子质量；u_r为气体分子与相关颗粒的相对速度，c_r是u_r的绝对值；k为Boltzmann常数；T_p为颗粒温度；Λ为气体转动自由度数；e_rot为单个分子转动能。

求和网格内所有仿真分子对颗粒的作用，可以得到时间步长内颗粒的动量和能量变化。

第二步，考虑固体颗粒对周围气体的影响。首先要确定在每个时间步长内哪个仿真分子将与颗粒进行碰撞。对Bird的非时间计数方法进行修正，来确定与所选的颗粒可能发生碰撞的计算分子数n_s。

式中：N_p为1个仿真颗粒所表示的实际固体颗粒的数量；n_g为与固体颗粒在同一网格里的气体仿真分子的数量；Δt为时间步长；(c_r)_max为网格内采样到的分子-颗粒对碰撞前最大相对速度。

对于这个颗粒，给定一个气体仿真分子与其发生碰撞，要么为以概率等于颗粒热适应系数τ的等温漫反射碰撞，要么为概率为1-τ的镜面反射。如果发生镜面反射，则相对速度c_r在碰撞中不发生改变，碰撞后的相对速度u_r^*可通过c_r与单位矢量n相乘得到。如果发生漫反射，碰撞后相对速度u_r^*围绕初始相对速度u_r的方位角ε在[0,2π]上等概率分布。在漫反射碰撞中，碰撞后相对速度c_r^*不能假定为等于初始相对速度c_r，而是需要通过使用“取舍”法，从如下分布函数来确定c_r^*的值

式中：β为气体在颗粒温度处最可几速度的倒数，β=[m/(2kT_p)]^1/2。

对于漫反射多原子分子气体，碰撞后转动能e_rot也必须改变。漫反射双原子气体分子的转动能可计算如下

式中：R_f为(0,1)之间的一个随机数。

整体坐标系下u_r的分量为u_r、v_r和w_r。采用Bird二元弹性碰撞，相对速度u_r^*的各分量可以由u_r、v_r、w_r、c_r和c_r^*以及角δ和ε计算得到。有

式中：δ为碰撞偏转角，定义为u_r与碰撞后的相对速度矢量u_r^*(u_r^*=u_m^*u_p^*)之间的夹角，这里u_m^*和u_p^*分别为下面碰撞时气体分子和颗粒的绝对速度。

镜面反射碰撞的偏转角分布为

漫反射碰撞的偏转角分布为

可以通过DSMC方法中常用的“取舍”法得到δ角。

碰撞后气体分子的绝对速度为

2   计算结果与讨论

本文仿真计算了3种不同稀薄程度的流场，流场条件如表 1所示。假设流场中纳米粒子在初始流场中分布均匀，并且跟随性很好，即纳米粒子与气流速度、温度一致。通过仿真考察两相流场在与模型相互作用时，颗粒分布与气相流场分布能否达到一致。

表  1  流场计算条件

Table  1  Flow conditions to be simulated

状态	1	2	3
Ma	3.0	5.0	10.0
数密度n/(1024·m-3)	4.4140	0.2806	0.0652
Knd	5.86	92.26	397.05
τp/μs	1.33	44.75	424.78
St	4.07	159	2215

下载: 导出CSV 

| 显示表格

流场气体为N₂，纳米粒子材料为试验中采用的TiO₂，纳米粒子的直径为50nm。纳米材料的密度为3840kg/m³。播撒在流场中纳米粒子数密度均为4×10¹¹/m³。

2.1   球锥绕流

本文首先计算了文献11]中的球锥绕流流场。计算条件为表 1中的状态2。球锥底部直径为80mm，长度为180mm。

图 1是计算结果与试验流场照片的对比。可以看出，计算结果的流场结构与试验照片一致。计算结果中，气相流场与纳米颗粒分布也基本一致。

图  1  钝锥绕流流场比较

Fig.  1  Comparison of bluntcone flowfield

下载: 全尺寸图片 幻灯片

图 2是x=165mm处气相和纳米粒子的密度分布比较，从图中可以看出，两者符合很好，表明在这种流场条件下，纳米粒子的跟随性很好，其分布能够反应流场的特性。

图  2  x=165mm截面上密度分布比较

Fig.  2  Comparison of number density distribution at x=165mm

下载: 全尺寸图片 幻灯片

2.2   压缩拐角流动

为了考察在复杂流场情况下纳米粒子的跟随性，本文计算了压缩拐角的流动。计算模型为二维带压缩拐角的平板(见图 3)。模型总长为150mm，其中水平部分为100mm，拐角为15°。计算中气体分子在模型表面进行漫反射，纳米颗粒进行镜面反射。

图  3  模型及尺寸(mm)

Fig.  3  Model configuration and scales (mm)

下载: 全尺寸图片 幻灯片

图 4是状态1流场数密度分布。流场在平板前缘形成贴体斜激波，在压缩拐角处，气体受到较强压缩，形成较大的高密度区域。对比气相和纳米粒子分布看出二者基本一致，特别是拐角处涡流结构(见图 5)也吻合得很好。

图  4  流场密度分布(状态1)

Fig.  4  Number density distributions (case 1)

下载: 全尺寸图片 幻灯片

图  5  拐角区域流线(状态1)

Fig.  5  Vortex structures at corner (case 1)

下载: 全尺寸图片 幻灯片

图 6是状态1在x=100mm处流场参数沿y向的分布比较。可以看出，气相和纳米粒子密度和速度的分布基本完全重合，这样，通过纳米粒子散射得到的图像能够准确地反映气体流场结构。

图  6  x=100mm截面上参数分布比较(状态1)

Fig.  6  Comparison of parameters distribution at x=100mm (case 1)

下载: 全尺寸图片 幻灯片

图 7是状态2的流场数密度分布比较。在这种条件下，纳米粒子与气相流场分布有明显差异。主要体现在压缩拐角区域。该区域内气体密度相对较低，而纳米粒子密度相对较高，这样纳米粒子的分布并不能反映气体流场的分布。原因分析如下：在压缩平板壁面附近边界层内，气体密度较低，纳米粒子在这个区域不能完全跟随气体的运动，在壁面附近聚集了较多的纳米粒子；而在压缩拐角附近，气相流场有涡产生，部分纳米粒子会随涡运动到拐角区域(见图 8)，拐角区域气体数密度也较低，纳米粒子会聚集在该区域，而不能完全随气体运动。

图  7  流场密度分布(状态2)

Fig.  7  Number density distributions (case 2)

下载: 全尺寸图片 幻灯片

图  8  拐角区域流线(状态2)

Fig.  8  Vortex structures at corner (case 2)

下载: 全尺寸图片 幻灯片

图 9是状态2在x=100mm处流场参数沿y向的分布比较。图中显示，2者数密度差别较大，由前文分析可知，在拐角壁面附近纳米粒子数密度很高。轴向速度比较接近，在拐角区域纳米粒子的速度比气体速度稍高，说明纳米粒子速度的改变落后于气相速度的改变。纳米粒子与气相法向速度分布规律类似，但纳米粒子速度较低，不能完全跟随气体的运动。这些比较表明，在状态2条件下纳米粒子的跟随性不好。

图  9  x=100mm截面上参数分布比较(状态2)

Fig.  9  Comparison of parameters distribution at x=100mm (case 2)

下载: 全尺寸图片 幻灯片

图 10是状态3的流场数密度分布。在这种条件下，纳米粒子与气相流场差别很明显。完全不能反映气相流场的流场结构。图 11中x=100mm处的流场参数对比显示，在此处纳米粒子的密度与气相的密度分布已经没有可比性。纳米粒子轴向的速度比气相高出很多，表明纳米粒子完全跟不上气相速度的改变。纳米粒子的法向速度低于气相速度，而且分布规律也不相同。在这个状态下，纳米粒子不能跟随气相流动的变化而变化，那么前文假设来流条件中纳米粒子与气相速度一致的条件有可能不成立。

图  10  流场密度分布(状态3)

Fig.  10  Number density distributions (case 3)

下载: 全尺寸图片 幻灯片

图  11  x=100mm截面上参数分布比较(状态3)

Fig.  11  Comparison of parameters distribution at x=100mm (case 3)

下载: 全尺寸图片 幻灯片

上述算例表明，在气相流场中，纳米粒子的跟随性与流场稀薄条件关系很大。表 1中St数中的特征时间是以网格尺度和来流速度为基础计算的。在St≤4时，纳米粒子跟随性很好，纳米粒子的分布基本上反映了气相流场的结构。在更大的St数条件下，纳米粒子的跟随性会随着St数的增加变得越来越差，导致纳米粒子在局部或整体流场中与气相流场的分布不一致，不能反映气相流场的结构。NPLS测量技术在这种条件下得到的激光散射信号不能反映气相流场结构。

由公式(1)可知，为了降低颗粒的响应时间，提高颗粒的跟随性，拓展NPLS技术在高超声速低密度流场中的应用，需要采用更小尺寸的纳米粒子。但降低粒子尺寸，会降低粒子对激光的散射信号。

另一种方法是采用密度较小的粒子，也可以提高颗粒的跟随性。由公式(2)可以推出，纳米粒子在气相流场中受力后加速度与粒子的半径和密度的乘积成反比(a~1/(R_pρ_p))。减小粒子密度，可以增加粒子加速度，加快粒子反应。从公式(1)中也可以看出，粒子松弛时间与粒子密度成正比。在颗粒尺寸受限的情况下，可以降低粒子的密度，以缩短粒子松弛时间，提高粒子跟随性。

3   结束语

本文采用基于DSMC方法的稀薄两相流的双向耦合算法，对NPLS测量技术在高超声速流场测量中TiO₂纳米粒子的跟随性进行了数值研究。研究表明，在高超声速流场条件下，伴随着流场稀薄度的增加，纳米粒子在流场的跟随性会降低，当流场中出现复杂流场结构时，TiO₂纳米粒子不能适应气相流场，建立与之相应的分布，这样限制了NPLS技术在高超声速低密度流场测量中的应用。