計算輻射噪聲的面聲源和點聲源結(jié)合方法

劉寶， 胡金華， 程廣利

(海軍工程大學(xué) 電子工程學(xué)院， 湖北 武漢 430033)

0 引言

求解結(jié)構(gòu)的聲輻射特性是聲學(xué)界普遍關(guān)注的問題，在解析法的基礎(chǔ)上，對不規(guī)則結(jié)構(gòu)通常采用有限元法(FEM)[1-2]，邊界元法(BEM)[3-4]計算。FEM在分析結(jié)構(gòu)聲振特性時，需要對結(jié)構(gòu)和三維流體介質(zhì)進(jìn)行單元離散和變量插值，計算量龐大。相對于FEM，BEM具有自動降低求解維度和適用于無限域的優(yōu)點。但是采用BEM計算時，需要處理奇異性積分問題和對應(yīng)內(nèi)部Dirichlet特征頻率處解的唯一性問題，雖然理論上可以采用CHIEF法[5]、Burton-Miller法[6]等獲得全波數(shù)范圍內(nèi)唯一的解，但無疑增加了計算的復(fù)雜性和存儲空間，降低了求解效率。

為此，19世紀(jì)90年代，Koopmann等[7-9]提出了采用簡單源匹配的虛擬聲源方法。作為邊界元的一種有效的替代方法，將結(jié)構(gòu)產(chǎn)生的聲場看作是由置于該輻射體內(nèi)部有限多個虛擬點聲源產(chǎn)生的聲場疊加，虛擬聲源的源強可以通過匹配結(jié)構(gòu)表面上若干點的法向速度獲得。進(jìn)一步利用源強可以直接求解聲場任一點聲壓、輻射聲功率等聲學(xué)參數(shù)。后來，Koopmann等開發(fā)了采用虛擬聲源計算的軟件POWER，并出版了相關(guān)專著[10]。Benthien等[11]指出：虛擬聲源所在區(qū)域與結(jié)構(gòu)表面不重合，因而不需要處理奇異性問題，并且可提高計算效率。實際應(yīng)用中為了計算上簡單方便，一般將點聲源布置在一個虛擬表面上。然而Jeans等[12]和Wilton等[13]的研究表明，正是由于采用了一個封閉曲面作為虛擬聲源的配置區(qū)域，導(dǎo)致了在關(guān)于虛擬面內(nèi)部Dirichlet問題的特征頻率處，解非唯一。為此他們提出了采用復(fù)數(shù)形式的Burton-Miller型組合層勢法，該類方法雖然解決了特征頻率處解的唯一性問題，但是相對于單獨采用單層勢或雙層勢，其計算時間增加了50%左右。

為了保證全波數(shù)域內(nèi)解唯一，Xiang等[14]提出了基于復(fù)數(shù)矢徑的虛擬聲源方法，他們將計算半徑改為復(fù)數(shù)，相當(dāng)于在結(jié)構(gòu)動力系統(tǒng)中加入阻尼，則在整個波數(shù)范圍內(nèi)保證了解的唯一性。運算該技術(shù)，向宇等[15]求解了空腔、空穴等二維結(jié)構(gòu)的聲振耦合問題，并進(jìn)一步求解了三維軸對稱空穴的聲輻射問題，結(jié)果表明：基于復(fù)數(shù)矢徑的虛擬聲源方法具有很高的計算精度。夏雪寶等[16]根據(jù)波疊加原理，研究了聲壓、聲源強度和表面法向振速之間的關(guān)系，推導(dǎo)了結(jié)構(gòu)聲輻射阻矩陣，并進(jìn)一步提出了附加源虛擬聲源方法，在單極子虛擬聲源法的基礎(chǔ)上，在虛擬面內(nèi)部添加一定數(shù)量的附加源，獲得聲場全波數(shù)域內(nèi)聲輻射特性的唯一解，案例表明其計算效率和精度都優(yōu)于Burton-Miller型組合層勢法的結(jié)果[17]。文獻(xiàn)[18]則采用面聲源代替點聲源構(gòu)造結(jié)構(gòu)外聲場，以二維脈動和擺動空穴為例，驗證了采用面聲源作為虛擬聲源計算的精確性，同時面聲源對位置的靈敏度小于點聲源。

在應(yīng)用虛擬聲源方法求解聲學(xué)參量時，虛擬聲源數(shù)目和位置的選擇對于求解的精度和效率有很大影響。采用虛擬聲源方法計算時，在聲輻射阻抗矩陣非對角線上元素占優(yōu)的情況下[19]，容易造成聲輻射阻抗矩陣病態(tài)，所計算的結(jié)果存在較大的誤差。根據(jù)靈敏度分析，當(dāng)結(jié)構(gòu)內(nèi)部虛擬聲源位置不變時，結(jié)構(gòu)表面振速分布是聲輻射參數(shù)的敏感變量，即振速分布的微小擾動可能使得虛擬聲源方法的計算結(jié)果產(chǎn)生較大誤差。文獻(xiàn)[8]針對不同表面振速分布的柱狀聲輻射體，采用虛擬聲源方法進(jìn)行了數(shù)值計算后指出，虛擬聲源方法的誤差與振速分布、波數(shù)、表面節(jié)點、虛擬聲源位置和個數(shù)都有關(guān)系。因此，為了保證求解的精度，需要針對不同結(jié)構(gòu)的表面振速分布，調(diào)整虛擬聲源的個數(shù)和位置。

為了使得聲輻射阻抗矩陣對角線上元素占優(yōu)，Koopmann等[9]由Borgiotti等[20]采用奇異值分解方法，建立結(jié)構(gòu)表面振速和遠(yuǎn)場聲壓，將聲輻射阻抗矩陣接近于0的奇異值截斷，形成新的聲輻射阻抗矩陣，結(jié)果表明，采用新的聲輻射阻抗矩陣計算的結(jié)果在精度上高于原聲輻射阻抗矩陣的結(jié)果。Zellers[21]則是利用無網(wǎng)格技術(shù)，將單極子、偶極子、三極子置于結(jié)構(gòu)離散單元的中心處，獲得了自輻射阻抗的解析表達(dá)式，用來替代聲輻射阻抗對角線上的元素，計算結(jié)果表明，該方法可以保證對角線上的元素占優(yōu)，且在低頻時可以取得較好的計算精度。上述工作針對聲輻射阻抗矩陣，進(jìn)行了不同程度的改進(jìn)，在一定程度上減小了計算誤差，但沒有考慮同樣會導(dǎo)致計算誤差的虛擬聲源位置因素。Hwang[22]的研究表明，虛擬聲源所在面應(yīng)該盡量與結(jié)構(gòu)表面共形，并且在理論上給出了虛擬面與結(jié)構(gòu)實際表面的最小距離，但沒有給出最大距離的限制。因此當(dāng)虛擬面與結(jié)構(gòu)實際表面的距離較大時，計算的結(jié)果仍可能不精確。Gounot等[23]、Pavic[24]則采用遺傳算法搜索一定范圍的虛擬空間，將虛擬聲源在結(jié)構(gòu)表面的振速與結(jié)構(gòu)表面實際的振速相比較確定虛擬聲源的位置。當(dāng)虛擬聲源數(shù)量龐大時，計算效率低，且搜索的虛擬聲源位置有可能為局部最優(yōu)位置。陳鴻洋等[25]以水下帶帽圓柱殼為對象，由聲場匹配搜索等效聲源分布，獲得最小二乘意義下的最優(yōu)等效源位置，所得結(jié)果在較寬頻段內(nèi)具有較好的適應(yīng)性，且適用于半無限大聲場。Wu等[26]提出了根據(jù)體積速度匹配原則確定虛擬聲源位置的方法，研究表明當(dāng)虛擬聲源位于優(yōu)化位置時，在中高頻上也可以保證良好的計算精度。商德江等[27]將波疊加原理、簡正波理論和多物理場有限元方法相結(jié)合，提出了一種淺海信道下彈性結(jié)構(gòu)聲輻射快速預(yù)報的聯(lián)合算法，并從數(shù)值法結(jié)果和解析解結(jié)果兩方面驗證了聯(lián)合算法的準(zhǔn)確性和可靠性。

本文根據(jù)波疊加原理，利用面聲源和點聲源作為非共面的虛擬聲源，搜索與邊界條件相匹配的真實聲壓與結(jié)構(gòu)表面振速，提出了一種保證全頻域解唯一的虛擬聲源方法。該方法在結(jié)構(gòu)內(nèi)部同時布置虛擬面聲源和點聲源離散聲場，面聲源位于結(jié)構(gòu)內(nèi)部一個與結(jié)構(gòu)邊界共形的表面，用其匹配結(jié)構(gòu)表面的聲壓與振速，點聲源則位于虛擬表面內(nèi)部，由于與虛擬表面不共面，因而可以保證解在特征頻率處的唯一性。文中以脈動球源和橫向振動球源為例，首先比較了點聲源和面聲源作為虛擬聲源計算的精確性，進(jìn)一步將本文所提方法與解析法、虛擬點聲源法、虛擬面聲源法相比較，說明該方法的優(yōu)點。

1 波疊加原理

結(jié)構(gòu)產(chǎn)生的聲場可以看作是由結(jié)構(gòu)內(nèi)部空間中聲源輻射聲場的疊加，此即波疊加原理。波疊加原理示意圖如圖1所示。

圖1 波疊加原理示意圖Fig.1 Schematic diagram of wave superposition principle

圖1中：SR為以R為半徑的球面；S為結(jié)構(gòu)表面積；V為結(jié)構(gòu)外空間；Vs為結(jié)構(gòu)所占體積；rs為結(jié)構(gòu)中心到表面一點處距離；r0為結(jié)構(gòu)中心到虛擬聲源所在位置處距離；n為結(jié)構(gòu)外法線方向；q(r0)為r0位置處虛擬聲源強度；r為結(jié)構(gòu)中心到聲場內(nèi)一點距離；p(r)為r位置處聲壓。

如圖1所示，場點r處的聲壓p(r)可以看做是由所有虛擬聲源輻射聲場的線性疊加，從而表示為

(1)

式中：j為虛數(shù)，j2=-1;k表示波數(shù)；ρf為介質(zhì)密度；c為聲速；q(r0)為虛擬聲源分布在r0處的聲源強度；Gk(r,r0)為自由空間的格林函數(shù)，滿足

(2)

理論上可以證明(1)式與Helmholtz邊界積分方程等價。根據(jù)(1)式，可以推導(dǎo)出虛擬點聲源法和虛擬面聲源法。

2 虛擬點聲源法

根據(jù)波疊加原理，將結(jié)構(gòu)產(chǎn)生的聲場看作是由結(jié)構(gòu)內(nèi)部空間中點聲源輻射產(chǎn)生的聲場疊加。即當(dāng)點聲源作為虛擬聲源時，q(r0)可以表示為內(nèi)部有限個點聲源的組合，即

(3)

式中：m為虛擬聲源數(shù)目，m=1,2,3,…,M；s(rm)為對應(yīng)rm位置處的虛擬聲源強度。

(3)式代入(1)式可得

(4)

(4)式為單層勢形式的聲壓計算公式。在此基礎(chǔ)上，可以給出雙層勢和混合勢形式的聲壓計算公式，分別如(5)式和(6)式所示。

(5)

(6)

式中：γ=j/k；nm表示對應(yīng)于rm位置處的結(jié)構(gòu)外法線方向。(4)式、(5)式、(6)式表示成矩陣形式為

p=Hs，

(7)

(8)

(9)

(10)

(8)式～(10)式表示成矩陣形式為

u=Us，

(11)

s=U-1u，

(12)

(12)式代入到(7)式中，即可求解結(jié)構(gòu)表面聲壓為

p=HU-1u.

(13)

除了采用(13)式計算結(jié)構(gòu)的表面聲壓外，還可以采用最小二乘法建立虛擬聲源強度求解方程，進(jìn)一步求解表面聲壓。以單層勢為例進(jìn)行計算，由最小二乘法原理，(8)式建立泛函形式的方差表達(dá)式如(14)式所示。

(14)

式中：F1表示泛函的方差；nrn為對應(yīng)位置rn處的法線方向；vm為第m個結(jié)構(gòu)表面節(jié)點位置處的振速。展開(14)式可得

(15)

(16)

(15)式代入(16)式，可得

(17)

利用(17)式即可求得最小二乘意義下的虛擬聲源強度，進(jìn)一步求解表面振速。

3 虛擬面聲源法

當(dāng)虛擬聲源連續(xù)分布于結(jié)構(gòu)內(nèi)部的某個虛擬邊界面時，由波疊加基本原理，(1)式可以寫為

(18)

式中:α、β分別為結(jié)構(gòu)表面上一點位置和虛擬聲源所在位置的符號表示。由于采用面單元作為虛擬聲源，因此需要對虛擬聲源強度q(β)進(jìn)行變量插值，在此以二次等參單元為例，將虛擬邊界面劃分為Ne個單元，每個單元含有8個節(jié)點，共含有N個節(jié)點。于是單層勢形式的聲壓計算式可以表示為

(19)

(20)

(21)

對同一個節(jié)點上的積分函數(shù)進(jìn)行合并，采用Gauss積分計算，最后表示成矩陣形式如下：

p=Gq，

(22)

(23)

(24)

(25)

(23)式、(24)式、(25)式最終表示成矩陣形式為

u=Qq，

(26)

p=GQ-1u.

(27)

同樣地，可以采用最小二乘法建立虛擬聲源強度的求解方程。以單層勢為例進(jìn)行計算，由最小二乘法原理，(26)式建立泛函形式的方差表達(dá)式如(28)式所示。

(28)

對(28)式進(jìn)行展開，可得

(29)

(30)

由(30)式即可計算出虛擬聲源強度，進(jìn)而求解結(jié)構(gòu)的表面聲壓。在引入基虛擬點聲源法和虛擬面聲源法后，本文繼續(xù)對比上述兩種虛擬聲源法的精度和計算效率。

4 不同聲源類型的計算結(jié)果比較

為了比較虛擬點聲源法和虛擬面聲源法各自的特點，以平面圓形空穴為例，其半徑為r，分別進(jìn)行脈動和橫向振動，其邊界線速度分別為v0和v0cosθ(θ是x軸方向與橫向振動方向的夾角)。以空穴中心為圓心，半徑為0.7r的圓作為虛擬邊界線，在虛擬邊界線上等間距地選擇8個節(jié)點作為點聲源。同樣地將該虛擬邊界線離散為4個二次等參單元作為面聲源，共含有8個節(jié)點，如圖2所示。設(shè)速度v0=1，r=1，采用上述兩種虛擬聲源匹配圓形邊界線上的速度，并將所得結(jié)果進(jìn)行比較，圖3、圖4、圖5分別為在頻率kr=π rad/s，kr=2π rad/s，kr=3π rad/s時的結(jié)果。

圖2 虛擬聲源布置模型Fig.2 Layout model of virtual sound sources

圖3 kr=π rad/s時速度匹配圖Fig.3 Velocity matching diagram for kr=π rad/s

圖4 kr=2π rad/s時速度匹配圖Fig.4 Velocity matching diagram for kr=2π rad/s

圖5 kr=3π rad/s時速度匹配圖Fig.5 Velocity matching diagram for kr=3π rad/s

由圖3～圖5可知，當(dāng)頻率從kr=π rad/s上升到頻率kr=3π rad/s時，采用虛擬點聲源法計算的誤差逐漸增大，不能較好地匹配結(jié)構(gòu)表面的振速。同時可以明顯看出，采用虛擬面聲源法的計算結(jié)果在頻率kr=π rad/s、kr=2π rad/s、kr=3π rad/s時匹配的結(jié)果都很好，沒有產(chǎn)生較大的誤差。因此，虛擬面聲源法計算的精度要好于虛擬點聲源法的結(jié)果。當(dāng)然，在計算效率上，虛擬點聲源法的計算效率要比虛擬面聲源法快很多，這是由于不同于虛擬面聲源法需要進(jìn)行積分運算，虛擬點聲源法不需要積分就可以直接求解聲阻抗矩陣，因此節(jié)約了大量時間。

5 面聲源和點聲源聯(lián)合的虛擬聲源法

在匹配聲場中的聲壓與振速方面，由于虛擬面聲源法的精度要高于虛擬點聲源法，故采用面聲源匹配聲場可以獲得更好地求解精度。由于面聲源作為虛擬聲源時通常位于同一個面上，因此在特征頻率處，會出現(xiàn)解的唯一性問題。為此在采用虛擬面聲源的基礎(chǔ)上，在虛擬面內(nèi)部布置若干虛擬點聲源，用來保證解唯一。則單層勢形式的聲壓可以表示為

(31)

該方法的基本思想是采用虛擬面聲源匹配結(jié)構(gòu)表面聲壓與振速，而點聲源不與面聲源共面，以保證解的唯一性，因此將該方法稱為面聲源和點聲源相結(jié)合的虛擬聲源法?？紤]到點聲源的存在是為了保證特征頻率處的解具有唯一性，但過多虛擬點聲源的存在會降低求解效率且影響計算精度。因此，要求點聲源的數(shù)量盡可能少，同時在結(jié)構(gòu)內(nèi)具有位置合適，這也是該方法的難點所在。即如何有效地確定虛擬點聲源的數(shù)目與位置。經(jīng)過實際運算，當(dāng)虛擬點聲源數(shù)目為虛擬面聲源的1/10左右且點聲源的位置均勻分布在結(jié)構(gòu)內(nèi)部時可以較好地保證解在全波數(shù)內(nèi)唯一，并且具有較高的計算精度。(31)式表示成矩陣形式為

p=[G?H][q?s].

(32)

由(31)式可以很自然地推導(dǎo)出振速的計算公式為

(33)

(33)式表示成矩陣形式為

v=[U?Q][q?s].

(34)

(34)式計算出的虛擬聲源強度代入(32)式中，即可計算出結(jié)構(gòu)的表面聲壓。

6 案例分析

算例1：針對無限大介質(zhì)中的脈動球源，采用面聲源和點聲源相結(jié)合的方法計算其全頻段的聲輻射特性參數(shù)。設(shè)球源半徑為r，結(jié)構(gòu)表面法向振速為v0，介質(zhì)聲速為c，介質(zhì)密度為ρf. 球體表面采用104個單元的8節(jié)點二次等參單元進(jìn)行離散插值。以球心為原點，采用半徑為0.7r的球面作為虛擬面，即放縮系數(shù)為0.7，應(yīng)用相同數(shù)目的二次等參單元對該虛擬面進(jìn)行離散插值，這些單元即為虛擬面源。在虛擬表面的內(nèi)部均勻分布10個點作為虛擬點聲源。采用Matlab編程進(jìn)行計算?？紤]到脈動球源表面各點理論上聲壓相同，以與z軸呈φ=45°(φ為z軸方向與橫向振動方向的夾角)處計算的無量綱化聲壓為例，將本文方法結(jié)果和虛擬面聲源法、虛擬點聲源法的結(jié)果、解析解進(jìn)行對比說明，如圖6所示。

圖6 脈動球源φ=45°處無量綱化聲壓比較圖Fig.6 Comparison of dimensionless sound pressures at pulsating spherical source for φ=45°

從圖6看出，點聲源計算的誤差隨著頻率的升高逐漸增大較多，在計算精度上本文方法和面聲源法的結(jié)果均優(yōu)于點聲源法的結(jié)果。當(dāng)頻率接近關(guān)于虛擬面內(nèi)部Dirichlet問題處的特征頻率時，即kr=4.48， 8.97時，采用點聲源法和面聲源法的計算結(jié)果偏離解析解較大，說明虛擬點聲源法和面聲源法均不能保證對應(yīng)特征頻率處聲輻射解的唯一性。采用面聲源和點聲源相結(jié)合的方法計算時，所得結(jié)果的實部和虛部與解析解在全部頻域內(nèi)吻合得很好，整體上與解析解一致，不會產(chǎn)生較大的偏差，說明該方法具有很好的計算精度。

算例2：針對無限大介質(zhì)中的橫向振動球源，采用面聲源和點聲源相結(jié)合的虛擬聲源法計算其全頻段的聲輻射特性參數(shù)。設(shè)球源半徑為r，表面振速可以表示為v0cosφ. 橫向振動球體和虛擬聲源表面離散單元的數(shù)量與離散形式和脈動球源相同。以與z軸呈φ=45°、φ=90°、φ=135°處計算的無量綱化聲壓和虛擬面聲源法、虛擬點源法的結(jié)果、解析解進(jìn)行對比說明，如圖7～圖9所示。

圖7 橫向振動球源φ=45°處無量綱化聲壓比較圖Fig.7 Comparison of dimensionless sound pressures at laterally oscillating spherical source for φ=45°

圖8 橫向振動球源φ=90°處無量綱化聲壓比較圖Fig.8 Comparison of dimensionless sound pressures at laterally oscillating spherical source for φ=90°

圖9 橫向振動球源φ=135°處無量綱化聲壓比較圖Fig.9 Comparison of dimensionless sound pressures at laterally oscillating spherical source for φ=135°

從圖7～圖9可知，本文方法和面聲源法的計算精度要高于點聲源發(fā)的結(jié)果，點聲源計算的誤差隨著頻率的升高逐漸增大較多。當(dāng)頻率接近特征頻率時，采用虛擬面聲源法和點源法的結(jié)果偏離解析解較大，不能保證特征頻率處解的唯一性。而采用面聲源和點聲源相結(jié)合的虛擬聲源法所得結(jié)果在全頻域內(nèi)與解析解吻合得很好，不會在特征頻率處發(fā)生偏差。

在橫向振動球源的計算效率上，采用本文方法、虛擬面聲源法和虛擬點聲源法的單頻計算時間分別為170.4 s、170.5 s和1.5 s，從中可以看出，點聲源法的計算效率要遠(yuǎn)高于采用面聲源計算聲輻射參數(shù)的方法，這是由于采用面聲源運算需要進(jìn)行矩陣組裝和積分運算，因而消耗了大量時間。

7 結(jié)論

根據(jù)波疊加原理，利用面聲源和點聲源作為非共面的虛擬聲源，搜索與邊界條件相匹配的真實聲壓與結(jié)構(gòu)表面振速，本文提出了一種保證全頻域解唯一的虛擬聲源計算方法。該方法在結(jié)構(gòu)內(nèi)部同時布置虛擬面聲源和點聲源離散聲場：面聲源位于結(jié)構(gòu)內(nèi)部一個與結(jié)構(gòu)邊界共形的表面，用其匹配結(jié)構(gòu)表面的聲壓與振速；點聲源則位于虛擬表面內(nèi)部，由于與虛擬表面不共面，因而可以保證解在特征頻率處的唯一性。案例分析結(jié)果表明：

1)采用面聲源匹配結(jié)構(gòu)表面振速的精度要優(yōu)于點聲源的結(jié)果；在頻率較高時，采用虛擬點聲源計算的誤差過大，因而虛擬點聲源法不適合處理高頻情形；為保證精度，可以采用虛擬面聲源法處理高頻情形，但相對于虛擬點聲源法運算量成倍增加，這是由于采用面聲源計算時需要進(jìn)行積分運算，因此采用點聲源的計算精度差但效率高。

2)與BEM相比，面聲源和點聲源相結(jié)合的虛擬聲源法中，由于虛擬源位于結(jié)構(gòu)內(nèi)部，因此不需要處理奇異性問題，簡化了計算過程，提高了計算效率；相對于單獨采用點聲源、面聲源作為虛擬聲源方法，聯(lián)合虛擬聲源方法可以保證解在全波數(shù)范圍內(nèi)的唯一性，且采用較少的單元就可以獲得很好的計算精度。

參考文獻(xiàn)(References)

[1] SUN F X, LIU C, CHEN Y M. An improved interpolating element-free Galerkin method based on nonsingular weight functions[J]. Mathematical Problems in Engineering, 2014(8)：249-269.

[2] 夏百戰(zhàn)，于德介，姚凌云. 二維多流體域耦合聲場的光滑有限元解法[J]. 聲學(xué)學(xué)報，2012，37(6)：601-609.

XIA B Z，YU D J，YAO L Y. A smoothed finite element method for two-dimensional coupling acoustic fields in multi-fluid domain[J]. Chinese Journal of Acoustics，2012，37(6)：601-609. (in Chinese)

[3] LIU J Y, RUDOLPHI T J. New identities for fundamental solutions and their applications to non-singular boundary element formulations[J]. Computational Mechanics，1999，24(4)：286-292.

[4] ATALLA N，BEMHARD R J．Review of numerical solutions for low-frequency structural-acoustic problems[J]．Applied Acoustics，1994，43(3)：271-294.

[5] COYETTE J P，F(xiàn)YFE K R．An improved formulation for acoustic eigenmode extraction from boundary element methods[J]．Journal of Vibration and Acoustics，1990，112(3)：392-397.

[6] BURTON A J，MILLER G F．The application of integral equation methods to the numerical solution of some exterior boundary-value problems[J]．Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences ，1971，323(1553)：201-210.

[7] KOOPMANN G H，SONG L，F(xiàn)AHNLINE J B．A method for computing acoustic fields based on the principle of wave superposition[J]．Journal of the Acoustical Society of America，1989，86 (6)：2433-2438．

[8] SONG L，KOOPMANN G H，F(xiàn)AHNLINEJ B．Numerical errors associated with the method of superposition for computing acoustic fields[J]．Journal of the Acoustical Society of America，1991，89(6)：2625-2633．

[9] KOOPMANN G H，BELEGUNDU A D．Tuning a wine glass via material tailoring an application of a method for optimal acoustic design[J]．Journal of Sound and Vibration，2001，239 (4)：665-678．

[10] KOOPMANN G H，F(xiàn)AHNLINE J B．Designing quiet structures: a sound power minimization approach [M]．San Diego, CA ,US: Academic Press， 1997．

[11] BENTHIEN G W，SCHENCK H A．Nonexistence and nonuniqueness problems associated with integral equation methods in acoustic[J]．Computers & Structures，1997，65(3)：295-305．

[12] JEANS R A，MATHEWS I C．The wave superposition method as a robust technique for computing acoustic fields[J]．Journal of the Acoustical Society of America，1992，92(2)：1156-1166．

[13] WILTON D T，MATHEWS I C，JEANS R A．A clarification of nonexistence problems with superposition method[J].Journal of the Acoustical Society of America，1993，94(3)：1676-1680．

[14] XIANG Y，LU J，HUANG Y Y．A fast wave superposition spectral method with complex radius vector combined with two-dimensional fast Fourier transform algorithm for acoustic radiation of axisymmetric bodies[J]．Journal of Sound & Vibration，2012，331(6)：1441-1454．

[15] 向宇，黃玉盈，馬小強．求解二維結(jié)構(gòu)—聲耦合問題的一種直接方法[J]．振動與沖擊，2003，22(4)：40-44,31．

XIANG Y，HUANG Y Y，MA X Q．Direct solution to 2-D sound-structure interaction problems[J]．Journal of Vibration and Shock，2003，22(4)：40-44,31．(in Chinese)

[16] 夏雪寶，向陽，王校青，等．基于波疊加法的聲輻射阻計算研究[J]．船舶力學(xué)，2015，19(1/2)：206-214．

XIA X B，XIANG Y，WANG X Q，et al．Research on the acoustic radiation resistance calculation based on wave superposition[J]．Journal of Ship Mechanics，2015，19(1/2)：206-214．(in Chinese)

[17] 夏雪寶，向陽．基于附加源波疊加法的聲輻射計算研究[J]．振動與沖擊，2015，34(1)：104-109．

XIA X B，XIANG Y．Acoustic radiation calculation based on additional sources wave superposition method[J]．Journal of Vibration and Shock，2015，34(1)：104-109．(in Chinese)

[18] 向宇, 黃玉盈. 伸縮虛擬邊界元法解二維Helmholtz外問題[J]. 力學(xué)學(xué)報，2003，35(3)：272-279．

XIANG Y，HUANG Y Y．The expanding-contracting virtual boundary element method for 2D-Helmholtz exterior problems[J]．Acta Mechanica Sinica，2003，35(3)：272-279．(in Chinese)

[19] GOLUB G H，VAN LOAN C F．Matrix computations[M]．Baltimore, MD, US: Johns Hopkins University Press，1983:11-29．

[20] BORGIOTTI G V，SARKISSIAN A，Williams E G，et al．Conformal generalized near-field acoustic holography for axisymmetric geometries[J]．Journal of the Acoustical Society of America，1990，88(1)：199-209.

[21] ZELLERS B C. An acoustic superposition method for computing structural radiation in spatially digitized domains[D]．State College, PA, US：The Pennsylvania State University，2006.

[22] HWANG J Y，CHANG S C. A retracted boundary integral equation for exterior acoustic problem with unique solution for all wave numbers[J]. Journal of the Acoustical Society of America，1991，90(2)：167-180.

[23] GOUNOTY J R，MUSAFIR R E. Genetic algorithms: a global search tool to find optimal equivalent source sets[J]. Journal of Sound & Vibration，2009，322(1/2)：282-298.

[24] PAVIC G. A technique for the computation of sound radiation by vibrating bodies using multipole substitute sources[J]. Acta Acustica United with Acustica，2009，92(1)：112-126.

[25] 陳鴻洋，商德江，李琪，等. 聲場匹配波疊加法的水下結(jié)構(gòu)聲輻射預(yù)報[J]. 聲學(xué)學(xué)報，2013，38(2)：137-146.

CHEN H Y，SHANG D J，LI Q，et al．Sound radiation prediction for underwater structure by field-matching wave superposition method[J]. Chinese Journal of Acoustics，2013，38(2)：137-146. (in Chinese)

[26] WU S W，XIANG Y．Location optimization of monopole equivalent sources in wave superposition method[J]．International Journal of Acoustics and Vibration，2018，23(2)：254-263.

[27] 商德江，錢治文，何元安，等．基于聯(lián)合波疊加法的淺海信道下圓柱殼聲輻射研究[J]．物理學(xué)報，2018，67(8)：125-138.

SHANG D J，QIAN Z W，HE Y A，et al．Sound radiation of cylinder in shallow water investigated by combined wave superposition method[J]．Acta Physica Sinica，2018，67(8)：125-138.(in Chinese)