亞音速軸流風(fēng)扇靜轉(zhuǎn)子相互作用線譜噪聲預(yù)報*

武星宇，魏應(yīng)三，靳栓寶，王 東，祝 昊，胡鵬飛，孫方旭

(海軍工程大學(xué) 艦船綜合電力技術(shù)國防科技重點試驗室， 湖北 武漢 430033)

軸流風(fēng)扇被廣泛應(yīng)用于民用和軍用航空領(lǐng)域，例如渦扇發(fā)動機，其低噪聲設(shè)計一直以來都受到國內(nèi)外許多學(xué)者的關(guān)注和研究[1]。從頻譜上看，風(fēng)扇氣動噪聲包含線譜噪聲和寬帶噪聲。其中，寬帶噪聲源于風(fēng)扇靜轉(zhuǎn)子葉片與隨機湍流脈動的相互作用，線譜噪聲主要源于周期性轉(zhuǎn)子尾流與下游靜子的相互作用[2-3]。對于大涵道比航空發(fā)動機，線譜噪聲為其主要噪聲源，因此深入研究轉(zhuǎn)子尾流與下游靜子相互作用產(chǎn)生的線譜噪聲對抑制風(fēng)扇噪聲具有重要意義。

隨著計算機技術(shù)的快速發(fā)展，數(shù)值分析法開始廣泛用于軸流機械非定常氣動特性和噪聲機理的研究，Verdon等[4]通過數(shù)值法求解三維線性歐拉方程，得到軸流渦輪機械靜子葉柵葉頻和倍葉頻下葉片表面非定常脈動壓力的分布。Elhadidi等[5]采用線性的歐拉方程和片條理論，計算分析了葉片彎掠設(shè)計對風(fēng)扇線譜噪聲的影響。為考慮流體黏性和靜子對轉(zhuǎn)子尾流的影響，Rumsey等[6]通過數(shù)值法求解非定常雷諾平均N-S方程，計算分析了風(fēng)扇靜轉(zhuǎn)子相互作用噪聲。Grace等[7]通過三維非定常雷諾平均N-S方程模擬計算靜子表面的非定常壓力，并據(jù)此預(yù)測葉柵出口線譜輻射聲功率。通過數(shù)值法可較準確地預(yù)報軸流風(fēng)扇的氣動噪聲[8-9]，但數(shù)值法建模時間和計算時間較長，在風(fēng)扇初期選型與參數(shù)設(shè)計時不適宜采用數(shù)值法。相對于數(shù)值計算方法，解析法在計算時間上具有很大的優(yōu)勢，適用于風(fēng)扇初期的低噪聲設(shè)計。

通過升力面理論計算流場與風(fēng)扇葉片相互作用產(chǎn)生的表面壓力脈動，結(jié)合聲類比方程可進一步計算得到風(fēng)扇的輻射噪聲。張偉光等[10]利用三維升力面理論計算轉(zhuǎn)子尾流與靜子葉片相互作用產(chǎn)生的線譜噪聲，并分析葉片彎掠設(shè)計對風(fēng)扇線譜噪聲的影響，但是升力面理論不能考慮風(fēng)扇靜子的安裝角[11]。為此，Carazo等[12]采用數(shù)值法模擬轉(zhuǎn)子尾流，結(jié)合Amiet[13]葉片非定常響應(yīng)理論模型，計算得到對轉(zhuǎn)螺旋槳線譜聲壓值，計算結(jié)果與試驗結(jié)果相近，但是該模型未考慮葉片之間的相互作用。Posson等[14]在Glegg[15]葉柵響應(yīng)函數(shù)的基礎(chǔ)上，推導(dǎo)得到葉片表面壓力脈動分布解析表達式。在此基礎(chǔ)上，de Laborderie等[16-17]推導(dǎo)出風(fēng)扇線譜噪聲計算公式，結(jié)合試驗得到的轉(zhuǎn)子葉片尾流參數(shù)，得到的風(fēng)扇輻射聲功率級計算值與試驗值在一倍葉頻處相近。此外，Rozenberg[18]、Sanjosé[19]、Magne[20]等也采用解析的葉柵非定常響應(yīng)模型計算風(fēng)扇線譜噪聲，但上述研究需要先計算葉片的非定常響應(yīng)，然后通過聲類比方程計算風(fēng)扇輻射噪聲，且需要數(shù)值法或試驗法提供的轉(zhuǎn)子尾流場參數(shù)。Philbrick等[21]試驗測量了葉片尾流參數(shù)，擬合曲線得到葉片經(jīng)驗尾流模型，并通過該模型預(yù)報風(fēng)扇線譜噪聲，發(fā)現(xiàn)采用試驗尾流數(shù)據(jù)和采用經(jīng)驗尾流模型得到的風(fēng)扇線譜噪聲結(jié)果相近。

本文將針對風(fēng)扇靜轉(zhuǎn)子葉柵相互作用線譜輻射噪聲進行研究。在文獻[15]中的簡諧湍流波作用下葉柵散射場速度勢解析表達式的基礎(chǔ)上，考慮風(fēng)扇轉(zhuǎn)子尾流的作用，通過傅里葉分解將轉(zhuǎn)子尾流分解為若干簡諧湍流波的疊加[21]，從而推導(dǎo)得到轉(zhuǎn)子尾流作用下靜子葉柵散射場的速度勢計算表達式，通過聲強與速度勢以及輻射聲功率與聲強的關(guān)系進一步推導(dǎo)得到了靜子葉柵線譜噪聲計算公式，并通過與NASA風(fēng)扇試驗?zāi)Ｐ秃秃舐?、傾斜靜子升力面計算模型對比驗證了該公式的適用范圍。并在此基礎(chǔ)上分析靜轉(zhuǎn)子間距、靜子葉片參數(shù)等對風(fēng)扇線譜噪聲的影響。

1 葉柵線譜噪聲公式推導(dǎo)

1.1 轉(zhuǎn)子葉柵尾流模型

忽略葉片厚度，在半徑r處截取風(fēng)扇靜轉(zhuǎn)子葉片并展開，可得到如圖1所示的靜轉(zhuǎn)子葉柵。圖1中：Ω為轉(zhuǎn)子的轉(zhuǎn)速，Sr為轉(zhuǎn)子葉柵柵距，αr為轉(zhuǎn)子葉柵安裝角，L為轉(zhuǎn)子尾流半尾跡寬度，DS為轉(zhuǎn)子后緣流向距離，(x,y)為固定在轉(zhuǎn)子葉柵上的坐標(biāo)系；αs為靜子葉柵的安裝角，(xs,ys)為固定在靜子葉柵上的坐標(biāo)系。

圖1 靜轉(zhuǎn)子葉柵示意Fig.1 Schematic diagram of stator and rotor cascades

大量試驗表明，單個葉片尾流可用高斯分布描述，在圖1中的o-xy坐標(biāo)系下其表達式[21-23]為

(1)

式中，w0(x,y)為葉片尾流場速度梯度分布，wr為轉(zhuǎn)子中心尾跡區(qū)軸向速度虧損。

對于轉(zhuǎn)子葉柵，轉(zhuǎn)子尾流場為單個葉片尾流場的疊加，即

(2)

式中，wR(x,y)為轉(zhuǎn)子尾流場速度梯度分布，表示的是轉(zhuǎn)子尾流場中某位置的流速與平均流速的差值，Z為整數(shù)。

以T=Srcosαr為周期將式(2)傅里葉級數(shù)展開，得

(3)

式中，xr·nr=ycosαr+xsinαr，F(xiàn)為xr·nr的函數(shù)。

(4)

(5)

由于式(4)中的半尾跡寬度L與DS有關(guān)，為此F可改寫為F′(DS)。 結(jié)合式(3)得到

Philbrick等[21]給出轉(zhuǎn)子中心尾跡區(qū)軸向速度虧損wr和半尾跡寬L經(jīng)驗計算公式。

(6)

(7)

其中：W0為轉(zhuǎn)子尾流平均速度；cr為轉(zhuǎn)子葉片弦長；Cd為轉(zhuǎn)子葉片阻尼系數(shù)，根據(jù)經(jīng)驗[21]有

(8)

式中，κ為葉片稠度，α1為進口相對氣流角，α2為出口相對氣流角，λ為壓力損失系數(shù)。

考慮轉(zhuǎn)子與靜子坐標(biāo)變換，在y方向上有

ys=y+Ωrt

(9)

得到周期性轉(zhuǎn)子尾流垂直于靜子葉柵的速度分量為

(10)

式中，B為轉(zhuǎn)子葉片數(shù)。

1.2 葉柵線譜輻射聲功率

Glegg[15]針對圖2所示葉柵模型提出葉柵響應(yīng)函數(shù)理論，該理論通過Wiener-Holf方法求解葉柵散射場速度勢函數(shù)積分表達式，最后得到簡諧湍流脈動作用下葉柵散射場的速度勢函數(shù)解析表達式，推導(dǎo)得到的表達式可由式(11)表示。

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

ikzzd-iωt]

(11)

通過片條理論，在半徑r處截取環(huán)形風(fēng)扇中Δr段，假設(shè)在分段內(nèi)的流體入流參數(shù)一致[16]，考慮葉柵入流為轉(zhuǎn)子尾流，葉柵的速度勢函數(shù)解析表達式可改寫為

i(σ-2πq)yd/h+ikzzd-imBΩt]

(12)

聲強與速度勢的關(guān)系[15]為

(13)

式中，“*”表示共軛，I為聲強。

代入速度勢函數(shù)，可得到

(14)

輻射聲功率Π與聲強的關(guān)系[15]為

(15)

對于圖1所示的葉柵結(jié)構(gòu)，其積分曲面定義在曲線x-yd/h=const和r

(16)

式中，Re表示取實部，這是由于輻射聲功率為實數(shù)時才能輻射出來。 通過式(16)可得到靜子第m階輻射聲功率計算表達式，即

(17)

相鄰葉片間湍流相位差角σ=kxdd+kydh，其中kxd為Glegg坐標(biāo)系下的軸向波數(shù)，kyd為側(cè)向波數(shù)。不考慮葉片后掠和傾斜，波數(shù)在Glegg坐標(biāo)系od-xdyd與坐標(biāo)系os-xsys之間的轉(zhuǎn)換關(guān)系[16]為

(18)

1.3 葉片后掠和傾斜的影響

圖3為后掠和傾斜葉片示意圖。圖3(a)中φ為后掠角，可表示為沿葉片徑向的函數(shù)φ(z)。圖3(b)中ψ為后掠角，可表示為沿葉片徑向的函數(shù)ψ(z)。此時固定在靜子葉片前緣的坐標(biāo)系os-xsyszs與Glegg坐標(biāo)系od-xdydzd之間的轉(zhuǎn)換關(guān)系如式(19)所示。

(a) 后掠(a) Swept (b) 傾斜(b) Lean圖3 后掠與傾斜葉片示意Fig.3 Schematic diagram of swept and lean blade

(xd,yd,zd)T=Q×(xs,ys,zs)T

(19)

式中，Q為轉(zhuǎn)換矩陣，其定義[16]為

(20)

通過Q可得到坐標(biāo)系od-xdydzd下的葉柵入流參數(shù)。

(21)

其中，[U0,Uyd,Uzd]T為od-xdydzd坐標(biāo)下的入流速度，[Uxs,Uys,Uzs]T為os-xsyszs坐標(biāo)下的入流速度，kzd為od-xdydzd坐標(biāo)下z向入流波數(shù)，kzs為os-xsyszs坐標(biāo)下z向入流波數(shù)。靜子葉柵來流波數(shù)K=[kxs,kys,kzs]T與流速U=[Uxs,Uys,Uzs]T滿足如式(22)所示關(guān)系[16]。

K·U=mBΩ

(22)

考慮葉片在zs方向上波數(shù)的變化，在坐標(biāo)系os-xsyszs下的波數(shù)為

(23)

其中，p為整數(shù)，Rt為葉梢半徑，Rh為輪轂半徑。通過式(21)可得到Glegg定義坐標(biāo)系下的葉柵入流參數(shù)。此外，由于傾斜和后掠角的影響，靜子在Glegg坐標(biāo)系下的弦長c、垂向間距h和軸向間距d也發(fā)生了變化。

(24)

其中，cs為os-xsyszs坐標(biāo)系下的弦長，Ss為靜子葉柵柵距。

2 結(jié)果和討論

2.1 不考慮葉片后掠與傾斜

NASA[22]在半消聲室中測量了不同轉(zhuǎn)速和靜子葉片數(shù)下風(fēng)扇的線譜噪聲，其試驗?zāi)Ｐ腿鐖D4所示，主要參數(shù)如表1所示。試驗測量了轉(zhuǎn)子轉(zhuǎn)速分別為1 700 r/min，1 750 r/min，1 800 r/min，1 850 r/min和1 887 r/min時的風(fēng)扇線譜噪聲，在轉(zhuǎn)速為1 800 r/min時，轉(zhuǎn)子葉片葉梢馬赫數(shù)為0.34，軸向馬赫數(shù)為0.15。

圖4 NASA試驗風(fēng)扇結(jié)構(gòu)[22]Fig. 4 NASA′s fan test structure[22]

表1 NASA風(fēng)扇試驗參數(shù)

Sutliff等[24]給出了轉(zhuǎn)子流場參數(shù)試驗測量結(jié)果，通過轉(zhuǎn)子入流流速，圖5給出了轉(zhuǎn)速為1 850 r/min時的轉(zhuǎn)子尾流場馬赫數(shù)梯度經(jīng)驗公式計算云圖。通過經(jīng)驗公式計算得到，在葉片中部轉(zhuǎn)子尾流中心線速度虧損約為18 m/s，Sutliff試驗測量的結(jié)果約為20 m/s，經(jīng)驗結(jié)果與試驗結(jié)果相近。

圖5 周期性轉(zhuǎn)子經(jīng)驗尾流公式計算云圖Fig.5 Cloud diagram of periodic rotor wake by empirical

假設(shè)導(dǎo)管壁面為硬壁面，即不考慮導(dǎo)管的吸聲作用，通過式(17)即可得到風(fēng)扇靜子的線譜輻射聲功率。圖6、圖7分別給出13葉靜子和14葉靜子的線譜輻射聲功率級(參考聲壓級為10-12W)試驗值與計算值對比曲線圖，計算值包括本文公式計算值與Laborderie[16]模型計算值，式(18)中的m=1時可得到葉頻處的輻射聲功率，m=2時可得2倍葉頻處的輻射聲功率。圖6(a)與圖6(b)分別為葉頻處風(fēng)扇上游和下游線譜輻射聲功率級，圖6(c)與圖6(d)分別為風(fēng)扇上游和下游2倍葉頻處線譜輻射聲功率級。

(a) 風(fēng)扇上游輻射聲功率級(葉頻)(a) Sound power level of the fan′s upstream in BPF

通過圖6可得，當(dāng)靜子葉片數(shù)為13、轉(zhuǎn)子轉(zhuǎn)速為1 700～1 887 r/min范圍時，在葉頻處本文公式計算得到的風(fēng)扇線譜輻射聲功率級與試驗值相差2～4 dB，特別是在葉頻下游工況，本文公式計算結(jié)果與試驗值相近。

通過圖7對比14葉靜子的線譜輻射聲功率級可得，在葉頻下游本文公式計算得到的輻射聲功率級與試驗值相差3 dB以內(nèi)；但是在2倍葉頻處，計算值與試驗值相差較大，接近10 dB。總體而言，在上述計算工況，相比于Laborderie模型，本文公式在葉頻下游處的預(yù)報結(jié)果與試驗值更加接近，但是在2倍的葉頻處，本文的計算結(jié)果偏大，這是由于使用經(jīng)驗?zāi)Ｐ皖A(yù)報轉(zhuǎn)子尾流參數(shù)時時未考慮轉(zhuǎn)子尾流沿靜子弦長的變化，導(dǎo)致傅里葉分解后2倍葉頻處的速度梯度偏大；而在葉頻處，本身速度梯度基數(shù)較大，弦長對轉(zhuǎn)子尾流速度梯度的影響較小。

(a) 風(fēng)扇上游輻射聲功率級(葉頻)(a) Sound power level of the fan′s upstream in BPF

2.2 考慮葉片后掠

Schulten[25]采用升力面法分析了后掠葉片對風(fēng)扇線譜噪聲的影響，模型參數(shù)如表2所示。

表2 Schulten后掠靜子風(fēng)扇模型主要參數(shù)

表2中，歸一化葉頻為轉(zhuǎn)子葉梢馬赫數(shù)與轉(zhuǎn)子葉片數(shù)的乘積。此外，轉(zhuǎn)子阻尼系數(shù)為0.01，靜轉(zhuǎn)子間距在葉中位置保持2.6倍的靜子弦長。在采用經(jīng)驗公式計算轉(zhuǎn)子尾流參數(shù)時，考慮后掠角導(dǎo)致的靜轉(zhuǎn)子間距的變化。圖8給出靜子的線譜輻射聲功率級本文公式計算值與Schulten升力面理論計算值對比曲線圖，由于風(fēng)扇葉頻噪聲被截止，因此圖中計算頻率為2倍的葉片通過頻率。

由圖8可得，在當(dāng)后掠角小于8°時，本文公式計算得到的輻射聲功率與Schulten升力面計算模型計算得到的輻射聲功率相近，但是當(dāng)后掠角大于8°時，本文公式計算得到的輻射聲功率快速下降，與Schulten計算模型得到的結(jié)果相差較大。這是由于本文公式基于平面葉柵推導(dǎo)而來，在大后掠角下，軸向波數(shù)變小，葉柵部分模態(tài)的輻射噪聲被“截止”，但對于環(huán)形葉柵，在平面葉柵中被“截止”的輻射模態(tài)不一定會在環(huán)形葉柵中被“截止”，這導(dǎo)致本文計算得到的結(jié)果與Schulten計算得到的結(jié)果相差較大。

(a) 靜子上游輻射聲功率(a) Sound power level of the stator′s upstream

2.3 考慮葉片傾斜

張偉光等[10]采用升力面法分析了傾斜角對風(fēng)扇線譜噪聲的影響，計算參數(shù)如表3所示，其中弦長與靜轉(zhuǎn)子間距參數(shù)通過導(dǎo)管直徑做歸一化處理，導(dǎo)管直徑為0.4 m。

表3 Schulten傾斜靜子風(fēng)扇模型主要參數(shù)

圖9給出靜子線譜輻射聲功率級本文公式計算值與文獻[10]升力面理論計算值對比曲線圖，計算頻率為1倍的葉片通過頻率，其中縱坐標(biāo)為計算結(jié)果與上游0°傾斜角的輻射聲級計算結(jié)果相對值。

(a) 靜子上游輻射聲功率(a) Sound power level of the stator′s upstream

由圖9可得，當(dāng)傾斜角在0°～10°范圍內(nèi)變化時，本文公式計算值與升力面計算值相差5 dB以內(nèi)，且兩者發(fā)展趨勢一致。對于上述模型，葉片的傾斜設(shè)計反而使靜子的線譜噪聲變大。

3 參數(shù)影響分析

在2.1節(jié)NASA 13葉靜子風(fēng)扇試驗?zāi)Ｐ偷幕A(chǔ)上，圖10給出了靜轉(zhuǎn)子間距、靜子弦長、靜子傾斜角和靜子葉片數(shù)等參數(shù)變化對靜子輻射噪聲的影響。

(a) 靜轉(zhuǎn)子間距的影響(a) Influence of the stator-rotor space

考慮到下游葉頻噪聲為主要噪聲貢獻量，圖10計算的均是下游葉頻噪聲，圖10(a)為靜轉(zhuǎn)子間距變化對靜子線譜噪聲的影響曲線圖，由圖可得，靜轉(zhuǎn)子間距對靜子線譜噪聲影響較大，隨著靜轉(zhuǎn)子間距增大，風(fēng)扇線譜噪聲降低，但同時降低的幅值越來越小。

圖10(b)為弦長變化對靜子線譜噪聲的影響曲線圖，由圖可得，弦長變化對NASA試驗?zāi)Ｐ偷木€譜噪聲影響較小。

圖10(c)為傾斜角變化對靜子線譜噪聲的影響曲線圖，由圖可得，當(dāng)傾斜角小于4°時，風(fēng)扇線譜噪聲變化較小，當(dāng)傾斜角大于4°時，傾斜角越大，風(fēng)扇線譜噪聲越小。

圖10(d)為葉片數(shù)變化對靜子線譜噪聲的影響圖，由圖可得，15葉靜子的線譜噪聲最小，其次是17靜子，20葉靜子的線譜噪聲最大。此外，在圖10(d)中16葉靜子由于與轉(zhuǎn)子葉片數(shù)相等，導(dǎo)致計算中出現(xiàn)奇點；當(dāng)葉片數(shù)大于21時，葉頻噪聲被“截止”。

4 結(jié)論

本文在葉柵響應(yīng)函數(shù)的基礎(chǔ)上推導(dǎo)了風(fēng)扇靜轉(zhuǎn)子相互作用線譜輻射聲功率預(yù)報公式，通過轉(zhuǎn)子經(jīng)驗尾流模型為該公式提供入流參數(shù)，并考慮靜子葉片后掠和傾斜設(shè)計對線譜噪聲的影響。通過與NASA試驗?zāi)Ｐ秃蜕γ嬗嬎隳Ｐ蛯Ρ?，驗證了本文公式適用性和適用范圍，得到如下結(jié)論：

1)本文公式能夠較準確地預(yù)報風(fēng)扇線譜噪聲，但是本文公式不適用于后掠葉片線譜噪聲的計算，這是因為平面葉柵中被“截止”的聲輻射模態(tài)在環(huán)形葉柵中不一定會被“截止”。

2)葉片的傾斜設(shè)計并非都能抑制線譜噪聲，對于部分模型葉片的傾斜設(shè)計反而會增大線譜噪聲。

3)針對NASA風(fēng)扇試驗?zāi)Ｐ?，增大靜轉(zhuǎn)子間距和傾斜設(shè)計能夠抑制線譜噪聲，但是靜子葉片弦長對線譜噪聲影響很?。淮送猓~片數(shù)對風(fēng)扇的線譜噪聲影響較大，特別是可通過增加葉片數(shù)的方式來“截止”靜子葉頻噪聲。因此，合理的靜轉(zhuǎn)子葉片數(shù)匹配設(shè)計對風(fēng)扇噪聲被動抑制具有重要意義。

綜上所述，本文解析公式可預(yù)報風(fēng)扇靜轉(zhuǎn)子相互作用線譜噪聲，相對于常用的數(shù)值法，在計算時間上具有一定的優(yōu)勢，適用于風(fēng)扇葉片的初期選型和低噪聲設(shè)計。與升力面理論相比，本文公式直接通過速度勢函數(shù)計算得到風(fēng)扇線譜噪聲，無須先求解葉片表面的非定常壓力脈動，其物理意義更加清晰，計算過程更加簡單。此外，升力面法無法考慮安裝角的影響。本文公式則可考慮安裝角和傾斜角對風(fēng)扇線譜噪聲的影響，但是也存在一定的缺陷，對于后掠角過大的模型其預(yù)報結(jié)果不是很準確，還需要進一步優(yōu)化本文公式。

致謝

清華大學(xué)黃振衛(wèi)博士和海軍潛艇學(xué)院段嘉希博士在數(shù)據(jù)處理提供了幫助和指導(dǎo)，謹致謝意！