淺部地下洞穴的瑞利波散射數值分析

王啟淵，柴華友，丁 瑋

武漢工程大學資源與土木工程學院，湖北 武漢 430074

隨著我國經濟的飛速發展，基礎建設也隨之得到了巨大的發展，尤其是地鐵、隧道等地下工程。近年來，因地下巖溶導致施工過程中地面塌陷，基礎不均勻沉降，造成施工機械設備被埋、人員傷亡等事故時有發生。巖土工程勘查中對地下洞穴的探測通常是大范圍的，常規的地質勘察手段難以同時兼顧效率和成本。波動測試是近年來新發展的一種物探方法［1-2］，Stokoe等［3］采用表面波譜 分 析 方 法（spectral analysis of surface wave，SASW）對高速公路路面及路基進行了探測，為瞬態瑞利波法在工程中的廣泛應用奠定了基礎。

對表面震源，瑞利波能量主導淺部波場，瑞利波遇遭淺部地下異質體會發生散射［4］。表面震源在介質中會產生體波和瑞利波，當介質中有洞穴存在時，體波和瑞利波都會發生散射［5-9］，散射波會產生頻散和譜擾動現象，其擾動的程度及形態與散射波中各類型波的能量及傳播特性有關。近年來，一些學者嘗試利用瑞利波在彈性半空間介質中的傳播特性，將瑞利波法應用于探測淺部地下洞穴。卞鵬等［10］將通過研究傳感器位置和間距對地下障礙物探測的影響，分析含洞穴半空間介質頻散曲線的特征，為探測地下障礙物時頻散曲線參數的變化研究做工作。柴華友等［11］通過提取激發震源后的介質表面質點的速度響應，從能量變化的角度分析洞穴對相速度擾動機制。劉中憲等［12］使用間接邊界積分方法對入射瑞利波在襯砌隧道中的散射現象進行了分析，發現了襯砌剛度對波散射的影響，同時發現隨著隧道埋深的加大，波的散射現象逐漸減弱。邵廣周等［13］基于洞穴衍射現象的研究提出了利用廣義S轉換的方法探測洞穴的埋深。蔣輝等［14］利用散射波現象檢驗地下洞穴的注漿效果，以波速的變化為判斷依據，通過研究表明波速速度及變化幅度越大，注漿的效果越好。張獻民等［15］通過有限元法進行數值模擬，總結了瑞利波信號受地下管道的形狀和管徑變化影響的相應規律。柴華友等［16］通過數值分析的方法分析散射波對表面波相速度的影響，由相速度曲線的變化確定異質體的參數。但相速度擾動影響因素很多，對各種干擾較敏感，相速度擾動與埋深規律難以判斷。

基于以上分析，研究散射波能量及傳播特性與洞穴埋深及幾何參數間關系有助于分析波場擾動特征，推動瑞利波在淺部洞穴探測中的應用。研究洞穴的形態特征主要考慮波在傳播過程中遇到洞穴時產生的入射、反射現象，得到散射波場的特征。通過散射波理論可以對介質中的洞穴進行分析，由于解析解以積分形式給出，需迭代求解，求解過程繁瑣，不利于研究洞穴散射波基本特征，因此使用數值模擬方法建立二維模型，分析洞穴對瑞利波的散射特征，是十分必要的。將數值模擬二維模型結果進行分析，得到相關特征，用于三維現象的分析中。本文主要從二維模型的角度分析洞穴對瑞利波散射特征。針對這一目的，本文主要開展以下工作：（a）利用有限元模擬計算均勻半空間中存在不同埋深及幾何參數洞穴情況下的散射波場；（b）將頻率域轉換成波長域，在偏移距-波長域振幅譜圖中研究波場中的能量變化及傳播特性與擾動特征；（c）建立波長與地下洞穴埋深及幾何參數間的關系。

1 有限元模型的建立與驗證

1.1 模型類型

采用平面及軸對稱模型研究瑞利波的散射，復雜邊界形狀的洞穴可近似用一組小矩形邊界洞穴組合建模，因此，采用矩形洞穴模型具有代表性。為了便于分析表面波場，取介質模型長25 m，寬20 m。在波的傳播過程中，當傳遞至人工模型邊界時，波會產生反射現象，有可能疊加到行進波中，對分析結果產生誤差。為了減少模型邊界產生的散射波，在有限元模型邊界施加無反射邊界，取模型土介質性質參數為：密度 ρ＝1 800 kg/m3，泊松比υ＝0.35，剪切波速Cs=250 m/s，具體模型如圖1所示。

圖1 模型示意圖Fig.1 Schematic diagram of numerical model

1.2 震源

在地震波研究中，如果子波在時域上很窄，則很容易從波場中識別出散射波。因此，狹窄的子波在時域上具有很高的分辨率。在頻域中，如果能量集中在一個狹窄的范圍內，有利于分析散射波的頻率行為。Ricker子波函數具有這些特性，取瞬態沖擊荷載變化用Ricker表示，其時域及頻率域表達式：

式中t0為峰值對應時間，f=1/t0為中心頻率（圖中箭頭表示位置）。假設 t0=0.01 s，則 fm=100 Hz，Ricker源及其頻譜如圖2所示。可以看出，Ricker源是由一個正值域和兩個負值域組成。

圖 2 Ricker函數：（a）時域；（b）頻域Fig.2 Ricker function：（a）Time domain；（b）Frequency domain

1.3 模型合理性驗證

蘭姆給出彈性均勻半空間的瑞利波位移解析解，通過數值模擬Lamb源下響應與解析解比較可驗證有限元模型、網格劃分及無反射邊界有效性。Lamb點荷載在時域的表達式：

頻域表達式為：

式中參數P及q用于調整Lamb源的振幅及頻域能量分布。

Lamb源的瑞利波的水平（徑向）和垂直位移近似表示為：

式中φ=arctan[(t-r/cR)/q]，K和H是與橫波、縱波和瑞利波速相關的常數，cR是瑞利波速度，μ是剪切模量，r是距離源的水平（徑向）距離（即偏移量）。

峰值t0的Lamb源可寫為：

根據傅里葉變換的時間位移規律，對應于上述響應的頻譜為：

式（5）、式（6）中的 φ變為：

在有限元模型中，Lamb源參數設置為：P=1，q=1×10-3，t0=16ms。解析解與數值解在偏移距為0.01 m、0.05 m、1 m的對比如圖3所示。由圖3可知，這兩種解的誤差很小，表明本文模型、網格劃分及無反射邊界是合理、有效的。

圖3 不同偏移距質點速度響應的數值解與解析解比較Fig.3 Comparison of numerical and analytical solutions of particle velocity responses at different offsets

2 洞穴對瑞利波散射現象

激發表面震源在介質中會產生直達波和瑞利波（R波），介質中若沒有洞穴，波不會發生散射現象。圖4（a）為表面豎直向質點速度等值云圖，可以看出表面波場由瑞利波主導。圖4（b）為某一時刻截面質點速度矢量幅值云圖，P、S波在介質體內傳播，R波沿表面傳播，主要能量分布在淺部區域。

圖4 均勻半空間中質點速度響應：（a）豎直向質點速度時程圖；（b）波場質點速度幅值等值云圖Fig.4 Responses of particle velocity in homogenous half space：（a）Time history of vertical particle velocity；（b）Contour plot of particle velocityin wave filed

由于瑞利波能量分布較淺，淺部信號能量主要由瑞利波主導，當淺部出現洞穴時，瑞利波質點軌跡會發生改變，瑞利波發生散射，在時域及頻域，散射波對波場擾動，其擾動的程度及形態與散射波中各類型波的能量及傳播特性有關。瑞利波的能量通常在一個波長的長度內，且多數能量分布于半波長的長度內，圖5給出高頻（短波長）與低頻（長波長）瑞利波能量在深度方向分布以及能量分布深度與洞穴埋深及幾何參數對比。由圖5可以分析洞穴對不同波長瑞利波的散射程度。瑞利波與洞穴相遇后，波長小于洞穴埋深的波被洞穴散射的能量很少，大部分能量繼續以原有的速度在洞穴上方介質傳播。而對于波長大于埋深的波，與洞穴相互作用，一部分能量發生反向散射，一部分能量在洞穴上方透射，在洞穴上方，透射波不僅傳播能量會發生變化而且傳播速度、衰減規律也會發生變化。均勻半空間中有洞穴的模擬示意圖如圖6所示，通過與前部分無洞穴情況相對比，討論瑞利波在遭遇洞穴后產生的散射現象對波場的影響。由圖 6（a），（b），（c）不同時刻質點速度幅值云圖可以看出入射波在洞穴散射過程，對淺部洞穴，入射波主要是R波，在波場中R波會發生反射、波的模式轉換、透射等現象。部分R波經過波的模式轉換轉化為橫波和縱波，還有一部分R波會直接穿透洞穴，形成透射波。其表面豎直向質點速度幅值見圖7，圖7和圖4相比，在洞穴近界面和遠界面處有明顯的變化，無洞穴的介質中，瑞利波沿傳播方向傳播至模型底部被無反射邊界吸收，無明顯能量衰減，而在有洞穴介質中傳播，在洞穴邊界處有明顯的散射波出現，這說明使用散射瑞利波探測地下洞穴是可行的。

圖5 瑞利波能量分布深度與波長示意圖Fig.5 Illustration of energy distribution depth and wavelength of Rayleigh wave

圖6 不同時刻豎直向質點速度幅值云圖：（a）波傳播至空洞近界面；（b）波傳播至空洞上方；（c）波傳播至空洞遠界面Fig.6 Cloud contours of vertical particle velocity at different times：（a）Waves before the cavity；（b）Waves over the cavity；（c）Waves after the cavity

圖7 表面豎直向質點速度云圖Fig.7 Contour of surface vertical particle velocity

3 洞穴埋深及幾何參數對波散射影響

3.1 散射波場分析

對表面各點豎直向質點速度響應作一維傅里葉變換得到頻譜圖，為了便于與深度量綱一致及消除波幾何擴散導致波衰減，將每個質點響應的頻率域譜轉化為波長域譜并將幅值進行歸一化得到偏移距-波長域振幅譜圖，通過對振幅譜能量變化分析預估洞穴的形狀，當均勻半空間中存在洞穴時，由前文可知瑞利波在遭遇洞穴后會發生散射現象，散射波在不同的區域表現出不同的特征。

圖8（a）給出介質無洞穴情形下的表面波場譜圖，不同位置譜能量相似，波場譜均勻，無擾動。當介質出現洞穴時，表面波場譜圖如圖8（b）所示，在此情況下，波在遇到洞穴前傳播狀態與均勻介質類似，當波傳遞至洞穴近界面處（b圖中Ⅰ處）時，瑞利波產生反射、繞射、透射等現象，導致在洞穴的近界面有明顯的能量區分界線。反射、透射使能量產生損失，同時洞穴近界面和頂面的相接拐點處（b圖中Ⅱ處）會使能量產生二次散射，導致在拐點處能量紊亂形成不規則的波場。不僅是洞穴頂部的拐點，底部與遠近界面的拐點也會產生這種現象，由于底部拐點，大部分能量向下散射，底部拐點附近散射波對表面波場的影響很小，導致很難判斷洞穴的深度。當波行進至洞穴上方（b圖中Ⅲ處），波在介質表面和洞穴頂層來回反射，反射波相互干涉，能量或加強或抵消。當上方透射波到達遠交界的拐點（b圖中Ⅳ處）同樣會產生散射現象，波的散射和繞射使上方能量明顯減少。在遠界面附近（b圖中Ⅴ處），透射波隨著傳播距離增加，傳播特性趨于入射瑞利波，能量分布深度逐漸趨于未遭遇洞穴前。通過在振幅能量譜圖中的特點，可以判斷出洞穴的大致特征。

圖8 偏移距-波長域表面豎直向質點振動速度振幅譜圖：（a）無洞穴；（b）有洞穴Fig.8 Amplitude spectra of vertical particle velocity of surface in offset-wavelength domain ：（a）Without cavity；（b）With cavity

3.2 散射波對表面波場擾動分析

下面以洞穴埋深、位置及幾何尺寸不同的幾種算例來驗證以上分析，主頻 fm=100Hz時主要瑞利波能量分布在波長0.5 m～2.5 m范圍，見圖8（a），為了避免埋深較小或較大洞穴的微弱的散射波被數值振蕩和頻譜泄漏引起的噪聲掩沒，取洞穴的埋深1.5 m～3 m，洞穴具體參數如表1所示，根據不同算例繪制偏移距-波長域質點振動速度振幅譜圖，如圖9所示。

基于波可能產生的散射行為，將研究區域分為洞穴前方（震源至洞穴近界面）、洞穴上方、洞穴后方（洞穴遠界面至遠場）三個區域。通過前文所述的洞穴形態判斷方法對由洞穴造成波場擾動的三個區域波場，建立波長與洞穴埋深的關系。

1）洞穴前方：波未傳播至洞穴近界面時，其傳播方式與均勻半空間中相同。當波與洞穴近界面遭遇，發生散射現象，洞穴近界面區域的反射波與

前行波相互干涉。算例 a、b、c、f研究洞穴高度、水平位置相同時埋深參數對前方散射波的影響，對應譜圖中特征并結合前文對洞穴近界面洞穴的判斷方法發現，前方出現條紋狀衰減，產生形狀規則的能量分界面。研究表明，不同埋深的洞穴前仍存在能量擾動，出現類似于能量加強或遞減的變化，這表明，只要洞穴位于瑞利波的能量影響深度內，不同埋深都會產生散射現象。算例c和d為埋深、水平位置相同但洞穴高度發生變化時洞穴對瑞利波散射的影響，將對應譜圖進行對比分析發現，洞穴高度增加，但是高度的下半部分散射向表面的能量較小。同時，在底部與近界面拐點處，波大部分能量向介質體內散射，一部分經過繞射至洞穴底部，向表面波場散射較少。因此，高度的影響對表面波場相對較小。算例c和e為同埋深同高度僅洞穴的水平位置變化時洞穴對瑞利波散射的影響，通過對應譜圖中能量變化情況結合前文所述判斷方法，譜圖中能量的變化趨勢是類似的，發現隨著洞穴水平位置的變化，只有能量變化區域的水平位置產生變化。

表1 洞穴模型位置和預估位置參數Tab.1 Actual and estimated locations of cavities

圖9 偏移距-波長域質點振動速度振幅譜圖（矩形實線表示空洞實際大小，虛線表示預估形狀，Ⅰ，Ⅱ，Ⅲ分別表示頂面，近界面，遠界面）Fig.9 Amplitude spectra of particle velocity in offset-wavelength domain（solid lines：actual shapes of cavities；dotted lines：estimated shapes，Ⅰ，Ⅱ and Ⅲ represent the top，near and far boundaries respectively）

2）洞穴上方：瑞利波通過洞穴近界面區域后，到達洞穴上方，繼續發生散射現象，散射波在自由表面和頂面之間來回反射。通過算例a、b、c、f研究同洞穴尺寸、水平位置但埋深不同時，洞穴上方波散射現象。當空洞在介質淺部時（算例a，埋深1.5 m），由于洞穴的埋深與介質表面距離較近波長較大的散射波無法通過，瑞利波具有波長與頻率成反比的特性，故能通過的波具有較高的頻率。由前文震源選取部分的描述可知，Ricker子波函數在高頻部分的頻率較集中，因此反射波在埋深與介質表面間發生干涉現象，在譜圖a中表現為洞穴頂部后半部分的能量明顯強于前半部分。當洞穴埋深加深后，與介質表面距離增大，干涉現象不再發生，故在算例b、c、f中這種特殊現象消失，在譜圖中受洞穴影響發生的能量變化特征與前文所述判斷方法一致。將算例c分別與算例d和e對比，發現洞穴的高度和質中水平位置的變化對洞穴導致的散射現象影響不大。

3）洞穴后方：上方透射波在洞穴遠界面繼續發生散射，遠界面附近的表面波場同樣產生擾動。當波離開洞穴遠界面區域后，其傳播行為逐漸趨于入射瑞利波傳播特性。在埋深較低時（1.5 m處），洞穴遠界面與頂部的拐點處沒有出現前文中洞穴形態預估中拐點處能量衰減區域，而是出現條紋狀的能量加深區域。這是由于頂部與介質表面間距離過短，散射波在洞穴上方發生干涉后在拐點處繼續發生干涉，在拐點處積聚也會使拐點處的能量加強，能量在繞過拐點后部分繞射至洞穴底部，部分繼續傳播。隨著埋深增加，沒有干涉現象產生，因此在算例b、c、f中都沒有出現算例a中的現象。將算例c和d相對比，當洞穴的埋深和水平位置不變高度變化時，與洞穴前的現象類似，洞穴的高度增加后，對表面波場造成的擾動不敏感，波在傳播至遠界面與底部交匯拐點時，大部分能量散射進洞穴內，一部分向下反射，一部分繼續傳播，傳播模式逐漸變為均勻半空間傳播方式。算例c和e為洞穴埋深和高度相同，水平位置變化的情況，對應譜圖進行分析，能量變化區域隨水平位移變化改變，但能量變化特征與前文所述一致。

3.3波場擾動與洞穴埋深及幾何參數的關系

波場能量擾動與洞穴埋深、幾何參數有關，根據能量變化程度可確定臨界擾動，由臨界擾動對應的波長可預估埋深。根據能量幅度變化形成的近界面、頂面和遠界面，可以判斷洞穴的形狀，建立λ與h的關系。預估洞穴位置參數如表1所示。

從波散射引起的擾動和能量模式的變化，可定位洞穴的遠近位置。但是從譜分析圖中，對于洞穴的高度很難估計。部分情況洞穴的遠近界面的估計位置比實際值大，這種情況可能是由于前向散射波在到達某個距離后干擾波的傳播，使反向散射波相對于遠近界面的前向散射波較弱。

通過洞穴模型位置和預估位置參數表及偏移距-波長域振幅譜圖，發現在洞穴的實際埋深（h）與波長（λ）存在一定的對應關系，各算例中無量綱幾何參數h/λ如表2所示。

表2 數值模擬中無量綱參數h/λ表Tab.2 Dimensionless parameters h/λ in numerical simulations

在埋深1.5 m～3 m的范圍內，無量綱幾何參數h/λ的變化范圍為 1.367～1.657，算例 a、b、c、f的數據表明隨深度變化，h/λ存在一定幅度的增長，但仍以1.5為中心浮動，算例c、d、e表明洞穴在大小和位置等影響因素下h/λ的值基本沒有變化，可以判斷洞穴的大小和在介質中的位置對h/λ沒有影響。綜上所述，在介質淺部存在空洞穴時，其擾動無量綱幾何參數h/λ介于1.3～1.6之間，即洞穴的埋深與偏移距-波長域振幅譜圖中洞穴對應波長的關系以1.5倍為中心上下浮動。

4 結 語

對均勻半空間介質中瑞利波的產生及傳播特征進行了分析，運用有限元軟件建立含洞穴均勻半空間模型，通過與解析解比較，驗證了有限元模型可行性，通過對不同埋深、幾何參數算例分析，得到如下結論：

1）當介質中存在洞穴時，瑞利波傳播至洞穴位置會發生散射，散射波導致表面波場能量發生擾動。

2）通過數值模擬波場質點速度響應分析了瑞利波在洞穴前方、上方、后方的散射現象及散射波對波場造成的擾動，由此從偏移距-波長域振幅譜圖的能量擾動變化判斷洞穴的形態，建立了波長與洞穴埋深的關系。

3）在偏移距-波長域，由振幅譜擾動，根據埋深與擾動特征波長的倍數關系可以預測洞穴埋深。

［1］ 欒明龍，魏紅，林萬順.瞬態瑞利波技術在工程勘察中的應用［J］. 物探與化探，2012，36（5）：878-883.

［2］ 孫笑，趙明階，汪魁，等.波動理論在巖土工程測試中的應用研究進展［J］.重慶交通大學學報（自然科學版），2013，32（1）：58-62.

［3］ STOKOE K H，NAZARIAN S.Effectiveness of ground improvement from spectral analysis of surface waves［C］//Proceeding of the eighth European conference on soil mechanics and foundation engineering，Helsinki，Finland.Netherlauds：AA Balkema Publishers，1983：91-95.

［4］ 劉優平，龔敏，徐斌.半空間飽和土中輸水管道對瑞利波的散射［J］.四川大學學報（工程科學版），2012，44（5）：71-77

［5］ 卞鵬，王媛，王篤強.有障礙物地基瑞利波頻散曲線的特征研究［J］.山東大學學報（工學版），2012，42（10）：99-103.

［6］ CHAI H Y，GOH S H，PHOON K K，et al.Effects of source and cavity depths on wave-fields in layered media［J］.Journal of Applied Geophysics，2014，107（8）：163-170.

［7］ CHAI H Y，GOH S H，PHOON K K，et al.Effects of source and cavity depths on wave fields in homogeneous half spaces［J］.Journal of Applied Geophysics，2013，93（2）：52-59.

［8］ CHAI H Y，PHOON K K，GOH S H，et al.Some theoretical and numericalobservations on scattering of Rayleigh waves in media containing shallow rectangularcavities［J］.Journal of Applied Geophysics，2012，83（6）：107-119.

［9］ CHAI H Y，PHOON K K，WEI C F，et al.Analysis of effects of active sourceson observed phase velocity based on the thin-layer method［J］.Journal of Applied Geophysics，2011，73（1）49-58.

［10］ CHAI H Y，WEI C F.Effects of cavities on the spectra of Rayleigh waves numerical aspects［C/CD］//HUANG A B，MAYNE P W.Geotechnical and Geophysical Site Characterization.Proceeding of the third International Conference on Site Charactorization，Taipei，Taiwan，April，1-4，2008.London，UK：Taylor&Francis Group，2008：807-711.

［11］ 柴華友，盧應發，李祺，等.半無限體介質內異質體方位探測數值研究［J］. 巖土力學，2007，28（1）：188-192

［12］ 劉中憲，梁建文，張賀.彈性半空間中襯砌隧道對瑞利波的散射［J］.巖石力學與工程學報，2011，30（8）：1627-1637.

［13］ SHAO G Z，TSOFLIAS G P，LI C J.Detection of near-surface cavities by generalized S-transform of Rayleigh waves［J］.Journal of Applied Geophysics，2016，129：53-65.

［14］ 蔣輝，柳超，賈超，等.地震散射方法在注漿效果檢測中的應用［J］.工程地球物理學報，2015，12（2）：260-265.

［15］ 張獻民，支雁飛.地下管道瞬態瑞雷波信號響應分析［J］. 工程勘察，2015，43（6）：81-87.

［16］ 柴華友，崔玉軍，盧應發，等.地下洞穴對表面波表觀相速度影響的數值分析［J］.巖石力學與工程學報，2006，25（5）：956-962.