高压涡轮全环非定常流动数值模拟

邵 飞 ，綦 蕾 ，罗建桥 ，刘东华

（1.北京航空航天大学能源与动力工程学院，北京 100191；2.中航工业动力机械研究所，湖南株洲 412002）

0 引言

涡轮内部流动受动静干涉的影响，其流动本质上是周期性非定常的[1-3]，叶片排之间的相互作用对涡轮叶片通道内的流动产生不可忽视的影响，如尾迹、位势作用、二次流和激波[4-8]等之间的相互作用，以及对流场的非定常扰动。尤其在高负荷跨声涡轮中，在涡轮导叶出口存在很强的尾缘激波和反射波等各种波系，对流场中的非定常扰动更为明显。从20世纪末开始，国内外学者已经开始密切关注高负荷跨声涡轮中的非定常现象。Miller等[9]通过研究某跨声高压涡轮级下游的非定常流动，发现高压涡轮级和下游静叶之间具有复杂的3维相互作用；Davis等[10]对跨声速高压涡轮、过渡段和低压涡轮进行数值模拟，详细讨论了高压涡轮的动静干涉、激波及其反射波对流场的非定常效应，以及对叶片表面压力脉动的影响；邹正平等[11]使用数值模拟方法研究了单级轴流涡轮内部上游叶片排尾迹和位势作用对下游叶片表面负荷分布的影响，结果表明在涡轮级环境中，上游叶片排的尾迹对通道内部流动、边界层流动损失的产生、发展和输运以及下游叶片表面吸力面负荷分布均产生明显的影响。

为了减小计算量以及获得更加详尽的细节流动，在非定常模拟计算中往往采用区域缩放法[12]对计算模型进行简化，而模型中对各级叶片数的选取也仅从最大公约数的角度来考虑。Clack等[13]发现使用区域缩放法模拟出的2倍叶片通过频率的扰动量级与未模化时的结果明显不同，而后者与试验结果更吻合。Yao[14]也通过对使用区域缩放法的算例进行数值计算并与试验结果对比，证明了模化比例越接近于真实比例，计算结果越精确。因此，区域缩放法通过修改叶片几何形状导致非定常流动的预测精度受到污染，尤其是在各种扰动频率的预测方面尤为突出。然而，对于全环涡轮内部非定常流动的数值计算研究相对较少[15]。

在具备大型并行机群硬件支持的条件下，本文尝试利用CFD模拟真实情况下高负荷跨声涡轮全环非定常内部流场，深入分析涡轮内部叶片排之间的非定常扰动，为高性能涡轮设计方法以及流场中的非定常扰动的预测提供参考。

1 研究方法

以某型高压涡轮为研究对象，其几何和气动参数见表1。采用商用软件CFX12.0求解3维非定常黏性雷诺平均N-S方程进行数值模拟，运用时间追赶的有限体积法进行数值计算。空间、时间离散分别采用2阶迎风格式和2阶后差欧拉格式，并运用多重网格技术加速收敛。湍流模型选用2方程k-ε模型。

表1 高压涡轮几何和气动参数

计算整体和布局网格如图1所示。高压涡轮定常计算网格总数为2092万，包含高压导叶网格数873万（单通道网格数22.4万）、高压动叶1219万（单通道网格数20万）、导叶和动叶沿展向34个网格；高压动叶叶尖间隙为0.4mm，展向给了7个网格节点。

采用1条槽缝来代替1排气膜孔的方法来模拟冷气流动，其中1条槽缝的面积近似等于1排气膜孔面积的总和。冷气边界条件给定质量流量、冷气流动方向角、冷气总温、冷气湍流度和涡黏性比。在轴向间隙处，冷气流动方向与旋转轴正方向之间的夹角约为30°；进口湍流度为5%。

图1 计算整体和局部网格

在计算中，进口边界条件给定总温、总压和气流角。其中总温按照已知进口总温沿展向的分布，总压展向均匀分布，来流设为轴向,出口给定静压平均值。计算工质的设定根据油气比计算变比热模拟真实燃气，壁面边界条件设为固壁绝热和无滑移条件。导叶与动叶的交界面以及单通道与其余通道之间的交界面均采用直接差值的方法。

在具体计算时，设定动叶旋转1个非定常周期为800个物理时间步，动叶通过1个导叶通道约用20个时间步。每个虚拟时间步设为10。为方便判断非定常计算是否稳定，在动叶前缘和尾缘附近分别设置监测点（如图2所示），分别监测各点处流场的静压和轴向速度。在计算各状态点时，首先计算2400个时间步长，得到稳定的计算结果，然后再读取800个时间步长的瞬时结果并输出其统计平均结果。定常与非定常数值模拟计算在大型并行机群上进行，使用8×10个核心（内存160G），计算耗时约1个月。

图2 监测点位置（50%叶高）

2 结果分析

2.1 监测点变量随时间的变化

非定常计算结果稳定后2个监测点处静压和轴向速度随时间的变化如图3所示。横坐标表示时间，并以动叶旋转1周的时间T做无量纲处理。从图中可见，监测点的静压和轴向速度呈现出良好的周期性，说明非定常计算结果收敛较好。在动叶通过导叶通道时，监测点1处的静压和轴向速度基本呈现出周期性波动，在动叶旋转1周后共历时39个小周期性波动；监测点2处的静压和轴向速度同样也存在周期性，与监测点1处不同的是，其在1个周期内存在4个较小的周期性波动。下文将对此作详细分析。

图3 非定常计算监测点变量随时间的变化

2.2 导叶非定常流动

高压导叶叶根、叶中、叶尖表面等熵马赫数分布如图4所示。图中蓝线是时间平均值，黑、红线分别是1个周期内当地等熵马赫数的最大值和最小值。从图中可见，叶片表面的非定常负荷脉动，吸力面静压最大脉动位置在0.8～1.0轴向弦长范围内，而在其他位置叶片表面几乎没有压力脉动。说明导叶通道内的流场受非定常效应的影响很小，范围也仅限于导叶吸力面尾缘附近。从图4（d）中的马赫数云图（时均结果）中可见，在导叶尾缘存在较强的燕尾波，其内尾波分支撞击在相邻导叶吸力面并发生反射。内尾波的存在导致下游动叶非定常效应对上游的影响受限于吸力面下游附近。

图4 导叶表面等熵马赫数分布及50%叶高马赫数

导叶出口截面50%叶高处静压分布如图5所示。横、纵坐标分别表示导叶节距和时间，其中τ为动叶经过1个导叶通道的时间。

图5 导叶出口截面50%叶高处静压分布

从图中可见，在下游无动叶的情况下，静压分布不随时间变化，图中的等值线应该平行于竖直方向，而下游动叶排的周期性扫过导叶出口流场，使得在导叶出口激波位置发生周期性变化，同时也引起导叶吸力面尾缘附近的负荷脉动。在0.1节距位置处也出现静压周期性变化的现象，此处对应的位置为导叶尾缘激波的内尾波经相邻导叶吸力面后产生的反射波，其位置也随动叶的转动而呈周期性变化。导叶尾缘激波的非定常变化将直接影响动叶表面的负荷分布。

2.3 动叶非定常流动

在动叶中采用静压系数来描述叶片表面的负荷分布

式中：Pref、ρref分别为时均结果动叶进口截面平均值；U为动叶各叶高处的牵连速度。

动叶的叶根、叶中、叶尖表面静压系数分布如图6所示。从图中可见，静压的脉动幅度相对较大，叶片表面负荷分布形式的变化也很大，意味着动叶表面受到很强的非定常气动力作用。结合导叶出口处压力场的时空图可知，导致动叶表面很强的非定常效应主要有以下原因：（1）在上游导叶出口有很强的燕尾波，其外尾波分支周期性作用于动叶前缘，产生很强的非定常效应；（2）由于导叶负荷分布是后加载叶型，尾缘负荷很大，导致对下游的位势作用很强；（3）导叶出口尾迹的作用。

图6 动叶表面静压系数分布

2.4 上游尾迹和激波对下游流场的扰动

不同瞬时时刻50%叶高位置扰动速度矢量在通道中的分布如图7所示，其定义为瞬时与时均速度矢量之差。由尾迹形成的逆射流周期性地被下游叶片切断之后继续向下游传播，并对气流产生从压力面到吸力面的相对输运。在尾迹上、下游附近位置出现1对方向相反的涡，分别对尾迹下游和上游流场加、减速。对比t和t+0.25τ时刻，叶片最大厚度位置处在尾迹经过前后产生方向相反的扰动速度，会引起或者增加叶片表面压力脉动，进而引起叶片负荷分布的变化。随着尾迹的输运，扰动涡对影响的范围也逐渐扩大，但上、下游涡核心向下游输运的速度明显不同。在t+0.75τ时刻，尾迹下游的涡核心随着尾迹输运至动叶尾缘附近，而尾迹上游的涡核心则位于0.5倍的轴向位置处，意味着扰动涡对引起的叶片表面脉动与尾迹通过频率一致，但叶片表面脉动滞后于尾迹的输运。

图7 不同瞬时50%叶高位置扰动速度矢量及熵增分布

激波在叶片通道之间传播的相对位置如图8所示。激波A、C指在不同位置处的外尾波，B指内尾波及其反射波。

不同叶高动叶表面静压系数如图9所示。横、纵坐标分别为动叶经过1个导叶通道的周期数和动叶轴向弦长的无量纲数。从图中可见，在动叶表面已知叶高任一位置的静压随时间的变化过程。吸力面和压力面的分隔点位于动叶轴向位置的最左端和最右端。导叶尾缘激波外尾波分支首先作用于叶片吸力面顶点附近，并随时间变化逐渐向上游迁移直至动叶前缘，之后再作用于压力面0.6轴向距离。内尾波分支则经相邻导叶吸力面反射后也以类似方式作用于动叶吸力面。从激波相对动叶的运动轨迹可见，随着动叶的旋转，激波在吸力面是从下游向上游推进，反映在图中为从左上到右下方向的压力梯度对角线。吸力面在图中显示在0～0.4轴向位置由导叶尾缘激波引起较大的周期性变化的压力波动，尤其是在0.1轴向位置附近，压力随时间的变化最为剧烈。

图8 激波在叶片通道之间传播的相对位置

图9 不同叶高动叶表面静压系数(左为吸力面，右为压力面)

导叶尾迹通过动叶通道时，对其吸力面、压力面的影响均为从上游向下游发展，在图9中一致表现为从左下到右上方向的压力梯度对角线。从图9中可见，在压力面0.1～0.6轴向位置处的压力受尾迹作用影响，呈现周期性变化。在压力面0.7轴向位置的下游的表面压力受尾迹作用很小；在压力面0.1轴向位置的上游仍会受到外尾波和内尾波的反射波的影响，此处激波在压力面推进的方向与尾迹相同。

图8中C处外尾波经过动叶压力面后反射到相邻动叶吸力面上，其运动方式在图9中表现为从左上到右下方向的压力梯度对角线，而图9中显示的方向是与压力面相同方向的压力梯度对角线。表明在吸力面0.5轴向位置下游受C处外尾波影响较小，受导叶尾迹的影响较大，在图9中显示的方向与导叶尾迹的传播方向一致。

综上所述，在动叶吸力面0～0.4倍轴向长度区域，激波引起明显的高压力梯度，而在0.5～1.0倍轴向长度区域，尾迹对吸力面的压力影响较大，在10%叶高处可以较清晰地看到由内尾波的反射波引起的压力梯度，而在其他叶高位置处，此反射波影响不大。在压力面，除前缘附近受激波的影响外，叶片表面的压力脉动主要受上游尾迹作用影响，在0.7～1.0倍轴向长度区域则受尾迹影响很小。图中各叶高的时空图均表面，导叶尾缘激波引起的动叶表面非定常压力脉动幅度要大于导叶尾迹的，尾缘激波对下游动叶表面负荷分布的影响比导叶尾迹作用产生的影响更大。

从图9中还可见由导叶、动叶几何关系引起的不同叶高处激波相对于动叶的位置变化。对比10%、50%和90%叶高激波位置可见，在叶尖前缘位置最早受到导叶尾缘激波的影响，但激波强度相对不大，引起的压力脉动幅度较小。经过约0.15个动叶通过周期后，导叶叶中位置的激波移动到动叶前缘位置，此时激波引起的压力脉动幅度最大；最后受到激波影响的位置在动叶根部，激波强度与叶中相比有所减弱。

动叶进口距前缘约16%轴向长度位置沿各叶高处的静压分布如图10所示。横、纵坐标分别表示动叶节距的无量纲数和其经过上游1个导叶栅距的时间。在10%叶高位置，清楚可见在0.56轴向距离位置处，由激波和其反射波引起的压力脉动。在50%叶高位置膨胀波对进口静压分布的影响较为明显，膨胀波的位置是随动叶与导叶的相对位置变化而变化的，故在图中并没有呈现与激波作用明显相同（从左上到右下对角线）的趋势。A-A和B-B位置处的静压系数随时间的变化如图11所示，分别反映了激波、反射波和膨胀波对压力场分布的影响。

图11 不同位置静压系数随时间的变化

2.5 低频扰动

从如图3所示的监测点处各参数的时间变化可见，监测点1处的静压值和轴向速度由于动叶扫过导叶出口处的流场时，受到激波和尾迹的扰动而呈现出非常均匀地周期性变化。静压和轴向速度经过FFT变换后的频谱如图12所示。从监测点1的频谱分析可以看出主频（NPF）对应于转子通过1个导叶通道的频率，同时也是激波和尾迹作用于动叶通道内流场和叶片表面的频率。另外，叶尾缘激波和尾迹的位置也存在由动叶周期性扫过导叶出口流场引起的周期性变化。

图12 静压频谱分析

受上述诸多因素的影响，动叶通道内的静压和轴向速度除了存在上述扰动外，还存在低频扰动。从监测点2静压频谱分析中清楚可见这种低频信号，其频率约为主频的10%，即动叶完整旋转1周出现4次周期性的波动。

动叶表面50%叶高位置选取8个关键点处的静压（无量纲）随时间的变化如图13所示。其中点1位于叶片前缘，点2、3分别位于激波最远作用点的上、下游，点4、5分别位于动叶激波上、下游，点6～8分别位于压力面前缘、叶中和尾缘附近。各点处的静压波动幅值均有所不同，幅值变化最大和最小的位置分别位于点1、8位置处。对比发现从点2到3位置，叶片表面静压波动程度明显减小，而且在点3位置处清晰可见上述提到的低频扰动。另外，在压力面上点7位置处仍可以受到导叶尾缘激波的影响，其压力波动的幅度与点6处的相比变化不大。激波向下游传播时强度越来越弱，到点8位置处激波对叶片表面静压波动的影响消失，尾迹的影响也有所减弱，而低频扰动仍然存在。从动叶出口的静压波动情况来看，低频扰动对涡轮总体性能参数影响很小，而对局部细节的流场影响相对较大，如在点3位置处低频静压脉动量占动静干涉引起脉动量的8%。

图13 动叶表面各点的相对静压变化（50%叶高）

3 结论

（1）引起高压涡轮内部非定常流动的因素主要有激波、位势和尾迹作用等。计算结果显示由导叶激波引起的非定常效应最为明显，其作用在动叶叶片表面，引起强烈的压力波动。尾迹对动叶内部流场的影响也是不可忽略的，动叶压力面和吸力面0.5倍轴向弦长位置下游处的静压波动很大程度上是由于尾迹输运导致的。

（2）从监测点处静压变化可见，在动叶通道内静压波动最大区域对应于导叶激波通过时刻，波动相对较小的区域对应于内尾波的反射波通过时刻，而尾迹的影响则很小。

（3）从监测点的变量及动叶叶片表面各关键点的静压变化数据可知，除了由导叶激波和尾迹引起的动叶通道及叶片表面静压周期性脉动外，还存在1种低频扰动，其频率远小于动叶通过频率。从动叶出口的静压波动情况来看，低频信号对涡轮总体性能参数影响很小，而对局部流场影响相对较大。

（4）从全环非定常计算结果可知，采用区域缩放法按最大公约数简化模型（如在本算例中模化为2∶3），在预测动叶通道内流场中的各种非定常扰动时必然产生较大误差，模化后的流场抹消了流场中存在的低频扰动。

致谢

感谢北京航空航天大学能源与动力工程学院邹正平教授在本课题研究中给予指导和帮助！

[1] Sharma O P，Pickett G F，Ni R H. Assessment of unsteady flows in turbines [J]. Journal of Turbomachinery，1992，114（1）：79-90.

[2] Pollan G. Secondary flows and loss caused by blade row interaction in a turbine stage [J]. Journal of Turbomachinery,2006，128（3）：484-491.

[3]季路成，周盛.转/静干扰效应纳入设计体系的途径探索[J].工程热物理学报，1999，20（5）：580-583.

JI Lucheng, ZHOU Sheng. Exploration of the rotor-stator interaction applying to design system[J]. Journal of Engineering Thermophysics，1999，20（5）：580-583.（in Chinese）

[4] Casciaro C，Treiber M，Sell M. Unsteady transport mechanisms in an axial turbine [J]. Journal of Turbomachinery，2000，122（4）：604-612.

[5] LaGraff J E，Ashpis D E，Oldfield，et al. Unsteady flows in turbomachinery[R]. NASA-CP-2006-214484.

[6]綦蕾.涡轮端区非定常相互作用机理及流动控制技术探讨[D].北京：北京航空航天大学，2010.

QI Lei. Investigation of unsteady interaction mechanism and flow control in the turbine endwall regions[D]. Beijing: Beihang University，2010.（in Chinese）

[7] Kachel C E, Denton J D. Experimental and numerical investigation of the unsteady surface pressure in a three stage model of an axial high pressure turbine[J]. Journal of Turbo ma chinery，2004，128（2）：261-272.

[8] Giles M B. Stator rotor interaction in a transonic turbine [J].Journal of Propulsion and Power，1990，6（5）：621-627.

[9] Miller R J，Moss R W，Ainsworth R W，et al. The development of turbine exit flow in a swan-necked inter stage diffuser [R].ASME 2003-GT-38174.

[10] Davis R L，Yao J X，Clack J P，et al. Unsteady interaction between a transonic turbine stage and downstream components [J]. International Journal of Rotating Machinery，2004，10（6）：495-506.

[11]邹正平,叶建,张永新，等.非定常流动对叶片表面负荷分布影响的数值模拟研究[J].燃气涡轮试验与研究，2006，19（1）：21-26.

ZOU Zhengping, YE Jian, ZHANG Yongxin, et al. Numerical simulation of unsteady flow effects on turbine blade loading distributions [J]. Gas Turbine Experiment and Research，2006，19（1）：21-26 .（in Chinese）

[12] Rai M M. Multi airfoil Navier-Stokes simulations of turbine rotor stator interaction [J]. Journal of Turbomachinery，1990，112：377-384.

[13] Clack J P，Stetson G M. The effect of airfoil scaling on the predicted unsteady loading on the blade of a 1 and 1/2 stage transonic turbine and a comparison with experimental results [R]. ASME 2000-GT-0446.

[14] Yao J X, Davis R L. Massively parallel simulation of the unsteady flow in an axial turbine stage [J]. Journal of Propulsion and Power，2002，18（2）：465-471.

[15] Vahdati M，Sayma A，Imregun M. Prediction of high and low engine order forced responses for an LP turbine blade [R].AIAA-98-3719.