Pekeris波導環境中的水下運動聲場仿真研究

姚舟磊，于曉林，許偉杰

(1.中國科學院聲學研究所東海研究站，上海 201815；2.中國科學院大學，北京 100049)

0 引言

Pekeris 波導是Pekeris 提出的由一個均勻海水層和一個液體半無限空間組成的波導環境[1]。Pe‐keris波導可以看作是淺海環境的一個合理的簡化模型，研究該模型中的聲傳播問題具有重要的理論和實際意義[2-3]。

波數積分方法是聲場建模的主要方法之一，適用于水平分層的波導環境。Pekeris[1]在1948年首次將波數積分方法引入到水聲學中，利用波數積分方法來處理Pekeris波導中的聲傳播問題。駱文于等[4]通過合理的歸一化，提出了一種數值穩定的波數積分方法。

國內學者對水下靜止聲源的聲場已有較多的研究[5-7]，但對運動聲源的聲場研究較少。如今魚雷等水下武器裝備的運動速度能達到80 m?s-1[8]，這對水下目標跟蹤提出了更高的要求。對水下運動聲源的聲場研究能為聲吶測速和測距提供理論基礎，具有現實意義。國外對水下運動聲源的聲場特性研究已有不少成果，但多數是用簡正波方法[9-10]和射線方法[11]進行分析。Schmidt 等[12]給出了水平不變環境中聲源和接收器聯合運動時的波數積分和簡正波表達式，并用修改的SAFARI模型進行仿真。

本文對Pekeris 波導中的運動聲源聲場進行分析，給出運動聲源二維波數積分的一維近似表示，并利用文獻[4]提出的數值穩定波數積分方法對聲場進行求解。首先在自由空間驗證了近似的一維波數積分表示的準確性，再分析Pekeris 波導中運動聲源的聲場特性。

1 理論推導

1.1 Pekeris波導中運動聲源的波數積分解

Pekeris波導環境示意圖如圖1所示，海水的密度和聲速為ρ1,c1,海水深度為D，海底的密度和聲速為ρ2,c2。聲源為單頻點源，以速度vs勻速運動。由于聲源的運動會破壞柱面的對稱性，因此引入笛卡爾坐標系。

圖1 Pekeris波導環境示意圖Fig.1 Schematic diagram of Pekeris waveguide environment

為了計算方便，令x-y平面為海平面且x軸正方向為聲源運動方向，在t=0時刻聲源正處于z軸，位置記為(0,0,zs)。該聲源滿足如下的波動方程[9]：

式中：ψ表示t時刻的聲壓，Ω為聲源的圓頻率，δ函數為狄拉克函數。利用傅里葉變換對：

可以得到運動點源波動方程對應的亥姆霍茲(Helmholtz)方程：

式中：k=ω/c，是頻率ω對應的波數。

使用二維空間傅里葉變換可將Helmholtz 方程在空間-波數域上進行轉換，其中x域的傅里葉變換對為

對式(4)進行空間傅里葉變換可將x域轉換到kr波數域：

利用δ函數的性質可對式(7)進行化簡，可以得到：

再將式(8)從y域變換到ky域可得：

式(9)是三維笛卡爾坐標系下運動點源深度分離的波動方程，其解的表達式為

其中：G(kx,ky,z,ω)是頻率為ω時的深度格林函數。對式(10)進行傅里葉逆變換就能得到聲場的時域解：

將式(11)代入式(10)并化簡，可得：

式(12)為運動點源聲場的二維波數積分解。式(12)中的雙重積分一般只能通過數值方法求解，而在遠場的情況下，為了提高解的精度需要增大水平波數kx、ky的采樣數量，但會大大增加運算時間。

因此在遠場條件下，假設聲源與接收器都位于x-z平面，且聲源與接收器的水平距離遠大于垂直距離，可認為聲源與接收器距離是勻速變化的。利用穩相法[16]可將式(12)的二維波數積分解近似為遠場條件下的一維波數積分：

式中：Ω±kxvs表示聲源運動產生的多普勒頻移對深度格林函數造成的影響，當聲源靠近接收器時多普勒頻移大于0，符號為正；相反當聲源遠離接收器時取負號。式中R(t)=R0?vst為t時刻時聲源與接收器的水平距離，R0為t=0時聲源與接收器的初始水平距離，當聲源接近接收器時取負號，聲源遠離接收器時取正號。

1.2 數值穩定格林函數的表示

在之前的研究中，文獻[13]和文獻[14]給出了兩種Pekeris 波導中點源聲場的深度分離格林函數表達式：

文獻[4]指出，式(17)和式(18)理論上是正確的，但是在求解待定系數矩陣方程時，由于上下行波選取的參考點不合理，導致在實際仿真中系數矩陣容易產生數值溢出的情況，會使求解得到的深度格林函數不穩定。文獻[4]通過合理的歸一化，將深度格林函數表示為

代入邊界條件(14)～(16)，可得線性方程組：

通過求解式(20)中的方程組，可得到Pekeris波導中數值穩定的深度格林函數，在計算過程中再將頻率替換成由聲源運動所產生的多普勒頻移Ω±kxvs并代入式(13)，便可得到運動聲源數值穩定的一維時域波數積分近似解。

2 仿真算例

第1 節給出了運動點源聲場的二維波數積分解。但由于二維解的計算量較大，因此在遠場情況下引入一定的假設，將解近似為一維波數積分解。本節先模擬自由空間下運動點源的聲場，通過對比遠場條件下的一維波數積分解與解析解，來驗證式(13)的準確性；再根據式(13)和式(20)的方程組對Pekeris 波導環境下的水下運動聲源的聲場進行研究。

2.1 運動聲源聲場一維波數積分數值解的驗證

自由空間中單頻運動點源所產生的聲場存在解析解。假設點源頻率為f0，以速度v沿x軸勻速運動，t=0時刻位于原點，則聲場解析解為[15]

其中：

式中：M為聲源運動速度v與介質聲速c的比v/c。

模擬自由空間參數如下：空間的介質聲速為1 500 m?s-1，介質密度為1 000 kg?m-3,運動點源的頻率為100 Hz，從原點(0,0,0)開始沿x軸正方向、以100 m?s-1速度勻速運動，接收器坐標為(10 000,0,50)。

圖2(a)、圖2(b)、圖2(c)分別是在0、20、40 s時計算解析解和數值解得到的聲壓幅值，對應的是聲源與接收器距離10、8和6 km的時刻。圖2(d)是0～80 s接收器接收聲壓級的變化曲線。圖2(d)中解析解與數值解得到的聲壓級變化曲線吻合，誤差為0.3 dB左右。圖2中藍色實線是解析解，黃點是式(13)計算的數值解?？梢钥闯鲈谧杂煽臻g中，聲源與接收器相距較遠的遠場條件下，一維波數積分數值解和解析解的計算結果符合的比較好。

圖2 不同時刻的接收波形和接收聲壓級的解析解與數值解對比Fig.2 Comparison between analytical and numerical solutions of the received waveforms and sound pressure levels at different times

因此，在遠場條件下，一維波數積分的數值解的精度較高，能滿足運動聲源聲場的計算。在自由空間下運動聲源多普勒頻移的理論值為

式中：f0為聲源的頻率，θ(t)為t時刻聲源、接收器連線方向與聲源運動方向的夾角。當水平距離遠大于垂直距離時可忽略角度θ的影響，計算得到距離較遠時的多普勒頻移理論值為107.14 Hz。圖3對比了0～99 s接收信號解析解、近似解的瞬時頻率和上式計算得到的多普勒頻移理論值，可以看到接收信號的瞬時頻率與理論值吻合較好。

圖3 信號瞬時頻率解析解、近似解和理論值的對比Fig.3 Comparison of the analytical and approximate solutions of signal instantaneous frequency and the theoretical values

2.2 Pekeris波導運動聲源數值穩定解的仿真

將Pekeris 波導中數值穩定深度格林函數代入式(13)，并進行數值仿真。Pekeris波導仿真環境如圖4所示，其中海水聲速為1 500 m?s-1，密度為1 000 kg?m-3，海深為100 m，海底聲速為1 800 m?s-1，密度為1 800 kg?m-3，海底的吸收系數為0.2 dB?λ-1。在初始時刻聲源位于z軸，深度為36 m，坐標為(0,0,36)。聲源頻率為50 Hz，以100 m?s-1速度沿x軸正方向運動。接收器深度為46 m，與聲源水平距離10 000 m，坐標為(10 000,0,46)。

圖4 Pekeris波導環境的仿真Fig.4 Simulation of Pekeris waveguide environment

當聲源靜止時，直接求解式(20)中的方程組得到該波導環境下靜止聲源的深度格林函數的幅度，結果如圖5所示。用KRAKENC聲場仿真軟件仿真該波導環境，計算得到該波導環境共會產生4階簡正波。表1 是圖5 的4 個峰值幅度對應的波數kx與KRAKENC求解得到的各階簡正波水平波數值。表1顯示該波導環境下靜止聲源的深度格林函數峰值對應的波數與各階簡正波的水平波數相符合。

表1 各階簡正波水平波數與深度格林函數峰值對應波數對比Table 1 Comparison of the horizontal wavenumber of each order of normal waves and the corresponding peak wavenumber of depth Green's function

圖5 靜止聲源的深度格林函數幅度Fig.5 The amplitudes of the depth Green's function for the stationary sound source

當聲源運動時，深度格林函數需要加上多普勒頻移，即式(13)中的Ω±kxvs，當聲源靠近接收器時取負號，遠離接收器時取正號。圖6是該波導環境下聲源不同運動狀態時深度格林函數的幅值，其中實線是聲源靜止時的深度格林函數的幅值，虛線是聲源以100 m?s-1速度遠離接收器時深度格林函數的幅值，點劃線是以100 m?s-1速度接近接收器時深度格林函數的幅值。由圖6可以看出，當聲源接近接收器時多普勒頻移使深度格林函數峰值對應的波數增大；聲源遠離接收器時深度格林函數的峰值波數減小。

圖6 聲源不同運動狀態時深度格林函數的幅度Fig.6 Amplitudes of the depth Green's function under different motion states of a sound source

求解式(13)得到時間區間為0～80 s時接收聲壓信號的波形，結果如圖7所示。其中藍線是聲源靜止時的接收信號波形，紅線是聲源以100 m?s-1速度勻速靠近接收器時接收的信號。

圖7 在0～80 s時間內的接收波形Fig.7 The received waveform from 0 to 80 s

圖8(a)～8(c)分別對應0、20、40 s時刻接收器的聲壓時域波形。由圖8可以看出，由于聲源與接收器水平距離隨時間發生變化，接收聲壓的幅度會發生起伏，并隨著聲源的接近幅度增大。由于聲源向接收器運動產生多普勒頻移，接收信號的頻率大于聲源靜止時接收信號的頻率。

圖8 不同時刻淺海Pekeris波導接收信號的時域波形Fig.8 Received signal waveforms in the shallow-sea Pekeris waveguide environment at different times

在淺海Pekeris 波導中，聲源產生的聲波會沿不同的路徑到達接收器，除了直達聲波外，還存在經過海面和海底多次反射后的聲波。由于聲源的運動，不同傳播路徑的聲波會產生不同的多普勒頻移，2.1 節的多普勒頻移求解公式僅適用于對直達聲波的求解。

圖9 為計算得到聲源通過接收器正上方的50～150 s 時間區間接收信號的瞬時頻率與直達聲波的多普勒頻移理論值的比較，可以看到當聲源靠近接收器時接收信號頻率大于原始頻率，遠離時小于原始頻率。由于直達聲波具有最大的多普勒頻移，沿其他路徑傳播的聲波多普勒頻移小于直達聲波，因此當聲源接近接收器時仿真得到接收信號的瞬時頻率略小于直達聲波的多普勒頻移。由于在Pekeris波導中運動聲源與接收器的水平距離隨時間發生變化，接收信號的幅度會隨時間發生起伏。在幅度較小處的瞬時頻率會發生較大的變化，因此在圖8中可以看到計算得到的瞬時頻率存在突變。

圖9 一維波數積分解得到的信號瞬時頻率與理論值對比Fig.9 Comparison between the signal instantaneous frequency obtained by the one-dimensional wavenumber integration solution and the theoretical value

截取并放大70～80 s 時刻的時域波形和瞬時頻率，如圖10所示。由圖10中的仿真結果可以看到瞬時頻率突變的時刻與波形幅度極小值相吻合。

對0～80 s的接收信號進行傅里葉變換得到接收信號的頻譜，如圖11所示?？梢钥吹浇邮招盘柕念l譜有多個譜峰，這是由于Pekeris 波導中存在多途效應，聲波沿多條路徑到達接收器。通過KRAKENC模擬得到該波導環境下50 Hz的聲源會產生4 階簡正波，不同的簡正波有不同的相速度，因此也就有不同的多普勒頻移[9]，當聲源的運動軌跡通過接收器的正上方，各階簡正波的頻移不會發生改變，因此接收的信號為多個單頻信號，在圖11中顯示了多個譜峰，譜峰的位置是由波導條件決定的。

圖11 在0～80 s時間內接收信號頻譜Fig.11 The spectrum of the received signal from 0 to 80 s

改變接收器的深度，并計算接收信號的頻譜，得到的結果如圖12所示，圖中橙色虛線和紫色點劃線分別為接收器在66 m、86 m 深度時接收信號的頻譜?？梢钥吹阶V峰的高低發生了變化，但是對應的頻率并沒有改變。將聲源頻率改為35 Hz，接收信號的頻譜如圖13所示。利用Krakenc計算該波導環境的簡正波函數，可以得到3階簡正波，與圖13中的3個譜峰相對應。

圖12 不同深度接收信號的頻譜Fig.12 Spectrums of the received signals at different depths

圖13 聲源頻率為35 Hz時接收信號的頻譜Fig.13 The spectrum of the received signal at the sound source frequency of 35 Hz

目前水下航行器和水下武器的運行速度還難以達到100 m?s-1，因此在相同波導環境下對相對低速運動的聲源進行仿真。聲源聲速為50 m?s-1和20 m?s-1時仿真得到0～80 s 接收器接收信號的時域波形如圖14所示。與圖7 相比，聲源運動速度減小，接收信號幅度的變化速度也變緩。

圖14 聲源不同速度運動時接收器的時域波形Fig.14 The time domain waveforms of the received signals when the sound source moves at different speeds

圖15 為不同運動速度時接收信號的時頻圖，可以看到當聲源運動速度減小時，接收信號的多普勒頻移也會減小。從圖15 中的顏色深淺變化也能看出速度越快接收信號變化得越快。

圖15 聲源不同速度運動時接收信號的時頻圖Fig.15 Time-frequency diagrams of the received signals when the sound source moves at different speeds

3 結論

本文為求解Pekeris 波導中運動聲源的聲場提供一種新思路，將Pekeris 波導中靜止聲源的數值穩定波數積分解應用到運動聲源的聲場分析中。在仿真實驗中，驗證了二維波數積分的一維近似數值解的正確性，并仿真分析了Pekeris 波導下運動的單頻點源產生的聲場。波數積分方法是直接數值求積分的方法，本文只在遠場條件下引入較小的誤差。利用數值穩定的計算方法來求解運動聲源的深度格林函數，不存在數值溢出的問題。因此本文的方法可以作為研究Pekeris 波導中運動單頻點源的遠場特性的參考模型來應用。