萬豪豪，朱 魯，陳金昊，李文喜

（山東科技大學(xué) 地球科學(xué)與工程學(xué)院，山東 青島 266590）

直流高密度電阻率法在礦井中的應(yīng)用已漸為成熟[1-2]，但隨著煤礦開采深度的逐漸增加，地質(zhì)條件變得尤為復(fù)雜，其中巖體突水已成為危害礦井安全的重大災(zāi)害之一，為確保井下開采工作的安全高效進行，科研工作者進行了深入研究。2003年，底青云分別運用佐迪法和2D積分法進行了高密度電阻率成像方面的研究[3]。2005年，為解決工作面采空區(qū)滯后突水問題，翟培合開發(fā)了底板突水動態(tài)監(jiān)測系統(tǒng)，實現(xiàn)了工作面底板突水監(jiān)測過程[4]；2008年，施龍青等運用三維高密度電法探測了煤礦巷道底板富水性[5]；劉志新、徐志剛和王大慶等[6-7]介紹了電阻率成像技術(shù)的特殊問題并對電阻率成像技術(shù)進行了展望；2016年，高衛(wèi)富等使用ANSYS進行全空間直流電法異常體正演模擬，精準地探測出了異常體的位置[8]。

1 地電場三維有限元理論

有限元理論的數(shù)學(xué)基礎(chǔ)是變分原理和剖分插值，將研究區(qū)剖分成有限個小單元，在每個單元上，將函數(shù)簡化成線性函數(shù)、二次函數(shù)或高次函數(shù)。它適用于物理性質(zhì)及參數(shù)復(fù)雜區(qū)域，其求解過程規(guī)范化，適于解決復(fù)雜數(shù)學(xué)模型問題[9]。

在直流電法數(shù)值模擬時，由于點電源周圍電場會因存在諸多干擾因素使得無法精確計算，采用二次場方法解決數(shù)據(jù)的誤差問題，使用正方形的網(wǎng)格對其進行劃分，每個網(wǎng)格節(jié)點的潛在值被設(shè)置為1個未知數(shù)，未知數(shù)被代入拉普拉斯方程進行計算，即將模型的電導(dǎo)率σ細分為背景電導(dǎo)率σp和異常電導(dǎo)率σs，其分別產(chǎn)生一次場up和二次場us，則總場u為二者疊加：

三維電場電位的邊值問題與下列變分問題等價：

式中：σ為地下介質(zhì)的電導(dǎo)率；r為點電源到邊界點的距離；Г∞為區(qū)域Ω的下邊界；Г為2種介質(zhì)的分界面；n為邊界外法線方向的向量坐標(biāo)；u0為正常電位。

2 模型建立與完善

在礦井中進行三維直流電法勘探時，根據(jù)礦井環(huán)境的特點，井下巷道之間大多數(shù)都是通過聯(lián)絡(luò)巷連接，但是每條聯(lián)絡(luò)巷并不都是在1條直線上，無法沿用地面勘測的傳統(tǒng)布線方法，且巷道與聯(lián)絡(luò)巷兩幫的圍巖均可導(dǎo)電，使得頂板與底板之間有著電力聯(lián)系，故并沒有采用傳統(tǒng)模擬實驗中，將巷道、聯(lián)絡(luò)巷設(shè)定為類空氣的超高阻體，而采用井下全空間地電模型[10]，將巷道與聯(lián)絡(luò)巷相結(jié)合形成測量通道。

巷道設(shè)置于地表深度600 m處，巷道長度設(shè)置為300 m，位于x坐標(biāo)-100～200 m處，聯(lián)絡(luò)巷長度設(shè)置為100 m，位于y坐標(biāo)-50～50 m處，總長為600 m，模型所要測量的每條邊長250 m，寬100 m，無窮遠邊界長2 400 m，寬2 400 m，深1 500 m，地質(zhì)模型位置平面圖如圖1。正方形GHKM表示無窮遠邊界，ABCD表示模型所要測量的范圍，ABEF表示巷道和聯(lián)絡(luò)巷的范圍。

對所建立的模型進行網(wǎng)格剖分[11]，沿著x軸剖分的平面數(shù)110個，沿著y軸剖分的平面數(shù)為110個，沿著z軸剖分的平面數(shù)為100個，三維立體空間總共剖分區(qū)塊數(shù)約為1.2×106個。為了提高正演計算的精度，在網(wǎng)格剖分過程中對所研究的目標(biāo)區(qū)域進行加密剖分，巷道，聯(lián)絡(luò)巷網(wǎng)格單元剖分大小為5 m，研究區(qū)域到無窮遠邊界處網(wǎng)格剖分逐漸變大，網(wǎng)格剖分大小不均勻。

圖1 地質(zhì)模型位置平面圖

供電電極A與測量電極M分布在ABCD的巷道上，電極之間的距離選擇10 m，電極總數(shù)為61個，電極編排的順序為1～61，電極測量的順序為CBAD，采用固定斷面連續(xù)測量掃描方式，測量深度為A到M距離的一半，測量點O位于A到M距離的中點，二級裝置AM[12]在三維空間中測量時，測量點O的位置如圖2。二極裝置A-M三維測量散點圖如圖2。

圖2 二極裝置A-M三維測量散點圖

根據(jù)實際礦井突水隱患建立在層狀介質(zhì)中的地電模型，巷道和聯(lián)絡(luò)巷的電阻率為10 000 Ω·m，圍巖的電阻率為2 000 Ω·m。模型1設(shè)置塊狀低阻異常體尺寸為20 m×20 m×20 m，電阻率值為10 Ω·m，位于2條巷道之間的煤層下方，頂部距離煤層底板20 m；模型2設(shè)置低阻陷落柱異常體模型尺寸為80 m×80 m×20 m，電阻率值為 10 Ω·m，位于 2條巷道之間的煤層下方，頂部距離煤層底板30 m；模型3在模型2的基礎(chǔ)上改變陷落柱在底板的位置，其余設(shè)置同模型1；模型4設(shè)置正方體低阻異常體尺寸為 40 m×40 m×40 m，電阻率值為 1 Ω·m，位于 2條巷道之間的煤層下方，頂部距離煤層底板20 m。

3 數(shù)值模擬結(jié)果與分析

3.1 塊狀低阻異常體地電模型

塊狀低阻異常體地電模型經(jīng)過有限元數(shù)值模擬計算后在經(jīng)過反演得到xz切面圖（圖3）。相比于正常地層明顯出現(xiàn)了低阻異常區(qū)，視電阻率值在10～30 Ω·m，異常體的位置和大小已經(jīng)在圖中用正方形線圈標(biāo)出，對比地電模型中異常體的位置發(fā)生了偏移，大小形態(tài)發(fā)生了改變，但總體來說二極裝置（A M）對于本次的塊狀低阻異常體的捕捉能力較強，分辨低阻異常體的能力較高。

圖3 塊狀低阻異常體地電模型xz切面圖（y=0）

塊狀低阻異常體地電模型經(jīng)模擬計算后的三維切片圖如圖4。其中y=50 m和y=-50 m的xz切片在位于巷道底板處，有2條相互平行的條帶，且條帶所顯示的視電阻率較高，推測是由于巷道高阻區(qū)的影響所致，無異常體區(qū)域模擬結(jié)果良好，巷道底板下方綠色區(qū)域的視電阻率值在100 Ω·m左右，與地電模型相符。但異常體x的坐標(biāo)范圍與實際存在偏差，因此，二極裝置（AM）在對異常體測量過程中，在y方向上測量效果好，分辨率高。

圖4 塊狀低阻異常體地電模型xz方向三維切片圖

3.2 陷落柱低阻異常體地電模型

通過陷落柱低阻異常體地電模型xz切面圖（圖5），二極裝置對本次分辨陷落柱低阻異常體的能力較強，斷面圖可以明顯地反映低阻異常區(qū)的大小以及埋深，而且也能夠判定異常體的形態(tài)，異常體的位置和大小已經(jīng)在圖中用長方形線圈標(biāo)出，視電阻率值在10～30 Ω·m，其頂?shù)茁裆钆c實際模型基本吻合。

圖5 陷落柱低阻異常體地電模型xz切面圖（y=0）

由陷落柱低阻異常體地電模型xz方向三維切片圖（圖6）可知，異常體視電阻率值在100 Ω·m左右。從低阻區(qū)的視電阻率范圍可推測異常體的y坐標(biāo)范圍為0～15 m，與實際地電模型中y坐標(biāo)的范圍0～20 m基本一致，異常體x的坐標(biāo)范圍與實際模型一致。因此，二極裝置（A-M）在對異常體測量過程中，在y方向上測量效果好，分辨率高。中間切片的地層，除了相對低阻區(qū)外，與正常地層和塊狀異常體的視電阻率基本相同。

圖6 陷落柱低阻異常體地電模型xz方向三維切片圖

3.3 柱狀低阻異常體地電模型

柱狀低阻異常區(qū)的大小和埋深與塊狀低阻異常體地電模型保持一致，異常體的長度基本和地電模型中的一樣，其頂?shù)茁裆钆c實際模型基本吻合。柱狀低阻異常體地電模型xz切面圖如圖7。但低阻體的左邊有明顯的上移，視電阻率值在20～50 Ω·m，其相較于陷落柱低阻異常體視電阻率值偏大，出現(xiàn)這種現(xiàn)象的原因有可能是異常體的中心坐標(biāo)發(fā)生偏移，離巷道較近，受到巷道高阻體的影響。

從柱狀低阻異常體地電模型xz方向三維切片圖（圖8）觀察可知，低阻區(qū)的視電阻率范圍可推測異常體的y坐標(biāo)范圍為-25～15 m，與實際地電模型中y坐標(biāo)的范圍0～20 m基本一致。

圖7 柱狀低阻異常體地電模型xz切面圖（y=0）

圖8 柱狀低阻異常體地電模型xz方向三維切片圖

3.4 正方體低阻異常體地電模型

正方體低阻異常體地電模型xz切面圖如圖9。正方體低阻體模擬結(jié)果可以反映低阻異常區(qū)的大小和埋深，視電阻率值在1～30 Ω·m，圖中的空白區(qū)域視電阻率值小于10 Ω·m。受低阻異常體的影響，主要異常區(qū)域集中在異常體的上部，地面下660 m處是異常區(qū)域的分界面，但660 m下方異常反應(yīng)在實際成像中并沒有出現(xiàn)，主要是受探測深度的影響，探測的越深，數(shù)據(jù)量越少，分辨率越低，因此，在700 m以下的位置視電阻率值并不準確。

圖9 正方體低阻異常體地電模型xz切面圖（y=0）

從正方體低阻異常體地電模型xz方向三維切片圖（圖10）可知，從低阻區(qū)的視電阻率范圍可推測異常體的y坐標(biāo)范圍為-25～25 m，與實際地電模型中y坐標(biāo)的范圍-20～20 m基本一致，在y=0方向切片圖中低阻異常區(qū)的范圍大于實際地電模型，這可能是由于低阻異常體的背景電阻率值非常小，電阻率值為1 Ω·m所致。

圖10 正方體低阻異常體地電模型xz方向三維切片圖

4 工程應(yīng)用

某礦1302上工作面（3煤）位于山西組中下部，3煤層開采過程中深受頂板砂巖水害的威脅，為保證1302上工作面開采的順利、安全進行，對該礦的1302上工作面進行井下高密度三維電法勘探。要求頂板最大探測高度為頂板上150 m，構(gòu)建頂板砂巖地層電阻率三維數(shù)據(jù)體，應(yīng)用切片技術(shù)圈定頂板砂巖地層的富水區(qū)域，結(jié)合具體鉆探情況指導(dǎo)1302上工作面開采過程中的防治水工作。

本次井下高密度三維電法勘探數(shù)據(jù)采集使用的是WDJD-4型高密度電阻率系統(tǒng)，采用A-M二極裝置形式的三維電法采集方法，采集參數(shù)中儀器所使用的電極道數(shù)是60道，測量的剖面數(shù)為30個，采用滾動的測量方式，2個相鄰電極之間的道間距為10 m，供電電壓選擇180 V的直流電壓。

把RES3DINV反演軟件所得出的三維立體成像圖數(shù)據(jù)輸出，然后在用Slicer Dicer和Voxler對三維立體成像圖進行多方向的三維切片，即為最終的成像成果圖。1302上工作面頂板砂巖地層電阻率三維數(shù)據(jù)體原點位于運輸巷與1聯(lián)絡(luò)巷向交匯處，沿運輸巷向里為y、沿1聯(lián)絡(luò)巷向下為正方向x，高度方向為垂直地層向上，本次勘探建立的1302上工作面外側(cè)部分頂板砂巖地層電阻率三維數(shù)據(jù)體長度750 m、寬度140 m、高度280 m。

1302上工作面頂板地層電阻率三維數(shù)據(jù)體如圖11，1302上工作面頂板地層電阻率三維數(shù)據(jù)體切片技術(shù)如圖12，1302上工作面頂板地層電阻率三維數(shù)據(jù)體垂層切片如圖13。

根據(jù)數(shù)據(jù)整體（圖11）和數(shù)據(jù)切片（圖12）可知，低阻區(qū)域明顯，對比圖13，可圈定3個低阻區(qū)進行解釋。第1區(qū)域位于橫向起始測點到280 m，頂板向上至高度180 m；第2區(qū)域位于橫向320～400 m之間，頂板向上至高度90 m處；第3區(qū)域位于橫向400～750 m之間，頂板向上至高度120 m處，以上3個區(qū)域，低阻情況明顯，推測為含水區(qū)域。

圖11 1302上工作面頂板地層電阻率三維數(shù)據(jù)體

圖12 1302上工作面頂板地層電阻率三維數(shù)據(jù)體切片技術(shù)

圖13 1302上工作面頂板地層電阻率三維數(shù)據(jù)體垂層切片（平行切眼、50 m間距）

5 結(jié)論

1）使用Visual Fortran編程軟件，對一定深度范圍內(nèi)所建立的三維全空間地電模型進行正演數(shù)值模擬，反演后的結(jié)果能夠準確反應(yīng)出所建地電模型異常體區(qū)域。Surfer軟件繪成的ρs視電阻率斷面圖，在一定的深度范圍內(nèi)，能大體上反映異常體的形態(tài)、位置和大小，分辨率高，基本上與地電模型吻合。

2）三維直流電法在煤礦中的應(yīng)用實例，圈定了3煤層頂板砂巖地層富水區(qū)，后期經(jīng)過鉆探驗證了富水區(qū)的存在，針對性強，極大地提高了鉆探的效率，從生產(chǎn)實際角度論證了高密度電法的可行性，為指導(dǎo)礦井水防治補充了新的方法。