基于纵波方位AVO的裂缝流体识别方法

摘要

本发明涉及一种基于纵波方位AVO的裂缝流体识别方法。裂缝既是油气的储集空间，又是流体的运移通道，探测裂缝对寻找油气储层意义重大。裂缝的识别是目前的研究热点和难点，而对裂缝储层中流体类型的识别更是难上加难。本发明在各向异性裂缝反演的基础上，进一步研究裂缝中流体类型识别方法。通过数值模型正演模拟表明，世界著名的Ruger方程中各向同性梯度和各向异性梯度项对裂缝流体类型具有一定的识别能力，特别是能够区分干裂缝(充填气)和湿裂缝(充填油和水)。本发明优势在于可应用于实际资料，对于各向异性程度较强的地区能够识别裂缝中的流体类型。

说明书

技术领域

本发明属于石油勘探领域，涉及一种裂缝储层综合评价中裂缝流体类型识别的方法。 

背景技术

裂缝型储层大约占世界油气储量的三分之一。碳酸盐岩、泥岩、砂岩、火山岩和变质岩中均有裂缝发育。裂缝一方面是油气的储集空间，另一方面作为流体的运移通道，能够将储层中分散的孔洞储集空间连通起来。因此，预测裂缝的发育和分布对于寻找油气储层意义重大。裂缝的存在使储层具有各向异性特征。所谓各向异性是指地震波的物理性质(振幅、速度、频率等)随观测方位的变化而变化。对于垂直排列的平行裂缝，可以将其等效为HTI介质(Horizontal Transverse Isotropy，即具有水平对称轴的横向各向同性介质)模型来研究。前人对基于纵波方位AVO特征(HTI介质的振幅各向异性特征)预测裂缝密度和方位进行了大量的理论和应用研究(马中高，2003；Ruger，1998；Gray，2008；Sun等，2011a；2011b；Sun等，2012)。下面给出常用的基于纵波方位AVO的裂缝密度和方位反演方法：对于小到中角度的入射角，世界著名的Ruger纵波反射系数方程为(Ruger，1998)： 

其中 

$P=\frac{1}{2}\frac{\mathrm{\Delta Z}}{\bar{Z}}$

$\mathrm{Biso}=\frac{1}{2}[\frac{\mathrm{\Delta \alpha }}{\bar{\alpha }}-{\left(\frac{2\bar{\beta }}{\bar{\alpha }}\right)}^2\frac{\mathrm{\Delta G}}{\bar{G}}]$

$\mathrm{Bani}=\frac{1}{2}[\Delta {\delta }^V+2{\left(\frac{2\bar{\beta }}{\bar{\alpha }}\right)}^2\mathrm{\Delta \gamma }]$

式中：α、β分别为纵波速度和横波速度；Δε^V、Δδ^V、Δγ是Thomsen参数； 为采集方位角； 为裂缝走向；i为入射角；Z为纵波垂直入射时的波阻抗；G为横波的切向模量。在(1)式中的第二项(即与sin²i相乘的项)是AVO梯度，梯度项分为两部分：各向同性梯度Biso和各向异性梯度Bani(该参数能够指示裂缝引起的各向异性程度即裂缝密度)。 

把包含裂缝方位角的余弦项展开，可以把公式(1)转化为如下的形式： 

其中 

C₁＝P，C₂＝Biso+Bani/2， 假设共反射点(CRP)点的覆盖次数为N，C₁、C₂、C₃、C₄是四个未知量，则上式可以简化为： 

AX＝B                            (3) 

其中 

上式是超定线性方程组，裂缝反演问题就是求解上述线性方程组中未知量X的问题。求解出未知向量X，即可解出P，Biso，Bani及 

P＝C₁，Biso＝C₂-Bani/2， 

为了进行更深层次的油气勘探，识别裂缝内的流体类型值得做深入的研究。但相比于裂缝密度和方位的预测，对裂缝内流体的识别目前发展尚不成熟。Bakulin等(2000)在Schoenberg提出的线性滑动理论(Schoenberg，1983)和Hudson理论(Hudson，1981；1986)基础上，讨论了裂缝充填不同流体的情况下Bani参数的差异，能够为裂缝流体识别提供一定的理论依据，即： 

(1)对于干裂缝(干的或充填气体)： 

$\mathrm{Bain}=\frac{4(-8g^2+12g-3)}{3(3-2g)(1-g)}\zeta ---\left(5\right)$

(2)对于湿裂缝(充填油或水)： 

$\mathrm{Bain}=\frac{16g}{3-2g}\zeta ---\left(6\right)$

其中 即横纵波速度比的平方；ζ是Hudson理论中的裂缝体密度。 

此外Shaw和Sen(2006)推导了纵波方位AVO与裂缝柔量的关系，指出正向柔量与切向柔量之比可以作为裂缝流体类型识别的指示因子。然而，对于裂缝流体识别在实际资料中的应用目前人们的研究很少。 

发明内容

本发明提出一种基于纵波方位AVO的裂缝流体类型识别方法，该方法能够用于实际资料的裂缝流体识别中。通过本发明提供的方法，能够识别出裂缝中充填的流体类型，特别是区分干裂缝(干的或者充填气)和湿裂缝(充填油或水)，为更精细的裂缝储层描述提供依据。 

本发明实现上述目的的具体实施方案如下： 

步骤1：对地震资料进行分方位处理(包括保幅偏移处理)，抽取分方位共反射点(CRP)道集； 

步骤2：将分方位的CRP道集进行数据合并，并为下一步作准备； 

步骤3：用本发明所采用的算法进行裂缝反演，算出每个CRP点上每个时间采样点上的裂缝各向同性梯度Biso和裂缝各向异性梯度Bani(Biso和Bani两参数取自公式1中给出的Ruger反射系数方程)； 

步骤4：输出各向同性梯度Biso和各向异性梯度Bani数据体； 

步骤5：将各向同性梯度Biso和各向异性梯度Bani进行交会，结合正演信息和已知井资料勾勒出干裂缝(干的或者含气)和湿裂缝(充填油或水)的范围； 

步骤6：通过交会分析得到干裂缝(干的或者含气)和湿裂缝(充填油或水)在平面上的范围。 

步骤7：进行地质意义分析。 

本发明所提供的基于纵波方位AVO识别裂缝流体类型的基本原理如下： 

该方法在常用的基于纵波方位AVO的裂缝反演方法基础上继续识别裂缝内的流体类型，公式(4)中的各向异性反演结果对Bani参数取正值，为了应用裂缝反演的结果，在裂缝流体识别中对干裂缝和湿裂缝 的Bani参数均做取绝对值处理。利用公式(5)和公式(6)，当固定裂缝体密度ζ后，可研究Bani绝对值随背景岩石横纵波速度比的变化的规律。研究可得，当背景岩石横纵波速度比小于0.7时，充填气体时的裂缝各向异性程度Bani更高。以本发明的是实施例中研究区新疆某区块碳酸盐岩储层为例，其横纵波速度比在0.5左右，一般小于0.7，因此在该区域中当背景岩石相同、裂缝走向和密集程度等因素不变时，含气裂缝的Bani比含油和水时大。 

接下来讨论介质中充填不同流体情况下，Ruger公式中各向同性梯度Biso参数的相对变化关系。推导如下： 

$\mathrm{Biso}=\frac{1}{2}(\frac{\mathrm{\Delta \alpha }}{\bar{\alpha }}-{\left(\frac{2\bar{\beta }}{\bar{\alpha }}\right)}^2\frac{\mathrm{\Delta G}}{\bar{G}})=\frac{{\alpha }_2^2-{\alpha }_1^2-8{\bar{\beta }}^2\frac{\mathrm{\Delta G}}{\bar{G}}}{{({\alpha }_2+{\alpha }_1)}^2}---\left(7\right)$

对于固定的模型来讲，上层介质的速度α1也不变，在充填不同的流体时，横波速度不会变化，即 不会变化。唯一变化的是裂缝在含不同流体情况下介质的纵波速度α₂，因此可以通过求解Biso对α₂的导数来分析Biso的变化规律。经过求导可得： 

$\frac{\partial \mathrm{Biso}}{{\partial \alpha }_2}=\frac{2({\alpha }_1{\alpha }_2+{\alpha }_1^2+8{\beta }^2\frac{\mathrm{\Delta G}}{\bar{G}})}{{({\alpha }_1+{\alpha }_2)}^3}>0---\left(8\right)$

因此Biso随α₂的增大而增加。对于相同背景岩石，充填气体时纵波速度α₂比充填水或油时的低，因此Biso更低。这样可以理解为，在背景岩石速度、裂缝性质等因素不变的情况下，Biso在裂缝含气时比含油和水时更低，同时结合前面所讨论的含气时Bani高于含油和含水情况下的Bani值，可以对气和油水进行区分。 

本发明中裂缝探测的适用范围如下：对于无裂缝和裂缝走向呈现不同方向的井(多组裂缝)，各向异性程度弱，建立在各向异性分析基础上的流体预测方法无法预测出井中的流体类型；只有对于裂缝走向呈现某一优势方向(单组裂缝)的井，其各向异性程度较强时，才能预测其中的流体类型。 

附图说明

图1是裂缝型储层正演模型示意图。 

图2是饱含不同流体分方位道集。 

图3(a)是不加噪音、饱含流体道集反演的Bani和Biso交会图； 

图3(b)是加50％噪音、饱含流体道集反演的Bani和Biso交会图。 

图4是加50％噪音、不同含流体饱和度道集反演的Bani和Biso交会图。 

图5是基于纵波方位AVO的裂缝流体识别方法流程图。 

图6是新疆某区块Bani和Biso交会图。 

图7是新疆该区块C井反演参数交会。 

图8是新疆该区块鹰山组预测裂缝流体分布。 

图9是新疆该区块鹰山组顶部构造图。 

具体实施方式

以下结合实施例与附图说明本发明具体实施方式。 

图1～图4用数值模型正演模拟说明该裂缝流体识别方法的有效性。基于本发明实施例中新疆某区块裂缝型储层特殊地质条件进行正演模拟，所用模型如图1所示，模型共两层，上层为灰岩基质，下层为裂缝介质。背景岩石基质参数以及充填气、油和水的介质参数附表1所示。 

表1裂缝充填物参数 

  充填物及其性质>  背景介质>  气>  油>  水>  纵波速度(m/s)>  6000>  350>  1300>  1647>  横波速度(m/s)>  3100>  10>  10>  10>  密度(kg/m³)>  2670>  270>  880>  1100>

当裂缝中充填流体饱和度不同时，可以根据Wyllie时间平均方程计算等效流体的速度和密度参数。Wyllie时间平均方程如下(Wyllie，1958)： 

$\frac{1}{v}=\frac{\phi }{v_1}+\frac{1-\phi }{v_2}---\left(9\right)$

式中：v为等效流体的速度，v₁为第一种包体的速度，v₂为第二种包体的速度，φ为第一种包体占流体的(体积)百分比。根据附表1和以上公式，计算出不同流体组分情况下的纵波速度和密度，如附表2所示(注：o表示油，g表示气，w表示水)。 

表2裂缝内流体的物性参数 

   流体组分>  纵波速度(m/s)>  密度(kg/m³)>  1>  o∶g＝10∶0>  1297.0>  880.0>  2>  o∶g＝8∶2>  842.0>  758.0>  3>  o∶g＝6∶4>  623.0>  636.0>  4>  o∶g＝4∶6>  494.0>  514.0>  5>  o∶g＝2∶8>  410.0>  392.0>  6>  o∶g＝0∶10>  350.0>  270.0>  7>  g∶w＝8∶2>  415.0>  436.0>  8>  g∶w＝6∶4>  511.0>  602.0>  9>  g∶w＝4∶6>  664.0>  768.0>  10>  g∶w＝2∶8>  946.0>  934.0>  11>  g∶w＝0∶10>  1648.0>  1100.0>  12>  w∶o＝8∶2>  1563.0>  1056.0>  13>  w∶o＝6∶4>  1487.0>  1012.0>  14>  w∶o＝4∶6>  1418.0>  968.0>  15>  w∶o＝2∶8>  1355.0>  924.0>

图2是饱含不同流体情况下的合成地震记录分方位道集，其目的是观察不同流体在道集中的差别。其具体做法如下：在图1所用的数值模型基础上，根据附表1给出的饱含气、油和水的参数计算Hudson弹性参数(Hudson，1981；1986)，利用Hudson弹性参数计算Thomsen参数(Ruger和Tsvankin，1997)，将其带入Ruger公式正演得到合成地震记录分方位道集。为显示反演结果方便，分别抽取固定偏移距，方位角为0°～180°的道集，分析方位振幅变化规律。从图2可以看出，裂缝饱含气后的方位振幅变化比饱含油和水后的方位振幅变化大，说明在背景岩石纵横波速度、密度和裂缝性质不变的情况下，饱含气后的各向异性程度比饱含油和水时大。

图3是本发明的裂缝流体交会分析结果，以说明不同流体数据点的范围。利用上述正演得到的合成地震记录道集，进行反演，得到Biso和Bani参数。用不含噪音和加50％噪音后的道集分别进行反演得到Bani与Biso参数交会图，如图3(a)、(b)所示。图中浅灰色椭圆勾勒出气的范围，深灰色椭圆勾勒出油和水的范围，黑色椭圆代表不含裂缝的各向同性岩石及噪音扰动引起的气与油水难以区分的区域。从图3(a)可以看出，不加噪音时气与油水的参数点能够明显的区分；如图3(b)所示，加入50％噪音后，黑色椭圆的范围(噪音扰动引起的气与油水难以区分的区域)扩大，但整体上气与油水在该模型中能够区分。 

图4是用不同含流体饱和度的道集进行反演得到的Biso和Bani参数交会图。根据附表2给出的不同流体饱和度情况下的参数计算Hudson弹性参数(Hudson，1981；1986)，利用Hudson弹性参数计算Thomsen参数(Ruger和Tsvankin，1997)，将其带入Ruger公式正演得到合成地震记录分方位道集，进行反演得到Biso和Bani参数交会图。浅灰色椭圆指示含气饱和度大于60％的区域，深灰色椭圆指示含油或水大于60％的区域，黑色椭圆指示不含裂缝的各向同性岩石和噪音引起的气与油水难以区分的区域。从图中可以看出，整体规律是主要含气参数点落在下方，而主要含油水参数点落在上方。 

图1～图4说明利用Bani与Biso参数交会分析方法能够对干裂缝(干的或者含气)和湿裂缝(含油水)进行区分。 

图5是本发明中裂缝流体识别在实际资料应用中的流程图： 

步骤1：对地震资料进行分方位处理(包括保幅偏移处理)，抽取分方位共反射点(CRP)道集； 

步骤2：将分方位的CRP道集进行数据合并，并为下一步作准备； 

步骤3：用本发明所采用的算法进行裂缝反演，算出每个CRP点上每个时间采样点上的裂缝各向同性梯度Biso和裂缝各向异性梯度Bani； 

步骤4：输出各向同性梯度Biso和各向异性梯度Bani数据体； 

步骤5：将各向同性梯度Biso和各向异性梯度Bani进行交会，结合正演信息和已知井资料勾勒出干裂缝(干的或者含气)和湿裂缝(充填油或水)的范围； 

步骤6：通过交会分析得到干裂缝(干的或者含气)和湿裂缝(充填油或水)在平面上的范围。 

步骤7：进行地质意义分析。 

以下以该方法在新疆某区块实际资料上的应用为例子具体说明该方法的实际应用效果。 

图6是新疆该区块Bani和Biso交会图。根据井上已知信息，划定了气和油水的范围，分别如图中浅色和深色线框所示。在勾勒气与油水范围时，主要考虑到以下两点： 

(1)为避免噪音的影响，选取Bani值较大的区域。 

(2)裂缝的含流体饱和度等因素影响对气和油水的区分，为突出含气和含有油水的储层，只勾勒含气或油水饱和度较高的区域。 

图7是新疆该区块C井反演参数交会与该井生产情况对比。从试油结果上看，该井日产油10.37方，产气54505方。通过交会分析后，该井产气和产油水的位置如图7右侧黑色区域所示，上段产气，下段产油和水。由于气的密度比油和水小，在储层中分布时一般气位于油和水上方，该井流体预测结果具有一定的合理性。 

图8是新疆该区块鹰山组裂缝流体预测结果，其中底色表示反演的裂缝密度。对有成像测井资料的井，统计鹰山组深度范围内高导缝的走向和条数，并用玫瑰图的形式在对应的裂缝矢量图显示出来，玫瑰图旁注明统计的FMI资料裂缝线密度(单位长度上的裂缝条数)。反演的裂缝密度大小与裂缝的条数和裂缝走向分布密切相关。当在一定深度范围内，裂缝走向一致时(裂缝类型定义为单组)，裂缝线密度越大，地震波在裂缝介质中传播各向异性程度越强(如A井、B井、C井、D井、E井)；若裂缝走向各不相同(裂缝类型定义为多组)，地震波的各向异性程度弱，应用各向异性进行裂缝探测效果不明显(F井、I井、G井、H井)。统计流体预测结果的吻合情况。多组裂缝井和无裂缝的井各向异性程度弱(F井、I井、G井、H井)，建立在各向异性分析基础上的流体预测方法无法预测出井中的流体类型。只有对于单组裂缝的井，其各向异性程度较强时，才能预测其中的流体类型。单组裂缝的井(A井、B井、C井、D井、E井)中，A井没有预测出流体分布，除该井外，其他都与井上的试油资料吻合，单组裂缝的井流体预测吻合率80％。通过与鹰山组顶部构造图(图9)对比可以发现，预测的流体主要集中在局部的构造高部位，这种分布具有一定的合理性。整体上来说，裂缝流体预测结果合理，证明该方法具有一定的实际应用价值。 

参考文献 

马中高.P波资料反演裂缝方法及实例.石油地球物理勘探，2003，38(5)：517-521. 

Gray D.Fracture Detection Using 3D Seismic Azimuthal AVO.CSEG Recorder，2008，38-49. 

Sun S Z，Xiao X，Wang Z et al.P-wave fracture prediction algorithm using prestack data with limited azimuthal distribution.2011a，73^rd EAGE annual meeting，Expanded abstract. 

Sun S Z，Wang Z，Yang H et al.P-wave fracture prediction algorithm using pre-stack data with limited azimuth distribution：A case study in the TZ45 area，Tarim Basin，China.Petroleum Science，2011b，8(4)：422-432. 

Sun S Z，Xiao X，Chen L et al.P-wave fracture prediction algorithm using data with limited azimuthal distribution. 

The Leading Edge，2012，31(2)：198-205. 

Ruger A.Variation of P-wave reflectivity with offset and azimuth in anisotropic media.Geophysics，1998，63(3)：935-947. 

Wyllie M R J，Gregory A，Gardner G.An experimental investigation of factors affecting elastic wave velocities in porous media.Geophysics，1958. 

Bakulin A，Grechka V，Tsvankin I.Estimation of fracture parameters from reflection seismic data-Part I：HTI model due to a single fracture set.Geophysics，2000，65(6)：1788-1802. 

Schoenberg M A，1983，Reflection of elastic waves from periodically stratified media with interfacial slip：Geophysical prospecting，31：265-292. 

Shaw R K，Sen M K.Use of AVOA data to estimate fluid indicator in a vertically fractured medium.Geophysics，2006，71(3)：C15-C24. 

Hudson J A.Wave speeds and attenuation of elastic waves in material containing cracks.Geophysical Journal of the Royal Astronomical Society，1981，64：133-150. 

Hudson J A.A higher order approximation to the wave propagation constants for a cracked solid.Geophysical Journal of the Royal Astronomical Society，1986，87：265-274. 

Ruger A，Tsvankin I.Using AVO for fracture detection：Analytic basis and practical solutions.The Leading Edge，1997，1429-1434。