變聲速彈性沉積層下壓縮波與剪切波的耦合影響?

劉亞琴3) 楊士莪1)2)3) 張海剛1)2)3)? 王笑寒3)

1)(哈爾濱工程大學,水聲技術重點實驗室,哈爾濱 150001)

2)(哈爾濱工程大學,海洋信息獲取與安全工業和信息化部重點實驗室,哈爾濱 150001)

3)(哈爾濱工程大學水聲工程學院,哈爾濱 150001)

(2018年8月27日收到;2018年9月27日收到修改稿)

1 引 言

海洋中的聲傳播是一個極其復雜的過程,淺海中聲傳播相較深海來說更為復雜.因為在淺海波導中,聲波的傳播受到海底界面的影響更嚴重.在實際的海洋環境中,海底沉積層中的聲速(壓縮波聲速和剪切波聲速)一般是隨深度變化的[1].多數的文獻和實驗結果表明,疏松的海洋沉積層中剪切波聲速一般隨深度連續增加,并且聲速正比于zυ(z以海水與沉積層交界面為基準)[2?7].壓縮波與剪切波聲速隨深度變化的特征意味著聲波與海底相互作用的性質在低頻和高頻時是截然不同的.對于大入射角、遠距離的傳播問題,海底的高頻響應可以將海底當作流體來處理[8].Liu等[9,10]在有液態沉積層(密度隨深度為廣義指數變化,壓縮波聲速為常數、n2線性或平方反比形式)的環境下,研究了粗糙海底的平面波散射聲場的空間功率譜密度以及平面波反射系數,為海底聲學的研究提供了典型的環境模型.但是其研究沒有考慮沉積層中的剪切特性.如果能量穿透沉積層透射在沉積層-彈性基底界面上,那么流體模型就不再適用了[8].彈性體中能量的向上折射導致壓縮波與剪切波之間的相互耦合、轉換.由于聲速的變化而導致的壓縮波與剪切波之間的耦合將會對聲場產生怎樣的影響成為本文研究的內容.

Ewing等[11]給出了在各向同性的非均勻介質中的波動方程.在接下來的有關研究中,一些學者(如Karal,Hook,Scholte)[12?14]在研究非均勻彈性介質中的波動方程時,同時考慮了剪切波聲速、壓縮波聲速以及密度的變化,得到的結果極其復雜.為了能夠使問題簡化且得到更加清晰的物理意義,接下來的有關變參數彈性層的研究主要包括兩個方向:一是在某些特殊情況下給出非均勻彈性層中的解析解或近似解析解.Gupta[15]對密度為常數,拉梅常數的變化形式為λ/λ1= μ/μ1=(1+bz)2(這里λ1,μ1和b為常數)的沉積層環境進行了分析,得到了P波與S波的去耦合方程,并且得到此種情況下的反射系數;Hall等[7]提出了半解析的反射模型,考慮了由于剪切模量隨著深度變化(G(z)∝zp)引起的耦合效應,導出了勢函數P和勢函數S的弱耦合方程、利用Born近似求解勢函數P的方程;另一個主要的研究方向是忽略剪切波及壓縮波聲速梯度引起的耦合.Vidmar和Foreman[16]提出了用傳遞矩陣的方法計算平面波反射系數的模型,其中彈性沉積層的密度、剪切波聲速、壓縮波聲速以及衰減可以任意變化.Westwood等[17]基于簡正波理論提出了ORCA模型,可以計算具有多層彈性環境下的聲場,且層中的剪切波速度和壓縮波速度可以是常數也可以是n2線性的,但是在彈性層中,P波和S波的勢函數滿足的波動方程忽略了耦合.

上述研究在處理彈性沉積層聲速變化的環境時,一般選擇忽略聲速變化引起的耦合,雖然Hall考慮了剪切模量隨深度變化引起的耦合并給出了半解析的反射模型,但是其并沒有分析耦合對反射的影響.Fryer[8]用數值解法研究了由于壓縮波、剪切波聲速梯度引起的P-SV耦合的反射率.通過比較總反射率與忽略耦合得到的部分反射率,得到了P-SV耦合的反射率,并研究了這種耦合對反射率的影響.此外,通過平面波響應函數清楚地揭示了耦合的主要結果是剪切向壓縮運動的轉換.但是數值解法的計算時間會隨著沉積層厚度的增加而增加.考慮彈性沉積層聲速引起的耦合效應,本文求解彈性沉積層某種聲速變化形式下聲場勢函數的近似解析解,進而研究此種耦合對低頻聲場的影響以及聲學參數對耦合效應的影響程度.

本文的內容主要包括:第2部分是聲場建模,重點研究變參數彈性沉積層中壓縮波與剪切波本征函數的近似解析解;第3部分仿真驗證本文方法的有效性,研究剪切波、壓縮波的耦合對低頻聲場(水中聲壓場以及沉積層中質點位移場)的影響;第4部分研究彈性沉積層中聲學參數對耦合的影響;第5部分是實驗數據處理;最后一部分是結論.

2 聲場建模

考慮流體層中單頻點聲源激發的聲場問題,在柱坐標系下建立如圖1所示的環境模型.時間因子為exp(?jωt),為了方便,在以下的推導過程中省略時間因子.

圖1 環境模型Fig.1 .Environment model.

模型中,聲源位于深度為H0,密度與聲速分別為ρ0和c0的均勻流體層中,點聲源的坐標為(0,zs);彈性沉積層厚度為(H1?H0),且沉積層中密度、壓縮波、剪切波聲速分別為ρ1,c1(z),b1(z);半無限彈性海底中密度、壓縮波、剪切波聲速分別為ρ2,c2,b2.各介質層中壓縮波與剪切波勢函數分別為φ0;φ1,ψ1;φ2,ψ2.

文獻[18]中已經證明,在柱坐標系下,ψ1和ψ2只有沿θ方向的分量,記θ方向上的分量分別為ψ1和ψ2,這樣,矢量勢函數就可以轉化為標量勢函數,將這些標量勢函數用Fourier-Bessel積分表示.記:

參看文獻[18,19],Z0(z,ξ),Z2(z,ξ)和G2(z,ξ)的表達式為:

在變參數彈性沉積層中,壓縮波與剪切波聲速隨深度變化將會引起二者之間的耦合,具體表征為函數Z1(z,ξ)和G1(z,ξ)形成耦合方程組.此前沒有相關文獻求解此種環境,本文建立、求解耦合方程組過程如下.

在變參數彈性介質中,φ1,ψ1所滿足的波動方程為[20]:

其中

對波動方程(4)式分別取散度和旋度,進行化簡整理(具體過程詳見附錄A),可得

在變參數彈性介質中Z1(z,ξ)和G1(z,ξ)所滿足的方程為:

其中

觀察(7)和(8)式組成的耦合方程組,考慮特殊情況σ=0時,即橫波聲速不隨深度變化,此時方程組(7)和(8)式為:

此時,無論縱波聲速如何變化,橫縱波之間無耦合.但在一般情況下,б=0,橫縱波之間發生耦合,此前一般都是利用數值處理方法來處理(7)和(8)式,隨著沉積層厚度的增加,數值處理方法計算時間增加.本文提出一種近似解析方法求解耦合方程(7)和(8).考慮到實際情況中α,σ均很小、其二階小量要遠小于α,σ的值,在求解過程中忽略有關二階小量.

和(8)記為:

對方程(11)和(12)分別進行逐次微分,可得

至此,耦合方程組轉化成了非耦合方程(13)和(14).根據方程(13)和(14)可求得函數P(z,ξ)和F(z,ξ),進而可求得本征函數Z1(z,ξ)和G1(z,ξ)的表達式.

方程(13)和(14)均是齊次微分方程,可通過變數代換的方法將其轉化為可解情況[21].具體求解過程在附錄B中給出,方程的解為:

其中

將P(z,ξ)和F(z,ξ)的表達式代入方程(11)和(12),可得Z1(z,ξ)和G1(z,ξ)的表達式為:

其中

從Z1和G1的表達式可以看出,剪切波聲速的σ值將會影響Z1的值,而壓縮波聲速中α對G1的影響在G1表達式中沒有體現,是由于在求解G1表達式的過程中,忽略了(B1)方程(詳見附錄B)中α的二階小量.這里反而是c10的值影響G1的表達式,但是在考慮是否引起壓縮波與剪切波之間的耦合,究其原因是由于聲速的變化,所以c10的值的改變不會對耦合產生影響.

根據點源條件、液/固邊界條件和固/固邊界條件,可得諸待定系數所滿足的方程組為

其中

[aij]9×9的具體表達式見附錄C.求解(19)式可得諸系數的值,進而得到勢函數的具體積分表達式.

進一步求勢函數的值可以用快速傅里葉變換(FFT)方法.以水中勢函數φ0(r,z)為例,φ0(r,z)的積分表達式為[22]

利用FFT變換就可得到勢函數的近似表達式:

同樣的方法可以計算出沉積層和半無限彈性海底中勢函數的表達式.

3 壓縮波與剪切波的耦合影響

3.1 耦合對水中聲壓場的影響

為驗證本文方法(考慮壓縮波與剪切波之間的耦合)的有效性,以COMSOL軟件仿真得到的聲壓傳播損失作為參照.此外,為了分析耦合對水中聲壓場的影響,仿真非耦合情況下的聲壓傳播損失.在仿真驗證之前,先簡單介紹非耦合情況與耦合情況聲場計算的不同.

1)未考慮耦合的彈性沉積層中,Z1(z,ξ)和G1(z,ξ)所滿足的方程(7)和(8)式退化為:

2)非耦合彈性沉積層中的法向應力和切向應力的表達式相較于耦合情況分別缺失項.

在仿真過程中,以具體的環境(如表1所列)為例,取聲源頻率f=25 Hz,聲源深度Zs=100 m,接收深度Zr=150 m.

圖2(a)和圖2(b)是仿真環境I下COMSOL軟件和耦合以及非耦合情況下計算出來的聲壓傳播損失.

表1 環境I的聲場參數Table 1 .Sound f i eld parameters of the ocean environment I.

從圖2(a)中可以看出,本文方法計算出來的傳播損失與COMSOL軟件計算出來的傳播損失基本一致,只有在某幾處距離上有些偏差,產生這些偏差的原因是因為本文方法在推導過程中忽略了α,σ的二階小量,且快速場方法在計算的過程中距離是抽樣間距(抽樣間距?r=2π/(M?kr),其中?kr是波數抽樣間距,M是計算的點數).仿真結果表明了本文方法的有效性.非耦合與COMSOL軟件計算的傳播損失曲線在近距離處基本符合,但是隨著距離的增加,兩條曲線有較大的偏差,如圖2(b)所示.圖3為耦合和非耦合算法中積分核函數Z0(z,ξ)的幅度沿實水平波數軸的分布圖.

圖2 聲壓傳播損失對比圖 (a)耦合情況;(b)非耦合情況Fig.2 .Comparison chart of TL curves for sound pressure propagation loss:(a)In coupled case;(b)in uncoupled case.

由于沉積層中壓縮波、剪切波衰減均為0.1 dB/λ,故在實軸上不存在極點,但是相應的極點明顯地表現為極陡的峰(極陡的峰對應的水平波數可以理解為簡正波理論中的本征值).從圖3中的局部放大圖可以看出,兩種方法計算的陡峰對應的水平波數(即本征值)有些微變化,具體值如表2所列.從表2中可以看出,考慮耦合分別導致3處陡峰處0.21%,0.25%,0.14%的本征值變化,本征值的改變導致水中聲場的變化.

圖3 積分核函數Z0(z,ξ)的幅度沿實水平波數的分布圖Fig.3 .Distribution of amplitude of integration kernel Z0(z,ξ)along real horizontal wavenumber axis.

表2 耦合與非耦合下陡峰對應的水平波數Table 2 .Horizontal wave numbers corresponding to steep peaks under coupled and uncoupled conditions.

圖4(a)為耦合與非耦合情況下計算出來的傳播損失與COMSOL軟件仿真的傳播損失之差,圖4(b)和圖4(c)分別為在0—5 km和5—20 km距離范圍內,耦合、非耦合兩種情況下計算的傳播損失與COMSOL軟件仿真的傳播損失之差絕對值的直方圖.其中紅色代表耦合情況、藍色代表非耦合情況.在0—5 km距離范圍內,耦合與非耦合兩種算法的傳播損失偏差均偏小,但是考慮耦合,偏差絕對值主要集中在0—2 dB(95.75%),而非耦合的偏差絕對值主要集中在0—5 dB(94.32%);在5—20 km距離范圍內,考慮耦合時傳播損失偏差絕對值主要集中在0—3 dB(90.1%),最大差值也僅為11.73 dB,而非耦合情況下,偏差絕對值在3—12 dB范圍內的數據所占比例達到49.99%,在12 dB誤差以上的數據也達到了9.68%的比重.對比圖2、圖3以及圖4的仿真結果,表明在此種環境,壓縮波、剪切波之間的耦合影響水中聲場,若忽略耦合,在5—20 km距離范圍內,傳播損失計算偏差絕對值最大可達35.41 dB,聲場計算偏差較大.

在運行環境一致的情況下,本文方法計算時間要遠小于COMSOL軟件計算時間.以環境I為基準,改變沉積層厚度而其他參數不變時,本文方法計算時間與COMSOL軟件計算時間如表3所列.

圖4 (a)傳播損失之差;(b)0—5 km的直方圖;(c)5—20 km的直方圖Fig.4 .(a)Dif f erence of TLs;(b)histogram in the range of 0 to 5 km;(c)histogram in the range of 5 to 20 km.

表3 計算時間對比表Table 3 .Comparison table of calculation time.

3.2 耦合對沉積層中質點位移場的影響

沉積層中質點水平位移、垂直位移的表達式:

其中urc,uzc是壓縮波產生的水平位移與垂直位移;urs,uzs是剪切波產生的水平位移與垂直位移.其中:

(25)—(28)式表明,壓縮波與剪切波產生的水平位移與垂直位移是Z1(z,ξ),G1(z,ξ),(z,ξ)和(z,ξ)的積分函數.

仿真環境如表4所列,聲源頻率f=25 Hz,聲源深度Zs=100 m.圖5給出了距離為3811 m時,質點水平位移與垂直位移隨深度的變化.其中紅色實線代表耦合算法下壓縮波與剪切波共同作用下產生的水平位移,黑色實線是耦合算法下只計算剪切波產生的水平位移,藍色實線是耦合算法下只計算壓縮波產生的水平位移;相對應的虛線代表的是非耦合情況.圖5表明,壓縮波、剪切波以及二者共同作用產生的水平位移、垂直位移在耦合和非耦合情況下均不同,即耦合改變沉積層中質點的位移場.究其原因,是由于耦合與非耦合情況下本征函數Z1(z,ξ)和G1(z,ξ)及其導數(z,ξ)和(z,ξ)的不同.圖6給出了本征值對應下本征函數及其導數隨深度的變化.耦合與非耦合情況下,壓縮波、剪切波以及二者共同作用產生的水平位移、垂直位移由于本征函數及導數的不同而不同.

表4 環境II的聲場參數Table 4 .Parameters of the ocean environment II.

圖5 水平位移與垂直位移隨深度變化Fig.5 .Horizontal and vertical displacement according to depth changes.

圖6 耦合與非耦合情況下Z1(z,ξ),G1(z,ξ)及其導數隨深度的變化Fig.6 .Functions Z1(z,ξ),G1(z,ξ)and derivative functions vary with depth in coupled and uncoupled cases.

圖5表明,耦合和非耦合情況下,剪切波產生的水平位移、垂直位移曲線分別與共同作用產生的水平位移、垂直位移曲線類似,隨深度變化有兩個極值;壓縮波產生的水平位移、垂直位移曲線隨深度單調變化.可以看出,剪切波在水平位移和垂直位移的產生中起主導作用.從urc,uzc,urs,uzs的表達式(25)—(28)看出,壓縮波產生的位移是y1(z),y2(z),(z),(z)的積分函數,剪切波產生的位移是的積分函數.圖7給出不同水平波數kr(即ξ)時,這些被積函數隨深度的變化圖.在kr值相同時,y3(z),y4(z)的值要遠大于的值大于y1(z),y2(z)的值.所以合成的垂直位移以及水平位移是y3(z),y4(z)以及起主導作用,即剪切波起主導作用.從圖7中可見,隨深度的變化主要有兩種,一種是有兩個極值點,另一種是只有一個極值點,故合成的水平位移、垂直位移隨深度的變化也應該是這兩種形式之一.

圖7 位移被積函數隨深度的變化Fig.7 .Changes of integrand of displacement according to depth.

4 沉積層參數對耦合程度的影響

彈性沉積層中壓縮波與剪切波聲速隨深度的變化引起二者之間的耦合,進而對聲場產生影響,所以這里僅考慮參數壓縮波聲速平方倒數的梯度α和剪切波聲速平方的梯度σ對耦合程度的影響.這里的耦合程度表現為耦合與非耦合情況下傳播損失的差值大小.仿真環境如表4所列,改變參數α,σ的取值,其余參數不變.參數α,σ的取值如表5所列.仿真過程中聲源頻率f=25 Hz,聲源深度Zs=100 m,接收深度Zr=30 m.圖8給出α,σ的不同取值時耦合與非耦合情況下聲壓傳播損失曲線.圖9和圖10分別給出了α,σ不同取值時(表5中編號(a)—(c)和編號(d)—(f))耦合與非耦合情況下5階本征波數的結果.圖9表明,當σ值保持不變、α值增大時,除了α=0.1時非耦合的某階本征波數值大于耦合情況,耦合與非耦合算法計算出來的本征波數值相同.在α值的增大過程中,考慮耦合對本征波數值的改變非常小,這也解釋了圖8(a),圖8(b)和圖8(c)中耦合與非耦合算法計算出來的傳播損失曲線基本一致的原因.圖8(a),圖8(b)和圖8(c)以及圖9表明,彈性沉積層中壓縮波聲速中的α值對耦合影響較小.圖10表明,當α值保持不變、σ值增大的過程中,考慮耦合在一定程度上改變了本征波數的值,導致耦合與非耦合算法計算出來的傳播損失曲線差異增大,如圖8(d)—(f)所示.圖8(d),圖8(e)和圖8(f)以及圖10表明剪切波聲速中的σ值對耦合影響較大.圖8、圖9和圖10表明,彈性沉積層中剪切波聲速中σ值對耦合的影響程度大于壓縮波聲速中α值的影響程度.

表5 聲場參數α,σ的取值Table 5 .Value of sound f i eld parameter α,σ.

圖8 α,σ取值不同時耦合與非耦合情況下的傳播損失 (a)α=0.001,σ=0.001;(b)α=0.08,σ=0.001;(c)α=0.1,σ=0.001;(d)α=0.001,σ=0.006;(e)α=0.001,σ=0.008;(f)α=0.001,σ=0.01Fig.8 .Transmission loss of coupling and uncoupling when α,σ changes:(a)α =0.001,σ =0.001;(b)α =0.08,σ=0.001;(c)α=0.1,σ=0.001;(d)α=0.001,σ=0.006;(e)α=0.001,σ=0.008;(f)α=0.001,σ=0.01.

圖9 耦合與非耦合情況時5階本征波數結果 表5中編號(a),(b)和(c)三種情況時本征波數結果Fig.9 .Results of eigenwave number in three cases of Tab.5(a),(b)and(c).

圖10 耦合與非耦合情況時5階本征波數結果 表5中編號(d),(e)和(f)三種情況時本征波數結果Fig.10 .Results of eigenwave number in three cases of Tab.5(d),(e)and(f).

5 實驗數據驗證

2018年1月,在青島某海域進行了聲傳播特性的實驗研究.一艘漁業運輸船作為發射船,兩條垂直陣列作為接收設備.發射船從遠處向接收設備行駛,航行過程中的輻射噪聲作為本次實驗的目標聲源.發射船行駛過程中進行GPS記錄.實驗海區深度約25 m,海底平坦.查閱該區海圖可知沉積層類型為砂-粉砂-黏土.對于黏土、粉砂和砂質沉積物,壓縮波和剪切波聲速不是恒定不變[23,24].故在聲場建模中,假設海洋環境如圖1所示.參數的獲取首先根據實驗海區海底性質以及文獻[23,24]中給出的該類海底參數的可能范圍,然后進行參數反演得到結果.

目標船有多條線譜,選取326 Hz的線譜,對其進行傳播損失的計算.使用的數據來自于1號垂直陣列上海底地震波拾振器(ocean-bottom seismometer,OBS)設備的聲壓傳感器接收的數據.OBS被放置于海底附近處.

圖11(a)表明:實驗結果與耦合理論計算出來的傳播損失曲線趨勢基本一致.圖11(b)表明:未考慮耦合時計算出來的傳播損失與實驗結果在某些距離處有明顯的趨勢差異(用黑色橢圓標記).對比圖11(a)和圖11(b)的結果,在實際海洋環境建模,疏松沉積層中將剪切波聲速假設成隨深度變化并考慮由此引起的耦合能夠更好地符合實際情況.

圖11 聲壓傳播損失對比圖 (a)耦合情況;(b)非耦合情況Fig.11 .Comparison chart of TL curves for sound pressure:(a)In coupled case;(b)in uncoupled case.

6 結 論

針對彈性沉積層典型聲速變化環境,求解出壓縮波與剪切波本征函數的近似解析解,并分析了由于壓縮波、剪切波聲速隨深度變化引起的耦合效應對低頻聲場(水中聲壓場、沉積層質點位移場)的影響,主要結論如下.

1)壓縮波與剪切波之間的耦合導致本征值減小,耦合對水中近場聲場計算影響較小,但是對遠場聲場影響較大,考慮耦合可以提高遠場聲場預報精度.

2)彈性沉積層中剪切波聲速平方的梯度σ值是壓縮波與剪切波耦合的決定性因素,當剪切波聲速不隨深度變化(即σ=0),壓縮波與剪切波之間無耦合.仿真結果表明,σ值對耦合的影響程度隨σ值增大而變大,且σ值對耦合影響程度比壓縮波聲速平方的倒數的梯度α值大.

3)壓縮波與剪切波之間的耦合導致沉積層中本征函數Z1(z,ξ),G1(z,ξ)及其導數的變化,由于沉積層中壓縮波、剪切波以及二者共同作用產生的水平位移、垂直位移為本征函數及導數的Fourier-Bessel積分,所以耦合與非耦合情況下沉積層中質點位移場不同.在仿真環境下,沉積層中壓縮波與剪切波共同作用形成的水平位移與垂直位移隨深度的變化曲線與剪切波單獨產生的水平位移、垂直位移曲線相類似,剪切波在此過程中起主要作用.

附錄A

為了推導方便,記:

對波動方程(4)式分別取散度和旋度,可得:

因為:

附錄B

其中:

忽略α2項,則

則P(z)的表達式為

附錄C