复杂地质层中电磁波测井响应特性的数值研究

刘国胜，杨海东，汤健超

(1.广东工业大学 信息管理工程系，广东 广州，510520；2.广东工业大学 机电工程学院，广东 广州，510006；3.华南师范大学 计算机学院，广东 广州，510631)

随着全球能源供需关系日趋紧张，世界各国对潜在石油、天然气资源的争夺越来越激烈，这同时促使了钻井设备与技术不断发展更新[1-2]。目前普遍采用的钻进技术——地质导向钻进技术，它可根据地下情况随时调整探测方向，目前这项技术已得到普遍应用，它大大降低开发成本，提高勘探成功率[3]。在钻井工作环境中，电磁波测井工具主要提供周围地质层介质电导率(电阻)的信息，进而判断地质构成。目前对电磁波测井工具的性能分析主要来源于大量工程试验，如何对电磁波测井工具进行有效建模成为科学研究的一个热点问题[4-10]。钻井工具中的电磁波探测部分由一系列极子天线和环形天线组成，各天线在钻轴上的位置分布和朝向存在差异。典型随钻电磁波测井工具的数学模型由1个(环形)发射天线和2个(环形)接收天线组成，其通用尺寸见图1(一般LWD工具的尺寸采用in为单位，1 in=25.4 mm)。图1中的发射天线附加正弦电压，被激励的电磁波与周围地质层作用后在接收天线上形成响应电压，通过比较2个接收天线中电压振幅和相位差异来判断周围地质的电导率[11]。目前对 LWD工具的分析方法包括解析法、伪解析法和纯数值方法。解析法将图1中的环形天线近似为磁极子天线，进而用势函数积分对问题进行求解[9,12-13]，这种方法简单易用，便于定量分析，但其忽略了天线尺寸和钻轴的影响，可能导致结果不准确；伪解析法作为一种半解析半数值的方法，是分析 LWD工具非常有效的方法之一[5,7-8]，但其无法处理复杂地质介质和地质结构，适用范围受到限制；典型的数值方法包括时域有限差分法[4,10]、有限元法和有限体积法[6]等，这些方法可用于复杂钻井环境的仿真，如偏心几何结构和各向异性地质层，但其计算量较大，不利于工程中处理大尺寸问题和反复多次模拟。Lovell等[7-8,14-15]对随钻电磁波测井模型中的伪解析法进行大量研究，对钻井电磁波测井中的仿真模型做出了较大贡献。在实际环境中，由于地质层受到长期重力作用，地质层介质通常在垂直和水平方向表现出不同的电磁特性，其垂直和水平方向的电导率存在差异，形成单轴各向异性介质[10]。伪解析法具有精度高、速度快的优点，将伪解析法推广到单轴各向异性介质的情况无疑具有重要意义。在此，本文作者讨论单轴各向异性介质中的伪解析法，给出在此环境中处理偏心多层结构的计算公式，并分析在钻轴倾斜挖掘时，如何对钻探方向进行最优选择，模拟在复杂地质环境中LWD工具的数值特性，并分析天线垂直距离、钻轴偏心等因素对LWD工具测量结果的影响。

图1 LWD结构模型(单位：mm)Fig.1 Configuration of LWD

1 单轴各向异性介质中的伪解析法

在单轴各向异性介质中，Maxwell方程可写为如下向量形式[11]：

式中：E为电场强度；H为磁场强度；μ为磁导率；ε为复介电张量。在地质层中，岩层受到重力挤压，通常在垂直和水平方向上表现出不同的电磁特性，因此，讨论ε其中σs和σz分别表示地质层介质在轴向和水平截面的电导率分量)。为简单起见，记从E和H中分离出 eikzz，也就是在式(1)中令 ∂ ∂z= ikz，并且将方程中算子和场量进行分解，方程(1)可分解为：

将式(4)和(5)代入式(2)和(3)可求得Es和Hs的表达式，将结果代回式(2)和(5)得到z分量满足的方程如下：

方程(6)的解可写成如下一般形式：

Gν=Jν或分别是第一类 Bessel函数和第一类Hankel函数。将式(7)代入到式(2)和(3)可求得φ分量的解。

由于在LWD工具中(如图1所示)，环形天线在水平截面上的投影为圆形，因此，天线周围电磁场的径向分量Eρ和Hρ对天线中的电压无贡献。在地质层边界处，电磁场的切向分量(φ和z分量)是连续的。图2所示为边界处的反射和透射。图2中边界Ri上的透射和反射矩阵分别定义为和，其可通过边界处的连续性条件推导得到[5,11,14]。

图2 边界处的反射和透射Fig.2 Reflection and transmission on boundary

2 复杂地质结构的处理

当受机械牵引或重力作用时，钻轴可能会偏离井孔中心。图3所示为两层偏心结构示意图。考虑到图3中2个偏心圆(层)所构成的几何结构，其中每层由轴各项异性介质填充。内、外圆的圆心分别对应圆柱坐标系(ρ′，φ′)和坐标系(ρ，φ)的原点。坐标系(ρ，φ)可以看成是由坐标系(ρ′，φ′)平移(d，φE)得到。考察两圆弧之间任意一点处的场强，其谱分量Ez和Hz可根据式(7)在坐标系(ρ′，φ′)和(ρ，φ)中分别展开[11,14]。

现在引入1个虚拟边界(见图3中的虚线圆)，此圆与内圆具有相同的圆心。值得注意的是，此处选择的虚拟圆边界与文献[5，14]中不同，此处虚拟圆边界的引入将会极大简化处理多层偏心结构，会为多层偏心结构推导出1个一致性的迭代解。当由内层入射的外行波被虚拟边界“反射”后，不同模ν′的波会相互影响，其结果等效于被外层边界反射的结果，因此，虚拟边界处的等效反射关系可表示为：

图3 两层偏心结构示意图Fig.3 Schematic diagram of two-layer eccentric structure

在坐标系(ρ，φ)中，外层边界的反射关系可用各个整数模ν独立表示，如

将式(8)和(9)代入到式(7)中，可解出：

因此，在引入虚拟边界后，1个两层偏心结构被简化成同心圆的结构，这将为计算带来方便。多层偏心单轴各向异性几何结构的反射矩阵可通过两层结构由1个迭代过程得到。

在得到最内层边界的反射矩阵之后，钻井测井模型可等价简化成包含钻轴、天线、井孔和虚拟边界在内的简单几何结构。经过变换[5,11]，可得到井孔中电磁场系数的方程：

式中：为钻轴的反射矩阵。可以看到：直接对方程(11)求解需要处理矩阵的逆问题，而此矩阵在波数kz、偏移距离d或各向异性比率η较大时条件数会指数增加，变成病态矩阵，求解其逆矩阵将引入较大误差。因此有必要给出一种迭代格式求解方程。通过将式(11)中左边第 2项移到右边，可以得到一种迭代格式：

式中：l=0，1，2，…，表示迭代次数。这种迭代格式具有直接的物理含义：迭代次数l代表波在井孔内被反射的次数。由于波在井孔中传播有损耗，因此，迭代可以在有限次数内停止。同时，亦可根据精度要求设置一种自适应停止算法。

3 数值实验

3.1 有向测井实例

图4 LWD测井工具靠近岩石层的情况及其几何模型Fig.4 Situation and geometry when LWD detection tool close to rock layer

当钻轴旋转时，倾斜测井工具可提供周围地质层的在圆柱坐标中的位角方向信息。在倾斜环状天线旋转1周覆盖360°，测井工具可提供角向数据以帮助实时决定下一步钻探方向。当周围地质层中电阻率存在明显差异时，可被测井工具探测到。这里采用一个与文献[3]中相似的例子。假设有1个导电岩层位于测井工具135°角方向，如图4所示。岩层的边界近似为直线，下面用1个半径足够大的圆来模拟此岩层边界。当钻井工具工作时，钻轴将做360°旋转，测量结果将随φ呈周期性变化。考虑由1个发射天线和2个接收天线组成的LWD工具，其中2个接收天线位于距发射天线上方609.6 mm和762.0 mm处。在图4中，导电岩层与井孔边缘的距离为685.8 mm，也就是d0=685.8 mm。用一系列偏心圆Red来近似岩层边界(其中d表示圆的偏心距离)。井孔和地质层的电导率分别为 5×10-4S/m和0.5 S/m，岩石层为单轴各向异性介质，其电导率水平分量σh=10 S/m，垂直分量σv=2.5 S/m。钻轴、环形天线和井孔半径分别为101.6，114.3和127.0 mm。在此例中，假设钻轴与井孔同心。发射天线、下部接收天线和上部接收天线的倾角分别为θt=-45°，θr1=0°和θr2=45°，测井工具的操作频率为2 MHz。图5所示为测井工具中电压相位角(PD)的响应结果，其中旋转角度φ和近似圆的偏移距离d为变量，用于模拟岩石层边界的圆半径与其偏移距离呈线性关系。由图5可以看到：当d(或圆半径)增大时，结果收敛；φ在135°和315°时达到局部极值点，此时天线正对或背对导电岩层，实际操作中可根据此信息调整钻探方向。在前面给出了伪解析法的归一化形式，这为用半径较大的圆边界来模拟岩石层奠定了基础。

图5 LWD工具对不同半径圆模拟岩石层边界的仿真响应情况Fig.5 LWD response through rock boundary for different radius rounds

3.2 可探测的深度

钻轴旋转的过程中，LWD工具可探测井孔周围地质层的电磁信息，但随着电磁波在地质层中衰减，当目标储层距离过远时，工具将无法有效探测到目标物，此例用数值方法分析 LWD工具的可探测深度。仍然采用图4中的模型，假设岩层电导率水平分量为σh=10 S/m，垂直分量为σv=4 S/m，地质层电导率为σ=0.5 S/m。对下面 2种情况，用接收天线在岩石层方向及其相反方向的相位差ΔPD来衡量工具灵敏度。模拟中取半径足够大的圆来近似岩层边界。

(1) 首先考察常规 LWD工具可探测的岩石层深度。在常规 LWD工具中，两接收天线离发射天线的垂直距离分别为zr1=609.6 mm和zr2=762.0 mm，考察钻轴与岩层的距离发生改变时的情况，假设图4中井孔与岩石层边界d0∈[127，1 143]mm变化，则d0变化过程中ΔPD的响应情况如图6所示。由图6可以看出：当测井工具远离目标储物层时，工具响应呈指数衰减。

图6 LWD工具对不同深度岩石层的响应情况Fig.6 Response of LWD tools at lithosphere of different depths

(2) 假设图 4中井孔与岩层的距离保持不变(d0=685.8 mm)，考察天线间距发生变化时的测井工具的响应情况。固定两接收天线间距zr2-zr1=152.4 mm，用接收天线的中心位置z=(zr2 +zr1 )与发射天线zt的距离来刻画所使用测井工具的尺寸。图 7所示为|zr-zt|在228.6～939.8 mm时ΔPD的变化情况。由图7可以看出：测井工具响应随发射天线和接收天线距离增大而呈线性增加。一般来说，LWD工具的探测深度约为2|zr-zt|。

3.3 对导电岩层的灵敏性

当目标储物层与周围地质层电磁属性存在较大差异时，储物层才可能被测井工具发现。为了简便，假设图4中的岩石层为各向同性介质，地质层电导率为0.5 S/m，岩层电导率在0.5～10 S/m 之间变化。在此例中钻轴中心与井孔中心重叠，无偏心情况，在测井工具随钻轴转动1周内，用伪解析法计算接收天线在各个方位角的相位差，结果如图8所示。由图8可以看到：当岩层电导率逐渐增大时，相位角响应差异明显，特别是在岩石层所在方位角φ=135°和相反方向φ=315°。

图7 LWD工具灵敏度随天线距离的变化情况Fig.7 Changes of sensitivity of LWD tools with antenna distance

图8 LWD工具对岩石层不同电导率的响应情况Fig.8 Response of LWD tools toward different conductivities of rock

3.4 钻轴偏移产生的影响

在工程实际中，由于机械牵引或重力等因素的作用，钻轴往往会发生偏移，其中心位置偏离所在井孔的中心。前面给出了一种有效处理这种偏心(以及多层偏心)结构的方法，这里考察钻轴偏移对探测结果的影响。假设岩层电导率水平分量σh=10 S/m，垂直分量σv=4 S/m，地质层为同性介质电导率σ=0.5 S/m。仍然假设井孔与岩层边界的距离d0=685.8 mm，用半径足够大的圆模拟岩层边界。实验中选取的圆半径为1 701.8 mm。钻轴偏移的位移分别选取de=25.4，50.8，…，127.0 mm，偏移方位角φe=15°，30°，45°，…，360°。图9所示为接收电压在岩层所在方向及其相反方向处ΔPD的变化情况。图9中观测点与坐标原点所成的角代表钻轴偏移的方位角φe，5条闭合曲线分别代表偏移位移de=25.4，50.8，…，127.0 mm 5种情况，与原点的距离表示相应的ΔPD的计算值。由图9可以看出：当钻轴发生偏移时，接收信号在方位角方向的分布会随之移动，并在钻轴偏移的方向(及其相反方向)形成峰值。此实验结论可用于对钻轴进行精确定位和校正。

图9 钻轴偏移对LWD工具测量结果的影响Fig.9 Impact of drilling shaft migration in measurement results of LWD tools

4 结论

(1) 推导了伪解析法在单轴各向异性介质中的计算公式，并提出处理复杂介质和大尺寸仿真问题的有效办法。归一化Bessel函数不仅能大大增加可计算目标的尺寸，还能提高算法的数值精度，同时为加强算法的稳定性提出了迭代求解策略。

(2) 本文所给出的算法计算速度快、精度高，并能有效仿真 LWD工具在测井工程中的各项性能，是一种高效、实用的方法。LWD工具的测井结果与目标储物层的距离呈指数关系，与天线间垂直距离呈线性关系。

致谢

感谢美国俄亥俄州立大学F.L.Teixeira教授对本文的前期指导工作，感谢 Halliburtion能源服务公司B.Donderici博士提供部分对比数据，同时感谢俄亥俄超级计算中心为本文提供实验平台。

[1]Wait J R.Geo-electromagnetism[M].New York∶ Academic Press,1982∶ 107-154.

[2]Bittar M S.Electromagnetic wave resistivity tool having a tilted antenna for determining the horizontal and vertical resistivities and relative dip angle in anisotropic earth formations∶ US 7557579B2[P].2009-07-07.

[3]Chou L, Li Q, Darquin A, et al.Steering toward enhanced production[J].Oilfield Review, 2005, 17(3)∶ 54-63.

[4]Hue Y K, Teixeira F L, Martin L S, et al.Three-dimensional simulation of eccentric LWD tool response in boreholes through dipping formations[J].IEEE Transaction on Geoscience and Remote Sensor, 2005, 43(2)∶ 257-268.

[5]Hue Y K, Teixeira F L.Analysis of tilted-coil eccentric borehole antennas in cylindrical multilayered formations for Well-Logging applications[J].IEEE Transaction on Antennas and Propagator, 2006, 54(4)∶ 1058-1064.

[6]Novo M S, Silva L C, Teixeira F L.Finite volume modeling of borehole electromagnetic logging in 3-D anisotropic formations using coupled scalar-vector potentials[J].IEEE Antennas and Wireless Propagator Letters, 2007, 6(1)∶ 549-552.

[7]Lovell J R, Chew W C.Response of a point source in a multicylindrically layered medium[J].IEEE Transaction on Geoscience and Remote Sensor, 1987, GE-25(6)∶ 850-858.

[8]Lovell J R, Chew W C, Effect of tool eccentricity on some electrical well-logging tools[J].IEEE Transaction on Geoscience and Remote Sensor, 1990, 28(1)∶ 127-136.

[9]Wang H, So P, Yang S, et al.Numerical modeling of multicomponent induction well-logging tools in the cylindrically stratified anisotropic media[J].IEEE Transaction on Geoscience and Remote Sensor, 2008, 46(4)∶ 1134-1147.

[10]Lee H O, Teixeira F L, Cylindrical FDTD analysis of LWD tools through anisotropic dipping-layered earth media[J].IEEE Transaction on Geoscience and Remote Sensor, 2007, 45(2)∶383-388.

[11]Hue Y K.Analysis of electromagnetic well-logging tools[D].Columbus, OH∶ The Ohio State University.Department of Electrical and Computer Engineering, 2006∶ 3-4.

[12]Zhong L, Li J, Bhardwaj A, et al.Computation of triaxial induction logging tools in layered anisotropic dipping formations[J].IEEE Transaction on Geoscience and Remote Sensor, 2008, 46(4)∶ 1148-1163.

[13]Wang G L, Torres-Verdin C, Gianzero S.Fast simulation of triaxial borehole induction measurements acquired in axially symmetrical and transversely isotropic media[J].Geophysics,2009, 74(6)∶ e233-e249.

[14]Chew W C.Waves and fields in inhomogeneous media[M].Piscataway∶ IEEE Press, 1995∶ 219-270.

[15]Chew W C.The singularities of a Fourier-type integral in a multicylindrical layer problem[J].IEEE Transaction on Antennas and Propagator, 1983, AP-31(4)∶ 653-655.