随钻方位电阻率测井仪器响应数值模拟

李 辉, 鄢志丹, 刘长波, 姜一波

(1.常州工学院电气与光电工程学院,江苏常州 213032; 2.中国石油大学(华东)信息与控制工程学院,山东青岛 266580; 3.中国石油集团渤海钻探定向井技术服务公司,天津 300280)

大斜度井或水平井钻井作业过程中,为使井眼轨迹保持在期望的储层内,在储层边界处需要随时调整钻头,迫切需要一种能开展实时地质导向及复杂油气层钻进的随钻测井系统[1]。传统的随钻电阻率测井仪器由于采用轴向缠绕发射和接收天线,从而不具备方位探测能力,因此无法独立确定地层方位走向和电阻率各向异性,即使利用多条响应曲线也很难消除电阻率各向异性的多解性[2]。近钻头测量仪器可以为随钻电阻率测井提供地质导向信息,但是在旋转测量过程中容易出现“晕转”,难以准确测量工程参数和地质参数[3]。旋转导向钻井仪器虽然可以提高造斜率,缩短钻井时间,但是将其与随钻电阻率测井仪器相结合必然导致仪器串过长,容易造成卡钻,并且仪器使用费用较高、维修困难,至今难以普及[4]。解决上述问题的最佳方法就是开展对随钻方位电阻率测井仪器及其测井资料处理与解释方法的研究。近年来,国外测井巨头斯伦贝谢、哈里伯顿、贝克休斯相继成功研制出了具有自主知识产权的随钻方位电阻率测井仪器并已投入商业化应用。在1993年,斯伦贝谢公司就推出了IDEAL地质导向系统,它是第一代集地质导向与地层评价为一体的随钻测井仪器。2005年,斯伦贝谢公司研制出具有未钻地层预测能力的随钻方位电阻率仪器PeriScope[5];1994年,哈里伯顿公司推出了随钻电阻率测井仪器FEWD[6];随后在2007年又成功研制出具有地质导向能力的随钻方位电阻率仪器ADR;1999年,贝克休斯公司研制出随钻电阻率测井仪器OnTrak,在2006年又推出了具有方位成像能力的方位电磁波电阻率仪器AziTrak[7]。以上研究成果,不仅提高了电磁波传播随钻测量技术的方位探测能力和信息量,而且提高了该技术对未钻地层界面位置及方位的钻前预测能力并初步实现了对复杂地层的高精度对比评价[8]。国内对随钻方位电阻率测井仪器的研究尚处于起步阶段,主要围绕复杂介质中电磁场数值模拟问题。2010年,魏宝君等[9]利用并矢Green函数法,对水平层状各向异性介质中定向电磁波传播随钻测量的幅度衰减和相位移进行了计算和分析,讨论了天线倾角及地层电阻率对比度对仪器测井响应的影响。2013年,杨震等[10]利用三维有限差分法,对随钻电阻率仪器激发的电磁场随仪器方位变化的规律,以及方位电阻率成像原理进行了分析和探讨。2015年,杨震等[11]又模拟了方位电阻率仪器的界面响应特征及边界走向与仪器工具面角的关系,分析了不同地层及仪器参数对随钻方位电阻率仪器探测深度以及界面距离反演结果的影响。2016年,杨震等[12]对带有横向天线的方位电阻率仪器进行快速响应数值模拟并对补偿方式进行了分析,提出通过对称发射补偿可以消除电阻率各向异性的影响。同年,刘乃震等[13]提出将“交联天线”应用于随钻电磁波传播方位电阻率仪器,通过分析交联天线的电压在仪器旋转过程中的变化规律实现对仪器所在地层的电阻率的计算,并通过下井实验验证了所提方案的有效性。2016年,王浩森等[14]采用三维有限体积法考察了钻铤、线圈倾斜角度以及地层各向异性等参数对仪器响应的影响。以上研究成果对于提升中国随钻测井技术水平以及打破国外技术封锁和垄断都将具有重要的现实意义。 随钻方位电阻率测井资料的精确性直接影响到测井解释与反演评价的可靠性,同时依据测量曲线以及电阻率成像进行实时地质导向又对仪器的测量精度提出了新的要求。为了提高正演数值模拟的精确性,笔者利用目标导向自适应有限元建立数值模拟模型,对随钻方位电阻率测井仪器响应进行计算和分析。

1 原理和方法

1.1 测量原理

传统的随钻电阻率测井仪器因补偿方式的不同,主要分为对称型和非对称型结构。但无论是对称型或非对称型结构,其天线缠绕方式都与钻铤轴向垂直,因此测量结果与仪器所处地层方位无关。随钻方位电阻率测井仪器采用倾斜缠绕接收天线,使测量结果与仪器工具面产生一定的关联,当仪器轴向旋转时测量到的电磁波信号就会出现定向幅值衰减和相位漂移,再利用地质导向信号的峰值和极性以及方位电阻率测井曲线实现储层边界探测。

设z轴为储层界面法向,在各向异性介质中接收天线处的磁场可以表示为

(1)

式中,Mx、My、Mz为由发射天线产生的磁场在x、y、z三个方向上的分量;Hx、Hy、Hz为接收天线在x、y、z三个方向所接收到的磁场分量;Cxx、Cxy、Cxz为x方向单位磁偶极子产生的三个磁场分量;Cyx、Cyy、Cyz为y方向单位磁偶极子产生的三个磁场分量;Czx、Czy、Czz为z方向单位磁偶极子产生的三个磁场分量。为了计算天线倾斜时的仪器响应,需要考虑发射天线和接收天线之间的信号耦合[1]。假设接收天线磁矩方向与仪器轴向(z轴)之间的夹角为θ,旋转过程中仪器高边与正北方向的夹角为方位角φ,如图1所示。接收天线处的感应电势会随方位角φ的变化而改变,因此接收天线处的感应电势如下式所示:

[(Czx+Cxz)cosφ+(Czy+Cyz)sinφ]=A1(φ)+A2(φ)+A3.

(2)

其中

(3)

式中,V(φ)为感应电势;θT为发射天线磁矩方向与仪器轴向(z轴)之间的夹角;θR为接收天线磁矩方向与仪器轴向(z轴)之间的夹角;φ为方位角;A1和A2为余弦信号;A3为常值信号。

图1 仪器坐标系与方位角Fig.1 Coordinate system and azimuth of tool

式(3)中,A1(φ)和A2(φ)可以表示为

(4)

式中,Adr、Adi、Asr、Asi为信号幅值;φdr、φdi、φsr、φsi为信号相位。

其中Adr、Adi、Asr、Asi表示为

(5)

对于传统的随钻电阻率测井仪器,由于所有天线的磁矩都沿仪器轴向,测量到的电磁波信号与方位角无关,因而不具备方位识别能力。当接收天线发生倾斜,天线磁矩就与方位角发生关联,测量到的电阻率就会随方位角的变化而改变。在实际应用中,利用A/D转换器对电磁波信号进行采样和量化,再通过数字信号处理器进行处理可以快速实现电磁波信号幅度比的测量。相位差数字化测量方法如过零检测法和相关函数分析法等[15],容易受到谐波和噪声的干扰,测量误差较大,而且还容易受到信号幅度及仪器硬件电路的影响,因此本文中采用幅度比作为主要信号,对随钻方位电阻率测井仪器响应进行研究。

1.2 数值方法

接收天线感应电动势是由交变电流产生,因此求解的电磁场为时谐场。根据电磁场原理,麦克斯韦方程可表示为

(6)

式中,H和E分别为磁场和电场;Mimp为外加磁场;ε、μ和σ分别为介电常数、磁导率和介质电导率;ρ为电荷密度;ω(ω≠0)为发射源角频率;Jimp为外加电流。将式(6)中法拉第定律代入安培定律,可得时谐麦克斯韦方程在求解域Ω的波动方程,如下式所示:

(7)

在矢量场中,任意旋度场一定是无散场,存在

(8)

那么,

(9)

将式(6)中高斯电场定律和式(8)代入式(9)可得

(10)

则式(10)为电场连续方程。

高频情况下,在H(curl)空间对电场进行求解,能够使求解过程中不会出现伪解。H(curl)空间定义为

H(curl,Ω)={E∈[L2(Ω)]2;×E∈L2(Ω)}.

(11)

假设求解域为Ω,且求解域边界为理想导体边界[16],则

(E,F)H(curl)=(E,F)L2(Ω)+(×E,×F)L2(Ω).

(12)

式中,F为任意测试函数。由于求解域边界是理想导体边界,那么求解电场E所用的矢量求解空间H可定义为

H={E∈H(curl;Ω);n×E=0,Γ1}.

(13)

将封闭区域边界条件作为约束方程代入波动方程式(7),依据法拉第定律,并引入任意测试函数F,在Ω域积分可得

(14)

再依据安培定律,

(15)

将式(15)代入式(14)可得电场变分方程:

(16)

同时,依据安培定律,并引入任意测试函数F,在Ω域积分可得

(17)

再依据法拉第定律,

(18)

将式(18)代入式(17)可得关于磁场的有限元变分方程:

(19)

式中,F∈H={F∈H(curl,Ω):(n×F)|Γ1=0};Γ1为理想导体边界条件;Γ2为发射线圈边界条件,Γ3为非理想导体边界条件[16];n为单位法线矢量。

联立式(16)和式(19)可得自适应hp有限元双线性函数表达式:

(20)

利用伽辽金法可将求解空间离散为若干个有限维子空间,并对其分别进行求解组成各自的刚度方程,将各个刚度方程组集成总刚度矩阵。通过施加不同边界条件和初始条件,再利用自适应hp有限元求解器对总刚度矩阵进行求解,最终可以得到在求解域Ω内接收天线附近的电磁场强度。

目标导向自适应有限元可以根据数值模型误差指示自动判断高误差区域,进行自适应调整,提高数值计算结果的准确性。该算法计算精度高,收敛速度快,可以使形状函数呈指数收敛速率快速准确地逼近真值。与其他数值算法相比,可以减少几个数量级的计算量,在实际工程问题的分析和计算中具有十分重要的意义。篇幅所限,详细算法见参考文献[16]～[18]。

双网格剖分法是目标导向自适应有限元的核心,其目的是在粗网格上求解方程组,利用求解细网格上方程的方法来求解粗网格方程,从而实现多网格法求解。通过递归的实现这一步骤,使每层网格的尺寸都是前一层网格的两倍,直到网格上的未知量个数小到求解代数方程直接法的效率与迭代法的效率一样高为止。双网格剖分法可归纳为:

(1)对于初始单元,首先采用h型细化方式对其进行细化,得到子单元K。同时再利用p型细化方式对子单元K中形状函数的阶次进行提升。

(2)对K进行p型细化后再采用h型细化方式,将其细化为4个子单元K1、K2、K3和K4。 然后对子单元K1再进行p型细化,从而使各子单元中形状函数的阶次对应发生改变。

(3)每种细化方式都需要对参考解进行计算,然后根据计算结果来计算单元区域误差,计算所得误差最小时细化方案为最优。

按照双网格剖分法对目标导向自适应有限元精确性进行验证。在初始网格中,利用目标导向自适应有限元对目标函数f(x)=sin (7πx2)ex,x∈(0,1) 进行数值逼近,主要步骤如图2所示。

图2 自适应数值逼近过程Fig.2 Adaptive numerical approximation process

目标导向自适应有限元可以使初始网格中的形状函数自适应地逼近目标函数,目标导向特性表现为对于计算误差较大的区域,算法通过自适应地调整形状函数的阶次和网格加密程度来达到提高计算精度的目的。在同等计算精度条件下,其计算时间、迭代次数、自由度以及全局误差等性能指标明显优于标准有限元、高阶有限元和矢量有限元等数值算法[19]。

2 数值模拟

随钻方位电阻率测井仪器基本结构如图3所示。其中,T1、T2、T3、T4、T5和T6为发射天线,R1、R2和R3为接收天线。

图3 随钻方位电阻率测井仪器基本结构Fig.3 Basic structure of azimuthal resistivity LWD tool

本文中数值模拟模型如图4所示。其中地层模型为水平层状各项同性均匀介质,上围岩和下围岩为低阻页岩层,上、下围岩垂深Lv1=Lv3=10 m,上、下围岩电导率σ1=σ3=1 S/m;目的层为高阻砂岩层,垂深Lv2=10 m,电导率σ2=0.1 S/m;地层倾角为α;接收天线R3倾角θ=45°,天线组合T6-R3源距为L1=40.64 cm,天线组合T5-R3源距为L2=81.28 cm,天线组合T4-R3源距为L3=121.92 cm;发射天线和接收天线带有磁缓冲器,磁缓冲器电阻率为1×104Ω·m,钻铤电阻率为1×10-6Ω·m,井眼环空充满钻井液,钻井液电阻率为0.2 Ω·m。

图4 随钻方位电阻率测井数值模型Fig.4 Numerical model of azimuthal resistivity LWD tool

2.1 磁场测量

随钻方位电阻率测井仪器发射频率f=2 MHz,采用T4-R3天线组合,对仪器接收天线R3附近的磁场进行测量,如图5所示。

仪器旋转过程中,测量到的磁场与方位角之间的关系如图6所示。

图5 磁场测量Fig.5 Measurement process of magnetic field

随钻方位电阻率测井仪器采用倾斜缠绕接收天线,旋转测量过程中所测量到的磁场与方位角产生一定的关联,当仪器高边在工具面上对应方向为正北,则磁场与方位角之间的变化趋势如图6所示。在x方向,磁场分量Cxx和Cxy的实部与虚部在第一、三象限呈递增趋势,在第二、四象限呈递减趋势;磁场分量Cxz的实部表现为在第一、二象限呈递减趋势,在第三、四象限呈递增趋势;Cxz的虚部所呈现出的变化趋势与实部完全相反。在y、z方向上可以测量到6个磁场分量,其中磁场分量Czz为常值信号,与方位角无法产生关联。利用磁场与仪器工具面之间的定量关系,可以得到不同方位角所对应的地质导向信号。

2.2 储层边界识别

仪器发射频率f=2 MHz,地层无倾斜,采用T4-R3天线组合。假设仪器处于上围岩且即将进入目的地层,当仪器从上围岩进入目的地层并到达下围岩,再从下围岩进入目的地层并到达上围岩,数值计算得到的地质导向信号如图7所示,地层电阻率曲线如图8所示。

通过图7可以看出,在地层边界处仪器所测量到的地质导向信号出现明显的极化异常。当仪器远离目的地层上边界并穿过下边界时,地质导向信号出现极性翻转,从正极性变为负极性,表明仪器已从上围岩穿过目的地层进入下围岩。当仪器从下围岩穿过目的地层到达上围岩时,地质导向信号极性同样发生两次翻转。此外,地质导向信号极值越大表明仪器越接近地层边界,通过地质导向信号的极值和极性可以准确判断仪器所处地层位置以及钻头与目的地层边界的距离。

2.3 天线源距

仪器天线组合为T4-R3、T5-R3和T6-R3,地层无倾斜,发射频率f=2 MHz,考察不同天线源距对地质导向信号及方位电阻率的影响。

通过图9(a)可以看出,随着发射与接收天线源距的增大,仪器的探测深度不断增强。源距较小时,地质导向信号极化不明显,幅值衰减较慢,探测距离较浅,对地层边界不敏感。源距较大时,地质导向信号幅度衰减加快,极性翻转迅速,能快速识别地层边界。图9(b)给出了不同源距时数值模型计算出的方位电阻率曲线,计算结果表明仪器旋转测量过程中,方位电阻率随源距的增大而增大。

2.4 发射频率

仪器发射频率f分别为125 kHz、400 kHz和2 MHz,地层无倾斜,采用天线组合为T4-R3,考查发射频率改变对地质导向信号及方位电阻率的影响。

通过图10(a)可以看出,仪器发射频率越高,在地层边界处地质导向信号极值越大,反映出钻头离地层边界越近。随着发射频率的降低,地质导向信号幅值衰减变缓,分辨率降低。图10(b)给出了不同发射频率时数值模型计算出的方位电阻率曲线,计算结果表明仪器旋转测量过程中,方位电阻率随发射频率的增大而增大。

图6 磁场分量与方位角Fig.6 Magnetic field component and azimuthal angle

图7 地质导向信号Fig.7 Geo-signal response

图8 地层电阻率曲线Fig.8 Resistivity curve of formations

图9 不同天线组合时地质导向信号与方位电阻率Fig.9 Geo-signal and azimuthal resistivity with different coil groups

图10 不同发射频率时地质导向信号与方位电阻率Fig.10 Geo-signal and azimuthal resistivity with different transmitting frequencies

2.5 地层倾角

地层倾角α分别为15°、30°、45°和60°,发射频率f=2 MHz,采用T4-R3天线组合,考察地层倾角改变对地质导向信号以及方位电阻率的影响。

通过图11(a)可以看出,地层倾角越大,所测量到的地质导向信号极值越大,极性变化越快。随着地层倾角的减小,地质导向信号幅值衰减变慢,信号的分辨率逐渐降低,对地层边界的敏感性有所减弱,信号的极性和极值变化缓慢。图11(b)给出了不同地层倾角时数值模型计算出的方位电阻率曲线。计算结果表明仪器旋转测量过程中,方位电阻率随地层倾角的增大出现剧烈变化,当地层倾角为60°时,在正南方向(低边工具面)方位电阻率出现极小值。

图11 不同地层倾角时地质导向信号与方位电阻率Fig.11 Geo-signal and azimuthal resistivity with different tilt angles of formations

2.6 围岩电导率

当上围岩电导率σ1分别为1、0.3、0.2、0.15和0.1 S/m,采用T4-R3天线组合,发射频率f=2 MHz,考查围岩电导率改变对接收天线R3感应电势的影响。

通过图12可以看出,当接收天线与围岩边界距离不变时,随着围岩电导率的增大,接收天线感应电势逐渐增大。当接收天线与围岩边界距离逐渐增大时,接收天线感应电势逐渐减弱。围岩电导率与边界距离改变对接收天线所接收到的定向电磁波信号的影响如表1所示。当接收天线与围岩边界距离不变时,随着围岩电导率的增大,接收天线所接收到的电磁波信号实部逐渐增大而虚部逐渐减小并呈现出规律性变化,围岩电导率的改变会影响定向电磁波信号幅值比和相位差。

图12 接收天线感应电势与围岩电导率Fig.12 Relationship of induced potential and surrounding rock conductivity of receiving antenna

边界距离/mσ1=1 S/m实部/10-4虚部/10-6σ1=0.35 S/m实部/10-4虚部/10-6σ1=0.2 S/m实部/10-4虚部/10-6σ1=0.15 S/m实部/10-4虚部/10-6σ1=0.1 S/m实部/10-4虚部/10-7 0.52.03-4.561.93-4.061.82-2.511.72-1.691.62-8.43 1.02.02-3.581.95-3.081.84-1.831.74-1.211.64-5.96 1.51.99-3.271.92-2.731.81-1.611.71-1.031.61-5.00 2.01.95-3.521.87-3.021.76-1.781.66-1.151.56-5.46 2.51.90-4.421.81-3.921.71-2.371.61-1.571.53-7.72 3.01.82-5.991.75-5.531.64-3.431.54-2.301.45-11.50 3.51.71-7.511.65-6.921.54-4.981.44-3.381.32-18.90 4.01.54-95.601.48-8.861.40-6.631.29-4.871.17-29.90 4.51.41-12.301.32-11.301.22-8.481.12-6.780.999-45.00 5.01.21-16.101.11-14.101.04-11.000.928-8.960.820-69.70

2.7 接收天线倾角

接收天线倾角θ分别为30°、60°和90°,发射频率f=2 MHz,采用T4-R3天线组合,考察当目的地层存在地磁异常时,接收天线的测井响应。

通过图13可以看出,随着接收天线倾角的增大,接收天线所测量到的地质导向信号极值不断增大。如图13(a)所示,当地层介质中无地磁异常时,随着接收天线倾角的增加,地质导向信号的极值不断增大,利用地质导向信号可以准确地指示出目的地层边界。如图13(b)所示,当目的地层存在地磁异常时,随着接收天线倾角的增加,地质导向信号出现漂移。同时在目的地层边界处,地质导向信号的强度有所减弱,对地层边界的识别精度将会降低。

图13 接收天线不同倾角时地质导向信号响应Fig.13 Geo-signal response with different tilt angles of receiving coil

3 结 论

(1)随钻方位电阻率测井仪器在旋转测量过程中所测量到的磁场与仪器工具面产生一定的关联。利用磁场分量与仪器工具面之间的定量关系,能够得到与方位电阻率相关的地质导向信息。

(2)发射频率与源距变化会影响随钻地质导向信号的强弱及衰减速率,会影响方位电阻率曲线的分辨率及其对地层边界的敏感性。利用随钻地质导向信号的极值和极性,结合井周电阻率成像可以为实时钻井提供准确的地质导向信息。

(3)地层倾角与围岩电导率的变化会引起随钻地质导向信号产生畸变,同时方位电阻率随地层倾角的增大会出现不规律的变化,对准确判断地层边界及仪器所处地层方位造成影响。

(4)电阻率成像的准确性将影响储层边界与裂缝边缘及走向识别的准确性,而井周电阻率成像过程容易受到地层各向异性及环境噪声的干扰,因此电阻率成像过程中需要进行补偿和去噪。下一步将对井周电阻率成像去噪和细节特征提取方法进行研究。