王 智,潘和平

(中國(guó)地質(zhì)大學(xué)(武漢)地球物理與空間信息學(xué)院,湖北武漢430074)

井地電阻率法是在井中供電地面接收電磁場(chǎng)的一類(lèi)電法勘探方法,最早由前蘇聯(lián)科學(xué)家于1958年提出[1],用來(lái)圈定煤層或油氣藏水平邊界。通常情況下，井地電阻率法的結(jié)果以視電阻率等值線圖或異常電位等值線圖表示。國(guó)內(nèi)外對(duì)井地電阻率法的正演數(shù)值模擬做了很多的研究。Scriba[2]和Spitzer[3]利用邊界元實(shí)現(xiàn)了三維地電斷面的正演問(wèn)題;Pridmore等[4]用有限元實(shí)現(xiàn)了三維正演模擬的計(jì)算;徐凱軍[5]和王志剛等[6]利用有限差分完成了垂直線源井地電位法的三維正演模擬;柯敢攀等[7]利用有限元實(shí)現(xiàn)了垂直線源井地電位法的三維數(shù)值模擬;戴前偉等[8]利用有限元法完成了2.5維的井地電位法的正演問(wèn)題;賈正元等[9]利用有限差分實(shí)現(xiàn)了三維井地電阻率法的數(shù)值模擬,并對(duì)點(diǎn)源、線源、傾斜源三種情況下產(chǎn)生的異常特征進(jìn)行了詳細(xì)的總結(jié);劉海飛等[10]采用有限元實(shí)現(xiàn)了三維連續(xù)電性介質(zhì)的數(shù)值模擬。

在井地電阻率圈定異常體邊界的研究方面,張?zhí)靷怺11-12]用物理水槽實(shí)驗(yàn)分別進(jìn)行在井中無(wú)套管和井中有套管情況下確定油氣藏邊界的研究,證實(shí)了利用井地電阻率法確定異常體邊界的可行性;并在后來(lái)提出了采用非無(wú)窮遠(yuǎn)處三極剖面法確定油氣藏邊界[13],且應(yīng)用此方法在新疆等地成功發(fā)現(xiàn)了小塊油氣藏[14-15]。湯井田等[16]利用井地電阻率法歧離率確定高阻油藏邊界;戴前偉等[17]利用梯度算子與拉普拉斯算子對(duì)異常響應(yīng)進(jìn)行增強(qiáng),取得了一定的效果。戴前偉等[8]利用不同點(diǎn)源深度時(shí)的視電阻率異常斷面圖的中心深度位置來(lái)確定異常體的上、下邊界;湯井田等[16]利用歧離率確定高阻油氣藏的邊界取得一定效果,但只限于鉆孔剛好打到油氣藏上,因此歧離率法發(fā)現(xiàn)井旁盲礦的能力很弱。

有時(shí)井地電阻率法的異常響應(yīng)信號(hào)非常弱小,特征不明顯，特別是當(dāng)異常體的埋深較深時(shí)，這種現(xiàn)象更加突出。因此，在總結(jié)前人研究成果的基礎(chǔ)上，我們從井地電阻率法滿足的邊值問(wèn)題出發(fā)，采用點(diǎn)電流源，井中供電、地表接收的觀測(cè)方式，利用有限元法實(shí)現(xiàn)了三維井地電阻率法的數(shù)值模擬；計(jì)算了多種不同情況下三維異常體模型的異常電位水平導(dǎo)數(shù)，分析了異常電位水平導(dǎo)數(shù)的極值與異常體邊界對(duì)應(yīng)關(guān)系的影響因素；引入歸一化總水平導(dǎo)數(shù)法[18]，以研究提高深部異常響應(yīng)特征和增強(qiáng)異常響應(yīng)的邊界范圍。

1 三維井地電阻率有限元數(shù)值模擬算法

1.1 基本原理

在穩(wěn)恒電流源供電下的研究區(qū)域中，各點(diǎn)電位滿足如下的邊值問(wèn)題:

(1)

與(1)式對(duì)應(yīng)的變分問(wèn)題為

(2)

式中:A是與源項(xiàng)有關(guān)的坐標(biāo)向量;σ為介質(zhì)的電導(dǎo)率;ωA是地面源點(diǎn)對(duì)地下區(qū)域的立體角,若源點(diǎn)在地下ωA=4π,若源點(diǎn)在地面ωA=2π;I為電流源項(xiàng);u為總電位;Γs為地面邊界;Γ∞為無(wú)窮遠(yuǎn)邊界;Ω為Γs與?！藿M成的封閉區(qū)域;r為源點(diǎn)到測(cè)點(diǎn)的距離;cos(r,n)是矢徑r與外法向方向n的夾角余弦。

1.2 有限單元法

本文利用有限單元法求解上述變分問(wèn)題,分為4個(gè)步驟。

1.2.1 網(wǎng)格剖分

采用四面體剖分,如圖1所示。

圖1 四面體單元

1.2.2 線性插值

設(shè)ui(i=1,2,3,4)是單元中4個(gè)節(jié)點(diǎn)的電位值,則四面體單元e內(nèi)任一點(diǎn)p(x,y,z)電位用這4個(gè)角點(diǎn)的電位進(jìn)行線性插值近似得到

(3)

式中：Ni(i=1,2,3,4)是形函數(shù),它是x,y,z的線性函數(shù):

(4)

式中：V是四面體的體積;Vi是插值點(diǎn)p(x,y,z)與四面體其它3個(gè)角點(diǎn)所組成的四面體體積;ai,bi,ci,di(i=1,2,3,4)是與四面體單元頂點(diǎn)坐標(biāo)有關(guān)的常數(shù)。

1.2.3 單元積分

變分問(wèn)題(2)式中，第1項(xiàng)對(duì)應(yīng)的單元積分為

(5)

其中,

(6)

第2項(xiàng)為電源項(xiàng),只與電源點(diǎn)的位置有關(guān);第3項(xiàng)為邊界單元項(xiàng),若單元的一個(gè)面的1,2,3節(jié)點(diǎn)落在邊界上,則

(7)

(8)

式中:Δ為單元e在邊界上的三角形的面積。

1.2.4 方程組的求解

將單元積分結(jié)果擴(kuò)展成由全體節(jié)點(diǎn)組成的矩陣后,進(jìn)行全部單元相加得到如下的線性方程組:

(9)

式中：K為剛度系數(shù)矩陣;U為待求電位矩;P為與電源有關(guān)的矩陣。

目前解決這類(lèi)方程的求解方法有很多,直接法有GS法、奇異值分解、LDLT分解等;迭代法有Newton法、共軛梯度法及其改進(jìn)方法等。其中吳小平等[19]與柳建新等[20]分別用對(duì)稱(chēng)超松弛預(yù)條件共軛梯度法(SSOR-PCG)與不完全LU分解預(yù)處理的穩(wěn)定雙共軛梯度算法(ILU-BICGSTAB)實(shí)現(xiàn)了快速穩(wěn)定的計(jì)算。結(jié)合前人的研究,我們采用不完全LU分解預(yù)處理的穩(wěn)定雙共軛梯度算法求解線性方程組(9)[20],可得到各個(gè)節(jié)點(diǎn)的電位。為了消除源點(diǎn)附近的奇異性,我們采取異常電位法計(jì)算,對(duì)于異常電位的計(jì)算過(guò)程與總電位計(jì)算相似,最后歸結(jié)為求解大型方程組得到的異常電位,異常電位和正常電位相加得到總電位。

1.3 算法正確性的驗(yàn)證

基于上述理論推導(dǎo)過(guò)程,編制了相應(yīng)的程序。以下用兩個(gè)模型驗(yàn)證程序的正確性和有效性。

模型一為3層水平層狀模型(圖2),每層參數(shù)

為:第1層電阻率50Ω·m,厚度為10m;第2層電阻率為100Ω·m,厚度為10m;第3層電阻率為10Ω·m。圖2坐標(biāo)軸上A為源點(diǎn),M為測(cè)點(diǎn)。對(duì)該模型進(jìn)行正演計(jì)算,并將計(jì)算結(jié)果與數(shù)值濾波法計(jì)算的結(jié)果進(jìn)行對(duì)比,如圖3所示。從圖中看出二者吻合得非常好,最大誤差為0.03%。

圖2 水平層狀模型

圖3 水平層狀模型數(shù)值計(jì)算與濾波法結(jié)果對(duì)比

模型二為地下低阻球體模型(圖4),源點(diǎn)和測(cè)點(diǎn)均在地面,模型參數(shù)為:A為源點(diǎn),坐標(biāo)為(50m,50m,0)M為測(cè)點(diǎn);圍巖電阻率100Ω·m,球體電阻率為10Ω·m,球體半徑為5m。球體中心埋深10m,坐標(biāo)為(50m,50m,10m),觀測(cè)電極距2m。對(duì)該模型進(jìn)行正演計(jì)算,將計(jì)算結(jié)果與均勻半空間內(nèi)球體外的電位解析解結(jié)果進(jìn)行對(duì)比(圖5),從圖5中看出二者吻合得非常好,最大誤差為2.85%。

圖4 低阻球體模型a 三維模型; b 模型XOY平面

圖5 低阻球體模型數(shù)值計(jì)算與解析法結(jié)果對(duì)比

2 水平導(dǎo)數(shù)極值圈定異常體范圍的影響因素分析

對(duì)于邊緣識(shí)別方法在重磁位場(chǎng)中的應(yīng)用研究較多[21],其中數(shù)值類(lèi)邊緣識(shí)別方法均是利用極大值位置或者零值位置來(lái)確定地質(zhì)體的邊緣位置,其理論基礎(chǔ)是二度體鉛垂臺(tái)階模型的重力異常特征,在該模型邊緣處重力異常總水平導(dǎo)數(shù)和解析信號(hào)振幅達(dá)到極大值，而垂向?qū)?shù)達(dá)到零值,故可利用這些特征位置來(lái)確定二度體鉛垂臺(tái)階的邊緣位置[22]。在電法勘探中,通常在異常體的邊界范圍處,即異常體與圍巖橫向接觸處,地表的電流密度變化同樣很大[16]。由此為利用水平導(dǎo)數(shù)確定其邊界位置提供了可能性。為了解水平導(dǎo)數(shù)極值確定異常體邊界范圍的有效性,下面從源的埋深、異常體的埋深、圍巖與異常體的差異、異常體的橫向位置對(duì)異常電位的水平導(dǎo)數(shù)的極值影響進(jìn)行分析?？偹綄?dǎo)數(shù)的定義如下[23]:

THD(x,y,0)=

(10)

式中:THD表示原始異常f的總水平導(dǎo)數(shù);f(x,y,0)為某一平面上的數(shù)據(jù)體,在本文中為地面的異常電位數(shù)據(jù)。

2.1 總水平導(dǎo)數(shù)的求取方式

我們給出兩種求取水平導(dǎo)數(shù)的方式。

方法一:從推導(dǎo)單元分析的積分方程過(guò)程中可得到異常電位法水平導(dǎo)數(shù)的解析求法。由式(3)可得:

(11)

其中，

(12)

同理可得

(13)

由式(12)與式(13)可知,某一點(diǎn)處的異常電位水平導(dǎo)數(shù)只與該單元內(nèi)各個(gè)網(wǎng)格節(jié)點(diǎn)的坐標(biāo)和異常電位有關(guān)。所以，在單元內(nèi)電位u的導(dǎo)數(shù)是常數(shù)。相鄰單元的?u/?x,?u/?y是不連續(xù)的。相鄰單元的u通過(guò)公共邊界點(diǎn)上的u連接起來(lái),u在全區(qū)域內(nèi)是連續(xù)的,但導(dǎo)數(shù)是不連續(xù)的。對(duì)于相鄰單元的公共節(jié)點(diǎn),可用各單元求得的導(dǎo)數(shù)的平均值,作為該點(diǎn)的導(dǎo)數(shù)值[24]。在正演的過(guò)程中需要計(jì)算ai,bi,ci,di(i=1,2,3,4)的值,因此在程序中保留原有的值,正演結(jié)束后,將地面四面體單元體內(nèi)的各個(gè)節(jié)點(diǎn)的異常電位回代到方程(12)與(13)中,即可得到該單元內(nèi)某一點(diǎn)處的異常電位的水平導(dǎo)數(shù)。

方法二:采用中心差分方法代替微積分進(jìn)行近似計(jì)算。詳細(xì)的計(jì)算過(guò)程請(qǐng)參見(jiàn)文獻(xiàn)[17]。

地下低阻異常體模型的總水平導(dǎo)數(shù)試算結(jié)果表明，上述兩種方法計(jì)算的結(jié)果是一致的，由此驗(yàn)證了方法一的正確性。

2.2 源的埋深對(duì)異常電位X方向水平導(dǎo)數(shù)的影響

圖6為一個(gè)地下低阻異常體的模型,模型的參數(shù)如下:A1(100m,100m,30m),A2(100m,100m,50m),A3(100m,100m,100m),A4(100m,100m,130m)為4個(gè)不同埋深的源點(diǎn),M為測(cè)點(diǎn),B為井口位置。井口在模型的中心位置,坐標(biāo)為(100m,100m,0),設(shè)圍巖電阻率ρ0=100Ω·m,異常體電阻率ρn=10Ω·m。異常體頂面埋深5m,距井距離為10m,長(zhǎng)、寬、高均為20m,正方異常體中心坐標(biāo)為(80m,80m,15m)。

圖6 不同點(diǎn)源深度的低阻正方體模型a 三維模型; b 模型XOY平面

圖7是異常電位X方向水平導(dǎo)數(shù)隨點(diǎn)源埋深的變化規(guī)律,從圖中可以看到,在異常體的兩個(gè)邊界處水平導(dǎo)數(shù)的極值位置基本不隨點(diǎn)源的埋深變化,極值之間的寬度保持不變,水平導(dǎo)數(shù)的幅值隨點(diǎn)源的埋深增加而減小。由此可見(jiàn)，點(diǎn)源的埋深對(duì)異常電位的水平導(dǎo)數(shù)的極值相對(duì)位置幾乎沒(méi)有影響,即點(diǎn)源的埋深對(duì)異常體的水平范圍的確定沒(méi)有影響。

圖7 在y=80m處不同深度點(diǎn)源計(jì)算的X方向水平導(dǎo)數(shù)結(jié)果對(duì)比

2.3 圍巖與異常體的電阻率差異對(duì)異常電位X方向水平導(dǎo)數(shù)的影響

采用圖6的低阻異常體模型分別計(jì)算不同電阻率下的X方向水平導(dǎo)數(shù),其中源點(diǎn)坐標(biāo)為(100m,100m,50m)。設(shè)圍巖的電阻率ρ0=100Ω·m,各異常體的電阻率ρn分別等于5,10,20,50,90,200,500Ω·m。圖8為不同電阻率差異下的X方向水平導(dǎo)數(shù)結(jié)果對(duì)比。

從圖8中可以看到,在異常體的兩個(gè)邊界處水平導(dǎo)數(shù)的極值位置基本不隨圍巖與異常體電阻率差異不同而變化,極值之間的寬度保持不變,水平導(dǎo)數(shù)的幅值隨兩者的電阻率差異減小而變小。由此可見(jiàn)，圍巖與異常體的電阻率差異對(duì)異常電位的水平導(dǎo)數(shù)的極值相對(duì)位置幾乎沒(méi)有影響,即圍巖與異常體的電阻率差異對(duì)異常體的水平范圍的確定沒(méi)有影響。

圖8 在y=80m處不同電阻率差異下X方向水平導(dǎo)數(shù)結(jié)果對(duì)比

2.4 偏離鉆井不同位置時(shí)對(duì)異常電位X方向水平導(dǎo)數(shù)的影響

基于圖6的低阻異常體模型,分別計(jì)算異常體右邊界面距離井5,10,15,20,25m時(shí)異常電位X方向的水平導(dǎo)數(shù),其中源點(diǎn)坐標(biāo)為(100m,100m,50m)。設(shè)圍巖電阻率ρ0=100Ω·m,異常體電阻率ρn=10Ω·m,圖9為偏離鉆井不同距離時(shí)計(jì)算X方向水平導(dǎo)數(shù)結(jié)果對(duì)比。

由圖9可以看出,在異常體的兩個(gè)邊界處水平導(dǎo)數(shù)的極值位置隨異常體偏離鉆井位置變化而變化,而極值之間的寬度保持不變,水平導(dǎo)數(shù)的幅值隨著偏離鉆井的距離增大而變小。由此可見(jiàn)，異常體偏離鉆井的距離對(duì)異常電位的水平導(dǎo)數(shù)的極值相對(duì)位置幾乎沒(méi)有影響，即異常體偏離鉆井的距離對(duì)異常體的水平范圍的確定沒(méi)有影響。

2.5 不同埋深的異常體對(duì)異常電位X方向水平導(dǎo)數(shù)的影響

基于圖6的低阻異常體模型，分別計(jì)算異常體頂面埋深為1,5,15,20,50m處的異常電位X方向的水平導(dǎo)數(shù)。設(shè)圍巖電阻率ρ0=100Ω·m,異常體電阻率ρn=10Ω·m,圖10為不同埋深下的異常體的X方向水平導(dǎo)數(shù)結(jié)果對(duì)比,圖11為水平導(dǎo)數(shù)的極值寬度D隨異常體埋深Z的變化規(guī)律。

從圖10與圖11 中可以看到,隨著異常體的埋深增加,異常體X方向的水平導(dǎo)數(shù)的極值位置向外移動(dòng)且幅值減小,異常電位的水平導(dǎo)數(shù)的極值之間的寬度隨埋深增加而增大。這遵從電阻率法的基本規(guī)律,異常體的埋深越淺,異常體的兩個(gè)邊界與極值的位置對(duì)應(yīng)越精確。因此,異常體的埋深對(duì)異常電位的水平導(dǎo)數(shù)的影響較大。

通過(guò)大量的模型實(shí)驗(yàn)研究,得出如下結(jié)論:源的埋深、圍巖與異常體的差異、異常體的橫向位置對(duì)異常電位的水平導(dǎo)數(shù)的極值位置幾乎沒(méi)有影響,而異常體的埋深對(duì)極值的位置,以及極值之間的寬度產(chǎn)生較大的影響,隨著異常體埋深的增加,異常體的水平導(dǎo)數(shù)的極值位置向外移動(dòng)且幅值減小。當(dāng)異常體的埋深較淺時(shí),異常體的水平導(dǎo)數(shù)的極值位置與異常體的邊界有著很好的對(duì)應(yīng)關(guān)系。

圖10 在y=80m處不同埋深的異常體計(jì)算的X方向水平導(dǎo)數(shù)結(jié)果對(duì)比

圖11 水平導(dǎo)數(shù)的極值寬度隨異常體埋深的變化

3 歸一化總水平導(dǎo)數(shù)的定義及應(yīng)用

當(dāng)存在兩個(gè)不同埋深的異常體時(shí),往往深部的異常體的異常響應(yīng)很弱小,特征不明顯,為了使深部與淺部的異常體的響應(yīng)特征及邊界范圍都能得到有效的顯示,引入歸一化水平總導(dǎo)數(shù)進(jìn)行異常體的特征識(shí)別。

歸一化總水平導(dǎo)數(shù)法是利用計(jì)算點(diǎn)總水平導(dǎo)數(shù)與一定窗口范圍內(nèi)總水平導(dǎo)數(shù)最大值的比值進(jìn)行地質(zhì)體邊界的識(shí)別,其具體的表達(dá)式為[18]

(14)

式中:THD表示原始異常u的總水平導(dǎo)數(shù);NTHD(i,j)表示點(diǎn)(i,j)處歸一化總水平導(dǎo)數(shù)值;m,n分別表示窗口的尺寸。

為了對(duì)歸一化總水平導(dǎo)數(shù)對(duì)異常體的邊界識(shí)別的有效性進(jìn)行檢驗(yàn),設(shè)計(jì)了2個(gè)雙異常體模型,源點(diǎn)A的坐標(biāo)為(100m,100m,80m),井口B的坐標(biāo)為(100m,100m,0),2個(gè)異常體的電阻率均為10Ω·m,圍巖的電阻率為100Ω·m。異常體的大小均為20m×20m×20m。①第1組模型:第1個(gè)異常體和第2個(gè)異常體的頂面埋深均為1m;②第2組模型:第1個(gè)異常體頂面埋深為1m,第2個(gè)異常體的頂面埋深為30m。各項(xiàng)模型參數(shù)如圖12與圖13所示。

圖14與圖15分別為兩個(gè)模型的異常電位平面圖。當(dāng)埋深較淺時(shí),異常電位的梯度帶都能顯示出異常體的范圍,但是異常的邊界范圍較為模糊；當(dāng)存在深部異常體時(shí),其異常電位的響應(yīng)很弱小,梯度帶變得更加稀疏,邊界范圍已無(wú)法準(zhǔn)確地識(shí)別出來(lái)。圖16中,紅色虛線框代表真實(shí)異常體的邊界在地面的投影,紅色實(shí)線框?yàn)楸痉椒ㄗR(shí)別出的異常體的邊界范圍。由圖16a與圖16b可以看出,當(dāng)異常體的埋深均較淺時(shí),總水平導(dǎo)數(shù)法與歸一化總水平導(dǎo)數(shù)均能精確地反映出異常體的水平邊界及范圍,異常體的邊界得到了顯著的增強(qiáng),這對(duì)異常的響應(yīng)范圍的圈定很有利;隨著深度的增加,導(dǎo)數(shù)的幅值迅速衰減,總水平導(dǎo)數(shù)法已無(wú)法識(shí)別出深部異常體,如圖16c所示;圖16d為異常電位的歸一化總水平導(dǎo)數(shù)法識(shí)別結(jié)果,可以看到,該法能清晰地識(shí)別出不同深度的異常體,提高了異常體的識(shí)別能力;由于異常體埋深的增加,導(dǎo)致了極值的位置與異常體實(shí)際的邊界位置有一定差異,這與前面分析得出的結(jié)論一致。

圖12 第1組雙異常體模型示意a XOY平面; b XOZ垂直斷面

圖13 第2組雙異常體模型示意a XOY平面; b XOZ垂直斷面

圖14 第1組模型的異常電位平面等值線

圖15 第2組模型的異常電位平面等值線

圖16 雙異常體模型總水平導(dǎo)數(shù)與歸一化總水平導(dǎo)數(shù)的識(shí)別結(jié)果a 第1組模型的總水平導(dǎo)數(shù); b 第1組模型歸一化總水平導(dǎo)數(shù); c 第2組模型的總水平導(dǎo)數(shù); d 第2組模型歸一化總水平導(dǎo)數(shù)

4 結(jié)束語(yǔ)

1) 從點(diǎn)源場(chǎng)滿足的基本關(guān)系式出發(fā),通過(guò)有限元實(shí)現(xiàn)了三維井地電阻率法正演計(jì)算,設(shè)計(jì)了兩個(gè)理論模型計(jì)算出解析解,證明了程序的正確性；并給出了異常電位水平導(dǎo)數(shù)的2種計(jì)算方法。

2) 通過(guò)大量的模型試算,詳細(xì)分析了水平導(dǎo)數(shù)極值與異常體邊界對(duì)應(yīng)關(guān)系的影響因素,異常電位的水平導(dǎo)數(shù)極值不受源的埋深、異常體偏離鉆井的位置、異常體與圍巖的電阻率差異的影響;異常電位的水平導(dǎo)數(shù)幅值會(huì)隨著異常體的埋深增加而減小，且極值間的寬度會(huì)逐漸變寬。

3) 當(dāng)異常體埋深較淺時(shí)，總水平導(dǎo)數(shù)法與歸一化總水平導(dǎo)數(shù)法均能較好地識(shí)別出異常體的邊界；當(dāng)存在兩種不同埋深的異常體時(shí)，深部異常體的異常電位信號(hào)特征不明顯，受淺部異常體的影響，邊界范圍不能同時(shí)顯示出來(lái)。而應(yīng)用歸一化總水平導(dǎo)數(shù)法能夠在一定程度上增強(qiáng)深部異常體的響應(yīng)特征?？梢?jiàn)在一定的地質(zhì)-地球物理?xiàng)l件下，歸一化總水平導(dǎo)數(shù)法有望為準(zhǔn)確圈定異常體的范圍提供一種較為有效的手段，以提高井地電阻率法地質(zhì)解釋的正確率和精確度。

4) 歸一化總水平導(dǎo)數(shù)法尚處于理論研究和數(shù)值模擬階段,還有很多問(wèn)題有待深入研究,如復(fù)雜的不均勻、不規(guī)則的異常體對(duì)異常電位的水平導(dǎo)數(shù)極值的影響;井地電阻率法中的多種裝置的研究(三極、四極);異常電位水平導(dǎo)數(shù)極值隨埋深變化的距離的定量計(jì)算等。

參 考 文 獻(xiàn)

[1] 黃仲良.石油重磁電法勘探[M].北京:石油大學(xué)出版社,1999:1-386

Huang Z L.Grativity magnetism electricity method in petroleum prospecting[M].Beijing:Petroleum University Press,1999:1-386

[2] Scriba H.Computation of the electrical potential in the three-dimension structures[J].Geophysical Prospecting 1981,29(5):790-802

[3] Spitzer K.3-D finite-difference algorithm for DC resistivity modeling using conjugate methods[J].Geophysical Journal International,1995,123(3):903-914

[4] Pridmore D F,Hohmann G W,Ward S H,et al.An investigation of finite-element modeling for electrical and electromagnetic data in three dimensions[J].Geophysics,1981,46(7):1009-1024

[5] 徐凱軍,李桐林.垂直有限線源三維地電場(chǎng)有限差分正演研究[J].吉林大學(xué)學(xué)報(bào):地球科學(xué)版,2006,36(1):137-141

Xu K J,Li T L.The forward modeling of three dimensional geoelectric field of vertical finite line source by finite difference method[J].Journal of Jilin University(Earth Science Edition),2006,36(1):137-141

[6] 王志剛,何展翔,魏文博.井地直流電法三維數(shù)值模擬中若干問(wèn)題研究[J].物探化探計(jì)算技術(shù),2006,28(4):322-327

Wang Z G,He Z X,Wei W B.Study on some problems upon 3d modeling of dc borehole-ground method[J].Computing Techniques for Geophysical and Geochemical Exploration,2006,28(4):322-327

[7] 柯敢攀,黃清華.井地電法的三維正反演研究[J].北京大學(xué)學(xué)報(bào):自然科學(xué)版,2009,45(2):264-272

Ke G P,Huang Q H.3D Forward and inversion problems of borehole-to-surface electrical method[J].Acta Scientiarum Naturalium Universitatis Pekinensis,2009,45(2):264-272

[8] 戴前偉,侯智超,王洪華.井地電位2.5維有限元數(shù)值模擬異常分析[J].物探化探計(jì)算技術(shù),2013,35(4):457-463

Dai Q W,Hou Z C,Wang H H.Anomalous analysis of the 2.5D finite element numerical simulation of borehole-fround electrical method[J].Computing Techniques for Geophysical and Geochemical Exploration,2013,35(4):457-463

[9] 賈正元,張剛.三維井地電阻率法有限差分?jǐn)?shù)值模擬研究[J].現(xiàn)代地質(zhì),2012,26(6):1199-2005

Jia Z Y,Zhang G.Research on finite-difference numerical simulation for 3D borehole-to-ground resistivity survey[J].Geoscience,2012,26(6):1199-2005

[10] 劉海飛,陳德鵬,戴前偉,等.連續(xù)電性介質(zhì)線源井-地電位三維有限元數(shù)值模擬[J].桂林理工大學(xué)學(xué)報(bào),2011,31(1):28-38

Liu H F,Chen D P,Dai Q W,et al.3D FEM modeling of borehole-surface potential with Line current source in semi-underground space of continuous variation of conductivity[J].Journal of Gu ilin University of Technology,2011,31(1):28-38

[11] 張?zhí)靷?用直流電阻率法確定油氣藏邊界的初步試驗(yàn)[J].石油地球物理勘探,1990,25(5):584-593

Zhang T L.An experiment in locating reservoir boundary using direct-current resistivity method[J].Oil Geophysical Prospecting,1990,25(5):584-593

[12] 張?zhí)靷?井中有套管情況下用直流電電阻率法確定油氣藏邊界的試驗(yàn)與研究[J].石油地球物理勘探,1993,28(3):314-324

Zhang T L.Experimental research on reservoir boundary determination using DC resistivity method at a cased borehole[J].Oil Geophysical Prospecting,1993,28(3):314-324

[13] 張?zhí)靷?張白林,聶荔.采用鉆井套管作電極的非無(wú)窮遠(yuǎn)三極剖面法尋找剩余油[J].西南石油學(xué)院院報(bào),1999,21(1):29-34.

Zhang T L,Zhang B L,Nie L.Looking for residual oils with nonboundless trielectrode section method using the drilling casing as electrodes[J].Journal of Southwest Petroleum Institute,1999,21(1):29-34

[14] 張?zhí)靷?張白林,聶荔.用地-井工作方式的三極梯度法尋找小塊油氣藏[J].石油地球物理勘探,1997,32(1):520-531

Zhang T L,Zhang B L,Nie L.Small reservoir discovery using ground borehole trielectrode gradient method[J].Oil Geophysical Prospecting,1997,32(1):520-531

[15] 張?zhí)靷?張白林,聶荔,等.非無(wú)窮遠(yuǎn)三極剖面法在新疆F油田的試驗(yàn)效果[J].西南石油學(xué)院學(xué)報(bào),2001,23(1):5-10

Zhang T L,Zhang B L,Nie L,et al.The test result of nonboundless trielectrode profiling method in F field[J].Journal of Southwest Petroleum Institute,2001,23(1):5-10

[16] 湯井田,張繼峰,馮兵,等.井地電阻率法歧離率確定高阻油氣藏邊界[J].地球物理學(xué)報(bào),2007,50(3):926-931

Tang J T,Zhang J F,Feng B,et al .Determination of borders for resistive oil and gas reservoirs by deviation rate using the hole-to-surface resistivity method[J].Chinese Journal of Geophysics,2007,50(3):926-931

[17] 戴前偉,陳德鵬,陳勇雄.電法勘探中異常響應(yīng)特征的增強(qiáng)算法及其實(shí)現(xiàn)[J].煤田地質(zhì)與勘探,2013,41(3):75-78

Dai Q W,Chen D P,Chen Y X.The enhanced algorithms and its implementation for the abnormal response characteristics in electrical exploration[J].Coal Geology & Exploration,2013,41(3):75-78

[18] Ma G,Li L.Edge detection in potential fields with the normalized total horizontal derivative[J].Computers&Geosciences,2012,41(1):83-97

[19] 吳小平,汪彤彤.利用共軛梯度算法的電阻率三維有限元正演[J].地球物理學(xué)報(bào),2003,46(3):428-432

Wu X P,Wang T T.A 3-D finite-element resistivity forward modeling using conjugate gradient algorithm[J].Chinese Journal of Geophysics,2003,46(3):428-432

[20] 柳建新,蔣鵬飛,童曉忠,等.不完全LU分解預(yù)處理的BICGSTAB算法在大地電磁二維正演模擬中的應(yīng)用[J].中南大學(xué)學(xué)報(bào):自然科學(xué)版,2009,40(2):484-491

Liu J X,Jiang P F,Tong X Z.Application of BICGSTAB algorithm with incomplete LU decomposition preconditioning to two dimensional magnetotelluric forward modeling[J].Journal of Central South University(Science and Technology),2009,40(2):484-491

[21] 王萬(wàn)銀,邱之云,楊永,等.位場(chǎng)邊緣識(shí)別方法研究進(jìn)展[J].地球物理學(xué)進(jìn)展,2010,25(1):196-210

Wang W Y,Qiu Z Y,Yang Y,et al.Some advance in the edge recognition of the potential field[J].Progress in Geophysics,2010,25(1):196-210

[22] 王萬(wàn)銀.位場(chǎng)總水平導(dǎo)數(shù)極值位置位置空間變化規(guī)律研究[J].地球物理學(xué)報(bào),2010,53(9):2258-2270

Wang W Y.Spatial variation law of the extreme value positions of total horizontal derivative for potential field data[J].Chinese Journal of Geophysics,2010,53(9):2258-2270

[23] Cordell L.Gravimetric expression of graben faulting in Santa Fe Country and the Espanola Basin[C]∥New Mexico Geological Society.30thFiled Conference.New Mexico:[s.n.],1979:59-64

[24] 徐世浙.地球物理中的有限單元法[M].北京:科學(xué)出版社,1994:120-122

Xu S Z.The finite element method in geophysics[M].Beijing:Science Press,1994:120-122