2. 德州农工大学 机械工程系, College Station TX 77843
2. Department of Mechanical Engineering, Texas A & M University, College Station, TX 77843, United States
随着高性能计算的快速发展，基于计算流体力学(CFD)的优化设计技术在叶轮机械叶片气动设计中已经得到了应用。基于CFD的气动设计优化(ADO)不仅能够缩短研制周期、降低成本，相对于依靠设计者经验的设计方法还能提高设计结果精度。对于压气机而言，通过优化叶片气动外形增加压比、降低流动损失对于提高压气机气动性能具有重要意义。
目前得到应用的优化方法分为两类：梯度方法和非梯度方法。进化算法[1-3]和代理模型[4]是典型的非梯度方法。Benini[2]采用进化算法对跨声速压气机转子叶片NASA Rotor 37的中弧线和厚度分布进行了优化设计，Lian和Liou[3]采用进化算法和代理模型对跨声速压气机转子叶片NASA Rotor 67进行了多目标优化设计，Samad等[4]采用多层次代理模型对转子叶片的叶片掠、叶片倾斜及叶片弯进行了优化设计。上述非梯度方法虽然在叶轮机械ADO中已经得到应用，但是对于多设计参数的复杂ADO问题，尤其是需要通过求解Navier-Stokes方程计算流场时，将需要大量的计算资源，制约了这类方法的工程应用。
基于灵敏度分析的优化方法在叶轮机械ADO中得到了应用，直接差分法(FDM)和伴随方法是两种常用的灵敏度计算方法。采用FDM计算灵敏度时，所需的流场计算次数和设计参数线性相关。在20世纪80年代，Jameson[5]成功地将伴随方法应用于飞行器气动外形优化设计。采用伴随方法，每个气动函数的灵敏度计算只需要大约两次流场计算，计算时间和设计参数个数基本无关。之后，基于伴随方法的翼型、机翼以及翼身组合体ADO研究得到了迅速发展[5-9]。由于灵敏度计算的高效性和精确性，伴随方法在叶轮机械ADO中的应用也越来越广泛[10-19]。Luo等[17-18]开展了基于伴随方法的跨声速压气机叶片的多工况、多目标ADO研究。近年来，掺混面模型方法也应用于伴随方程的数值求解，推动了基于伴随方法的叶轮机械多排ADO研究的发展[13, 14, 16, 20]。
本文将发展一种基于掺混面模型的叶轮机械多排伴随方法，并应用于某型低速、低压比压气机ADO。首先采用自主开发的初步设计方法设计4.5级压气机叶片的气动外形，通过数值模拟发现：最后级静子叶片的端壁和吸力面角区内存在较大分离，严重影响了压气机效率。之后介绍适用于多排ADO的伴随方法，并在近失速工况下对最后级静子叶片的中弧线及安装角进行优化，优化目标是降低出口熵增，同时约束出口流量。最后开展多工况ADO，在两个不同转速下进行最后级静子叶片的外形优化设计，分析外形优化设计对设计工况、非设计工况气动性能的影响。
1 压气机初步设计该低速4.5级压气机由第1排进口导叶和后4级叶排构成，给定设计压比和等熵效率。在初步设计之前需要确定一些前提条件：单级多变效率和压气机相等，给定单级温比，压气机进口到导叶出口的面积收缩率等。采用基于经验修正[21]的4.5级压气机初步设计的基本流程为：
第1步 给定效率初值，确定压气机进、出口流道几何参数、速度三角形和气动参数。
第2步 根据给定的温升分布、静子出口气流角和流道参数等，进行逐级设计，迭代计算各级转子、静子出口截面的相关参数。
第3步 计算整机压比和效率，判断是否满足设计条件，如果不满足则用计算效率替代初始效率，重复设计。
第4步 采用简单径向平衡方程确定轮毂和机匣的相关参数。
第5步 确定压气机轴向尺寸，进行二维叶片造型和径向堆积。
图 1为该4.5级压气机轮毂、叶中和机匣处的叶型，图中：x表示轴向，Rθ表示周向。值得注意的是：初步设计时，导叶和静子叶片均采用直叶片设计。流场分析发现：直叶片设计使得在前三级的影响下，最后级静子叶片的端壁和吸力面角区流动分离严重；因而，该初步设计的4.5级压气机存在较大优化空间，可以通过叶片外形优化设计抑制流动分离，提高压气机气动性能。
|
| 图 1 3个不同径向位置4.5级压气机叶片叶型 Fig. 1 Blade profiles of the 4.5-stage axial compressor at three different spans |
ADO的流场计算采用流场求解器Turbo90[22-24]，给定目标函数及优化约束的数学建模以确定伴随方程的边界条件并实现数值求解。确定伴随灵敏度后，采用简单、易实现的最速下降法来进行优化设计，气动外形扰动步长根据经验给定。
2.1 流场数值模拟压气机最后级的流场通过求解雷诺平均Navier-Stokes (RANS)方程及Spalart-Allmaras (SA)湍流模型[25]确定。采用多重网格、隐式残值光顺、当地时间步长等技术加速计算收敛。最后级进口流动参数如总压、总温、流动角由4.5级压气机流场数值模拟结果给定，出口处给定轮毂压力并由径向平衡方程及其他流动参数确定压力径向分布，转-静交接面上采用掺混面模型[26]实现流场连续求解，并采用无反射边界条件[27]。本文所采用的Mixed-out掺混面模型能够保证交接面上流量、动量和能量的通量守恒，其详细数学形式参考附录A。
计算采用H-型结构网格，不考虑叶尖间隙影响，S1及S2流面网格如图 2所示。首先在3套不同的网格上进行流向和周向网格无关性研究，径向网格数固定为73，由此确定收敛的流向网格数；之后在3套不同网格上进行周向和径向网格无关性研究，叶面流向网格数固定，周向和径向网格增长比为2。网格无关性研究的具体实现流程可参考本文作者前期工作[19]。表 1为近失速工况下周向和径向网格无关性研究结果，包括：归一化流量、静子效率、级效率、压比。由表可知：相对于网格2与网格1的流场偏差，网格3与网格2的流场偏差明显减小，流量、级效率的相对偏差分别为0.31%和0.46%，压比偏差接近于零。后续研究将采用网格2，转子和静子的网格规模均为：叶面流向97、周向49、径向57个网格点。
|
| 图 2 最后级网格 Fig. 2 Grids for the last stage |
| 网格 | 流量 | 静子效率/% | 级效率/% | 级压比 |
| 1 | 0.926 1 | 80.659 | 89.921 | 1.068 4 |
| 2 | 1.000 0 | 91.844 | 94.586 | 1.070 6 |
| 3 | 1.003 1 | 92.267 | 94.987 | 1.070 8 |
实际上，为了进一步证实网格无关性解的有效性，研究中对优化后的最后级也采用了相同规模的密网格(网格3)进行流场数值模拟，其气动性能也得到了提升。
2.2 参数化建模中弧线、厚度分布、安装角是叶片的3个重要设计参数。本文的ADO研究中，只改变最后级静子叶片的中弧线和安装角径向分布，其厚度分布不变。中弧线的变化通过叠加Hicks-Henne型函数[28]实现，安装角变化在前期工作中已经介绍[15]。选取10%、15%、20%、30%、70%、80%、85%、90%这8个径向位置，每个径向位置的安装角都作为设计参数，叶型中弧线变化由8个型函数叠加确定；共有72个设计参数。
图 3为叶型中弧线变化和安装角变化的几何示意图，其中：δ表示变化量，下标camber和stagger分别表示中弧线和安装角。
|
| 图 3 设计参数变化 Fig. 3 Design parameter variations |
叶型变化后，在原始网格基础上采用插值网格扰动技术[7]重新生成计算网格：
| $ {x_{{\rm{ new }}}} = {x_{{\rm{ old }}}} + r({x_{{\rm{W, new }}}} - {x_{{\rm{W ,old }}}}) $ | （1） |
式中：x表示网格坐标；下标old和new分别表示原始和更新后的网格, W表示叶面网格点；r为归一化的空间距离，在叶面网格点r=1。
2.3 目标函数本文ADO的主要目标是降低流动损失，同时约束出口流量。叶轮机械的流动损失可以用熵增来衡量，因而目标函数定义为
| $ \mathit{\boldsymbol{I}} = {S_{{\rm{gen}}}} + \Lambda |\sigma - 1| $ | （2） |
式中：σ表示设计流量与原始流量的比值；Λ为约束系数；Sgen表示平均熵增，定义为
| $ {S_{{\rm{gen}}}} = \frac{{\int {\int \Delta } s\rho \mathit{\boldsymbol{u}} \cdot {\rm{d}}\mathit{\boldsymbol{s}}}}{{\int {\int \rho } \mathit{\boldsymbol{u}} \cdot {\rm{d}}\mathit{\boldsymbol{s}}}} $ | （3） |
其中：Δs表示每个单元上的熵增，其定义见文献[15]。
基于伴随方法的多工况ADO方法在前期工作中已经介绍[17]，定义加权目标函数为
| $ \mathit{\boldsymbol{I}} = \sum\limits_{i = 1}^M {{\lambda _i}} {\mathit{\boldsymbol{I}}_i} $ | （4） |
式中：I表示气动参数；M为工况数；λi表示权重且满∑λi=1。相应地，灵敏度计算表达式为
| $ \mathit{\boldsymbol{G}} = \sum\limits_{i = 1}^M {{\lambda _i}} {\mathit{\boldsymbol{G}}_i} $ | （5） |
式中：Gi表示不同工况下的气动参数灵敏度。
本文研究采用最速下降法作为优化方法，所需要的灵敏度由伴随方法计算。接下来将重点介绍伴随方法的基本原理。
3 伴随方法伴随方法基本原理较为成熟，这里只作简单介绍。定义I为气动参数，一般与流场w和气动外形相关，而气动外形与设计参数V相关，即
| $ \mathit{\boldsymbol{I}} = \mathit{\boldsymbol{I}}(\mathit{\boldsymbol{w}},\mathit{\boldsymbol{V}}) $ | （6） |
气动参数关于设计参数的灵敏度为
| $ \frac{{\delta \mathit{\boldsymbol{I}}}}{{\delta \mathit{\boldsymbol{V}}}} = \frac{{\partial \mathit{\boldsymbol{I}}}}{{\partial \mathit{\boldsymbol{w}}}} \cdot \frac{{\delta \mathit{\boldsymbol{w}}}}{{\delta \mathit{\boldsymbol{V}}}} + \frac{{\partial \mathit{\boldsymbol{I}}}}{{\partial \mathit{\boldsymbol{V}}}} $ | （7） |
式中：右端第1、第2项分别表示流场变分、气动外形变分对灵敏度的影响，灵敏度计算的关键是流场变分的计算。如果采用FDM计算灵敏度，虽然原理简单，但是计算效率偏低。值得注意的是，流场和设计参数满足流体力学方程R，即
| $ \mathit{\boldsymbol{R}}(\mathit{\boldsymbol{w}},\mathit{\boldsymbol{V}}) = {\bf{0}} $ | （8） |
类似于方程(7)，方程(8)可以写成
| $ \delta \mathit{\boldsymbol{R}} = \frac{{\partial \mathit{\boldsymbol{R}}}}{{\partial \mathit{\boldsymbol{w}}}}\delta \mathit{\boldsymbol{w}} + \frac{{\partial \mathit{\boldsymbol{R}}}}{{\partial \mathit{\boldsymbol{V}}}}\delta \mathit{\boldsymbol{V}} = {\bf{0}} $ | （9） |
伴随方法的主要目标是通过将流体力学方程作为约束引入到灵敏度计算中，消除流场变分对灵敏度计算的影响。定义伴随变量Ψ=[Ψ1, Ψ2, Ψ3, Ψ4, Ψ5]T，将伴随变量与方程(9)相乘后再与方程(7)相减，有
| $ \frac{{\delta \mathit{\boldsymbol{I}}}}{{\delta \mathit{\boldsymbol{V}}}} = (\frac{{\partial \mathit{\boldsymbol{I}}}}{{\partial \mathit{\boldsymbol{w}}}} - {\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}^{\rm{T}}}\frac{{\partial \mathit{\boldsymbol{R}}}}{{\partial \mathit{\boldsymbol{w}}}})\frac{{\delta \mathit{\boldsymbol{w}}}}{{\delta \mathit{\boldsymbol{V}}}} + (\frac{{\partial \mathit{\boldsymbol{I}}}}{{\partial \mathit{\boldsymbol{V}}}} - {\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}^{\rm{T}}}\frac{{\partial \mathit{\boldsymbol{R}}}}{{\partial \mathit{\boldsymbol{V}}}}) $ | （10） |
若方程(10)右端第1项的系数为零，即
| $ \frac{{\partial \mathit{\boldsymbol{I}}}}{{\partial \mathit{\boldsymbol{w}}}} - {\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}^{\rm{T}}}\frac{{\partial \mathit{\boldsymbol{R}}}}{{\partial \mathit{\boldsymbol{w}}}} = 0 $ | （11） |
则有
| $ \mathit{\boldsymbol{G}} = \frac{{\delta \mathit{\boldsymbol{I}}}}{{\delta \mathit{\boldsymbol{V}}}} = \frac{{\delta \mathit{\boldsymbol{I}}}}{{\delta \mathit{\boldsymbol{V}}}} - {\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}^{\rm{T}}}\frac{{\partial \mathit{\boldsymbol{R}}}}{{\partial \mathit{\boldsymbol{V}}}} $ | （12） |
方程(11)称为方程(8)的伴随方程，Ψ为伴随变量，方程(12)为灵敏度计算公式且与流场变分无关。由此可知：只需求解一次流体力学方程和一次伴随方程，即可由方程(12)计算气动参数I的灵敏度。此外，伴随方程是线性方程组，可以采用和RANS方程相同的数值格式进行求解[11, 29]。
Giles和Pierce[29]指出了伴随方程进、出口边界条件确定方法：伴随方程的特征线方向与流体力学方程相反。对于多排压气机叶片，转子出口、静子入口的伴随方程边界条件需要采用掺混面方法[13, 16]确定。本文采用基于Mixed-out掺混面模型的伴随掺混面方法，首先计算伴随通量Fa1、Fa2、Fa3、Fa4和Fa5：
| $ \left\{ {\begin{array}{*{20}{l}} {{F_{{\rm{a1}}}} = \int {\int {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_1}} } (\rho \mathit{\boldsymbol{u}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}}\\ {{F_{{\rm{a2}}}} = \int {\int {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_2}} } (\rho {u_x}\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}}\\ {{F_{{\rm{a3}}}} = \int {\int {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_3}} } (\rho {u_\theta }\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}}\\ {{F_{{\rm{a4}}}} = \int {\int {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_4}} } (\rho {u_r}\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}}\\ {{F_{{\rm{a5}}}} = \int {\int {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_5}} } (\rho H\mathit{\boldsymbol{u}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}} \end{array}} \right. $ | （13） |
式中：r、p、H分别表示密度、压力和总焓；ds为面积单元；u =[ux, uθ, ur]表示速度矢量；n为单位法向矢量。之后由方程(13)可以计算平均伴随变量Ψa，即
| $ \left\{ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_{{\rm{a1}}}} = {F_{{\rm{a1}}}}/\int {\int {\rho \mathit{\boldsymbol{u}} \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_{{\rm{a2}}}} = {F_{{\rm{a2}}}}/\int {\int {(\rho {u_x}\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_{{\rm{a3}}}} = {F_{{\rm{a3}}}}/\int {\int {(\rho {u_\theta }\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_{{\rm{a4}}}} = {F_{{\rm{a4}}}}/\int {\int {(\rho {u_r}\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{\mathit{\boldsymbol{ \boldsymbol{\varPsi} }}_{{\rm{a5}}}} = {F_{{\rm{a5}}}}/\int {\int {\rho H\mathit{\boldsymbol{u}} \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } } \end{array}} \right. $ | （14） |
转-静交接面上仍采用Giles和Pierce提出的无反射边界条件。具体实现流程为[10, 13]：在亚声速交接面下游(对应下游静子入口)，伴随特征变量Ψa1、Ψa2、Ψa3、Ψa4由内场确定，Ψa5由交接面上游确定；在亚声速交接面上游(对应上游转子出口)，Ψa1、Ψa2、Ψa3、Ψa4由交接面下游确定，Ψa5由内场外插确定。确定所有伴随特征变量后，伴随变量Ψa1、Ψa2、Ψa3、Ψa4、Ψa5可通过求解线性方程组确定。
4 结果分析在基于伴随方法的ADO研究中，首先需要对伴随灵敏度进行精度验证，将基于FDM的灵敏度收敛解作为精确解。图 4给出了30%叶高的灵敏度分布，前8个为中弧线设计参数灵敏度，第9个为安装角变化灵敏度；AD表示伴随方法。
|
| 图 4 灵敏度对比 Fig. 4 Comparisons of gradients |
整体上，从叶片前缘(LE)至叶片中部，伴随灵敏度与FDM灵敏度较为接近；从叶片中部至叶片尾缘(TE)，伴随灵敏度与FDM存在一定偏差，主要原因是没有考虑湍流模型方程的伴随方程求解，而是采用“定涡”假设(Frozen Eddy Viscosity)，即忽略湍流变量变分对灵敏度的影响。实际上，当流动分离较大时，忽略湍流影响会降低伴随灵敏度精度[30]。图 4中，后3个中弧线设计参数(6, 7, 8)正处于静子叶片分离区。虽然精度有所下降，但是伴随灵敏度的分布与FDM基本一致，本文将采用基于“定涡”假设的伴随方法进行ADO研究。
4.1 近失速点单点优化该4.5级压气机最后级静子角区存在流动分离，且在接近失速工况时分离区愈大。单点优化设计在设计转速近失速工况下进行，约束函数系数Λ= 0.4。图 5给出了等熵效率和流量的收敛曲线，流量最大相对偏差不超过0.1%，约束效果明显。表 2给出了优化前后的熵增、流量及静子和最后级的等熵效率，优化后静子和级效率分别增加了约0.29%、0.15%。
|
| 图 5 单点优化熵增和流量收敛曲线 Fig. 5 Convergence history of entropy and mass flow rate under one point optimization |
| 叶片 | 熵增 | 流量 | 静子效率/% | 级效率/% |
| 优化前 | 1.000 0 | 1.000 0 | 91.844 | 94.586 |
| 优化后 | 0.990 9 | 1.001 1 | 92.133 | 94.731 |
研究中采用相同规模的密网格(表 1中的网格3)对优化叶片流场进行数值模拟，表 3给出了优化前后密网格的计算结果。与原始叶片相比，优化后静子和级效率分别增加了0.28%、0.13%。这种较好的网格一致性结果进一步证实了本文网格无关性验证的有效性及优化设计结果的可靠性。
| 叶片 | 熵增 | 流量 | 静子效率/% | 级效率/% |
| 优化前 | 1.000 0 | 1.003 1 | 92.267 | 94.987 |
| 优化后 | 0.991 2 | 1.003 8 | 92.544 | 95.115 |
图 3为优化前后静子叶片的气动外形。由图 6(a)可知：优化设计后自前缘至20%弦长，叶片中弧线曲率增加；自20%弦长至90%弦长，叶片中弧线曲率降低；叶片气动外形的这种变化将使前缘附近流动加速、来流角与叶片前缘角偏差降低，有利于抑制叶背流动分离。
|
| 图 6 优化前后静子叶片气动外形 Fig. 6 Stator blade aerodynamic shape before and after optimization |
为了比较分析优化设计对最后级流场的影响，图 7给出了静子吸力面附近的流线，图 8给出了最后级10%叶高S1流面的马赫数云图。由图 7可知：原始静子轮毂、机匣角区均存在流动分离，轮毂角区尤为明显；优化设计使轮毂流动分离区朝尾缘和端壁移动，流动分离得到了较好抑制；机匣流动分离区也得到了较好抑制。由图 8可以更好地看出优化设计对静子叶片叶背流动分离的抑制作用。整体上，对于低速压气机，流动分离区的减小将是效率增加的主要原因。
|
| 图 7 静子吸力面附近流线 Fig. 7 Streamlines near stator suction side |
|
| 图 8 10%叶高S1流面马赫数云图 Fig. 8 Mach number contours at 10% span S1 stream surface |
上述优化设计能够改善最后级在指定工况下的气动性能，在压气机全工况范围内对气动性能的影响更值得关注。图 9给出了最后级、静子叶片的等熵效率η特征工作线。虽然在其他工况条件下气动性能改善并不明显，但是静子、级效率在全工况范围内都得到一定的提升。若对最后级的转子叶片也进行优化设计，甚至对整个4.5级压气机进行优化设计，通过调整上游叶排对下游流动的影响，将能进一步提升压气机气动性能。
|
| 图 9 优化前后效率工作线 Fig. 9 Adiabatic efficiency for both original and optimized blades |
值得注意的是：单点优化设计的工况点远离失速边界，叶背流动分离区较小，通过控制流动分离获得效率提升有限。实际上，在前期工作中采用基于本征正交分解的混合模型方法对该型号压气机最后级也进行了优化设计[31]：工况点距离失速边界较近，流动分离较剧烈，通过最后级静子安装角优化使等熵效率提高了1.47%。因此，虽然通过本文的优化设计获得的效率提升较为有限，但是由上述分析可知：优化设计结果是可靠的，在工程中仍具有应用价值。
4.2 多工况优化多工况优化设计研究中，将对最后级静子叶片在100%设计转速的近失速工况(P1)、140%设计转速的近失速工况(P2)进行中弧线和安装角的优化设计。在每个转速下分别计算其灵敏度，最后由方程(5)计算加权灵敏度，并由最速下降法实现优化设计, 权重λ1=λ2=0.5。
图 10给出了两个不同转速下(P1和P2)熵增和出口流量的收敛曲线。表 4给出了优化前后叶片气动性能。值得注意的是：表 4中P2转速下原始和优化叶片的熵增和流量都统一采用P1转速下原始叶片的熵增和流量进行归一化处理。多工况优化设计使静子效率、级效率在P1转速下分别提升了0.29%和0.15%，在P2转速下分别提升了0.24%和0.13%，而两个转速下流量偏差都低于0.15%，表明两个转速条件下流量都得了较好约束。此外，由表 4可知高转速条件下熵增增加，其主要原因是：高转速条件下压气机负荷增大，角区分离加剧，流动损失增大。
|
| 图 10 多工况优化熵增和流量收敛曲线 Fig. 10 Convergence history of entropy and mass flow rate under multi-points optimization |
| 转速 | 熵增 | 流量 | 静子效率/% | 级效率/% |
| 原始P1 | 1.000 0 | 1.000 0 | 91.844 | 94.586 |
| 原始P2 | 1.067 3 | 1.463 9 | 94.009 | 96.109 |
| 优化P1 | 0.992 6 | 1.001 2 | 92.139 | 94.733 |
| 优化P2 | 1.053 5 | 1.465 5 | 94.249 | 96.238 |
图 11给出了P2转速下最后级10%叶高S1流面的马赫数云图。值得注意的是：图 11的最大相对马赫数为0.45，远高于图 8的0.375。对比图 11和图 8发现：随着转速的提高，原始叶片叶背分离更严重，如图 8(a)和图 11(a)所示；高转速下，通过优化设计静子吸力面的分离区也得到了较好抑制；此外，相对于图 8，高转速下分离区的控制效果更加明显。
|
| 图 11 140%设计转速下10%叶高S1流面马赫数云图 Fig. 11 Mach number contours at 10% span S1 stream surface with 140% design rotation speed |
上述结果表明：通过单点、多工况气动外形优化设计能够提升近失速工况最后级的等熵效率。非设计转速下气动外形设计优化对设计转速下气动性能的影响需要进一步深入分析。图 12给出了单点、多工况优化后静子、最后级在非设计转速全工况范围内的效率变化。整体上，多工况优化设计能使静子、最后级的等熵效率在非设计转速全工况范围内得到提升。主要原因是：静子端壁和吸力面角区的流动分离在不同转速条件下都存在，通过气动外形设计优化抑制流动分离的基本原理是一致的：减小来流角和叶片前缘角的偏差、流动加速等。此外，由于ADO在近失速工况下进行，当工作点远离失速边界时，效率提升不断降低；当工作点向阻塞边界靠近时，效率提升不断增加，这一点有待深入研究。
|
| 图 12 单工况、多工况优化效率提升 Fig. 12 Improvements of adiabatic efficiency for single and multi-points optimization |
1) 采用伴随方法计算梯度优化方法所必需的灵敏度信息，即使对于多设计参数的ADO问题，优化效率仍非常高；对于强湍流问题，忽略湍流影响的伴随灵敏度精度会稍微下降，但仍可以用于气动外形设计优化工作。
2) 通过单/多工况优化设计，最后级静子叶片中弧线曲率在前缘附近增加，降低了相对来流角并且使叶片吸力面流动加速，对叶片端壁和吸力面角区的流动分离具有较好的抑制效果。
3) 最后级及静子叶片效率的提升主要来源于流动分离的抑制，在全工况范围内叶片等熵效率都得到了提升；此外，多工况设计优化能够在不同转速条件下实现叶片全工况范围内气动性能的提升。
本文的研究工作初步实现了基于黏性伴随方法的叶轮机械多叶排ADO，对于后续开展考虑湍流影响的伴随方法研究、基于伴随方法的多级ADO研究提供了理论和技术基础。
附录AMixed-out是叶轮机械多排流场数值模拟的一种主要的掺混面方法, 该方法的基本表达式为
| $ \left\{ {\begin{array}{*{20}{l}} {{F_{\rm{1}}} = \int {\int {\rho \mathit{\boldsymbol{u}} \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{F_{\rm{2}}} = \int {\int {(\rho {u_1}\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{F_{\rm{3}}} = \int {\int {(\rho {u_2}\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{F_{\rm{4}}} = \int {\int {(\rho {u_3}\mathit{\boldsymbol{u}} + p\mathit{\boldsymbol{n}}) \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } }\\ {{F_{\rm{5}}} = \int {\int {\rho H\mathit{\boldsymbol{u}} \cdot {\rm{d}}\mathit{\boldsymbol{s}}} } } \end{array}} \right. $ | （A1） |
式中：H表示总焓；p、ρ分别表示压力、密度；u、n、ds分别表示速度矢量、单位外法线矢量、单元面积矢量。
方程(A1)经过整理后，可以简化成一个关于平均压力的二次方程。确定平均压力之后，再由方程(A1)即可计算其他平均流动参数：
| $ \left\{ \begin{array}{l} \bar p = \frac{{a/b \pm \sqrt {{{(a/b)}^2} + ({\gamma ^2} - 1)(c - 2{F_1}{F_5})/b} }}{{\gamma + 1}}\\ \overline {{u_x}} = ({F_2} - \bar p\overline {{\rm{d}}{s_x}} )/{F_1}\\ \overline {{u_\theta }} = ({F_3} - \bar p\overline {{\rm{d}}{s_\theta }} )/{F_1}\\ \overline {{u_r}} = ({F_4} - \bar p\overline {{\rm{d}}{s_r}} )/{F_1}\\ \bar \rho = {F_1}/\left( {\overline {{u_x}} \overline {{\rm{d}}{s_x}} + \overline {{u_\theta }} \overline {{\rm{d}}{s_\theta }} + \overline {{u_r}} \overline {{\rm{d}}{s_r}} } \right) \end{array} \right. $ | （A2） |
式中：下标x、θ、r分别表示轴向、周向和径向；平均压力表达式中根号前面的“+”对应亚声速流动解，根号前面的“-”对应超声速流动解；参数a、b、c分别表示为
| $ \left\{ {\begin{array}{*{20}{l}} {a = {F_2}\overline {{\rm{d}}{s_x}} + {F_3}\overline {{\rm{d}}{s_\theta }} + {F_4}\overline {{\rm{d}}{s_r}} }\\ {b = {{\overline {{\rm{d}}{s_x}} }^2} + {{\overline {{\rm{d}}{s_\theta }} }^2} + {{\overline {{\rm{d}}{s_r}} }^2}}\\ {c = F_2^2 + F_3^2 + F_4^2} \end{array}} \right. $ | （A3） |
| [1] | OYAMA A, LIOU M, OBAYASHI S. Transonic axial-flow blade shape optimization using evolutionary algorithm and three-dimensional Navier-Stokes solver[J]. Journal of Propulsion and Power, 2004, 20(4): 302-310. |
| Click to display the text | |
| [2] | BENINI E. Three-dimensional multi-objective design optimization of a transonic compressor rotor[J]. Journal of Propulsion and Power, 2004, 20(3): 559-565. |
| Click to display the text | |
| [3] | LIAN Y, LIOU M. Multi-objective optimization of transonic compressor blade using evolutionary algorithm[J]. Journal of Propulsion and Power, 2005, 21(6): 979-987. |
| Click to display the text | |
| [4] | SAMAD A, KIM K, GOEL T, et al. Multi-objective optimization of transonic compressor blade using evolutionary algorithm[J]. Journal of Propulsion and Power, 2008, 24(2): 302-310. |
| Click to display the text | |
| [5] | JAMESON A. Aerodynamic design via control theory[J]. Journal of Scientific Computing, 1988, 3(3): 233-260. |
| Click to display the text | |
| [6] | JAMESON A. Optimum aerodynamic design using CFD and control theory: AIAA-1995-1729[R]. Reston: AIAA, 1995. |
| [7] | REUTHER J, JAMESON A. Aerodynamic shape optimization of wing and wing-body configurations using control theory: AIAA-1995-0123[R]. Reston: AIAA, 1995. |
| [8] | KIM S, ALONSO J J, JAMESON A. Two-dimensional high-lift aerodynamic optimization using the continuous adjoint method: AIAA-2000-4741[R]. Reston: AIAA, 2000. |
| [9] | KIM S, ALONSO J J, JAMESON A. Multi-element high-lift configuration design optimization using viscous continuous adjoint method[J]. Journal of Aircraft, 2004, 41(5): 1082-1098. |
| Click to display the text | |
| [10] | YANG S, WU H Y, LIU F, et al. Aerodynamic design of cascades by using an adjoint equation method: AIAA-2003-1068[R]. Reston: AIAA, 2000. |
| [11] | WU H Y, LIU F, TSAI H. Aerodynamic design of turbine blades using an adjoint equation method: AIAA-2005-1006[R]. Reston: AIAA, 2005. |
| [12] | LI H, SONG L, LI Y, et al. 2D viscous aerodynamic shape design optimization for turbine blades based on adjoint method[J]. Journal of Turbomachinery, 2011, 133(1): 031014. |
| Click to display the text | |
| [13] | WANG D, HE L. Adjoint aerodynamic design optimization for blades in multistage turbomachines Part Ⅰ:Methodology and verification[J]. Journal of Turbomachinery, 2010, 132(2): 021011. |
| Click to display the text | |
| [14] | WANG D, HE L. Adjoint aerodynamic design optimization for blades in multistage turbomachines-Part Ⅱ:Validation and applications[J]. Journal of Turbomachinery, 2010, 132(2): 021012. |
| Click to display the text | |
| [15] | LUO J, XIONG J, LIU F, et al. Three-dimensional aerodynamic design optimization of a turbine blades by using an adjoint method[J]. Journal of Turbomachinery, 2011, 133(1): 011026. |
| Click to display the text | |
| [16] | WALTHER B, NADARAJAH S. Constrained adjoint-based aerodynamic shape optimization of a single-stage transonic compressor[J]. Journal of Turbomachinery, 2013, 135(3): 021017. |
| Click to display the text | |
| [17] | LUO J, ZHOU C, LIU F. Multipoint design optimization of a transonic compressor blade by using an adjoint met-hod[J]. Journal of Turbomachinery, 2014, 136(5): 051005. |
| Click to display the text | |
| [18] | LUO J, LIU F. Multi-objective optimization of a transonic compressor blade by using an adjoint method[J]. AIAA Journal, 2015, 53(3): 797-801. |
| Click to display the text | |
| [19] | LUO J, LIU F, MCBEAN I. Turbine blade row optimization through endwall contouring by an adjoint method[J]. Journal of Propulsion and Power, 2015, 31(2): 505-518. |
| Click to display the text | |
| [20] | WALTHER B, NADARAJAH S. Optimum shape design for multirow turbomachinery configurations using a discrete adjoint approach and an efficient radial basis function deformation scheme for complex multiblock grids[J]. Journal of Turbomachinery, 2015, 137(3): 081006. |
| Click to display the text | |
| [21] | CUMPSTY N A. Compressor aerodynamics[M]. London: Longman Scientific & Technical, 1989: 47-61. |
| [22] | LIU F, JAMESON A. Multigrid Navier-Stokes calculations for three-dimensional cascades[J]. AIAA Journal, 1993, 31(10): 1785-1791. |
| Click to display the text | |
| [23] | ZHENG X, LIU F. Staggered upwind method for solving Navier-Stokes and k-ω turbulence model equations[J]. AIAA Journal, 1995, 33(6): 991-998. |
| Click to display the text | |
| [24] | LIU F, ZHENG X. A strongly-coupled time-marching method for solving the Navier-Stokes and k-ω turbulence model equations with multigrid[J]. Journal of Computational Physics, 1996, 128(2): 289-300. |
| Click to display the text | |
| [25] | SPALART P R, ALLMARAS S R. A one-equation turbulence model for aerodynamic flow: AIAA-1992-0439[R]. Reston: AIAA, 1992. |
| [26] | DENTON J D. The calculation of three-dimensional viscous flow through multistage turbomachines[J]. Journal of Turbomachinery, 1992, 114(1): 18-26. |
| Click to display the text | |
| [27] | GILES M. Nonreflecting boundary conditions for Euler equation calculations[J]. AIAA Journal, 1990, 28(12): 2050-2058. |
| Click to display the text | |
| [28] | HICKS R M, HENNE P A. Wing design by numerical optimization[J]. Journal of Aircraft, 1978, 15(7): 407-412. |
| Click to display the text | |
| [29] | GILES M, PIERCE N. An introduction to the adjoint approach to design[J]. Flow, Turbulence and Combustion, 2000, 65(10): 393-415. |
| Click to display the text | |
| [30] | ZYMARIS A, PAPADIMITRIOU D, GIANNAKOGLOU K, et al. Continuous adjoint approach to the Spalart-Allmaras turbulence model for incompressible flows[J]. Computers & Fluids, 2009, 38(8): 1528-1538. |
| Click to display the text | |
| [31] | LUO J, ZHU Y, TANG X, et al. Flow reconstruction and aerodynamic shape optimization of turbomachinery blades by POD-based hybrid models[J]. Science China Technological Sciences, 2017, 60(11): 1658-1673. |
| Click to display the text |