基于波疊加近場聲全息射線波函數的信息補償

陳巖豪 向宇 石梓玉

摘? 要：基于波疊加法的近場聲全息技術因其在適應性和數值計算上的優勢，近年來已經被廣泛應用于聲源識別、定位和聲場分析.在基于波疊加法的近場聲全息技術中，以射線波函數為波疊加法的積分核函數可有效改善系統矩陣的病態性，提高聲場計算的數值穩定性.但由于射線波函數在低波數處的指向性過強，聲場信息主要集中在射線波函數的主指向處，導致了其他方位的聲場信息缺失.針對以上問題，提出了一種射線波函數信息的補償方法，即在射線波函數中加入一定比例的單極子球波函數對重建聲場進行信息補償.數值算例結果表明：補償信息后的射線波函數不僅保證了聲場重建過程中數值計算的穩定性，而且彌補了單一射線波函數在低波數處的缺陷，進一步提高了聲場的重建精度.

關鍵詞：近場聲全息;波疊加法;射線波函數;信息補償 ;聲場重建

中圖分類號：TB52? ? ? ? ? ? DOI：10.16375/j.cnki.cn45-1395/t.2021.03.002

0? ? 引言

基于波疊加法的近場聲全息理論[1-2]因其算法簡單、計算效率較高，被廣泛應用于聲源識別、定位和聲場重建，其基本思想是：聲源產生的聲場可由該聲源內部一系列虛擬等效源產生的聲場疊加替代，而虛擬等效源的源強可以通過匹配全息面的測量聲學信息得到[3].波疊加法一般采用的是只與距離有關的球形單極子函數，其構成的系統矩陣通常為一病態矩陣，影響了數值計算的穩定性.對于弱病態的系統矩陣，通常采用正則化方法[4-5]過濾或截斷較小奇異值來抑制測量誤差的影響，提高數值計算的穩定性.為了獲取更高的計算精度，賀春東等[6-8]對近場聲全息的正則化方法作了進一步的優化.然而，單一的正則化方法無法從根本上改變系統矩陣的固有缺陷，要確保重建結果與實際聲場信息相吻合，有時還需要對等效源的位置作進一步選擇[9]和優化[10-11].無論何種優化方式，計算過程中都是通過改變系統矩陣形態來對重建結果施加影響.Song等[12]對系統矩陣研究后指出，要獲得較高的計算精度，系統矩陣需滿足對角占優且盡量對稱.

為了使系統矩陣滿足上述條件，文獻[13]提出以格林函數的[n（n=1， 2， 3， …）]階方向導數為積分核函數的射線波疊加法.與傳統波疊加法相比，射線波疊加法有效地改善了系統矩陣的病態性，并獲得了較為穩定的計算結果.但是，射線波疊加法在低波數處聲場重建效果并不理想.其原因在于射線波函數在低波數處的指向性過強，其疊加聲場信息將主要集中在各個射線波等效源的主指向處，而在非主指向處出現了信息缺失.

為了克服以上問題，本文在射線波函數中添加一定權重的球形單極子函數來彌補聲場缺失的信息.研究發現，經過信息補償后的射線波函數不僅保證了系統矩陣趨于良態，而且提高了聲場的重建精度.最后通過脈動球聲源以及雙脈動球干涉聲源驗證了本文方法的優越性.

1? ? 波疊加法及射線波函數的性質分析

將虛擬等效源布置在輻射聲源內部的封閉曲面[SE]上，可得輻射體外部區域任意一點聲壓[p（r）][3]：

[pr=SEqrEgr，rEdSE]? ? ? ? ? ? ? （1）

其中：[q（rE）]為待求等效虛擬源強，[g（r，rE）=eikR/（4πR）]，[i=-1]，[k]為波數，[R=r-rE].若將等效源面[SE]網格離散成N個節點，同時在節點處配置[N]個單極子聲源，得到離散后的積分方程：

[p（r）=i=1Nq（rEi）g（r，rEi）]? ? ? ? ? ? ? ?（2）

由式（2）可知，聲源外部的輻射聲場可由內部一系列虛擬聲源輻射聲場疊加得到.若全息面上存在[N]個全息測點，可將式（2）轉化為矩陣形式如下：

[PH=GHQ]? ? ? ? ? ? ? ? ? ? ? ? ? ?（3）

式中：[PH]為測點聲壓列向量，[GH]為全息面測點與等效源之間的系統矩陣，[Q]為待求等效源源強向量.通過對[GH]求逆即可計算得到[N]個虛擬等效源的源強，然后將求得的源強代入式（2）即可求解聲源面外部任意位置的聲壓[p（r）].但由于格林函數[g（r，rE）]是一個球波函數，其構成的系統矩陣[GH]往往是一個病態矩陣，數值求逆時會放大測量誤差的影響，可能導致重建結果失真.為了克服以上問題，文獻[13]構造了具有主指向的射線波函數來替代傳統波疊加法的單極子球波函數，極大地提高了聲場重建的數值穩定性.

2? ? 射線波函數的特點及缺陷

Ochmann[14-15]指出只要滿足Helmholtz方程和Sommerfield輻射條件的解析函數均可以作為波疊加法的積分核函數.基于這一研究思路，可以通過改變積分核函數的方式提高聲場重建的數值穩定性[16-17].在文獻[13]中構造了如下的射線波函數：

[Ki（r，rE）=?（i）g（r，rE）?l（i），i=1，2，3，…]? ? ? （4）

式中：[Ki（r，rE）]為具有指向性的射線波函數，[l]為其主指向，[i]為求導階數.當[i=1]且[l]為等效源面法向量n時，[?g（r，rE）?n]為常用的雙層勢函數.經過計算發現，隨著求導階數的增加，射線波函數的指向性逐漸增強，如圖1所示.

若將上述的N個射線波等效源的主指向與全息面測點一一對應，系統矩陣如下：

[Ki（r，rE）=Ki（r1，rE1）Ki（r1，rE2）…Ki（r1，rEn）Ki（r2，rE1）Ki（r2，rE2）…Ki（r2，rEn）???Ki（rN，rE1）Ki（rN，rE2）…Ki（rN，rEn）N×N]（5）

由于在[rEi]處的射線波等效源在其對應的全息測點[ri]處的聲波激勵最大，其構成的系統矩陣[Ki（r，rE）]呈現主對角占優且近似對稱.與傳統單極子等效源構成的系統矩陣相比，[Ki（r，rE）]具有較好的數值穩定性，不會因為很小的測量誤差導致計算結果失真[12].

經過進一步的研究發現，當射線波函數指向性過強時，系統矩陣會出現以下情況：

[Ki（rn，rEn）?Ki（rn，rEm），n， m=1， 2， 3， …且 m≠n]? ?（6）

即系統矩陣主對角元素遠大于非主對角元素，此時盡管系統矩陣的數值穩定性較好，但重建聲場與實際聲場信息不匹配.原因如圖2所示.

指向性過強的射線波等效源輻射的聲場信息主要集中于與之對應的全息面測點附近，而在遠離全息面測點的位置將會出現信息缺失.因此，選取的射線波求導階數不宜過高，一般為2～3階.其次，文獻[13]的數值算例指出，射線波函數的重建效果與波數也有一定的關系，其在低波數處重建效果不太理想.經分析發現，如圖3所示，在相同的求導階數下，射線波函數在低波數處的指向性較強，而在波數較大時，指向性強度相對較弱.說明了波數較小時，由于射線波等效源的指向性過強，導致了聲場重建信息缺失.

3? ? 射線波函數的信息補償

由于在低波數處射線波函數的聲場信息主要集中在各個射線波的主指向性處，導致其他位置的聲場信息缺失，因此，需要采用一定的策略對缺失的聲場信息進行補償.單極子球面波函數是一種全空間覆蓋且在各個方向均勻衰減的波函數，因此，本文考慮在射線波函數中添加一定權重的單極子波函數對聲場信息進行補償，表達式如下：

[K（i）a（r，rE）=ag（r，rE）+?（i）g（r，rE）?l（i），i=0， 1， 2， 3，…]（7）

其中：i為求導階數，[a]為組合權重系數.通過選取合適的權重系數就可以保證射線波函數構成的系統矩陣趨于良態的同時，提高低波數處聲場的重建精度.為了便于對比，給出不同權重系數時的射線波函數[K（3）a（r，rE）]指向性，如圖4所示.

由圖4可以看出隨著單極子波函數權重的增加，射線波函數的指向性逐漸減弱，因此，可以通過調整參數[a]的大小，補償合適的球面波信息，在保證改善低波數處聲場重建精度的同時，提高聲場重建的數值穩定性.

4? ? 數值算例與分析

為了驗證加入單極子球波函數后射線波函數的數值穩定性以及對低波數處的改善效果，選取文獻[13]中數值計算穩定較好的三階射線波函數[?（3）g（r，rE）?l（3）]進行組合對比驗證，其中重建聲源模型為脈動球聲源和雙脈動球干涉聲源.

4.1? ?脈動球聲源重建效果對比

首先采用脈動球聲源進行驗證，半徑為[Rs]的脈動球聲源在距離圓心為[r]處聲壓如下[18]：

[p（r）ρ0cv0=RsrikRs1+ikRse-ik（r-Rs）]? ? ? ? ? ? （8）

其中：[ρ0]為介質的平均密度，[v0]為脈動球表面質點法向振速幅值，[c]為聲波在介質中的傳播速度，[i=-1]，? [k]為波數.設定[ρ0=1.21] [kg/m?]，[v0=1] [m/s]，[c=340] [m/s].脈動球聲源半徑[Rs]=1.0 m.全息面半徑為1.5 m的球面，按照緯度間隔[π/11]，經度間隔[2π/11]布置全息測點，除去兩極，共100個節點.等效源面半徑? ? [rE=0.2] [m]，其等效源點經緯度分布方式與全息測點? 相同.

根據上述模型布置，分別采用常用的單極子波函數[g（r，rE）]、雙層勢波函數[?g（r，rE）?n]、三階射線波函數[?（3）g（r，rE）?l（3）]以及信息補償后的射線波函數[K（3）a（r，rE）]進行聲場重建，其中經過大量算例驗證[a=1～10]時系統矩陣較為穩定，而且改善效果較為明顯.本文選取[a=8].在數值計算中，系統矩陣條件數大小可以衡量矩陣的病態程度，所以首先對比在不同波數下各個波函數構成的系統矩陣條件數如圖5所示.

圖5表明，單極子波函數構成的系統矩陣條件數較大，采用雙層勢波函數[?g（r，rE）?n]可以降低系統矩陣條件數，而采用三階射線波函數[?（3）g（r，rE）?l（3）]和[K（3）a（r，rE）]降低效果更加明顯，說明射線波函數和組合后的射線波函數均能夠有效改善系統矩陣的病態性，提高聲場計算的數值穩定性.

為了驗證以上各個波函數在不同波數的重建精度和穩定性，在全息測點聲壓中添加40 dB的噪聲，并利用以上各個波函數對脈動球表面聲壓實部和虛部進行重建.不同波函數的重建效果如圖6所示.

圖6（a）可以看出，單極子波函數[g（r，rE）]重建數值解完全偏離了解析解，圖6（b）中雙層勢波函數[?g（r，rE）?n]的重建效果有所改善，但是部分波數處吻合程度較差;采用三階射線波函數[?（3）g（r，rE）?l（3）]在低波數處重建數值解偏離了解析解，在波數較大時重建效果較好;信息補償后的射線波函數[K（3）8（r，rE）]在各個波數處均有較好的重建效果.綜合條件數對比和重建結果來看，補償信息后的射線波函數既改善了系統矩陣的病態性，又彌補了單一射線波函數在低波數處的缺陷.

4.2? ?雙脈動球干涉聲源節點聲壓重建

在實際工程中，聲場的分布通常較為復雜且無解析解可循.本文根據替代法[19]的原理，利用雙脈動球形成的干涉聲場代替實際聲場，其中2個脈動球中心坐標分別為（0.2， 0， 0）和（-0.2， 0， 0），半徑均為0.4 m，振速為[v=1] [m/s].全息面為半徑[RH=1.0] m的球形，節點分布緯度間隔為[π12]，經度間隔為[2π12]，除去兩極節點，共121個節點.等效源布置面半徑為[r=0.2] m，等效源點分布與全息測點相同.聲場重建節點均勻分布在一半徑為0.8 m的圓形之上，共121個重建節點.其二維示意圖如圖7? ? ?所示.

依舊采用上述單極子波函數[g（r，rE）]、雙層勢波函數[?g（r，rE）?n]、三階射線波函數[?（3）g（r，rE）?l（3）]以及信息補償后的射線波函數[K（3）8（r，rE）]對聲場進行重建.為了驗證補償后射線波函數在低波數處的改善效果，設定波數為[k=0.5]，并對全息測點聲壓數據添加40 dB的噪聲.對比重建節點解析解聲壓以及重建數值解聲壓如圖8所示.

根據圖8可以看出單極子波函數[g（r，rE）]、雙層勢波函數[?g（r，rE）?n]、三階射線波函數[?（3）g（r，rE）?l（3）]在各個節點的重建聲壓均偏離了解析解，而補償信息后的射線波函數[K（3）8（r，rE）]在節點處重建結果與解析解的吻合程度較好，說明低波數處補償信息后的射線波函數在空間不同位置處均有較好的重建效果.根據表1中各個系統矩陣條件數對比可以發現，單極子波函數[g（r，rE）]和雙層勢波函數[?g（r，rE）?n]由于條件數較大，在利用MATLAB計算過程中已經出現提示求逆計算結果可能不準確.而相比單極子波函數[g（r，rE）]和雙層勢波函數[?g（r，rE）?n]，雖然三階射線波函數[?（3）g（r，rE）?l（3）]條件數最小，但是重建效果依舊不太理想，結合射線波函數缺陷的分析，說明在低波數處利用射線波函數[?（3）g（r，rE）?l（3）]重建聲場時出現了信息缺失.而補償信息后的射線波函數[K（3）8（r，rE）]不僅降低了系統矩陣條件數，改善了系統矩陣的病態性，而且保證各個節點處聲場的重建精度，因此，足以說明本文方法的優越性.

5? ? 結論

本文綜合考慮了系統矩陣以及各個波函數的特性，提出了一種射線波函數信息補償法，即在射線波函數中添加一定權重的單極子球波函數對聲場信息進行補償.數值算例結果表明：信息補償后的射線波函數顯著降低了系統矩陣的條件數，既保證了聲場重建的數值穩定性，又彌補了單一射線波函數在低波數的缺陷，提高了聲場重建的計算精度.在實際工程應用當中，結構體振動輻射的聲場更加復雜，經驗性地選取組合參數可能無法滿足重建精度的要求，因此，如何穩定選取最優的組合參數使射線波函數能夠適應于任意聲場的重建是下一步的研究重點.

參考文獻

[1]? ?KOOPMANN G H，SONG L M，FAHNLINE J B.A method for computing acoustic fields based on the principle of wave superposition[J].The Journal of the Acoustical Society of America，1989，86（6）：2433-2438.

[2]? ? LEE S. Review：the use of equivalent source method in computational acoustics[J].The Journal of Computational Acoustics，2017，25（1）：1630001.

[3]? ?JEANS R，MATHEWS I C.The wave superposition method as a robust technique for computing acoustic fields[J].The Journal of the Acoustical Society of America，1992，92（2）：1156-1166.

[4]? ?WILLIAMS E G.Regularization methods for near-field acoustical holography[J].The Journal of the Acoustical Society of America，2001，110（4）：1976-1988.

[5]? ? 陳達亮， 舒歌群， 衛海橋. 近場聲全息正則化方法比較[J]. 天津大學學報， 2008，41（6）： 696-702.

[6]? ? 賀春東， 張永斌， 畢傳興，等. 分部優化正則化在近場聲全息技術中的應用[J]. 聲學學報， 2011，36（4）： 419-426.

[7]? ?CHELLIAH K，RAMAN G G，MUEHLEISEN R T. Enhanced nearfield acoustic holography for larger distances of reconstructions using fixed parameter Tikhonov regularization[J].The Journal of the Acoustical Society of America，2016，140（1）：114-120.

[8]? ? XIAO Y. A new method for determining optimal regularization parameter in near-field acoustic holography[J].Shock and Vibration，2018，2018： 1-13.

[9]? ? BAI M R，CHEN C C，LIN J H.On optimal retreat distance for the equivalent source method-based nearfield acoustical holography[J].The Journal of the Acoustical Society of America，2011，129（3）：1407-1416.

[10] 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.

[11] ? GOUNOT Y J R，MUSAFIR R E.Genetic algorithms：a global search tool to find optimal equivalent source sets[J].Journal of Sound and Vibration，2009，322（1-2）：282-298.

[12] ? SONG L M，KOOPMANN G H，FAHNLINE J B. Numerical errors associated with the method of superposition for computing acoustic field[J].The Journal of the Acoustical Society of America，1991，89（6）：2625-2633.

[13] ? 向宇，石梓玉，陸靜，等.基于波疊加法的非共形近場聲全息波函數的構造與選擇[J].振動與沖擊，2020，39（15）：183-192.

[14] ? OCHMANN M A.The source simulation technique for acoustic radiation problems[J].Acta Acustica united with Acustica，1995，81（6）：512-527.

[15] ? OCHMANN M A.The full-field equations for acoustic? radiation and scattering[J].Journal of the Acoustical Society of America，1999，105（5）：2574-2584.

[16] ? 石梓玉，向宇，陸靜，等.一種提高聲場重構穩定性的射線等效源法[J].廣西科技大學學報，2019，30（3）：1-7，21.

[17] ? 張陽，向宇，石梓玉.二維聲全息的δ函數約束型射線波疊加法[J].廣西科技大學學報，2020，31（4）：19-26，48.

[18] ? CHIEN C C，RAJIYAH H，ATLURI S N.An effective method for solving the hyper-singular integral equations in 3-D acoustics[J].The Journal of the Acoustical Society of? America，1998，88（2）：918-937.

[19]? KIRKUP S M，HENWOOD D J.Methods for speeding up the boundary element solution of acoustic radiation problems[J].Journal of Vibration and Acoustics，1992，114（3）：374-380.

Information compensation of ray wave function based on wave

superposition near field acoustic holography

CHEN Yanhao1，2，XIANG Yu*2，SHI Ziyu2

（1.School of Mechanical and Traffic Engineering， Guangxi University of Science and Technology，

Liuzhou 545006， China; 2.Guangxi Key Laboratory of Automobile Component and Vehicle Technology （Guangxi University of Science and Technology）， Liuzhou 545006， China）

Abstract： Near field acoustic holography based on wave superposition has been widely used in sound source identification， location and field analysis in recent years because of its advantages in adaptability and numerical calculation. In the near-field acoustic holography technology based on the wave? ? ? ? ? ? superposition method， the ray wave function as the integral kernel function of the wave superposition method can effectively improve the ill-conditionedness of the system matrix and improve the numerical stability of the sound field calculation. However， due to the strong directivity of the ray wave function at low wave numbers， the sound field information is mainly concentrated at the main direction of the ray wave function， resulting in the lack of sound field information in other directions. To solve the above problems， a ray wave function information compensation method is proposed， that is， a certain proportion of monopole spherical wave function is added to the ray wave function to compensate the? reconstructed sound field. Numerical example results show that the ray wave function after? ? ? ? ? ? ? ?supplementary information not only guarantees the stability of the numerical calculation in the sound field reconstruction process， but also makes up for the defect of a single ray wave function at low? ? ?wavenumbers， and further improves the accuracy of sound field reconstruction.

Key words： near-field acoustic holography; wave superposition method; ray wave function; information compensation; sound field reconstruction.

（責任編輯：黎? 婭）