淺海水平變化波導下低頻聲能量傳輸特性

祝捍皓，肖 瑞，朱 軍，唐 駿

(1.浙江海洋大學 海洋科學與技術學院，浙江 舟山 316022；2.浙江大學 浙江省海洋觀測-成像試驗區重點實驗室，浙江 舟山 316021；3.中國科學院聲學研究所 聲場聲信息國家重點實驗室，北京 100190；4.浙江海洋大學 船舶與機電工程學院，浙江 舟山 316022；5.挪威科技大學 聲學組，挪威 特隆赫姆 7491)

淺海中的聲傳播問題是水聲領域研究的熱點之一，是理解、預測及應用各種淺海聲學現象的基礎[1-3].隨著潛艇隱身技術的發展和對各類物理海洋現象的監測需求，各類淺海水下聲學探測手段逐漸轉向低頻(小于1 kHz)，這顯然與低頻信號在淺海中的傳播問題緊密聯系，因而對淺海低頻聲傳播問題的研究獲得了越來越廣泛的關注[4-6].但目前在研究淺海環境下的聲傳播問題時，多數只關注海水層中聲場的計算建模和傳播規律，且通常假設海底為水平分層的液態介質，而有關考慮海底為彈性介質下地形變化復雜的淺海低頻聲場建模的研究則相對較少[7].

淺海海底結構一般可視為由近孔隙介質的沉積層和近彈性介質的基巖層構成，沉積層厚度一般在10 m左右[8-10].由于低頻水聲信號波長長且有較強的穿透力，其在淺海環境中傳播時能穿透沉積層到達基巖層，所以，在研究淺海低頻聲傳播問題時將海底視為彈性介質更為合理.此外，實際淺海海底受內陸河流沖刷及大陸板塊擠壓等因素影響，其海底地形、海底底質均存在較大差異.以海底地形為例，有的海域較平整，有的則成斜坡狀，有的甚至為陡峭的山狀強切割地形.若仍采用傳統水平分層的淺海波導模型將無法準確描述低頻聲信號在淺海中的傳播特性，更無法為淺海水聲工程應用提供合理的理論參考.因此，在對淺海低頻聲傳播問題的研究中，亟待開展在淺海水平變換波導下的討論.

任何聲傳播問題的研究基礎是合理的聲場計算方法，經過近五、六十年的發展，目前在水聲傳播研究領域已建立起多種聲場計算方法，如簡正波方法、射線方法、拋物方程方法、快速場方法以及一系列衍生計算方法[11].但上述方法都是在對波動方程和海洋環境作出一些假設和近似的基礎上推導得到，其通用性均受到了一定的限制，尤其是對淺海水平變化波導下低頻聲傳播問題的計算[12].隨著對淺海水平變換波導下低頻聲傳播問題的日益關注，積極建立通用性更強的淺海聲場計算方法日益重要.有限元法(FEM)[13]把環境劃分為離散的單元進行求解，因而更適用于描述聲場中的各類變化，但由于其計算量大，以往除用于提供基準解外，一般在海洋聲傳播問題中很少使用.近年來，計算機硬件水平的快速提高為FEM在海洋聲傳播問題中的應用提供了有效的支持.

基于上述研究現狀，針對淺海水平波導下低頻聲信號的傳播特性的研究需求，本文采用COMSOL Multiphysics(以下簡稱COMSOL)[14]有限元軟件平臺建立了淺海低頻聲傳播的有限元計算模型.針對當前淺海聲場建模研究中多數只關注海水層中聲場的計算建模和傳播規律且通常假設海底為液態介質的研究現狀，以聲能流為研究對象[15]，同步討論低頻聲波在淺海水體/彈性海底中的傳播變化規律，以期為低頻聲信號在淺海中的探測應用提供理論基礎.

1 理論概述

1.1 聲場的有限元形式

針對本文所關注的淺海環境下的低頻聲傳播問題，在三維柱坐標系下建立了符合淺海環境特點的波導模型.如圖1(a)所示，模型中深度為zs的簡諧點聲源位于柱坐標對稱軸上，海水層與海底層分別近似為均勻各向同性的流體介質和彈性體介質；ρw、cw分別為流體層中密度和聲速；ρb、cp及cs分別為海底層中的密度、縱波聲速和橫波聲速；αp、αs分別為海底層中的縱波聲速衰減和橫波聲速衰減；海水層深度H(r)隨水平距離r變化；Sw為海水層，在模型中表示為均勻流體；Sb為海底層，在模型中表示為均勻彈性介質.由于柱坐標系的軸對稱性，可以將三維問題轉化為二維(r,z)平面上求解[14].

圖1 淺海聲場模型Fig.1 Acoustic field model in shallow water

有限元法求解過程中的定義域如圖1(b)所示.圖中：L1邊界為柱坐標中心對稱軸邊界；L2為海水層與淺海海表交互作用邊界；L3為海水層與彈性海底層的交互作用邊界；模型四周為完美匹配層(Perfectly Matched Layer，PML)[13].

對Helmholtz方程進行權重積分并結合高斯理論，Sw中的聲學有限元方程以及Sb中的有限元振動方程可表示為[13-14]

(Kw+iωCw-ω2Mw)pi=[Fwi]

(1)

(Kb+iωCb-ω2Mb)si=[Fbi]

(2)

式中：Kw、Cw及Mw分別為剛度矩陣、阻尼矩陣及質量矩陣；ω為模型中簡諧點聲源發出聲波的角頻率；pi為海水層節點i處的聲壓；Fwi為聲學激勵；Kb、Cb及Mb分別為彈性海底層中上沒有受到邊界約束部分的剛度矩陣、阻尼矩陣及質量矩陣；si為海底層節點i處的位移；Fbi為結構體上的激勵載荷.

在計算時，設定L1為Dirichlet邊界條件，即邊界L1上聲壓p滿足恒等于0的物理假設，其表達式為

p=0,在L1上

(3)

在海水層與彈性海底層的相互耦合邊界L4上建立聲-彈耦合方程，即滿足法向位移、法向應力連續，切向應力為0，其表達式為

(4)

(5)

(6)

為模擬聲信號在無限大的海洋環境中的傳播，在模型四周增加完美匹配層(Perfectly Matched Layer，PML)[13]，通過在有限元方程中增加吸收系數轉化為PML方程.

通過耦合上述有限元法建立的物理場方程、流-彈耦合方程和PML方程，便可實現對淺海環境下海水層中聲壓p和海底層中位移場s的計算，通過p、s與各層其他物理場間的對應關系，可進一步得到模型中如振速場v、應力場T等其它物理場數值[15].在本文研究中，上述計算過程將在有限元軟件COMSOL平臺上實現.

1.2 波導中聲能流的計算

在傳統標量聲場的建模中，對聲場能量的研究多是基于傳播過程中聲壓的變化特性展開.相比聲壓，聲能流更適合揭示聲波能量傳播中的一般性“流動”特征且更適合用于對水體/海底中聲能量特性的統一討論.因此本文研究中，也將借助聲能流密度這一矢量對低頻聲信號在淺海環境下的傳播特性展開討論.在本文所定義的(r,z)平面內，頻域上海水層與彈性海底層中的聲能流Iw、Ib定義為

(7)

(8)

(9)

(10)

下文以In為研究對象，對淺海水平變化波導下的低頻聲信號傳播特性進行研究.

2 有限元法聲場計算結果驗證

為驗證FEM對海洋聲傳播問題計算結果的正確性，本節首先使用FEM仿真計算了水平[15]及楔形上坡[11]兩類海底地形下的聲壓傳播特性，并與已有聲場計算程序所得結果進行對比.所設定兩類類海洋環境模型與參數示意圖如圖2所示.

兩類模型中，均設定水體中cw=1 500 m/s，ρw=1 000 kg/m3.圖2(a)為一類水平分層海底模型，水深在r=0至r0=4 km (r0為模型最遠距離)處恒為100 m.圖2(b)上坡海底海洋模型設定為美國聲學學會(Acousticcal Society of America，ASA)提出的楔形海域檢驗模型，即在水深在r=0至r0=4 km處由200 m線性減小為0，上坡傾斜角度約為α1=2.68°.兩類模型中均考慮海底為彈性介質，圖2(a)中ρb=1 500 kg/m3，cp=2 000 m/s，cs=1 000 m/s，αpλ=αsλ=0.1 dB (λ為波長)；圖2(b)中ρb=1 500 kg/m3，cp=1 700 m/s，cs=800 m/s，αpλ=λαs=0.5 dB.圖2(a)中zs=20 m，頻率100 Hz.圖2(b)中zs=100 m，頻率25 Hz.

圖3所示為分別采用FEM以及現有聲場計算方法模擬得到的圖2中兩類海洋環境下的單接收點(z=30 m)聲壓傳播損失(Transmission Loss，TL)對比圖[13]，TL計算公式為

(11)

式中：pref│r=1 m為距離聲源1 m處的參考聲壓[1].

圖3中虛線為FEM計算結果，實線對應現有聲場計算方法計算結果.計算中，根據環境模型與現有聲場計算方法適用性，分別選擇快速場程序(Scotter)和拋物方程程序(RAMs)進行求解[16].

可以看出，FEM與兩種聲場計算程序仿真計算得到的傳播損失曲線完全一致，驗證了本文所研究FEM方法在計算水聲場結果的準確性.從仿真結果也可以看到，在聲源頻率及聲源深度一致的前提下，海底地形、海底聲學參數的變化將對波導中聲能量的傳播產生重要影響，而海底地形在水平方向上起伏、變化又廣泛存在于實際淺海環境中.因此，深入研究海底地形變化對低頻聲信號的影響規律，才能準確將低頻聲傳播特性用于水聲工程實踐應用中.

針對上述目標，在本節研究基礎上，下一節將通過數值模擬進一步討論海底地形變化對低頻聲傳播的影響規律，以期為低頻聲信號的應用提供理論基礎.考慮FEM在使用中對波動方程不做任何假設和近似，可以處理各種幾何形狀下的物理場問題，因此更適用于對淺海實際水平復雜變化波導環境下的計算，下文研究中在COMSOL平臺上利用FEM方法開展仿真研究.

圖2 兩類海底地形下海洋環境模型示意圖Fig.2 Schematic diagram for marine environment models of two types of seabeds

圖3 接收深度z=30 m時兩類海洋環境下TL曲線對比圖Fig.3 TL versus r at z=30 m in two marine environments

3 海底地形變化對聲信號傳播特性的影響

3.1 楔形海底

圖4～6分別為當海底聲學參數不變，海底地形分別為水平、上坡及下坡時水下20 m處100 Hz點聲源激發起聲能流在波導中的傳播特征，聲能流強度用聲強級SIL表示，其計算式為

Iref=6.76×10-19(W/m2)

(12)

式中：Iref為參考聲強[1].

以下圖中各點聲能流的大小和傳播方向分別用箭頭的長度和方向來表示.除地形變化外，圖4～6其余仿真條件均與圖2(a)中一致，ΔH1為水平方向4 km處海面與海底的垂直距離.

圖4 水平海底仿真環境下聲能流傳播圖(ΔH1=100 m)Fig.4 Distribution of acoustic energy flux in simulation of horizontal seabed (ΔH1=100 m)

圖5 楔形上坡海底仿真環境下聲能流傳播圖Fig.5 Distribution of acoustic energy flux in simulation of up-sloping wedge-shaped seabed

圖6 楔形下坡海底仿真環境下聲能流傳播圖Fig.6 Distribution of acoustic energy flux in simulation of down-sloping wedge-shaped seabed

在圖2(a)仿真條件下，100 Hz聲波在其水體將激發起共9階波導簡正波，其對應在水平海底界面上的掠射角θ∈[4°,40°][1].對比分析圖4和5的仿真結果可以看到，相比水平海底，由于受α1的影響，水體層中i′次反射后聲波與海底接觸時的掠射角將增加(2i′-1)α1，原本小掠射角對應的低階簡正波傳播特性不斷向大掠射角、易衰減的高階簡正波耦合，聲能量向海底的“泄漏”效應顯著增強且起伏跨度縮短.α1越大，低階簡正波向高階簡正波的耦合越多，高階聲能量向海底“泄漏”效應越強，對應水體層中總聲能量相應的衰減也越快，傳輸至遠場后所保留的簡正波階數也越少，干涉特性越簡單.在圖5仿真條件下，隨著ΔH1從100 m變化至0，水層中100 Hz聲能量傳輸至r=4 km處衰減為20 dB左右.

對比分析圖4和6可以看到，由于有α2的存在，聲源位于淺水處，隨著水平傳播距離的增大水深不斷增加，水體層中i′次反射后聲波與海底接觸時的掠射角將減小(2i′-1)α2，與楔形上坡海底的影響規律相反，原本大掠射角對應的高階簡正波傳播特性不斷向小掠射角的低階簡正波耦合，聲能量向海底的“泄漏”效應明顯減弱且起伏跨度增大.α2越大，高階簡正波向低階簡正波的耦合越多，聲能量向遠處傳播時的衰減越慢，聲場的干涉特性也越簡單，即楔形下坡海底的存在對聲場形成了增強效應.

3.2 小型海底山，小型海溝

圖7 海底山、海溝地形下海洋環境模型示意圖Fig.7 Schematic diagram of marine environment model of seamount and trench

圖8和9分別為當海底聲學參數不變，海底存在小型海底山與小型海溝時100 Hz點聲源激發起聲能流在波導中的傳播特征，仿真中設定海底山與海溝中心位置位于水平方向2 km處.除地形變化外，圖8和9中其余仿真條件均與圖2(a)中一致.

圖8中設定在海底水平距離1.5～2.5 km處存在一小型海底山.

圖8 海底山仿真環境下聲能流傳播圖Fig.8 Distribution of acoustic energy flux in simulation of seamount

對比分析圖4和8的仿真結果可以看到，相比圖4水平海底環境下的聲能量傳播特性，由于小型海底山的“阻擋”，整個仿真傳播過程中聲能量的衰減增加，且海底山越高，聲能量的衰減越大.從圖8的對比中可以看出，當海底存在海底山影響時，海底山前、后聲能量的傳播規律不相同.由于仿真條件下的海底山地形具有底邊短、變化快的起伏特點，低頻聲波受其影響主要集中于海山前后，在點聲源激發起的聲能量到達海底山前受海底山地形的影響并不明顯.隨著聲波不斷臨近海底山，受海底山前沿的反射作用，使部分原沿水平方向傳播的聲能量向海底山上方水層反射，聲能量在海底山上方升高，即對應z<ΔH2處的聲能量受海底山的“阻擋”作用并不顯著.而同時在此過程中，部分聲能量由海底山上坡處透射進入海底，該部分聲能量主要對應掠射角較大，易衰減的中、高階簡正波，因此當聲能量在“翻越”海底山的過程中其包含的中、高階簡正波能量更易快速衰減，在“翻越”海底山后只余少量低階簡正波存在，因此海底山后聲能量起伏穩定且衰減較慢；而當對接收點深度z>ΔH2時，海底山對聲能量傳播產生的明顯的“阻擋”作用，聲能量衰減明顯.

與圖8仿真結果相對應，圖9所示為當海底存在海溝時的聲波傳播特性的仿真圖.仿真中設定在海底r=1.5～2.5 km處存在一海溝.

圖9 海溝仿真環境下聲能流傳播圖Fig.9 Distribution of acoustic energy flux in simulation of trenches on seabed

4 結論

本文基于有限元法，結合具體仿真算例討論了淺海水平變化波導下聲能量的傳播變化特性及機理，以期為淺海水聲工程設備的研發，測試及使用等提供一些理論參考.對比淺海水平分層波導下的聲能量傳播規律，結果表明：

(1)利用有限元法可以準確實現對各類復雜海洋環境下聲場分布的求解.

(2)無論是楔形上坡或是下坡海底，二者均會對淺海中聲能量的傳播產生重要影響.相比水平海底，上坡海底傾角α1越大，高階聲能量向海底“泄漏”效應越強，聲能量向遠處傳播時的衰減越快；而楔形下坡海底地形的影響則正好相反，傾斜角度α2越大，聲能量向遠處傳播時的衰減越慢.正波的轉移越多，聲能量向遠處傳播時的衰減越慢，即楔形下坡海底的存在對聲場形成了增強效應.

(3)當海底存在小型海底山影響時，由于海底山上坡面的反射作用，點聲源激發起的聲能量在“翻越”海底山的過程中，將在海底山上方產生提升，海底山后聲能量起伏穩定且衰減較慢.

(4)與海底山的影響類似，在點聲源激發起的聲能量到達海溝前其受地形的影響并不明顯，當海溝深度較淺，聲能量在“跨越”海溝過程中將產生部分聲能量的額外衰減.隨著海溝深度的增加，當存在夾角大于波導中部分低階簡正波的掠射角時，聲能量在海溝上方水層產生提升，隨著海溝深度的繼續增加，“跨越”后聲能量的起伏特性與水平海底基本一致.

致謝在此感謝哈爾濱工程大學水聲工程學院樸勝春教授課題組為本論文的撰寫提供有限元軟件支持.