消聲器傳遞損失預測的邊界元數(shù)值配點混合方法

楊　亮， 季振林， Wu T W

(1.哈爾濱工程大學 動力與能源工程學院，哈爾濱　150001；2.肯塔基大學 機械工程學院，列克星敦　40506)

?

消聲器傳遞損失預測的邊界元數(shù)值配點混合方法

楊亮1， 季振林1， Wu T W2

(1.哈爾濱工程大學 動力與能源工程學院，哈爾濱150001；2.肯塔基大學 機械工程學院，列克星敦40506)

摘要：將邊界元法與數(shù)值配點法結合形成混合方法用于計算任意截面形狀消聲器的傳遞損失。消聲器劃分為若干子結構，用邊界元法計算具有非規(guī)則形狀的子結構阻抗矩陣，用二維有限元法提取等截面子結構特征值及特征向量，用配點法獲得阻抗矩陣；將每個子結構阻抗矩陣連接用于傳遞損失計算。為減少計算時間提出簡化方法計算消聲器傳遞損失。結果表明，混合法在保證計算精度前提下可節(jié)省計算時間。

關鍵詞：消聲器；傳遞損失；邊界元方法；數(shù)值配點方法

解析方法、數(shù)值方法廣泛用于消聲器傳遞損失計算。前者雖計算速度快、精度高，但只能用于規(guī)則形狀消聲器計算，Selamet等[1-3]用解析模態(tài)匹配法計算膨脹腔、內插管膨脹腔及回流腔等消聲器傳遞損失。Selamet等[4]用二維解析模態(tài)匹配法計算填充吸聲材料的阻性消聲器聲學性能。數(shù)值方法適用性強，理論上可用于任意復雜結構消聲器計算，Peat等[5]用有限元、Ji[6]用邊界元、Wu等[7-8]用直接混體邊界元法對不同類型消聲器進行廣泛研究。其中，直接混體邊界元法更具靈活性，尤其適用消聲器的聲學性能計算。對大尺度問題，數(shù)值方法計算時間長、內存需求大。因此，考慮將解析法與數(shù)值法結合用于消聲器聲學性能研究具有實際意義。Kirby等[9-12]用二維有限元法獲得消聲器或空調管道截面的特征值及特征函數(shù)，用配點法獲得相關聲學性能評價參數(shù)，并用數(shù)值模態(tài)匹配法對阻性消聲器進行研究。Fang等[13-14]用數(shù)值模態(tài)匹配法對抗性消聲器及有流穿孔阻性消聲器進行計算。

無論數(shù)值配點法或數(shù)值模態(tài)匹配法，對沿軸向方向消聲器或消聲管道截面形狀變化情況均不適用，因此，不同計算方法的耦合及優(yōu)化需進一步研究。

本文擬用邊界元-數(shù)值配點混合法進行消聲器傳遞損失計算。將消聲器劃分為若干子結構，對沿軸向管道截面變化或具有復雜內部結構的子結構用直接混體邊界元法計算其阻抗矩陣，考慮三維邊界元法與二維有限元法網(wǎng)格劃分的一致性，對截面均勻一致的子結構考慮用二維有限元方法進行特征值分析，獲得橫向波數(shù)，進而用配點法計算獲得阻抗矩陣，將不同子結構阻抗矩陣連接獲得消聲器的傳遞損失。為減少計算時間，用簡化方法壓縮子結構阻抗矩陣，并通過算例驗證混合方法的有效性。

1直接混體邊界元方法

該方法源于傳統(tǒng)子結構邊界元方法，即將具有復雜內部結構的消聲器分成若干具有明確邊界的幾部分，對每部分用邊界積分方程描述，據(jù)各部分交界面的聲學參數(shù)連續(xù)條件獲得整體控制方程。交界面有兩個未知量時，用超奇異積分方程提供額外方程。本文涉及的邊界積分方程(完整積分見文獻[8])為

(1)

(2)

圖1　消聲器子結構劃分Fig.1 Substructure division of silencer

源點、場點距離較近時，對公式中涉及的奇異積分及超奇異積分需妥善處理[15-16]。不同于文獻[7-8]求解矩陣的直接方法，本文用帶預處理的BI-CG求解器。為描述問題，將消聲器分三部分(圖1)，對子結構A、C，考慮用邊界元方法獲得阻抗矩陣；對子結構A，阻抗矩陣可表示為

(3)

式中：PAI，vAI分別為進口聲壓及質點振速向量；PAO，vAO分別為出口聲壓及質點振速向量。

阻抗矩陣作為結構聲學特性描述，是唯一確定的。對式(3)需在不同邊界條件下用BEM計算。如，為獲得阻抗矩陣第一列數(shù)值，對進口第一單元設v=1，進出口其它單元設v=0，進出口聲壓值即為第一列阻抗矩陣值。盡管需多次運行BEM，但每次求解的矩陣相同。對結構C，得阻抗矩陣為

(4)

2數(shù)值配點法

圖1的子結構B沿軸向等截面，但由于截面形狀任意，無法用解析法獲得橫向波數(shù)，故用數(shù)值方法進行特征值分析，考慮二維有限元及三維邊界元網(wǎng)格劃分的一致性，為便于各子結構阻抗矩陣連接，用二維有限元計算特征值及特征向量。以阻性穿孔管道為例，圖1中Ⅰ的內部為空氣，Ⅱ為填充吸聲材料，Ⅰ與Ⅱ之間由穿孔管分隔。區(qū)域Ⅰ、Ⅱ內橫向聲壓控制方程可表示為

(5)

(6)

式中：pxy1，pxy2分別為區(qū)域Ⅰ、Ⅱ的橫向聲壓分量；kxy1，kxy2分別為空氣、吸聲材料中橫向波數(shù)，兩區(qū)域軸向波數(shù)kz相同，與橫向波數(shù)滿足關系為

(7)

(8)

穿孔管內外滿足邊界條件為

(9)

用伽遼金加權余量法并代入邊界條件獲得阻性穿孔管道橫截面的本征方程[17]為

(10)

(11)

(12)

求解式(10)可獲得軸向波數(shù)及本征向量Φn(x,y)。對穿孔管抗性消聲器，式(10)簡化[18]為

(13)

對剛性壁直管道，橫截面本征方程表達式為

(14)

管道內聲壓、質點振速表示為

(15)

(16)

(17)

在B的進出口各選N個配點，則有N個模態(tài)，據(jù)式(16)，獲得B進出口配點振速為

(18)

同理，據(jù)式(15)得進出口聲壓表達式為

(19)

3消聲器傳遞損失計算

若使子結構B進出口配點分別與子結構A的出口截面、子結構C的進口截面節(jié)點一致，可連接二者阻抗矩陣獲得整體阻抗矩陣用于傳遞損失計算。對圖1消聲器，按結構連接順序，建立前兩子結構A、B進出口阻抗矩陣關系，即

(20)

式中：

用相同方法，將式(20)中阻抗矩陣與子結構C的阻抗矩陣連接可得整個結構進出口的阻抗矩陣關系，即

(21)

對含更多子結構的消聲器，阻抗矩陣計算方法相同。在消聲器進出口，聲波可視為平面波，此時式(21)被壓縮為2×2的矩陣。消聲器四極參數(shù)表示為

(22)

傳遞損失計算公式可表示為

(23)

式中：SI，SO分別為消聲器進出口面積。

4算例驗證

以一方形膨脹腔消聲器為例驗證本文方法的正確性，見圖2。將消聲器劃分為3個子結構，lA=lC=0.05 m，lB=0.125 m，進出口管道長寬均b=0.05 m，膨脹腔長寬均bc=0.15 m。子結構A、C的阻抗矩陣由BEM獲得，子結構B的阻抗矩陣由配點法計算，各子結構交界面均生成144個單元，即可考慮144個模態(tài)，計算結果見圖3。由圖3看出，混合方法計算結果與三維邊界元吻合較好，該方法的正確性獲得驗證。用數(shù)值配點方法計算阻抗矩陣，只需在子結構進出口劃分網(wǎng)格，可節(jié)省計算時間，該時間與子結構長度無關，適用于具有較長等截面結構的大尺度聲學問題。需注意的是，由于橫向波數(shù)及本征向量由二維有限元方法獲得，較解析方法相比存在一定誤差，且隨模態(tài)階數(shù)提高誤差有增大趨勢。據(jù)文獻[9],本算例認為最高到第29個模態(tài)結果較準確，其后高階模態(tài)對應的特征值、特征向量與準確值誤差較大。因高階模態(tài)對結果貢獻量小，故傳遞損失計算結果相對準確，解析法與數(shù)值法對前5個高階模態(tài)計算結果比較見表1。

圖2　方形膨脹腔消聲器Fig.2 Square expansion chamber silencer

圖3　方形膨脹腔消聲器傳遞損失Fig.3 TL of square expansion chamber silencer

模態(tài)解析解數(shù)值解相對誤差/%(1,0)20.94421.0040.29(1,1)29.61929.7040.29(2,0)41.88842.3681.15(1,2)46.83247.2880.97(2,2)59.23859.9171.15

若模型尺寸較大，為保證與子結構A、C截面網(wǎng)格一致，B截面會劃分較多網(wǎng)格，對式(18)、(19)，求解模態(tài)幅值系數(shù)會費較多計算時間，實際應用中，考慮過多模態(tài)無意義。本算例計算頻率最高到2 500 Hz，(3,0)模態(tài)或(0,3)模態(tài)激發(fā)頻率約3 440 Hz，實際計算中考慮最高到第二階模態(tài)即可滿足需要。

因此，如何保證計算精度的同時控制模態(tài)數(shù)量尚待進一步研究。本文考慮壓縮子結構A、C的阻抗矩陣，在A的出口截面及C的進口截面劃分若干子區(qū)域，用區(qū)域內部節(jié)點平均聲學變量近似表示整個區(qū)域聲學變量數(shù)值，對子結構B的進出口截面，選擇與子結構A、C劃分子區(qū)域中心對應的配點，進而用數(shù)值配點法計算阻抗矩陣。子結構A進口及C出口認為平面波傳播，只劃分1個子區(qū)域。

以方形穿孔管抗性消聲器為例，驗證簡化方法的正確性。圖4的結構尺寸與圖2相同，區(qū)別為膨脹腔內部添加穿孔管，穿孔率φ=0.08，穿孔孔徑dh=0.002 49 m，穿孔壁厚tw=0.000 9 m。子結構交界面區(qū)域簡化劃分見圖5，各子結構交界面劃分為9個子區(qū)域，bs= 1/3bc=0.05 m。簡化后式(3)、(4)中矩陣規(guī)模減小為10×10。對子結構B，選擇與子區(qū)域1～9中心一致的配點，阻抗矩陣計算同前，式(17)中矩陣規(guī)模減小為18×18，大大減少計算時間。傳遞損失計算結果見圖6，可見與三維邊界元吻合良好，驗證所提方法實際應用的有效性。

圖4　方形穿孔管消聲器結構尺寸Fig.4 Structure dimension of square perforated tube silencer

圖5　子結構A、C截面劃分Fig.5 Cross-section division of substructures A and C

圖6　方形穿孔管消聲器傳遞損失Fig.6 TL of square perforated tube silencer

通過兩個算例驗證所提方法的有效性，并與傳統(tǒng)子結構方法計算時間比較，說明本文方法的優(yōu)越性。考慮穿孔管阻性消聲器，模型結構尺寸同圖4，區(qū)別為膨脹腔內填充吸聲材料，流阻率R=4 896 Rayls/m。子結構劃分形式與交界面子區(qū)域劃分同圖5，簡化的混合方法與三維邊界元計算結果見圖7，可見吻合良好。

圖7　方形穿孔管阻性消聲器傳遞損失Fig.7 TL of square dissipative perforated tube silencer

橢圓截面消聲器亦廣泛應用，故以橢圓膨脹腔消聲器為例說明簡化混合方法的應用。模型尺寸見圖8，lA=lC=0.05 m，lB=0.182 3 m，橢圓長軸b1=0.2 m，短軸b2=0.14 m，進出口插管長l1=l2=0.025 m，消聲器進出口管道截面為圓形，半徑r=0.025 m。各子結構交界面子區(qū)域劃分見圖9，傳遞損失計算結果見圖10，可見計算結果與三維邊界元吻合良好。

圖8　橢圓膨脹腔消聲器結構尺寸Fig.8 Structure dimension of elliptical expansion chamber silencer

圖9　橢圓截面子區(qū)域劃分Fig.9 Cross-section division of elliptical interface

圖10　橢圓膨脹腔消聲器傳遞損失Fig.10 TL of elliptical expansion chamber silencer

簡化混合方法與傳統(tǒng)子結構邊界元方法在頻率1 000 Hz時計算時間比較見表2。以方形穿孔管阻性消聲器與橢圓膨脹腔消聲器為例，所有計算均在Intel i7 3632QM, 2.2GHz計算機上進行。由于用簡化混合方法計算B的阻抗矩陣，整個子結構無需劃分網(wǎng)格，只選適當?shù)呐潼c，每次計算求解的矩陣規(guī)模遠小于傳統(tǒng)邊界元方法，子結構B的阻抗矩陣計算時間基本可忽略不計，因此簡化混合方法能節(jié)省大量機時。實際應用中，對尺寸較大消聲器，劃分子結構進行傳遞損失計算時，混合方法優(yōu)勢較明顯。

表2　消聲器傳遞損失計算時間比較

5結論

(1) 對含若干等截面、非等截面子結構消聲器，可分別考慮用配點法、邊界元方法計算各自阻抗矩陣，并將其連接進行傳遞損失計算。

(2) 壓縮阻抗矩陣可節(jié)省機時。簡化混合方法在保證計算精度同時能有效減少計算時間，尤其適用具有較長等截面子結構大尺度聲學問題計算。

參 考 文 獻

[1] Selamet A, Ji Z L. Acoustic attenuation performance of circular expansion chambers with offset inlet/outlet: I.analytical approach[J]. Journal of Sound and Vibration, 1998, 213(4): 601-617.

[2] Selamet A, Ji Z L. Acoustic attenuation performance of circular expansion chambers with extended inlet/outlet[J]. Journal of Sound and Vibration, 1999, 223(2): 197-212.

[3] Selamet A, Ji Z L. Acoustic attenuation performance of circular flow-reversing chambers[J]. Journal of the Acoustical Society of America, 1998, 104(5): 2867-2877.

[4] Selamet A, Xu M B, Lee I J, et al. Analytical approach for sound attenuation in perforated dissipative silencers with inlet/outlet extensions[J]. Journal of the Acoustical Society of America, 2005, 117(4): 2078-2089.

[ 5] Peat K S, Rathi K L. A finite element analysis of the convected acoustic wave motion in dissipative silencers[J]. Journal of Sound and Vibration, 1995, 184(3): 529-545.

[6] Ji Z L. Boundary element acoustic analysis of hybrid expansion chamber silencers with perforated facing[J]. Engineering Analysis with Boundary Elements, 2010,34(7): 690-696.

[7] Wu T W, Cheng C Y R, Zhang P. A direct mixed-body boundary element method for packed silencers[J]. Journal ofthe Acoustical Society of America, 2002, 111(6): 2566-2572.

[8] Lou G, Wu T W, Cheng C Y R. Boundary element analysis of packed silencers with a substructuring technique[J]. Engineering Analysis with Boundary Elements, 2003,27(7): 643-653.

[9] Kirby R, Lawrie J B. A point collocation approach to modelling large dissipative silencers[J]. Journal of Sound and Vibration, 2005, 286(1): 313-339.

[10] Kirby R, Williams P T, Hill J. A three dimensional investigation into the acoustic performance of dissipative splitter silencers[J]. The Journal of the Acoustical Society of America, 2014, 135(5): 2727-2737.

[11] Kirby R, Amott K, Williams P T, et al. On the acoustic performance of rectangular splitter silencers in the presence of mean flow[J]. Journal of Sound and Vibration, 2014, 333(24):6295-6311.

[12] Kirby R. A comparison between analytic and numerical methods for modelling automotive dissipative silencers with mean flow[J]. Journal of Sound and Vibration, 2009, 325(3):565-582.

[13] Fang Z, Ji Z L. Acoustic attenuation analysis of expansion chambers with extended inlet/outlet[J]. Noise Control Engineering Journal, 2013, 61(2):240-249.

[14] Fang Z, Ji Z L. Numerical mode matching approach for acoustic attenuation predictions of double-chamber perforated tube dissipative silencers with mean flow[J]. Journal of Computational Acoustics, 2014, 22(2):1450004-1-1450004-15.

[15] Kirkup S M. Fortran codes for computing the discrete Helmholtz integral operators[J]. Advances in Computational Mathematics, 1998, 9(3/4): 391-409.

[16] Krishnasamy G, Schmerr L W, Rudolphi T J, et al. Hypersingular boundary integral equations: some applications in acoustic and elastic wave scattering[J]. Journal of Applied Mechanics, 1990, 57(2): 404-414.

[17] 方智,季振林. 穿孔管阻性消聲器橫向模態(tài)和聲學特性計算與分析[J]. 振動與沖擊, 2014, 33 (7):138-146.

FANG Zhi, JI Zhen-lin. Transversal modes and acoustic attenuation performance of a perforated tube dissipative silencer[J]. Journal of Vibration and Shock, 2014, 33(7): 138-146.

[18] 方智,季振林. 穿孔管消聲器橫截面模態(tài)及消聲特性的有限元分析[J]. 振動與沖擊, 2012, 31(17): 190-194.

FANG Zhi, JI Zhen-lin. Finite element analysis of cross-section modes and acoustic attenuation characteristics of a perforated tube silencer[J]. Journal of Vibration and Shock, 2012, 31(17):190-194.

基金項目：國家自然科學基金(11174065)

收稿日期：2014-11-12修改稿收到日期：2015-01-20

通信作者季振林 男，教授，博士生導師，1965年生

中圖分類號:TB535

文獻標志碼：A

DOI:10.13465/j.cnki.jvs.2016.02.025

Transmission loss prediction for silencers by using the combination of boundary element method and numerical collocation approach

YANG Liang1, JI Zhen-lin1, WU T W2

(1. School of Power and Energy, Harbin Engineering University, Harbin 150001, China; 2. College of Mechanical Engineering, University of Kentucky, Lexington 40506, Lexington)

Abstract:A technique combining the boundary element method (BEM) and numerical collocation approach was proposed to calculate the transmission loss of silencers with arbitrary cross-section. The silencer was divided into several substructures. The BEM and the numerical collocation approach were employed combinedly to evaluate the impedance matrices of substructures with non-regular cross-section and axially uniform cross-section respectively, and then the impedance matrices for all substructures were achieved. The transmission loss predictions of several silencers agree well with the results by BEM. The computational efficiency was analyzed and compared with that of the traditional substructure method, the combined technique may save computational time efficiently.

Key words:silencer; transmission loss; boundary element method; numerical collocation approach

第一作者 楊亮 男，博士生，1989年9月生