核電超長末葉氣動穩定性的流固耦合數值研究

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

(1.西安交通大學航天航空學院, 710049, 西安; 2.西安交通大學機械結構強度與振動國家重點實驗室, 710049, 西安; 3.東方電氣集團東方汽輪機有限公司, 618000, 四川德陽)

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

目前,許多研究文獻雖利用上述簡化耦合模型研究了流場流動和結構振動特征,然而大多數文獻集中于研究機翼或壓氣機葉片的顫振穩定性[11-14],對核電汽輪機末級葉片,尤其是本文所研究的超長核電末級葉片的氣動特征鮮有涉及。

因此,本文采用三維流固耦合時域響應算法,研究了新開發的核電汽輪機超長末級葉片在設計背壓3 600 Pa、高背壓7 200 Pa和10 800 Pa工況下的氣動穩定性問題。計算結果表明,在所計算的背壓工況內葉片維持氣動穩定,未發現顫振現象,但穩定性逐漸降低。

1 計算原理和方法

1.1 流體控制方程

為完全模擬核電葉片不同葉高的流動特征,采用三維Navier-Stokes(N-S)方程作為控制方程,該方程組表示流體運動所應滿足的質量守恒,動量守恒和能量守恒三大物理學定律,可寫為

(1)

式中:Q為瞬態矢量;F、G、H分別為3個方向上的通量;S為源性通量。

雖然瞬時N-S方程可以用來描述湍流,但其非線性特征不易求解,工程上常引入雷諾應力的方法將N-S方程中瞬態的變化量以時均化的形式表現,引入湍流控制方程把湍流的脈動值與時均值聯系起來,方程變為雷諾平均納維斯托克斯方程(RANS)。本文采用的k-ε湍流控制方程為

(2)

(3)

1.2 固體求解算法

固體域的動力學方程可寫為

(4)

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

1.3 流固耦合的方法及計算流程

流固耦合計算的實質是流體和固體的交界面上發生相互作用,進而產生耦合并相互影響。由于流場網格與固體網格并不一一對應,所以在流固耦合交界面上,需運用樣條方法,將流固交界面節點位移做轉換

Δxa=GΔxs

(5)

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

qs=GTqa

(6)

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

·(Γdispδ)=0

(7)

(8)

圖1 時域響應流固耦合算法的流程圖

2 末級汽輪機有限元網格模型

本文研究的對象為新開發的核電汽輪機,末級葉片長度超過2 m,靜葉數為66,動葉數為76,蒸汽從末級葉片的靜葉入口流入,隨后經動葉出口流出,幾何示意圖如圖2所示。

圖2 末級葉片通道幾何示意圖

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

圖3 設計背壓工況下的動葉中部背弧面流體壓力分布

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

圖4給出了核電末級葉片的網格,其中葉片固體網格含有13 041個單元。動靜葉流體計算域的網格均為結構化網格,對葉片周圍的流動區域設置了邊界層,內部計算域網格環繞邊界層分布直至充滿整個計算域。動、靜葉網格的平均第一層邊界層網格高度約為3.1×10-5m,對第一層網格高度做無量綱化處理,可得流場網格的y+≈4<5。

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

(c)靜葉葉頂流道 (d)動葉葉頂流道圖4 核電末級葉片的計算網格

3 計算結果與討論

3.1 流場特征分析

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

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

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

3.2 葉片氣動穩定性分析

(a)3 600 Pa響應曲線 (b)3 600 Pa頻譜圖

(c)7 200 Pa響應曲線 (d)7 200 Pa頻譜圖

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

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

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

表1 不同背壓下葉片振動至最大、最小位移時的流場軸向力

表2 不同背壓下葉片振動至最大、最小位移時的流場周向力

對數阻尼比的計算公式為

(9)

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

圖8 氣動阻尼比與背壓的關系

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

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

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

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

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

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

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

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

(a)最大位移時 (b)最小位移時圖12 背壓為10 800 Pa時葉片不同振動時刻的瞬態流線分布

4 結 論

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