倾斜层理地层随钻电磁波测井响应特征

范宜仁，胡旭飞，邓少贵，袁习勇，李海涛

（1.中国石油大学（华东）地球科学与技术学院，山东青岛 266580；2.深层油气地质与勘探教育部重点实验室，山东青岛266580；3.海洋国家实验室海洋矿产资源评价与探测功能实验室，山东青岛 266071；4.中国石油大学CNPC测井重点实验室，山东青岛 266580；5.中国石油大学（北京）克拉玛依校区石油学院，新疆克拉玛依 834000）

1 问题提出

随钻电磁波电阻率测井数据能够用于流体识别和流体饱和度计算等储集层评价研究。为了实现快速反演，实时获取地层的真实地层参数，通常需要将复杂井眼和地层环境（三维问题）简化为一系列一维（1D）地层模型[1-5]，这种方法在随钻测井中是可行的，原因主要有两点：①随钻测井过程中，井眼及钻井液侵入影响较小，可以忽略；②测井尺度下，复杂的井眼轨迹可通过逐步开窗的方式进行约束，简化为一系列形状单一的井眼轨迹（即一系列1D地层模型）。传统的1D快速随钻电磁波数值模拟主要基于含垂直对称轴的横向各向同性（VTI）地层模型。低能沉积环境条件下，如深海和深湖中沉积形成的地层，可呈现出标准的VTI介质特征（见图1a），VTI介质属于单轴各向异性介质的特例，表现为任意水平方向上，介质具有相同电阻率Rh；垂直方向上，介质具有不同的电阻率Rv，一般情况下Rh≤Rv。该类介质的电磁场计算研究较为深入，可通过引入赫兹位函数，将电磁场波动方程简化为两个互相不耦合的标量波动方程，最终得到的电磁场表达式为含0阶和1阶Bessel函数的一重积分[5-7]。

图1 不同沉积环境下的地层特征

随着石油勘探开发的不断深入，各种复杂地层情况下的油气藏成为研究重点，这使得VTI地层模型应用受到一定限制，并非所有的地层都能简化为VTI模型。高能沉积环境（河流和沙漠环境）条件下形成的地层，常见倾斜层理（见图1b），此时的地层为含倾斜对称轴的横向各向同性（TTI）地层，具体表现为：沿层理方向电阻率相同，垂直层理方向的电阻率不同，如果仍然将其简化为VTI地层模型，则反演得到的电阻率等参数就不再可靠，若将其用于定性评价和定量计算，则解释的符合率会下降。明确TTI介质条件下的电磁波类测井响应特征，能够更加全面有效地进行测井定性及定量评价。

基于TTI介质的电磁场理论和计算均较复杂，目前国内在这方面未见相关文献，斯伦贝谢Anderson等人[8-9]简单讨论了由倾斜层理引起的地层各向异性的感应测井响应特征；Wang等人[10-12]给出了倾斜层理引起的地层各向异性的多分量感应测井正反演数值模拟算法，并讨论分析了不同的倾斜层理情况下多分量感应测井的响应特征。Anderson和Wang等人均指出[8-12]，倾斜层理对测井响应特征的影响不容忽视，尤其在实际资料运用中，要考虑实际地层条件，建立合适的地层模型。笔者首先从矢量格林函数出发，采用各向异性角和各向异性方位表征含倾斜层理TTI介质的各向异性特性，推导出TTI介质的电磁场求解方法；然后采用数值模拟方法，以半空间介质为例，分析VTI介质和TTI介质的电磁波反射和透射特征，指出了传统的赫兹位函数方法以及广义反射系数等方法的局限性，并验证了算法的可靠性，得到不同井斜、各向异性角和各向异性方位情况下的随钻电磁波响应特征。

2 TTI地层正演数值模拟

2.1 基于并矢格林函数的各向异性介质电磁场

无限厚各向异性介质中，电场的波动方程可以表示为：

VTI介质中，各向异性地层的电导率张量可表示为：

TTI介质中，可定义各向异性角Ψ和各向异性方位χ表征其各向异性特性（见图2），则TTI介质中电导率张量可表述为：

图2 单轴各向异性地层电导率张量旋转示意图

利用傅里叶变换，空间域电场表达式为：

（4）式可进一步改写为：

其中

行列式|W|是如下关于kz的一个4阶多项式：

式中a，b，c，d，e为关于kx，ky的表达式，进一步可表示为：

当|W|=0时，kz存在4个解，表示两种类型的波，为Ⅰ型波和Ⅱ型波向上的z方向的波数，为Ⅰ型波和Ⅱ型波向下的z方向的波数[13-14]。

将（7）式代入（5）式得：

式中R=r-r′，dk=dkxdkydkz。

利用留数定理，当接收线圈位置z在源位置z′的上方时，则电场E(R)为：

其中：

当接收线圈位置z在源位置z'的下方时，电场E(R)为：

其中：

N层介质中，设信号源在第S层，接收线圈在第L层，接收线圈位置电场为：

δLS为Kronecker delta函数，当L=S时，δLS=1；当L≠S时，δLS=0；Ed（R）为直达场，Es（R）为散射场（包括反射和透射场）。因此，只要推导出直达场和散射场，就可以得到各层介质中的电场。相应的磁场分量可以通过如下表达式，由电场分量获取：

多层TTI介质中的电磁场推导可参考Wang等人文献[10-12]，限于篇幅，本文不再详述。

2.2 传统随钻电磁波电阻率测井原理

基于上节推导可得到TTI介质中接收线圈处的电磁场，进而计算得到接收线圈的电动势。以传统随钻电磁波测井的单发双收线圈结构为例，设两个接收线圈处的电动势分别为V1和V2，则幅度比EATT和相位差Δφ定义如下：

利用相位差和幅度比的刻度图版即可将相位差和幅度比转换为幅度比电阻率Rad和相位差电阻率Rps。

3 算法验证及半空间反射和透射分析

首先利用相位差和幅度比电阻率的刻度图版来验证算法的正确性，图3为传统随钻电磁波频率为2 MHz、发射线圈到仪器记录点距离为40.64 cm（16英寸）的相位差电阻率和幅度比电阻率刻度图版，可见基于并矢格林函数的解析方法得到的计算结果与传统的利用Hertz位函数的解析方法得到的计算结果完全吻合。其次，采用两个5层地层模型进行验证（井斜分别为0和30°，层厚均为2.0 m，电阻率参数见表1），如图4a和4b所示，可以看到，本文提出的算法在多层介质中同样适用。

为研究不同介质界面处的电磁波反射和透射特征，笔者给出半空间的反射和透射系数。图5为不同情况下的半空间介质模型，表2为对应不同半空间介质模型下的参数。σ1u和σ2u表示上层电导率，σ1d和σ2d表示下层电导率。上层地层如为各向同性地层，则σ1u=σ2u；如为各向异性地层，则σ1u≠σ2u，下层地层类似。

各向同性地层和VTI地层中，Ψ=0，χ=0；TTI地层中，Ψ≠0，限于篇幅，设定χ等于0（当χ≠0时，反射和透射情况相同，此处不再详述）。图6为对应不同模型界面处的Fresnel反射和透射系数随入射角的变化情况。Ri,j和Ti,j分别表示Fresnel反射和透射系数，i和j为电磁波类型（Ⅰ型或Ⅱ型波）。从图6a到图6c可以看到，当界面两侧的介质为各向同性或VTI介质时，即Ψ=0，χ=0的情况下，RⅠ,Ⅱ=RⅡ,Ⅰ=TⅠ,Ⅱ=TⅡ,Ⅰ=0，仅存在RⅠ,Ⅰ，RⅡ,Ⅱ，TⅠ,Ⅰ和TⅡ,Ⅱ，说明不存在交叉反射波和透射波，即当入射波为Ⅰ型时，会产生反射Ⅰ型和透射Ⅰ型波；当入射波为Ⅱ型时，情况类似。因此，界面两侧为各向同性或VTI介质时，Ⅰ型波和Ⅱ型波在界面处的反射和透射不耦合，可以分别考虑及单独计算Ⅰ型波和Ⅱ型波（也称为TE波即横电波和TM波即横磁波）的电磁场，这就是传统基于VTI介质计算TE波和TM波电磁场的理论基础。从图6d到6f可以看出，当界面一侧存在TTI地层时，各个交叉分量（RⅠ,Ⅱ，RⅡ,Ⅰ，TⅠ,Ⅱ，TⅡ,Ⅰ）均不再为0，说明经过反射和透射后，Ⅰ型波和Ⅱ型波耦合，其电磁场为两种类型波的耦合场，此时传统基于VTI介质的电磁场计算方法已不再适用，可以用本文提出的基于并矢格林函数的方法进行计算。

图3 幅度比和相位差刻度图版验证

表1 两个5层地层模型的电阻率参数

图4 5层地层模型算法验证（模型1和模型2）

图5 不同半空间介质模型

表2 不同半空间介质模型地层参数

图6 不同半空间介质界面处的Fresnel系数

4 数值算例

设各向异性地层为无限厚，一个方向的电阻率R1=2.0 Ω·m，另一个方向的电阻率R2=20.0 Ω·m，以传统随钻电磁波频率为2 MHz、发射线圈到仪器记录点距离为40.64 cm（16英寸）的相位差为例，分析其受各向异性角和各向异性方位影响。图7为传统随钻电磁波相位差随井斜角θ的变化情况，当ψ=0时，相位差随井斜角θ增大呈单调递减趋势，各向异性方位χ对相位差随井斜角θ的变化没有影响；当ψ≠0时，如果χ=0，则存在一个临界井斜角θc=90-ψ，即当ψ分别为30°、45°和60°时，θc分别为60°、45°和30°（见图7a）。当井斜角θ≤θc时，随着井斜角的增大，相位差递减；当θ＞θc时，随着井斜角的增大，相位差增大。如果χ≠0，仍存在临界井斜角θc，规律与图7a所示相同（见图7b），但θc≠90-ψ。

图7 相位差随井斜角θ的变化

图8为传统随钻电磁波相位差随各向异性角ψ的变化。首先可以看到，当θ=0时，相位差随ψ增大呈单调递减趋势，χ对相位差随ψ的变化没有影响；当θ≠0，若χ=0，则存在一个临界各向异性角ψc=90-θ，即当θ分别为30°、45°和60°时，ψc分别为60°、45°和30°（见图8a），当ψ≤ψc时，随着各向异性角的增大，相位差递减；当ψ＞ψc时，随着各向异性角的增大，相位差增大。如果χ≠0，仍存在临界各向异性角ψc（见图8b），变化规律与图8a所示相同，但ψc≠90-θ。

图8 相位差随各向异性角ψ的变化

值得注意的是，通常认为，随钻电磁波在界面处产生的“犄角”主要由井斜和电阻率对比度等因素造成（本质为界面两侧的电场不连续[15]），但实际上各向异性角和各向异性方位也是影响“犄角”产生的重要因素。图9a为两层地层模型及其参数，上部地层设为各向同性地层（χ0=0，ψ0=0），下部为单轴各向异性地层。

首先设下部地层的各向异性方位角χ1=0，考察下部地层不同ψ1和θ情况下，传统随钻电磁波测井在界面处的响应。当下部地层为VTI地层，即ψ1=0和χ1=0（见图9b和图9c），可见随着井斜角的增大，下部地层幅度比电阻率Rad和相位差电阻率Rps均增大，且当井斜角θ小于60°时，界面处不会产生犄角；当下部为TTI地层，设各向异性角ψ1=30°时（见图10a和10b），当θ分别为30°、45°和60°时，与图9b和9c相比，下层的Rad和Rps电阻率幅度发生明显变化，此外，Rad在界面处产生明显“犄角”；当θ=60°时，Rps在界面处产生明显犄角。因此，与VTI地层相比，相同井斜角情况下，各向异性角的存在可能会使随钻电磁波测井在界面处产生“犄角”。

图9 两层介质、不同井斜角传统随钻电磁波测井响应（ψ1=0，χ1 =0）

图10 两层介质、不同井斜角传统随钻电磁波测井响应（ψ1=30°，χ1=0）

基于图9a的两层地层模型及参数，考查下部地层不同各向异性方位角χ1≠0时传统随钻电磁波测井响应在界面处的响应。由图11中可以看到，Rad和Rps在界面处的“犄角”幅度会受各向异性方位的变化而变化，因此各向异性方位也是影响“犄角”产生的原因之一。

图11 不同各向异性方位χ情况下传统随钻电磁波测井响应（θ=30°，ψ1=30°）

图12为两个模拟实例，均为5层地层，井斜为45°，具体地层参数见表3。其中图12a中各层不存在倾斜层理，图12b中第2和第4层均存在倾斜层理，但其各向异性角和各向异性方位不同。图13为发射频率2.0 MHz、发射线圈到仪器记录点距离为40.64 cm（16 in）的传统随钻电磁波测井响应。由图13可见，尽管两个地层模型的电阻率参数相同，但其各向异性的特性不同，随钻电磁波测井幅值及其在界面处的响应特征均会发生变化，因此，如果将地层始终简化为VTI地层而不考虑其内部层理变化，则无论其测井响应还是反演得到的地层真电阻率均不能表征真正的地层电性特征。

图12 5层地层模型实例

表3 两个5层地层模型的地层参数

图13 5层地层模型传统随钻电磁波测井响应特征

5 结论

为满足快速电阻率反演需求，传统随钻电磁波测井正演主要基于VTI地层模型，但实际地层更加复杂，VTI地层仅仅是各向异性地层的一个特例，笔者从并矢格林函数出发，详细推导了一套多层TTI地层的随钻电磁波测井正演计算方法。基于数值算例对算法进行了验证，结果表明笔者提出的算法不仅适用于VTI地层，也适用于TTI地层，具有较强的普适性。

通过半空间电磁波反射和透射特征数值模拟，可以得到当界面两侧为各向同性或者VTI地层时，两种类型的电磁波（Ⅰ型波和Ⅱ型波）产生的电磁场不耦合，随钻电磁波的正演数值模拟可采用传统的赫兹位函数方法；当界面的某侧为TTI地层时，两种类型的电磁波（Ⅰ型波和Ⅱ型波）产生的电磁场会耦合，随钻电磁波的正演数值模拟可采用本文算法。

各向异性角ψ和各向异性方位χ对传统随钻电磁波测井的响应不容忽视。TTI地层中，存在临界井斜角θc和临界各向异性角ψc，当θ≤θc时，随着井斜角的增大，相位差递减；当θ＞θc时，随着井斜角的增大，相位差增大；当ψ≤ψc时，随着各向异性角的增大，相位差递减；当ψ＞ψc时，随着各向异性角的增大，相位差增大。此外，当χ=0时，θc=90-ψ，ψc=90-θ。

随钻电磁波在界面处产生的“犄角”，不仅与井斜、电阻率对比度等因素有关，而且也受各向异性角和各向异性方位的影响，如含倾斜层理地层中，在低井斜角情况下，也有可能产生“犄角”。

符号注释：

E——电场，V/m；Ed（R）——直达场，V；Es（R）——散射场，V；——电磁场波数域表达式；EATT——幅度比，dB；f——频率，Hz；H——磁场，A/m；i，j——电磁波类型（Ⅰ型波和Ⅱ型波）；k——波数，无因次；kx——x方向波数，rad/m；ky——y方向波数，rad/m；kz——z方向波数，rad/m；——Ⅰ型波和Ⅱ型波向上的z方向的波数，rad/m；——Ⅰ型波和Ⅱ型波向下的z方向的波数，rad/m；M——磁流源，A·m2；M0——磁流源强度，A·m2；r——接收位置（x，y，z），m；r′——源位置（x′，y′，z′），m；Rad——幅度比电阻率，Ω·m；Rh——水平电阻率，Ω·m；Ri,j——Fresnel反射系数；Rps——相位差电阻率，Ω·m；Rv——垂直电阻率，Ω·m；——各向异性角旋转矩阵，（°）；——各向异性方位旋转矩阵，（°）；Ti,j——Fresnel透射系数；V1和V2——两个接收线圈的感应电动势，V；|V1|和|V2|——感应电动势V1和V2的幅度值，V；μ——自由空间的磁导率，H/m；ω——角频率,rad/s；σh——水平电导率，S/m；σν——垂直电导率，S/m；Ψ——各向异性角，（°）；χ——各向异性方位，（°）；→cσ——VTI介质中的电导率张量，S/m；——TTI介质中的电导率张量，S/m；δLS——Kronecker delta函数，当L=S时，δLS=1，当L≠S时，δLS=0；Δφ——相位差，（°）；φ1，φ2——感应电动势V1和V2的相位角，（°）；θ——井斜角，（°）；θc——临界井斜角，（°）；ψc——临界各向异性角，（°）；σ1u，σ2u——两层地层模型中的上部地层两个方向的电导率，S/m；σ1d，σ2d——两层地层模型中的下部地层两个方向的电导率，S/m；∇——哈密顿算子。