謝關(guān)寶，楊震，趙文杰

（1．中石化石油工程技術(shù)研究院，北京100101；2．中國(guó)石油大學(xué)（北京），北京102249；3．中石化勝利油田博士后工作站，山東 東營(yíng)257017）

0 引 言

傳統(tǒng)隨鉆電磁波電阻率測(cè)井儀采用同軸發(fā)射－接收線圈，這種儀器結(jié)構(gòu)和方法不能反映地層方位信息，在地質(zhì)導(dǎo)向應(yīng)用過(guò)程中存在很大局限性。新一代隨鉆電磁波電阻率測(cè)井儀裝有傾斜（或水平）發(fā)射／接收線圈，可提供地層方位信息［1］。2005年斯倫貝謝公司推出第一代隨鉆方位電磁波電阻率測(cè)井儀PeriScope，2006年貝克休斯公司推出隨鉆方位電磁波電阻率測(cè)井儀APR，2007年哈里伯頓公司推出隨鉆方位電磁波電阻率測(cè)井儀 ADR［2－4］，其特點(diǎn)是采用復(fù)合線圈系、多種頻率測(cè)量，提供多種探測(cè)深度的幅度比和相位差。本文以三線圈系為基礎(chǔ)，模擬帶有傾斜接收線圈的方位電磁波電阻率儀器的響應(yīng)規(guī)律，為儀器研發(fā)和地質(zhì)導(dǎo)向應(yīng)用提供理論依據(jù)。

1 數(shù)學(xué)建模

1．1 數(shù)學(xué)物理模型

計(jì)算模型為單發(fā)雙收線圈系，發(fā)射線圈T與儀器軸垂直（見(jiàn)圖1），接收線圈相對(duì)于儀器軸存在一定夾角。假設(shè)地層坐標(biāo)系為xyz，儀器軸相對(duì)于水平地層的傾角為α，方位角為β（相對(duì)于x軸），接收線圈磁矩與儀器軸線所在垂面夾角為γ。儀器的參考坐標(biāo)系為x′y′z′，地層坐標(biāo)系與儀器坐標(biāo)系之間的磁矩轉(zhuǎn)換關(guān)系為

圖1 三線圈系坐標(biāo)示意圖

電磁場(chǎng)數(shù)值模擬中的頻率域麥克斯韋方程［5］可表示為

式中，H、E分別為磁場(chǎng)強(qiáng)度和電場(chǎng)強(qiáng)度；Ji為感應(yīng)電流密度；Js為源電流密度；ω為角頻率；μ0＝4π×10－7（H／m）；σ′為復(fù)電導(dǎo)率張量，σ′＝σ＋iωε；ε為介電常數(shù)。為了得到場(chǎng)源附近精確的電磁場(chǎng)分布，需要精細(xì)的網(wǎng)格劃分，因此把式（2）中的總電場(chǎng)分解為背景場(chǎng)Eb和散射場(chǎng)E′，即

由式（2）至式（5）得

式中，Js＝［（σ′－σ′0I）＋iω（ε－ε0）I］E0；σ′0為背景電導(dǎo)率；ε和ε0分別為介電常數(shù)和背景介電常數(shù)，由于ε和ε0的數(shù)量級(jí)很小，故Js≈（σ′－σ′0I）E0。

在直角坐標(biāo)系中，對(duì)式（6）進(jìn)行有限差分離散。離散時(shí)采用KANE S Yee提出的交錯(cuò)網(wǎng)格［6］（見(jiàn)圖2）。在交錯(cuò)網(wǎng)格中，電場(chǎng)E在單元棱邊上采樣，磁場(chǎng)H在單元面中心采樣。該方法的優(yōu)點(diǎn)在于微分Maxwell方程組成的梯度、旋度、散度關(guān)系式對(duì)其差分模擬依然成立，從而保證了在無(wú)源區(qū)電流守恒，即

使用交錯(cuò)網(wǎng)格對(duì)微分方程進(jìn)行有限差分離散不會(huì)產(chǎn)生贗解［7］。

圖2 交錯(cuò)網(wǎng)格示意圖

考慮地層為有耗媒質(zhì)，電磁波衰減較快，在數(shù)值模擬過(guò)程中須采用有限空間作為求解區(qū)域。在區(qū)域邊界上，電場(chǎng)E滿足Dirichlet邊界條件0，其中?Ω為求解區(qū)域Ω的邊界。

散射電磁場(chǎng)可分解為

式中，σxx（r）、σyy（r）和σzz（r）分別為x、y和z方向的電導(dǎo)率；σb為背景電導(dǎo)率。

對(duì)每一個(gè)節(jié)點(diǎn)的x、y、z分量分別按式（9）至式（11）進(jìn)行差分，所有節(jié)點(diǎn)形成等效方程組

式中，A為超大型正定對(duì)稱稀疏復(fù)數(shù)矩陣，完全存儲(chǔ)需較大內(nèi)存，對(duì)方程進(jìn)行差分計(jì)算發(fā)現(xiàn)，矩陣A每行最多只有13個(gè)非0元素（見(jiàn)圖3），這是應(yīng)用ICCG方法的基礎(chǔ)。在安裝矩陣A的過(guò)程中，只需安裝內(nèi)點(diǎn)，對(duì)于和Dirichlet邊界相鄰的內(nèi)點(diǎn)，在其對(duì)應(yīng)方程中，略去和邊界點(diǎn)對(duì)應(yīng)的項(xiàng)即可。建立地層模型在x、y、z方向上分別有Nx、Ny、Nz個(gè)節(jié)點(diǎn)，施加Dirichlet邊界條件后，共有3NxNyNz－5（NxNy＋NyNz＋NxNz）＋8（Nx＋Ny＋Nz）－12個(gè)未知數(shù)，即形成3NxNyNz－5（NxNy＋NyNz＋NxNz）＋8（Nx＋Ny＋Nz）－12個(gè)方程。

圖3 稀疏矩陣非0元素分布

得到地層電磁場(chǎng)分布后，則可進(jìn)一步得到2個(gè)接收線圈上的感應(yīng)電動(dòng)勢(shì)VR1和VR2，則實(shí)際測(cè)井記錄的幅度比和相位差為［8］

式中，arg表示取相位角。根據(jù)幅度衰減和相位差轉(zhuǎn)換就可以得到隨鉆電磁波電阻率。

1．2 改進(jìn)的ICCG算法

共軛梯度法理論上n步之內(nèi)可以收斂，但一般情況下n次迭代難以得到滿意精度的解，需要反復(fù)迭代，大大限制了共軛梯度的應(yīng)用，需要引入預(yù)條件技術(shù)，預(yù)條件的思路是尋找A的一個(gè)近似分解

盡管ICCG方法實(shí)現(xiàn)了對(duì)稱正定矩陣高效分解，但在計(jì)算過(guò)程中并不能保證不出現(xiàn)對(duì)角元為負(fù)或幅值很小的情況。當(dāng)對(duì)角元為負(fù)時(shí)，會(huì)影響預(yù)處理矩陣的正定性，與矩陣A的近似程度也會(huì)下降；當(dāng)對(duì)角元幅值很小時(shí)，必然引起舍入誤差增大，從而影響分解質(zhì)量和求解效率。本文采用張永杰［9］等提出的對(duì)角元修正方法提高ICCG方法的分解質(zhì)量和求解效率。

令T為預(yù)處理矩陣M＝LLT第i行的1范數(shù)，即

式中，lij、lji分別為下三角矩陣L中的元素。當(dāng)分解過(guò)程中對(duì)角元dii＜T時(shí)，取dii＝T；對(duì)角元dii＞T時(shí)，不進(jìn)行修正。

圖4為分別用改進(jìn)的ICCG和ICCG計(jì)算同一方程的收斂過(guò)程比較，很明顯，改進(jìn)的ICCG的殘差隨迭代次數(shù)的增加迅速下降，而普通ICCG方法中，殘差的震蕩性比較強(qiáng)，需迭代更多的次數(shù)［10］。

圖4 改進(jìn)的ICCG和ICCG收斂過(guò)程比較

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

計(jì)算模型中2個(gè)接收線圈與發(fā)射線圈之間的距離分別為0．6m和0．76m，儀器工作頻率為2MHz，地層的相對(duì)介電系數(shù)εr＝10。地層模型為三層各向同性地層，電阻率分別為1、10Ω·m和1Ω·m，井眼于地層的相對(duì)傾角為80°（見(jiàn)圖5），地層界面位置為z0＝1．5m，z1＝3．5m。圖6、圖7分別模擬了在鉆進(jìn)過(guò)程中，接收線圈處于不同方位角時(shí)的幅度比和相位差（接收線圈與儀器軸的夾角為45°）。該計(jì)算方法的正確性和準(zhǔn)確性在以前發(fā)表的文章［11］中與一維解析解進(jìn)行了對(duì)比，結(jié)果相吻合，證明了該計(jì)算方法準(zhǔn)確可靠。

圖5 井眼軌跡示意圖

圖6 方位電磁波電阻率儀器隨接收線圈方位變化幅度比模擬結(jié)果（θR＝45°）

圖7 方位電磁波電阻率儀器隨接收線圈方位變化相位差模擬結(jié)果（θR＝45°）

通過(guò)模擬結(jié)果可以看出，當(dāng)方位電磁波儀器以較大相對(duì)傾角（80°）穿過(guò)地層界面時(shí)，儀器響應(yīng)同樣出現(xiàn)極化角［12］，該極化角對(duì)地層界面同樣有良好的識(shí)別作用。但與傳統(tǒng)的電磁波電阻率儀器響應(yīng)不同，極化角幅度和位置與接收線圈的方位有關(guān)系，正是這種關(guān)系的存在，采用該種結(jié)構(gòu)的儀器可以實(shí)現(xiàn)對(duì)方位的識(shí)別。接收線圈在界面處的感應(yīng)電動(dòng)勢(shì)隨方位角γ周期性變化。這種響應(yīng)的變化只是出現(xiàn)在地層界面處，而離開(kāi)地層界面一定距離后，儀器的響應(yīng)不再隨接收線圈方位角的變化而變化，這就使儀器在能夠區(qū)別地層方位進(jìn)行精確地質(zhì)導(dǎo)向的時(shí)候同時(shí)具備了地層評(píng)價(jià)的能力。

在儀器旋轉(zhuǎn)過(guò)程中，當(dāng)接收線圈的磁矩與地層界面垂直時(shí)（γ＝0°和γ＝180°），極化角幅度最大，此時(shí)對(duì)地層定位能力最強(qiáng)［12］。當(dāng)接收線圈磁矩與地層界面平行時(shí)（γ＝90°和γ＝270°），儀器響應(yīng)與常規(guī)電磁波電阻率儀器相同。

圖8、圖9為模擬接收線圈與儀器軸夾角為60°的情況下的儀器響應(yīng)。隨著接收線圈傾斜角增加，幅度比和相位差的極化角在地層界面處幅度有所增加，變化更加明顯。這是由于隨著相對(duì)傾斜角的增加，在地層界面附近會(huì)導(dǎo)致磁場(chǎng)的ZX分量增加，這對(duì)于定向應(yīng)該是有利的；但隨著傾斜角的增加，接收線圈處的感應(yīng)電動(dòng)勢(shì)中能反映地層電阻率的電磁場(chǎng)ZZ分量比重降低，導(dǎo)致幅度比和相位差的誤差增大，不利于地層評(píng)價(jià)。因此目前國(guó)外公司商業(yè)應(yīng)用的方位電磁波儀器普遍選擇θR＝45°的方案。

圖8 方位電磁波電阻率儀器隨接收線圈方位變化幅度比模擬結(jié)果（θR＝60°）

圖9 方位電磁波電阻率儀器隨接收線圈方位變化相位差模擬結(jié)果（θR＝60°）

在圖10地層模型中，儀器在鉆進(jìn)過(guò)程中可能會(huì)遇到10Ω·m的目的層，分別位于儀器的上方或下方，常規(guī)隨鉆電磁波電阻率儀器響應(yīng)是相同的，反映不出地層模型的差異。但在實(shí)際鉆井過(guò)程中，如果遇到類似地質(zhì)情況，需要調(diào)整鉆頭鉆進(jìn)方向，以使鉆頭進(jìn)入可能含油的高電阻率地層，常規(guī)隨鉆電磁波電阻率儀器就無(wú)能為力了，但通過(guò)方位電磁波電阻率儀器的響應(yīng)（見(jiàn)圖11），可以明顯區(qū)別判斷2種地層模型，實(shí)現(xiàn)地質(zhì)導(dǎo)向。

圖10 地層模型示意圖

圖11 方位電磁波電阻率響應(yīng)模擬結(jié)果

3 結(jié) 論

（1）隨鉆方位電磁波電阻率測(cè)井儀在鉆進(jìn)過(guò)程中，儀器響應(yīng)會(huì)出現(xiàn)極化角，該現(xiàn)象只出現(xiàn)在地層界面處，且在離開(kāi)界面時(shí)不隨接收線圈方位角的變化而變化，從而結(jié)合方位信息，可以指導(dǎo)鉆頭鉆進(jìn)方向，實(shí)現(xiàn)地質(zhì)導(dǎo)向。

（2）在儀器旋轉(zhuǎn)過(guò)程中當(dāng)接收線圈的磁矩與地層界面垂直時(shí)（γ＝0°和γ＝180°），極化角的幅度最大，此時(shí)對(duì)地層的定位能力最強(qiáng)；當(dāng)接收線圈的磁矩與地層界面平行時(shí)（γ＝90°和γ＝270°），儀器的響應(yīng)與常規(guī)的電磁波電阻率儀器相同。

（3）接收線圈與儀器軸的夾角影響極化角的幅度，夾角太小會(huì)導(dǎo)致界面指示作用減弱，夾角太大則會(huì)影響儀器地層評(píng)價(jià)的精度。

［1］Prensky S．Recent Advances in LWD／MWD and Formation Evaluation［J］．World Oil，2006，227（3）：69－74．

［2］Li Q，Omeragic D，Chou L，et al．New Directional Electromagnetic Tool for Proactive Geosteering and Accurate Formation Evaluation While Drilling［C］∥SPWLA 46th Annual Logging Symposium，Paper UU，June 2005．

［3］Bittar M S，Klein J，Beste R，et al．A New Azimuthal Deep Reading Resistivity Tool for Geosteering and Advanced Formation Evaluation［C］∥SPE Annual Technical Conference and Exhibition Transactions，SPE Pa－per 109971，November 2007．

［4］Zhang Zhiyi，Christophe Gonguet，Vinit Rajan，et al．Diectional LWD Resistivity Tools and Their Business Impacts［C］∥ SPWLA 49th Annual Logging Symposium，Paper FFFF，May 2008．

［5］梁燦彬，秦光戎，梁竹?。姶艑W(xué)［M］．北京：高等教育出版社，1985．

［6］KANE S Yee．Numerical Ssolution of Initial Boundary Problem Involving Maxwell’s Equation in Isotropic Media［C］∥IEEE Transactions Antennas Propagation，1966，14（3）：302－307．

［7］楊震，范宜仁，文藝，等．三維頻率域隨鉆電磁波測(cè)井?dāng)?shù)值模擬［J］．地球物理學(xué)進(jìn)展，2009，24（5）：1833－1838．

［8］楊錦舟，魏寶君，林楠．傾斜線圈隨鉆電磁波電阻率測(cè)量?jī)x器基本原理及其在地質(zhì)導(dǎo)向中的應(yīng)用［J］．中國(guó)石油大學(xué)學(xué)報(bào)：自然科學(xué)版，2009，33（1）：44－49．

［9］張永杰，孫秦，李江海．大型稀疏線性方程組的改進(jìn)ICCG方法［J］．計(jì)算物理，2007，24（5）：581－584．

［10］楊震．非均勻復(fù)雜地層隨鉆電磁波測(cè)井響應(yīng)研究［D］．東營(yíng)：中國(guó)石油大學(xué)（華東），2009．

［11］楊震，劉慶成，岳步江．隨鉆電磁波測(cè)井資料中“極化角”的形成機(jī)理及影響因素分析［J］．測(cè)井技術(shù)，2010，34（3）：210－214．

［12］魏寶君，田坤，張旭．定向電磁波傳播隨鉆測(cè)量基本理論及其在地層界面探測(cè)中的應(yīng)用［J］．地球物理學(xué)報(bào)，2010，53（10）：2507－2515．