李亞，劉忠族 ，許影博 ，張楠

1船舶振動噪聲重點實驗室，江蘇無錫214082

2中國船舶科學研究中心，江蘇無錫214082

3中國船舶科學研究中心水動力學重點實驗室，江蘇無錫214082

0 引 言

離心風機管道中的噪聲一般來源于葉輪輻射聲、蝸殼和管道壁面邊界層噪聲、尾流場的湍流噪聲、蝸殼和管道振動引起的噪聲。噪聲在離心風機管道中傳播時，還會受到管道噪聲模態(tài)的影響。另外，離心風機運行時，電機噪聲和環(huán)境噪聲也會傳入管道。葉輪輻射聲和壁面邊界層噪聲屬于偶極子類型，輻射效率高于四極子類型的湍流噪聲，離心風機振動噪聲相比氣動聲小很多［1］。

將計算流體動力學（CFD）與計算聲學相結(jié)合，可以進行氣動噪聲預(yù)報。Kato等［2］運用大渦模擬（Large Eddy Simulation，LES）方法得到葉片表面的脈動壓力，采用Lighthill聲類比理論計算遠場輻射聲。Younsi等［3］采用滑移網(wǎng)格方式及SST k-ω湍流模型計算了離心風機內(nèi)部流場，利用非定常流場參數(shù)作為FW-H方程［4］輸入項預(yù)報了輻射聲。陳敏等［5］運用多塊混合網(wǎng)格和動網(wǎng)格技術(shù)對螺旋槳的水動力性能進行大渦模擬，并結(jié)合FW-H方程預(yù)報了螺旋槳的低頻流噪聲。在上述噪聲預(yù)報中，采用的是自由空間格林函數(shù)，其必須假設(shè)離心風機的幾何尺寸滿足聲學緊致結(jié)構(gòu)特征的要求。實際上，由于風機蝸殼和葉輪壁面的存在，機殼內(nèi)部聲源激發(fā)的噪聲會在葉輪和蝸殼內(nèi)壁面上發(fā)生多次反射及散射，最終通過風機進、出口向外輻射。

早在 1974 年，Moreland［6］就通過實驗證實了蝸殼腔體結(jié)構(gòu)對聲波反射及散射會導(dǎo)致部分頻段向外輻射的聲功率大大增強。劉厚林等［7］采用大渦模擬方法獲得了離心泵內(nèi)部瞬態(tài)流場，在計算葉片表面偶極子輻射聲場時，將葉片分為10個部分，分別積分得到3個方向的時域力，以此作為葉片旋轉(zhuǎn)偶極子聲源；以進、出口為吸聲邊界，其他壁面作為全反射壁面條件，聲場計算結(jié)果與實驗結(jié)果吻合良好。

離心風機形狀復(fù)雜，過去流場模擬多采用非結(jié)構(gòu)化網(wǎng)格，本文擬采用ICEM-CFD軟件繪制結(jié)構(gòu)化網(wǎng)格，使用Fluent軟件模擬流場，通過得到的聲源預(yù)報葉輪輻射聲和壁面邊界層這2種途徑下的管道噪聲，并將兩者疊加結(jié)果與實際結(jié)果進行比較分析。

1 基于大渦模擬的風機流場模擬

1.1 模擬對象

計算對象為某型離心風機，如圖1所示。離心風機重約117 kg，高約1.25 m；葉輪為后向葉輪，包含10個葉片，葉輪直徑為503 mm；轉(zhuǎn)速為1 470 r/min。采用管道進口、管道出口方式［8］對離心風機進行實驗。進口管道與出口管道直徑均為267 mm，進口管道長約12.3 m，出口管道長約11.6 m。測量情況如圖2所示，其中傳聲器安裝鼻錐，迎著來流，氣體從鼻錐光順地流過，可減少湍流脈動的影響。

圖1 離心風機與葉輪Fig.1 Centrifugal fan and impeller

圖2 離心風機噪聲實驗Fig.2 Noise measurement of the centrifugal fan

1.2 網(wǎng)格劃分

根據(jù)實驗中離心風機管道的連接情況，建立了整個流域的幾何模型。流域含旋轉(zhuǎn)葉輪，因此采用滑移網(wǎng)格方法模擬風機葉片的旋轉(zhuǎn)運動：在風機葉片所在區(qū)域劃分出一個小圓柱體區(qū)域（圖3），部分管道伸入葉輪內(nèi)部，導(dǎo)致交界面有階梯。小圓柱體內(nèi)包含葉片，此區(qū)域作為內(nèi)域，為轉(zhuǎn)動部分；其他區(qū)域作為外域，為靜止部分。靜止部分和轉(zhuǎn)動部分設(shè)置交界面并進行數(shù)據(jù)交換，以此來保證各物理量守恒。

采用結(jié)構(gòu)化網(wǎng)格形式劃分網(wǎng)格，內(nèi)域網(wǎng)格如圖4所示；外域網(wǎng)格存在由方變圓段，無法保證結(jié)構(gòu)化網(wǎng)格質(zhì)量，因此也分為2部分，網(wǎng)格劃分如圖5所示。

圖3 交界面劃分示意圖Fig.3 Schematic diagram of interface division

圖4 風機內(nèi)域結(jié)構(gòu)化網(wǎng)格Fig.4 Structured grid of internal zone of the centrifugal fan

圖5 風機外域結(jié)構(gòu)化網(wǎng)格Fig.5 Structured grid of exterior zone of the centrifugal fan

在大渦模擬之前，采用劃分方式相同而網(wǎng)格數(shù)量不同的多套網(wǎng)格（第1層網(wǎng)格厚度保持不變）進行穩(wěn)態(tài)流場模擬，以出口流量為判據(jù)進行網(wǎng)格收斂性分析，考慮到大渦模擬網(wǎng)格數(shù)量要求以及單機計算能力，選取了其中1套網(wǎng)格：內(nèi)域網(wǎng)格數(shù)為 5.19×106，蝸殼與管道網(wǎng)格數(shù)為 2.34×106，出口管道網(wǎng)格數(shù)為2.7×105；蝸殼第1層網(wǎng)格尺度設(shè)置為0.6 mm，葉片為0.3 mm，出口段為0.3 mm。在運算穩(wěn)定后，無因次壁面法向高度 y+的范圍為5～35，滿足計算要求。

1.3 大渦模擬計算風機內(nèi)部流場

湍流數(shù)值模擬方法可以分為直接數(shù)值模擬（DNS）方法和非直接數(shù)值模擬方法。大渦模擬既有直接模擬又有非直接模擬特征。它將流動分解為大、小2種尺度的湍流，對大尺度湍流進行直接模擬，而對小尺度湍流做近似處理。由于小尺度渦不直接依賴邊界條件，而且一般具有各向同性性質(zhì)，可以采用合適的湍流模型模擬。

大渦模擬計算對初始條件要求較高。首先采用RNG k-ε湍流模式下的結(jié)果作為初始條件，然后采用壁面自適應(yīng)局部渦粘性（Wall-Adapting Local Eddy-viscosity，WALE）亞格子模型，進行大渦模擬計算，其中內(nèi)域采用運動參考系。在計算穩(wěn)定后，內(nèi)域采用網(wǎng)格運動形式進行模擬，在運行一段時間后進行數(shù)據(jù)采集。數(shù)據(jù)采集的起始步數(shù)為第4 200步，對應(yīng)時刻為0.105 s，采樣間隔時間為2.5×10-5s。

葉輪旋轉(zhuǎn)1/10周（即一個葉片旋轉(zhuǎn)到上一個葉片位置）需要163步，1/10周再分為3個部分，分別繪制截面的壓力分布云圖，如圖6所示。從圖中可以看出，在計算過程中，葉片在不斷旋轉(zhuǎn)。在第4 212步與第4 374步，下一個葉片來到了同樣位置，但葉輪附近的壓力分布并不一致，這說明受渦流影響，葉輪所受壓力并不是完全重復(fù)的。在第4 212步時，葉片離渦舌部位很近，氣體受到擠壓導(dǎo)致壓強增大，葉片轉(zhuǎn)過之后，壓力下降。

圖6 葉輪旋轉(zhuǎn)1/10周時流場絕對壓力的變化Fig.6 Absolute pressure changes of flow field in the 1/10 cycle of impeller

2 離心風機管中噪聲預(yù)報

離心風機管中噪聲主要來源于旋轉(zhuǎn)葉片的氣動噪聲，另外蝸殼表面壓力脈動也會輻射噪聲，這2個激勵源均由大渦模擬計算得到。噪聲在管中傳播時還應(yīng)考慮管道的影響。

2.1 葉輪輻射聲計算

2.1.1 旋轉(zhuǎn)噪聲理論

首先，分析旋轉(zhuǎn)機械在自由場的輻射聲。葉輪工作時，葉片面元上會產(chǎn)生脈動壓力，特別是在葉輪出口附近存在渦脫落時，壓力脈動更明顯，這些壓力脈動會產(chǎn)生噪聲。而且葉輪在不斷旋轉(zhuǎn)，因此葉輪噪聲相當于旋轉(zhuǎn)偶極子輻射聲。

對于其中一塊面元，F(xiàn)={Fx，F(xiàn)y，F(xiàn)z} ，可以得到各分力產(chǎn)生的總聲壓［9］：

2.1.2 聲學有限元計算原理

葉輪在蝸殼中旋轉(zhuǎn)時，其噪聲傳播會受到管道聲模態(tài)的影響。對于有限空間中的輻射聲問題，一般采用邊界元法或有限元法。邊界元法盡管可以將維數(shù)降低，但求解速度慢［10］。聲學有限元法的一般思路是：聲場滿足微分方程（含邊界條件），可以構(gòu)建與之等價的泛函，整個聲場泛函可以用網(wǎng)格單元的形函數(shù)表達，對其求極值，即其變分為0，可以構(gòu)建整體的矩陣方程。

在均勻介質(zhì)、非粘性和絕熱狀態(tài)下，流體內(nèi)的聲學波動方程為

式中：ρ0為靜態(tài)密度；c=ω/k。令聲壓p'=p(x，y，z)ejωt，聲源q'=q(x，y，z)ejωt，代入式（2），得到Helmholtz方程：

采用加權(quán)余量法，可以得到聲場滿足的矩陣方程［11］：

式中：K為聲學剛度矩陣；C為聲學阻尼矩陣；Μ為質(zhì)量矩陣；p為聲壓向量；Q為聲源向量；V為速度向量；P為輸入聲壓向量。這些量均由計算域單元的形函數(shù)表達。

2.1.3 葉輪輻射聲預(yù)報

采用大渦模擬計算葉輪旋轉(zhuǎn)1周（從4200步到5831步）時的流場，輸出了葉輪表面在每一步的受力情況，文件格式為計算流體力學通用符號系統(tǒng)數(shù)據(jù)（CFD General Notation System Data，CGNS）。在Virtual.Lab軟件中導(dǎo)入相應(yīng)文件，由于采用有限體積法計算得到的脈動壓力結(jié)果存儲在網(wǎng)格形心上，需要將形心結(jié)果傳遞到網(wǎng)格各節(jié)點［10］。經(jīng)過傳遞，在某一時刻的壓力云圖如圖7所示。

圖7 第4 265步時葉片上的壓力云圖Fig.7 Pressure contours on impeller at the 4 265thstep

在Virtual.Lab軟件中進行扇聲源計算時，扇聲源被當作一類專門的計算對象。將其葉輪等效為一個個旋轉(zhuǎn)偶極子，即將葉片表面進行劃分，得到的每個小區(qū)域相當于1個聲源。聲源數(shù)與計算的最高頻率有關(guān)。圖8（a）所示的彩色色塊為在設(shè)置計算頻率后葉輪表面區(qū)域自動劃分的結(jié)果。然后，添加扇聲源的邊界條件，旋轉(zhuǎn)速度設(shè)置為1 470 r/min，傅里葉變換加窗選擇漢寧（Hanning）窗。圖8（b）為葉輪旋轉(zhuǎn)參數(shù)設(shè)置結(jié)果，由圖可以看出旋轉(zhuǎn)方向以及各個面積塊的離散力分布。

圖8 在Virtual.Lab軟件中的葉輪設(shè)置Fig.8 Impeller setting in software Virtual.Lab

通常在劃分聲學有限元網(wǎng)格時，1個波長不少于6個單元，但由于外域尺度大，網(wǎng)格數(shù)量較多，計算時間長且對計算機硬件要求高，所以采用自適應(yīng)階次有限元（FEM Adaptive Order，F(xiàn)EMAO）算法，其主要原理是采用高階多項式表達單元內(nèi)的參數(shù)，比較粗糙的網(wǎng)格也能得到同樣精度的計算結(jié)果［12］。

原來蝸殼入口處有部分伸到葉輪中，如圖9（a）所示，這是由于蝸殼厚度為2 mm，而且有卷邊。單個網(wǎng)格尺寸比較大，難以擬合其邊界，導(dǎo)致網(wǎng)格生成質(zhì)量不高，所以將卷起部分進行了簡化，如圖 9（b）所示。圖9（c）為劃分好的網(wǎng)格，其中網(wǎng)格均為四面體，單元最大邊長為64 mm，網(wǎng)格數(shù)為12.9×104。

圖9 聲場區(qū)域網(wǎng)格劃分過程Fig.9 Meshing process of sound field

進風管口與出風管口暴露在空氣中，所以管口可作為無反射邊界條件。在進、出風管口處的單元組上定義聲阻抗來模擬無反射邊界，將阻抗設(shè)置為416.5 N·s/m3。最大計算頻率為6 kHz，頻率間隔為24.5 Hz。求解場點坐標分別為：進風管（0 mm，0 mm，8 860 mm），出風管（-8 150 mm，368.7 mm，175.5 mm），這2點與實際測量點位置一樣。圖10為1/3倍頻程噪聲預(yù)報結(jié)果。由圖可看出，噪聲最大值為85 dB，進、出風管測點的噪聲在高頻時比較接近。

圖10 葉輪輻射聲預(yù)報結(jié)果Fig.10 Prediction results of impeller noise

2.2 蝸殼脈動壓力輻射聲計算

2.2.1 有固體壁面時湍流輻射聲

當湍流區(qū)內(nèi)存在剛性表面時，湍流區(qū)內(nèi)的密度變化量可用Curle方程表示［13］：

2.2.2 蝸殼壁面脈動壓力輻射聲預(yù)報

在大渦模擬時，蝸殼表面力無法同時輸出，可以采用另外一種軟件操作方法：即輸出計算過程中整體域的解，文件格式仍為CGNS。由于輸出文件包含了整個流體域的結(jié)果，自然也包括蝸殼與管道界面的結(jié)果。采用該方法后，每個文件均較大，因此步長增加1倍，為5×10-5s。數(shù)值計算從第4 200步開始，共計算816步，即葉輪旋轉(zhuǎn)1周。

管道內(nèi)表面是光順的，湍流強度低，僅產(chǎn)生邊界層的噪聲，且該噪聲的量級很低；而蝸殼內(nèi)表面有大量的渦分離，所以本次只對蝸殼內(nèi)表面部分進行積分計算。

數(shù)據(jù)存儲于體單元中，實際計算只需要蝸殼表面網(wǎng)格上的數(shù)據(jù)，所以要進行數(shù)據(jù)傳遞。流場計算時的網(wǎng)格較密，而殼體單元組的網(wǎng)格較稀疏，本文采用無能量損失算法（Conservative Maximum Distance）將力從流體網(wǎng)格映射到蝸殼殼體網(wǎng)格。在CFD時域計算時，考慮到聲學分析需要，對其進行快速傅里葉變換，圖11為150 Hz時的脈動壓力云圖。

圖11 150 Hz時的蝸殼壁面脈動壓力云圖Fig.11 Pressure fluctuation contours of volute wall at 150 Hz

同樣，在入口單元和出口單元組分別定義無反射邊界條件，將完成數(shù)據(jù)傳遞后的結(jié)果作為邊界條件，物理類型為表面偶極子（Surface Dipole）。最大計算頻率6 kHz，頻率間隔為25 Hz。1/3倍頻程噪聲計算結(jié)果如圖12所示。由圖12的噪聲預(yù)報結(jié)果可以看出，不同位置處的輻射聲在低頻時有差別、在高頻時趨于一致。

圖12 蝸殼壁面脈動壓力輻射聲預(yù)報結(jié)果Fig.12 Prediction results of volute wall pressure fluctuation noise

3 風機管中噪聲計算結(jié)果

3.1 進風管中葉輪輻射聲與蝸殼壁面脈動壓力輻射聲對比

圖13為前述進風管中葉輪輻射聲與蝸殼壁面脈動壓力輻射聲的對比結(jié)果。由圖可以看出，進風管內(nèi)蝸殼流噪聲在低頻部分高于葉輪輻射聲，而在高頻時又低于葉輪輻射聲。這表明葉輪和蝸殼壁面的非定常力是離心風機最主要的氣動噪聲源，證實了文獻［14］中的結(jié)論。

圖13 進風管中葉輪輻射聲與蝸殼脈動壓力輻射聲的對比結(jié)果Fig.13 Contrast between impeller noise and volute pressure fluctuation noise

一般葉片噪聲較大，但本例中離心風機葉輪與蝸殼距離非常近，在蝸舌部位距離僅有28.5 mm，表明分離的渦（脫離渦將造成較大的脈動，從而產(chǎn)生較大輻射聲）尾流還會直接作用于蝸殼表面。另外，從聲源作用面積來講，10個葉片的總面積（含葉面和葉背2面）為0.36 m2，而蝸殼內(nèi)表面面積為2.2 m2，后者遠大于前者。這些因素導(dǎo)致蝸殼內(nèi)流噪聲較大。但高頻部分葉輪輻射聲又增加了，這部分跟直接渦發(fā)放有關(guān)。在大于臨界雷諾數(shù)（Re＞2×106）情況下，葉片渦的發(fā)放頻率滿足以下公式［15］：

式中：Yf為靠近葉片出氣邊處尾流的剪切層厚度，實際相當于渦街的橫向?qū)挾龋籙s為出氣邊處的切向速度，St為斯特勞哈爾數(shù)。由于葉片為厚2.5 mm的直平頭，Yf厚度為2.25 mm；葉輪出氣邊速度由流場得到，約為38 m/s，由此得到渦發(fā)放頻率為3 378 Hz。一般在渦發(fā)放頻率附近的噪聲會明顯增高，由圖13可知，在3 150～4 000 Hz中心頻段處噪聲級較大，跟3 378 Hz這一頻率吻合。

3.2 計算結(jié)果與實驗驗證

葉輪與蝸殼表面脈動力雖然從同一時刻開始計算，但步長不一樣，所以將他們當作無規(guī)相位聲波進行疊加，并對比疊加結(jié)果與實驗結(jié)果（圖14），由圖可以看出，與試驗結(jié)果相比，數(shù)值計算結(jié)果在中低頻處偏低，在高頻處有點偏高，最大偏差約14 dB，兩者整體趨勢一致。

圖14 實驗結(jié)果與數(shù)值計算結(jié)果對比Fig.14 Comparison between experimental results and predicted results

數(shù)值計算結(jié)果的精度與整個預(yù)報過程的每個環(huán)節(jié)都有關(guān)系：在流域建模、流場預(yù)報、輻射聲預(yù)報都會存在影響。其中流場預(yù)報時與網(wǎng)格個數(shù)、交界面設(shè)置（由于內(nèi)域旋轉(zhuǎn)需設(shè)置此項）、大渦模擬計算本身精度有關(guān)；聲學計算只考慮力源項的結(jié)果，采用軟件Virtual.Lab求解時，將葉輪表面劃分為上千個集中面源，這也降低了精度。另外，由圖14還可以看出，預(yù)報結(jié)果出現(xiàn)了較多的尖峰，這與采樣時間短、在計算時沒有將采樣結(jié)果當作隨機信號有關(guān)。而且葉輪對聲場的反射無法計算，流速對聲傳播的影響也未考慮。

4 結(jié) 論

本文采用結(jié)構(gòu)化網(wǎng)格，使用大渦模擬方法獲得了離心風機流場多個時刻的壓力分布，分別計算了風機管中葉輪輻射聲和蝸殼脈動壓力輻射聲，得出以下結(jié)論：

1）葉輪輻射聲和蝸殼脈動壓力輻射噪聲相比，在低頻時前者明顯低于后者，這與蝸殼輻射聲面積較大有關(guān)；在高頻時前者又會高于后者，這與葉片較強的渦發(fā)放有關(guān)。

2）與噪聲實驗結(jié)果相比，噪聲數(shù)值計算結(jié)果在低頻時偏低，在高頻時偏高，兩者整體趨勢基本一致，能夠滿足工程分析需要。

3）數(shù)值計算結(jié)果與實驗結(jié)果較吻合，驗證了葉輪輻射聲和蝸殼壓力脈動輻射噪聲預(yù)報方法的有效性，也從側(cè)面反映出兩者是管中噪聲主要的噪聲源。