胡書(shū)凡， 趙永輝* ， 吳健生2,， 葛雙成

1 同濟(jì)大學(xué)海洋與地球科學(xué)學(xué)院， 上海 200092 2 九三學(xué)社上海市委員會(huì)， 上海 200041 3 浙江水利水電學(xué)院， 杭州 310018

0 引言

瑞雷波是由縱波和橫波垂直分量在自由表面發(fā)生相長(zhǎng)干涉而形成的一種沿自由表面?zhèn)鞑サ牡卣鸩?在垂向非均勻介質(zhì)中，瑞雷波具有幾何頻散特性，即不同頻率分量以不同的相速度傳播，進(jìn)而可用來(lái)探測(cè)地下介質(zhì)的速度結(jié)構(gòu)信息，其應(yīng)用領(lǐng)域包括淺地表地球物理(曾求等，2020；宋政宏等，2020；林融冰等，2020)、勘探地球物理(孟小紅和郭良輝，2007；吳曲波等，2019)及區(qū)域與全球地震學(xué)(顧勤平等，2020；付媛媛和肖卓，2020)等.在淺地表面波勘探中，早期采用的面波譜分析法(Spectral Analysis of Surface Waves, SASW)通過(guò)測(cè)量?jī)傻澜邮招盘?hào)的相位差來(lái)獲取頻散曲線(Nazarian et al.，1986)，但其計(jì)算結(jié)果存在易受隨機(jī)噪聲、體波以及多模式疊加等影響的問(wèn)題(Neducza，2007；Hashemi et al.，2019).而之后發(fā)展的面波多道分析法(Multichannel Analysis of Surface Waves, MASW)將接收到的時(shí)空域多道面波信號(hào)變換到另一個(gè)域中以此來(lái)構(gòu)建頻散圖像(Park et al.，1999；Xia et al.，1999)，根據(jù)面波的強(qiáng)能量特征可從此圖像中提取得到高質(zhì)量的頻散數(shù)據(jù)(Socco et al.，2010)，常用的變換方法有：頻率-波數(shù)變換(Gabriels et al.，1987；Rix et al.，2002)、傾斜疊加變換(McMechan and Yedlin，1981；Xia et al.，2007)以及相移法(Park et al.，1998)等.當(dāng)前，面波多道分析方法以其非侵入性、無(wú)損、高效及經(jīng)濟(jì)的特點(diǎn)而越來(lái)越受到學(xué)術(shù)界的重視(夏江海等，2015).

瑞雷波頻散曲線的正演是其反演的基礎(chǔ)，并可為實(shí)際工作調(diào)查的設(shè)計(jì)提供理論指導(dǎo)依據(jù).當(dāng)前，面波多道分析法中采用的理論頻散曲線計(jì)算方法大多基于水平層狀模型假設(shè)，主要有Thomson-Haskell法(Thomson，1950；Haskell，1953)、Schwab-Knopoff法(Schwab and Knopoff，1970，1972)、δ矩陣法(Watson，1970；Buchen and Ben-Hador；1996)、Abo-Zena法(Abo-Zena，1979)、反射透射系數(shù)法(Kennett and Kerry，1979)以及廣義反射透射系數(shù)法(Chen，1993；何耀鋒等，2006)等.上述方法通過(guò)對(duì)頻散方程進(jìn)行搜根可快速實(shí)現(xiàn)瑞雷波理論相速度的計(jì)算，但實(shí)際淺地表介質(zhì)中往往呈現(xiàn)較明顯的橫向不均勻性，接收排列下方是否能用水平層狀模型來(lái)表達(dá)還有待商榷.而基于波動(dòng)方程的全波場(chǎng)模擬理論上可考慮瑞雷波在任意橫向變化介質(zhì)中的傳播，模擬過(guò)程通常采用有限元(劉雪明等；2009)、有限差分(周竹生等，2007；袁士川等；2018)等數(shù)值計(jì)算方法來(lái)實(shí)現(xiàn)，其中又以交錯(cuò)網(wǎng)格有限差分算法(Virieux，1986)應(yīng)用最為廣泛.在交錯(cuò)網(wǎng)格有限差分算法中，自由表面邊界的處理將對(duì)瑞雷波模擬精度產(chǎn)生較大的影響，早期常采用的真空近似法(Boore，1972)存在數(shù)值穩(wěn)定性差和計(jì)算精度低的問(wèn)題(Graves，1996).之后，通過(guò)許多學(xué)者的研究，分別提出了應(yīng)力鏡像法(Levander，1988)、橫向各向同性介質(zhì)替換法(Mittet，2002)、聲學(xué)-彈性邊界近似法(Xu et al.，2007)和改進(jìn)的真空處理公式(Zeng et al.，2012)等自由表面邊界的處理方法，有效地提高了瑞雷波的模擬精度.在模擬得到全波場(chǎng)記錄后，可與實(shí)際處理一樣采用面波多道分析法來(lái)提取瑞雷波頻散曲線.然而，在實(shí)際正演模擬中往往需要對(duì)多組模型及其參數(shù)變化的影響進(jìn)行分析，由于計(jì)算精度的要求以及穩(wěn)定性條件的限制導(dǎo)致該類(lèi)方法需要較長(zhǎng)的計(jì)算時(shí)間.

鑒于現(xiàn)有頻散曲線計(jì)算方法受限于水平層狀模型假設(shè)或全波場(chǎng)模擬計(jì)算量較大等問(wèn)題，本文從面波多道分析法中聚束分析的計(jì)算公式入手，通過(guò)推導(dǎo)振幅歸一化后聚束分析輸出功率譜中相速度與局部相速度之間的關(guān)系，采用類(lèi)似于從頻散圖像中估計(jì)多道瑞雷波相速度的方式，利用黃金分割極值搜索算法求得理論相速度值，以實(shí)現(xiàn)橫向非均勻介質(zhì)中多道瑞雷波頻散曲線的快速正演.并針對(duì)兩類(lèi)典型的橫向不均勻介質(zhì)模型，分析了局部頻散曲線、平均頻散曲線以及本文方法計(jì)算得到的理論頻散曲線與二維彈性波場(chǎng)模擬分析得到的頻散曲線之間的差異.

1 振幅歸一化聚束分析

(1)

其中，j2=-1；ω為角頻率，單位為rad/s；k為波數(shù)，單位為rad/m；wl為每道的權(quán)重.

上式也可寫(xiě)成矩陣形式，表示為

Z(ω,k)=eH(k)WS(ω),

(2)

在頻散分析過(guò)程中，通常計(jì)算的是聚束分析的輸出功率譜(Rix et al.，2002；Zywicki and Rix，2005)，表示為

P(ω,k)=Z(ω,k)Z(ω,k)H

=eH(k)WS(ω)SH(ω)WHe(k)

=eH(k)WR(ω)WHe(k),

(3)

其中，R(ω)為空間譜相關(guān)矩陣，表示為

(4)

(5)

其中，Δφm n(ω)為第m與第n道之間的相位差.由式(3)可知，利用式(5)進(jìn)行頻散分析得到的相速度與利用相移法估計(jì)得到的相速度仍然是一致的，只是該式計(jì)算的是功率譜，而相移法生成的頻散圖像所表征的是幅度譜.

2 多道瑞雷波頻散曲線正演

(6)

其中，dm n=xn-xm為兩道之間的距離；km n(ω)為兩道之間的空間波數(shù)，確定了瑞雷波在兩接收點(diǎn)之間的傳播；Am n(ω)為與頻率相關(guān)的衰減函數(shù).

在正問(wèn)題中可將地下模型離散為一系列的網(wǎng)格單元，網(wǎng)格單元的厚度可根據(jù)需求設(shè)置相同或不同(如圖1a所示)，若只考慮彈性介質(zhì)，則每個(gè)網(wǎng)格單元的彈性參數(shù)由橫波速度、縱波速度以及密度表征.在此模型離散方式下，(6)式可表示為如下形式

(7)

其中，kq(ω)為局部空間波數(shù)，由m與n道之間第q列網(wǎng)格單元的彈性參數(shù)確定；Δxq為m與n道之間第q列網(wǎng)格單元的寬度.

根據(jù)式(7)，第m與第n道之間的相位差可表示為

(8)

將(8)式代入(5)式，并根據(jù)關(guān)系v=ω/k，可得到振幅歸一化聚束分析輸出功率譜中相速度與局部相速度之間的關(guān)系，表示為

(9)

其中，vq(ω)為局部相速度.

在實(shí)際多道瑞雷波數(shù)據(jù)的頻散分析過(guò)程中，頻散曲線通常是根據(jù)頻散圖像中的譜峰值來(lái)確定的.同樣地，多道瑞雷波數(shù)據(jù)的理論頻散曲線也可通過(guò)求取式(9)的局部極大值來(lái)求得，其具體實(shí)現(xiàn)方法：

(1)利用傳遞矩陣算法求得由每一列模型網(wǎng)格單元構(gòu)成的層狀模型的局部相速度vq(ω).

(2)將局部相速度vq(ω)轉(zhuǎn)換為局部空間波數(shù)kq(ω)，并根據(jù)(8)式計(jì)算每接收對(duì)之間的相位差Δφm n(ω).

(3)利用黃金分割極值搜索算法求得式(9)的局部極大值v(ω)，計(jì)算完所有頻率點(diǎn)后即可得到多道瑞雷波數(shù)據(jù)的理論頻散曲線.

正演算法流程如圖1所示.對(duì)于復(fù)雜地質(zhì)模型式(9)可能存在多個(gè)局部極大值，因此在極值搜索過(guò)程中應(yīng)充分考慮頻散曲線的連續(xù)性，即將上一個(gè)頻點(diǎn)的相速度值作為初值進(jìn)行搜索，從而得到一條與實(shí)際頻散曲線相匹配的理論頻散曲線.

圖1 正演算法流程圖(a) 模型網(wǎng)格剖分示例，黑色表示接收點(diǎn)位置； (b) 局部頻散曲線； (c) 接收對(duì)之間的相位差； (d) 振幅歸一化聚束分析的輸出功率譜，黑色點(diǎn)線為通過(guò)黃金分割極值搜索得到的多道瑞雷波頻散曲線.Fig.1 Workflow for the forward algorithm(a) Example of model discretization. The black indicates the position of receivers; (b) Local dispersion curves; (c) Phase difference between a receiver pair; (d) Output power of the amplitude-normalized beamforming. The black dotted line represents the multichannel Rayleigh wave dispersion curve obtained by the Golden Section search.

3 理論模型試驗(yàn)

為驗(yàn)證本文所提出的基于聚束分析的多道瑞雷波頻散曲線正演算法的可靠性，對(duì)典型的橫向不均勻介質(zhì)模型(斷層模型、含低速異常模型)進(jìn)行了試驗(yàn).首先通過(guò)二維交錯(cuò)網(wǎng)格時(shí)間域有限差分算法(Virieux，1986)，對(duì)各向同性彈性模型進(jìn)行全波場(chǎng)模擬.全波場(chǎng)數(shù)值模擬中采用時(shí)間二階、空間十二階的差分精度，自由邊界條件采用應(yīng)力鏡像法(Levander，1988)進(jìn)行處理；模型大小為60 m(寬)×15 m(深)，橫縱向網(wǎng)格剖分為0.1 m×0.1 m；激發(fā)震源采用中心頻率為20 Hz的雷克子波，時(shí)間采樣間隔取0.25 ms；接收排列共48道，最小偏移距10 m，道間距1 m.在得到了瑞雷波的波形記錄后，采用式(5)對(duì)多道瑞雷波數(shù)據(jù)進(jìn)行頻散分析并根據(jù)譜峰值提取相應(yīng)的頻散曲線.之后，采用基于聚束分析的多道瑞雷波頻散曲線正演算法對(duì)相同模型進(jìn)行計(jì)算得到理論頻散曲線，并分析該理論頻散曲線與波場(chǎng)模擬分析得到的頻散曲線之間的誤差.

由于在一定偏移距范圍內(nèi)瑞雷波以柱面波方式傳播且此范圍內(nèi)體波能量較強(qiáng)，從而對(duì)頻散分析造成影響并導(dǎo)致低頻分量的相速度值被低估或高估，又稱(chēng)之為近場(chǎng)效應(yīng)(Zywicki，1999)，因此本文采用歸一化排列中心距(Normalized Array Center Distance，NACD)和歸一化瑞雷波速度(Normalized Rayleigh Wave Velocity，NRWV)來(lái)進(jìn)行誤差分析(Yoon and Rix，2009)，其計(jì)算公式分別為

(10)

(11)

3.1 斷層模型

如圖2a所示，在該模型中間區(qū)域存在地層錯(cuò)斷，模型中泊松比為常數(shù)0.333，每層密度均為1900 kg·m-3，在橫向大于40 m的區(qū)域內(nèi)地層橫波速度由上到下依次為100、120、170 m·s-1和270 m·s-1.圖中紅色圓圈為激發(fā)震源位置，黑色三角形為接收點(diǎn)位置，黑色虛線代表接收排列中點(diǎn).圖2b為利用全波場(chǎng)模擬得到的共炮點(diǎn)道集，可以看出在偏移距為20 m左右處瑞雷波記錄的視速度發(fā)生了明顯變化，意味著地下存在較強(qiáng)的橫向不均勻性.

圖2 斷層模型及相應(yīng)的單炮記錄(a) 斷層模型，紅色○和黑色分別表示激發(fā)震源和接收點(diǎn)的位置，黑色虛線代表接收排列中點(diǎn)； (b) 模擬的共炮點(diǎn)道集.Fig.2 Fault model and the corresponding single-shot record(a) Fault model. The red ○ and black represent the positions of source and receivers, respectively; (b) The simulated common-shot gather.

圖3a為利用式(9)計(jì)算得到的振幅歸一化聚束分析的理論輸出功率譜, 圖中黑色虛線代表波數(shù)k=2π/L，其中L為接收排列長(zhǎng)度.在模型離散過(guò)程中，每個(gè)網(wǎng)格單元的寬度固定為0.1 m、厚度根據(jù)模型設(shè)定.需要說(shuō)明的是，在理論頻散曲線的正演過(guò)程中并不需要像生成頻散圖像那樣計(jì)算每個(gè)相速度值對(duì)應(yīng)的理論功率譜，而只需進(jìn)行極大值搜索即可.由圖可見(jiàn)，在橫向不均勻模型中理論功率譜將存在多個(gè)局部極大值，而利用考慮頻散曲線連續(xù)性的策略可較好地搜索得到一條光滑變化的理論頻散曲線(圖3a中白色點(diǎn)線).圖3b為利用式(5)對(duì)模擬單炮記錄進(jìn)行波場(chǎng)變換得到的頻散圖像，圖中黑色虛線代表波數(shù)k=2π/L.相較于圖3a，此頻散圖像中的能量分布在一定程度上受到了體波及散射面波的影響.在頻率小于5 Hz的范圍內(nèi)，由于低頻分量的缺失以及泄露誤差的影響，已難以準(zhǔn)確地獲取瑞雷波相速度，同時(shí)在高頻段相速度趨于常速度值，因此在接下來(lái)的對(duì)比中只對(duì)5～50 Hz頻段范圍內(nèi)的分量進(jìn)行分析，從頻散圖像中提取的頻散曲線如圖3b中的黑色點(diǎn)線所示.采用基于聚束分析的多道瑞雷波頻散曲線正演算法計(jì)算得到的理論相速度值如圖3b中的白色點(diǎn)劃線所示，可以看出在整個(gè)頻段范圍內(nèi)該理論頻散曲線與提取的頻散曲線吻合較好，只是在5～6 Hz頻段內(nèi)相速度值要低于提取值，這主要是由于低頻瑞雷波不發(fā)育以及近場(chǎng)效應(yīng)所導(dǎo)致的，理論頻散曲線與提取頻散曲線之間的均方根誤差為4.02 m·s-1，大于6 Hz頻段內(nèi)的均方根誤差為0.94 m·s-1.接收排列中點(diǎn)處的局部頻散曲線如圖3b中的黑色點(diǎn)劃線所示，在所有頻段內(nèi)局部頻散曲線的相速度值都要低于提取得到的相速度值，兩者之間的均方根誤差為6.83 m·s-1，大于6 Hz頻段內(nèi)的均方根誤差為5.54 m·s-1，說(shuō)明當(dāng)?shù)叵麓嬖谳^強(qiáng)的橫向不均勻性時(shí)利用面波多道分析法提取的頻散曲線并不滿(mǎn)足中點(diǎn)假設(shè)，將提取得到的頻散曲線置于排列中點(diǎn)會(huì)給之后的反演解釋帶入較大的誤差.接收排列下方所有局部頻散曲線的平均值如圖3b中的藍(lán)色點(diǎn)劃線所示，在大于6 Hz頻段范圍內(nèi)平均相速度值皆高于提取得到的相速度值，兩者之間的均方根誤差為6.96 m·s-1，大于6 Hz頻段內(nèi)的均方根誤差為6.16 m·s-1，說(shuō)明面波多道分析法得到的頻散曲線并不是排列下方局部頻散曲線的簡(jiǎn)單算術(shù)平均.

圖3 斷層模型頻散分析結(jié)果(a) 基于振幅歸一化聚束分析得到的理論功率譜，白色點(diǎn)線為極值搜索得到的理論多道瑞雷波頻散曲線，黑色虛線代表波數(shù)k=2π/L； (b) 模擬的單炮記錄進(jìn)行波場(chǎng)變換得到的頻散圖像，黑色點(diǎn)線為提取得到的頻散曲線，白色點(diǎn)劃線為正演得到的理論頻散曲線，黑色點(diǎn)劃線為排列中點(diǎn)處的局部頻散曲線，藍(lán)色點(diǎn)劃線為排列下方的平均頻散曲線，黑色虛線代表波數(shù)k=2π/L； (c) NACD-NRWV曲線圖，黑色點(diǎn)劃線、○標(biāo)注線和×標(biāo)注線分別為以理論相速度、局部相速度、以及平均相速度作為參考速度計(jì)算得到.Fig.3 Results of dispersion analysis for fault model(a) Theoretical out power of the amplitude-normalized beamforming. The white dotted line represents the multichannel Rayleigh wave dispersion curve obtained by the local maximum search. The black dashed line indicates wavenumber k=2π/L; (b) Dispersion image generated from the simulated common-shot gather. The black dotted line indicates the extracted dispersion curve, the white dot-solid line represents the theoretical dispersion curve, the black dot-solid line is the local dispersion curve at the midpoint of receiver array, and the blue dot-solid line depicts the average phase velocity beneath the entire receiver array. The black dashed line indicates wavenumber k=2π/L; (c) The plot of normalized Rayleigh wave velocity versus normalized array center distance. The black dot-solid line, the ○ marked line, and the × marked line are calculated by the theoretical phase velocity, local phase velocity, and average phase velocity.

圖3c為將理論頻散曲線、局部頻散曲線和平均頻散曲線作為參考相速度計(jì)算得到的歸一化排列中心距-歸一化瑞雷波速度曲線圖(NACD-NRWV曲線圖).根據(jù)水平層狀模型中近場(chǎng)效應(yīng)的研究(Yoon and Rix，2009；Bodet et al.，2009；Roy and Jakka，2017)，當(dāng)?shù)貙又袡M波速度呈速度遞增型時(shí)，歸一化排列中心距大于1的范圍內(nèi)相對(duì)誤差小于10%；而對(duì)于所有類(lèi)型地層，歸一化排列中心距大于2的范圍內(nèi)相對(duì)誤差小于5%.利用理論頻散曲線作為參考值得到的曲線如圖3c中的黑色點(diǎn)劃線所示，在歸一化排列中心距大于1的范圍內(nèi)最大相對(duì)誤差為7.8%，在歸一化排列中心距大于2的范圍內(nèi)最大相對(duì)誤差為2.1%，這與水平層狀模型中的研究結(jié)果一致，進(jìn)一步說(shuō)明了本文正演方法的可靠性.利用排列中心點(diǎn)處局部頻散曲線作為參考值得到的曲線如圖3c中的圓號(hào)標(biāo)注線所示，在歸一化排列中心距大于1和大于2的范圍內(nèi)最大相對(duì)誤差皆為16.4%，可知在此情況下排列中點(diǎn)假設(shè)并不成立.利用排列下方局部頻散曲線的平均值作為參考相速度得到的曲線如圖3c中的叉號(hào)標(biāo)注線所示，在歸一化排列中心距大于1和大于2的范圍內(nèi)最大相對(duì)誤差皆為11.0%，相比于理論頻散曲線的結(jié)果誤差值較大.

本文所有算法的計(jì)算部分代碼都是在Intel Core i3-4150處理器、Windows操作系統(tǒng)下利用Intel Fortran編譯器進(jìn)行編譯執(zhí)行的.對(duì)于二維彈性波時(shí)間域有限差分算法，由于需要對(duì)人工截?cái)噙吔邕M(jìn)行處理并考慮瑞雷波的影響深度(實(shí)際模型大小為80 m×60 m)，此算例中的計(jì)算耗時(shí)為390 s，若將網(wǎng)格剖分為0.05 m×0.05 m并同時(shí)滿(mǎn)足相同采樣間隔的要求，則計(jì)算耗時(shí)上升至3192 s.而對(duì)于圖2a所示的層狀模型，基于聚束分析的多道瑞雷波頻散曲線正演算法無(wú)需對(duì)模型進(jìn)行深度方向上的網(wǎng)格細(xì)分，在此算例中計(jì)算耗時(shí)僅為0.5 s，若將網(wǎng)格單元的寬度設(shè)置為0.05 m，計(jì)算耗時(shí)也僅為0.7 s，能夠滿(mǎn)足多組模型以及模型中參數(shù)變化影響的快速計(jì)算分析.

3.2 含低速異常模型

如圖4a所示，在該模型20～40 m范圍內(nèi)存在一橢圓形低速異常體，該橢圓的長(zhǎng)軸和短軸分別為10 m和1 m，模型中泊松比為常數(shù)0.333、密度均勻?yàn)?900 kg·m-3，地層中橫波速度由上到下依次為120 m·s-1、170 m·s-1和270 m·s-1，低速異常體的橫波速度為100 m·s-1.圖4b為利用全波場(chǎng)模擬得到的共炮點(diǎn)道集，由于此模型在橫向上的不均勻性相較于上個(gè)模型要更弱，因此在單炮記錄中瑞雷波視速度的變化也更緩一些.

圖4 含低速異常模型及相應(yīng)的單炮記錄(a) 含低速異常模型，紅色○和黑色分別表示激發(fā)震源和接收點(diǎn)的位置，黑色虛線代表接收排列中點(diǎn)； (b) 模擬的共炮點(diǎn)道集.Fig.4 Low-velocity anomaly model and the corresponding single-shot record(a) Low-velocity anomaly model. The red ○ and black represent the positions of source and receivers, respectively; (b) The simulated common-shot gather.

圖5a為利用式(9)計(jì)算得到的理論輸出功率譜以及利用黃金分割極值搜索得到的理論多道瑞雷波頻散曲線(白色點(diǎn)線), 圖中黑色虛線代表波數(shù)k=2π/L.在模型離散過(guò)程中，每個(gè)網(wǎng)格單元的寬度固定為0.1 m、厚度根據(jù)模型設(shè)定.圖5b為利用式(5)對(duì)模擬單炮記錄進(jìn)行波場(chǎng)變換得到的頻散圖像，圖中黑色虛線代表波數(shù)k=2π/L.同樣地，這里只對(duì)5～50 Hz頻段范圍內(nèi)的分量進(jìn)行分析，從頻散圖像中提取的頻散曲線如圖5b中的黑色點(diǎn)線所示.基于振幅歸一化聚束分析的多道瑞雷波頻散曲線正演算法計(jì)算得到的理論相速度值如圖5b中的白色點(diǎn)劃線所示，在整個(gè)頻段范圍內(nèi)理論頻散曲線與提取的頻散曲線一致性較高，并且在此模型中頻散圖像受近場(chǎng)效應(yīng)的影響要更弱，差異主要集中在8.5～10 Hz頻段范圍內(nèi)，這是由于低速異常體引起的散射瑞雷波導(dǎo)致的，理論頻散曲線與提取頻散曲線之間的均方根誤差為1.51 m·s-1.接收排列中點(diǎn)處的局部頻散曲線如圖5b中的黑色點(diǎn)劃線所示，由于模型淺層的橫波速度無(wú)任何橫向變化，局部頻散曲線與提取頻散曲線的差異主要集中在低于18 Hz的頻段范圍內(nèi)，此范圍內(nèi)相速度值要低于從頻散圖像中提取得到的相速度值，兩者之間的均方根誤差為7.12 m·s-1.接收排列下方所有局部頻散曲線的平均值如圖5b中的藍(lán)色點(diǎn)劃線所示，平均頻散曲線與提取頻散曲線的差異也主要是在頻率低于18 Hz的范圍內(nèi)，兩者之間的均方根誤差為3.03 m·s-1.

圖5 含低速異常模型頻散分析結(jié)果(a) 基于振幅歸一化聚束分析得到的理論功率譜，白色點(diǎn)線為極值搜索得到的理論多道瑞雷波頻散曲線，黑色虛線代表波數(shù)k=2π/L； (b) 模擬的單炮記錄進(jìn)行波場(chǎng)變換得到的頻散圖像，黑色點(diǎn)線為提取得到的頻散曲線，白色點(diǎn)劃線為正演得到的理論頻散曲線，黑色點(diǎn)劃線為排列中點(diǎn)處的局部頻散曲線，藍(lán)色點(diǎn)劃線為排列下方的平均頻散曲線，黑色虛線代表波數(shù)k=2π/L； (c) NACD-NRWV曲線圖，黑色點(diǎn)劃線、○標(biāo)注線和×標(biāo)注線分別為以理論相速度、局部相速度、以及平均相速度作為參考速度計(jì)算得到； (d) 利用最小偏移距為1 m、5 m和10 m的單炮記錄得到的疊加頻散圖像，圖中線條表征與(b)一致.Fig.5 Results of dispersion analysis for low-velocity anomaly model(a) Theoretical output power of the amplitude-normalized beamforming. The white dotted line represents the multichannel Rayleigh wave dispersion curve obtained by the local maximum search. The black dashed line indicates wavenumber k=2π/L; (b) Dispersion image generated from the simulated common-shot gather. The black dotted line indicates the extracted dispersion curve, the white dot-solid line represents the theoretical dispersion curve, the black dot-solid line is the local dispersion curve at the midpoint of receiver array, and the blue dot-solid line represents the average phase velocity beneath the entire receiver array. The black dashed line indicates wavenumber k=2π/L; (c) The plot of normalized Rayleigh wave velocity versus normalized array center distance. The black dot-solid line, the ○ marked line, and the × marked line are calculated by the theoretical phase velocity, local phase velocity, and average phase velocity. (d) Stacked dispersion image generated by common-shot gathers with the nearest offset of 1 m, 5 m and 10 m. The representation of plotted lines is the same as (b).

圖5c為將理論頻散曲線、局部頻散曲線和平均頻散曲線作為參考相速度值計(jì)算得到的NACD-NRWV曲線圖.利用理論頻散曲線作為參考值得到的曲線如圖5c中的黑色點(diǎn)劃線所示，在歸一化排列中心距大于1和大于2的范圍內(nèi)最大相對(duì)誤差都為4%，該誤差對(duì)應(yīng)的相速度點(diǎn)位于受散射瑞雷波影響的頻段內(nèi).利用排列中心點(diǎn)處局部頻散曲線作為參考值得到的曲線如圖5c中的圓號(hào)標(biāo)注線所示，在歸一化排列中心距大于1的范圍內(nèi)最大相對(duì)誤差為15.9%，在大于2的范圍內(nèi)最大相對(duì)誤差為14.8%，可知在此情況下排列中點(diǎn)假設(shè)仍然是不成立的.利用排列下方局部頻散曲線的平均值作為參考相速度值得到的曲線如圖5c中的叉號(hào)標(biāo)注線所示，在歸一化排列中心距大于1的范圍內(nèi)最大相對(duì)誤差為6.9%，在大于2的范圍內(nèi)最大相對(duì)誤差為4.6%.

地下存在橫向不均勻性產(chǎn)生的散射瑞雷波將嚴(yán)重影響頻散圖像的質(zhì)量(Yilmaz and Kocaoglu，2012；Mi et al.，2017)，而通過(guò)疊加不同偏移距數(shù)據(jù)的頻散圖像可以較好地壓制非一致性噪聲(Neducza，2007；Socco et al.，2009；Pasquet and Bodet，2017).對(duì)于相同模型，保持接收排列固定不變，將最小偏移距分別設(shè)置為1 m和5 m并進(jìn)行全波場(chǎng)模擬，然后將1 m、5 m和10 m三種偏移距單炮記錄得到的頻散圖像進(jìn)行疊加，疊加結(jié)果如圖5d所示，圖中黑色點(diǎn)線為從疊加頻散圖像上提取的頻散曲線，稱(chēng)之為疊加頻散曲線.由圖可見(jiàn)，疊加頻散圖像受散射瑞雷波的影響較弱，基于振幅歸一化聚束分析的多道瑞雷波頻散曲線正演算法計(jì)算的理論相速度值與提取得到的相速度值之間的一致性有所提升，兩者之間的均方根誤差為0.95 m·s-1，理論頻散曲線與疊加頻散曲線之間的平均相對(duì)誤差為0.4%，最大相對(duì)誤差為2.5%.在頻率低于18 Hz范圍內(nèi)，排列中點(diǎn)局部頻散曲線的相速度值仍低于提取得到的相速度值，兩者之間的均方根誤差為6.60 m·s-1，局部頻散曲線與疊加頻散曲線之間的平均相對(duì)誤差為2.0%，最大相對(duì)誤差為12.2%.平均相速度值與提取得到的相速度值之間的均方根誤差有所上升，為3.38 m·s-1，兩者之間的平均相對(duì)誤差為1.1%，最大相對(duì)誤差為7.6%.可以看出，基于振幅歸一化聚束分析的多道瑞雷波頻散曲線正演算法計(jì)算得到的相速度值與波場(chǎng)模擬分析得到的相速度值之間的誤差較小，具有較高的可靠性.

4 結(jié)論

本文基于振幅歸一化聚束分析，提出了一種適用于橫向不均勻介質(zhì)模型情況下的多道瑞雷波頻散曲線快速正演算法，并通過(guò)與二維彈性波時(shí)間域有限差分模擬結(jié)果對(duì)比驗(yàn)證了該算法的可靠性.理論模型分析揭示，對(duì)于橫向不均勻地下介質(zhì)，利用面波多道分析法提取得到的頻散曲線是接收排列下方介質(zhì)的綜合反映，但不是接收排列下方局部頻散曲線的簡(jiǎn)單算術(shù)平均，更不能認(rèn)為是接收排列中點(diǎn)處的介質(zhì)響應(yīng).本文算法研究可為面波多道分析法的觀測(cè)系統(tǒng)優(yōu)化提供技術(shù)支撐，并且將該方法引入頻散曲線的反演解釋有望克服當(dāng)前面波多道分析法在橫向不均勻介質(zhì)探測(cè)中水平分辨率的限制，提升反演結(jié)果的可靠性.從理論上講，該算法也適用于其他類(lèi)型導(dǎo)波(如勒夫波、Scholte波等)在橫向不均勻介質(zhì)中的多道頻散曲線計(jì)算，而如何考慮多模式面波的傳播，以及相應(yīng)的多道頻散曲線正演算法，將是更具理論意義及實(shí)際應(yīng)用價(jià)值的新探索.