目前，壁面放电等离子体作为一种流动控制技术受到广泛的关注. 它可以用来提高超音速喷管和超音速近壁面流动的气动性能，例如减小阻力^{[1, 2, 3]}，控制分离^{[4, 5]}，改变进气道内的激波结构^{[6, 7, 8, 9]}以及提高燃烧室内燃气的混合度^{[10, 11, 12]}.基于壁面放电等离子体的流动控制机理分为两大类：一类是基于磁流体动力学的基本观点，认为其控制效应是非热效应；另一类是基于放电的焦耳热效应的基本观点，认为控制机理是热效应.Webb等^[13]通过实验研究了壁面电弧放电对激波边界层干扰而引起边界层的分离控制机理，结果显示电弧放电对激波和边界层分离控制机理主要为放电热效应，边界层内放电加热导致边界层特性改变，进而影响整个流场结构. Leonov和Yarantsev^[14]针对壁面放电等离子体在超音速流动中的控制机理做了实验方面的研究，结果表明：在无外磁场作用的情况下，电弧等离子体流动控制的效应为电流的热效应，具体表现为在放电区下游区域形成热电弧层，该区域内气流温度显著升高，气流密度显著降低，边界层厚度增大并产生诱导激波，最终改变核心流场的流动参数.

壁面放电等离子体在流动控制方面的实验研究已经取得了长足的发展，但通过实验研究等离子体在流动控制中的应用也存在以下几方面的不足：所需的成本太高，实验台的搭建需要大量的财力和时间，而且每次实验需要大量的能源；由于实际条件限制无法得到流场中任意点上的任意参数值，并且由于试验时间非常短，实验结果记录有一定难度，并且对记录数据的仪器要求很高.因此，通过数值模拟来研究放电等离子体的流动控制是一种不可缺少的方法^{[15, 16, 17]}.

高超音速进气道是新一代冲压发动机的重要部件之一，其利用超音速气流经过斜压缩面对气流增压，因此没有多余压气部件，结构简单.但超音速进气道内含有复杂的激波系结构，在非设计点工作时，进气道内的激波系会造成流量系数和总压恢复系数的下降^{[18, 19, 20]}.

Yan^[21]在2马赫 超音速进气道内，通过数值模拟的方式对壁面放电等离子体的流动控制原理进行了研究. 结果显示，壁面放电加热气流可以诱导产生斜激波，从而改变放电区域下游流场中的激波结构，并且放电能量越大，诱导激波越强，流场激波结构的改变越明显. 但是，随着放电能量的增大，放电加热区域局部的气流温度将变得非常高.根据Leonov等^[22]的研究，对于仅仅在2 kW功率的放电加热下，局部气流的温度就高达5 000 K，以上这样的温度足以将空气中的氧气分子离解(O$_{2} \to $2O)，在更高的能量输入下，空气中氮气也将会被离解(N$_{2} \to $2N).当温度更高时，气体将发生电离过程(N$ \to $N$^{ + }$+e$^{ - }$,O$ \to $O$^{ +}$+e$^{-})$，伴随着各种粒子之间非常复杂的电离和复合过程^{[23, 24]}.数值模拟电离过程是存在很大困难的，需要考虑不同粒子种类以及其反应平衡模型.对于空气，为了得到较准确的结果，通常需要考虑11种成份(O$_{2}$，N$_{2}$,O，N,NO,O$^{ + }$,N$^{ + }$,O$_{2}^{ + }$,N$_{2}^{ + }$,NO$^{ + }$,e$^{ - })$,这样就大大增加了计算量，对于工程问题而言，这样的方法几乎是无法进行的.

根据Leonov和Yarantsev^[14]研究超音速流场中壁面放电等离子体的实验结果，显示当放电输入能量小于2 kW时，电弧内气流最高温度在3 000 $\sim $5 000 K之间，本文中考虑的放电功率分别为0.9 kW和1.8 kW，计算采用定等压比热不太合理了，所以必须要考虑等压比热随温度非线性变化的影响.表1中可以直观的看出，当气体温度小于9 000 K时，空气的电离度很低，只有万分之一的量级，因此，假设空气组分没有发生改变是合理的，即可以认为气体常数$R$为定值，依旧可以使用理想气体模型^{[26, 27]}.

表 1(Table 1)

表 1 粒子数目与气体温度的对应关系[27]Table 1 Number densities as a function of gas temperature[27]

表 1 粒子数目与气体温度的对应关系[27] Table 1 Number densities as a function of gas temperature[27]

作为气体动力学中的一个重要影响参数比热比$\kappa $，由于流场中的局部范围内气流的温度很高，比热比$\kappa$不再是一个 常数，而比热比与等压比热存在固定的换算关系，当$R$保持不变时，比热比随着等压比热发生变化.

本文数值模拟了5马赫的超音速进气道内壁面放电对激波的控制过程，分别采用定等压比热和变等压比热进行计算，其中变等压比热采用了Balakrishnan^[28]得到的等压比热$C_p$随温度的拟合关系式.

根据Leonov和Yarantsev^[14]研究超音速流场中壁面放电等离子体的实验结果，超音速气流内的等离子电弧结构如图1所示，电弧沿流向为长条状，电弧内气流被加热并膨胀.随着输入电功率增大，电弧能量密度增强，而电弧的高度和长度基本保持在一定范围内.在数值模拟壁面放电流动控制过程中，电弧加热气流过程模拟的是否准确是关注的重点.实验结果显示在输入能量小于2 kW时，电弧内气流最高温度在3 000 K$\sim$5 000 K之间，并且随着输入功率的增大，电弧内气流温度逐渐增加.理想状态下电弧内气流输入能量与温升的关系为

其中，${\rm{d}} W$为单位时间内的输入能量，$m$是放电加热层的气流总质量，${\rm{d}} T$为单位时间内放电加热层气体的温升. 可以看出：(1)当不考虑等压比热$C_p$的变化时，放电加热区域的温升与输入能量呈线性关系；(2)当考虑$C_p$随温度的变化，放电加热区域的温升与输入能量呈非线性关系. 因此，这两种等压比热数值模拟的结果会存在一定的差异.

图 1(Fig. 1)

图 1 超音速气流中壁面放电图[14] Fig. 1 Surface discharge in supersonic flow from literature[14]

根据$C_p$与温度的拟合关系式^[28]

其中,参数$a,b$和$c$是分段给出的，具体值见表2. 图2绘出了$C_p$随温度的变化曲线，并且对比了变$C_p$与定$C_p$ (该值为1.006 kJ/(kg$\cdot$K))的区别.

表 2(Table 2)

表 2 拟合多项式的参数Table 2 The parameters of the fitting polynomial

表 2 拟合多项式的参数 Table 2 The parameters of the fitting polynomial

图 2(Fig. 2)

图 2 两种等压比热随温度变化关系图，其中变等压比热数据来自文献[28] Fig. 2 Constant Cp and varied Cp from literature[28]

5马赫 的二元超音速进气道结构如图3所示，其中$A(0,0)$点为坐标原点,$HA=0.08$ m,$HG=0.06$ m,$DE =0.08$ m. 其中$HG$为入口边界，$DE$为出口边界，$EFG$和$HABCD$分别为上下壁面边界.进气道下壁面AB和BC与水平面的夹角分别为$7^ \circ $和$14^ \circ $.虚线$AF$和$BF$分别是由两个斜坡面压缩气流而产生的斜激波.

图 3(Fig. 3)

图 3 5 马赫进气道结构简图 Fig. 3 The Mach 5 inlet structure

计算域为三维结构，展向长度为0.009 m. 热激励器为一个长方体结构，热激励器的展向长度为0.003 m，热激励器处于展向中心位置. $Z=0$ m和$Z=0.009$ m分别为计算域的两个侧面，这两个面采用对称边界条件，这样可有效的减少计算资源.

采用在能量方程中添加热源项来模拟放电等离子体的加热效应，热源项大小与放电功率相等. 三维非定常Navier-Stokes方程如下

气体状态方程

其中，$e = C_p T -RT $,$R$为空气的气体常数，$C_p$为等压比热容，$S$为热激励器输入的能量项，其表达式为

其中，$E$为壁面放电功率(单位:W)，$x_1$,$x_2$,$y_1$,$y_2$,$z_1$和$z_2$是长方体热激励器尺寸的参数，分别是$-40$ mm,$-20$ mm,0,3 mm,1.5 mm和4.5 mm. 文中分别考虑了3个大小的$E$分别为0,0.9 kW和1.8 kW，对应的$S$分别为 0,5 GW/m$^3$,和10 GW/m$^3$. 其中，$E=0$表示无热激励作用.

根据飞行走廊图,马赫数为5的吸气式飞行器可以在15 km左右的高度飞行. 计算中采用15 km左右高空的气流参数:为了计算方便起见，静压选用$p_\infty = 1.2\times 10^4$ Pa,静温采用$T_\infty = 216$ K.入口给定总压、总温以及马赫数，出口边界上的值由内场一阶插值得到.上下壁面均为采用绝热无滑移壁面条件，展向使用对称边界条件. 采用适用于计算强逆压流以及分离流的K-WSST湍流模型^[29]，采用Sutherland插值公式来考虑黏性随温度的变化^[30].时间精度采用二阶隐格式，时间步长选择为$5\times 10^{ - 8}$ s，满足CFL数为0.5，空间上采用三阶MUSCL (monotonicupwind scheme for conservation laws)格式，通量项离散采用 AUSM (advection upstream splittingmethod)格式^[31]. 由于流场中加入稳态热源项，因此我们将对最终的准稳态流场进行分析.

对加热区域和边界层进行了网格加密. 如图4所示为计算域的网格平面图.展向的网格尺寸保持不变，$\Delta z$恒为0.000 25 m.

图 4(Fig. 4)

图 4 计算域的平面网格图 Fig. 4 Two dimensional grid

为了消除网格对于计算结果的影响，在热源$E=0.9$ kW的情况下进行了网格相关性的研究. 5种网格中均保持$x$和$z$方向的网格数分别为388和18，$y$方向在壁面处加密，拉伸比为1.02. 表3给出了这5种网格的一些特征量的值.这些值反映了边界层和热源处网格的精度，其中$\Delta y$是下壁面的第一层网格的高度，$\Delta x$为加热区域$X$方向最小网格尺寸. 5种网格的变化趋势是：网格1,2和3保持$\Delta x$不变，逐渐减小$\Delta y$的值；网格4相比于网格3，$\Delta x$减小到一半而$\Delta y$减小到1/4；网格5相比网格4，保持$\Delta x$不变，$\Delta y$减小一半.

表 3(Table 3)

表 3 5 种网格的特征尺寸量Table 3 The feature sizes of the five grids

表 3 5 种网格的特征尺寸量 Table 3 The feature sizes of the five grids

下面给出了5种网格的计算结果，图5为流场中加热区域下壁面的压力分布图(压力值经过来流静压无量纲化). 可以看出，在加热区域前的压力分布十分依赖于$\Delta y$的值，随着$\Delta y$的逐渐变小，压力逐渐趋于均匀.这5种网格均可以扑捉到加热区域压力的变化，并且网格精度较高的4和5两种网格的压力分布已经非常接近，相比于前面3种比较粗糙的网格，网格4和5可以较好地捕捉到激波的产生位置和加热区域附近的膨胀波系.

图 5(Fig. 5)

图 5 进气道前段下壁面压力的分布图 Fig. 5 The bottom wall pressure distribution in the front of the inlet

图6是垂直壁面方向上第一个$\Delta Y^{ + }$的分布图，根据K-WSST湍流模型的要求，$\Delta Y^{ + }$的值要在1的量级，由图6看出只有网格5满足该条件. 因此，本文采用网格5进行数值计算.

图 6(Fig. 6)

图 6 5 种网格下壁面ΔY+ 值的分布图 Fig. 6 The bottom wall ΔY+ value distribution under the five grids

图7为中心面上温度云图，温度以来流静温为特征量进行无量纲处理. 热源项为0表示原始流场，热源项不为0表示有热激励的流场，可以看出热源项的值越大，下壁面附近的高温区域就越大，流场的最高温度也越高.高温层仅仅处于下壁面附近，主流区域的温度改变不明显.

图 7(Fig. 7)

图 7 温度云图 Fig. 7 Temperature contours

图8为中心面上压力云图，压力以来流静压为特征量进行了无量纲处理.结果显示，在热源的前端，产生了一道诱导激波，使原 始激波的激波角增大，激波强度降低.随着热源项值的增大，诱导激波强度增大，原始激波的改变越明显.激波结构的改变可以归结为以下两方面原因：一方面，气流经过诱导激波后，马赫数会降低；另一方面，压缩角处的边界层由于高温而变厚，降低了压缩角的压缩效果. 在这两个因素的共同作用下，原来由于压缩坡面产生的斜激波的强度降低，激波角增大.

图 8(Fig. 8)

图 8 压力云图 Fig. 8 Pressure contours

图9以压力等值面的形式给出了流场展向结构，压力等值面可以被看作为热源诱导激波的激波面，在热源附近压力等值面受热源形状影响非常明显，而远离热源压力等值面趋于平面.加入热源后，压力等值面发生明显的变化，但是流场的展向结构不会随着热源的能量大小发生明显的变化，并且等压比热值对流场展向结构的影响很小.

图 9(Fig. 9)

图 9 p = 20 kPa 的压力等值面 Fig. 9 Pressure iso-suface，p = 20 kPa

图10为中心面上马赫数云图，可以清晰地看到激波和膨胀波的位置以及边界层区域的流动.随着热源项的加入，加热区域下游边界层附近马赫数迅速降低，这主要因为热源项的能量大部分转化为气流的内能，使气流温度显著升高，气流速度的变化相比于温度而言很小，故马赫数的变化主要受温度的影响，随着温度的升高而降低.由于边界层内的温度很高，所以边界层内马赫数非常低，并且热源项越大这种现象越明显.但是这样的影响仅仅发生在近壁区域，而远离壁面的区域，气流马赫数没有发生明显的变化.

图 10(Fig. 10)

图 10 马赫数云图 Fig. 10 Mach number contours

通过前面的分析发现，在流场中加入热源来改变流场结构，是完全可以实现的，并且其效果也是相当明显的.在流场中加入热源会使局部温度急剧升高，为了能更真实地模拟热激励流动控制过程，必须要考虑到高温下气体分子离解和电离的过程，宏观上表现为输运系数和热力学参数是温度的函数. 本文重点研究了热力学参数等压比热对计算结果的影响.下面将详细地分析不同热源能量下，定等压比热与变等压比热的计算结果.图11分别给出了不同热源能量下，2种等压比热得到的下壁面中心线上温度分布曲线.结果显示：变等压比热的温度始终低于定等压比热的温度，这是由于在相同的温度下，变等压比热始终大于等于定等压比热的值，因此变等压比热的温升低于定等压比热的温升.当热源能量越高时，这种差异就越明显(见图11(b)和11(c))，由图2可知，温度越高，变等压比热与定等压比热之间的差异越明显，当温度大于10 000 K后，两者的差别在2到3倍. 这样的差异会使气流最高温度的差别非常大.

图 11(Fig. 11)

图 11 下壁面温度分布图 Fig. 11 The temperature distribution of bottom wall

图12分别给出了不同热源能量下，下壁面中心线上压力分布曲线，其中对比了两种等压比热的计算结果.图11(a)为热源能量为0时，两种等压比热下的压力分布，由于流场中的温度较低，变等压比热与定等压比热的值非常接近，故压力变化非常接近.图12(b)为热源能量较低时，2种等压比热下的壁面压力分布，可以看出在加热区域($-0.04$ m$ < x < -0.02$ m)压力值有一定的差别，而其他区域的压力分布几乎一致.对比图12(c)与12(b)发现，在高热源能量下，变等压比热与定等压比热的结果之间的差异较为明显.原因在于，高的热源能量会得到高的气流温度，而高气流温度导致高等压比热值，最终会导致变等压比热与定等压比热之间的差异增大.因此，随着热源能量的增大，两种等压比热下的结果差异越来越明显.

图 12(Fig. 12)

图 12 下壁面压力分布图 Fig. 12 The pressure distribution of bottom wall

图13为中心面$y=0.04$ m上的温度分布曲线图，图13(a)表明，气流经过原始斜激波后，虽然温度升高，但是温度上升的幅度不足以使$C_p$发生明显的改变. 因此变等压比热与定等压比热的结果相差很小.图13(b)和13(c)表明随着热源能量的加入，激波位置和强度都发生了改变，采用定等压比热得到的温度高于变等压比热的结果，也就是采用定等压比热会过高的估计诱导激波的强度以及流场的改变幅度.

图 13(Fig. 13)

图 13 直线y = 0:04m 上的温度分布图 Fig. 13 The temperature distribution along the straight line y = 0:04m

图14为中心面$y=0.04$ m上压力分布曲线，图14(a)表明，无能量加入时，两种等压比热得到的压力分布一致.图14(b)和14(c)显示热源能量不为零时，使用定等压比热会过高地估计诱导激波的强度，不能准确确定原始激波的改变幅度.随着热源能量的增大，流场中诱导激波增强，原始激波结构的改变会愈明显，同时2种等压比热的计算结果之间的差别也会越来越大.

图 14(Fig. 14)

图 14 直线y = 0:04m 上压力分布图 Fig. 14 The pressure distribution along the straight line y = 0:04m

图15为中心面$y=0.04$ m上总压恢复系数分布曲线，其中图15(a)表明无能量加入时，两种等压比热下的结果均显示激波后总压下降，在激波前(低温区)总压恢复系数的分布一致，而在激波后(温度较高区域)两者的结果出现差异.图15(b)和15(c)表示加入热量后，总压恢复系数下降.对比发现，使用定等压比热会低估总压恢复系数的值，主要原因是：定等压比热的结果会过高估计流场的温度.一般而言，温度增高，总压下降. 随着热源能量的增加，总压恢复系数在降低.尽管高的热源能量可以带来强有力的激波变化，但激烈的热激励也带来了气流总压大幅降低的弊端.因此，应该权衡利弊，谨慎选择壁面放电功率的大小.

图 15(Fig. 15)

图 15 直线y = 0:04m 上总压恢复系数分布图 Fig. 15 The total pressure recovery coefficient distribution along the straight line y = 0:04m

出口截面的总压恢复系数是衡量气流品质的一个重要参数.表4为不同热源能量下，两种等压比热下的出口截面平均总压恢 复系数.可以看出(1)总压恢复系数低于无热源的情况，说明热源在增加了系统对激波控制能力的同时，对总压恢复系数具有一定负面影响；(2)出口截面的平均总压恢复系数随热源能量的增加而非单调变化，其主要原因：加入热源可以改变激波的强度，这里主要用来降低流场中的激波强度(使总压恢复系数增大)，但是加热又会产生诱导激波并且使气流总温增大(使总压恢复系数降低)，这就意味着加入热源既有增大总压恢复系数的效果，也有减小总压恢复系数的效果.总的效果如何，主要取决于这两种效果相互竞争的结果.因此才会出现非单调变化的趋势；(3)使用变等压比热得到的出口平均总压恢复系大于使用定等压比热的结果.

表 4(Table 4)

表 4 出口截面平均总压恢复系数Table 4 Surface averaged total pressure recovery coefficient at outlet

表 4 出口截面平均总压恢复系数 Table 4 Surface averaged total pressure recovery coefficient at outlet

本文采用数值方法研究了马赫数为5的超音速进气道内壁面放电等离子体在激波控制中的作用，重点研究了等压比热对模拟结果的影响. 具体过程是：采用在壁面附近处加入一长方体均匀热源模拟壁面放电的加热效应.通过数值计算模拟了壁面放电对超音速气流的热激励效应，发现热激励对于流场中激波的改变非常明显.同时发现加入热源后流场的局部温度升高很多，甚至高达50 000 K以上，这样的结果明显是不合理的.因此，文中进一步采用根据实验数据拟合的等压比热公式来代替定等压比热进行计算，结果发现，流场的局部高温下降到了一个可以接受的范围，同时由于温度的大幅下降，流场中的压力分布也发生相应的变化.

对比了在不同热源能量下，定等压比热和变等压比热的计算结果，发现：(1)相比于变等压比热的计算结果，定等压比热的预测结果中，流场温度偏高，同时流场中激波的变化幅度偏大.(2)随着热源能量的增大，定等压比热与变等压比热得到的结果之间的差异更加明显，特别的是，对于温度以及和温度关联比较密切的物理量.

文中采用的等压比热随温度的变化关系式是通过大量实验得到的拟合结果，因此在特定的条件下，可以反应真实的物理过程，得到比较准确的温度分布，并能够较为精确地捕捉到流场的变化，对准确理解和认识热激励流动控制的机理有着重要的影响.

Yan Hong, Wang Song