董 毅，李靜蕊，喬世海，仝景陽

（1.中煤航測遙感集團(tuán)有限公司，陜西 西安 710100；2.渭南陜煤啟辰科技有限公司，陜西 渭南 714026）

1 引言

早期煤礦開采留下的采空區(qū)形成了大型儲(chǔ)水空間，一旦形成裂隙，導(dǎo)通含水層，將直接引發(fā)突水事故［1］。目前煤礦采空區(qū)探測常用的探測方法包括鉆探法和地球物理方法［2］。鉆探法直接準(zhǔn)確，但僅為“一孔之見”，無法識(shí)別采空區(qū)空間分布及延伸方向。地球物理方法中如高密度直流電法、瞬變電磁法［3］、等值反磁通瞬變電磁法［4］、地空瞬變電磁法［5］和可控源音頻大地電磁法［6］，是通過分析地下介質(zhì)電性特征來推斷采空區(qū)分布。然而，高密度直流電法探測深度較淺，為接觸式探測方法。瞬變電磁法一般使用非接觸式探測，即利用不接地回線向地下發(fā)射一次性場，關(guān)斷后由接收回線采集地下介質(zhì)受激發(fā)產(chǎn)生的渦流場。采空區(qū)的存在會(huì)影響地下電磁場的傳播特征，具體為采空區(qū)積水時(shí)為相對(duì)低阻，電磁場向下傳播較慢；空區(qū)不積水時(shí)為相對(duì)高阻，電磁向下傳播較快。瞬變電磁法是通過分析二次渦流場的衰減規(guī)律，從而推斷地下介質(zhì)分布特征［7］。該方法對(duì)低阻體靈敏，探測精度高，野外工作便捷，在煤礦采空區(qū)探測中應(yīng)用廣泛［8，9］。瞬變電磁法對(duì)低阻體探測靈敏的特性決定了其對(duì)含水體探測的優(yōu)勢(shì)。然而，當(dāng)?shù)叵聹\層存在低阻層時(shí)，會(huì)延緩電磁場向下傳播的速度，相當(dāng)于形成淺層“屏蔽層”，嚴(yán)重影響對(duì)下方含水體的探測精度。2005年，石顯新等［10］證明了低阻層對(duì)瞬變電磁法探測的屏蔽作用，指出探測同樣深度需要更長的探測時(shí)間，從而導(dǎo)致晚期信噪比降低，不利于對(duì)下方低阻體的探測精度。2009年，石顯新等［11］分析了華北煤田地層屬于A型地層，在進(jìn)行瞬變電磁法探測時(shí)需要充分考慮低阻屏蔽層問題，提出需要增加探測時(shí)間和發(fā)射功率，保證足夠的探測深度和信噪比強(qiáng)度。侯彥威等［12，13］在大回線瞬變電磁法迭代過程中，分別利用控制步長的約束反演方法和擬大地電磁法（Magnetotelluric Method，MT）反演技術(shù)，實(shí)現(xiàn)了在低阻覆蓋層下淺埋積水采空區(qū)的精細(xì)探測。程輝等［14］利用數(shù)值模擬研究了雙層采空區(qū)組合下的瞬變電磁場響應(yīng)特征，指出下層采空區(qū)電磁響應(yīng)容易受上層低阻影響，瞬變電磁法只能分辨一定層間距的高阻—低阻模型。Yan Shu等［15］分析了多層采空區(qū)的電磁響應(yīng)特征，認(rèn)為當(dāng)保證足夠的觀測時(shí)間和下方采空區(qū)異常響應(yīng)幅度大于觀測誤差條件時(shí)，利用中心回線瞬變電磁法探測多層充水采空區(qū)是可行的。Dong Yi等［16］采用地震與可控源音頻大地電磁聯(lián)合反演方法，通過數(shù)值仿真和工程應(yīng)用，驗(yàn)證了該方法能夠準(zhǔn)確探測砂巖低阻層下方煤層頂板巖層的富水性。Wang Peng等［17，18］分析了大定源回線瞬變電磁法探測不同尺度采空區(qū)時(shí)的電磁響應(yīng)特征，為瞬變電磁探測采空區(qū)資料解釋提供了可靠依據(jù)。

以上研究表明，淺部低阻層對(duì)下方采空區(qū)存在電磁屏蔽作用，但沒有深入研究電磁場在三維采空區(qū)模型空間中的真實(shí)傳播特征，且未分析不同參數(shù)對(duì)低阻層下方采空區(qū)瞬變電磁法探測的影響規(guī)律。本次研究以采空區(qū)瞬變電磁法探測為例，建立低阻層覆蓋下采空區(qū)三維地質(zhì)模型，利用數(shù)值模擬技術(shù)，分析低阻覆蓋層下采空區(qū)瞬變電磁場響應(yīng)特征，揭示采空區(qū)富水程度、埋深以及低阻層厚度對(duì)電磁場傳播的影響規(guī)律，指導(dǎo)該條件下采空區(qū)瞬變電磁法探測資料的合理解釋及野外數(shù)據(jù)采集時(shí)的參數(shù)設(shè)置。

2 瞬變電磁法三維有限差分正演理論

本次研究采用李飛等提出的基于全空間初始場源的半空間瞬變電磁法三維有限差分正演方法［19］，其正演精度已被驗(yàn)證準(zhǔn)確。均勻全空間條件下垂直磁偶極子源在直角坐標(biāo)系下展開公式為：

式中，H為磁場強(qiáng)度分量，A／m；E為電場強(qiáng)度，V／m；（x，y，z）為空間位置，r是收發(fā)距，r＝，單位m；φ（u）＝為概率積分，u＝，ρ是均勻全空間介質(zhì)電阻率，單位Ω·m。M是發(fā)射磁矩，M＝IS；其中I為發(fā)射電流，單位A；S為發(fā)射面積，單位m2。

初始時(shí)刻和時(shí)間迭代步長的選擇參考Wang等提出的方法［20］，均勻半空間模型上半空間設(shè)為空氣層，避免網(wǎng)格向上延拓帶來的計(jì)算誤差。該方法另一優(yōu)勢(shì)為最大和最小網(wǎng)格比例不受限制，很適合第一類邊界條件，具有更快的計(jì)算速度。

3 采空區(qū)三維電磁響應(yīng)特征分析

選擇地面采集的垂直磁場分量轉(zhuǎn)換的感應(yīng)電動(dòng)勢(shì)參數(shù)為研究對(duì)象［21，22］，利用數(shù)值模擬方法，分析不同三維地質(zhì)模型在相同時(shí)刻的地下空間電磁場分布特征，總結(jié)電磁場傳播規(guī)律。以陜北侏羅系含煤地層為例，對(duì)煤系地層進(jìn)行簡化，建立的三維地質(zhì)—地球物理模型如圖1所示。地層由下至上分為延安組煤系地層、直羅組砂巖層和上層覆蓋層。其中煤系地層中包含采空區(qū)，各地層及采空區(qū)參數(shù)如表1所示。圖1中采空區(qū)（包含裂隙區(qū)）體積為100m×100m×50m（長×寬×高），頂部距地面200m。其充水時(shí)的電阻率為10Ω·m，未充水時(shí)電阻率為1000Ω·m。數(shù)值模擬時(shí)，發(fā)射線框平鋪于地面，邊長600m，在線框中心1／3范圍內(nèi)共布設(shè)13個(gè)接收點(diǎn)，點(diǎn)距10m。發(fā)射電流10A，接收面積1m2。

表1 地層劃分及電阻率賦值Table 1 Strata division and resistivity assignment

圖1 采空區(qū)三維地質(zhì)—地球物理模型Fig.1 Three－dimensional geological－geophysical model of goaf

3.1 采空區(qū)不同充水程度

建立采空區(qū)不同充水程度三維地質(zhì)—地球物理模型，沿著地面接收點(diǎn)繪制二維地質(zhì)截面，結(jié)果如圖2所示，其中圖2（a）為不充水采空區(qū)，圖2（b）為半充水采空區(qū)，圖2（c）為全充水采空區(qū)。圖2中地層和采空區(qū)的厚度及電性參數(shù)保持不變，僅改變采空區(qū)的充水程度。當(dāng)采空區(qū)半充水時(shí)，上部采空區(qū)電阻率1000Ω·m，下部采空區(qū)電阻率10Ω·m，厚度均為25m。

圖2 不同充水程度采空區(qū)地質(zhì)—地球物理模型二維斷面Fig.2 Two－dimensional section diagram of geologicalgeophysical model of goaf with different content of water filling

圖3 為圖2在2ms時(shí)在接收點(diǎn)范圍內(nèi)各地質(zhì)—地球物理斷面對(duì)應(yīng)的地下空間電磁場分布特征，以200m×500m的xoz電磁場剖面進(jìn)行成圖分析。圖中藍(lán)色虛線為地質(zhì)斷面中砂巖低阻層位置，紅色虛線為采空區(qū)位置（下同）。由圖3可以看出，采空區(qū)不同充水程度時(shí)的地下空間電磁場分布特征差別較大。地層由上到下，感應(yīng)電動(dòng)勢(shì)幅值整體為“低—高—低”分布特征，其中低阻層范圍內(nèi)的感應(yīng)電動(dòng)勢(shì)幅度最高，煤系地層感應(yīng)電動(dòng)勢(shì)幅值最低。分析是由于砂巖層電阻率最低，電磁場衰減較慢，相同時(shí)刻的感應(yīng)電動(dòng)勢(shì)幅值相對(duì)最大。煤系地層的電阻率最高，感應(yīng)電動(dòng)勢(shì)衰減速度最快，相同時(shí)刻的感應(yīng)電動(dòng)勢(shì)幅值相對(duì)最低，證明電磁場響應(yīng)特征與地質(zhì)模型的電性特征吻合。此外，當(dāng)采空區(qū)不充水及半充水時(shí)，地下空間電磁場分布基本一致，分析為采空區(qū)富水體積較小，電磁場傳播對(duì)其敏感性低導(dǎo)致的。當(dāng)采空區(qū)全充水時(shí)，在低阻電性層下方的感應(yīng)電動(dòng)勢(shì)等值線呈明顯的下凹特征，且采空區(qū)左右兩側(cè)感應(yīng)電動(dòng)勢(shì)等值線彎曲變形明顯，據(jù)此可識(shí)別采空區(qū)的底及左右界面。然而，采空區(qū)內(nèi)部感應(yīng)電動(dòng)勢(shì)幅值略低于上方含水層，分析為受上方砂巖含水層影響，其對(duì)下方含水采空區(qū)有一定屏蔽作用，并且砂巖含水層和全充水采空區(qū)頂部的電磁場響應(yīng)等值線相融，導(dǎo)致無法識(shí)別積水采空區(qū)的頂部埋深位置，不利于資料的準(zhǔn)確解釋。

圖3 不同充水程度采空區(qū)地下空間電磁場分布Fig.3 Electromagnetic field distribution map of underground pace in goaf with different content of water filling

3.2 不同低阻層厚度

建立不同低阻層厚度的三維地質(zhì)—地球物理模型，沿著地面接收點(diǎn)繪制的二維地質(zhì)斷面如圖4所示，其中圖4（a）中低阻層厚度為10m，圖4（b）中低阻層厚度為100m。數(shù)值模擬時(shí)僅改變上方低阻層的厚度，其余地層及采空區(qū)參數(shù)不變，并且低阻層底界面埋深保持150m位置。

圖4 不同低阻層厚度的地質(zhì)—地球物理模型二維斷面Fig.4 Two－dimensional cross section diagram of geological－geophysical model with different thickness of low resistivity layer

圖5 為圖4在2ms時(shí)在接收點(diǎn)范圍內(nèi)各地質(zhì)—地球物理斷面對(duì)應(yīng)的地下空間電磁場分布。聯(lián)合圖3（c）進(jìn)行對(duì)比可以看出，不同低阻層厚度時(shí)的地下空間電磁場分布特征差別較大。當(dāng)?shù)妥鑼右来巫兒駮r(shí)，在深度小于250m范圍內(nèi)的電磁場響應(yīng)幅值依次增大，分析是由于低阻層的存在，使相同時(shí)刻的電磁場衰減速度減慢導(dǎo)致的。以電磁場等值線彎曲變形作為其傳播深度的識(shí)別標(biāo)志，在2ms時(shí)，圖5（a）中電磁場已傳播至350m深度，圖3（c）中電磁場傳播至300m深度，圖5（c）中電磁場傳播至約270m深度。這表明淺部低阻體的存在將會(huì)延緩電磁場向下擴(kuò)散的速度。若使電磁場傳播至相同深度，當(dāng)上方低阻層越厚時(shí)，則需要更長的觀測時(shí)間。此外，當(dāng)?shù)妥鑼虞^薄時(shí)，下方采空區(qū)的電磁場響應(yīng)與上方低阻層能夠完全分離，其向下傳播速度也較快，有助于準(zhǔn)確識(shí)別采空區(qū)的頂部位置。但隨著低阻層厚度增大，受上方低阻層屏蔽影響，積水采空區(qū)頂部的電磁場等值線與低阻層逐漸難以分離，無法準(zhǔn)確識(shí)別采空區(qū)的頂部埋深，僅能據(jù)此圈定積水采空區(qū)的底及左右邊界位置。

圖5 不同低阻層厚度地下空間電磁場分布Fig.5 Electromagnetic field distribution of underground space with different thickness of low resistance layer

3.3 不同采空區(qū)埋深

建立不同采空區(qū)埋深的三維地質(zhì)—地球物理模型，沿著地面接收點(diǎn)繪制的二維地質(zhì)斷面如圖6所示，其中圖6（a）中采空區(qū)頂界面埋深170m，距上方低阻層底界面20m，圖6（b）中采空區(qū)頂界面埋深250m，距上方低阻層底界面100m。圖中僅改變積水采空區(qū)頂界面埋深，其余地層及采空區(qū)參數(shù)不變。

圖6 不同采空區(qū)埋深地質(zhì)—地球物理模型二維斷面Fig.6 Two－dimensional section diagram of geological－geophysical model with different buried depth of goaf

圖7為圖6在2ms時(shí)在接收點(diǎn)范圍內(nèi)各地質(zhì)—地球物理斷面對(duì)應(yīng)的地下空間電磁場分布。

圖7 不同采空區(qū)埋深的地下空間電磁場分布Fig.7 Electromagnetic field distribution map of underground space with different buried depth of goaf

聯(lián)合圖3（c）進(jìn)行對(duì)比可以看出，隨著積水采空區(qū)頂部埋深逐漸增大，相同時(shí)刻的電磁場向下傳播也越深，積水采空區(qū)處均表現(xiàn)為電磁場等值線出現(xiàn)明顯下凹特征。當(dāng)積水采空區(qū)頂部埋深距低阻層底界面20m時(shí)，采空區(qū)電磁場與低阻層未能分離，難以識(shí)別采空區(qū)頂界面位置；當(dāng)積水采空區(qū)頂部埋深距低阻層底界面100m時(shí)，積水采空區(qū)范圍內(nèi)的電磁場與低阻層已經(jīng)能夠分離，根據(jù)感應(yīng)電動(dòng)勢(shì)等值線下凹彎曲特征，能夠圈定采空區(qū)的空間位置；隨著采空區(qū)埋深減小，其范圍內(nèi)的感應(yīng)電動(dòng)勢(shì)幅度依次增大，分析為相同時(shí)刻砂巖低阻層使電磁場衰減較慢導(dǎo)致的，不利于對(duì)采空區(qū)位置的準(zhǔn)確圈定。

3.4 感應(yīng)電動(dòng)勢(shì)曲線對(duì)比

為進(jìn)一步分析不同采空區(qū)三維地質(zhì)模型的電磁場響應(yīng)特征，取不同地質(zhì)模型在采空區(qū)正上方7號(hào)測點(diǎn)的感應(yīng)電動(dòng)勢(shì)衰減曲線進(jìn)行對(duì)比，總結(jié)其隨時(shí)間的衰減規(guī)律。感應(yīng)電動(dòng)勢(shì)衰減曲線如圖8所示，圖8（a）對(duì)應(yīng)不同富水程度采空區(qū)地質(zhì)模型，圖8（b）對(duì)應(yīng)不同低阻層厚度地質(zhì)模型，圖8（c）對(duì)應(yīng)不同采空區(qū)埋深地質(zhì)模型。

可以看出，圖8（a）中采空區(qū)不含水和半含水地質(zhì)模型在7號(hào)點(diǎn)的感應(yīng)電動(dòng)勢(shì)衰減曲線基本重合，而相同時(shí)刻全充水采空區(qū)在7號(hào)點(diǎn)的感應(yīng)電動(dòng)勢(shì)幅值略高，分析是由于低阻采空區(qū)延遲電磁場衰減速度導(dǎo)致的，結(jié)果與圖3中吻合。圖8（b）中各曲線間差別較大，在0.3ms之前，低阻層厚度10m時(shí)的感應(yīng)電動(dòng)勢(shì)幅值最大，其次是低阻層厚度50m，低阻層厚度100m時(shí)感應(yīng)電動(dòng)勢(shì)值最低。0.3ms之后，低阻層厚度10m時(shí)感應(yīng)電動(dòng)勢(shì)衰減速度最快，其次是低阻層厚度50m，而低阻層厚度100m時(shí)7號(hào)點(diǎn)的感應(yīng)電動(dòng)勢(shì)衰減速度最慢，幅值相對(duì)最高。分析為0.3ms之前受地表覆蓋層厚度的影響較大，0.3ms之后受砂巖含水層厚度影響較大導(dǎo)致的。圖8（c）中7號(hào)測點(diǎn)感應(yīng)電動(dòng)勢(shì)曲線的衰減規(guī)律不同，隨著采空區(qū)頂界面埋深增加，相同時(shí)刻感應(yīng)電動(dòng)勢(shì)衰減曲線越慢，幅值越大；隨著采空區(qū)頂界面埋深減小，感應(yīng)電動(dòng)勢(shì)衰減曲線越快，幅值越小，這是受砂巖層與采空區(qū)之間煤系地層的影響。由此可以得出，采空區(qū)富水性、低阻層厚度和采空區(qū)埋深均是影響電磁法探測低阻層下方采空區(qū)分辨率的重要因素。

4 結(jié)論

利用數(shù)值模擬方法，對(duì)比分析了低阻層覆蓋下不同采空區(qū)地質(zhì)模型的三維瞬變電磁場響應(yīng)特征，得出以下結(jié)論：

1）淺部低阻層等效電磁“屏蔽層”會(huì)延緩電磁場向地下傳播的速度，若想獲得較深處地層的電性響應(yīng)，需要延長地面數(shù)據(jù)采集時(shí)間。低阻層越厚，地面接收的電磁場響應(yīng)幅度越大，然而，積水采空區(qū)引起的電磁場響應(yīng)在頂部位置越難與其分離。

2）低阻層下積水采空區(qū)會(huì)引起電磁場等值線出現(xiàn)下凹特征，采空區(qū)左右兩側(cè)電磁場等值線出現(xiàn)明顯彎曲變形。隨著采空區(qū)含水程度減弱、埋深增大，采空區(qū)引起的電磁場響應(yīng)幅度下降明顯，加上淺部低阻層的存在，不利于電磁法對(duì)下方采空區(qū)的精細(xì)探測。

3）數(shù)值模擬分析結(jié)果可為相同地質(zhì)情況下，瞬變電磁法準(zhǔn)確探測采空區(qū)時(shí)的參數(shù)設(shè)置及資料合理解釋提供依據(jù)。