大粗糙度起伏下一維海面聲散射建模與分析

周 江，王 斌，王文歡，范 軍

(1.上海交通大學(xué)，上海 200240；2.海洋智能裝備與系統(tǒng)教育部重點(diǎn)實(shí)驗(yàn)室，上海 200240)

0 引 言

起伏海面的聲散射特性研究在水聲通信、主動(dòng)聲吶探測(cè)和水聲對(duì)抗等領(lǐng)域有著廣泛的應(yīng)用[1-5]。對(duì)于中遠(yuǎn)距離的聲場(chǎng)傳播問(wèn)題，大入射角甚至掠入射條件下的起伏海面聲散射特性一直是研究的重點(diǎn)之一。然而，對(duì)于水面目標(biāo)或水面附近目標(biāo)如水面艦船、水下航行器等的聲吶探測(cè)和預(yù)警問(wèn)題，考慮到海面與目標(biāo)之間復(fù)雜的相互作用，需要對(duì)起伏海面的聲散射進(jìn)行建模與分析，尤其是小入射角下大粗糙度起伏海面的多次聲散射不容忽視。近年來(lái)國(guó)內(nèi)外對(duì)粗糙海面的聲散射發(fā)展了多種精確和近似的解析方法，主要有：可以近似計(jì)算絕對(duì)軟邊界條件下粗糙海面散射聲場(chǎng)的積分方程法[6-7]，其優(yōu)點(diǎn)是考慮了多次散射和影區(qū)的貢獻(xiàn)，但因其計(jì)算量較大無(wú)法應(yīng)用于工程中；當(dāng)隨機(jī)起伏方差和相關(guān)長(zhǎng)度遠(yuǎn)小于入射波波長(zhǎng)時(shí)，適用于低頻近似條件的小尺度起伏粗糙面的微擾法[8-10]；當(dāng)粗糙面平均曲率半徑遠(yuǎn)大于入射波波長(zhǎng)、隨機(jī)起伏方差和相關(guān)長(zhǎng)度遠(yuǎn)大于入射波波長(zhǎng)時(shí)，適用于高頻近似條件下的大尺度起伏粗糙面的基爾霍夫近似(Kirchhoff Approximation,KA)方法[11-15]等，其優(yōu)點(diǎn)是形式簡(jiǎn)單，但是不能計(jì)算小入射角下大粗糙度起伏海面的多次聲散射；可以近似計(jì)算均方根斜率較小的粗糙面散射聲場(chǎng)的小斜率近似法[16-17]。隨著計(jì)算機(jī)技術(shù)的快速發(fā)展，一系列數(shù)值計(jì)算方法如有限元方法[18]、邊界元方法[19]等也可用于起伏海面聲散射計(jì)算，但是由于網(wǎng)格劃分工作量巨大，僅適用于低頻計(jì)算。針對(duì)KA方法的缺點(diǎn)，研究人員也做出了相應(yīng)的修正，比如Mc-Cammon等[20-21]利用修正的KA方法計(jì)算了絕對(duì)軟邊界條件下正弦表面的散射場(chǎng)，同時(shí)也考慮影區(qū)、多次散射和表面曲率的影響。Welton[22]將KA方法與菲涅爾(Fresnel)相位近似相結(jié)合，研究了幾何影區(qū)、二階散射和發(fā)射(接收)指向函數(shù)對(duì)絕對(duì)軟邊界條件下正弦表面散射聲場(chǎng)的影響。盡管文獻(xiàn)[20-22]已經(jīng)將海面多次散射引入到KA方法中，但其應(yīng)用局限于確定性正弦海面聲散射情況。Richards等[23]將KA方法與射線理論相結(jié)合，根據(jù)本征聲線計(jì)算海面散射聲場(chǎng)，但是該方法不適用于焦散區(qū)和聲影區(qū)，難以計(jì)算大粗糙度起伏海面的多次散射聲場(chǎng)。

本文在傳統(tǒng)KA方法的基礎(chǔ)上考慮了海面的陰影區(qū)和亮區(qū)之間的多次散射，推導(dǎo)了大粗糙度起伏海面多次聲散射的迭代求解公式，從而提出了一種適用于起伏海面高頻散射聲場(chǎng)計(jì)算的改進(jìn)KA方法。該方法拓寬了傳統(tǒng)KA方法的適用范圍，對(duì)海面的形式、均方根斜率和起伏高度以及聲波的入射角度沒(méi)有限制。本文以一維余弦和高斯譜海面的有限元結(jié)果為標(biāo)準(zhǔn)解，驗(yàn)證了本文方法的準(zhǔn)確性和適用性。以一維余弦海面為例，本文分析了不同入射聲波頻率和角度下海面散射強(qiáng)度的分布情況，并解釋了Bragg散射的產(chǎn)生機(jī)理。以一維高斯譜海面為例，本文討論了不同海面均方根高度和相關(guān)長(zhǎng)度情況下高斯譜海面散射強(qiáng)度的變化規(guī)律。

1 改進(jìn)的Kirchhoff近似方法

針對(duì)聲波入射到一維起伏海面上的情形，本文采用了如圖1所示的聲散射模型，海面S起伏高度為z=η(x)，聲波入射角和散射角分別為θi和θs，海面水平長(zhǎng)度L取決于發(fā)射、接收指向性以及距離。利用改進(jìn)的KA方法求解該模型的中高頻散射聲場(chǎng)，如圖2所示該方法將其聲散射分為兩部分，(1)入射聲波照射海面亮區(qū)的直接散射，忽略了海面影區(qū)對(duì)散射聲場(chǎng)的貢獻(xiàn)。(2)海面的多次聲散射，考慮了海面亮區(qū)內(nèi)某點(diǎn)處的第一次聲散射對(duì)海面上其他點(diǎn)(包括亮區(qū)和影區(qū))的作用，從而引起第二次聲散射，同樣也引起了其他高階散射。

圖1 起伏海面聲散射幾何模型Fig.1 Geometric model of acoustic scattering from undulating sea surface

圖2 起伏海面多次聲散射機(jī)理Fig.2 The mechanism of multiple acoustic scattering from the undulating sea surface

根據(jù)惠更斯第二積分公式[24]，可得滿足自由邊界條件的起伏海面散射聲壓(略去時(shí)間因子e-iωt)為

式中：是單位表面法向量，r2是接收點(diǎn)的位置矢徑，rs是起伏海面的位置矢徑。為了簡(jiǎn)化書(shū)寫(xiě)形式，將向量表示為rα-rβ?Rαβ，其模|rα-rβ|?Rαβ，單位向量RαβRαβ?，下角標(biāo)α和β對(duì)應(yīng)不同的位置矢徑。格林函數(shù)及其法向量可表示為

式中：cosθs=·(rs)，H0(1)和H1(1)分別為第一類零階和一階漢克爾函數(shù)。文獻(xiàn)[25]也指出，海面散射相當(dāng)于一系列偶極子聲源的輻射，此時(shí)ptot(rs)是該等效偶極子模型的總聲源密度[25]，于是有：

式中：pq，q=1，2，…表示第q階偶極子聲源密度，式(3)可以通過(guò)迭代方法進(jìn)行求解。

考慮到海面的第一次聲散射是由入射聲波激勵(lì)的，且海面表面滿足自由邊界條件，根據(jù)Kirchhoff近似可得

式中：入射聲壓pinc(rs)=H0(1)(kR1s)，r1是發(fā)射點(diǎn)矢徑；可見(jiàn)性函數(shù)γ1(R1s)定義為

考慮海面的第二次聲散射，此時(shí)海面上再輻射聲源是由第一階偶極子聲源激發(fā)的，于是有：

式中：rs，r's∈S分別為場(chǎng)點(diǎn)和源點(diǎn)，cosθs'=·(rs')；互可見(jiàn)性函數(shù)γ2(Rss')定義為

依次類推，該等效模型的第q(q≥2)階偶極子聲源密度可歸納為

為了對(duì)式(1)中的積分進(jìn)行數(shù)值求解，將總長(zhǎng)度為L(zhǎng)的海面離散成N個(gè)單元，則該起伏海面的散射聲場(chǎng)近似等于這N個(gè)單元各階散射聲場(chǎng)相干疊加之和。以如圖3所示的線單元對(duì)為例，c(m)和c(n)是它們的幾何中心點(diǎn)，其矢徑分別為和，A(m)和A(n)是單元上任意點(diǎn)，其矢徑分別為和。為了簡(jiǎn)化書(shū)寫(xiě)，將向量表示為其模單位向量上角標(biāo)(m)和(n)分別對(duì)應(yīng)第m和n個(gè)線單元(m≠n)。式(1)可以離散表示為

圖3 線單元對(duì)等效偶極子模型Fig.3 Equivalent dipole models of line element pairs

第n個(gè)單元的第q(q≥2)階偶極子聲源密度為

此時(shí)，對(duì)于式(1)的求解就轉(zhuǎn)化成對(duì)N個(gè)線單元的積分計(jì)算，而線單元的積分式可以進(jìn)一步轉(zhuǎn)化成與線單元坐標(biāo)相關(guān)的求和解析式。利用一階泰勒公式進(jìn)行近似，可以得到：

于是有：

該迭代過(guò)程的終止條件為

式中：ε為收斂限定的最小正數(shù)閾值。得到起伏海面的單位長(zhǎng)度散射強(qiáng)度為

在上述的推導(dǎo)過(guò)程中，考慮了海面的影區(qū)遮擋和多次聲散射效應(yīng)，故式(9)可應(yīng)用于中高頻率、大粗糙度起伏海面的單、雙站聲散射計(jì)算。而經(jīng)典的物理聲學(xué)方法通常假設(shè)幾何影區(qū)對(duì)散射聲場(chǎng)的貢獻(xiàn)可以忽略，故只限于式(9)中q=1時(shí)的一階散射聲場(chǎng)計(jì)算，式(9)僅討論了一維情況。將式(9)中二維格林函數(shù)改為三維格林函數(shù)，即可適用于二維海面的計(jì)算，從而更加符合工程實(shí)際應(yīng)用情況。另外，下文中一維海面相關(guān)結(jié)論均可推廣到二維情況。由于本文側(cè)重點(diǎn)在于方法的提出和驗(yàn)證，以及海面一般規(guī)律總結(jié)，所以對(duì)于不同參數(shù)下任意形狀二維海面的相關(guān)散射計(jì)算不再贅述。

2 方法驗(yàn)證

起伏海面可以看作是由無(wú)數(shù)個(gè)不同波長(zhǎng)、波高和相位的諧波海面疊加形成的。基于兩種不同類型的起伏海面模型對(duì)本文所提方法進(jìn)行驗(yàn)證，對(duì)于規(guī)則的余弦波海面，其波面方程為

式中：-L2≤x≤L2，H、Λ和φ分別為波高、波長(zhǎng)和相位。對(duì)于隨機(jī)起伏海面，選擇典型的高斯譜描述粗糙海面的結(jié)構(gòu)和能量分布，其一維表達(dá)式為[11]

式中：δ為海面的均方根高度，決定了粗糙海面的起伏程度；l為相關(guān)長(zhǎng)度，它決定了粗糙海面的橫向變化頻率。再利用蒙特卡洛(Monte Carlo)方法[26]，結(jié)合高斯譜進(jìn)行一維海面幾何建模。圖4(a)給出了一維余弦海面高度曲線η(x)，其中H=0.4 m，Λ=1 m，φ=π，L=10 m。圖4(b)給出了基于高斯譜生成的一維海面高度曲線η(x)，其中δ=0.4 m，l=1 m，L=10 m。

圖4 一維余弦和高斯譜海面模型Fig.4 One-dimensional cosine and Gaussian sea surface models

下面分別利用改進(jìn)的KA方法、KA方法和有限元方法計(jì)算如圖4所示的兩種一維海面模型的散射強(qiáng)度。圖5比較了θi=10°和f=10 kHz時(shí)三種方法計(jì)算的散射強(qiáng)度隨散射角變化結(jié)果，其中圖5(a)為余弦海面模型，圖5(b)為高斯譜海面模型。圖6比較了θi=10°和θs=10°下的三種方法計(jì)算的海面鏡向散射強(qiáng)度隨頻率變化結(jié)果，其中圖6(a)為余弦海面模型，圖6(b)為高斯譜海面模型。

圖5 一維海面模型三種方法散射強(qiáng)度隨散射角變化計(jì)算結(jié)果Fig.5 Variations of scattering intensity with scattering angle calculated by three different methods for the two one-dimensional sea surface models

圖6 一維海面模型三種方法散射強(qiáng)度隨頻率變化計(jì)算結(jié)果Fig.6 Variations of scattering intensity with frequency calculated by three different methods for the two one-dimensional sea surface models

由圖5和圖6可以看出，對(duì)于大傾角的余弦和大粗糙度的隨機(jī)起伏海面，當(dāng)入射角和散射角均較小時(shí)，傳統(tǒng)KA方法的計(jì)算結(jié)果與有限元結(jié)果誤差較大，且誤差隨著頻率的增大呈增大趨勢(shì)，而本文所提出的改進(jìn)的KA方法很好地解決了傳統(tǒng)KA方法存在的問(wèn)題，其計(jì)算結(jié)果與有限元結(jié)果具有很好的一致性，這為改進(jìn)的KA方法計(jì)算起伏海面聲散射的準(zhǔn)確性和適用性提供了理論依據(jù)。

3 計(jì)算結(jié)果與分析

3.1 不同入射波發(fā)射參數(shù)對(duì)余弦海面散射強(qiáng)度的影響

為了清晰描述不同參數(shù)對(duì)海面聲散射特性的影響，在以下的計(jì)算中均利用改進(jìn)的KA方法計(jì)算式(14)定義的一維余弦海面的分置散射強(qiáng)度。若無(wú)特殊說(shuō)明，海面參數(shù)選取H=0.4 mm，Λ=1 m，L=10 m。

圖7給出了f=10 kHz下四個(gè)典型聲波入射角度下海面的總體、一階和二階的散射強(qiáng)度對(duì)比結(jié)果，圖7(a)～7(e)分別為θi=0°，30°，45°和60°下的對(duì)比，圖7(f)是上述4個(gè)典型入射角度下總體、一階和二階的散射強(qiáng)度的概率分布結(jié)果。

圖7 不同入射角下余弦海面的一階、二階和總散射強(qiáng)度隨散射角的變化及其概率分布Fig.7 Variations of the first order, second order and total scattering intensities with scattering angle and their probability distributions under cosine sea surface and different incidence angles

從圖7(a)和7(c)可以看出，對(duì)該最大波傾角為72°的余弦海面，當(dāng)入射角為0°或30°時(shí)，在部分散射角度下海面的二階散射強(qiáng)度比一階散射強(qiáng)度高出5～15 dB(在θi=0°時(shí)補(bǔ)充了放大圖)，在散射強(qiáng)度[0，15]dB區(qū)間內(nèi)海面的總散射強(qiáng)度概率要大于一階散射強(qiáng)度概率，這說(shuō)明海面的二次散射對(duì)總散射聲場(chǎng)的貢獻(xiàn)不可忽略；當(dāng)入射角為45°和60°時(shí)，只有少數(shù)散射角度下海面的二階散射起作用；而當(dāng)入射角為60°時(shí)，在散射強(qiáng)度[0，15]dB區(qū)間內(nèi)海面的二階散射強(qiáng)度概率幾乎為0，這說(shuō)明海面的高階散射貢獻(xiàn)可以忽略。進(jìn)一步分析圖7可以看出，當(dāng)入射角小于45°，且散射角|θs|小于45°時(shí)，余弦海面的二階散射聲場(chǎng)對(duì)總散射聲場(chǎng)的貢獻(xiàn)較大，而隨著入射角和散射角|θs|的減小，二階散射的影響是逐漸增大的。這說(shuō)明在小入射角和小散射角條件下，海面的多次散射效應(yīng)予以考慮，而且越接近垂直入射，海面的多次散射效應(yīng)越明顯。

為了更清晰地觀察聲波入射角度對(duì)散射強(qiáng)度的影響，圖8進(jìn)一步給出了f=10 kHz下海面散射強(qiáng)度的入射角度-散射角度譜圖。圖8(a)～8(c)分別對(duì)應(yīng)總體、一階和二階的散射強(qiáng)度。圖9給出了θi=0°時(shí)余弦海面散射強(qiáng)度的頻率-散射角度分布圖，圖9(a)～9(c)分別對(duì)應(yīng)總體、一階和二階的散射強(qiáng)度。

圖8 余弦海面散射強(qiáng)度隨入射角度和散射角度變化云圖Fig.8 Variation nephograms of the first order, second order and total scattering intensities with incidence angle and scattering angle under cosine sea surface

圖9 余弦海面散射強(qiáng)度隨頻率和散射角度變化云圖(θi=10°)Fig.9 Variation nephograms of the first order, second order and total scattering intensities with frequency and scattering angle under cosine sea surface and θi=10°

圖8(b)和圖9(b)中的亮線可以用Bragg散射來(lái)解釋。Bragg散射的根本原因是海面波浪間散射回波的相干疊加。Bragg散射示意圖如圖10所示。α和β是相鄰的兩個(gè)波谷，入射聲波通過(guò)點(diǎn)α的路徑長(zhǎng)度為A+D，通過(guò)點(diǎn)β的路徑長(zhǎng)度為D+E。只有通過(guò)兩點(diǎn)α和β的聲程差是聲波波長(zhǎng)λ的整數(shù)倍，即滿足條件(D-E)-(B-A)=±nλ時(shí)，Bragg散射才會(huì)發(fā)生。根據(jù)幾何關(guān)系：(D-E)=Λsinθs和(B-A)=Λsinθi，Bragg散射的條件[27-28]為

圖10 Bragg散射示意圖Fig.10 Schematic diagram of Bragg scattering

圖8、9中黑色虛線均是根據(jù)式(16)畫(huà)出的Bragg散射線。由圖8(b)和圖9(b)可以看出，黑色線的位置與亮線吻合得非常好，且亮線的數(shù)量隨著頻率的提高而增多。這說(shuō)明在中高頻段，當(dāng)聲波波長(zhǎng)遠(yuǎn)小于海面波長(zhǎng)時(shí)，Bragg散射對(duì)海面一階散射聲場(chǎng)的影響很大，而且隨著頻率的提高，Bragg散射效應(yīng)更加明顯。在圖8(c)和圖9(c)中，另一個(gè)發(fā)現(xiàn)是，海面的二階聲散射亮線和一階聲散射亮線出現(xiàn)在相同的位置。

3.2 不同海面參數(shù)對(duì)高斯譜海面散射強(qiáng)度的影響

3.2.1 不同均方根高度下的高斯譜海面散射強(qiáng)度

圖11給出了不同波高下海面的散射強(qiáng)度隨散射角變化結(jié)果，其中θi=10°和f=10 kHz。保持高斯譜海面相關(guān)長(zhǎng)度不變，改變均方根高度δ=0.4和1.2。隨著高斯譜海面均方根高度的增大，海面粗糙度增大，海面二階散射對(duì)總散射聲場(chǎng)的貢獻(xiàn)越來(lái)越大。當(dāng)聲波以10°入射到均方根高度δ=1.2的高斯譜海面時(shí)，海面的多次散射無(wú)法忽略。在小入射角條件下，海面的均方根高度越大，海面的多次散射效應(yīng)越明顯。

圖11 不同海面均方根高度對(duì)高斯譜海面散射強(qiáng)度的影響Fig.11 Influence of different root mean square heights of sea surface on the scattering intensity from Gaussian sea surface

3.2.2 不同相關(guān)長(zhǎng)度下的高斯譜海面散射強(qiáng)度

圖12給出了不同波高下海面的散射強(qiáng)度隨散射角變化結(jié)果，其中θi=10°和f=10 kHz。保持海面均方根高度不變，改變相關(guān)長(zhǎng)度l分別為0.5和1.5。隨著高斯譜海面相關(guān)長(zhǎng)度的增大，海面二階散射對(duì)總散射聲場(chǎng)的貢獻(xiàn)越來(lái)越小，海面的一階散射對(duì)總散射聲場(chǎng)起主要貢獻(xiàn)。當(dāng)聲波以10°入射到相關(guān)長(zhǎng)度l=1.5的高斯譜海面時(shí)，海面的多次散射可以忽略。在小入射角條件下，海面的相關(guān)長(zhǎng)度越小，海面的多次散射效應(yīng)越明顯。

圖12 不同海面相關(guān)長(zhǎng)度對(duì)高斯譜海面散射強(qiáng)度的影響Fig.12 Influence of different correlation lengths of sea surface on the scattering intensity from Gaussian sea surface

4 結(jié) 論

本文針對(duì)起伏海面高頻聲散射計(jì)算問(wèn)題，提出了一種改進(jìn)的物理聲學(xué)方法。與經(jīng)典KA方法的不同之處在于，這種方法考慮了海面的陰影區(qū)和亮區(qū)之間的多次聲散射，較大地提高了計(jì)算精度。本文以一維余弦和一維高斯譜海面為例，該方法與有限元方法計(jì)算結(jié)果的比較，表明該方法適應(yīng)于中高頻、小入射角下大傾角起伏海面的聲散射計(jì)算。根據(jù)一維余弦海面散射強(qiáng)度的計(jì)算結(jié)果，分析了不同入射聲波角度和頻率對(duì)海面散射強(qiáng)度的影響規(guī)律，解釋了Bragg散射的產(chǎn)生機(jī)理，同時(shí)討論了不同海面均方根高度和相關(guān)長(zhǎng)度情況下高斯譜海面散射強(qiáng)度的變化規(guī)律。結(jié)果表明：

(1) 在小入射角和小散射角條件下，均方根高度較大的高斯起伏海面的多次散射效應(yīng)要予以考慮。隨著入射角的減小和海面均方根高度的增大，海面的多次散射對(duì)總散射強(qiáng)度的貢獻(xiàn)越來(lái)越大；

(2) 在小入射角和小散射角條件下，相關(guān)長(zhǎng)度較小的高斯起伏海面的多次散射效應(yīng)要予以考慮。隨著入射角的減小和海面相關(guān)長(zhǎng)度的減小，海面的多次散射對(duì)總散射強(qiáng)度的貢獻(xiàn)越來(lái)越大；

(3) 當(dāng)聲波波長(zhǎng)遠(yuǎn)小于余弦海面波長(zhǎng)時(shí)，Bragg散射對(duì)海面聲散射的影響較大。隨著余弦海面波長(zhǎng)的增大和聲波頻率的提高，Bragg散射條紋增多，海面的Bragg散射效應(yīng)更加明顯。