計算輻射噪聲的面聲源和點聲源結合方法

劉寶， 胡金華， 程廣利

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

0 引言

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

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

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

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

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

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

1 波疊加原理

結構產生的聲場可以看作是由結構內部空間中聲源輻射聲場的疊加，此即波疊加原理。波疊加原理示意圖如圖1所示。

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

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

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

(1)

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

(2)

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

2 虛擬點聲源法

根據波疊加原理，將結構產生的聲場看作是由結構內部空間中點聲源輻射產生的聲場疊加。即當點聲源作為虛擬聲源時，q(r0)可以表示為內部有限個點聲源的組合，即

(3)

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

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

(4)

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

(5)

(6)

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

p=Hs，

(7)

(8)

(9)

(10)

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

u=Us，

(11)

s=U-1u，

(12)

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

p=HU-1u.

(13)

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

(14)

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

(15)

(16)

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

(17)

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

3 虛擬面聲源法

當虛擬聲源連續分布于結構內部的某個虛擬邊界面時，由波疊加基本原理，(1)式可以寫為

(18)

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

(19)

(20)

(21)

對同一個節點上的積分函數進行合并，采用Gauss積分計算，最后表示成矩陣形式如下：

p=Gq，

(22)

(23)

(24)

(25)

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

u=Qq，

(26)

p=GQ-1u.

(27)

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

(28)

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

(29)

(30)

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

4 不同聲源類型的計算結果比較

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

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

5 面聲源和點聲源聯合的虛擬聲源法

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

(31)

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

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

(32)

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

(33)

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

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

(34)

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

6 案例分析

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

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

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

算例2：針對無限大介質中的橫向振動球源，采用面聲源和點聲源相結合的虛擬聲源法計算其全頻段的聲輻射特性參數。設球源半徑為r，表面振速可以表示為v0cosφ. 橫向振動球體和虛擬聲源表面離散單元的數量與離散形式和脈動球源相同。以與z軸呈φ=45°、φ=90°、φ=135°處計算的無量綱化聲壓和虛擬面聲源法、虛擬點源法的結果、解析解進行對比說明，如圖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可知，本文方法和面聲源法的計算精度要高于點聲源發的結果，點聲源計算的誤差隨著頻率的升高逐漸增大較多。當頻率接近特征頻率時，采用虛擬面聲源法和點源法的結果偏離解析解較大，不能保證特征頻率處解的唯一性。而采用面聲源和點聲源相結合的虛擬聲源法所得結果在全頻域內與解析解吻合得很好，不會在特征頻率處發生偏差。

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

7 結論

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

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

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

參考文獻(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] 夏百戰，于德介，姚凌云. 二維多流體域耦合聲場的光滑有限元解法[J]. 聲學學報，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，FYFE 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，FAHNLINE 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，FAHNLINEJ 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，FAHNLINE 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] 向宇，黃玉盈，馬小強．求解二維結構—聲耦合問題的一種直接方法[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]．船舶力學，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]. 力學學報，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] 陳鴻洋，商德江，李琪，等. 聲場匹配波疊加法的水下結構聲輻射預報[J]. 聲學學報，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] 商德江，錢治文，何元安，等．基于聯合波疊加法的淺海信道下圓柱殼聲輻射研究[J]．物理學報，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)