地面瞬變電磁勘探中源外地形的影響

饒麗婷，武 欣，黨 博，黨瑞榮，李 勇

(1. 西安石油大學(xué) 電子工程學(xué)院,陜西 西安 710065; 2. 中國科學(xué)院地質(zhì)與地球物理研究所中國科學(xué)院礦產(chǎn)資源研究重點(diǎn)實(shí)驗(yàn)室，北京 100029; 3. 北京中科地垣科技有限公司，北京 100190)

0 引 言

瞬變電磁法是近年來地球物理領(lǐng)域發(fā)展較快的一種重要方法，具有探測深度大、靈敏度高、方便快捷、成本低的特點(diǎn)，廣泛應(yīng)用于金屬礦、石油、天然氣、煤炭、地?zé)?、地下水等探測中。邱志忠等研究表明，地形是影響電磁場分布的主要因素之一。中國山地、高原、丘陵眾多，這些起伏地形的國土面積約占總國土面積的2/3，復(fù)雜地形為中國礦產(chǎn)資源的發(fā)現(xiàn)和勘查造成了極大困難。目前，地面瞬變電磁勘探數(shù)據(jù)的處理與解釋通?；谄教沟匦?；在復(fù)雜地形區(qū)域?qū)嵤┛碧焦ぷ鲿r(shí)，若直接基于平坦地形進(jìn)行數(shù)據(jù)處理與解釋，容易將地形影響錯(cuò)誤地解釋為地下異常體，并掩蓋地下真實(shí)目標(biāo)體的響應(yīng)，造成數(shù)據(jù)解釋的較大誤差。近年來，隨著國內(nèi)地形地質(zhì)條件良好的礦產(chǎn)資源已逐步勘探完畢，為了進(jìn)一步加大礦產(chǎn)資源的勘查力度，前人逐漸開展了考慮地形條件的電磁勘探法相關(guān)研究。2012年邱衛(wèi)忠通過測道圖判斷不同尺寸回線源的地形影響特征，并針對不同情況給出了簡單的地形校正方法；殷長春等分別在2015～2017年先后開展了起伏地形上航空電磁法的數(shù)值模擬研究，并探討了山谷、山脊等簡單地形條件對航空電磁場的影響規(guī)律；2016年薛國強(qiáng)等通過對地形影響基本規(guī)律、經(jīng)典比值校正原理的分析，提出了一種高效且實(shí)用的地形影響校正算法；2019年曹鳳鳳基于有限元法計(jì)算了起伏地形條件下地空瞬變電磁響應(yīng)，分析了地空瞬變電磁系統(tǒng)中地形效應(yīng)在斜坡、山谷、山脊等典型地形條件下的影響特征；2021年齊彥福等通過結(jié)合非結(jié)構(gòu)時(shí)間域有限元法和L-BFGS反演算法，成功實(shí)現(xiàn)了考慮關(guān)斷時(shí)間的地面瞬變電磁三維帶地形反演；2021年霍軍鵬研究了傾斜地形對地面回線源的影響，提出傾斜地形校正方法，并將其成功應(yīng)用于存在地形影響的陜北侏羅系煤田采空區(qū)探測工程中。

發(fā)射源下方存在地形起伏的情況，發(fā)射源、地形等多種變化參數(shù)的影響相互耦合，難以分離各種影響因素的作用，數(shù)據(jù)解釋也異常困難。為了避免地形造成復(fù)雜影響，實(shí)際勘探中鋪設(shè)發(fā)射源和接收點(diǎn)會盡量避開起伏地形，此時(shí)不需要考慮發(fā)射源形變，只需要考慮源外地形的影響。目前，對地形影響的分析大多從測道曲線圖或者視電阻率斷面圖進(jìn)行定性判斷，較少從整個(gè)觀測時(shí)間段內(nèi)量化評估地形影響的程度。為了量化評估地形影響的程度，需要借助三維數(shù)值模擬方法分別獲得帶地形、不帶地形的瞬變電磁響應(yīng)。帶地形的瞬變電磁響應(yīng)可以利用各種數(shù)值方法進(jìn)行模擬，如積分方程(IE)法、有限差分法、有限元(FE)法、有限體積法等。其中，有限元法憑借對復(fù)雜地電模型適用性強(qiáng)、形成的總體系數(shù)矩陣稀疏且對稱、易于并行等優(yōu)點(diǎn)，已成為解決電磁法三維正演問題應(yīng)用最廣泛的方法。本文在不考慮發(fā)射源形變基礎(chǔ)上，基于有限元法計(jì)算了帶地形的地面瞬變電磁響應(yīng)。首先，以均勻半空間模型解析解和前人塊狀高導(dǎo)異常體模型的數(shù)值結(jié)果為標(biāo)準(zhǔn)，驗(yàn)證了瞬變電磁法三維正演建模與仿真的正確性；然后，在矩形回線源外構(gòu)造了山脊地形、山谷地形，以均勻半空間的下階躍響應(yīng)為參考，將地形與發(fā)射源的距離、大地電導(dǎo)率、異常體電導(dǎo)率作為變化參數(shù)，給出了整個(gè)觀測時(shí)間段內(nèi)在不同接收點(diǎn)處地形造成的相對誤差分布，分析并總結(jié)了源外地形對電磁響應(yīng)的影響規(guī)律。

1 瞬變電磁三維正演方法

本文將利用COMSOL AC/DC模塊的磁場接口進(jìn)行瞬變電磁三維正演。COMSOL Multiphy-sics是一款以有限元法為基礎(chǔ)，求解偏微分方程(單場)或偏微分方程組(多場)來實(shí)現(xiàn)真實(shí)物理現(xiàn)象的仿真軟件，廣泛應(yīng)用在聲學(xué)、化學(xué)反應(yīng)、電磁學(xué)、流體動力學(xué)、熱傳導(dǎo)、微波工程、光學(xué)等諸多領(lǐng)域。其中，COMSOL AC/DC模塊可以對電場、磁場和電磁場在穩(wěn)態(tài)、瞬態(tài)、頻域下的應(yīng)用進(jìn)行建模與仿真。該模塊的基礎(chǔ)理論是麥克斯韋方程組和邊界條件，通過矢量有限元法求解給定幾何模型的介質(zhì)內(nèi)電磁場分布。

為了簡化麥克斯韋方程組的求解過程，COMSOL AC/DC模塊引入了標(biāo)量位函數(shù)和矢量位函數(shù)替代電場()和磁場()。其表達(dá)式為

(1)

(2)

式中：為矢量磁位；為標(biāo)量電位函數(shù)；為時(shí)間。

在準(zhǔn)靜態(tài)(不考慮位移電流)條件下，有限元法的時(shí)間域控制方程組為

(3)

(4)

式中：為外加場源電流密度；為地電模型電導(dǎo)率；為磁導(dǎo)率，通常取真空中磁導(dǎo)率。

令Γ為計(jì)算區(qū)域外邊界，為外邊界單位法向量，采用自然Dirichlet邊界條件，即

(×)=0

(5)

=0

(6)

根據(jù)矢量有限元法，利用弱公式化，對計(jì)算區(qū)域內(nèi)模型的控制方程進(jìn)行離散化，形成相應(yīng)的數(shù)值模型方程組，然后通過迭代法或直接法求解各個(gè)網(wǎng)格單元處矢量磁位的近似解。具體離散過程可參考文獻(xiàn)[22]。瞬變電磁法屬于時(shí)間域電磁法，本文直接在時(shí)間域獲得電磁場響應(yīng)數(shù)值解，在COMSOL AC/DC模塊的磁場接口下選擇瞬態(tài)研究，時(shí)間步進(jìn)方法選擇后向差分公式。

電磁勘探法的頻率一般較低，需要足夠大的計(jì)算區(qū)域才能避免較為嚴(yán)重的截?cái)嘈?yīng)，從而獲得較為精確的數(shù)值解。為了兼顧精度和效率，在計(jì)算區(qū)域外圍可以引入無限元域。無限元表示沿特定坐標(biāo)軸拉伸的區(qū)域，其作用是近似形成無限大的域，有效避免截?cái)嘈?yīng)。在每個(gè)坐標(biāo)方向上，無限元域內(nèi)實(shí)際建模坐標(biāo)被歸一化為0～1的無量綱坐標(biāo)()；拉伸尺度可以定義為的函數(shù)(())，通過拉伸尺度函數(shù)得到經(jīng)過延伸的坐標(biāo)，可以認(rèn)為是在給定方向上的新位置；在每個(gè)方向上分別拉伸，最終形成總拉伸區(qū)域。拉伸尺度函數(shù)表達(dá)式為

(7)

式中：為整個(gè)計(jì)算區(qū)域的對角線長度；為大于1的無量綱數(shù)字，由(+)/計(jì)算，其中為預(yù)設(shè)無限元域的拉伸厚度，應(yīng)設(shè)置為遠(yuǎn)大于，以獲得較大的拉伸值，一般可設(shè)置為1 000。

由式(7)可知，拉伸尺度函數(shù)具有非線性特征，靠近內(nèi)邊界的無限元域單元拉伸尺度較小，靠近外邊界的無限元域單元拉伸尺度較大。

2 三維正演驗(yàn)證

本節(jié)將驗(yàn)證在COMSOL AC/DC模塊的磁場接口下進(jìn)行瞬變電磁三維正演的正確性。瞬變電磁法屬于時(shí)間域電磁法，其觀測時(shí)間范圍廣、電磁響應(yīng)數(shù)值量級跨度大，對建模的細(xì)節(jié)參數(shù)(計(jì)算區(qū)域尺度、網(wǎng)格剖分、時(shí)間步長設(shè)置)非常敏感，不合理的建模參數(shù)會導(dǎo)致較大的誤差。為了保證帶地形的瞬變電磁響應(yīng)計(jì)算以及后續(xù)分析的有效性，需要通過解析解以及前人已有數(shù)值解來驗(yàn)證瞬變電磁三維正演的正確性，從而獲得合理的建模參數(shù)。本文的計(jì)算環(huán)境為Win10系統(tǒng)下的COMSOL Multiphysics 5.3版本，硬件為Intel(R) Core(TM) i7-8550U處理器、8G內(nèi)存的個(gè)人筆記本電腦。

首先，基于均勻半空間模型進(jìn)行驗(yàn)證。將地電模型設(shè)置為0.01 S·m的均勻半空間模型，矩形發(fā)射回線尺度為100 m×100 m，電流為1 A，接收點(diǎn)位于回線中心，這種情況下可以采用等面積的圓形回線中心點(diǎn)的解析解來檢驗(yàn)其數(shù)值精度。均勻半空間三維正演建模如圖1(a)所示。計(jì)算區(qū)域?yàn)榘霃? km的圓球，球體被分成兩個(gè)半球，其中一個(gè)半球作為大地區(qū)域，另一個(gè)作為空氣區(qū)域。矩形回線源布設(shè)在大地半球與空氣半球的交界平面上，回線中心與坐標(biāo)原點(diǎn)重合。由于計(jì)算區(qū)域有限，為了避免截?cái)嘈?yīng)，將圓球的外圍劃分出300 m厚的區(qū)域作為無限元域。經(jīng)過網(wǎng)格剖分后，共計(jì)54 158個(gè)網(wǎng)格域單元。其中，內(nèi)圓球部分選擇自由四面體網(wǎng)格；外圓球部分為無限元域。無限元域的網(wǎng)格剖分選擇掃掠方式，網(wǎng)格設(shè)置為棱柱，這些棱柱單元以內(nèi)圓球的三角面為輻射源面。為了保證精度，在發(fā)射源和接收點(diǎn)部分進(jìn)行了網(wǎng)格加密，如圖1(b)所示。這種局部加密是保證精度和提高計(jì)算效率的一種平衡實(shí)施方式。發(fā)射電流為下階躍電流，時(shí)間步選擇非均勻?qū)?shù)定義方法，觀測時(shí)間為(1～10 000)×10s。其中，在(1～100)×10s內(nèi)每十倍頻率步數(shù)為10，在(1～10)×10s內(nèi)每十倍頻率步數(shù)為50，在(1～10)×10s內(nèi)每十倍頻率步數(shù)為100。這種非均勻的時(shí)間步長設(shè)置能夠保證數(shù)值精度和計(jì)算效率。時(shí)間步進(jìn)方法采用向后差分方式，案例總耗時(shí)約為11 min。

圖1 均勻半空間三維正演模型和網(wǎng)格剖分示意圖Fig.1 Schematic Views of Homogenous Half-space 3D Forward Model and Mesh Generation

圖2(a)為均勻半空間模型上矩形回線源中心點(diǎn)處磁場分量()時(shí)域有限元解與解析解對比。

圖2 均勻半空間模型上矩形回線源磁場z分量響應(yīng)Fig.2 Magnetic Responses of z Component of Rectangular Loop Source over Homogenous Half-space Model

從圖2(a)可以看出，有限元解與解析解完全重合。圖2(b)為有限元解與解析解的相對誤差，整個(gè)觀測時(shí)間范圍內(nèi)相對誤差均小于2%，由于在(1～100)×10s內(nèi)有充足的時(shí)間采樣，其相對誤差小于0.5%。

以含塊狀高導(dǎo)異常體的三維大地模型進(jìn)一步驗(yàn)證瞬變電磁三維正演的有效性。如圖3(a)所示，該地電模型將電導(dǎo)率為2.0 S·m的塊狀高導(dǎo)異常體埋于電導(dǎo)率為0.1 S·m的均勻半空間模型中，塊狀體頂部距離地面30 m，其長、寬、高分別為100、40、30 m，空氣電導(dǎo)率設(shè)置為1.0×10S·m。矩形發(fā)射回線尺寸為100 m×100 m，測點(diǎn)位于回線中心，電流為1 A。該模型最早由Newman等利用積分方程法實(shí)施三維正演時(shí)引入。計(jì)算區(qū)域和網(wǎng)格剖分設(shè)置與均勻半空間模型案例一致。球形計(jì)算區(qū)域外圍也包含了300 m厚的無限元域，整個(gè)計(jì)算區(qū)域被劃分為62 696個(gè)網(wǎng)格域單元，其中發(fā)射回線與接收點(diǎn)處進(jìn)行了網(wǎng)格加密，對塊狀高導(dǎo)異常體也進(jìn)行了局部網(wǎng)格加密。圖3(b)對比了本文的有限元解、Newman等的積分方程解、Wang等的時(shí)域有限差分(FDTD)解、Li等的有限元解。在(1～10)×10s內(nèi)，本文的有限元解與積分方程解、其他有限元解吻合較好，比時(shí)域有限差分解略??；在(1～10)×10s內(nèi)，本文的有限元解與時(shí)域有限差分解、其他有限元解吻合較好，比積分方程解略大。需要說明的是，由于積分方程解、時(shí)域有限差分解、其他有限元解都是直接從這些研究學(xué)者的文章中讀取，其精度低于原始計(jì)算的數(shù)值解，所以沒有計(jì)算本研究與這些數(shù)值解的相對誤差。

圖3 含塊狀高導(dǎo)異常體大地模型及其響應(yīng)曲線Fig.3 Geodetic Model with Massive High Conductive Anomaly and Its Response Curve

3 源外不同地形對電磁響應(yīng)的影響特征分析

本節(jié)考察矩形回線源外不同地形對電磁響應(yīng)的影響規(guī)律。所有案例中，回線源均保持位置、形狀、發(fā)射電流不變。案例總體上分成兩種情況：不含異常體和含異常體。其中，不含異常體情況下，地形存在山脊、山谷兩種起伏地形，將地形設(shè)置在回線源外的不同距離處；在大地電導(dǎo)率不變的條件下，分別考察不同距離時(shí)山脊、山谷地形對電磁響應(yīng)的影響規(guī)律；在距離相同的條件下，考察不同大地電導(dǎo)率時(shí)山脊地形對電磁響應(yīng)的影響規(guī)律。含異常體情況下，起伏地形為山脊地形，與發(fā)射源保持固定距離，發(fā)射源下方一定埋深處存在異常體，考察異常體電導(dǎo)率不同時(shí)山脊地形對電磁響應(yīng)的影響規(guī)律。

3.1 不含異常體情況

3.1.1 不同參數(shù)下山脊地形對電磁響應(yīng)的影響

圖4為回線源外存在山脊地形的模型示意圖。該模型采用100 m×100 m矩形回線源，源中心與坐標(biāo)原點(diǎn)重合，回線源中心點(diǎn)在坐標(biāo)原點(diǎn)，在回線源外采用六面體模擬出山脊地形；該六面體上、下表面分別為200 m×200 m和400 m×400 m的矩形，其中心在=0平面上投影重合在軸上，上表面與地面的距離為100 m；空氣和大地電導(dǎo)率分別設(shè)置為1.0×10S·m和0.01 S·m。如圖4(b)所示，將發(fā)射回線右棱邊與山脊地形最左邊棱邊的距離設(shè)置為可變，考察距離分別為10、50、100、200 m處山脊地形對電磁響應(yīng)的影響。整個(gè)計(jì)算區(qū)域?yàn)榘霃? km的圓球，在地面上沿著軸從點(diǎn)(0，-25，0)到點(diǎn)(0，25，0)均勻鋪設(shè)52個(gè)接收點(diǎn)，時(shí)域響應(yīng)的觀測范圍為(1～10 000)×10s，分別計(jì)算出均勻半空間模型、山脊地形模型的接收點(diǎn)處時(shí)域磁場分量響應(yīng)。以均勻半空間模型的數(shù)值解為參考，計(jì)算存在山脊地形模型磁場響應(yīng)的相對誤差。

圖4 回線源外存在山脊地形的模型示意圖Fig.4 Schematic Views of Model Containing a Ridge Topography Outside the Rectangular Loop Source

圖5為不同距離處存在山脊地形的時(shí)域數(shù)值(TDEM)解關(guān)于均勻半空間模型的相對誤差分布?？v軸為已轉(zhuǎn)換成log域的時(shí)間()點(diǎn)，實(shí)際范圍為(1～10 000)×10s，橫軸為軸上接收點(diǎn)的距離。從圖5可以看出，誤差較大區(qū)域集中在某個(gè)時(shí)間段內(nèi)，接收點(diǎn)越靠近山脊地形，其相對誤差越大，隨著距山脊地形的距離不斷增大，整體誤差水平降低，誤差較大區(qū)域逐漸移至更晚時(shí)間上。這種誤差分布的形態(tài)主要與電磁波信號的穿透特性有關(guān)，發(fā)射源發(fā)射下階躍電流后，形成非常寬頻帶的激勵，引起地下介質(zhì)的感應(yīng)渦流。頻率較高的電磁信號穿透距離較淺；頻率較低的電磁信號穿透距離更遠(yuǎn)。頻率較高的電磁信號引起近距離地下介質(zhì)的感應(yīng)渦流，該感應(yīng)渦流形成的二次場在較早時(shí)間上回到觀測點(diǎn)；頻率較低的電磁信號引起遠(yuǎn)距離地下介質(zhì)的感應(yīng)渦流，該感應(yīng)渦流形成的二次場在較晚時(shí)間上回到觀測點(diǎn)。接收時(shí)間越早，包含的頻率成分越高；接收時(shí)間越晚，包含的頻率成分越低。當(dāng)山脊地形處于不同距離時(shí)，引起感應(yīng)渦流的電磁波頻率成分不同，因此，距離越遠(yuǎn)，誤差較大區(qū)域在時(shí)間上更晚。同時(shí)，山脊地形相對頻率較高電磁波的電尺寸更大，因而近距離時(shí)山脊地形造成的影響也更大一些。從整體誤差水平來看，與發(fā)射源距離10 m的山脊地形引起的相對誤差高達(dá)14%左右，對實(shí)際勘探的影響較大；隨著山脊地形與發(fā)射源的距離增大，誤差逐漸減小，距離200 m的山脊地形引起的相對誤差小于2%，對實(shí)際勘探的影響可以忽略。因此，在實(shí)際勘探中，將發(fā)射源盡量遠(yuǎn)離山脊能夠減小地形影響。

圖5 不同距離處存在山脊地形的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布Fig.5 Distributions of Relative Errors of TDEM Response of the Model Containing Ridge Topography at Different Distances with Respect to Homogeneous Half-space Model

圖6為不同大地電導(dǎo)率時(shí)存在山脊地形的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布。該山脊地形固定在距離發(fā)射源100 m處，整個(gè)大地電導(dǎo)率依次設(shè)置為0.001、0.005、0.050、0.100 S·m。從圖6可以看出，不同大地電導(dǎo)率的整體誤差水平基本一致，誤差較大區(qū)域集中在某個(gè)時(shí)間段內(nèi)，接收點(diǎn)越靠近山脊地形，其相對誤差越大，隨著大地電導(dǎo)率的增大，誤差中心逐漸移至更晚時(shí)間上。這種誤差分布形態(tài)也是與電磁波信號的穿透特性有關(guān)，電磁波在導(dǎo)電媒質(zhì)中的傳播距離與頻率、電導(dǎo)率有關(guān)；根據(jù)趨膚深度計(jì)算公式，當(dāng)電磁信號頻率固定不變時(shí)，媒質(zhì)的電導(dǎo)率越高，電磁能量衰減越快，穿透深度也越淺。因此，對于某一固定距離的山脊地形，當(dāng)大地電導(dǎo)率增大時(shí)，頻率更低的電磁信號才能傳播到此處，進(jìn)而引起山脊地形的感應(yīng)渦流，感應(yīng)渦流形成的二次場將在更晚時(shí)間上回到觀測點(diǎn)。從整體誤差水平來看，不同大地電導(dǎo)率在不同時(shí)間段上引起的誤差水平相近，傳播到山脊地形的頻率成分所覆蓋的區(qū)域大致相當(dāng)。因此，在接收點(diǎn)處，山脊區(qū)域的貢獻(xiàn)也是相對固定的。綜上所述，對于固定距離的山脊地形，大地電導(dǎo)率的改變主要影響山脊地形的電磁響應(yīng)時(shí)間；相對均勻半空間模型，山脊地形引起的誤差水平基本一致；大地電導(dǎo)率越高，山脊地形引起的誤差分布在時(shí)間上越晚。

圖6 不同大地電導(dǎo)率時(shí)存在山脊地形(距離為100 m)的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布Fig.6 Distributions of Relative Errors of TDEM Response of the Model Containing Ridge Topography (Distance of 100 m) at Different Conductivities of Earth with Respect to Homogeneous Half-space Model

3.1.2 不同參數(shù)下山谷地形對電磁響應(yīng)的影響

圖7為回線源外存在山谷地形的模型示意圖。該模型采用100 m×100 m矩形回線源，源中心與坐標(biāo)原點(diǎn)重合，回線源中心點(diǎn)在坐標(biāo)原點(diǎn)，在回線源外采用六面體模擬出山谷地形；該六面體上、下表面分別為400 m×400 m和200 m×200 m的矩形，其中心在=0平面上投影重合在軸上，下表面與地面的距離為100 m；空氣和大地電導(dǎo)率分別設(shè)置為1.0×10、0.01 S·m。如圖7(b)所示，將發(fā)射回線右棱邊與山谷地形最左邊棱邊的距離設(shè)置為可變，考察距離分別為10、50、100、200 m處存在山谷地形對電磁響應(yīng)的影響。與圖4存在山脊地形的模型示意圖相比，山谷地形與山脊地形的幾何結(jié)構(gòu)是關(guān)于地平面對稱的，其他參數(shù)設(shè)置都是相同的，這樣的模型設(shè)定是為了更好地對比不同類型(凹陷、凸起)地形的影響。

圖7 回線源外存在山谷地形的模型示意圖Fig.7 Schematic Views of Model Containing a Valley Topography Outside the Rectangular Loop Source

圖8為不同距離處存在山谷地形的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布。縱軸為已轉(zhuǎn)換成log域的時(shí)間點(diǎn)，實(shí)際范圍為(1～10 000)×10s，橫軸為軸上接收點(diǎn)的距離。從圖8可以看出，較大誤差形成的區(qū)域也是集中在一定時(shí)間段內(nèi)，接收點(diǎn)越靠近山谷地型，其相對誤差越大，隨著距山谷地形的距離不斷增大，整體誤差水平降低，較大誤差區(qū)域逐漸移至更晚時(shí)間上，距離上靠近山谷地形的誤差區(qū)域呈現(xiàn)減小趨勢。與山脊地形的時(shí)域數(shù)值解相對誤差分布(圖5)相比，山谷地形的時(shí)域數(shù)值解相對誤差分布形態(tài)與山脊地形呈現(xiàn)相似的特點(diǎn)，引起這種誤差分布形態(tài)的原因也是一樣的。不同的是，山谷地形引起的誤差水平更高，與發(fā)射源距離10 m的山谷地形引起的相對誤差高達(dá)25%左右，而同等距離山脊地形引起的相對誤差為14%左右，可見山谷地形比山脊地形對實(shí)際勘探的影響更大；隨著山谷地形與發(fā)射源的距離增大，相對誤差逐漸減小，距離為200 m的山谷地形引起的相對誤差小于4%。因此，為了減小地形對電磁響應(yīng)的影響，有效措施依然是盡量遠(yuǎn)離山谷地形。

圖8 不同距離處存在山谷地形的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布Fig.8 Distributions of Relative Errors of TDEM Response of the Model Containing Valley Topography at Different Distances with Respect to Homogeneous Half-space Model

3.1.3 不同距離下兩種地形的時(shí)間域電磁響應(yīng)

以均勻半空間模型的電磁響應(yīng)作為參考，考察不同距離下山脊、山谷地形的時(shí)間域電磁響應(yīng)(圖9)。從圖9(a)可以看出，山脊地形的貢獻(xiàn)是增大一定時(shí)間段的電磁響應(yīng)。根據(jù)前文分析，發(fā)射源引起地下導(dǎo)電介質(zhì)的感應(yīng)渦流，感應(yīng)渦流產(chǎn)生二次場，然后觀測點(diǎn)處接收該二次場；與均勻半空間模型相比，山脊地形相當(dāng)于延展了一部分導(dǎo)電介質(zhì)區(qū)域，增加了一定的感應(yīng)渦流，導(dǎo)致產(chǎn)生了更強(qiáng)的二次場。從圖9(b)可以看出，山谷地形的影響是減小一定時(shí)間段的電磁響應(yīng)。與均勻半空間模型相比，山谷地形相當(dāng)于減小了一部分導(dǎo)電介質(zhì)區(qū)域，缺少了部分區(qū)域的感應(yīng)渦流，導(dǎo)致二次場比均勻半空間模型更弱。

圖9 兩種地形在中心點(diǎn)處的時(shí)域數(shù)值解Fig.9 TDEM Responses at the Central Receiving Point for Two Kinds of Topography

3.2 含異常體情況

圖10為含異常體與山脊地形的模型示意圖。矩形回線源的參數(shù)和山脊地形的幾何尺寸與3.1.1節(jié)一致。如圖10(b)所示，將發(fā)射回線右棱邊與山脊地形最左邊棱邊的距離設(shè)置為100 m。異常體的頂部埋深為50 m，長為100 m，寬為200 m，高為40 m，其中心在=0平面上投影為坐標(biāo)原點(diǎn)。在本節(jié)所有案例中，異常體的幾何信息不變，其電導(dǎo)率是可變的，以考察不同電導(dǎo)率異常體對電磁響應(yīng)的影響。計(jì)算區(qū)域與網(wǎng)格剖分的參數(shù)與前面案例一致。在地面上沿著軸從點(diǎn)(0，-25，0)到點(diǎn)(0，25，0)均勻鋪設(shè)52個(gè)接收點(diǎn)，分別計(jì)算出均勻半空間模型、異常體(不含地形)模型、異常體(含地形)模型這3種情況下接收點(diǎn)處時(shí)域磁場分量響應(yīng)。

圖10 含異常體與山脊地形的模型示意圖Fig.10 Schematic Views of Model Containing Ridge Topography and Anomaly

3.2.1 異常體(不含地形)模型

以均勻半空間模型的數(shù)值解為參考，計(jì)算異常體(不含地形)模型關(guān)于均勻半空間模型磁場響應(yīng)的相對誤差，以考察不同電導(dǎo)率異常體對電磁響應(yīng)的影響。異常體電導(dǎo)率依次設(shè)定為1.00、0.50、0.10、0.05 S·m。圖11為不同電導(dǎo)率異常體(不含地形)模型的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布。從圖11可以看出，不同電導(dǎo)率異常體對電磁響應(yīng)的影響主要集中在一定的時(shí)間段內(nèi)，隨著異常體電導(dǎo)率的減小，異常體引起的誤差水平逐漸降低，同時(shí)誤差中心區(qū)域在時(shí)間上逐漸移至更早期。如3.1節(jié)的案例分析一樣，異常體引起的誤差分布形態(tài)主要與電磁波信號的穿透特性、感應(yīng)渦流有關(guān)。異常體電導(dǎo)率越高，不僅能夠產(chǎn)生更強(qiáng)的二次場響應(yīng)，同時(shí)二次場在接收時(shí)間上的影響持續(xù)時(shí)間越長，其最大影響中心在時(shí)間上也越晚。

圖11 不同電導(dǎo)率異常體(不含地形)模的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布Fig.11 Distributions of Relative Errors of TDEM Response of the Model Containing Anomaly (Without Topography) at Different Conductivities with Respect to Homogeneous Half-space Model

為了更好地理解不同電導(dǎo)率異常體對電磁響應(yīng)的影響，本文給出了在回線源中心處不同電導(dǎo)率異常體(不含地形)模型的時(shí)間域磁場分量響應(yīng)(圖12)。從圖12可以看出：不同電導(dǎo)率異常體磁場響應(yīng)逐漸增強(qiáng)的起始時(shí)間是一致的，這是由于異常體的位置是固定的，電磁信號傳播到不同電導(dǎo)率異常體所經(jīng)歷的導(dǎo)電媒質(zhì)和距離是不變的；同時(shí)，可以觀察到異常體的貢獻(xiàn)是增大一定時(shí)間段內(nèi)的電磁響應(yīng)，電導(dǎo)率越高，其二次場響應(yīng)增大的幅度越高，持續(xù)時(shí)間越長。

圖12 不同電導(dǎo)率異常體(不含地形)模型在中心點(diǎn)處時(shí)域數(shù)值解Fig.12 TDEM Responses at the Central Receiving Point for Models with Anomaly (Without Topography) at Different Conductivities

3.2.2 異常體(含地形)模型

本文計(jì)算了同時(shí)包含異常體和地形模型的時(shí)域數(shù)值解，然后分別以均勻半空間模型、異常體(不含地形)模型的數(shù)值解為參考，計(jì)算了異常體(含地形)模型分別關(guān)于這兩種模型的磁場響應(yīng)的相對誤差。其中，異常體電導(dǎo)率依次設(shè)定為0.10、0.05 S·m。圖13為不同電導(dǎo)率異常體(含地形)模型的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布。與圖11(c)、(d)不含地形的情況相比，圖13整體誤差水平與分布形態(tài)基本一致，難以分辨地形的影響特征。因此，為了單獨(dú)觀察地形對電磁響應(yīng)的影響，圖14給出了不同電導(dǎo)率異常體(含地形)模型的時(shí)域數(shù)值解關(guān)于異常體(不含地形)模型的相對誤差分布。當(dāng)異常體電導(dǎo)率不同時(shí)，圖14的誤差分布形態(tài)和水平基本一致，這也與圖5(c)的結(jié)論基本一致。這說明地形、異常體對電磁響應(yīng)的影響可以分離。對于一個(gè)包含地形效應(yīng)的實(shí)測數(shù)據(jù)，如果異常體產(chǎn)生的二次場比較強(qiáng)，地形效應(yīng)的影響相對非常小，可以忽略，將地面當(dāng)作理想平面來處理數(shù)據(jù)；如果異常體產(chǎn)生的二次場較弱，則不能忽略地形效應(yīng)的影響，此時(shí)如將地面當(dāng)作理想平面來處理數(shù)據(jù)，則異常體電導(dǎo)率的反演值與真實(shí)值差距較大；如果地形效應(yīng)比較強(qiáng)，例如起伏地形比較靠近發(fā)射源，不考慮地形去反演數(shù)據(jù)，也容易將地形效應(yīng)產(chǎn)生的二次場差異當(dāng)作是異常體。

圖13 不同電導(dǎo)率異常體(含地形)模型的時(shí)域數(shù)值解關(guān)于均勻半空間模型的相對誤差分布Fig.13 Distributions of Relative Errors of TDEM Response of the Model Containing Anomaly (with Topography) at Different Conductivities with Respect to Homogeneous Half-space Model

圖14 不同電導(dǎo)率異常體(含地形)模型的時(shí)域數(shù)值解關(guān)于異常體(不含地形)模型的相對誤差分布Fig.14 Distributions of Relative Errors of TDEM Response of the Model Containing Anomaly (with Topography) at Different Conductivities with Respect to the Model Containing Anomaly (Without Topography)

4 結(jié) 語

本文針對地面電磁勘探法，不考慮發(fā)射源形變，構(gòu)建了矩形回線源外存在不同地形的多種勘探場景，基于有限元法計(jì)算了不同場景下的電磁響應(yīng)，考察了地形、電性等參數(shù)變化對瞬變電磁響應(yīng)的影響規(guī)律。

(1)源外地形對電磁響應(yīng)的影響主要集中在某個(gè)時(shí)間段內(nèi)，越靠近地形的觀測點(diǎn)，受影響程度越嚴(yán)重；隨著與發(fā)射源距離增大，地形作用減弱，受影響區(qū)域也推移到較晚期時(shí)間段上。

(2)當(dāng)?shù)匦闻c發(fā)射源的距離不變時(shí)，隨著大地電導(dǎo)率增大，地形主要影響區(qū)域趨向于集中在更晚期時(shí)間上，但是受影響的程度基本一致。

(3)相同幾何尺寸的山谷地形比山脊地形對電磁響應(yīng)的影響更大。山谷地形對電磁響應(yīng)的影響是減小某個(gè)時(shí)間段內(nèi)二次場響應(yīng)，而山脊地形則是增強(qiáng)一定時(shí)間段內(nèi)二次場響應(yīng)。

(4)地形、異常體對電磁響應(yīng)的影響可以分離，異常體的變化并不會增強(qiáng)或減弱地形對電磁響應(yīng)的貢獻(xiàn)。對于包含地形效應(yīng)的勘探數(shù)據(jù)，如果異常體產(chǎn)生的二次場比較強(qiáng)，地形效應(yīng)的影響相對較小，可以忽略地形；如果異常體產(chǎn)生的二次場較弱，則不能忽略地形效應(yīng)的影響。

謹(jǐn)以拙作祝賀西安地質(zhì)調(diào)查中心組建六十周年！西安地質(zhì)調(diào)查中心一直承擔(dān)著國家在西北地區(qū)的基礎(chǔ)性、公益性地質(zhì)調(diào)查和戰(zhàn)略性礦產(chǎn)勘查任務(wù)及相關(guān)綜合研究工作，取得了豐碩的研究成果和榮譽(yù)。西安石油大學(xué)與西安地質(zhì)調(diào)查中心多年來在油氣成藏地質(zhì)學(xué)、非常規(guī)油氣地質(zhì)與勘探、深部地質(zhì)構(gòu)造、含油氣盆地構(gòu)造學(xué)、油氣田開發(fā)地質(zhì)學(xué)及地震勘探等多個(gè)領(lǐng)域進(jìn)行了深入合作，形成了良好的合作關(guān)系。作為西安石油大學(xué)的一份子，期待未來參與到西安石油大學(xué)和西安地質(zhì)調(diào)查中心進(jìn)一步的合作中，攜手為國家能源安全和可持續(xù)發(fā)展做出新的貢獻(xiàn)！