劉穎超，陳 榴，戴 韌，楊其國

(上海理工大學(xué) 能源與動(dòng)力工程學(xué)院，上海 200093)

大功率凝汽式汽輪機(jī)排汽壓力低，低壓末級(jí)在高蒸汽比容和速度下工作。百萬千瓦級(jí)汽輪機(jī)的末級(jí)全速鋼制葉片長度已達(dá)1 300 mm，葉片頂部圓周速度約600 m/s。如果軸向絕對(duì)出口馬赫數(shù)接近1，則動(dòng)葉頂部截面葉型相對(duì)出口馬赫數(shù)可達(dá)1.7左右[1]。末級(jí)動(dòng)葉的氣動(dòng)設(shè)計(jì)是個(gè)難題，末級(jí)的低輪轂比和超長葉片導(dǎo)致末級(jí)從葉根到葉頂?shù)牧鲃?dòng)速度、角度和葉型變化劇烈，尤其是葉片頂部截面葉型在滿足跨音速流動(dòng)氣動(dòng)設(shè)計(jì)時(shí)，還必須滿足嚴(yán)苛的許可應(yīng)力要求。

跨音速葉柵流動(dòng)中不僅有激波和膨脹波，還有波系之間的相互作用，激波與邊界層之間相互干涉，流動(dòng)復(fù)雜而難以準(zhǔn)確計(jì)算葉柵出流參數(shù)。早期計(jì)算方法是二元特征線法，其通過假設(shè)兩相鄰葉片喉部處音速線是一條直線來求解渦輪葉柵跨音速流動(dòng)。Deich[2]通過分析得出特征線法在葉柵出口馬赫數(shù)低于1.4、低節(jié)弦比和薄直背葉型的葉柵流動(dòng)中相當(dāng)準(zhǔn)確。而對(duì)于高節(jié)弦比、高反動(dòng)度的葉柵頂部薄葉型，F(xiàn)orster[3]的實(shí)驗(yàn)結(jié)果說明葉柵喉部的直線音速線假設(shè)已不成立，計(jì)算方法失效。

葉柵跨音速流動(dòng)的實(shí)驗(yàn)測量也十分困難。平面葉柵實(shí)驗(yàn)中，葉片尾緣的激波在葉柵兩側(cè)自由流邊界反射，形成“反射激波”，影響葉柵出口流動(dòng)參數(shù),在大安裝角和葉片間距大的葉柵中很容易發(fā)生這種現(xiàn)象。另一方面，在超音速氣流中，氣動(dòng)探針難以精確測量當(dāng)?shù)乜倝海捎眉y影法可以獲得激波形態(tài)，但是不能量化評(píng)定葉柵效率[4]。

基于時(shí)間推進(jìn)法的CFD求解方法為分析跨音速葉柵流動(dòng)提供了有力工具。分析激波相關(guān)流動(dòng)，需要各種高精度格式，如總變差減小的差分格式(TVD)[5]、Roe迎風(fēng)格式[6]、非自由參數(shù)耗散有限差分格式(NND)[7]等，目前對(duì)包含中等激波強(qiáng)度的跨音速葉柵流動(dòng)計(jì)算基本沒有困難，大渦模擬(LES)[8]計(jì)算也是可行的。具有2階迎風(fēng)格式的計(jì)算方法基本可以滿足葉柵跨音速流動(dòng)分析的需要，結(jié)合適當(dāng)?shù)膬?yōu)化算法，可以完成跨音速葉柵的優(yōu)化設(shè)計(jì)[9-11]。

系列實(shí)驗(yàn)[12-13]證明，跨音速葉柵流動(dòng)損失與葉型幾何形狀密切相關(guān)。葉型背弧的后半段曲率是重要的影響因素，大曲率出氣邊的葉型損失在超音速時(shí)是亞音速時(shí)的4倍，其中激波誘導(dǎo)邊界層分離是損失增加的主要原因，這時(shí)采用直線背弧是更好的選擇。前緣氣流加速主要由葉片形狀決定，并在某種程度上與壓比無關(guān)，而尾緣對(duì)流場的影響主要取決于壓比，尾緣形狀的影響較小。各種經(jīng)典葉頂截面葉柵跨音速流動(dòng)[14]的流動(dòng)性能大致相當(dāng)，主要損失是葉柵尾緣激波和吸力面上的反射波損失。

筆者應(yīng)用CFD軟件分析某汽輪機(jī)末級(jí)動(dòng)葉頂部截面葉型流動(dòng)損失機(jī)理，結(jié)合葉柵氣動(dòng)設(shè)計(jì)方法，通過參數(shù)優(yōu)化，調(diào)控葉柵流動(dòng)中的激波強(qiáng)度和反射位置，降低葉柵內(nèi)與激波有關(guān)的流動(dòng)損失，提高葉柵流動(dòng)效率。

1 葉柵流動(dòng)分析的數(shù)值方法

1.1 驗(yàn)證模型與邊界條件

流動(dòng)分析選用CFX軟件，采用有限體積法求解雷諾平均N-S方程，守恒方程中的對(duì)流項(xiàng)和湍流輸運(yùn)方程均選擇2階差分格式。選取類似末級(jí)動(dòng)葉頂部截面的雙圓弧葉柵實(shí)驗(yàn)[15]作為參照，驗(yàn)證數(shù)值方法對(duì)跨音速葉柵流動(dòng)特性計(jì)算的可靠性。該葉型相對(duì)柵距為1.5，安裝角為30°，最大相對(duì)厚度為0.1。

參考實(shí)驗(yàn)工況，計(jì)算域進(jìn)口和出口分別給定總壓p*和靜壓p2，壓比p2/p*=0.23，進(jìn)口氣流角為35°，進(jìn)口總溫T*=280 K。周向采用周期性邊界條件，葉片表面為絕熱邊界條件。

1.2 網(wǎng)格驗(yàn)證

計(jì)算網(wǎng)格采用結(jié)構(gòu)網(wǎng)格，計(jì)算域進(jìn)口距葉型前緣距離為C(其中C為弦長)，出口距葉型尾緣距離為2C，邊界層伸長比小于1.1。如表1所示,采用3套不同質(zhì)量網(wǎng)格進(jìn)行流動(dòng)計(jì)算，其中Cp為葉型表面壓力系數(shù),湍流模型選用SST模型，湍流度為5%，收斂后求解殘差小于10-6。

表1 網(wǎng)格無關(guān)性比較Tab.1 Comparison of the grid independence

由表1可以看出，對(duì)出口馬赫數(shù)、激波前Cp峰值及相應(yīng)位置的預(yù)測，3套網(wǎng)格之間略有差別，中等網(wǎng)格與細(xì)網(wǎng)格的預(yù)測效果基本一樣，表明中等網(wǎng)格的求解精度滿足性能分析要求。

1.3 湍流模型驗(yàn)證

分別采用標(biāo)準(zhǔn)k-ε、k-ω和SST湍流模型進(jìn)行定常計(jì)算，得到的葉型表面壓力分布Cp如圖1所示，其中p表示當(dāng)?shù)貕毫?X表示葉型的橫坐標(biāo)。與實(shí)驗(yàn)對(duì)比發(fā)現(xiàn)，3種湍流模型預(yù)測葉型壓力面的壓力分布與實(shí)驗(yàn)值均較吻合，差異不大，但對(duì)吸力面的預(yù)測效果存在偏差。標(biāo)準(zhǔn)k-ε模型對(duì)近壁區(qū)求解精度不高，在0.8C以后壓力一直降低，無法準(zhǔn)確預(yù)測激波。k-ω和SST模型預(yù)測所得的吸力面激波前壓力高于實(shí)驗(yàn)值，激波位置與實(shí)驗(yàn)值吻合。在吸力面前緣附近，SST模型與實(shí)驗(yàn)值更加接近，因此筆者選擇SST模型。

圖1 采用不同湍流模型計(jì)算所得葉型表面壓力分布Fig.1 Blade surface pressure distribution of different turbulent flow models

2 跨音速葉柵的流動(dòng)分析

2.1 計(jì)算模型及邊界條件

以某汽輪機(jī)末級(jí)動(dòng)葉跨音速頂部截面葉型為原葉型，其葉柵結(jié)構(gòu)如圖2所示，葉片弦長為150 mm，相對(duì)柵距為0.926，最大厚度為6.09 mm，其他相關(guān)幾何參數(shù)如表2所示。

圖2 葉柵結(jié)構(gòu)示意圖Fig.2 Schematic diagram of cascade structure

表2 葉型幾何參數(shù)Tab.2 Geometric parameters of the profile

計(jì)算域進(jìn)口總壓p*為17 021 Pa，出口靜壓p2為4 845.4 Pa，進(jìn)口總溫T*為308.1 K。不考慮末級(jí)流道在子午面的流線翹曲和子午加速率。湍流模型采用SST模型，湍流度為5%。計(jì)算網(wǎng)格仍采用上述結(jié)構(gòu)網(wǎng)格。

2.2 葉柵流動(dòng)的激波結(jié)構(gòu)

圖3為原葉型表面和沿葉柵流道中線MN(見圖4中)的等熵馬赫數(shù)分布。圖4給出了葉柵流道內(nèi)馬赫數(shù)分布云圖，由圖4能清晰看到膨脹波與激波結(jié)構(gòu)。在吸力面前緣，氣流迅速達(dá)到超音速后再減速(AB1段)，產(chǎn)生了圖4中AB位置的激波。然后,超音速氣流在葉柵流道內(nèi)流經(jīng)相鄰尾緣壓力面附近，遇膨脹波K1加速，并產(chǎn)生反射膨脹波K2，加速至吸力面53%C處達(dá)到第2次速度峰值，即C點(diǎn)，等熵馬赫數(shù)接近1.8。之后遇到尾緣壓力面DC處斜激波減速(CD1段)，并在CE位置產(chǎn)生反射波(C2E1段減速)。由于出口馬赫數(shù)高，在尾緣吸力面FG處也形成一道激波。

圖3 葉柵內(nèi)等熵馬赫數(shù)分布Fig.3 Isentropic Mach number distribution in the cascade channel

圖4 葉柵流道內(nèi)馬赫數(shù)分布Fig.4 Mach number distribution in the cascade channel

尾緣處吸力面和壓力面激波通常稱為“燕尾波”。各種波的存在、傳播和相互作用，尾緣吸力面激波、反射波與尾跡(K3)干涉，以及激波誘導(dǎo)邊界層分離，使得葉柵流動(dòng)更加復(fù)雜。

2.3 葉柵氣動(dòng)損失分析

跨音速葉柵氣動(dòng)損失來源包括激波損失、邊界層損失和尾緣損失(取決于尾緣厚度和基壓)。衡量和評(píng)價(jià)不同損失對(duì)提高葉柵氣動(dòng)效率十分重要。

從熱力學(xué)分析，流場中不可逆熵增是損失來源，可以用熵增來衡量葉柵內(nèi)各種損失?？缫羲偃~柵中復(fù)雜激波造成的氣動(dòng)損失是熵增來源之一。利用激波前后的質(zhì)量、動(dòng)量和能量守恒，能計(jì)算激波引起的熵增Δsshock，表達(dá)式[16]如下：

(1)

式中：cv為比定容熱容；k為絕熱指數(shù)；Manormal為激波上游法向馬赫數(shù)，是由圖4中激波位置確定的平均變量。

以進(jìn)口溫度T0與進(jìn)口速度u0對(duì)激波熵增無量綱化,得到無量綱激波熵增Δsshock-normalized：

(2)

表3統(tǒng)計(jì)了葉柵內(nèi)不同位置的激波損失。由表3可知,尾緣處激波引起的損失最大，占比58.5%，其次是反射波，這是因?yàn)闇u輪的高壓比導(dǎo)致尾緣后基壓低，當(dāng)流體在尾緣周圍膨脹到這種低壓，在吸力面和壓力面氣流相遇的地方被強(qiáng)激波重新壓縮[17]，使熵值增加。而反射波位于喉口后，是高馬赫數(shù)集中區(qū)域，損失增大。

表3 葉柵不同位置激波損失Tab.3 Shock loss at different positions of the cascade

圖5給出了葉型邊界層厚度分布，在經(jīng)過相鄰尾緣壓力面激波后，吸力面位移厚度突增形成波峰，表明激波干涉邊界層，造成局部邊界層分離，邊界層增厚，氣動(dòng)損失增加。葉片邊界層內(nèi)黏性效應(yīng)導(dǎo)致能量耗散，是不可逆過程，邊界層損失用熵增Δsbl計(jì)算，定義[16]如下：

圖5 葉型邊界層厚度分布Fig.5 Boundary layer thickness distribution of the blade

(3)

由進(jìn)口溫度與進(jìn)口速度對(duì)邊界層引起的無量綱熵增為：

(4)

通過上式計(jì)算，邊界層無量綱熵增損失為0.157 8，其中還有壓力側(cè)激波與吸力面邊界層干涉而附加的損失，而無量綱激波熵增損失為0.255 4，能量損失更大，故削弱激波強(qiáng)度、減小激波相關(guān)損失對(duì)提升葉柵氣動(dòng)性能是首要的。通過調(diào)整優(yōu)化葉型幾何形狀，改變曲率分布，控制尾緣激波及其反射波，以降低葉柵內(nèi)激波損失。

3 葉型氣動(dòng)優(yōu)化方法

3.1 葉型的參數(shù)化

采用中弧線疊加厚度分布生成葉型[18]。中弧線采用3階Bezier曲線表示，如圖6所示。n階Bezier曲線表達(dá)式如下：

圖6 采用3階Bezier曲線表示的葉型中弧線Fig.6 Camber line of the blade represented by the third-order Bezier curve

(5)

式中：p(t)表示構(gòu)成Bezier曲線上的點(diǎn);Bi,n(t)為Bernstein函數(shù)；Pi為直角坐標(biāo)系下控制點(diǎn)的位置矢量。

葉型中弧線首、末兩點(diǎn)由葉型前緣點(diǎn)、尾緣點(diǎn)和安裝角確定。前后緣的斜率設(shè)定為進(jìn)、出氣流角的約束，控制點(diǎn)P2、P3只能分別沿直線P0P2和P3P4移動(dòng)，故中弧線形狀由2個(gè)參數(shù)XP2、XP3控制。

葉型厚度分布如圖7所示，以最大厚度控制點(diǎn)T3為界，左側(cè)用2階Bezier曲線表示，右側(cè)用3階Bezier曲線表示。通過葉型前緣厚度和尾緣厚度確定首、末2點(diǎn)T1、T6。控制點(diǎn)T2、T3、T4的厚度均為最大厚度，葉型厚度分布由控制點(diǎn)T3、T4和T5確定，故葉型厚度分布由5個(gè)參數(shù)XT3、YT3、XT4、XT5、YT5控制。

圖7 Bezier曲線表示葉型厚度分布Fig.7 Thickness distribution of the blade represented by Bezier curve

3.2 基于Kriging代理模型的優(yōu)化

采用基于克里金法(Kriging)的代理模型，對(duì)動(dòng)葉頂部截面葉型進(jìn)行氣動(dòng)優(yōu)化。Kriging代理模型[19]具有精確的全局近似設(shè)計(jì)空間，能準(zhǔn)確地預(yù)測線性和非線性函數(shù)，其表達(dá)式如下：

(6)

試驗(yàn)設(shè)計(jì)(DOE)是構(gòu)建Kriging代理模型的第一步，建立設(shè)計(jì)空間生成樣本。采用最優(yōu)拉丁超立方設(shè)計(jì)在設(shè)計(jì)變量范圍(見表4)內(nèi)生成30個(gè)樣本，作為Kriging代理模型的訓(xùn)練樣本集。

表4 設(shè)計(jì)變量取值范圍Tab.4 Value range of design variables

優(yōu)化的目標(biāo)函數(shù)選擇能評(píng)價(jià)葉柵氣動(dòng)損失的出口總壓損失系數(shù)Cp,t，其具體定義為：

(7)

式中：p1為進(jìn)口靜壓；p02為出口總壓。

葉型前緣和尾緣直徑分別為2.5 mm和1.4 mm。設(shè)定的約束條件是葉型截面面積的變化相對(duì)于原葉型不超過±10%。在Isight軟件上，基于DOE生成樣本，選擇Kriging代理模型建立目標(biāo)函數(shù)的近似模型，并在設(shè)計(jì)變量范圍內(nèi)，采用多島遺傳算法尋優(yōu)，得到最低Cp,t所對(duì)應(yīng)的設(shè)計(jì)變量數(shù)值,然后用CFX校核最優(yōu)葉型的性能。

3.3 優(yōu)化結(jié)果分析

經(jīng)過優(yōu)化計(jì)算，最終獲得優(yōu)化后葉型，圖8給出了優(yōu)化前后葉型的幾何構(gòu)型對(duì)比,圖中采用葉型弦長對(duì)葉型軸向和周向坐標(biāo)進(jìn)行無量綱化。優(yōu)化后葉型截面面積比原葉型小9.6%，最大厚度為5.78 mm，位置基本不變。在葉型80%C～90%C處，優(yōu)化葉型壓力面曲率明顯減小，型線為直線，有利于減小壓力梯度，削弱尾緣壓力面激波強(qiáng)度。通過Kriging模型預(yù)測得到Cp,t為0.215 6，CFX實(shí)際計(jì)算所得Cp,t為0.216 7，代理模型不確定度為0.51%。

圖8 優(yōu)化前后葉型對(duì)比Fig.8 Comparison between optimized and baseline profiles

葉片表面連續(xù)曲率影響葉片性能[20-21]。圖9比較了優(yōu)化前后葉型吸力面曲率的分布，優(yōu)化后葉型吸力面曲率整體更加連續(xù)平滑，在20%C～75%C處曲率均低于原葉型，75%C后優(yōu)化后葉型曲率略大于原葉型，從而決定葉型表面等熵馬赫數(shù)的分布形狀(見圖10)。由圖10可知，優(yōu)化后葉型前緣吸力面的速度減小，吸力面上的第2次速度峰值降低并向后移動(dòng)，延遲了壓力峰值點(diǎn)，尾緣壓力面激波入射位置后移。從圖11優(yōu)化前后沿葉柵流道中線MN的等熵馬赫數(shù)分布也可以看到，反射波位置后移，激波強(qiáng)度削弱。同時(shí)，優(yōu)化后葉型接近尾緣的壓力面維持了較高壓力，起到類似后加載的作用，有利于降低葉型氣動(dòng)損失。圖12給出了優(yōu)化后葉柵流道內(nèi)馬赫數(shù)分布。由圖12可以看出,優(yōu)化后葉柵喉口后高馬赫數(shù)集中區(qū)域縮小，吸力面反射波區(qū)域縮小，強(qiáng)度減弱。

圖9 優(yōu)化前后葉型吸力面曲率分布對(duì)比Fig.9 Comparison of curvature distribution of blade suction surface before and after optimization

圖10 優(yōu)化前后葉型表面等熵馬赫數(shù)對(duì)比Fig.10 Comparison of isentropic Mach number distribution of the blade surface before and after optimization

圖11 優(yōu)化前后沿葉柵流道中線MN的等熵馬赫數(shù)分布對(duì)比Fig.11 Comparison of isentropic Mach number distribution along the centerline MN of the cascade channel before and after optimization

圖12 優(yōu)化后葉柵流道內(nèi)馬赫數(shù)分布Fig.12 Contour distribution of Mach number distribution in cascade channel of optimized blade

圖13比較了優(yōu)化前后葉型不同位置的激波損失，優(yōu)化后葉型的各處與激波相關(guān)的損失均下降，前緣AB處激波基本消失，反射波損失減少26%，尾緣處激波損失減少34%，總的激波損失相比原葉型下降約32%。同時(shí)，激波強(qiáng)度減弱能減少對(duì)邊界層的干擾，邊界層損失降低。優(yōu)化后葉型的邊界層損失Δsbl-normalized為0.130 8，相比原葉型降低了17%。

圖13 優(yōu)化前后葉型不同位置激波損失對(duì)比Fig.13 Comparison of shock loss at different positions of the blade profile before and after optimization

圖14和圖15給出了優(yōu)化前后葉柵出口1.2倍軸向弦長Cax處尾跡速度和尾跡氣流角沿相對(duì)柵距的分布,圖中采用葉柵出口截面平均速度對(duì)尾跡速度進(jìn)行無量綱化。由圖14和圖15可以看出，優(yōu)化葉型低速幅值降低,尾跡寬度有所減小，尾跡氣流角幅值有所降低，有利于減小尾跡摻混損失。優(yōu)化后，吸力面?zhèn)任槽E氣流角的變化更加緩慢，壓力面?zhèn)任槽E速度降低，尾跡氣流角的變化相對(duì)平緩，表明尾緣燕尾波減弱，尾流干涉減弱。

圖14 優(yōu)化前后葉型尾跡速度分布Fig.14 Circumferential velocity distribution of the blade before and after optimization

圖15 優(yōu)化前后葉型尾跡氣流角分布Fig.15 Flow angle circumferential distribution of the blade before and after optimization

圖16比較了優(yōu)化前后葉型沿流向方向的總壓損失系數(shù)。由圖16可以看出，優(yōu)化后葉型的總壓損失降低，特別是在中間吸力面后，相鄰尾緣壓力面激波干擾減弱，降幅增大，表明優(yōu)化后葉型的氣動(dòng)性能提升。優(yōu)化后葉型出口總壓損失比原葉型降低約24%，質(zhì)量流量增加4.8%。

圖16 優(yōu)化前后葉型流向方向的總壓損失系數(shù)對(duì)比Fig.16 Comparison of total pressure loss in the blade flow direction before and after optimization

4 結(jié) 論

(1)末級(jí)動(dòng)葉跨音速頂部截面的葉柵流動(dòng)損失主要來自激波，激波帶來的熵增損失比邊界層產(chǎn)生的熵增損失多61%，則跨音速葉柵主要通過削弱激波來減少氣動(dòng)損失。

(2)優(yōu)化葉型曲率分布是降低激波損失的關(guān)鍵。葉型尾緣的厚度決定了尾緣燕尾波的強(qiáng)度，吸力面的曲率決定了反射波的位置和強(qiáng)度。優(yōu)化后葉型激波產(chǎn)生的熵增損失相比原葉型下降約32%，葉型出口總壓損失系數(shù)比原葉型減少約24%，優(yōu)化后葉型的氣動(dòng)性能明顯提升。