基于有限元法離心風機氣動噪聲預報研究

閆 瑩，肖友洪，盧華兵，許順利

（哈爾濱工程大學 動力與能源工程學院，哈爾濱 150001）

0 引言

近年來離心風機等旋轉機械的噪聲問題越發(fā)引人注目，隨著計算機技術的發(fā)展，氣動噪聲的問題亟待解決。離心風機在運轉時的噪聲主要分為氣動噪聲和振動噪聲，對于轉速不高的離心風機來說，其氣動噪聲是主要噪聲。氣動噪聲又分為離散噪聲和寬頻噪聲［1］。其中離散噪聲又稱為葉片通過頻率（BPF）噪聲，是葉片通過帶動周圍氣體波動從而引發(fā)的噪聲［2-4］。寬頻噪聲則是在流動的固體表面的湍流引起渦流而產(chǎn)生的噪聲。通常情況下離散噪聲的峰值更加明顯，因此對于風機的氣動噪聲的研究主要是在離散噪聲上。

對于氣動噪聲研究過程，由于聲音和流動波動都是N-S方程的計算結果，因此氣動聲學有兩種不同的解法。其中直接數(shù)值模擬求解N-S方程，通過捕捉各個尺度的渦，從而精準地得到聲音的產(chǎn)生和傳播。這種求解方式通常要求高階、低耗散的方法用于精確捕捉和傳播聲波。然而其計算成本大，不適宜解決工程問題。對于工程問題，采用Lighthill聲類比方法更為普遍。

目前對Lighthill的求解方法主要是積分方法如Curle和FW-H方程［5］。外流場問題可以利用自由空間格林函數(shù)求解。而對于幾何邊界復雜不對稱的內(nèi)流問題，難以用格林函數(shù)求解，對于遠場聲輻射需要做進一步的簡化，導致計算誤差較大。另一種方法是基于有限元求解，采用Lighthill變分形式，將聲源轉化為體聲源和面聲源項。這種方法將聲源在空間上進行離散，考慮了聲與結構的相互作用，對于幾何復雜的問題有一定優(yōu)勢。

隨著空氣動力學，計算機和氣動聲學的發(fā)展，對氣動噪聲數(shù)值研究逐漸增多。LIN等［6］對前向離心式風機進行研究，采用數(shù)值和試驗相結合的方式，得到試驗值與數(shù)值模擬的誤差小于5.4%，證明了數(shù)值模擬對于離心風機內(nèi)流模擬的可行性。羅凱等［7］針對多翼離心風機進行風機流噪聲計算，得到風機的蝸舌半徑、葉片數(shù)、葉片進口安放角對風機的風量和流噪聲均有影響的結論。趙征等［8］以葉片數(shù)作為變量，細致研究了其對高比轉速離心風機氣動性能和噪聲的影響，得到了增加葉片數(shù)可以提高氣動性能降低噪聲的結論，并計算得到最佳葉片數(shù)。柳琦等［9］對某車用交流發(fā)電機進行噪聲預測，采用大渦模擬和Lighthill聲類比理論，考慮了聲傳播對氣動噪聲的影響。胡四兵等［10］針對小型高速離心風機進行試驗研究得到輻射噪聲為風機主要噪聲類型?？祻姷龋?1］采用有限元方法對車用旋渦風機進行氣動噪聲數(shù)值模擬，對遠場氣動噪聲進行計算，驗證有限元法計算氣動噪聲的準確性，為降噪提供依據(jù)。

本文首先對離心風機進行噪聲試驗，然后對三維流場數(shù)值模擬，將計算得到的流場數(shù)據(jù)作為聲源信息插值到聲學網(wǎng)格上，通過有限元方法計算得到風機空間聲場分布曲線，并與試驗值進行對比，從而驗證ACTRAN有限元法預測氣動噪聲的準確性。根據(jù)上述仿真流程，進一步探究葉輪轉速對聲壓級的影響及原因，以驗證仿真實用性。

1 離心風機試驗研究

此研究的離心風機葉輪為后向葉輪，試驗開啟轉速為1 788 r/min。測點布置在風機蝸殼出口中心平面左、右兩側以邊長為1m的正方形的對角線長度45°方向上，具體測點布置見圖1。

圖1 試驗及測點布置Fig.1 Layout of experimental and measuring points

采用型號為B&K4957麥克風，及相關的數(shù)據(jù)采集軟件采集信號。試驗時自然進氣，在確認離心風機的運轉已經(jīng)平穩(wěn)時采集時域數(shù)據(jù)，經(jīng)快速傅里葉變換（FFT）得到2 000 Hz以下的不同頻率的聲壓級。從圖2測試結果可以看到，風機基頻和二階諧頻的離散噪聲明顯，其基頻的聲壓級為89 dB。左、右兩側的噪聲頻譜存在差異，主要是離心風機分布的不對稱性，在進口段附近也存在噪聲，因此測點2噪聲頻譜略高于測點1。

圖2 試驗監(jiān)測點聲壓頻譜Fig.2 Spectrum diagram of sound pressure at experimental monitoring points

2 離心風機數(shù)值模擬

2.1 流場模擬

試驗研究的離心風機的三維結構如圖3所示。其主要部分包括葉輪、蝸殼、集流器，相應的結構參數(shù)見表1。

圖3 三維結構Fig.3 Three-dimensional structure diagram

表1 離心風機結構參數(shù)Tab.1 Structural parameters of centrifugal fan

當葉輪旋轉時，流體從進風口吸入進入集流器導流到葉輪內(nèi)，風向由軸向變成徑向，氣體由蝸殼收集，最終排出風機完成工作。為了保證在風機仿真中流場更容易收斂，增加在進出口各延長一段，減少風機的回流現(xiàn)象［12］。對上述的風機結構抽取三維流域并作適當簡化得到仿真計算的流域，如圖4所示。

圖4 各部分流域Fig.4 Fluid domain diagram of each part

為了避免仿真中離心風機的回流現(xiàn)象，網(wǎng)格劃分過程分成5個部分，分別為進口延長段、集流器、葉輪、蝸殼和出口延長段，對尺寸較小的區(qū)域進行加密，各部分的網(wǎng)格數(shù)量見表2。

表2 各部分網(wǎng)格數(shù)Tab.2 Grid number of each part

網(wǎng)格采用多面體-六面體核心的非結構化網(wǎng)格。另外，考慮聲學分辨率的要求，最終總網(wǎng)格數(shù)達到730萬左右。其網(wǎng)格劃分如圖5所示。

圖5 網(wǎng)格細節(jié)Fig.5 Grid detail diagram

流場計算定常部分采用Realizable k-ε湍流模型，壓力-速度耦合采用SIMPLE方法，初始邊界條件設置為總壓入口，靜壓出口。模型采用了無滑移壁面多重參考系MRF模型。將穩(wěn)態(tài)結果作為瞬態(tài)的初值，采用滑移網(wǎng)格（SMM）模型，模擬真實的葉輪旋轉。關于動靜交界的部分，設置了interface面，將蝸殼與葉輪區(qū)域和葉輪與進口區(qū)域進行動靜結合，實現(xiàn)流場數(shù)據(jù)的交互。

計算中風機轉速設置為1 788 r/min，為滿足聲學和流體分辨率，設置時間步長為Δt=5×10-5s，最大迭代步數(shù)為25步，則根據(jù)Δt=60/（nz）可知，葉輪旋轉1r所需要的步數(shù)為672步。由于開始時葉輪旋轉不具有穩(wěn)定性，在瞬態(tài)計算時先讓葉輪旋轉3～6周，當監(jiān)測點呈現(xiàn)出周期性時認為瞬態(tài)收斂開始輸出計算結果。根據(jù)奈奎斯特定理，最大頻率為10 000 Hz。本文在葉輪旋轉5周監(jiān)測點呈現(xiàn)周期性，因此從葉輪第6周開始輸出聲源信息用于聲學計算。

建立離心風機氣動聲學的模型如圖6所示。由于試驗進口布置管道，因此在聲學建模時將進口部分建立管道模態(tài)。將葉輪和蝸殼的交界面作為面聲源，將旋轉域的聲源結果插值到面聲源上，體聲源則用蝸殼靜域來表示。在蝸殼出口處，建立直徑為1 m的半圓區(qū)域，來模擬蝸殼出口的聲傳播的有限元區(qū)域。半圓表面定義為無限元組件，用來模擬無反射自由空間的邊界條件。

圖6 聲學模型Fig.6 Diagram of acoustic model

聲學網(wǎng)格如圖7所示，為了滿足ACTRAN計算精度的要求，通常假設1個波長內(nèi)至少有6個網(wǎng)格［13-15］。網(wǎng)格單元的長度應滿足式（1）。

式中 L ——聲學網(wǎng)格長度，mm

c ——聲音的傳播速度，m/s；

fm——聲學計算最大頻率，Hz。

本文所關注的最大頻率在2 000 Hz左右，因而需要的最大網(wǎng)格尺度需要小于28.3 mm。則最終的聲學網(wǎng)格數(shù)在230萬左右。

3 結果與討論

3.1 流場分析

圖8，9分別示出z=0和x=0時的壓力云圖。流體經(jīng)過葉輪后從徑向流出到蝸殼區(qū)域，壓力是逐漸增高的，主要是葉輪的機械能轉換流體的動能。流體流經(jīng)蝸殼后，其壓力逐步減小，蝸殼將流體的壓力勢能轉化為動能。從蝸殼進口到出口，流體的壓力逐漸減小，且在蝸舌的附近壓力達到了最低值，根據(jù)伯努利方程，其速度最高。從x=0中的云圖可以看到，由于蝸殼的結構不對稱，蝸殼右側的壓力高于左側。集流器到葉輪部分先縮小再擴張的結構使得葉輪和集流器的交界處出現(xiàn)壓力波動。

圖8 z=0時壓力云圖Fig.8 Pressure nephogram at z=0

圖9 x=0時的壓力云圖Fig.9 Pressure nephogram at x=0

葉片吸力面是流體由于壓力減小而沖擊葉片的面，壓力面是葉輪施加于流體的面。圖10示出葉片表面壓力云圖，可以看到，葉片外部為壓力面，內(nèi)部為吸力面。可以發(fā)現(xiàn)，在徑向方向上，靜壓值隨著徑向距離的增大而增大。葉片尾緣靜壓值與葉片前緣相比更大，葉片壓力側靜壓值也比葉片吸力側的要大。葉尾根部的靜壓值明顯偏大，從等值線上看葉片表面的壓力分布并不均勻。主要是由于葉片在旋轉過程中，是動能轉化為壓力能的過程。而蝸殼壓力分布不均是因為其不對稱性，從圖11中可以看出，在蝸殼上部出現(xiàn)了渦，存在二次流，流動不均勻。因此葉片上部的壓力大于下部。

圖10 葉片表面壓力云圖Fig.10 Pressure nephogram of blade surface

圖11 x=0的速度矢量Fig.11 The velocity vector at x=0

在流場計算的過程中，會通過設置一些特定位置的監(jiān)測點來反映葉輪旋轉過程中的壓力脈動情況，故在蝸殼和葉輪附近設置監(jiān)測點，位置如圖12所示。

圖12 監(jiān)測點位置示意Fig.12 Schematic diagram of locations of monitoring points

圖13示出蝸殼各點的壓力脈動。從波動的周期性來看，各點壓力脈動都呈現(xiàn)出比較規(guī)律的周期性，體現(xiàn)了瞬態(tài)求解已經(jīng)達到了穩(wěn)定。另外從壓力波動的幅度圖來看，蝸舌處的監(jiān)測點w1的壓力脈動最為劇烈，這是因為在蝸舌的附近的渦流較大，所以蝸舌部分成為風機的主要噪聲源。另外，從幅值占比圖14來看監(jiān)測點w1的壓力波動遠大于 w2，w3，w4的壓力脈動。

圖13 監(jiān)測點的壓力脈動Fig.13 Pressure pulsation at the monitoring points

圖14 各監(jiān)測點壓力幅值占比Fig.14 Ratio of pressure amplitude at the monitoring points

因此，蝸殼上的尾跡減弱是因為w2，w3，w4距離蝸舌較遠，而且葉輪出口和蝸殼之間的徑向間隙隨著距離蝸舌的增加而增大。因此3個監(jiān)測點呈現(xiàn)較小的壓力脈動，壓力幅值也逐漸減小，實際上，蝸殼的壓力波動和近壁面的流動區(qū)域網(wǎng)格數(shù)以及所選湍流模型密切相關，因此采用精度更高的湍流模型，如大渦模擬等會達到更好的效果。從壓力脈動來看，主要噪聲源是葉輪出口及蝸舌附近。為了驗證這一點，對穩(wěn)態(tài)部分進行了寬頻噪聲的仿真，使用寬頻噪聲源模型來提取關于噪聲源的有用診斷信息，以確定流體的哪一部分對噪聲的貢獻最大。

圖15～17分別示出風機轉速為1 788 r/min的寬頻噪聲云圖，從圖中可以看出，蝸殼表面蝸舌附近的聲功率級較大，且其最大部分靠近葉輪區(qū)域。而由于葉輪流道出口壓力較大，氣流分布比較紊亂無序。高速轉動的葉輪甩出的氣流直接撞擊蝸舌及其附近區(qū)域，致使蝸舌表面附近的氣體密度變化很大，產(chǎn)生幅度較大的壓力波動，這種壓力波動引發(fā)出離心風機主要的噪聲。本文中的離心風機的葉輪葉片后緣離蝸舌很近，被葉片甩出的氣流來不及減速撞擊蝸舌，更加重了蝸舌附近的壓力波動。葉輪上噪聲值比較大的點集中在葉片后部，前蓋后部，以及靠近蝸舌附近的出口。在這些位置壓力波動劇烈，而聲音是由波動產(chǎn)生的，所以這部分噪聲值的貢獻也較大。

圖15 蝸殼表面聲功率級云圖Fig.15 Nephogram of sound power level on volute surface

圖16 葉輪表面聲功率級云圖Fig.16 Nephogram of acoustic power level on impeller surface

圖17 z=0寬頻噪聲Fig.17 Broadband noise graph z=0

3.2 噪聲分析

從圖18中可以看出，仿真對于基頻和諧頻噪聲的捕捉也較為明顯，其中基頻噪聲監(jiān)測點1試驗值達到89 dB，而仿真值在86.69 dB左右，相差2.31 dB，基頻仿真誤差在4%左右，誤差較小，且仿真與試驗值寬頻噪聲的趨勢也比較對應，仿真與試驗的總聲壓級分別為94.38 dB和93.804 dB，誤差為0.6%，仿真誤差在容許范圍內(nèi)，可以用作噪聲的預測，減小試驗成本。監(jiān)測點2的仿真值在86 dB左右，而試驗值在95 dB左右，試驗與仿真的誤差相對較大。其原因主要是試驗存在進口噪聲，且測點2 附近的空間相對狹窄，因此測點2的試驗值是略高于測點1的，與仿真值存在差異。

圖18 試驗與仿真聲壓頻譜對比Fig.18 Comparison of experimental and simulated sound pressure spectra

從聲壓云圖19來看，隨著頻率的升高，倍頻的最大聲壓值逐漸降低，這與試驗測點的聲壓頻譜呈現(xiàn)的規(guī)律基本一致。在云圖中聲壓最大的地方主要集中在蝸舌附近區(qū)域，這與之前噪聲源產(chǎn)生位置的分析是相吻合的。

圖19 1-3階頻率的聲壓云圖Fig.19 Sound pressure nephogram of frequencies of orders 1～3

關于仿真的噪聲與試驗值存在的差異，主要是由于在仿真中出口設置了無反射邊界，而實際試驗并不是在消聲室完成，試驗空間有限存在著聲波的反射。另外，在機械運行過程中還伴隨著其他噪聲，均會增加試驗與仿真誤差。此外，仿真的總周期用64核處理器在一周內(nèi)基本完成，對于試驗周期來說，從模型成型到反復調(diào)試得到最終聲壓結果需要2個月左右的時間。因此仿真預測工作可以大幅降低試驗加工，運輸，調(diào)試的時間和成本，縮短周期。

3.3 進一步討論

根據(jù)上述分析，可以驗證ACTRAN在噪聲預測中的實用性。因此進行進一步探究葉輪轉速對聲壓級的影響及原因。改變?nèi)~輪轉速1 200 r/min，2 400 r/min，對比不同轉速下聲壓級的變化，最終得到不同轉速下的聲壓頻譜。

從圖20中可以看出，在不同轉速下，其離散噪聲變化較為明顯，寬頻噪聲并沒有大幅度的變化。葉輪轉速對離散噪聲的影響比寬頻噪聲大得多。其主要原因和2種噪聲不同的產(chǎn)生機理有關。

圖20 不同轉速下的聲壓級曲線Fig.20 SPL curves at different speeds

不同轉速下的基頻噪聲值，其基頻噪聲隨著轉速升高提升了19.29 dB，見表3。轉速升高，葉輪周期性拍打流體而產(chǎn)生的壓力脈動越明顯，因此其離散噪聲也越高。

表3 不同轉速下的基頻噪聲Tab.3 Fundamental frequency noise values at different speeds

4 結論

（1）采用Realizable k-ε模型對離心風機非定常流場進行模擬，并結合壓力云圖和壓力波動，以及寬頻噪聲模型，判斷得出主要噪聲源位置在蝸舌和葉輪處，為進一步進行遠場噪聲求解奠定基礎。

（2）采用有限元法對遠場聲輻射進行求解并與試驗結果進行對比。結果表明試驗值與仿真值的離散噪聲誤差在4%左右，仿真值與試驗值寬頻噪聲趨勢基本一致，其總聲壓級的誤差為0.6%，證明了有限元方法對于離心風機氣動噪聲預測的可行性。且從仿真預測時間與試驗周期來看，仿真過程減少了試驗模型反復鑄造以及安裝調(diào)試的重復時間和成本，可以一定程度上替代實驗模型反復調(diào)試工作時間和成本。

（3）對不同轉速下離心風機遠場噪聲的變化情況進行討論。結果得到，轉速從1 200 r/min升高到2 400 r/min，引起基頻噪聲升高了19.29 dB。葉輪轉速對離散噪聲的影響比寬頻噪聲大得多。轉速越高，引起離散噪聲的變化越明顯，合理控制轉速范圍可以從源頭抑制噪聲產(chǎn)生。

（4）有限元方法雖然克服了格林函數(shù)對于復雜邊界的求解問題，但其假設條件仍然沒有考慮到聲音對于流場的耦合作用，研究可以朝流聲耦合的仿真方向進行，以提高仿真的準確性。