基于模態(tài)展開法的二維消聲器傳遞損失分析

林 力，張 旭，李 云，程光旭

（西安交通大學(xué) 能源與動力工程學(xué)院，西安 710049）

壓縮機(jī)進(jìn)排氣系統(tǒng)產(chǎn)生的壓力脈動不僅影響壓縮機(jī)的熱力性能而且會誘發(fā)管路振動。由于壓縮氣體所含的雜質(zhì)會阻塞阻性消聲器的微孔和多孔材料，抗性消聲器更為廣泛地用于壓縮機(jī)系統(tǒng)以便減弱氣流脈動［1］。對于小型全封閉制冷壓縮機(jī)，薄層的抗性消聲器因為能更有效地利用壓縮機(jī)有限的內(nèi)部空間，而得到了更多的應(yīng)用。

Munjal［2］推薦傳遞矩陣法用于分析抗性消聲器傳遞損失。這種方法在一維平面波范圍內(nèi)計算方便，效率高但無法反映消聲器的高次波效應(yīng)。Kim和Soedel［3］采用模態(tài)展開法分析了具有高次波效應(yīng)的三維聲腔的壓力脈動。Lai和 Soedel［4－7］采用了相同的方法研究了二維消聲器的壓力脈動。但是，Soedel等人［3－7］之前的研究都采用了點(diǎn)聲源模型計算消聲器的壓力響應(yīng)函數(shù)，而Zhou和Kim［8］指出點(diǎn)聲源模型會在聲源區(qū)域附近產(chǎn)生數(shù)值收斂的問題。為了解決這個問題，Zhou和Kim［8］建議采用面聲源模型計算壓力響應(yīng)函數(shù)。Kadam和Kim［9］采用了面聲源模型計算了三維矩形消聲器的四極參數(shù)。

雖然二維消聲器元件的性能在文獻(xiàn)［4，5］中已經(jīng)討論過，但是由于其采用的點(diǎn)聲源模型在聲源附近存在數(shù)值收斂問題，得到壓力響應(yīng)的值不一定準(zhǔn)確。另外，面聲源模型雖然也已經(jīng)在消聲器性能的分析中得到了應(yīng)用［8，9］，但是這些分析的對象都局限于矩形消聲器。根據(jù)文獻(xiàn)查閱，目前并沒有有關(guān)二維消聲器傳遞損失分析方面的文獻(xiàn)。因此，本文運(yùn)用模態(tài)展開法分析了不同形狀的二維抗性消聲器的傳遞損失。為了克服點(diǎn)聲源模型在源點(diǎn)區(qū)域的數(shù)值收斂問題，本文采用面聲源推導(dǎo)了矩形和圓柱形二維消聲器的傳遞損失計算模型，并將計算結(jié)果與三維有限元法進(jìn)行了比較、驗證。此外，還將一個二維矩形消聲器傳遞損失的結(jié)果與文獻(xiàn)［10］中提到格林函數(shù)法的計算結(jié)果進(jìn)行了比較和分析，討論了兩種方法的差異。

1 傳遞損失數(shù)學(xué)模型

在推導(dǎo)不同形狀的二維消聲器的傳遞損失之前，首先對運(yùn)用模態(tài)展開法求解壓力響應(yīng)函數(shù)的方法進(jìn)行簡要的回顧。傳遞矩陣定義了頻域范圍內(nèi)聲學(xué)系統(tǒng)變量進(jìn)出口的關(guān)系。

圖1 消聲器聲學(xué)系統(tǒng)示意圖Fig.1 Schematic representation of acoustic system of a muffler

對于圖1所示的聲學(xué)系統(tǒng)，其傳遞函數(shù)可以表示為：

下標(biāo)1和2分別代表入口和出口。P和Q是聲壓和體積流率的諧振幅。A，B，C，D代表四極參數(shù)。它們可以用壓力響應(yīng)函數(shù)表示如下［3］：

fij表示系統(tǒng)在j點(diǎn)受到單位容積流率的激勵并同時阻塞另一點(diǎn)i的壓力響應(yīng)。從方程（2）可知，如果已知聲學(xué)系統(tǒng)的壓力響應(yīng)函數(shù)，那么其相應(yīng)的傳遞矩陣也可以得到。

根據(jù)文獻(xiàn)［2］，使用四極參數(shù)推導(dǎo)的傳遞損失的公式可以如下表示：

Z1和Z2分別是系統(tǒng)入口和出口的特征阻抗。Z1=ρc/S1，Z2=ρc/S2，其中S1和S2分別代表入口和出口管道的截面積。ρ和c分別代表流體介質(zhì)的密度和音速。

正交曲面坐標(biāo)系下的二維非齊次聲波方程可以表示為：

式中，p代表聲壓是質(zhì)量流率，▽2是二維曲面坐標(biāo)系下的拉普拉斯算子［4］：

式中，α1，α2代表一組正交的曲線坐標(biāo)，A1和A2表示拉梅系數(shù)。

1.1 二維矩形消聲器

對于如圖2所示的薄矩形消聲器，采用矩形坐標(biāo)系：α1=x，α2=y，A1=1，A2=1，進(jìn)行分析。波動方程表示為：

圖2 采用面源模型的二維矩形消聲器Fig.2 Rectangular thin muffler element with surface source model

使用正方形面聲源模型，質(zhì)量流率可以表示為：

其中δ（·）是狄拉克函數(shù)，H（·）表示單位階躍函數(shù)，L是正方形聲源的邊長。（x1，y1）是聲源的中心坐標(biāo)。h是該消聲器的厚度。

使用模態(tài)展開法［3］，可以得到源點(diǎn)中心的壓力響應(yīng)函數(shù)：

其中：

式中：Pmn和ωmn分別代表（m，n）階復(fù)模態(tài)振幅和固有圓周頻率。Lx和Ly分別代表矩形消聲器的長度和寬度。Nmn由下式給出：

由于在面源上的壓力分布不是常數(shù)，壓力響應(yīng)函數(shù)是由計算面源上壓力響應(yīng)的平均值得到［8］：

令方程（12）中（x，y）=（x1，y1）并且Q=1，可以得到壓力響應(yīng)函數(shù)f11。通過相同的方法，也可以得到f21，f22，f12。通過計算得到的壓力響應(yīng)函數(shù)可以由方程（2）和（3）計算不同尺寸，不同進(jìn)出口位置的二維矩形消聲器的傳遞損失。

1.2 二維圓柱形消聲器

圖3為常用于小型全封閉制冷壓縮機(jī)的二維圓柱形消聲器。方程（4）采用圓柱坐標(biāo)系：α1=θ，α2=z，A1=r和A2=1，進(jìn)行計算。其波動方程為：

圖3 采用面源模型的二維圓柱形消聲器Fig.3 Cylindrical thin muffler element with surface source model

采用面聲源模型，質(zhì)量流率表示為：

其中Θ代表面源的弧度，Z代表面源的高度。（θ1，z1）是面源中心點(diǎn)的坐標(biāo)。R代表此圓柱的半徑。采用模態(tài)展開法，可以得到面源中點(diǎn)的壓力響應(yīng)函數(shù)：

其中：

Pmnl和ωmnl分別代表（m，n，l）階復(fù)模態(tài)振幅和固有圓周頻率。R代表圓柱體的半徑。Lz是圓柱的高。Nmnl由下式給出：

由前述分析可知。由于矩形消聲器在面源上的壓力分布并不均勻，因此平均的壓力響應(yīng)函數(shù)可以定義為：

通過方程（19）計算f11，f21，f22和f12，就可以得到圓柱形消聲器的四極參數(shù)和傳遞損失。

2 有限元法驗證

為了驗證模態(tài)展開法推導(dǎo)出的二維薄型消聲器的傳遞損失模型，本節(jié)首先采用聲學(xué)有限元法建立了三維有限元模型。分析采用商用有限元軟件COMSOL［11］。對1.2節(jié)中提到的圓柱薄形消聲器建立完整的三維有限元模型。其尺寸為：Lz=0.185 m，R=0.07 m，h=0.01 m，Θ=0.142，Z=0.01 m。入口、出口中心點(diǎn)坐標(biāo)分別為：（0.092 5，0）和（0.092 5，π）。分析選用密度ρ=6.04 kg/m3聲速c=162.9 m/s的冷媒作為介質(zhì)。傳遞損失計算通過文獻(xiàn)［2］和文獻(xiàn)［10］中給出的進(jìn)出口的速度和聲壓得到。傳遞損失可以表示為：

其中pin是入口的平均聲壓，uin和uout分別為進(jìn)出口的平均速度值。對此消聲器建立由8 668個四面體單元構(gòu)成的網(wǎng)格。該網(wǎng)格所對應(yīng)的最高頻率為2 000 Hz的單位波長，其節(jié)點(diǎn)數(shù)為9，滿足聲學(xué)有限元計算中單位波長至少為6的基本要求。經(jīng)過驗證，進(jìn)一步加密網(wǎng)格并不能得到更為精確的傳遞損失計算結(jié)果。為了抑制進(jìn)口處的壓力反射，進(jìn)口采用平面波輻射邊界條件；為了模擬無回聲的出口邊界條件，在出口處采用阻抗邊界條件。傳遞損失的計算范圍為10 Hz到2 000 Hz，分辨率為10 Hz。為了比較三維有限元法和模態(tài)展開法的計算效率，計算都采用相同配置的計算機(jī)（CPU：Intel? CoreTM2 2.33GHz，ROM：2.96GB）進(jìn)行。對于模態(tài)展開法，采用面聲源模型對于200個頻率步長的傳遞損失計算耗費(fèi)的CPU時間為91 s，而三維有限元法耗費(fèi)的CPU時間為243 s，該時間不包括前后處理的時間。

計算結(jié)果如圖4中所示，模態(tài)展開法得到的計算結(jié)果與三維有限元法得到的結(jié)果在低頻范圍內(nèi)的結(jié)果相當(dāng)一致。模態(tài)展開法得到的結(jié)果反映了消聲器的高次波效應(yīng)。由于模態(tài)展開建立的計算模型是二維的，而有限元法建立的模型本質(zhì)上是三維的，因此兩種方法得到的結(jié)果在1 500 Hz～2 000 Hz范圍內(nèi)有微小的偏差。另外，對于模態(tài)展開法而言，由于參與計算的模態(tài)數(shù)是有限的（本例中為40 000），因此會產(chǎn)生一定截斷誤差，這也是兩種方法計算結(jié)果差異的原因之一。

圖4 傳遞損失計算結(jié)果比較Fig.4 Sound transmission loss of a cylindrical thin muffler

圖5 二維矩形消聲器傳遞損失Fig.5 Transmission loss of a two-dimensional rectangular muffler

圖6 平均相對誤差Fig.6 Averaged relative errors

3 與格林函數(shù)法的比較討論

文獻(xiàn)［10］使用格林函數(shù)法推導(dǎo)了三維矩形消聲器的四極參數(shù)和傳遞損失。該方法同樣可以用于分析計算二維矩形消聲器的傳遞損失。為了將本文的方法與格林函數(shù)法進(jìn)行比較，選用一尺寸為Lx=0.1 m，Ly=0.08 m，h=0.01 m，進(jìn)出口位置為：（0.02，0.04）和（0.08，0.04）的二維矩形消聲器進(jìn)行計算。圖5比較了采用兩種方法得到的該消聲器傳遞損失的計算結(jié)果。從圖5可以看出，當(dāng)選用的面聲源的尺寸控制在較小范圍內(nèi)（本例為L=0.005 m）時，采用兩種方法得到的計算結(jié)果并沒有顯著的差異。為了比較兩種方法的計算誤差，對此消聲器建立了聲學(xué)三維有限元模型進(jìn)行計算，并將有限元方法得到的傳遞損失作為標(biāo)準(zhǔn)，比較模態(tài)展開法和格林函數(shù)法的計算結(jié)果的誤差。誤差定義式由下式給出：

從圖6中可以看出，模態(tài)展開法和格林函數(shù)法計算結(jié)果的誤差都隨著面源尺寸的不斷增加而增加，而且格林函數(shù)法的誤差相對于模態(tài)展開法的誤差隨面源邊長增加得更快。對于格林函數(shù)法，此方法是建立在面源上均勻速度分布［10］的假設(shè)上的，該假設(shè)只能在進(jìn)出口的尺寸比相關(guān)波長尺寸小得多的條件下才成立。但是隨著面聲源的尺寸不斷加大，這一條件并不能得到滿足，因此計算的誤差也隨之增大。而對于模態(tài)展開法，在推導(dǎo)其壓力響應(yīng)函數(shù)時，主要的過程是對質(zhì)量流率點(diǎn)源在進(jìn)出口面上對壓力響應(yīng)進(jìn)行積分，然后通過求解壓力響應(yīng)的平均值得到相應(yīng)的四極參數(shù)，在此過程中并不需要假設(shè)面源上速度的均勻分布。因此，當(dāng)面聲源尺寸相對較大時，采用模態(tài)展開法進(jìn)行計算的精度較格林函數(shù)法更高。另一方面，從計算效率來看，由于模態(tài)展開法相對于格林函數(shù)法在計算四極參數(shù)時需要多一個步驟，即求解壓力響應(yīng)在面源上的平均值［8］，因此耗費(fèi)更多的計算時間。在本例的計算中，格林函數(shù)法耗費(fèi)的CPU時間為13 s而模態(tài)展開法為75 s。因此在計算二維消聲器的傳遞損失時到底采用哪種方法還需要根據(jù)具體的要求進(jìn)行選擇。如果面源的尺寸比較小，選用格林函數(shù)法具有較高的計算效率。但是當(dāng)面源尺寸增加時，模態(tài)展開法的計算精度更高。

4 結(jié)論

采用模態(tài)展開法和面源模型建立了矩形和圓柱形二維消聲器的傳遞損失計算模型，并通過有限元法驗證了該方法的正確性。通過與三維有限元法分析比較，模態(tài)展開法在計算二維矩形和二維圓柱形消聲器的傳遞損失時具有更高的計算效率。另外，還與格林函數(shù)法進(jìn)行了比較。如果面源的尺寸較小，選用格林函數(shù)法需要的計算時間較少，但是當(dāng)面源尺寸不能滿足速度均勻分布的假設(shè)時，格林函數(shù)法與模態(tài)展開法相比誤差較大。

［1］Soedel W.Sound and vibrations of positive displacement compressors［M］.CRC Press，Boca Raton，F(xiàn)L，2007.

［2］Munjal M L.Acoustics of ducts and mufflers with application to exhaust and ventilation system design［M］.John Wiley＆Sons，INC，1987.

［3］Kim J， SoedelW. Generalformulation offourpole parameters forthree-dimensionalcavitiesutilizing modal expansion，with special attention to the annular cylinder［J］.Journal of Sound and Vibration，1989，129：237 －254.

［4］Lai P C C，Soedel W.Two dimensional analysis of thin，shell or plate like muffler elements［J］.Journal of Sound and Vibration，1996，194：137 －171.

［5］Lai P C C，Soedel W.Two dimensional analysis of thin，shell or plate like muffler elements of non-uniform thickness［J］.Journal of Sound and Vibration，1996，195：445 －475.

［6］Lai P C C，Soedel W.Free gas pulsations in acoustic system composed of two thin，curved or flat，two dimensional gas cavities which share a common open boundary［J］.Journal of Sound and Vibration，1996，198：225 －248.

［7］Lai P C C，Soedel W.Gas pulsation in thin，curved or flat cavities due to multiple mass flow sources with special attention to multi-cylinder compressor［J］.Journal of Sound and Vibration，1996，197：45 －66.

［8］Zhou W，Kim J.Formulation of four poles of threedimensional acoustic systems from pressure response functions with special attention to source modeling［J］.Journal of Sound and Vibration，1999，219：89 －103.

［9］ Kadam P，Kim J.Experimental formulation of four poles of three-dimensional cavities and its application［J］.Journal of Sound and Vibration，2007，307：578 －590.

［10］ Venkatesham B，Tiwari M，Munjal M L.Transmission loss analysis of rectangular expansion chamber with arbitrary location of inlet/outlet by means of Green's functions［J］.Journal of Sound and Vibration，2009，323：1032 －1044.

［11］ COMSOL Multiphysics rev.3.5a Acoustic Module User's Guide［M］.COMSOL Inc，USA，2008.