唐明智 ,金東海 ,2,郭 昕 ,桂幸民 ,2
(1.北京航空航天大學能源與動力工程學院,北京100191;2.先進航空發動機協同創新中心,北京100191;3.中國航發四川燃氣渦輪研究院,成都610500)
現階段3D數值模擬在風扇/壓氣機設計中得到了廣泛使用。然而在設計的初級階段,通流模型由于一方面可以快速提供風扇/壓氣機的性能及內部流動特性,另一方面也易于將經驗參數納入到數值模擬中[1-2],因此仍具有重要作用[3-4]。目前應用最廣泛的通流方法是流線曲率法,而基于周向平均Euler[5]或Navier-Stokes(N-S)方程[6-8]的通流模型也同樣被廣泛研究,并應用于壓氣機的特性預測和流場模擬。這些通流模型適用于跨聲壓氣機的流場模擬和性能預測[9]。
在對N-S方程進行周向平均的過程中,不可避免地產生一些附加項,其中與表面力相關的項被模化為無黏葉片力項和葉片黏性項,而對于由于流場的周向不均勻性及控制方程的非線性所引發的高階項,即周向脈動應力項,早期研究認為其對流場的影響可以忽略,即假設葉片通道進口前流動參數沿周向均勻分布,這主要是由于在早期研究中,相對于控制方程中的其他項,反映周向不均勻性的項的量級很小;另一個原因則是設計者缺乏模型來預測這些項的值[10]。然而,隨著當前壓氣機負荷的提高,以及彎掠葉片技術的廣泛應用,周向不均勻性的影響正被逐步認識到,有助于更為精確地應用彎掠葉片技術。事實上,這些項的作用可能會比黏性項的作用更大,并且隨著負荷的增大,其作用也會增強[11]。此外,周向不均勻性會影響流動參數的軸向和展向的分布[9,12-13],并能反映角區失速和徑向摻混等現象。對于掠葉片,葉片通道進口的流動平衡也會被周向不均勻性所影響,并會誘導流動參數的重新分配,而引發周向不均勻性的1個主要來源就是無黏葉片力[14]。
為了模化周向不均勻性,通常的方法是計算周向脈動應力項,其中1種方法是開展S1流面上的計算,然后提取出周向脈動應力項。吳仲華采用中心流線法來簡化對S1流面的計算[15],并在之后得到了一定的發展[16-17]。由于周向脈動應力項具有諧波特征[18],因此非線性諧波平衡法也被用來模化周向不均勻性[19],不過這種方法并未獲得廣泛應用。
文獻[20]通過應力輸運模型對周向脈動應力項進行建模,分析了具有不同掠角的掠葉柵在0°迎角下流場中周向不均勻性的影響,結果表明周向不均勻性會改變掠葉片進口流動參數,并引發進口新的徑向平衡。本文在此基礎上應用該應力輸運模型分析迎角改變時掠葉柵進口流場中周向不均勻性的作用。
在相對柱坐標系中對N-S方程進行周向平均[10],并將周向脈動應力項整合為周向脈動源項P,可得到通流模型的主控方程

方程中各項定義為

式中:U為守恒量;F和G為對流(無黏)通量;Fv和Gv為擴散(黏性)通量;S為N-S方程組在相對柱坐標下導出的源項;FB為無黏葉片力;FF為黏性葉片力。
定義由于葉片切向厚度產生的堵塞系數b

式中:Δφ=φs-φp,為葉片通道周向寬度;N為葉片數。
在葉片區b<1,在非葉片區b=1。周向脈動源項P中的高階項即為周向脈動應力項,這些項是由于葉輪機的周向不均勻性及N-S方程的非線性所導致的,其與周向脈動源項均能夠反映周向不均勻性的大小。
與文獻[20]類似,為分析葉片通道進口徑向平衡的改變,對方程(1)中的徑向動量方程進行推導及簡化,可以得到無黏形式的完全徑向平衡方程

式中:等式左邊為徑向壓力梯度(RGP);等式右邊分別為周向速度引發的離心加速度項(CENT_W)、由于子午流線的曲率而引發的離心加速度的徑向分量(CENT_M)、子午速度變化所引發的加速度的徑向分量(AC_M)、周向脈動源項的徑向分量 P'r,P'r與 Pr的區別在于二者的比值為密度,而下面將仍以Pr表示周向脈動源項的徑向分量,并不再重復說明;FBr為無黏葉片力的徑向分量,由于主要分析周向不均勻性對葉片通道進口流動平衡的影響,該處FBr=0,因此后面將不給出FBr的分布情況。
類似雷諾平均后所產生的雷諾應力項,式(1)中的周向脈動應力項代表著周向不均勻性所帶來的影響,本節將簡要介紹周向脈動應力項的應力輸運模型,詳細建模過程可參考文獻[20]。類比湍流模型中的應力輸運模型,可以推導出對于周向脈動應力項的應力輸運方程,以張量形式表示的無黏形式的應力輸運方程為

式中:方程左側為對流項(CON);右側分別為生成項(PRO)、輸運項(TRA)、速度 -壓力關聯項(PRE,簡稱壓力項)以及速度-源項(V-S)關聯項,其中生成項表征應力輸運方程與平均后的軸對稱流場之間的關聯,速度-壓力關聯項表征由于壓力的脈動性而做的功,輸運項、壓力項及速度-源項關聯項均無法直接獲得,需要進行建模,而在這幾項中,壓力項為關鍵項,其余幾項則可忽略,因此下面簡要介紹壓力項的建模方法,詳細過程可參考文獻[20]。
假設速度和壓力可沿周向進行傅里葉級數展開,并保留到1階,則壓力項可進一步推導為

考慮到密度加權平均的定義,φ0和φ1之間滿足關系

從式中可見,當確定3個速度和壓力的周向偏導數后,可解出壓力項。進一步引入S2流面無黏形式的周向動量方程以及絕對坐標系下流動無旋假設,并結合連續方程和無黏形式下的能量方程,可推導出壓力和3個速度分量的周向偏導數,其形式為

由此可實現對壓力項的求解,從而實現應力輸運方程的封閉。另外從式(7)、(8)中可見,引發速度和壓力的周向不均勻性的1個主要因素是環量沿軸向和徑向的偏導數,即負荷的分配情況,而文獻[21]中也基于S2流面上周向動量方程分析了周向不均勻性產生的原因,認為對于定常問題周向不均勻性來源于黏性力的存在和氣流的環量沿流線的導數不為零,而在未發生大范圍分離時,相比黏性力,環量改變所反映的無黏葉片力是最主要的源項,因此環量改變對周向不均勻性的產生應起主要作用。
環量的分布特征在設計問題中是給定量,在分析問題中則主要是源于無黏葉片力的作用,因此該模型同時適用于設計問題和分析問題。
采用1個前掠20°葉柵作為研究對象,該葉柵基元的厚度分布采用文獻 [22]中的設計進口馬赫數為0.62的可控擴散葉型,中弧線采用任意造型的方式獲得。葉柵的基本參數見表1。由于缺乏試驗數據作為驗證,本文采用3D數值模擬的結果作為基準,并應用了NUMECA FINE/Turbo3D數值模擬軟件對各葉柵算例進行數值模擬研究。通過對某單級跨聲速風扇、某風扇/增壓級以及某單級離心壓氣機的數值計算結果與試驗結果的對比校核,NUMECA FINE/Turbo的計算精度能夠滿足壓氣機設計過程中的流場驗算和分析[23]。邊界條件的設置為:在計算域的進口處給定總溫、總壓以及氣流方向,通過改變進口氣流角來獲得不同的名義迎角,這里的名義迎角是指遠前方均勻來流氣流角與葉片幾何進口角的夾角,后面簡稱為迎角,給定的名義迎角為0°、±2°和±4°;采用當地迎角的概念,即葉片通道進口處氣流角與葉片幾何進口角的夾角。在出口給定流量以控制進口馬赫數,對于各個迎角,進口馬赫數均為0.62;環壁和葉片表面采用絕熱無滑移邊界。湍流模型采用Spalart-Allmaras(S-A)方程模型,網格劃分保證固壁面第1層網格單元y≤10。3D數值模擬計算網格采用AUTOGRID模塊自動生成,如圖1所示。對于不同迎角算例,均采用同一套網格,總網格數約為162萬。

表1 葉柵基本參數

圖1 3D數值模擬網格
另外,為分析周向不均勻性的影響,開展了以下2種類型的周向平均通流計算:
(1)CAM:未考慮周向不均勻性的周向平均通流計算;
(2)CAM+MODEL:加入了以應力輸運模型計算的周向脈動應力項的周向平均通流計算。
周向平均通流模型的計算網格如圖2所示,總網格數為3185。通流模型采用基于時間推進的有限體積方法來求解控制方程,時間離散格式為顯式格式為Jameson4步Runge-Kutta法。對無黏通量的空間離散采用Edwards[24]提出的低耗散通量分裂格式(LDFSS)。對于邊界條件的設置,進出口邊界條件與3D數值模擬相同,環壁采用對稱邊界條件。

圖2 通流計算網格
本文主要關注的是周向不均勻性對葉片通道進口流場的影響,此處黏性作用不明顯,雖然無黏葉片力在此處為0,但壓力勢已經會以周向不均勻性的形式造成一定影響,下面將分別從葉片通道進口氣流角和徑向平衡的角度進行分析。
為考察迎角改變時周向不均勻性對葉片通道進口流場的影響,首先分析在引入周向不均勻性后葉片通道進口氣流角的變化,本節將靠近葉片前緣處的氣流角與計算域進口氣流角做差以便于分析這一變化,結果如圖3所示。從圖3(a)中可見,對于20%展高處,隨著迎角的增大,3D和CAM預測的差值有一定的增長趨勢,而CAM+MODEL則呈現出下降趨勢,這主要是由于隨著迎角增大,Pu的作用增強(圖4(a)),導致周向速度變化更大,但總體來說CAM+MODEL所預測的差值CAM要更接近3D的,在名義迎角從-4°增大到4°的過程中,CAM+MODEL與CAM所預測的氣流角差值間的偏差量從約0.1°增加到了約0.7°,預測精度提高了25%以上;對于50%和80%展高處,在加入周向不均勻性的影響后,對氣流角的預測精度同樣提高,并且隨著名義迎角的增大,3D和CAM+MODEL預測出的αLE-αinlet的值均呈現出下降趨勢;與20%展高處的結果類似,加入周向不均勻性模型后通流模型對氣流角改變量的預測精度也同樣有所提高,提高50%以上。

圖3 掠葉柵不同展高處αLE-αinlet隨名義迎角的變化
由于葉片通道進口處無黏葉片力為0,因此上述氣流角的改變是由于周向不均勻性所誘發的,周向脈動源項的軸向和周向分量Px、Pu會分別影響葉片通道進口軸向和周向速度,從而造成氣流角的改變,并且隨著迎角的增大和負荷的提高,周向不均勻性的影響也在增大。
由于Pu會對環量造成影響,因此對其預測是否準確很重要。葉柵不同展高截面葉片上游距前緣約4%弦長處無量綱Pu隨迎角的變化情況如圖4所示,其中Pu是由計算域進口對應流動參數進行的無量綱化。對于3個不同展高處,Pu的分布呈現出相近的特征,都是隨著迎角的增大Pu的絕對值增加,說明其作用將增強。此外,對于不同迎角,應力輸運模型計算得到的Pu與3D結果的偏差在不同截面都很小,都在20%以內,說明對于非設計工況,應力輸運模型能夠很好地預測周向脈動源項。
在加入應力輸運模型后,通流模型的計算時間比未加入應力輸運模型時的約增加40%,這是由于在主控方程中加入6個微分方程,但計算時間仍遠遠少于3D數值模擬。各算法的計算時間見表2。


圖4 周向分量Pu隨名義迎角的變化掠葉柵不同展高處周向脈動源項
上述分析表明,在加入周向不均勻性后,葉片通道進口氣流角發生了改變,即當地迎角發生了改變。對于掠葉柵,各展向截面當地迎角的變化會改變葉片負荷沿展向的分布情況,進而改變進口徑向流動平衡。為進一步考察名義迎角和當地迎角改變對徑向平衡的影響,給出了不同名義迎角下葉柵上游距前緣約4%弦長處徑向平衡方程(3)中各項沿展向的分布情況(如圖5所示),并以流道進口的流動參數對這些項進行了無量綱化。圖中的命名規則以“-2_3D”為例,其中-2表示名義迎角為-2°,3D指3維計算的結果。
對比圖 5(a)~(e)中各項,從量上來看,與 RGP處于同一量級的項為CENT_M和Pr,而CENT_W和AC_M則要比其余3項小1~2個量級,表明周向不均勻性明顯影響葉片通道進口徑向平衡及徑向壓力梯度。

表2 葉柵基本參數


圖5 掠葉柵進口徑向平衡參數沿展向的分布隨名義迎角的變化
對于RGP隨名義迎角的變化情況(圖5(a)),對于 3D 結果,迎角從 -2°~0°增大時,RGP 有所減小,但并不明顯;而從0°~2°增大時,RGP的值則明顯減小。這說明名義迎角的增大會使徑向壓力梯度得到加強,而CAM的結果中RGP隨著名義迎角的增大則呈現出相反的趨勢,對徑向壓力梯度變化的描述不足。在加入周向不均勻性的影響之后,CAM+MODEL計算得到的各名義迎角下RGP的分布都更接近3D結果,隨迎角的變化趨勢也要更接近3D結果。結合圖5(b)~(e)可見,在加入周向不均勻性的影響后,RGP的變化主要來源于CENT_M和Pr,并且Pr對于RGP的變化起主要作用。
對于CENT_W(圖5(b)),3D、CAM+MODEL和CAM之間的偏差反映出周向分速度的變化,即在靠近葉片通道進口處環量已經由于周向不均勻性的作用而發生改變。CAM+MODEL相比CAM偏離3D更多,說明加入周向不均勻性模型后通流計算得到的vu所受到的影響要大于3D數值模擬,但從圖中橫坐標的數值可見,CAM+MODEL與CAM和3D之間的實際區別并不大。而從趨勢上來看,CAM與3D之間確實存在偏差,而且在加入周向不均勻性模型后,CAM+MODEL的結果向著3D的方向移動,表現出了一定的一致性,說明在葉片通道進口vu確實會由于周向不均勻性的作用而發生改變;此外,當迎角從-2°~2°增大時,CAM與3D之間的偏差增大,反映出了周向不均勻性的作用在增強,與此同時,CAM與CAM+MODEL之間的偏差也隨著迎角的變化呈現出了相似變化,說明加入周向不均勻性模型后通流模型能夠與3D數值模擬一樣反映出負荷變化對葉片通道進口流場的影響。
CENT_M(圖 5(c))在名義迎角從 -2°~2°增大后,其值呈現出減小趨勢,說明迎角增大后子午流線曲率增大。而CAM對CENT_M的預測要高于3D的結果,趨勢也相反;在引入周向不均勻性的影響后,CAM+MODEL所預測的CENT_M的精度有所提高,可達60%以上,趨勢也與3D一致,說明加入周向不均勻性影響后通流模型能夠更好地刻畫子午流線曲率。
AC_M(圖5(d))在名義迎角發生改變時變化并不是很明顯,其值略有增加。在引入周向不均勻性的影響后,CAM+MODEL的分布相比CAM更接近3D。
從圖5(e)中可見,隨著名義迎角的增大,周向脈動源項的徑向分量Pr的絕對值增大,反映出負荷沿展向有增大的趨勢。此外,應力輸運模型計算得到的Pr的分布與3D結果基本一致,最大偏差在50%以內。
綜上所述,CAM+MODEL能夠得到與3D相接近的分布,說明應力輸運模型在迎角改變的情況下有能力反映出與3D計算相近的周向不均勻性的影響,提高通流計算的預測精度。同時,加入周向不均勻性后,通流模型能夠更好地反映出掠葉柵在不同迎角下流場中周向不均勻性的影響。
研究表明,周向不均勻性會對葉輪機流場造成影響,因此有必要在通流模型中計入周向不均勻性的影響,以使其更好地計入真實3D流動效應,提高設計、分析通流模型的預測精度。為模化并分析迎角改變時掠葉柵進口流動周向不均勻性的影響,本文將1種模化周向平均通流模型中周向脈動應力項的應力輸運模型結合到周向平均通流模型中,針對掠葉柵在不同迎角下的結果對比分析了該模型的可行性及周向不均勻性的作用,得出了以下結論:
(1)隨著迎角的增大,周向脈動源項的絕對值及其對葉片通道進口流場的影響都隨之提高,在本文算例中受周向不均勻性的影響葉片通道進口氣流角的改變量可達0.7°。
(2)在加入周向不均勻性模型后,在迎角改變的條件下通流計算預測的進口氣流角均更接近3D結果,對葉片通道進口氣流角的改變量的預測精度普遍提升30%以上,對葉片通道進口流動徑向平衡的描述精度提高60%以上。引入周向不均勻性有助于提升通流分析軟件對采用掠設計的壓氣機流場的計算精度。
(3)對于處于不同迎角下的葉柵流場,該應力輸運模型都能較好的預測出周向脈動源項,在前緣前模型預測結果與3D結果的相對偏差普遍在20%之內,最大不超過50%。
(4)由于應力輸運模型在主控方程中增加了6個微分方程,因此對通流模型的計算時間影響較大,加入該模型后計算時間在原通流模型的基礎上增加了約40%,但仍遠遠少于3D數值模擬所用時間。該模型將周向不均勻性與環量分布關聯,可以為設計或分析問題提供周向脈動應力項的預測。