庾明達(dá), 徐自力, 范志飛, 曹守洪, 范小平, 方宇

(1.西安交通大學(xué)航天航空學(xué)院, 710049, 西安; 2.西安交通大學(xué)機(jī)械結(jié)構(gòu)強(qiáng)度與振動(dòng)國(guó)家重點(diǎn)實(shí)驗(yàn)室, 710049, 西安; 3.東方電氣集團(tuán)東方汽輪機(jī)有限公司, 618000, 四川德陽(yáng))

汽輪機(jī)中的顫振是一種持續(xù)、不衰減的自激振動(dòng)現(xiàn)象,是由汽輪機(jī)流道內(nèi)的氣流力與葉片發(fā)生相互作用引起,所以研究葉片的氣動(dòng)穩(wěn)定性需考慮流固耦合的作用。目前,研究顫振穩(wěn)定性應(yīng)用最多的流固耦合方法為能量法,仲繼澤等通過(guò)能量法預(yù)測(cè)了相位角、流量等參數(shù)對(duì)壓氣機(jī)葉片的顫振穩(wěn)定性的影響[1-2]。能量法思路清晰、簡(jiǎn)便,但由于假設(shè)過(guò)多,與實(shí)際情況不很相符[3]。另有很多學(xué)者采用流固耦合方法來(lái)模擬葉片的實(shí)際振動(dòng)響應(yīng)特征,Vahdati等用強(qiáng)耦合CFD方法研究了寬弦風(fēng)扇的失速顫振現(xiàn)象,但此方法求解難度大,占用資源高[4];Dowell等基于帶副翼操縱面的機(jī)翼3DOF模型,運(yùn)用ROM降階技術(shù)研究了臨界音速條件下馬赫數(shù)對(duì)機(jī)翼顫振速度和穩(wěn)定性的影響[5],這種氣動(dòng)力降階模型目前多局限于機(jī)翼等具有簡(jiǎn)單幾何外形結(jié)構(gòu)的流固耦合計(jì)算;文獻(xiàn)[6-7]利用模態(tài)疊加法將葉片響應(yīng)形式分解為單階模態(tài)的線(xiàn)性振動(dòng),研究了跨音速壓氣機(jī)的顫振邊界。還有不少學(xué)者使用二維模型模擬壓氣機(jī)的流場(chǎng)特征[8-10],但二維模型難以完全模擬高負(fù)載葉輪機(jī)械的流動(dòng)狀況。

目前,許多研究文獻(xiàn)雖利用上述簡(jiǎn)化耦合模型研究了流場(chǎng)流動(dòng)和結(jié)構(gòu)振動(dòng)特征,然而大多數(shù)文獻(xiàn)集中于研究機(jī)翼或壓氣機(jī)葉片的顫振穩(wěn)定性[11-14],對(duì)核電汽輪機(jī)末級(jí)葉片,尤其是本文所研究的超長(zhǎng)核電末級(jí)葉片的氣動(dòng)特征鮮有涉及。

因此,本文采用三維流固耦合時(shí)域響應(yīng)算法,研究了新開(kāi)發(fā)的核電汽輪機(jī)超長(zhǎng)末級(jí)葉片在設(shè)計(jì)背壓3 600 Pa、高背壓7 200 Pa和10 800 Pa工況下的氣動(dòng)穩(wěn)定性問(wèn)題。計(jì)算結(jié)果表明,在所計(jì)算的背壓工況內(nèi)葉片維持氣動(dòng)穩(wěn)定,未發(fā)現(xiàn)顫振現(xiàn)象,但穩(wěn)定性逐漸降低。

1 計(jì)算原理和方法

1.1 流體控制方程

為完全模擬核電葉片不同葉高的流動(dòng)特征,采用三維Navier-Stokes(N-S)方程作為控制方程,該方程組表示流體運(yùn)動(dòng)所應(yīng)滿(mǎn)足的質(zhì)量守恒,動(dòng)量守恒和能量守恒三大物理學(xué)定律,可寫(xiě)為

(1)

式中:Q為瞬態(tài)矢量;F、G、H分別為3個(gè)方向上的通量;S為源性通量。

雖然瞬時(shí)N-S方程可以用來(lái)描述湍流,但其非線(xiàn)性特征不易求解,工程上常引入雷諾應(yīng)力的方法將N-S方程中瞬態(tài)的變化量以時(shí)均化的形式表現(xiàn),引入湍流控制方程把湍流的脈動(dòng)值與時(shí)均值聯(lián)系起來(lái),方程變?yōu)槔字Z平均納維斯托克斯方程(RANS)。本文采用的k-ε湍流控制方程為

(2)

(3)

1.2 固體求解算法

固體域的動(dòng)力學(xué)方程可寫(xiě)為

(4)

式中:M、C、K分別為質(zhì)量、阻尼和剛度矩陣;q(t)為位移;p(t)為作用在流固耦合交界面上的氣動(dòng)壓力。采用直接積分法求解葉片的振動(dòng)方程得到葉片的有限元節(jié)點(diǎn)振動(dòng)位移。

1.3 流固耦合的方法及計(jì)算流程

流固耦合計(jì)算的實(shí)質(zhì)是流體和固體的交界面上發(fā)生相互作用,進(jìn)而產(chǎn)生耦合并相互影響。由于流場(chǎng)網(wǎng)格與固體網(wǎng)格并不一一對(duì)應(yīng),所以在流固耦合交界面上,需運(yùn)用樣條方法,將流固交界面節(jié)點(diǎn)位移做轉(zhuǎn)換

Δxa=GΔxs

(5)

式中:xa為流場(chǎng)位移;xs為固體位移;G為樣條矩陣。反之,流場(chǎng)氣動(dòng)力需要施加到固體交界面上,根據(jù)做功等效原理有

qs=GTqa

(6)

式中:qa、qs為施加在流固交界面上的氣動(dòng)力。在流固耦合計(jì)算過(guò)程中,流場(chǎng)完成內(nèi)部迭代計(jì)算后,按照式(6),與固體交界面交換壓力等變量,之后固體產(chǎn)生位移,此時(shí)流場(chǎng)網(wǎng)格變形,以便進(jìn)行下一次的迭代計(jì)算。在流固耦合交界面處,流場(chǎng)網(wǎng)格節(jié)點(diǎn)直接進(jìn)行變形,其他較遠(yuǎn)處的網(wǎng)格變形遵循以下位移擴(kuò)散方程

·(Γdispδ)=0

(7)

(8)

圖1 時(shí)域響應(yīng)流固耦合算法的流程圖

2 末級(jí)汽輪機(jī)有限元網(wǎng)格模型

本文研究的對(duì)象為新開(kāi)發(fā)的核電汽輪機(jī),末級(jí)葉片長(zhǎng)度超過(guò)2 m,靜葉數(shù)為66,動(dòng)葉數(shù)為76,蒸汽從末級(jí)葉片的靜葉入口流入,隨后經(jīng)動(dòng)葉出口流出,幾何示意圖如圖2所示。

圖2 末級(jí)葉片通道幾何示意圖

本文采用三維計(jì)算模型和循環(huán)對(duì)稱(chēng)邊界的單通道流道作為計(jì)算對(duì)象,用以計(jì)算末級(jí)流場(chǎng),這樣可有效減少計(jì)算規(guī)模,提高計(jì)算效率。為了選取合適的網(wǎng)格,本文采用由疏到密具有不同單元數(shù)量、不同網(wǎng)格尺寸的4套流場(chǎng)網(wǎng)格,對(duì)設(shè)計(jì)背壓工況的定常流動(dòng)進(jìn)行了計(jì)算,得到了流體施加在動(dòng)葉片中部截面背弧面的壓力分布,如圖3所示。圖3中,X為葉片中部展弦方向任意一點(diǎn)的位置坐標(biāo),C為葉片中部的弦長(zhǎng)。

圖3 設(shè)計(jì)背壓工況下的動(dòng)葉中部背弧面流體壓力分布

由圖3可知,隨著網(wǎng)格單元數(shù)的不斷增多,計(jì)算得到的背弧面靜壓力分布趨于定值,在流場(chǎng)網(wǎng)格單元數(shù)為707 100時(shí),繼續(xù)增加網(wǎng)格單元數(shù),背弧面的靜壓力計(jì)算結(jié)果也基本保持不變。因此,對(duì)該核電末級(jí)葉片,單元數(shù)選取707 100。

圖4給出了核電末級(jí)葉片的網(wǎng)格,其中葉片固體網(wǎng)格含有13 041個(gè)單元。動(dòng)靜葉流體計(jì)算域的網(wǎng)格均為結(jié)構(gòu)化網(wǎng)格,對(duì)葉片周?chē)牧鲃?dòng)區(qū)域設(shè)置了邊界層,內(nèi)部計(jì)算域網(wǎng)格環(huán)繞邊界層分布直至充滿(mǎn)整個(gè)計(jì)算域。動(dòng)、靜葉網(wǎng)格的平均第一層邊界層網(wǎng)格高度約為3.1×10-5m,對(duì)第一層網(wǎng)格高度做無(wú)量綱化處理,可得流場(chǎng)網(wǎng)格的y+≈4<5。

(a)末級(jí)流道整體 (b)葉片固體

(c)靜葉葉頂流道 (d)動(dòng)葉葉頂流道圖4 核電末級(jí)葉片的計(jì)算網(wǎng)格

3 計(jì)算結(jié)果與討論

3.1 流場(chǎng)特征分析

采用單通道模型,在流場(chǎng)計(jì)算域中,動(dòng)靜葉流道兩側(cè)均施加循環(huán)對(duì)稱(chēng)邊界,將靜葉與動(dòng)葉流場(chǎng)域統(tǒng)一視為內(nèi)部區(qū)域,以此模擬動(dòng)葉與靜葉之間的相對(duì)旋轉(zhuǎn)運(yùn)動(dòng),葉片表面、葉頂、葉轂為無(wú)滑移壁面,靜葉進(jìn)口設(shè)置蒸汽的流量和溫度,在動(dòng)葉出口設(shè)置壓力值。對(duì)設(shè)計(jì)背壓3 600 Pa和7 200、10 800 Pa兩個(gè)高背壓工況下的定常流動(dòng)進(jìn)行了計(jì)算,得到了3 600、7 200 Pa工況下末級(jí)流道截面流線(xiàn)分布,如圖5所示,可知在設(shè)計(jì)背壓時(shí),流體從動(dòng)葉流出后,流動(dòng)保持平穩(wěn),而在提高背壓后,流動(dòng)未發(fā)生分離或是方向改變。而當(dāng)背壓提高至10 800 Pa時(shí),動(dòng)葉尾部流道內(nèi)發(fā)生了旋渦現(xiàn)象,此時(shí)流線(xiàn)分布及壓力分布如圖6所示。由圖6可知,流體是順著壓力梯度向下游流動(dòng),其流動(dòng)趨勢(shì)與壓力梯度保持一致,然而在該背壓工況下動(dòng)葉流道中部的壓力梯度發(fā)生了明顯的波動(dòng),出現(xiàn)幾個(gè)局部低壓凸口,流體順著壓力梯度梯度流經(jīng)這些區(qū)域時(shí),由于動(dòng)能已經(jīng)減小,不足以克服壓差導(dǎo)致流向改變,少量流體在動(dòng)葉流道近出口處形成了旋渦。

(a)背壓為3 600 Pa (b)背壓為7 200 Pa圖5 子午面流線(xiàn)分布

(a)流線(xiàn)分布 (b)壓力分布圖6 背壓為10 800 Pa時(shí)子午面流線(xiàn)分布及壓力分布

3.2 葉片氣動(dòng)穩(wěn)定性分析

(a)3 600 Pa響應(yīng)曲線(xiàn) (b)3 600 Pa頻譜圖

(c)7 200 Pa響應(yīng)曲線(xiàn) (d)7 200 Pa頻譜圖

(e)10 800 Pa響應(yīng)曲線(xiàn) (f)10 800 Pa頻譜圖圖7 各背壓工況下的尾緣響應(yīng)曲線(xiàn)及其頻譜圖

在不同背壓工況下動(dòng)葉頂尾緣點(diǎn)上的時(shí)間歷程響應(yīng)曲線(xiàn)及其頻譜圖如圖7所示。由圖7可知,3個(gè)背壓工況下的響應(yīng)曲線(xiàn)全部現(xiàn)為衰減趨勢(shì),響應(yīng)曲線(xiàn)中的主頻成分與動(dòng)葉的第1階固有動(dòng)頻ω1的比值接近1,即葉片是以動(dòng)葉的第1階固有動(dòng)頻為主進(jìn)行自由衰減響應(yīng),對(duì)葉片進(jìn)行的模態(tài)分析結(jié)果顯示,葉片第2,3,4階模態(tài)與第1階模態(tài)的固有頻率之比分別為1.61,2.57,3.39。由圖7d可知,頻譜圖中含有葉片2階固有頻率附近的振動(dòng)成分,表明在背壓為7 200 Pa時(shí),葉片表現(xiàn)為1階、2階模態(tài)的疊加振動(dòng)。隨著背壓的逐步增大,葉片的幅值均分別呈現(xiàn)出逐漸減小的趨勢(shì),即葉片在這3個(gè)背壓工況下都將保持氣動(dòng)穩(wěn)定,不會(huì)發(fā)生顫振。

表1、表2給出了在3種工況下,葉片振動(dòng)至最大位移和最小位移時(shí)流場(chǎng)的軸向力和周向力。由表1、表2可知,在3個(gè)背壓條件下,當(dāng)葉片從最大位移處振動(dòng)至最小位移處時(shí),流場(chǎng)的軸向力和周向力均發(fā)生了變化,即固體的振動(dòng)影響了流體氣動(dòng)力的分布。

表1 不同背壓下葉片振動(dòng)至最大、最小位移時(shí)的流場(chǎng)軸向力

表2 不同背壓下葉片振動(dòng)至最大、最小位移時(shí)的流場(chǎng)周向力

對(duì)數(shù)阻尼比的計(jì)算公式為

(9)

式中:m為周期間隔;A為振動(dòng)幅值。計(jì)算時(shí),由于該核電葉片的材料阻尼很小,式(7)的阻尼比主要為氣動(dòng)阻尼比。氣動(dòng)阻尼比與背壓的關(guān)系如圖8所示,可知隨著背壓的提高,3個(gè)工況下氣動(dòng)阻尼均為正,但提高背壓后阻尼比隨之減小,葉片的氣動(dòng)穩(wěn)定性隨之減弱。

圖8 氣動(dòng)阻尼比與背壓的關(guān)系

背壓工況下振動(dòng)歷程的相平面圖如圖9所示。由圖9可知:在背壓為3 600、7 200 Pa時(shí)葉片振動(dòng)的速度和位移均呈現(xiàn)出雜亂的變化趨勢(shì),相平面圖呈現(xiàn)出多個(gè)頻率和幅值的疊加運(yùn)動(dòng),氣動(dòng)阻尼對(duì)結(jié)構(gòu)有顯著影響,相平面呈現(xiàn)不規(guī)則的收斂趨勢(shì),體現(xiàn)了強(qiáng)烈的非線(xiàn)性特征;背壓升至10 800 Pa時(shí),頻率成分變得單一,此時(shí)的阻尼比相當(dāng)小,有限時(shí)間段內(nèi)葉片在相平面圖上收斂趨勢(shì)變得更加緩慢,隨著時(shí)間的進(jìn)一步推進(jìn),葉片最終收斂到靜平衡位置。

(a)背壓為3 600 Pa (b)背壓為7 200 Pa

(c)背壓為10 800 Pa圖9 各背壓下振動(dòng)響應(yīng)的相平面圖

80%葉高截面流場(chǎng)馬赫數(shù)云圖和截面展弦向葉片壓力分布如圖10、11所示。由圖10、11可知,葉柵流道內(nèi)出現(xiàn)了激波,激波不斷朝葉片前緣移動(dòng),在背壓為3 600、7 200 Pa時(shí),流場(chǎng)的背弧面均可觀察到穩(wěn)定的激波;到背壓為10 800 Pa時(shí),葉片背弧面上出現(xiàn)了低能流團(tuán),同時(shí)邊界層出現(xiàn)了分離,并促使激波結(jié)構(gòu)發(fā)生了變化,此時(shí)的激波分為兩支,一支為前緣正激波,一支與流道平行使整個(gè)激波呈“人”字形分布。

(a)背壓為3 600 Pa (b)背壓為7 200 Pa

(c)背壓為10 800 Pa圖10 80%葉高截面流場(chǎng)馬赫數(shù)云圖

圖11 80%葉高截面展弦向葉片壓力分布

文獻(xiàn)[15-16]指出,流道內(nèi)激波的位置變化是導(dǎo)致葉片氣動(dòng)彈性穩(wěn)定性差異的主要誘因,當(dāng)激波從葉尖前緣脫體后,葉片從氣流中吸收能量導(dǎo)致振動(dòng)發(fā)散并迅速喪失其氣動(dòng)穩(wěn)定性。本文研究中,背壓為10 800 Pa工況時(shí)動(dòng)葉流道內(nèi)雖已出現(xiàn)靠近前緣的分岔激波,但并沒(méi)有脫離動(dòng)葉前緣葉柵流道,背壓的提高使動(dòng)葉柵內(nèi)的激波向前移動(dòng),但激波沒(méi)有發(fā)生脫體,這在一定程度上解釋了3個(gè)工況下核電葉片保持氣動(dòng)穩(wěn)定的原因;在葉片振動(dòng)相對(duì)劇烈的葉高中上部流道內(nèi),流動(dòng)保持穩(wěn)定,沒(méi)有紊流產(chǎn)生。當(dāng)背壓為10 800 Pa及葉片開(kāi)始振動(dòng)后的最大、最小位移時(shí),提取葉片0.5、0.7葉高截面和靠近葉尖振動(dòng)觀察點(diǎn)的0.95葉高截面的流動(dòng)分布情況,如圖12所示。由圖12可知,葉片最大、最小位移時(shí)葉片3個(gè)截面的動(dòng)葉葉柵的氣流分別沿著葉片內(nèi)、背弧面向下游流動(dòng),未發(fā)生附面層分離,葉片保持氣動(dòng)穩(wěn)定,同時(shí)流場(chǎng)在各個(gè)截面處的流動(dòng)平穩(wěn)。

(a)最大位移時(shí) (b)最小位移時(shí)圖12 背壓為10 800 Pa時(shí)葉片不同振動(dòng)時(shí)刻的瞬態(tài)流線(xiàn)分布

4 結(jié) 論

本文采用三維時(shí)域響應(yīng)流固耦合算法,研究了新開(kāi)發(fā)核電汽輪機(jī)末級(jí)葉片從設(shè)計(jì)背壓為3 600 Pa逐漸提高到高背壓為7 200、10 800 Pa條件下的氣動(dòng)穩(wěn)定性問(wèn)題。計(jì)算結(jié)果表明,隨著背壓的升高,在葉柵流道內(nèi)的近出口處引起了流體壓力的擾動(dòng),在壓差的影響下逐漸形成了旋渦現(xiàn)象。在從設(shè)計(jì)背壓提高至兩個(gè)高背壓工況下,核電汽輪機(jī)葉片的時(shí)域響應(yīng)均以1階動(dòng)頻為主分別進(jìn)行自由衰減,表明該核電超長(zhǎng)葉片在所計(jì)算的背壓范圍內(nèi)將保持氣動(dòng)穩(wěn)定。背壓提高后,動(dòng)葉柵和流道內(nèi)的激波隨著背壓的增大不斷向前緣移動(dòng),但并未發(fā)生脫體,葉片振動(dòng)至不同時(shí)刻時(shí)在葉高各截面流道內(nèi)流動(dòng)沒(méi)有發(fā)生分離,這解釋了葉片保持氣動(dòng)穩(wěn)定的原因。但是,氣動(dòng)阻尼比隨背壓的提高逐漸消減,氣動(dòng)穩(wěn)定性隨背壓的提高而逐漸減弱。