基于K-空間積分格林函數(shù)的近場聲全息技術(shù)

邵光輝，牛悅嬌，馬佳男

（1.中船重工第715研究所第十研究室，浙江 杭州 310012；2.河北聯(lián)合大學(xué)，河北 唐山 063009；3.哈爾濱工程大學(xué)，黑龍江 哈爾濱 150001）

基于K-空間積分格林函數(shù)的近場聲全息技術(shù)

邵光輝1，牛悅嬌2，馬佳男3

（1.中船重工第715研究所第十研究室，浙江 杭州 310012；2.河北聯(lián)合大學(xué)，河北 唐山 063009；3.哈爾濱工程大學(xué)，黑龍江 哈爾濱 150001）

目前基于Neumann邊界條件空間聲場變換技術(shù)格林函數(shù)的主要有限離散化算法是k-空間抽樣格林函數(shù)法.然而，該條件下的k-空間積分格林函數(shù)法卻很少應(yīng)用在相關(guān)的課題研究中.所以對于該算法的重構(gòu)特性以及這兩種算法在空間聲場變換中的優(yōu)劣是人們一直普遍關(guān)注的問題.本文通過計(jì)算修正了該算法在Neumann邊界條件下的計(jì)算參數(shù)，并通過數(shù)值仿真分析了k-空間積分格林函數(shù)的相關(guān)特性，最后通過仿真計(jì)算給出了這兩種方法在不同重構(gòu)參數(shù)下的誤差分析.這些分析結(jié)果可為進(jìn)一步的工程實(shí)踐提供參考.

空間聲場變換；平面近場聲全息；格林函數(shù)；法向質(zhì)點(diǎn)振速

近年來，近場聲全息技術(shù)得到了廣泛的應(yīng)用和發(fā)展，其中以基于空間聲場變換（STSF）[1]~[3]、邊界元（BEM）和最小二乘法這三種算法最為主流.目前這三種方法以基于空間Fourier變換和逆變換方法的NAH技術(shù)最為成熟，工程實(shí)現(xiàn)最容易，應(yīng)用最廣泛；BEM法和HELS方法雖然能適應(yīng)各種形狀的聲源，卻都有其自身缺陷，且工程實(shí)現(xiàn)存在不同程度困難.隨著矢量水聽器的發(fā)展，可以通過測量質(zhì)點(diǎn)振速來對聲場進(jìn)行重建和預(yù)測，但現(xiàn)有的重建算法大部分都采用了K-空間抽樣格林函數(shù)法對聲場進(jìn)行重建的，但是由于Neumann邊界條件下k-空間振速-聲壓格林函數(shù)在輻射圓周上存在奇異性[4]，這種奇異性會使得格林函數(shù)的幅值在輻射圓周上具有很大的躍變，從而影響重構(gòu)精度.

本文將給出基于Neumann邊界條件下的另一種格林函數(shù)重建算法--K-空間積分格林函數(shù)法，對聲場進(jìn)行重建，通過仿真計(jì)算驗(yàn)證這種格林函數(shù)離散方法的正確性和可行性，并用這種方法與K-空間抽樣格林函數(shù)法進(jìn)行比較，并給出相應(yīng)的結(jié)論.

1 Neumann邊界條件下的格林函數(shù)

理想流體介質(zhì)中小振幅聲波傳播的波動方程[5]：

對于單頻聲波速度勢Φ=覬e-jωt,ω為聲波的角頻率，覬是空間分布函數(shù)，它滿足Helmholtz方程：

其中，c為聲速，k為聲波波數(shù)，故波數(shù)空間又稱為k-空間.通過Helmholtz方程，可以推導(dǎo)出Helmholtz公式.此公式用聲場邊界函數(shù)值表示聲場（穩(wěn)態(tài)單頻波動聲場）的積分形式解.當(dāng)點(diǎn)源都集中在某一封閉曲面s內(nèi)時，Helmholtz公式表示為

由式（3）可見，Helmholtz公式用覬和鄣覬/鄣n邊界值的面積分來確定聲場中任意一點(diǎn)的速度勢函數(shù)值，因此當(dāng)已知邊界質(zhì)點(diǎn)振速的分布和聲壓的分布值時，就可以用Helmholtz積分求出場中任意點(diǎn)的速度勢函數(shù)值.

格林函數(shù)表示一定邊界條件下點(diǎn)源的場，與邊界條件一一對應(yīng).單頻聲場中的格林函數(shù)滿足下面的方程[6]：

其中，r軆'代表聲源的位置，r軆代表場點(diǎn)的位置，δ函數(shù)表示點(diǎn)源，時間因子取e-jωt，解得：

式（5）表示的是自由場中的格林函數(shù)，稱為格林函數(shù)的基本解.該式為具有1/r奇點(diǎn)形式的函數(shù)且滿足Helmholtz方程，可以成為式（3）的輔助函數(shù).利用格林函數(shù)這種可選性可以選擇適當(dāng)?shù)母窳趾瘮?shù)形式來簡化Helmholtz公式.

Helmholtz方程與第二類邊界條件構(gòu)成的定解問題叫做第二邊值問題或Neumann問題.對于式（2），第二類邊界條件是指鄣覬/鄣n在區(qū)域邊界上為給定函數(shù).相應(yīng)地，該邊界條件下滿足式（5）和Neumann邊界條件的解稱為Neumann格林函數(shù).

根據(jù)“虛源法”，平面邊界下格林函數(shù)的數(shù)學(xué)表達(dá)式為：

式（3）中，封閉曲面s外一點(diǎn)的速度勢可以看作為s上次級元波在場點(diǎn)o的速度勢迭加之總和，相當(dāng)于源點(diǎn)位于s上.當(dāng)s為平面時，有由式（7）可算出Helmholtz公式中的輔助函數(shù)項(xiàng)：

其中uz(x,y,z)，p(x,y,z)分別為空間點(diǎn)(x,y,z)處為法向質(zhì)點(diǎn)振速、聲壓.

通過式（8）可以得到Neumann邊界條件下，實(shí)空間域下法向質(zhì)點(diǎn)振速-聲壓的格林函數(shù):

gN即為Neumann格林函數(shù).由尤拉公式：

對式（3）進(jìn)行二維空間Fourier變換并利用尤拉公式，整理得到k-空間法向質(zhì)點(diǎn)振速-聲壓的格林函數(shù)：

2 K-空間積分格林函數(shù)法

K-空間積分格林函數(shù)法是由W.A.Veronesi等人提出的，哈爾濱工程大學(xué)的金莉萍修正了該格林函數(shù)法在輻射圓周上計(jì)算公式，本文將通過仿真計(jì)算分析該格林函數(shù)對重建結(jié)果的影響，并通過對雙點(diǎn)源的仿真進(jìn)一步驗(yàn)證該方法對重建分辨率和重建精度的影響.

由式（2）可以看出在基于Neumann邊界條件下的格林函數(shù)在輻射圓周上具有奇異性，這種奇異性會使得格林函數(shù)的幅值在輻射圓周上具有很大的跳變，從而影響到重建的精度.K-空間積分格林函數(shù)法是通過格林函數(shù)在K-空間的積分值來改善函數(shù)在輻射圓周上的奇異性.其原理圖如圖1所示.

圖1 K-空間積分原理圖

在波數(shù)域的點(diǎn)(kx,ky)附近的環(huán)形區(qū)域帶k2r2≤k2r≤k2r1上進(jìn)行積分求格林函數(shù)G的平均值，以克服輻射圓周上的奇異性.記kr=(kx2+ky2)1/2，積分環(huán)帶內(nèi)徑 kr1=(kx2+ky2)1/2-△kr，外徑kr2=(kx2+ky2)1/2+△kr，其中，為環(huán)帶寬度的一半.積分分為三個部分，即積分在輻射圓內(nèi)小于kr2的傳播波區(qū)域、倏逝波區(qū)域即大于kr2的區(qū)域，以及傳播波和倏逝波混合的區(qū)域.

于是通過計(jì)算可得，在Neumann邊界條件下的K-空間積分格林函數(shù)：

3 數(shù)值仿真分析

圖2給出的是當(dāng)重建距離增大時，K-空間積分格林函數(shù)的幅值隨kr/k變化曲線.我們知道源面上的聲波向全息面?zhèn)鞑ミ^程中，低波數(shù)的傳播波成分幅值不發(fā)生變化；而高波數(shù)的倏逝波成分幅度將會按指數(shù)規(guī)律衰減.所以通過圖2我們可以看到不同波數(shù)成分對格林函數(shù)的影響.

圖2 距離不同時|GN|與kr/k間的關(guān)系

通過對比可以看出：（1）基于K-空間積分格林函數(shù)法的格林函數(shù)對重建距離更敏感，相較而言K-空間抽樣格林函數(shù)對重建距離并不敏感；（2）當(dāng)kr/k>1時（即高波數(shù)成分），K-空間抽樣格林函數(shù)的幅值衰減迅速在kr/k=1.5時，幅值已衰減殆盡，而基于K-空間積分的格林函數(shù)幅值衰減速度遠(yuǎn)低于K-空間抽樣格林函數(shù)，重建距離越小，該函數(shù)衰減速度越慢；（3）當(dāng)kr/k=1時（即輻射圓周），可以看到兩種格林函數(shù)都有幅值上的躍變，但K-空間積分格林函數(shù)在圓周上的幅值明顯遠(yuǎn)小于K-空間抽樣格林函數(shù)，說明K-空間積分格林函數(shù)法明顯改善了式（2)在輻射圓周上的奇異性；（3）通過對兩種格林函數(shù)衰減速度的對比，我們可以知道K-空間積分格林函數(shù)由于其衰減特性可以獲得更多的高空間波數(shù)成分，尤其當(dāng)重建距離較小時，可能會得到更好的重建精度.

圖3 頻率不同時|GN|與kr/k間的關(guān)系

圖4 不同測量距離下，兩種格林函數(shù)法的聲場重建

由圖3可以看出不同頻率對格林函數(shù)衰減特性的影響，隨著頻率的減小，K-空間積分格林函數(shù)可以得到更多的高空間波數(shù)成分，而由于K-空間抽樣格林函數(shù)倏逝波衰減的十分迅速，使得在低頻條件下，很難捕捉到更多的高空間成分，這樣會對重建精度產(chǎn)生相應(yīng)的影響；隨著頻率的增大，兩種格林函數(shù)在衰減特性上都沒有明顯的變化.

下面通過仿真驗(yàn)證，來進(jìn)一步了解K-空間積分格林函數(shù)法對聲場重建精度的影響.

首先，令全息面為4m×4m，重構(gòu)頻率f=1500hz，兩個點(diǎn)源相距0.2m，采樣點(diǎn)數(shù)N=64，聲速c=1500m/s，dz=0.01m對聲場進(jìn)行振速-聲壓的聲場逆向重構(gòu).

圖4給出了在相同重構(gòu)距離下，不同測量距離對聲場重建精度和分辨率的影響.當(dāng)測量距離Zh=0.02時，兩種格林函數(shù)算法下的聲場重建都獲得了較好的重建分辨率，但是在基于k-空間抽樣格林?jǐn)?shù)算法下重建的聲壓在孔徑邊緣存在較小的起伏，這是由于“卷繞誤差”和其格林函數(shù)本身存在奇異性引起的.將Zh增大到0.06m并保持重構(gòu)距離dz不變，可以看到：隨著測量距離的增大，兩種格林函數(shù)下的重建分辨率都存在減小的現(xiàn)象，但是基于k-空間積分格林函數(shù)算法的重建結(jié)果仍然可以看出輻射聲源的特性；但是隨著Zh的增大，基于k-空間抽樣格林函數(shù)算法下的重建聲壓幅值在孔徑邊緣的起伏增大，即使利用漢寧窗和k域?yàn)V波，都不能抑制卷繞誤差對重構(gòu)結(jié)果的影響，從而影響對聲源特性的判斷.

4 結(jié)論

通過數(shù)值仿真結(jié)果可以得出如下結(jié)論：當(dāng)通過測量法向質(zhì)點(diǎn)振速對聲壓場進(jìn)行重構(gòu)時，相對于k-空間抽樣法，k-積分格林函數(shù)法有效地抑制了重構(gòu)中常存在的“卷繞誤差”并改善了k-空間抽樣格林函數(shù)在輻射圓周的奇異性.當(dāng)增大測量距離保持重構(gòu)距離不變時，基于k-空間積分格林函數(shù)法相對于k-空間抽樣格林函數(shù)算法的聲場重建，具有相對較好的重構(gòu)精度和重建分辨率，可以預(yù)見該算法在復(fù)雜聲源存在的聲場條件下會有較好的應(yīng)用前景，其應(yīng)用特點(diǎn)也可為進(jìn)一步的工程實(shí)踐提供參考.

〔1〕J.D. Maynard, E.G. Williams, Y.Lee, Near field Acoustic Holography: I. Theory of Generalized Holography and the Development of NAH,J.Acoust. Soc.m.1985,78(4):1395-1413.

〔2〕陳曉東.近場平面聲全息的測量和重構(gòu)誤差研究[D].合肥工業(yè)大學(xué)，2004

〔3〕VeronesiW A,Maynard JD.Nearfield acoustic holography （NAH）II. Holographic reconstruction algorithms and computer implementation,J.Acoust.Soc. Am.1987,81(5):1307-1322.

〔4〕金莉萍.基于格林函數(shù)的典型聲場反演技術(shù)[D].哈爾濱工程大學(xué)，2006.

〔5〕何祚鏞.聲學(xué)逆問題一聲全息變換技術(shù)及源特性判別[J].物理學(xué)進(jìn)展，1996（3）.

〔6〕W.A.Veronesi and J.D.Maynard.Digital holographic reconstruction of source with arbitrarily shaped surfaces. J.Acoust.Soc.Am.,1989,85(2):588~598.

O438.1

A

1673-260X（2012）05-0008-04