王 帥，王 斌，范 軍

（上海交通大學(xué) 海洋工程國家重點實驗室，上海 200240）

0 引 言

作為主動聲納的重要判據(jù)，海洋環(huán)境中復(fù)雜目標聲散射特性一直是海洋計算聲學(xué)的研究熱點與挑戰(zhàn)，尤其是在淺海環(huán)境中，目標聲散射特性會變得異常復(fù)雜。針對海洋環(huán)境中復(fù)雜目標聲散射建模問題，國內(nèi)外開展了深入而廣泛的研究，涉及到了幾乎所有的海洋計算聲學(xué)方法，例如：有限元、邊界元、物理聲學(xué)、拋物方程、簡正波方法、波數(shù)積分以及射線聲學(xué)等。

Ingenito[1]采用積分方程方法研究理想波導(dǎo)中目標的聲散射，為單體散射近似奠定了理論基礎(chǔ)；Sarkissian[2]在此基礎(chǔ)上用替代法改進了表面積分方程；Jensen等[3]結(jié)合Ingenito的方法忽略了界面的多次散射，通過用物體的平面波散射函數(shù)來表示該物體，需要計算格林函數(shù)，而任何涉及其中散射積分的計算都是相當復(fù)雜的，即使是簡單的球體或圓柱；Schmidt[4]和Abawi[5-6]提出了一種疊加方法，將散射場近似等于目標內(nèi)表面上大量點源產(chǎn)生的場，但為了獲得點源強度，表面格林函數(shù)的梯度必須精確計算，這在通常情況下是非常困難和低效的；Belibassakis等[7]提出了逐級耦合模態(tài)法，用于計算目標必須為軸對稱的波導(dǎo)中目標的散射場，而克服上述障礙的有限元方法不能直接解決低速大存儲的問題；基于單散射模型和射線追蹤技術(shù)，Chen 和Fan[8]提出了一種幾何聲學(xué)方法來計算淺水目標回波，但是由于Kirchhoff 近似在高頻和聲影區(qū)的限制以及射線理論只適用于高頻，該方法不適用于涉及目標散射聲影區(qū)域的波導(dǎo)中前向散射情況；Fan[9]結(jié)合了邊界元方法和簡正波聲傳播模型建立了聯(lián)合模型，可以計算淺海波導(dǎo)在聲速剖面隨深度變化、與距離無關(guān)的情況下，復(fù)雜形狀目標的散射聲場。另外結(jié)合了簡正波聲傳播模型以及目標散射的Kirchhoff近似方法建立了混合模型，一定程度上解決了大尺寸、復(fù)雜形狀目標在淺海波導(dǎo)中散射聲場的工程預(yù)報問題，僅適用于當目標遠離波導(dǎo)界面情況，即目標和界面兩者之間多次散射可以近似忽略。同時，因為物理聲學(xué)方法不能計算目標聲影區(qū)的散射，所以方法也不適用于淺海波導(dǎo)中的目標前向散射的情況。Jiang[10]提出了基于格林函數(shù)-拋物方程的淺海環(huán)境中目標聲輻射計算方法，不可避免地引入拋物方程誤差。

當目標與邊界的距離相對于聲波波長可以比擬時，不能忽略多次散射效應(yīng)，基于單體散射近似的方法誤差很大。Sarkissian[2]發(fā)現(xiàn)即使在較低頻率，例如目標離界面為10倍目標半徑，多次散射的影響已不能忽略；Makris[12]發(fā)現(xiàn)對比自由空間目標的散射，淺海中目標散射最重要的區(qū)別是，當目標靠近界面時，目標和界面之間的多次散射將起重要作用；Fan[8]論證了淺海信道中，散射聲場中經(jīng)界面反射等多途散射聲的貢獻不可忽略，否則在散射聲壓幅值、散射方向性、頻率響應(yīng)特性等多方面的誤差都較大，這種區(qū)別在目標靠近界面時尤其明顯；Qian[11]也研究了目標靠近海面和海底時，其耦合頻率和自由場差別較大，只有當目標達到一定深度時，才逐漸和自由場結(jié)果基本趨于一致。

目前，考慮了目標與界面之間多次散射的計算方法主要以分離變量或者有限元為主，但局限于理想海洋環(huán)境、簡單目標或者計算效率低下。本文提出了一種基于虛擬柱島、有限元和簡正波混合數(shù)值方法。這種方法利用分離變量或者有限元的優(yōu)勢，即利用有限元方法精確計算目標附近的散射聲壓情況，并根據(jù)簡正波快速計算遠場散射聲場。該方法的優(yōu)勢有：（1）考慮了目標和界面的多次散射，可以計算目標靠近海面或者海底的情況；（2）計算速度相對于有限元方法較快，計算效率更高。同時，該方法適用于低頻散射，尤其是極低頻的條件，可用于復(fù)雜目標在淺海信道中前向散射問題。

1 三維淺海信道聲散射建模方法

1.1 理論模型

對于圖1 所示淺海信道中復(fù)雜目標聲散射問題，本文提出的虛擬柱島、有限元和簡正波混合方法分為三個步驟：（1）利用簡正波方法計算出無目標情況下海洋信道中任意一點的入射聲場；（2）用一個上界面是水面、下界面是海底的虛擬柱島將目標及其所在區(qū)域從淺海信道中分離出來，利用有限元方法計算該區(qū)域內(nèi)散射聲場空間分布，其中背景聲場為步驟（1）計算的入射聲場、虛擬柱島表面為無反射邊界；（3）利用周向諧波、深度方向簡正波正交性，根據(jù)虛擬柱島表面散射聲場確定散射聲場展開系數(shù)以及外部區(qū)域散射聲場。

圖1 淺海信道目標聲散射模型Fig.1 Model of acoustic scattering from objects in shallow waters

分別對上述三個步驟展開敘述，時間因子為e-jωt。

1.1.1 海洋信道入射聲場

假設(shè)水層聲速為c1、密度為ρ1，聲源深度為z0，海底聲速為c2、密度為ρ2、深度為h，在深度H處近似滿足聲壓自由條件。以聲源為中心建立圓柱坐標系，其中z= 0為海面、z軸正方向指向海底。對于單極子點聲源，其輻射聲場具有圓周對稱性，根據(jù)非齊次Helmholtz方程[13-14]、分離變量法以及邊界條件，可以得到入射聲場遠場的表達式為

1.1.2 虛擬柱島內(nèi)散射聲場

目標表面和外部流體相接的耦合面上法向振速連續(xù)、切向應(yīng)力連續(xù)以及切向應(yīng)力為零，由此可寫出目標與流體的耦合方程為

其中，K、C、M分別代表剛度矩陣、阻尼矩陣以及質(zhì)量矩陣，下標w、s和τ分別表示聲學(xué)、力學(xué)以及耦合，u代表結(jié)構(gòu)的位移，ω代表角頻率，ps、pr分別代表彈性和剛性散射聲場，兩者之和為總散射聲場。

在虛擬柱島表面上，采用完美匹配層（PML）技術(shù)[11]模擬無反射邊界，該技術(shù)通過增加吸收系數(shù)將波動方程轉(zhuǎn)換為吸收層控制方程，即

式中，σi為吸收系數(shù)，p為匹配層中的聲壓。

1.1.3 虛擬柱島外散射聲場

相對于單極子源而言，水下目標散射聲場空間分布更為復(fù)雜，但其本質(zhì)仍是聲輻射問題，只是聲源更為復(fù)雜而已。因此，散射聲場也可以表示為滿足遠場輻射、邊界條件和連續(xù)條件的柱坐標系下波動方程形式解[17]，即

根據(jù)步驟（2）計算得到的虛擬柱島表面散射聲壓，利用周向諧波、深度方向簡正波正交性即可確定Anm，即

其中，r0為虛擬柱島半徑。

1.2 算例和驗證

1.2.1 典型算例

計算參數(shù)如下：海水聲速為1500 m/s、密度為1000 kg/m3；海底深度為100 m、聲速為1700 m/s、密度為1500 kg/m3、吸收系數(shù)為0.5 dB/λ[16]、截止深度為300 m；單位幅度單極子聲源頻率為75 Hz，剛性球半徑為10 m，兩者深度均為50 m、水平相距2500 m。圖2給出了無目標情況下聲源所在垂直平面內(nèi)入射聲場空間分布。

圖2 淺海信道中入射聲場空間分布Fig.2 Spatial distribution of incident sound field in shallow waveguide

根據(jù)仿真結(jié)果可以看到，聲波能量大部分集中在海水層中，而海底層中聲波能量隨著水平距離增加迅速衰減。圖3 給出了以剛性球為中心、半徑為20 m 的虛擬柱島內(nèi)剛性球球心與聲源所在垂直平面內(nèi)散射聲場的空間分布。其中，水平向右為x軸正方向，即前向散射，PML厚度為5 m。

圖3 全海深虛擬柱島內(nèi)散射聲壓級Fig.3 Scattering pressure level in virtual cylindrical column in double layers

可以觀察到，散射聲場聲能量也主要集中在海水層，透入海底層的聲波能量會隨著深度的增加而快速衰減、在PML近似無反射。為了更加直觀地觀察目標散射聲場空間分布規(guī)律，圖4給出了海水層內(nèi)虛擬柱島表面散射聲場，可直觀地看出前向散射聲場幅度更大，符合理論規(guī)律。

圖4 海水層虛擬柱島表面散射聲壓Fig.4 Scattering pressure on the surface of virtual cylindrical column in water layer

根據(jù)式（7）計算待定系數(shù)Anm，代入式（6）計算海洋信道中目標的散射聲場，圖5給出了虛擬柱島外剛性球球心與聲源所在垂直平面內(nèi)散射聲場的空間分布，其中周向、深度方向截斷階次分別為10和28。

圖5 淺海信道中散射聲場空間分布Fig.5 Spatial distribution of scattered sound field in shallow waveguide

可以觀察到，散射聲場隨著水平距離增加呈現(xiàn)衰減趨勢；由于高階次簡正波衰減速度快，深度方向散射聲場起伏度逐漸削弱。散射聲波每階簡正波可以分解為一對上行波和下行波，經(jīng)上、下界面多次反射向外傳播。至此，完成了三個步驟的計算，并得到了淺海信道中目標散射聲場。

1.2.2 仿真驗證

由于1.2.1 小節(jié)仿真算例沒有解析解，為此通過對比有限元計算結(jié)果驗證本文方法的正確性，其中本文方法截斷域半徑為15 m、有限元方法截斷域半徑分別為15 m、25 m 和35 m，其它計算參數(shù)和1.2.1小節(jié)相同。圖6對比了不同水平距離上前向散射聲壓級隨深度的變化規(guī)律。

圖6 不同距離上前向散射聲壓級隨深度變化規(guī)律Fig.6 Forward scattering sound pressure levels with respect to depth at different distances

可以看出，兩種方法計算結(jié)果吻合較好，驗證了本文方法計算淺海信道中目標散射聲場的準確性。另外，從計算效率和計算時間上來看，三個截斷域半徑下有限元方法分別需要15 min、70 min 和200 min，而本文方法僅需要15 min，其計算速度、效率優(yōu)勢隨著水平距離的增加而迅速增大。

2 計算參數(shù)選取原則

考慮到界面與目標之間多次散射，本文采用有限元方法計算虛擬柱島內(nèi)散射聲場，為此虛擬柱島半徑與深度直接決定了本文方法的計算速度與效率。

2.1 虛擬柱島半徑

對于有限元方法而言，計算域截斷越小、計算速度越高，剛好覆蓋目標的虛擬柱島最為理想。然而，虛擬柱島表面PML對目標散射場瞬失波吸收效果很差，計算誤差隨著虛擬柱島半徑的減小迅速增大。為此，計算精度和速度是矛盾的，需要結(jié)合實際要求進行平衡和取舍。

圖7 給出了不同虛擬柱島半徑情況下剛性球球心與聲源所在垂直平面散射聲場，其半徑分別為10 m和60 m，即二分之一聲波波長、三倍聲波波長，其它計算參數(shù)與圖3相同。

圖7 不同半徑虛擬柱島情況下目標附近散射聲壓級Fig.7 Scattering pressure levels in virtual column with different radii

對比圖3和圖7可以明顯看出，虛擬柱島半徑為一倍或三倍聲波波長時目標附近散射聲場計算結(jié)果吻合較好，相比之下虛擬柱島半徑為二分之一波長時計算偏差較大。

圖8 給出了根據(jù)圖3、圖7 虛擬柱島表面散射聲壓計算得到的目標前向散射聲場隨深度的變化規(guī)律，其中，水平距離為35 m。

圖8 不同半徑虛擬柱島情況下目標前向散射聲場隨深度變化規(guī)律Fig.8 Forward scattering fields in virtual columns of different radii with respect to depth

可以看到，隨著虛擬柱島半徑增加，淺海信道中目標散射聲場計算結(jié)果逐漸收斂，例如半徑為二分之一聲波波長時計算結(jié)果相對于三倍聲波波長偏差接近15 dB，半徑擴大至一倍聲波波長時計算結(jié)果偏差縮小至2 dB。誤差主要來源于PML，它只能吸收向外傳播的聲波，而目標散射聲場近場及其受上、下界面影響向內(nèi)傳播的聲波則被PML 放大，導(dǎo)致計算誤差增大，因此虛擬柱島半徑不可過小，建議在一倍聲波波長以上。

在仿真的計算時間上，虛擬柱島半徑10 m、20 m 和60 m 情況下計算時間分別為5 min、15 min 和34 h，即虛擬柱島半徑越大，計算時間越長，但計算精度越高。因此，需要結(jié)合實際情況，權(quán)衡計算時間和計算精度。

2.2 海底深度選取

由于聲波在海底沿深度方向傳播過程中能量快速衰減，假設(shè)特定深度方向上聲波近似衰減為0是合理的，那么，海底截斷深度越大，計算精度就越高，而計算速度卻越慢。為此，海底截斷深度是該方法需要重點考慮的另一項重要參數(shù)。

圖9給出了海底截取深度分別為20 m和400 m，即1倍聲波波長和20倍聲波波長情況下剛性球球心與聲源所在垂直平面散射聲場，其他計算參數(shù)與圖3相同。

對比圖9和圖3，虛擬柱島海底深度為1倍聲波波長時目標附近散射聲場計算結(jié)果與20倍聲波波長間的偏差較大，當海底深度增加至10倍聲波波長時與20倍聲波波長之間的偏差已經(jīng)難以察覺。

圖10給出了根據(jù)圖3、圖9虛擬柱島表面散射聲壓計算得到的目標前向散射聲場隨深度的變化規(guī)律，其中，水平距離為15 m。

圖9 不同截止深度虛擬柱島情況下目標附近散射聲壓級Fig.9 Scattering pressure levels in virtual column of different cut-off depths

圖10 不同截止深度虛擬柱島情況下目標前向散射聲場隨深度變化規(guī)律Fig.10 Forward scattering fields in virtual columns of different cut-off depths with respect to depth

可以看到，當海底截止深度取得過小，例如一倍聲波波長時，計算結(jié)果未收斂、偏差較大；當海底截止深度需要足夠大時，例如10 倍聲波波長，即衰減幅度閾值為5 dB，計算結(jié)果已經(jīng)收斂。

在仿真的計算時間上，海底截止深度為20 m、200 m 和400 m 時分別需要2 min、15 min 和20 min，海底截止深度越大，則計算時間越長，計算精度越高。在實際應(yīng)用中，可以根據(jù)需要權(quán)衡計算精度和計算速度，通常情況下海底截止深度衰減幅度閾值為5 dB即可。

3 淺海信道中近海面目標聲散射特征仿真

3.1 剛性球

構(gòu)建的模型如下：海水深度為30 m，海底采用吸收系數(shù)為2 dB/λ 的泥沙；目標是半徑為5 m 的剛性球、深度為10 m；聲源是100 Hz 的單位幅度單極子聲源，位于距離目標水平距離400 m、深度20 m。根據(jù)第2 章分析，虛擬柱島半徑、海底截止深度分別設(shè)置為30 m、37.5 m，即兩倍聲波波長和幅度衰減5 dB所對應(yīng)的深度。

圖11給出了虛擬柱島內(nèi)剛性球球心與聲源所在的垂直平面、目標表面以及虛擬柱島表面散射聲壓級。

圖11 淺海信道中近海面剛性球散射聲場Fig.11 Sound scattered by a rigid sphere near the sea surface in shallow waveguide

結(jié)合示意圖和偽彩圖可見，對于目標而言，目標背向散射能量集中在斜下方，可以判斷出入射聲波是從左下方斜入射的，同時目標前向散射能量也集中在斜下方，這是由于海面反射造成的。相較于上一個案例中水平入射，目標散射聲場能量明顯向海底方向偏轉(zhuǎn)。由于目標靠近界面，受海面多次散射的影響，目標的散射聲壓由于疊加而增強。

3.2 圓柱組合體

除了計算標準體如球之外，本節(jié)以簡化的Benchmark模型，來說明該計算方法對于較為復(fù)雜目標的應(yīng)用。

目標的結(jié)構(gòu)主體是長25 m、半徑2.5 m 的圓柱體，兩端是半徑為2.5 m 的兩個半球，組合而成的圓柱組合體模型，總長度為30 m，深為10 m。目標所在的環(huán)境參數(shù)如下：海水深度為30 m，海底取吸收系數(shù)為2 dB/λ 的泥沙。100 Hz的單位幅度單極子聲源位于水下5 m，目標和聲源水平距離為400 m，目標頭部指向聲源。將上述結(jié)構(gòu)用計算模型構(gòu)建可以得到圖12。

圖12 淺海信道中圓柱組合體模型Fig.12 Cylindrical composite model in shallow waveguide

圖13 淺海信道中近海面組合體散射聲場Fig.13 Sound scattered by a cylindrical composite model near the sea surface in shallow waveguide

通過仿真計算，可以得到虛擬柱島內(nèi)目標中心與聲源所在垂直平面、目標表面以及虛擬柱島表面散射聲壓級，其中虛擬柱島半徑、深度參數(shù)與3.1節(jié)的相同。

可以看出，目標的前向散射聲場能量更強，表面亮點集中在目標主體中后段、且強度是非均勻的。此外，由于目標靠近界面，受界面的多次散射影響，目標上、下方散射空間分布是不對稱的。

4 結(jié) 語

本文提出了一種基于虛擬柱島、有限元法和簡正波的淺海信道中目標低頻聲散射的混合數(shù)值方法，該方法相對于其他方法，具有的特征為：（1）考慮了目標和界面的多次散射，可以計算的目標和界面之間的距離較近，即目標靠近海面或海底的散射情況；（2）有限元和簡正波混合方法比傳統(tǒng)有限元方法計算速度快、效率高。另外，該混合數(shù)值方法適用于低頻散射的情況，尤其是極低頻的條件，且適用于目標在淺海信道中前向散射的情況。為了同時滿足計算精度與速度的需求，虛擬柱島半徑取一倍聲波波長，海底截止深度對應(yīng)幅度衰減5 dB。仿真計算了近海面典型目標在淺海信道中的聲散射特征，對散射聲場方位偏轉(zhuǎn)以及表面亮點分布進行了分析。該方法可應(yīng)用于淺海信道中復(fù)雜目標回聲特征的仿真和評估。