基于黏彈性介質(zhì)的隧道地震二維正演模擬

凌 飛，肖宏躍，朱夏樂(lè)，雷 宛，蔣亞?wèn)|

(成都理工大學(xué) a.地球物理學(xué)院; b.環(huán)境與土木工程學(xué)院，成都 610059)

基于黏彈性介質(zhì)的隧道地震二維正演模擬

凌 飛a，肖宏躍a，朱夏樂(lè)a，雷 宛b，蔣亞?wèn)|a

(成都理工大學(xué) a.地球物理學(xué)院; b.環(huán)境與土木工程學(xué)院，成都 610059)

針對(duì)隧道中各種不良地質(zhì)體的超前預(yù)報(bào)問(wèn)題，以黏彈性介質(zhì)為基礎(chǔ)，在空氣隧道周圍添加一定厚度的開(kāi)挖損傷帶，建立符合實(shí)際情況的觀測(cè)系統(tǒng)，通過(guò)交錯(cuò)網(wǎng)格有限差分正演，模擬地震波在斷層、軟弱夾層、溶洞等地質(zhì)災(zāi)害體模型中的傳播，根據(jù)不同時(shí)刻的波場(chǎng)快照來(lái)研究相應(yīng)地質(zhì)模型地震記錄的響應(yīng)特征。通過(guò)對(duì)隧道地震二維波場(chǎng)模擬結(jié)果進(jìn)行分析，研究結(jié)果表明：對(duì)一階速度-應(yīng)力黏彈性波動(dòng)方程，采用交錯(cuò)網(wǎng)格有限差分，并且對(duì)隧道添加約束條件，可以快速、有效地實(shí)現(xiàn)隧道地震的二維正演模擬，并且波場(chǎng)快照能夠清楚地反映地震波在不同地質(zhì)模型中的傳播特征以及波的轉(zhuǎn)換特征，地震記錄能夠反映不同地質(zhì)模型的地震響應(yīng)特征，可以為實(shí)際隧道超前預(yù)報(bào)提供理論依據(jù)。

黏彈性介質(zhì)；隧道地震；交錯(cuò)網(wǎng)格；開(kāi)挖損傷帶；波場(chǎng)模擬

1 研究背景

地下隧道施工過(guò)程中，由于地質(zhì)條件的變化，常易發(fā)生突發(fā)性的意外情況。據(jù)國(guó)內(nèi)隧道施工的不完全統(tǒng)計(jì)，施工過(guò)程中因塌方、涌水、突泥、巖爆、高瓦斯等地質(zhì)災(zāi)害事故造成的停工時(shí)間大約占總工期的30%。通過(guò)隧道地質(zhì)超前預(yù)報(bào)，及時(shí)預(yù)報(bào)隧道掌子面前方的不良地質(zhì)體及其性質(zhì)、位置、規(guī)模、產(chǎn)狀與成災(zāi)可能性，預(yù)防涌水、突泥、坍塌等災(zāi)害性事故的發(fā)生，為正確選擇開(kāi)挖斷面、支護(hù)設(shè)計(jì)參數(shù)和優(yōu)化施工方案提供依據(jù)。以工程地質(zhì)分析法為基礎(chǔ)，以地球物理探測(cè)方法為手段的隧道綜合地質(zhì)超前預(yù)報(bào)方法已成為隧道工程界的共識(shí)，并得到廣泛的應(yīng)用。在隧道超前預(yù)報(bào)中，常用的地震波法有：負(fù)視速度法、水平剖面法(HSP)、隧道地震反射波超前預(yù)測(cè)法(TSP)、真反射層析成像法(TRT)、隧道地震波層析成像法(TST)、隧道地質(zhì)超前預(yù)報(bào)(TGP)等[1]。由于隧道超前地質(zhì)預(yù)報(bào)的特點(diǎn)決定了現(xiàn)場(chǎng)工作必須在狹小的隧道施工現(xiàn)場(chǎng)實(shí)時(shí)完成，因此，超前探測(cè)技術(shù)必須滿足狹小場(chǎng)地條件下進(jìn)行長(zhǎng)距離、高精度、高效率探測(cè)要求。大量的工程實(shí)踐表明，隧道開(kāi)挖過(guò)程中最常見(jiàn)的不良地質(zhì)體是斷層破碎帶、巖溶、軟弱夾層、采空區(qū)等。有不少隧道處在復(fù)雜巖性、構(gòu)造及高地應(yīng)力地區(qū)，尤其是在復(fù)雜的巖溶地區(qū)，因受巖性、地層組合、構(gòu)造、地貌等諸多因數(shù)的影響控制，巖溶發(fā)育復(fù)雜多變、大小不一、形態(tài)各異，給勘察工作帶來(lái)更大難度[2]。研究和了解這些地質(zhì)體的地質(zhì)地球物理特征將對(duì)于準(zhǔn)確的探測(cè)與預(yù)報(bào)十分關(guān)鍵。

近年來(lái)，許多學(xué)者在隧道地震的數(shù)值模擬方面做了相關(guān)的研究[3-7]。杜立志[3]對(duì)不同傾角的地層進(jìn)行了波場(chǎng)模擬，但并未將隧道作為空氣介質(zhì)考慮；魯光銀等[5]將隧道空間作為自由邊界，劉江平等[4]對(duì)隧道添加吸收邊界條件，兩者均對(duì)不同的地質(zhì)模型進(jìn)行了正演模擬，但均未進(jìn)行黏彈性介質(zhì)的數(shù)值模擬。

為了模擬結(jié)果更接近實(shí)際地質(zhì)情況，本文以開(kāi)爾芬黏彈性介質(zhì)為基礎(chǔ)，從一階速度-應(yīng)力黏彈性波動(dòng)方程出發(fā)，采用交錯(cuò)網(wǎng)格有限差分[8]對(duì)含開(kāi)挖擾動(dòng)帶的復(fù)雜隧道介質(zhì)模型進(jìn)行二維正演模擬。對(duì)模擬結(jié)果進(jìn)行對(duì)比分析，為實(shí)際隧道超前預(yù)報(bào)提供理論依據(jù)及指導(dǎo)。

2 黏彈性介質(zhì)波動(dòng)方程及其有限差分格式

在非均勻各向同性黏彈性介質(zhì)中，假設(shè)不考慮外力項(xiàng)，二維開(kāi)爾芬(Kelvin)黏彈性波動(dòng)方程的一階速度-應(yīng)力方程形式如下[9]：

(1)

式中：B是密度ρ=ρ(x,z)的倒數(shù)，或者叫做疏度；vx(x,z,t)和vz(x,z,t)分別是x和z方向的速度分量；σxx(x,z,t)，τxz(x,z,t)和σzz(x,z,t)是應(yīng)力分量；λ和μ為均勻彈性各向同性介質(zhì)的拉梅系數(shù);λ′和μ′為黏性介質(zhì)的拉梅系數(shù)。式(1)中拉梅系數(shù)為：

(2)

式中：vp和vs分別為介質(zhì)的縱波和橫波速度；ω為圓頻率；Qp，Qs分別為縱波和橫波的品質(zhì)因子。

(3)

(4)

(5)

(6)

(7)

(8)

(9)

3 隧道約束條件及邊界條件

由于超前預(yù)報(bào)是基于空氣隧道介質(zhì)模型來(lái)進(jìn)行，所以正演模擬時(shí)，需要對(duì)隧道添加約束條件。首先是對(duì)隧道中空氣密度的處理。當(dāng)數(shù)值計(jì)算至空氣隧道壁時(shí)，因?yàn)榧羟心Ａ客蝗幌陆抵?，并且壓力會(huì)下降一個(gè)數(shù)量級(jí)，這將導(dǎo)致數(shù)值模擬發(fā)生錯(cuò)誤。為了避免這種錯(cuò)誤的發(fā)生，需要將隧道中空氣的密度設(shè)置成一個(gè)較大值，本文將其密度設(shè)置成150 kg/m3。

其次是在隧道周圍添加開(kāi)挖損傷帶(EDZ)。由于隧道在開(kāi)挖過(guò)程中，會(huì)對(duì)隧道壁周圍的巖體造成損傷。通常情況下，越靠近隧道壁的圍巖，其損傷程度越大，而檢波器則需要安置在未受損傷的圍巖中。圍巖的損傷程度取決于許多因素，如開(kāi)挖方法、巖體力學(xué)特性、初始地應(yīng)力場(chǎng)、裂隙分布等[10]。本文將隧道損傷帶的寬度設(shè)置為0.5 m。

最后，對(duì)隧道里的數(shù)值模擬采用時(shí)間二階、空間四階的彈性波有限差分公式[5]進(jìn)行計(jì)算，并且進(jìn)行吸收邊界處理[4]。

對(duì)于模型邊界條件的處理，本文采用了PML完全匹配層吸收邊界[11]，可以較好地吸收傳播到模型邊界的地震波。

4 觀測(cè)系統(tǒng)及地質(zhì)模型的建立

為了使模擬結(jié)果更接近實(shí)際情況，本文建立了與實(shí)際TSP隧道超前預(yù)報(bào)相應(yīng)的觀測(cè)系統(tǒng)以及各種地質(zhì)模型。正演模擬采用的觀測(cè)系統(tǒng)如圖1所示，隧道長(zhǎng)50 m，寬8 m，震源位于隧道壁一側(cè)，沿隧道壁一共安置了24個(gè)檢波點(diǎn)，最小偏移距為10 m，道間距為1.5 m，震源埋深2 m，檢波點(diǎn)埋深1.5 m。

圖1 觀測(cè)系統(tǒng)示意圖

本文建立了3類地質(zhì)模型，即斷層、軟弱夾層以及溶洞模型，如圖2所示。通過(guò)建立不同類型的模型來(lái)了解并研究地震波在不同地質(zhì)體中的傳播規(guī)律，以及地震記錄的響應(yīng)特征；建立不同傾角(模型分界面與隧道中心軸線的夾角)的斷層和軟弱夾層模型來(lái)研究隨界面傾角變化波場(chǎng)以及地震記錄的變化規(guī)律。本文建立了界面傾角為90°和60°的2種地質(zhì)模型，模型的分界面距離掌子面80 m，軟弱夾層的寬度為10 m。通過(guò)建立不同物理特性參數(shù)的溶洞模型來(lái)研究不同充填物質(zhì)對(duì)地震記錄的影響。溶洞位于隧道掌子面正前方，溶洞中心距離掌子面83 m，其半徑為3 m。模型具體參數(shù)見(jiàn)表1。

圖2 地質(zhì)模型

模型材料Vp/(m·s-1)Vs/(m·s-1)ρ/(kg·m-3)QpQs隧道3400．1150圍巖30001732．020005040EDZ22001270．016004030斷層25001443．018004535軟弱夾層20001154．014003022溶洞11500866．010003020溶洞23400．1150105

圖3 直立斷層模型波場(chǎng)快照

5 波場(chǎng)數(shù)值模擬

通過(guò)添加隧道約束條件，運(yùn)用黏彈性介質(zhì)中的差分格式，針對(duì)建立的各種不良地質(zhì)體實(shí)現(xiàn)二維波場(chǎng)模擬。本文模型大小為250 m×100 m，網(wǎng)格間隔Δx=Δz=0.5 m，采樣間隔Δt為0.05 ms，地震記錄長(zhǎng)度為175 ms，震源采用200 Hz的Ricker子波。

5.1 斷層模型模擬結(jié)果及分析

對(duì)直立斷層模型不同時(shí)刻的波場(chǎng)快照(圖3所示)分析，當(dāng)t=35 ms時(shí)，可以觀察到由震源激發(fā)產(chǎn)生的直達(dá)縱波(P)、橫波(S)、沿隧道壁傳播的瑞利面波(R)和繞射波。瑞利面波在此刻到達(dá)隧道掌子面，其部分能量轉(zhuǎn)換成橫波(RS)繼續(xù)向前傳播，也有部分能量的瑞利面波會(huì)沿隧道壁往后傳播。t=55 ms時(shí)，直達(dá)縱波(P)在斷層分界面處產(chǎn)生了反射波和透射波，此時(shí)可以清楚地看到反射縱波(PP)和透射縱波(TP)。t=75 ms時(shí)，由于縱橫波存在波速差異，導(dǎo)致縱橫波之間的距離逐漸增加，此時(shí)可以明顯看到反射轉(zhuǎn)換橫波(PS)和透射轉(zhuǎn)換橫波(TPS)。t=95 ms時(shí)，橫波(S)已到達(dá)斷層分界面并產(chǎn)生了反射橫波(SS)。

在分析地震記錄響應(yīng)特征的過(guò)程中，先確定不同波的同相軸，再比較相同道數(shù)記錄的峰值大小，從而判斷不同波之間的能量強(qiáng)弱關(guān)系。圖4和圖5分別是在彈性和黏彈性介質(zhì)中對(duì)直立斷層模型做正演模擬得到的地震記錄。與彈性介質(zhì)相比，黏彈性介質(zhì)中的反射波能量明顯減弱，原因是黏彈性介質(zhì)考慮了地層對(duì)地震波能量的吸收作用，使模擬更接近實(shí)際地質(zhì)情況。

圖4 直立斷層模型地震記錄(彈性介質(zhì))

圖5 直立斷層模型地震記錄(黏彈性介質(zhì))

圖6是傾角為60°斷層模型的地震記錄，觀察黏彈性介質(zhì)中2種傾角斷層模型的地震記錄發(fā)現(xiàn)，t≤60 ms時(shí)，地震記錄主要由以下幾種波構(gòu)成：直達(dá)縱波(P)、由隧道壁以及開(kāi)挖損傷帶產(chǎn)生的反射波以及面波(R)，這幾種波的能量都很強(qiáng)。t>60 ms時(shí)，地震記錄主要反映來(lái)自斷層分界面的反射波。對(duì)比圖5中直立斷層模型的x和z分量地震記錄，發(fā)現(xiàn)x分量記錄中反射縱波(PP)能量強(qiáng)于z分量，而z分量記錄中反射轉(zhuǎn)換橫波(PS)和反射橫波(SS)的能量均強(qiáng)于x分量。對(duì)比60°傾角斷層模型的x和z分量地震記錄，發(fā)現(xiàn)z分量記錄中反射轉(zhuǎn)換橫波(PS)和反射橫波(SS)的能量依然強(qiáng)于x分量。分別對(duì)比2種斷層模型的x分量和z分量記錄，較直立斷層模型，60°傾角斷層模型的x分量反射縱波(PP)能量增強(qiáng)，但是反射橫波(SS)能量減弱；其z分量記錄反射縱波(PP)和反射轉(zhuǎn)換橫波(PS)能量均增強(qiáng)，但是反射橫波(SS)能量減弱。

圖6 傾角60°斷層模型地震記錄

5.2 軟弱夾層模型模擬結(jié)果及分析

不同傾角的軟弱夾層模型地震記錄如圖7和圖8所示，x分量和z分量記錄均可以看出地震波到達(dá)軟弱夾層2個(gè)分界面產(chǎn)生的反射波，x分量記錄的反射縱波能量要強(qiáng)于z分量記錄，而z分量記錄的反射橫波能量要強(qiáng)于x分量記錄。此外，2個(gè)分量記錄中，2條反射縱波同相軸的波形反向(由PP1和PP2表示)，這也驗(yàn)證了在堅(jiān)硬圍巖中間存在著一軟弱夾層。對(duì)比2個(gè)軟弱夾層模型的x分量和z分量記錄，對(duì)于傾角為60°的軟弱夾層模型而言，其x分量和z分量記錄中，反射縱波和反射轉(zhuǎn)換縱波能量增強(qiáng)，但是反射橫波及反射轉(zhuǎn)換橫波能量減弱。

圖7 直立軟弱夾層模型地震記錄

圖8 傾角60°軟弱夾層模型地震記錄

圖9 溶洞模型1(水)地震記錄

圖10 溶洞模型2(空氣)地震記錄

5.3 溶洞模型模擬結(jié)果及分析

溶洞模型1全部充填水，溶洞模型2全部充填空氣，空氣溶洞的物理特性參數(shù)與隧道一樣。圖9和圖10為2種不同填充物質(zhì)溶洞的x和z分量地震記錄。對(duì)于2種溶洞模型的x分量記錄而言，反射縱波(PP)能量明顯強(qiáng)于反射橫波(SS)能量，富水溶洞的反射縱波(PP)和反射橫波(SS)能量均弱于空氣溶洞，并且富水溶洞其反射波的疊加現(xiàn)象比空氣溶洞明顯。對(duì)于z分量記錄而言，反射橫波(SS)能量突出，反射縱波(PP)能量減弱。富水溶洞的反射縱波(PP)和反射橫波(SS)能量依然弱于空氣溶洞，但觀察初至的反射橫波(SS)之后的反射波，富水溶洞的反射波能量強(qiáng)于空氣溶洞。

6 結(jié) 論

(1) 對(duì)一階速度-應(yīng)力黏彈性波動(dòng)方程，采用交錯(cuò)網(wǎng)格有限差分，模型外部邊界采用PML完全匹配層吸收邊界條件，對(duì)隧道添加約束條件，能夠快速、準(zhǔn)確以及有效地實(shí)現(xiàn)隧道地震的二維正演模擬；波場(chǎng)快照能夠清楚地反映地震波在不同地質(zhì)模型中的傳播特征以及波的轉(zhuǎn)換特征，地震記錄能夠反映不同地質(zhì)模型的地震響應(yīng)特征，可以為實(shí)際隧道超前預(yù)報(bào)提供理論依據(jù)。

(2) 由于在黏彈性介質(zhì)中地層對(duì)彈性波的吸收衰減作用，使地震波能量減弱，更符合實(shí)際地質(zhì)的情況。通過(guò)對(duì)隧道地震模擬結(jié)果進(jìn)行分析，可以更好地了解地震波在地層中的傳播和衰減規(guī)律。

(3) 對(duì)本文建立的坐標(biāo)系而言，x分量記錄中，反射縱波能量突出，有利于對(duì)P波的提取與研究；z分量記錄中，反射橫波能量突出，有利于對(duì)S波的提取與研究。

(4) 對(duì)于不同傾角的地質(zhì)模型，當(dāng)?shù)刭|(zhì)界面傾角減小時(shí)，反射縱波能量增強(qiáng)，反射橫波能量減弱；對(duì)于軟弱夾層，夾層第1個(gè)分界面的反射波表現(xiàn)為正視速度特征，而第2個(gè)界面的反射波表現(xiàn)為負(fù)視速度特征；對(duì)于溶洞，富水溶洞反射波能量弱于空氣溶洞，其反射橫波疊加現(xiàn)象比空氣溶洞明顯。

[1] 趙永貴. 國(guó)內(nèi)外隧道超前預(yù)報(bào)技術(shù)評(píng)析與推介[J]. 地球物理學(xué)進(jìn)展，2007，22(4)：1344-1352.(ZHAO Yong-gui. Analysis and Recommendation of Tunnel Prediction Techniques at Home and Abroad[J]. Progress in Geophysics，2007，22(4)：1344-1352. (in Chinese))

[2] 宋先海,顧漢明,肖柏勛. 我國(guó)隧道地質(zhì)超前預(yù)報(bào)技術(shù)述評(píng)[J]. 地球物理學(xué)進(jìn)展，2006，21(2)：605-613. (SONG Xian-hai, GU Han-ming, XIAO Bo-xun. Overview of Tunnel Geological Advanced Prediction in China[J]. Progress in Geophysics，2006，21(2)：605-613.(in Chinese))

[3] 杜立志. 隧道施工地質(zhì)地震波法超前探測(cè)技術(shù)研究[D]. 長(zhǎng)春：吉林大學(xué)，2008. (DU Li-zhi. Study of Tunnel Seismic Geological Prediction Technology in the Process of Construction[D]. Changchun：Jilin University，2008. (in Chinese))

[4] 劉江平，程 飛，范承余，等. 基于隧道空間全波場(chǎng)二維數(shù)值模擬與特征分析[J]. 巖土工程學(xué)報(bào)，2012，34(9)：1705-1711.(LIU Jiang-ping, CHENG Fei, FAN Cheng-yu,etal. Two-dimensional Numerical Simulation of Tunnel-based Seismic Full-wave Fields[J]. Chinese Journal of Geotechnical Engineering，2012，34(9)：1705-1711.(in Chinese))

[5] 魯光銀，熊 瑛，朱自強(qiáng). 隧道反射波超前探測(cè)有限差分正演模擬與偏移處理[J]. 中南大學(xué)學(xué)報(bào)(自然科學(xué)版)，2011，42(1)：136-141. (LU Guang-yin, XIONG Yin, ZHU Zi-qiang. Detection Simulation ahead of Tunnel Face and Reverse-time Migration with Reflection Wave Method[J]. Journal of Central South University (Science and Technology)，2011，42(1)：136-141. (in Chinese))

[6] KNEIB G, LEYKAM A. Finite-difference Modelling for Tunnel Seismology[J]. Near Surface Geophysics，2004，2(2)：71-93.

[7] BOHLEN T, LORANG U, RABBEL W,etal. Rayleigh-to-shear Wave Conversion at the Tunnel Face: From 3D-FD Modeling to Ahead-of-drill Exploration[J]. Geophysics，2007，72(6):67-79.

[8] ROBERTSSON J O A, BLANCH J O, SYMES W W. Viscoelastic Finite-difference Modeling[J]. Geophysics, 1994, 59(9):1444-1456.

[9] 奚 先,姚 姚. 二維粘彈性隨機(jī)介質(zhì)中的波場(chǎng)特征分析[J]. 地球物理學(xué)進(jìn)展，2004，19(3)：608-615. (XI Xian，YAO Yao. The Analysis of the Wave Field Characteristics in 2-D Viscoelastic Random Medium[J]. Progress in Geophysics，2004，19(3)：608-615. (in Chinese))

[10]吉小明. 隧道開(kāi)挖的圍巖損傷擾動(dòng)帶分析[J]. 巖石力學(xué)與工程學(xué)報(bào)，2005，24(10)：1697-1702. (JI Xiao-ming. Study on Mechanical and Hydraulic Behavior of Tunnel Surrounding Rock Masses in Excavation-Disturbed Zone[J]. Chinese Journal of Rock Mechanics and Engineering，2005，24(10)：1697-1702. (in Chinese))

[11]COLLINO F, TSOGKA C. Application of the Perfectly Matched Absorbing Layer Model to the Linear Elastodynamic Problem in Anisotropic Heterogeneous Media[J]. Geophysics，2001，66(1):294-307.

(編輯：黃 玲)

Two-dimensional Forward Seismic Modeling for TunnelBased on Viscoelastic Medium

LING Fei1，XIAO Hong-yue1，ZHU Xia-le1，LEI Wan2，JIANG Ya-dong1

(1.College of Geophysics，Chengdu University of Technology，Chengdu 610059，China； 2.College of Environment and Civil Engineering, Chengdu University of Technology, Chengdu 610059,China)

In view of the problems in advanced detection of defective geological structures in tunnel, a geological tunnel model was established based on viscoelastic medium to simulate the seismic wave propagation in fault, weak interlayer, karst cave and other geological hazard body. The tunnel is surrounded by the excavation damaged zone. A practical observation system and the staggered grid finite-difference method were employed for the simulation. According to the wave field snapshots at different instants, the response of seismic record acquired from the corresponding forward modelling was researched. Result shows that the two-dimensional forward simulation based on viscoelastic medium can be achieved by using a staggered grid finite-difference algorithm and adding tunnel constraints. The wave field snapshots could clearly reflect the seismic wave propagation and transformation in different geological models. The seismic record reflects the seismic response characteristics of different models, and offers a theoretical basis for the seismic prediction of actual tunnels.

viscoelastic medium; tunnel seismic; staggered grid; excavation damaged zone; wave field simulation

2014-07-07；

2014-08-18

凌 飛(1989-)，男，四川成都人，碩士，主要研究方向?yàn)楣こ膛c環(huán)境地球物理勘探，(電話) 028-84073510(電子信箱)825729075@qq.com。

10.3969/j.issn.1001-5485.2015.05.023

2015,32(05):121-126

P631.4

A

1001-5485(2015)05-0121-06