井外溶洞储层的偶极横波反射特征研究∗

杨 岩 张 波 李 超 张晋言 许孝凯 陈 浩 王秀明

(1 中国科学院声学研究所 声场声信息国家重点实验室 北京 100190)

(2 中国科学院大学 北京 100049)

(3 北京市海洋深部钻探研究中心 北京 100190)

(4 中石化胜利石油工程有限公司测井公司 东营 257061)

0 引言

碳酸盐岩缝洞型油气藏是我国重要的油气藏类型之一，具有广阔的勘探和开发前景。由于其非均质性强，埋藏较深，尺度分布不均，给储层识别带来了挑战。其中溶洞体在地震资料上常表现为多个强能量团的串珠状响应特征[1−2]，但受到地震勘探分辨率的限制，一些尺度较小的溶洞往往难以识别。

利用阵列声波测井中记录的反射波来进行井旁构造成像的远探测声反射成像测井技术是近年来声波测井领域的重要进展之一，它可以探测裂缝、溶洞等储层在井外的扩展情况，以及水平井的储层边界，还可以识别地震勘探无法探测到的小构造如小断层和岩丘等[3]，目前逐渐成为复杂油气藏勘探的有力手段。由于低频偶极源有较好的辐射特征及对反射方位的敏感性，利用正交偶极声波测井仪测量的反射横波来进行反射体成像及方位识别的技术近年来发展迅速，在碳酸盐岩等非均质储层勘探中发挥着越来越重要的作用。Tang 等[4]于2009年提出采用偶极声源激发的横波进行反射声波测井，通过对比分析，指出辐射的SH 横波在井孔所在的垂直平面中有较好的辐射覆盖及较高的反射灵敏度，可以得到有方位性的更加清晰的声成像结果。2012年，唐晓明等[5]论述了偶极反射声波测井技术的基本原理与方法，主要指出快地层中SH 型横波具有较宽的辐射覆盖与反射灵敏度，慢地层中偶极辐射的纵波具有较好的探测特性。偶极横波远探测技术的发展使得反射声波测井的探测范围从几米延伸到了几十米甚至上百米。Li等[6]研究了二维情况下偶极远探测中井外缝、洞的反射波响应特征及逆时偏移成像方法。

随着计算机的发展，三维有限差分(Finitedifference time-domain,FDTD)技术成为了研究偶极反射声波测井的主要工具。何峰江[7]综合运用二维与三维数值模拟方法，系统研究了反射波的幅度随声源频率、源距和反射界面位置变化的关系。Wei 等[8]采用三维有限差分对偶极声源的辐射特性、SH 反射波与SV 反射波的幅度特性、反射波幅度与界面方位的关系等问题进行了系统研究。前人在研究远探测响应特征一般都采用无限延伸的平面裂缝模型，或者在二维坐标系下的简化模型，而针对井外溶洞型储集体的声波测井响应研究还较少，尤其是随深度连续变化的偶极横波，这是由于考虑了井孔在内，井周三维大尺度范围内的远探测模拟计算量较大。本文基于三维有限差分及并行计算方法，研究了不同尺寸以及离井眼不同距离的溶洞体的偶极横波反射特征，尤其是SH 反射横波的特征。研究结果为实际的现场资料解释工作提供了必要的理论基础。

1 方法

在直角坐标系中，对于各向同性介质，一阶应力与速度方程可以写为

其中，

式(1)∼(2)中，vi(i=x,y,z)、τij(i,j=x,y,z)和gij(i,j=x,y,z)分别为速度分量、应力分量和体积源分量。式(3)中，λ和µ是介质的拉梅系数，与纵横波速度Vp和Vs的关系为对于液体介质，比如井孔流体和储层流体，被认为是一种具有零剪切速度的特殊各向同性介质。速度-应力方程采用交错网格在空间-时间域上进行离散。各速度分量和应力分量交错分布在空间网格上，相邻分量之间相差半个网格；在时间上各分量也是交错分布的，相邻的速度分量和应力分量之间相差半个时间步长。以vx和txx为例，其离散方程为[9]

其中，δi(i=x,y,z)是差分算子，∆t是时间步长。为了减弱数值频散，空间和时间步长应该满足：

其中，Cmin和Cmax分别是计算模型的最小速度和最大速度，fmax是声源的最高频率，∆x、∆y和∆z为x、y和z方向的空间步长。本文中采用非分裂完全匹配层作为吸收边界条件[10]，无需分解波场的各个变量，而是直接进行复数扩展代换，虽然在模拟中引入了卷积运算，但是各个物理量的完整性得到了保证。

由于三维计算模型较大，整个有限差分计算采用非均匀交错网格并行算法来提高计算效率。在井孔附近对网格进行精细划分，再在一般地层介质中则使用大网格以节约计算资源[11]。采用消息传递和共享存储的方式，进程间采用消息传递并行编程MPI，各进程内部采用共享内存的OpenMP 并行，这样的设计可以使得共享了进程内存地址空间的多个线程共同执行相应进程的指令。基于天河二号并行计算集群，结合上述并行计算及网格划分方法，包含仪器提升的偶极三维远探测波场模拟的大计算量问题得到了有效解决。

2 数值计算

2.1 溶洞模型

首先建立包含井孔的三维溶洞模型。图1 为三维数值模拟模型的xz方向切面示意图，模型大小为x=16 m、y=8 m、z=16 m，井孔半径为0.1 m，采用非均匀网格，在井孔附近网格大小为0.01 m，地层中网格大小为0.02 m，时间步长为16.0×10−3s。如图1 所示，溶洞为一球体，球内充满水，地层为碳酸盐岩地层，地层参数如表1所示，球心距离井轴为L，球体直径为d。

图1 溶洞模型示意图Fig.1 Schematic diagram of karst cave model

表1 地层参数Table 1 Formation Parameters

图2 是距离井轴4 m、直径为1.2 m 的溶洞SH反射波场，其中图2(a)、图2(b)和图2(c)为共偏移距道集(COG)，D为源距；图2(d)、图2(e)和图2(f)为共源道集(CSG)，P为声源提升高度。从反射波到时分析可以得出，红色方框内能量较弱的为来自溶洞近井壁界面的SH 反射波，后续能量较强的为来自溶洞远井壁界面的反射波，初至SH 反射波能量较弱，由于溶洞尺寸较小，跟后续振幅较强波场较难区分(主要是指球内聚焦波，见下文)，这也会给后续的成像及解释工作带来困难，可能会出现解释错误。图3为距离井轴10 m、直径为4 m 的溶洞SH反射波场，可以看到随着溶洞与井孔距离的变大，溶洞的反射到时以球心为中心成圆弧状的这种特点越来越明显，且反射波的在接收源距上的分布也表现出了特点，即在小源距接收器范围内会出现较强的后续反射波，随着源距的增大，该反射波逐渐变弱甚至消失。图4 为距离井轴20 m、直径为20 m 的溶洞SH 反射波场，与图3 比较发现，随着溶洞直径的进一步增大，来自溶洞远井壁界面的后续反射波在小源距接收器范围内仍然具有很强的振幅，而且其会随着仪器的提升而变化，在整个接收器源距的分布范围变大。

为了研究后续反射波随着源距变化的分布特征，以离井轴距离为10 m、直径为4 m 的溶洞模型为例，提取波场快照，如图5 所示。图5(a)和图5(b)分别为声波在地层中向外传播和经过溶洞前壁时形成的反射波，可以看出溶洞前壁形成的反射波是一个能量均匀的球面波，因此在波形图中不管是大源距和小源距的接收器中都能接收到；图5(c)为经过溶洞后壁反射但还未传播出溶洞的反射波；图5(d)为经过溶洞后壁反射在地层中传播的反射波，可以清晰地看到经过溶洞后壁会聚的声波不再是能量均匀的球面波，而是有方向性的声波，从波场快照中可以明显看到，溶洞后壁相当于一个凹面镜，对声波有聚焦效果，而随着溶洞直径的变大，该聚焦效果变差，聚焦范围随之变大，这也解释了当溶洞直径较小时反射波场后续的接收波形只能在小源距的接收器中接收到，而随着溶洞直径的变大该反射波不仅能在小源距的接收器中接收到，随着仪器在井中的提升在其他源距接收器中也可以接收到后续能量较强的反射波形，该反射与前部SH 反射波的时间间隔与溶洞大小以及离井轴距离有关。

因此，当溶洞体尺寸较大时，在反射波场中可以明显地观察到分别来自溶洞前后两个界面的反射波，对比两种不同的反射波，其特点是在小源距接收器中总能够接收到来自溶洞后壁(远离井壁一侧)能量较强的反射波，而在大源距接收器中可以接收到来自溶洞前壁(即靠近井壁一侧)能量较强的反射波，利用该特征可以对前后两个界面分别成像，分析溶洞体的尺寸。研究结果也对之后的仪器设计提供了理论支持，这就要求对于溶洞储层的探测要采用变源距的声波测井仪器，提高接收器的覆盖范围。

图2 不同数据道集模拟结果(L=4 m, d=1.2 m)Fig.2 Simulation results of different data trace sets (L=4 m, d=1.2 m)

图3 不同数据道集模拟结果(L=10 m, d=4 m)Fig.3 Simulation results of different data trace sets (L=10 m, d=4 m)

图4 不同数据道集模拟结果(L=20 m, d=20 m)Fig.4 Simulation results of different data trace sets (L=20 m, d=20 m)

图5 不同时刻的波场快照Fig.5 Wave field snapshot at different time

2.2 缝洞结合体模型

在实际地层中，裂缝和溶洞往往是交错存在的，在溶洞模型基础上，建立缝洞结合体模型，如图6所示，裂缝穿过溶洞，倾角为70◦，裂缝厚度为0.04 m，L=10 m，d=4 m。图7为缝洞模型的共偏移距道集，在反射波场中，裂缝SH 反射波清晰可见，后续在小源距范围内的接收器中同样观察到了来自溶洞后壁会聚的反射波(红色方框内)，由于裂缝的干扰反射波形会出现断断续续的现象，并且其能量要弱于来自裂缝的反射波。不同于之前单溶洞模型的是，来自溶洞前壁的反射波湮没在裂缝反射中，无法识别。因此，在缝洞结合体的反射波场中，溶洞的反射波场要弱于裂缝的反射波场，且只能识别出来自溶洞后界面的反射波，这会对后续成像以及估计溶洞尺寸带来困难。

图6 缝洞体模型示意图Fig.6 Schematic diagram of crack-cave model

2.3 溶洞反射体的方位识别

将图3 对应的模型(L=10 m,d=4 m)得到的四分量数据进行旋转，可以获得SH 波和SV 波不同方位的波场切片[4]：

图7 缝洞体COG 道集模拟结果(L=10 m, d=4 m)Fig.7 Simulation results of COG trace set of crack-cave model (L=10 m, d=4 m)

图8 溶洞体方位分析Fig.8 Orientation analysis of karst cave model

其中，角度θ为走向角，对于远场而言，一般SH反射波的幅度大于SV 波，SH 反射能量最大的角度即为反射界面走向角度，即偶极源极化方向和反射界面平行时，产生的SH 波反射幅度最大。四分量的共偏移距道集波形及方位分析结果见图8，图8(e)显示SH 波最大幅度的对应的走向角在每个深度位置都是90◦，并且在中间(球心位置)这种对应更明显，所以可以判断球体的走向角度为90◦，由于球体中心轴在x轴方向(0◦)，所以方位应该是0◦或者180◦，与模型中球心在0◦方位一致，说明传统用于确定界面走向的偶极四分量的分析方法仍然适用于溶洞。但是由于偶极源的辐射特性，以上方法求得的走向仍然存在180◦的不确定性，即无法区分在井外的倾向[12]。

2.4 成像

图9 为模型中L=10 m、d=4 m 时的不同源距阵列(图3(a)和图(b))的反射波成像图，图9(a)为最小源距为0 m 的阵列波场成像结果，图9(b)为最小源距为3 m的结果。通过成像可以进一步证实，采用小源距阵列成像，可以获得球体的第二个界面信息，不同源距有两部分界面，前后距离大于模型本身距离(球直径)，这是由于偏移时采用横波速度，而后续的波形有一部分路径是在球内传播的(如图9所示)，速度较慢，采用较高的横波速度会使距离加大。

根据波的旅行时以及速度，易得以下估算球直径的公式：

其中，s为两种阵列成像结果溶洞距离的差，vs为横波速度，vf为溶洞内介质速度，这里为水，实际中可以根据储层参数估算。由图9 可以估算出溶洞直径约4 m，这与模型一致。

图9 溶洞不同源距反射波成像(L=10 m, d=4 m)Fig.9 Reflection wave imaging of karst cave with different source spacing (L=10 m, d=4 m)

3 结论

本文通过大规模时域有限差分的方法数值计算了三维溶洞储层的不同测量位置、不同源距的连续反射波场，得出了以下结论：

(1)溶洞体的偶极源反射，除了首先到达的前界面的SH反射波外，后续还有能量较强的波列，应为球内的凹界面聚焦作用的结果。因而不同源距成像结果中，会有两部分界面，前后距离大于模型本身尺度，这是由于后续的波形有一部分路径并不是直线且不以横波速度传播，对于小源距后部反射强，大源距范围靠近井壁一侧更强，该特征可用于分析溶洞体的尺寸。本研究可为偶极横波远探测的资料解释以及仪器设计提供依据。

(2)缝洞结合体的偶极源反射波场，在小源距范围内仍然可以观察到能量较强的后续反射波，但总体来说裂缝反射强度要大于溶洞的反射强度，且来自溶洞前壁(靠近井壁一侧)的反射波会湮没在裂缝反射波中，这会对成像结果带来一定的干扰。本文模型是高角度缝洞模型，下一步的研究应该围绕缝洞结合体倾角对反射波的影响来展开讨论。

(3)对溶洞的偶极四分量反射波进行方位分析，得到的方位与溶洞中心轴的一致，说明传统用于分析界面走向的方法仍然适用于溶洞，只不过得到的方位是溶洞中心轴的。

(4)本文分析虽然是假定溶洞为一理想球体，实际中溶洞体可能形状不规则，但是只要溶洞存在凹界面，可以近似为球体，本文的结论都是适用的。