魏寶君, 任 臣, 吳康康, 黨 峰

(1.中國石油大學(華東)理學院,山東青島 266580; 2.中國石油集團測井有限公司,陜西西安 710077)

井間電磁成像技術(shù)是油氣藏探測的重要手段之一,該技術(shù)將低頻可移動電磁信號源置于發(fā)射井中并在周圍一口或幾口井中移動探測器測量感應(yīng)信號,然后采用電磁成像技術(shù)將測量信號進行反演以獲得井間地層電導(dǎo)率的分布。該技術(shù)作為地球物理探測的新方法,其探測深度更深、范圍更廣,是對傳統(tǒng)單井電磁測井技術(shù)的突破,大大提高了對油氣藏描述的準確性[1-2]。由于井間實際地層的電參數(shù)分布是三維的,因此最理想的井間電磁成像技術(shù)應(yīng)是通過成像獲得井間地層電導(dǎo)率的三維分布,但由于測量數(shù)據(jù)信息量的限制要想準確獲得完全三維的電導(dǎo)率分布圖像是不可能的。一般在對實際測量數(shù)據(jù)進行反演時采用軸對稱二維模型進行成像,即假設(shè)電導(dǎo)率異常體相對于發(fā)射井呈軸對稱分布,從而獲得電導(dǎo)率在發(fā)射井—接收井剖面的二維近似分布[3-6]。由于成像模型與實際地層存在差別,因而有必要研究采用軸對稱二維井間地層模型對實際地層電磁響應(yīng)進行成像的精確性。筆者采用所開發(fā)的二維積分方程[7]和體積分方程[8-10]分別模擬軸對稱二維井間地層和三維井間地層模型中異常體的響應(yīng),用所開發(fā)的Born迭代反演方法[7,11]通過反演不同模型的模擬數(shù)據(jù)獲得發(fā)射—接收剖面的二維電導(dǎo)率成像,由于金屬套管對井間電磁響應(yīng)的影響可以消除[12-13],在模擬時假設(shè)發(fā)射井和接收井中均無金屬套管;基于模擬結(jié)果分析采用軸對稱二維井間地層模型實現(xiàn)高質(zhì)量成像的可行性及影響因素。

1 基本理論

在進行數(shù)值模擬時忽略位移電流的影響并假設(shè)地層是非磁性的,其磁導(dǎo)率取真空中的值。水平層狀背景地層的電導(dǎo)率設(shè)為σb,不同背景層的σb值可不同。將發(fā)射線圈視為磁偶極子并假設(shè)發(fā)射源隨時間的變化關(guān)系為exp(iωt),其中ω為角頻率。圖1給出了軸對稱二維井間地層模型和三維井間地層模型的簡圖,可采用積分方程進行數(shù)值模擬。

1.1 模擬軸對稱二維井間地層模型響應(yīng)的積分方程

所采用的軸對稱二維井間地層模型見圖1(a),該模型電導(dǎo)率異常體相對于發(fā)射井呈軸對稱分布。在模擬時采用圓柱坐標系,則所有場量均只是徑向坐標ρ和軸向坐標z的函數(shù),與方位角坐標φ無關(guān),計算區(qū)域為發(fā)射—接收剖面。接收井中接收線圈處磁場強度的軸向分量可表示為如下積分方程的形式[7]:

(1)

圖1 井間地層模型簡圖Fig.1 Schematic configuration of cross-hole formation model

在積分區(qū)域D內(nèi),總電場E(r,rT)滿足如下積分方程:

E(r,rT)=Eb(r,rT)+

(2)

其中

Eb(r,rT)=-iωμ0NTITΓ(r,rT).

式中,Eb為電場強度的背景值;Γ(r,r′)為縱向成層背景地層磁矢勢Green函數(shù)。式(2)稱為目標方程,該式是計算積分區(qū)域內(nèi)總電場分布的第二類Fredholm積分方程。式(2)可采用穩(wěn)定型雙共軛梯度(BCGS)方法進行迭代求解,從而得到總電場在積分區(qū)域D內(nèi)的分布。在采用反演成像技術(shù)對井間電磁實際測量數(shù)據(jù)進行處理時一般采用軸對稱二維模型進行成像,從而獲得電導(dǎo)率在發(fā)射井—接收井剖面的二維近似分布。在該成像技術(shù)的每一步迭代過程中均將上述積分方程作為其正演算法。

1.2 模擬三維井間地層模型響應(yīng)的體積分方程

對于實際地層的電導(dǎo)率分布情況,可采用圖1(b)所示的三維井間地層模型。在模擬這種模型時采用直角坐標系,接收井中接收線圈處的磁場強度矢量可表示為如下體積分方程的形式[9]:

E(r′,rT)dr′.

(3)

式中,rT=(xT,yT,zT)、rR=(xR,yR,zR)分別為發(fā)射線圈和接收線圈的位置坐標;V為包含三維電導(dǎo)率異常體的積分區(qū)域;GHJ(rR,r′)為層狀背景地層中r′處的單位電流元在接收點rR處的磁型并矢Green函數(shù);E(r′,rT)為發(fā)射線圈在積分區(qū)域內(nèi)產(chǎn)生的總電場強度矢量;Hb(rR,rT)為發(fā)射源在層狀背景地層中產(chǎn)生的磁場強度矢量。若積分區(qū)域內(nèi)總電場強度矢量E(r′,rT)已知,由該式可獲得接收線圈處的磁場強度的3個分量。

在積分區(qū)域V內(nèi),總電場強度矢量E(r,rT)滿足如下體積分方程:

E(r,rT)=Eb(r,rT)+(-iωμb+

(4)

式中,GAJ(r,r′)為層狀背景地層中r′處的單位電流元在場點r處的磁矢勢并矢Green函數(shù);Eb(r,rT)為發(fā)射源在層狀背景地層中產(chǎn)生的電場強度矢量。式(4)亦可采用穩(wěn)定型雙共軛梯度(BCGS)方法進行迭代求解,從而得到總電場強度矢量在積分區(qū)域V內(nèi)的分布。

設(shè)發(fā)射線圈的磁偶極矩為MT,則式(3)和式(4)中發(fā)射源在背景地層中產(chǎn)生的電場強度矢量和磁場強度矢量可表示為

Eb(r,rT)=GEM(r,rT)·MT,

(5)

Hb(rR,rT)=GHM(rR,rT)·MT.

(6)

式中,GEM和GHM分別為水平層狀背景地層中單位磁偶極子源產(chǎn)生的電型和磁型并矢Green函數(shù)。

1.3 軸對稱二維井間地層模型的成像方法

在成像過程中,假設(shè)所采集到的數(shù)據(jù)點的總數(shù)目為M,成像區(qū)域D被離散為N個相同單元,每個單元內(nèi)的總電場強度和電導(dǎo)率均恒定,則式(1)經(jīng)離散后可表示為如下形式:

f=M·x.

(7)

式中,f為M維復(fù)列矢量,其元素由測量到的散射磁場的軸向分量構(gòu)成,本文中用于成像的測量數(shù)據(jù)既可以根據(jù)三維井間地層模型由體積分方程模擬得到,也可以根據(jù)軸對稱二維井間地層模型由二維積分方程模擬得到;x為由未知電導(dǎo)率差值構(gòu)成的N維實列矢量;M為M×N維復(fù)矩陣,其各元素的具體表達式可由式(1)的離散形式得到。

(8)

式中,iT為發(fā)射源位置序號;iR為接收器位置序號,發(fā)射源位置總數(shù)與接收器位置總數(shù)的乘積即為數(shù)據(jù)點總數(shù)目M;j為成像單元的序號,j=1,…,N。

由于式(8)中總電場強度依賴于x,因而矩陣M亦依賴于x,故式(7)的反演是一個非線性問題;又由于在一般情況下測量數(shù)據(jù)的信息量相對于未知量而言是非常有限的,故式(7)的反演結(jié)果一般是不唯一的。反演式(7)可通過采用Born迭代方法(BIM)實現(xiàn),該方法具有反演過程穩(wěn)定的優(yōu)點[11]。在BIM的第n+1步迭代中可通過使下列正則化成本函數(shù)最小化獲得x的近似值xn+1:

(9)

(10)

式中,上標“+”表示復(fù)共軛轉(zhuǎn)置。求解式(10)采用正則化最小二乘法實現(xiàn)。

2 正演數(shù)值模擬對比

無論是軸對稱二維井間地層模型,還是三維井間地層模型,其背景地層是相同的,發(fā)射源在背景地層中產(chǎn)生的磁場強度的軸向分量也相同,區(qū)別只在于電導(dǎo)率異常體產(chǎn)生的散射場差異有多大,故在通過數(shù)值模擬對這兩種模型的散射場進行對比時假設(shè)背景地層為均勻介質(zhì)。

2.1 異常體沿y方向厚度變化時的影響

圖2為井間地層模型1,所采用的發(fā)射頻率為1 kHz,在井間區(qū)域有兩個相對低電導(dǎo)率異常體,其電導(dǎo)率自上而下分別為0.05和0.01 S·m-1。對軸對稱二維模型,異常體相對于發(fā)射井呈圓周對稱;對三維模型,異常體沿垂直于紙面方向(y向)居中放置,即相對于y=0平面保持對稱,但沿y方向的厚度Δy可不同。在模擬時,固定發(fā)射源位于軸向坐標zT=0處,接收器在接收井中沿井軸方向在-100 m至100 m間移動。圖3給出了該模型二維和三維情況下散射磁感應(yīng)強度Bs的軸向分量隨接收器軸向坐標zR變化時的對比結(jié)果,并給出了不同厚度三維異常體散射場與二維異常體散射場的平均誤差,其中三維異常體沿y向的厚度Δy共取9種數(shù)值。

圖2 地層模型1Fig.2 Formation model 1

圖3 地層模型1數(shù)值模擬結(jié)果Fig.3 Numerical simulation results for formation model 1

由圖3(a)和圖3(e)可以看出,當三維異常體沿y方向的厚度較小時,其產(chǎn)生的散射場的幅度與二維異常體產(chǎn)生的場相比要小得多,平均誤差也較大,這主要是三維異常體體積太小導(dǎo)致的。隨著Δy的增加,三維異常體產(chǎn)生的散射場的幅度逐漸增加,越來越接近軸對稱二維異常體產(chǎn)生的散射場,平均誤差亦逐漸減小,但在該模型中其散射場最大值并沒有達到后者的散射場。當Δy增加到一定數(shù)值時(對本模型約是120 m),三維異常體產(chǎn)生的散射場幅度達到最大值,其后再增加Δy其散射場幅度無明顯改變,平均誤差亦無明顯變化。由圖3(b)和圖3(e)可以看出,該模型兩類異常體散射場相位的差異小于幅度的差異,平均誤差最小。當Δy在很小值(5 m)的基礎(chǔ)上增加時,首先出現(xiàn)三維異常體散射場相位逐漸偏離軸對稱二維異常體散射場相位的情況,平均誤差亦相應(yīng)增加。其后二者再逐漸接近,平均誤差亦隨Δy的增加而逐漸減小。其后再增加Δy,其散射場相位雖稍有增加,但無明顯改變,平均誤差亦無明顯變化。圖3(c)中散射場實部隨Δy的變化規(guī)律與圖3(a)的變化規(guī)律類似,其平均誤差隨Δy的變化規(guī)律亦相似。由圖3(d)和圖3(e)可以看出,兩類異常體散射場虛部的差異較大,平均誤差亦較大,平均誤差隨Δy變化的規(guī)律類似于散射場相位的變化規(guī)律。由于在接收器軸向坐標zR連續(xù)改變時散射場虛部出現(xiàn)符號變化,且其強度遠小于實部分量,故該分量相對誤差較大。

模型2在井間區(qū)域有兩個相對高電導(dǎo)率異常體,其電導(dǎo)率自上而下分別為20 和10 S·m-1,其余參數(shù)及測量方式同模型1。圖4為該模型二維(2D)和三維情況下散射磁感應(yīng)強度軸向分量隨接收器軸向坐標zR變化時的對比結(jié)果,并給出了不同厚度三維異常體散射場與二維異常體散射場的平均誤差,三維異常體沿y向的厚度Δy亦取上述9種數(shù)值。

由圖4(a)可得到與圖3(a)相同的結(jié)論,只是平均誤差最小值對應(yīng)的Δy約為160 m。由圖4(b)和(e)可知,當Δy較小時,三維異常體散射場的相位與二維異常體散射場的相位相比差別較大,平均誤差也較大。隨著Δy的增加,三維異常體散射場的相位逐漸減小,越來越接近二維異常體散射場的相位,平均誤差亦逐漸減小,到Δy為160 m時達到最小,其后變化較緩慢。由圖4(c)和(e)可看出,當Δy較小時,兩類異常體散射場實部的差異較大,平均誤差亦較大。當Δy在很小值(5 m)的基礎(chǔ)上增加時,首先出現(xiàn)三維異常體散射場實部逐漸偏離軸對稱二維異常體散射場實部的情況,平均誤差亦相應(yīng)增加。其后二者再逐漸接近,平均誤差亦隨Δy的增加而逐漸減小。其后再增加Δy其散射場實部雖稍有增加但改變較小。該模型由于在接收器軸向坐標zR連續(xù)改變時散射場實部出現(xiàn)符號變化,故相對誤差較大一些。圖4(d)中散射場虛部隨Δy的變化規(guī)律與圖4(a)的變化規(guī)律更為接近,其平均誤差隨Δy的變化規(guī)律亦相似。

綜合上述兩種模型的計算結(jié)果可以看出,在將三維地層環(huán)境中得到的井間測量數(shù)據(jù)在發(fā)射—接收剖面采用軸對稱二維模型進行成像時其成像結(jié)果的質(zhì)量與異常體沿y方向的厚度Δy有關(guān)。若異常體沿y方向厚度較小,其產(chǎn)生的散射場與二維異常體產(chǎn)生的散射場相比要小,平均誤差較大,成像效果與實際模型相比差別較大。三維異常體沿y方向尺寸越大,誤差越小,成像效果越接近于真實情況,但當Δy增加到一定數(shù)值后再增加厚度時成像效果改善不大。上述結(jié)果也說明,接收器所接收到的電磁信號主要來自于發(fā)射—接收剖面附近地層的貢獻,越偏離該剖面,地層對信號的貢獻越小。

圖4 地層模型2數(shù)值模擬結(jié)果Fig.4 Numerical simulation results for formation model 2

2.2 異常體沿y方向平移的影響

模型3為模型1取三維電導(dǎo)率異常體的厚度Δy為120 m時得到,并假設(shè)三維異常體沿y方向可以平移,其中心點坐標yc可變化。圖5為該模型異常體為二維軸對稱(2D)和三維分布情況下散射磁感應(yīng)強度軸向分量隨接收器軸向坐標zR變化時的對比結(jié)果,并給出了不同的中心點坐標yc對應(yīng)的三維異常體散射場與二維異常體散射場的平均誤差,其中三維異常體中心點坐標yc共取11種數(shù)值。

由圖5(a)和(e)可以看出,當三維電導(dǎo)率異常體沿y方向偏心程度較低即中心點坐標yc較小時(對本例yc可到20 m),其產(chǎn)生的散射場的幅度與異常體居中時散射場的幅度無明顯差別,與二維異常體產(chǎn)生的散射場幅度較為接近,平均誤差也較小。其后隨著偏心程度的增加,三維異常體產(chǎn)生的散射場的幅度逐漸減小,越來越偏離軸對稱二維異常體產(chǎn)生的散射場,平均誤差亦隨之增大。由圖5(b)和(e)可以看出,三維異常體沿y方向偏心程度對其散射場相位的影響遠小于對幅度的影響,在yc<60 m的范圍內(nèi)散射場的相位無明顯改變,平均誤差也較小。只有當yc達到70 m之后相位的變化才逐漸明顯,相位逐漸落后,平均誤差亦隨之增大。圖5(c)中散射場實部隨yc的變化規(guī)律與圖5(a)的變化規(guī)律類似,其平均誤差隨yc的變化規(guī)律亦相似。由圖5(d)和(e)可看出,由于在接收器軸向坐標zR連續(xù)改變時散射場虛部出現(xiàn)符號變化且其強度遠小于實部分量,兩類異常體散射場虛部的平均誤差較大,平均誤差隨yc變化的規(guī)律類似于散射場相位的變化規(guī)律。在yc<40 m的范圍內(nèi)散射場的虛部改變不明顯。

模型4為模型2取三維電導(dǎo)率異常體的厚度Δy為160 m時得到,并假設(shè)三維異常體沿y方向可以平移,其中心點坐標yc可變化。圖6為該模型異常體為二維軸對稱和三維分布情況下散射磁感應(yīng)強度軸向分量隨接收器軸向坐標zR變化時的對比結(jié)果,并給出了不同的中心點坐標yc對應(yīng)的三維異常體散射場與二維異常體散射場的平均誤差,其中三維異常體中心點坐標yc亦取11種數(shù)值。

由圖6(a)和(e)可得到與圖5(a)相似的結(jié)論,對本例當yc<40 m時三維異常體產(chǎn)生的散射場的幅度與二維異常體產(chǎn)生的散射場幅度較為接近。由圖6(b)和(e),該模型三維異常體沿y方向偏心程度對其散射場相位的影響既小于對幅度的影響也小于對模型3中相位的影響。由圖6(c)和圖6(e)可以看出,該模型由于在接收器軸向坐標zR連續(xù)改變時散射場實部出現(xiàn)符號變化,在yc較大時兩類異常體散射場實部的平均誤差較大,平均誤差隨yc變化的規(guī)律類似于散射場相位的變化規(guī)律。在yc小于50 m的范圍內(nèi)散射場的實部改變不明顯。圖5(d)中散射場虛部隨yc的變化規(guī)律與圖5(a)的變化規(guī)律類似,其平均誤差隨yc的變化規(guī)律亦相似。

圖6 地層模型4數(shù)值模擬結(jié)果Fig.6 Numerical simulation results for formation model 4

綜上所述兩種模型的計算結(jié)果可以看出,在將三維地層環(huán)境中得到的井間測量數(shù)據(jù)在發(fā)射—接收剖面采用軸對稱二維模型進行成像時,其成像結(jié)果的質(zhì)量與異常體沿y方向的偏心程度有關(guān)。在厚度Δy足夠大且不變的前提下,若異常體沿y方向偏心程度較低,其產(chǎn)生的散射場與二維異常體產(chǎn)生的散射場相比較為接近,平均誤差較小,成像效果越接近于真實情況。三維異常體沿y方向偏心程度越大,相對誤差越大,成像效果越偏離真實情況。

2.3 異常體沿x方向平移的影響

模型5中井間區(qū)域有兩個相對高電導(dǎo)率異常體,其電導(dǎo)率自上而下分別為20 和10 S·m-1,其余參數(shù)及測量方式同模型1。三維異常體沿y方向的厚度Δy為較小值20 m、中心點坐標yc=0 m。無論是軸對稱二維異常體還是三維異常體均可在發(fā)射—接收剖面沿x方向移動,其中心點坐標xc共取8種數(shù)值。圖7給出了當異常體沿x方向中心點坐標xc取不同的值時軸對稱異常體(2D)和三維異常體(3D)的散射場隨接收器軸向坐標zR的變化關(guān)系。

對比圖7(a)和(b)可以看出,由于三維異常體沿y方向的厚度較小,其散射場幅度的極大值要小于軸對稱異常體產(chǎn)生的散射場幅度的極大值。隨著xc的變化,軸對稱異常體產(chǎn)生的散射場幅度的極大值變化不大,而三維異常體產(chǎn)生的散射場幅度的極大值出現(xiàn)較大變化。隨著xc的增加,三維異常體產(chǎn)生的散射場幅度的極大值逐漸減小,在xc=100 m即井間距的一半時達到最小,其后又逐漸增大。這說明在厚度Δy較小時,中心點坐標xc越接近發(fā)射—接收井中間位置,三維異常體產(chǎn)生的散射場的幅度與二維異常體產(chǎn)生的磁場相比越小,平均誤差也越大。對比圖7(c)和(d)可以看出,由于xc改變導(dǎo)致的軸對稱異常體和三維異常體散射場相位的變化規(guī)律相似,但二者存在一定的相位差。對比圖7(e)和(f)中散射場的實部、圖7(g)和(h)中散射場的虛部,可以得到與散射場幅度相似的結(jié)論。

圖7 地層模型5數(shù)值模擬結(jié)果Fig.7 Numerical simulation results for formation model 5

3 成像結(jié)果對比

基于上述數(shù)值模擬得到的分析結(jié)果,對實際成像效果進行對比。在圖8所示的井間地層模型中有兩個相對高電導(dǎo)率異常體,所采用發(fā)射頻率為1 kHz。在模擬成像所需的測量數(shù)據(jù)時設(shè)發(fā)射源在發(fā)射井中的垂向移動范圍為-50～50 m,移動間隔為5 m,則發(fā)射源位置共為21個。接收器在接收井中的垂向移動范圍、間隔與發(fā)射源相同,接收器位置共計也為21個,故數(shù)據(jù)點的總數(shù)目即信息量數(shù)目為M=441。在模擬時可分別假定待成像的電導(dǎo)率異常體相對于發(fā)射井軸呈軸對稱分布、沿y方向以不同厚度居中或偏心分布,而在利用這些模擬的測量數(shù)據(jù)進行成像時均采用軸對稱模型,即在成像時無論異常體的實際形狀如何均假設(shè)其相對于發(fā)射井呈軸對稱分布，從而實現(xiàn)在發(fā)射—接收剖面的二維成像。在該模型中井間成像范圍設(shè)定為水平方向0～100 m、垂向-50～50 m,每個成像單元的尺寸為5 m×5 m,則成像單元的總數(shù)目為N=400。

圖8 成像模型Fig.8 Imaging model

圖9(a)為實際異常體相對于發(fā)射井軸呈軸對稱時的剖面成像結(jié)果,由于異常體的實際分布方式與成像所采用的軸對稱模型是相同的,故成像結(jié)果與原模型最為接近,無論是分布范圍還是異常體的電導(dǎo)率數(shù)值均具有較好的成像效果。圖9(b)為實際異常體為三維異常體、沿y方向居中且厚度較大(yc=0 m、Δy=80 m)時的成像結(jié)果。根據(jù)前面分析,此時三維異常體產(chǎn)生的散射場與軸對稱二維異常體產(chǎn)生的散射場雖有差別但較為接近,故從該圖也能得到較為理想的成像結(jié)果,但無論是分布范圍還是異常體的電導(dǎo)率數(shù)值均不如圖9(a)精確。圖9(c)為三維異常體沿y方向居中且厚度較小(yc=0 m、Δy=40 m)時的成像結(jié)果,圖9(d)為三維異常體沿y方向偏心(yc=40 m、Δy=80 m)時的成像結(jié)果。這兩種情況下三維異常體產(chǎn)生的散射場與軸對稱二維異常體產(chǎn)生的散射場差別較大,成像質(zhì)量均不如圖9(b),與前面分析結(jié)論也是一致的。由圖9的成像結(jié)果可以看出,只要異常體沿y方向厚度足夠大、其中心點距離發(fā)射—接收剖面較近,是可以采用軸對稱二維模型實現(xiàn)較高質(zhì)量成像的。異常體沿y方向厚度越小、偏心程度越大,成像效果越差。

圖9 成像結(jié)果對比Fig.9 Comparison of imaging results

4 結(jié) 論

(1)將三維地層環(huán)境中得到的井間測量數(shù)據(jù)在發(fā)射—接收剖面采用軸對稱二維模型進行成像時，其成像質(zhì)量與異常體沿垂直于發(fā)射—接收剖面方向的厚度有關(guān)。若異常體厚度較小,其產(chǎn)生的散射場與二維異常體產(chǎn)生的散射場相比要小,成像效果與實際模型相差較大。厚度越大,成像效果越接近于真實情況,但當厚度增加到一定數(shù)值后再增加厚度時成像效果變化不大。

(2)將三維地層環(huán)境中得到的井間測量數(shù)據(jù)在發(fā)射—接收剖面采用軸對稱二維模型進行成像時,其成像質(zhì)量與異常體沿垂直于發(fā)射—接收剖面方向的偏心程度有關(guān)。在厚度足夠大的前提下,若異常體偏心程度較低,其產(chǎn)生的散射場與二維異常體產(chǎn)生的散射場相比較為接近,成像效果越接近于真實情況。偏心程度越大,成像效果越偏離真實情況。

(3)三維異常體在發(fā)射—接收剖面水平方向中心點坐標越接近于發(fā)射—接收井中間位置,其產(chǎn)生的散射場的信號強度越小。

(4)只要異常體沿垂直于發(fā)射—接收剖面方向的厚度足夠大,其中心點距離發(fā)射—接收剖面較近,是可以采用軸對稱二維模型實現(xiàn)高質(zhì)量成像的。異常體厚度越小,偏心程度越高,成像效果越差。