叶片顶端跨音速叶型的气动性能分析与优化

2021-10-23刘颖超杨其国

动力工程学报 2021年10期

刘颖超，陈榴，戴韧，杨其国

(上海理工大学能源与动力工程学院，上海 200093)

大功率凝汽式汽轮机排汽压力低，低压末级在高蒸汽比容和速度下工作。百万千瓦级汽轮机的末级全速钢制叶片长度已达1 300 mm，叶片顶部圆周速度约600 m/s。如果轴向绝对出口马赫数接近1，则动叶顶部截面叶型相对出口马赫数可达1.7左右[1]。末级动叶的气动设计是个难题，末级的低轮毂比和超长叶片导致末级从叶根到叶顶的流动速度、角度和叶型变化剧烈，尤其是叶片顶部截面叶型在满足跨音速流动气动设计时，还必须满足严苛的许可应力要求。

跨音速叶栅流动中不仅有激波和膨胀波，还有波系之间的相互作用，激波与边界层之间相互干涉，流动复杂而难以准确计算叶栅出流参数。早期计算方法是二元特征线法，其通过假设两相邻叶片喉部处音速线是一条直线来求解涡轮叶栅跨音速流动。Deich[2]通过分析得出特征线法在叶栅出口马赫数低于1.4、低节弦比和薄直背叶型的叶栅流动中相当准确。而对于高节弦比、高反动度的叶栅顶部薄叶型，Forster[3]的实验结果说明叶栅喉部的直线音速线假设已不成立，计算方法失效。

叶栅跨音速流动的实验测量也十分困难。平面叶栅实验中，叶片尾缘的激波在叶栅两侧自由流边界反射，形成“反射激波”，影响叶栅出口流动参数,在大安装角和叶片间距大的叶栅中很容易发生这种现象。另一方面，在超音速气流中，气动探针难以精确测量当地总压，采用纹影法可以获得激波形态，但是不能量化评定叶栅效率[4]。

基于时间推进法的CFD求解方法为分析跨音速叶栅流动提供了有力工具。分析激波相关流动，需要各种高精度格式，如总变差减小的差分格式(TVD)[5]、Roe迎风格式[6]、非自由参数耗散有限差分格式(NND)[7]等，目前对包含中等激波强度的跨音速叶栅流动计算基本没有困难，大涡模拟(LES)[8]计算也是可行的。具有2阶迎风格式的计算方法基本可以满足叶栅跨音速流动分析的需要，结合适当的优化算法，可以完成跨音速叶栅的优化设计[9-11]。

系列实验[12-13]证明，跨音速叶栅流动损失与叶型几何形状密切相关。叶型背弧的后半段曲率是重要的影响因素，大曲率出气边的叶型损失在超音速时是亚音速时的4倍，其中激波诱导边界层分离是损失增加的主要原因，这时采用直线背弧是更好的选择。前缘气流加速主要由叶片形状决定，并在某种程度上与压比无关，而尾缘对流场的影响主要取决于压比，尾缘形状的影响较小。各种经典叶顶截面叶栅跨音速流动[14]的流动性能大致相当，主要损失是叶栅尾缘激波和吸力面上的反射波损失。

笔者应用CFD软件分析某汽轮机末级动叶顶部截面叶型流动损失机理，结合叶栅气动设计方法，通过参数优化，调控叶栅流动中的激波强度和反射位置，降低叶栅内与激波有关的流动损失，提高叶栅流动效率。

1 叶栅流动分析的数值方法

1.1 验证模型与边界条件

流动分析选用CFX软件，采用有限体积法求解雷诺平均N-S方程，守恒方程中的对流项和湍流输运方程均选择2阶差分格式。选取类似末级动叶顶部截面的双圆弧叶栅实验[15]作为参照，验证数值方法对跨音速叶栅流动特性计算的可靠性。该叶型相对栅距为1.5，安装角为30°，最大相对厚度为0.1。

参考实验工况，计算域进口和出口分别给定总压p*和静压p2，压比p2/p*=0.23，进口气流角为35°，进口总温T*=280 K。周向采用周期性边界条件，叶片表面为绝热边界条件。

1.2 网格验证

计算网格采用结构网格，计算域进口距叶型前缘距离为C(其中C为弦长)，出口距叶型尾缘距离为2C，边界层伸长比小于1.1。如表1所示,采用3套不同质量网格进行流动计算，其中Cp为叶型表面压力系数,湍流模型选用SST模型，湍流度为5%，收敛后求解残差小于10-6。

表1 网格无关性比较Tab.1 Comparison of the grid independence

由表1可以看出，对出口马赫数、激波前Cp峰值及相应位置的预测，3套网格之间略有差别，中等网格与细网格的预测效果基本一样，表明中等网格的求解精度满足性能分析要求。

1.3 湍流模型验证

分别采用标准k-ε、k-ω和SST湍流模型进行定常计算，得到的叶型表面压力分布Cp如图1所示，其中p表示当地压力,X表示叶型的横坐标。与实验对比发现，3种湍流模型预测叶型压力面的压力分布与实验值均较吻合，差异不大，但对吸力面的预测效果存在偏差。标准k-ε模型对近壁区求解精度不高，在0.8C以后压力一直降低，无法准确预测激波。k-ω和SST模型预测所得的吸力面激波前压力高于实验值，激波位置与实验值吻合。在吸力面前缘附近，SST模型与实验值更加接近，因此笔者选择SST模型。

图1 采用不同湍流模型计算所得叶型表面压力分布Fig.1 Blade surface pressure distribution of different turbulent flow models

2 跨音速叶栅的流动分析

2.1 计算模型及边界条件

以某汽轮机末级动叶跨音速顶部截面叶型为原叶型，其叶栅结构如图2所示，叶片弦长为150 mm，相对栅距为0.926，最大厚度为6.09 mm，其他相关几何参数如表2所示。

图2 叶栅结构示意图Fig.2 Schematic diagram of cascade structure

表2 叶型几何参数Tab.2 Geometric parameters of the profile

计算域进口总压p*为17 021 Pa，出口静压p2为4 845.4 Pa，进口总温T*为308.1 K。不考虑末级流道在子午面的流线翘曲和子午加速率。湍流模型采用SST模型，湍流度为5%。计算网格仍采用上述结构网格。

2.2 叶栅流动的激波结构

图3为原叶型表面和沿叶栅流道中线MN(见图4中)的等熵马赫数分布。图4给出了叶栅流道内马赫数分布云图，由图4能清晰看到膨胀波与激波结构。在吸力面前缘，气流迅速达到超音速后再减速(AB1段)，产生了图4中AB位置的激波。然后,超音速气流在叶栅流道内流经相邻尾缘压力面附近，遇膨胀波K1加速，并产生反射膨胀波K2，加速至吸力面53%C处达到第2次速度峰值，即C点，等熵马赫数接近1.8。之后遇到尾缘压力面DC处斜激波减速(CD1段)，并在CE位置产生反射波(C2E1段减速)。由于出口马赫数高，在尾缘吸力面FG处也形成一道激波。

图3 叶栅内等熵马赫数分布Fig.3 Isentropic Mach number distribution in the cascade channel

图4 叶栅流道内马赫数分布Fig.4 Mach number distribution in the cascade channel

尾缘处吸力面和压力面激波通常称为“燕尾波”。各种波的存在、传播和相互作用，尾缘吸力面激波、反射波与尾迹(K3)干涉，以及激波诱导边界层分离，使得叶栅流动更加复杂。

2.3 叶栅气动损失分析

跨音速叶栅气动损失来源包括激波损失、边界层损失和尾缘损失(取决于尾缘厚度和基压)。衡量和评价不同损失对提高叶栅气动效率十分重要。

从热力学分析，流场中不可逆熵增是损失来源，可以用熵增来衡量叶栅内各种损失。跨音速叶栅中复杂激波造成的气动损失是熵增来源之一。利用激波前后的质量、动量和能量守恒，能计算激波引起的熵增Δsshock，表达式[16]如下：

(1)

式中：cv为比定容热容；k为绝热指数；Manormal为激波上游法向马赫数，是由图4中激波位置确定的平均变量。

以进口温度T0与进口速度u0对激波熵增无量纲化,得到无量纲激波熵增Δsshock-normalized：

(2)

表3统计了叶栅内不同位置的激波损失。由表3可知,尾缘处激波引起的损失最大，占比58.5%，其次是反射波，这是因为涡轮的高压比导致尾缘后基压低，当流体在尾缘周围膨胀到这种低压，在吸力面和压力面气流相遇的地方被强激波重新压缩[17]，使熵值增加。而反射波位于喉口后，是高马赫数集中区域，损失增大。

表3 叶栅不同位置激波损失Tab.3 Shock loss at different positions of the cascade

图5给出了叶型边界层厚度分布，在经过相邻尾缘压力面激波后，吸力面位移厚度突增形成波峰，表明激波干涉边界层，造成局部边界层分离，边界层增厚，气动损失增加。叶片边界层内黏性效应导致能量耗散，是不可逆过程，边界层损失用熵增Δsbl计算，定义[16]如下：

图5 叶型边界层厚度分布Fig.5 Boundary layer thickness distribution of the blade

(3)

由进口温度与进口速度对边界层引起的无量纲熵增为：

(4)

通过上式计算，边界层无量纲熵增损失为0.157 8，其中还有压力侧激波与吸力面边界层干涉而附加的损失，而无量纲激波熵增损失为0.255 4，能量损失更大，故削弱激波强度、减小激波相关损失对提升叶栅气动性能是首要的。通过调整优化叶型几何形状，改变曲率分布，控制尾缘激波及其反射波，以降低叶栅内激波损失。

3 叶型气动优化方法

3.1 叶型的参数化

采用中弧线叠加厚度分布生成叶型[18]。中弧线采用3阶Bezier曲线表示，如图6所示。n阶Bezier曲线表达式如下：

图6 采用3阶Bezier曲线表示的叶型中弧线Fig.6 Camber line of the blade represented by the third-order Bezier curve

(5)

式中：p(t)表示构成Bezier曲线上的点;Bi,n(t)为Bernstein函数；Pi为直角坐标系下控制点的位置矢量。

叶型中弧线首、末两点由叶型前缘点、尾缘点和安装角确定。前后缘的斜率设定为进、出气流角的约束，控制点P2、P3只能分别沿直线P0P2和P3P4移动，故中弧线形状由2个参数XP2、XP3控制。

叶型厚度分布如图7所示，以最大厚度控制点T3为界，左侧用2阶Bezier曲线表示，右侧用3阶Bezier曲线表示。通过叶型前缘厚度和尾缘厚度确定首、末2点T1、T6。控制点T2、T3、T4的厚度均为最大厚度，叶型厚度分布由控制点T3、T4和T5确定，故叶型厚度分布由5个参数XT3、YT3、XT4、XT5、YT5控制。

图7 Bezier曲线表示叶型厚度分布Fig.7 Thickness distribution of the blade represented by Bezier curve

3.2 基于Kriging代理模型的优化

采用基于克里金法(Kriging)的代理模型，对动叶顶部截面叶型进行气动优化。Kriging代理模型[19]具有精确的全局近似设计空间，能准确地预测线性和非线性函数，其表达式如下：

(6)

试验设计(DOE)是构建Kriging代理模型的第一步，建立设计空间生成样本。采用最优拉丁超立方设计在设计变量范围(见表4)内生成30个样本，作为Kriging代理模型的训练样本集。

表4 设计变量取值范围Tab.4 Value range of design variables

优化的目标函数选择能评价叶栅气动损失的出口总压损失系数Cp,t，其具体定义为：

(7)

式中：p1为进口静压；p02为出口总压。

叶型前缘和尾缘直径分别为2.5 mm和1.4 mm。设定的约束条件是叶型截面面积的变化相对于原叶型不超过±10%。在Isight软件上，基于DOE生成样本，选择Kriging代理模型建立目标函数的近似模型，并在设计变量范围内，采用多岛遗传算法寻优，得到最低Cp,t所对应的设计变量数值,然后用CFX校核最优叶型的性能。

3.3 优化结果分析

经过优化计算，最终获得优化后叶型，图8给出了优化前后叶型的几何构型对比,图中采用叶型弦长对叶型轴向和周向坐标进行无量纲化。优化后叶型截面面积比原叶型小9.6%，最大厚度为5.78 mm，位置基本不变。在叶型80%C～90%C处，优化叶型压力面曲率明显减小，型线为直线，有利于减小压力梯度，削弱尾缘压力面激波强度。通过Kriging模型预测得到Cp,t为0.215 6，CFX实际计算所得Cp,t为0.216 7，代理模型不确定度为0.51%。

图8 优化前后叶型对比Fig.8 Comparison between optimized and baseline profiles

叶片表面连续曲率影响叶片性能[20-21]。图9比较了优化前后叶型吸力面曲率的分布，优化后叶型吸力面曲率整体更加连续平滑，在20%C～75%C处曲率均低于原叶型，75%C后优化后叶型曲率略大于原叶型，从而决定叶型表面等熵马赫数的分布形状(见图10)。由图10可知，优化后叶型前缘吸力面的速度减小，吸力面上的第2次速度峰值降低并向后移动，延迟了压力峰值点，尾缘压力面激波入射位置后移。从图11优化前后沿叶栅流道中线MN的等熵马赫数分布也可以看到，反射波位置后移，激波强度削弱。同时，优化后叶型接近尾缘的压力面维持了较高压力，起到类似后加载的作用，有利于降低叶型气动损失。图12给出了优化后叶栅流道内马赫数分布。由图12可以看出,优化后叶栅喉口后高马赫数集中区域缩小，吸力面反射波区域缩小，强度减弱。

图9 优化前后叶型吸力面曲率分布对比Fig.9 Comparison of curvature distribution of blade suction surface before and after optimization

图10 优化前后叶型表面等熵马赫数对比Fig.10 Comparison of isentropic Mach number distribution of the blade surface before and after optimization

图11 优化前后沿叶栅流道中线MN的等熵马赫数分布对比Fig.11 Comparison of isentropic Mach number distribution along the centerline MN of the cascade channel before and after optimization

图12 优化后叶栅流道内马赫数分布Fig.12 Contour distribution of Mach number distribution in cascade channel of optimized blade

图13比较了优化前后叶型不同位置的激波损失，优化后叶型的各处与激波相关的损失均下降，前缘AB处激波基本消失，反射波损失减少26%，尾缘处激波损失减少34%，总的激波损失相比原叶型下降约32%。同时，激波强度减弱能减少对边界层的干扰，边界层损失降低。优化后叶型的边界层损失Δsbl-normalized为0.130 8，相比原叶型降低了17%。

图13 优化前后叶型不同位置激波损失对比Fig.13 Comparison of shock loss at different positions of the blade profile before and after optimization

图14和图15给出了优化前后叶栅出口1.2倍轴向弦长Cax处尾迹速度和尾迹气流角沿相对栅距的分布,图中采用叶栅出口截面平均速度对尾迹速度进行无量纲化。由图14和图15可以看出，优化叶型低速幅值降低,尾迹宽度有所减小，尾迹气流角幅值有所降低，有利于减小尾迹掺混损失。优化后，吸力面侧尾迹气流角的变化更加缓慢，压力面侧尾迹速度降低，尾迹气流角的变化相对平缓，表明尾缘燕尾波减弱，尾流干涉减弱。