魯光銀，羅帥，朱自強(qiáng)，石克亮，夏成志

基于柱坐標(biāo)系的隧道空間全波場(chǎng)數(shù)值模擬與分析

魯光銀，羅帥，朱自強(qiáng)，石克亮，夏成志

(中南大學(xué) 地球科學(xué)與信息物理學(xué)院，湖南 長(zhǎng)沙 410083)

不同于地面地震，隧道空間的特殊結(jié)構(gòu)使得其波場(chǎng)十分復(fù)雜，基于實(shí)際隧道模型的全波場(chǎng)數(shù)值模擬是有效認(rèn)識(shí)隧道空間內(nèi)波場(chǎng)傳播規(guī)律的重要手段。關(guān)于這方面的研究以往多是基于直角坐標(biāo)系，而隧道作為類似空氣柱的特殊地質(zhì)模型，采用柱坐標(biāo)系下的全波場(chǎng)數(shù)值模擬對(duì)認(rèn)識(shí)隧道空間的波場(chǎng)特征更有意義。基于柱坐標(biāo)系下的一階彈性波速度-應(yīng)力方程，推導(dǎo)任意偶數(shù)階交錯(cuò)網(wǎng)格差分格式，設(shè)計(jì)了人工截?cái)噙吔绾退淼雷杂蛇吔?，?duì)巖性分界面、斷層破碎帶、溶洞等隧道典型不良地質(zhì)模型進(jìn)行數(shù)值模擬與波場(chǎng)特征分析。研究結(jié)果表明：基于柱坐標(biāo)系下的一階彈性波速度-應(yīng)力方程，采用交錯(cuò)網(wǎng)格有限差分算法結(jié)合FCT技術(shù)，可實(shí)現(xiàn)隧道空間全波場(chǎng)高精度數(shù)值模擬，波場(chǎng)特征符合波的運(yùn)動(dòng)學(xué)與動(dòng)力學(xué)特征，為隧道空間復(fù)雜波場(chǎng)特征的認(rèn)識(shí)和資料解譯提供有效的依據(jù)。

超前預(yù)報(bào)；隧道空間；柱坐標(biāo)系；全波場(chǎng)；數(shù)值模擬

隨著交通運(yùn)輸?shù)目焖侔l(fā)展，越來(lái)越多的隧道被建造。然而，隧道作為一項(xiàng)隱蔽性地下工程，在勘察設(shè)計(jì)階段很難精確的發(fā)現(xiàn)不良地質(zhì)體。在隧道開(kāi)挖過(guò)程中，有可能遇到斷層、溶洞，暗河等不良地質(zhì)現(xiàn)象，極易誘發(fā)突水、涌泥、巖爆等地質(zhì)災(zāi)害，進(jìn)而造成工期延誤，經(jīng)濟(jì)損失，甚至人員傷亡等重大安全事故[1?2]。地震法作為隧道超前預(yù)報(bào)的主要方法之一，包括TSP，TST，HSP，TRT和陸地聲納法等[3?8]，以其長(zhǎng)距離、相對(duì)高精度而得以廣泛應(yīng)用。但是隧道空間內(nèi)的波場(chǎng)是十分復(fù)雜的，接收記錄上往往存在直達(dá)波、反射波、面波、縱波、橫波轉(zhuǎn)換波[9]，這無(wú)疑為數(shù)據(jù)處理與成果解譯增加了難度?；谒淼揽臻g的全波場(chǎng)數(shù)值模擬是認(rèn)識(shí)隧道空間內(nèi)復(fù)雜波場(chǎng)傳播規(guī)律的基礎(chǔ)，同時(shí)也是觀測(cè)系統(tǒng)設(shè)計(jì)，數(shù)據(jù)采集、處理、解釋的重要依據(jù)。近年來(lái)有不少學(xué)者進(jìn)行了相關(guān)研究，其中絕大部分是基于直角坐標(biāo)系，魯光銀等[10]采用直接法處理隧道自由邊界條件，并進(jìn)行了地質(zhì)界面模型的正演模擬與偏移成像；劉江平，程飛等[9, 11]采用吸收邊界條件處理隧道自由邊界，并進(jìn)行了不同觀測(cè)系統(tǒng)和不同傾角地質(zhì)界面對(duì)波場(chǎng)特征和偏移成像的影響研究；凌飛等[12]基于黏彈性介質(zhì)并在隧道周圍引入一開(kāi)挖損傷帶對(duì)隧道典型不良地質(zhì)現(xiàn)象進(jìn)行了全波場(chǎng)數(shù)值模擬；宋杰[13]針對(duì)隧道工程建設(shè)中可能遇到的典型不良地質(zhì)體，建立了巖性分界面、斷層以及巖溶體等多個(gè)三維地質(zhì)模型，進(jìn)行了三維地震波探測(cè)正演模擬，分析總結(jié)了隧道典型不良地質(zhì)體三維波場(chǎng)響應(yīng)和地震成像特征；王京[14]基于二階彈性波位移波動(dòng)方程采用有限元方法對(duì)隧道典型不良地質(zhì)體進(jìn)行了數(shù)值模擬，并重點(diǎn)探討了有限元方法實(shí)現(xiàn)自由邊界條件的處理方法；Jetschny等[15]考慮了隧道空間的影響，利用面波和轉(zhuǎn)換橫波對(duì)隧道前方的地質(zhì)構(gòu)造進(jìn)行預(yù)測(cè)，并且討論了這2種波的轉(zhuǎn)換機(jī)理；YANG等[16]基于三維彈性波方程交錯(cuò)網(wǎng)格有限差分算法，實(shí)現(xiàn)了煤礦巷道模型的三維槽波勘探數(shù)值模擬；Harmankaya等[17]實(shí)現(xiàn)了利用隧道掘進(jìn)過(guò)程中產(chǎn)生的噪聲作為震源進(jìn)行超前探測(cè)的數(shù)值模擬，并驗(yàn)證了該方法對(duì)散射體的定位能力；LIU等[18]進(jìn)行了極坐標(biāo)系下的全波場(chǎng)數(shù)值模擬，并處理了復(fù)雜的邊界條件；查欣杰等[19]通過(guò)構(gòu)建二維含低速異常的隧道介質(zhì)模型，研究了隧道彈性波場(chǎng)傳播規(guī)律和異常體邊界成像準(zhǔn)確性。模擬隧道空間三維波場(chǎng)最經(jīng)濟(jì)的方法是把構(gòu)造模型包括震源近似為軸對(duì)稱模型。而相較于直角坐標(biāo)系采用柱坐標(biāo)系進(jìn)行數(shù)值模擬具有以下2點(diǎn)優(yōu)勢(shì)[20]：1) 直角坐標(biāo)系中的有限差分算法的公式比較簡(jiǎn)單, 但是處理軸對(duì)稱地質(zhì)構(gòu)造需要比較多的計(jì)算成本；2) 對(duì)于具有軸對(duì)稱性的波場(chǎng)問(wèn)題, 可以轉(zhuǎn)化為二維問(wèn)題, 大大簡(jiǎn)化計(jì)算。論文基于柱坐標(biāo)系下的一階彈性波速度?應(yīng)力方程，將隧道與圍巖邊界作自由邊界處理，同時(shí)引入柱坐標(biāo)系下的分裂式完全匹配層，采用時(shí)間二階，空間十階的高階交錯(cuò)網(wǎng)格有限差分法并結(jié)合FCT技術(shù)，實(shí)現(xiàn)了柱坐標(biāo)系下隧道空間的高精度全波場(chǎng)數(shù)值模擬，通過(guò)計(jì)算巖性分界面、斷層破碎帶、溶洞等復(fù)雜地質(zhì)模型，分析了隧道典型不良地質(zhì)現(xiàn)象的波場(chǎng)特征，為隧道空間復(fù)雜波場(chǎng)特征的認(rèn)識(shí)和資料解譯提供了依據(jù)。

1 波動(dòng)方程及其有限差分離散

考慮到隧道空間的近似軸對(duì)稱結(jié)構(gòu)以及三維數(shù)值模擬的計(jì)算成本，本文將隧道空間作軸對(duì)稱處理。對(duì)于具有軸對(duì)稱的完全彈性介質(zhì)，介質(zhì)參數(shù)僅是和的函數(shù)，三維模型可以簡(jiǎn)化為二維。此時(shí)柱坐標(biāo)系下的一階彈性波速度?應(yīng)力方程為：

式中：為介質(zhì)密度；和為拉梅常數(shù)；v和v分別為軸向(方向)和徑向(方向)的速度；和分別為軸向和徑向的正應(yīng)力；為切應(yīng)力。

對(duì)上述方程采用交錯(cuò)網(wǎng)格[21]進(jìn)行離散，得到任意偶數(shù)階差分格式：

2 穩(wěn)定性條件與邊界條件

正演模擬過(guò)程中為保證離散數(shù)值計(jì)算的穩(wěn)定性需要設(shè)置穩(wěn)定性條件，與直角坐標(biāo)系不同, 柱坐標(biāo)系要求時(shí)間步長(zhǎng)更小,穩(wěn)定性條件為[22]：

式中：min為最小波長(zhǎng)；min(Δ,Δ)為空間最小網(wǎng)格步長(zhǎng)；max為最大波速。

若考慮為深埋隧道，軸對(duì)稱柱坐標(biāo)系下隧道空間全波場(chǎng)數(shù)值模擬需要處理2種邊界條件，即計(jì)算區(qū)域外的人工截?cái)噙吔绾退淼雷杂蛇吔鐥l件(圖1)。對(duì)于人工截?cái)噙吔绮捎肔IU[23]提出的柱坐標(biāo)系下的分裂式完全匹配層處理；而隧道自由邊界條件則采用橫向各項(xiàng)同性介質(zhì)替換法處理[24]，即用橫向各向同性介質(zhì)近似代替自由邊界，直接令自由邊界處的正應(yīng)力為0，切應(yīng)力為0的條件則通過(guò)對(duì)自由邊界上物性參數(shù)的設(shè)定，在波動(dòng)方程交錯(cuò)網(wǎng)格有限差分迭代求解中實(shí)現(xiàn)。二維情況下的處理方式為：

式中：τxx和τrr分別為隧道掌子面和側(cè)壁的正應(yīng)力；ρ0，λ0和μ0為自由邊界上的密度和拉梅常數(shù)。ρ，λ和μ為彈性介質(zhì)的密度和拉梅常數(shù)。

3 觀測(cè)系統(tǒng)與地質(zhì)模型

為了獲取較好的波場(chǎng)記錄，采用如圖2所示的觀測(cè)系統(tǒng)。隧道長(zhǎng)70 m，寬10 m，與實(shí)際隧道相仿。炮點(diǎn)置于隧道掌子面與右邊墻的交點(diǎn)處，埋深1 m；在隧道左右邊墻各布置了47檢波器，偏移距15 m，道間距0.5 m，埋深1.5 m。

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

為驗(yàn)證本文算法的正確性，以及分析隧道典型不良地質(zhì)現(xiàn)象的波場(chǎng)特征設(shè)計(jì)了如圖3所示的3個(gè)地質(zhì)模型。模型Ⅰ為垂直巖性分界面，界面距掌子面40 m；模型Ⅱ?yàn)閿鄬悠扑閹?，破碎帶?0 m，傾角45°；模型Ⅲ為溶洞，溶洞半徑5 m，各模型的彈性參數(shù)見(jiàn)表1。

(a) ModelⅠ；(b) ModelⅡ；(c) Model Ⅲ

表1 模型彈性參數(shù)表

4 數(shù)值模擬結(jié)果及分析

針對(duì)上述理論模型采用高階差分格式進(jìn)行數(shù)值模擬，其中模型大小為，200 m×200 m，Δ=Δ=0.25 m，Δ=0.025 ms，震源選擇中心頻率為100 Hz的雷克子波。

4.1 巖性分界面模型

如圖4為垂直界面模型徑向分量20，35和50 ms的波場(chǎng)快照。由圖4可知，在=20 ms時(shí)，由于沒(méi)有遇到巖性分界面，在波場(chǎng)快照中可以看到直達(dá)縱波(P波)、直達(dá)橫波(S波)、隧道自由表面上產(chǎn)生的面波(R波)；同時(shí)可以看到，面波延隧道側(cè)壁傳播并在掌子面與側(cè)壁交點(diǎn)處轉(zhuǎn)換為繞射橫波(RS波)[25]，繞射橫波傳至側(cè)壁時(shí)又重新轉(zhuǎn)換成面波；在=35 ms時(shí)，在巖性分界面處處產(chǎn)生的反射縱波(PP波)和透射縱波(TP波)，以及由P波產(chǎn)生的反射轉(zhuǎn)換橫波(PS波)以及透射轉(zhuǎn)換橫波(TPS波)；在=50 ms時(shí)SS波出現(xiàn)，但由橫波產(chǎn)生的轉(zhuǎn)換縱波(SP波)能量較弱，不太明顯。從數(shù)值模擬結(jié)果可以看出，波場(chǎng)傳播特征符合波的運(yùn)動(dòng)學(xué)和動(dòng)力學(xué)規(guī)律，證明了本文算法的正確性。

圖5為垂直界面模型模擬記錄。由圖5可知，徑向分量中的SS波、PS波能量強(qiáng)于軸向分量，而軸向分量中的PP波能量則強(qiáng)于徑向分量；軸向分量中炮檢同側(cè)接收的PP波能量弱于炮檢異側(cè)，而徑向分量中炮檢同側(cè)接收的PS波、SS波弱炮檢異側(cè)；界面一次反射縱波到達(dá)掌子面時(shí)發(fā)生反射，該反射波向前傳至界面處再次反射形成多次反射波(MV波)，當(dāng)其到達(dá)檢波器的時(shí)間和界面一次反射波相近時(shí)就會(huì)與界面一次反射波相互疊加，若該多次反射波能量較強(qiáng)，就會(huì)嚴(yán)重影響有效波的識(shí)別，這在軸向分量的記錄中體現(xiàn)的十分明顯；同時(shí)，由單炮記錄可估算出P波和R波的視速度分別約為1 879 m/s和1 000 m/s，而PP波和SS波的旅行時(shí)分別約為48.5 ms和83.5 ms這與模型設(shè)計(jì)相符，進(jìn)一步證明了數(shù)值模擬的正確性。

(a) t=20 ms；(b) t=35 ms；(c) t=50 ms

(a) X分量炮檢同側(cè)接收；(b) X分量炮檢異側(cè)接收；(c) R分量炮檢同側(cè)接收；(d) R分量炮檢異側(cè)接收

4.2 斷層破碎帶模型

圖6為斷層破碎帶模型模擬記錄。如圖6所示，與巖性分界面模型相比斷層破碎帶模型仍然具有軸向分量PP波能量較強(qiáng)，徑向分量PS波、SS波能量較強(qiáng)的特征；炮檢同側(cè)接收的記錄中頂、底界面的反射縱波(PP1，PP2)能量較強(qiáng)且同相軸較為清晰，可較好的區(qū)分破碎帶頂、底界面，但在一定程度上受波場(chǎng)疊加干擾，而炮檢波異側(cè)接收記錄波場(chǎng)疊加干擾嚴(yán)重，影響了有效波的識(shí)別；同時(shí)頂、底界面的反射橫波(SS1，SS2)受波場(chǎng)疊加干擾較小，利于頂、底界面的的識(shí)別與橫波成像；不同于單反射界面模型，斷層破碎帶模型在其頂?shù)捉缑婢鶗?huì)產(chǎn)生反射及透射，并且在頂?shù)捉缑嬷g還會(huì)發(fā)生多次反射，因此模擬記錄上波場(chǎng)比較復(fù)雜。

4.3 溶洞模型

圖7為溶洞模型模擬記錄。如圖7所示，對(duì)于上述溶洞作模型，溶洞頂部的反射縱波(PP1)在軸向分量及徑向分量的炮檢同側(cè)的接收記錄中有比較清晰的體現(xiàn)，但溶洞底部的反射縱波(PP2)由于受波場(chǎng)疊加影響，對(duì)溶洞底部的識(shí)別會(huì)有一定程度的干擾，徑向分量炮檢異側(cè)接收的反射縱波能量十分微弱，在記錄中幾乎無(wú)法體現(xiàn)；溶洞頂、底的反射橫波(SS1，SS2)在徑向分量的記錄中有較好的體現(xiàn)，利于橫波成像；由于溶洞內(nèi)部的低速填充物對(duì)地震波能量具有聚集作用[26]，所以溶洞底部的反射波PP2能量較溶洞頂部反射波PP1能量更強(qiáng)；不同于界面反射，溶洞作為局部地質(zhì)體，產(chǎn)生的繞射波會(huì)嚴(yán)重干擾有效波的識(shí)別。

(a) X分量炮檢同側(cè)接收；(b) X分量炮檢異側(cè)接收；(c) R分量炮檢同側(cè)接收；(d) R分量炮檢異側(cè)接收

(a) X分量炮檢同側(cè)接收；(b) X分量炮檢異側(cè)接收；(c) R分量炮檢同側(cè)接收；(d) R分量炮檢異側(cè)接收

5 結(jié)論

1) 考慮隧道空間的情況下，基于柱坐標(biāo)系下的一階彈性波速度?應(yīng)力方程，采用高階交錯(cuò)網(wǎng)格算法結(jié)合FCT技術(shù)，實(shí)現(xiàn)了柱坐標(biāo)系下隧道空間典型不良地質(zhì)模型的高精度全波場(chǎng)數(shù)值模擬，系統(tǒng)分析了其傳播機(jī)理與波場(chǎng)特征，為隧道空間復(fù)雜波場(chǎng)的認(rèn)識(shí)以及數(shù)據(jù)處理和解譯提供了有效的依據(jù)。

2) 巖性分界面模型，當(dāng)界面傾角為90°時(shí)，掌子面及地質(zhì)界面間的多次反射波對(duì)有效波的識(shí)別會(huì)產(chǎn)生干擾，在數(shù)據(jù)處理和解譯時(shí)應(yīng)當(dāng)加以區(qū)分。

3) 斷層破碎帶模型，炮檢同側(cè)接收的記錄中，頂、底界面的反射縱波能量較強(qiáng)且同相軸較為連續(xù)，可較好的區(qū)分頂、底界面，但多次波以及波場(chǎng)疊加干擾會(huì)影響有效波的識(shí)別。

4) 溶洞頂?shù)走吔绲姆瓷淇v波以及反射橫波在記錄中均有體現(xiàn)，但是溶洞底部反射縱波受波場(chǎng)疊加的影響，會(huì)對(duì)溶洞頂?shù)走吔绲淖R(shí)別造成一定程度的影響，而反射橫波能夠有效的區(qū)分溶洞頂?shù)走吔?，利于橫波成像。

[1] 張楊, 楊君, 周黎明, 等. TSP在隧道工程施工中的常見(jiàn)干擾和對(duì)巖體裂隙水及軟弱夾層等的預(yù)報(bào)研究[J]. 地球物理學(xué)進(jìn)展, 2018, 33(2): 892?899. ZHANG Yang, YANG Jun, ZHOU Liming, et al. Common interference and the prediction of rock fissure water and weak interlayer in tunnel construction using TSP[J]. Progress in Geophysics, 2018, 33(2): 892?899.

[2] 周輪, 李術(shù)才, 許振浩, 等. 隧道綜合超前地質(zhì)預(yù)報(bào)技術(shù)及其工程應(yīng)用[J]. 山東大學(xué)學(xué)報(bào)(工學(xué)版), 2017, 47(2): 55?62. ZHOU Lun, LI Shucai, XU Zhenhao, et al. Integrated advanced geological prediction technology of tunnel and its engineering application[J]. Journal of Shandong University (Engineering Science), 2017, 47(2): 55?62.

[3] 李術(shù)才, 劉斌, 孫懷鳳, 等. 隧道施工超前地質(zhì)預(yù)報(bào)研究現(xiàn)狀及發(fā)展趨勢(shì)[J]. 巖石力學(xué)與工程學(xué)報(bào), 2014, 33(6): 1090–1113. LI Shucai, LIU Bin, SUN Huaifeng, et al. State of art and trends of advanced geological prediction in tunnel construction[J]. Chinese Journal of Rock Mechanics and Engineering, 2014, 33(6): 1090?1113.

[4] 趙永貴. 國(guó)內(nèi)外隧道超前預(yù)報(bào)技術(shù)評(píng)析與推介[J]. 地球物理學(xué)進(jìn)展, 2007, 22(4): 1344?1352. ZHAO Yonggui. Analysis and recommendation of tunnel prediction techniques at home and abroad[J]. Progress in Geophysics, 2007, 22(4): 1344?1352.

[5] ZHAO Yonggui, JIANG Hui, ZHAO Xiaopeng. Tunnel seismic tomography method for geological prediction and its application[J]. Applied Geophysics, 2006, 3(2): 69? 74.

[6] Dickmann T, Sander B. Drivage concurrent tunnel seismic prediction[J]. Felsbau-Rock and Soil Engineering, 1996(14): 406?411.

[7] Inazaki T, Isahai H, Kawamura S, et al. Stepwise application of horizontal seismic profiling for tunnel prediction ahead of the face[J]. The Leading Edge, 1999, 18(12): 1429?1431.

[8] Neil D W, Haramy K Y, Hanson D H, et al. Tomography to evaluate site condition during tunneling[J]. Geotechnical Speccial Publicaion, 1999(90): 71?281.

[9] 劉江平, 程飛, 范承余, 等. 基于隧道空間全波場(chǎng)二維數(shù)值模擬與特征分析[J]. 巖土工程學(xué)報(bào), 2012, 34(9): 1705?1711. LIU Jiangping, CHENG Fei, FAN Chengyu, et al. Two- dimensional numerical simulation of tunnel-based seismic full-wave fields[J]. Chinese Journal of Geotechnical Engineering, 2012, 34(9): 1705?1711.

[10] 魯光銀, 熊瑛, 朱自強(qiáng). 隧道反射波超前探測(cè)有限差分正演模擬與偏移處理[J]. 中南大學(xué)學(xué)報(bào)(自然科學(xué)版), 2011, 42(1): 136?141. LU Guangyin, XIONG Ying, ZHU Ziqiang. 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.

[11] CHENG Fei, LIU Jiangping, QU Niannian, et al. Two-dimensional pre-stack reverse time imaging based on tunnel space[J]. Journal of Applied Geophysics, 2014 (104): 106?113. DOI:10.1016/j.jappgeo.2014.02. 013.

[12] 凌飛, 肖宏躍, 朱夏樂(lè), 等. 基于黏彈性介質(zhì)的隧道地震二維正演模擬[J]. 長(zhǎng)江科學(xué)院院報(bào), 2015, 32(5): 121?126. LING Fei, XIAO Hongyue, ZHU Xiale, et al. Two- dimensional forward seismic modeling for tunnel based on viscoelastic medium[J]. Journal of Yangtze River Scientific Research Institute, 2015, 32(5): 121?126.

[13] 宋杰. 隧道施工不良地質(zhì)三維地震波超前探測(cè)方法及其工程應(yīng)用[D]. 濟(jì)南: 山東大學(xué), 2016. SONG Jie. The three-dimensional seismic ahead prospecting method and its application for adverse geology in tunnel construction[D]. Jinan: Shandong University, 2016.

[14] 王京. 隧道空間有限元地震波場(chǎng)模擬與逆時(shí)偏移成像[D]. 武漢: 中國(guó)地質(zhì)大學(xué), 2017. WANG Jing. Finite element seismic wave field simulation and reverse time migration based on tunnel space[D]. Wuhan: China University of Geosciences, 2017.

[15] Jetschny S, Bohlen T, Kurzmann A. Seismic prediction of geological structures ahead of the tunnel using tunnel surface waves[J]. Geophysical Prospecting, 2011, 59(5): 934?946.

[16] YANG Sitong, WEI Jiuchuan, CHENG Jiulong, et al. Numerical simulations of full-wave fields and analysis of channel wave characteristics in 3-D coal mine roadway models[J]. Applied Geophysics, 2016, 13(4): 621?630.

[17] Harmankaya U, Kaslilar A, Wapenaar K, et al. Locating scatterers while drilling using seismic noise due to tunnel boring machine[J]. Journal of Applied Geophysics, 2018, 152: 86?99.

[18] LIU Yutao, LIU Jiangping, CHENG Fei, et al. Full wave field finite difference modeling of tunnel space in polar coordinates[C]// International Geophysical Conference, Qingdao, China, 17-20 April 2017, Qingdao, China. Society of Exploration Geophysicists and Chinese Petroleum Society, 2017.

[19] 查欣潔, 高星, 王偉, 等. 隧道工程勘察中的超前預(yù)報(bào)成像方法研究[J]. 地球物理學(xué)報(bào), 2018, 61(3): 1151? 1156. ZHA Xinjie, GAG Xing, WANG Wei, et al. Advanced prediction migration method research in tunnel engineering investigation[J]. Chinese J Geophys, 2018, 61(3): 1151?1156.

[20] Hiroshi Takenaka, Hiroki Tanaka, Taro Okamoto, et al. Quasi-cylindrical 2.5D wave modeling for large-scale seismic surveys[J]. Geophysical Research Letters, 2003, 30(21):2086.DOI:10.1029/2003GL018068.

[21] Graves R W. Simulating seismic wave propagation in 3D elastic media using staggered-grid finite difference[J]. Bulletin of the Seismological Society of America, 1996, 86(4): 1091?1106.

[22] LIU Qinghuo, Sinha B K. A 3D cylindrical PML/FDTD method for elastic waves in fluid-filled pressurized boreholes in triaxially stressed formations[J]. Geophysics, 2003, 68(5): 1731?1743.

[23] LIU Q H. Perfectly matched layers for elastic waves in cylindrical and spherical coordinates[J]. The Journal of the Acoustical Society of America, 1999, 105(4): 2075? 2084.

[24] Mittet R. Free-surface boundary conditions for elastic staggered-grid modeling schemes[J]. Geophysics, 2002, 67(5): 1616?1623.

[25] Bohlen T, Lorang U, Rabbel W, et al. Rayleigh-to-shear wave conversion at the tunnel face: From 3D-FD modeling to ahead-of-drill exploration[J]. Geophysics, 2007, 72(6): T67?T79.

[26] 董良國(guó), 黃超, 劉玉柱, 等. 溶洞地震反射波特征數(shù)值模擬研究[J]. 石油物探, 2010, 49(2): 121?124, 15. DONG Liangguo, HUANG Chao, LIU Yuzhu, et al. Numerical simulation of seismic wave propagation in cave carbonate reservoir[J]. Geophysical Prospecting for Petroleum, 2010, 49(2): 121?124, 15.

Full wave filed numerical simulation and analysis of tunnel space based on cylindrical coordinate system

LU Guangyin, LUO Shuai, ZHU Ziqiang, SHI Keliang, XIA Chengzhi

(School of Geosciences and Info-physic, Central South University, Changsha 410083, China)

Differing from ground seismic exploration, the special structure of tunnel space makes its wave field very complex. Full-wave field numerical simulation based on actual tunnel model is an important method to effectively understand the propagation law of wave field in tunnel space. Previous studies on this field are mostly based on Cartesian coordinate system. As a special geological model similar to air column, numerical simulation of full wave field in cylindrical coordinate system is more meaningful for understanding the characteristics of wave field in tunnel space. An arbitrary even-order staggered grid difference scheme was deduced, artificial truncated boundary and tunnel free boundary were designed, and numerical simulation and wave field characteristics analysis were performed for typical problematic geological models of tunnels such as lithologic interface, fault fracture zone and karst cave based on the first-order elastic wave velocity-stress equation in cylindrical coordinates. The results show that the staggered grid finite difference method based on the first-order elastic wave velocity-stress equation in cylindrical coordinates combined with FCT technology can achieve high-precision numerical simulation of the full-wave field in tunnel space. The wave field characteristics are in accordance with the kinematic and dynamic characteristics of the wave, which provides an effective basis for understanding the complex wave field characteristics in tunnel space and interpreting the data.

advance prediction; tunnel space; cylindrical coordinates; full wave field; numerical simulation

P631.4

A

1672 ? 7029(2020)02 ? 0388 ? 08

10.19713/j.cnki.43?1423/u.T20190267

2019?04?08

國(guó)家自然科學(xué)基金資助項(xiàng)目(41974148)；湖南省安全生產(chǎn)監(jiān)督管理局資助項(xiàng)目(201907)

魯光銀(1976?)，男，湖北宜昌人，教授，博士，從事工程地球物理勘探研究；E?mail：13975894898@139.com

(編輯 蔣學(xué)東)