1 引言

由于页岩气储层的复杂性，目前还不存在能够广泛适用的能精确描述其主要储层性质及速度之间关系的岩石物理模型.自洽模型通过计算将包含物置于无限大背景材料产生的形变所引起的应变能来估计岩石的体积模量与剪切模量(Budiansky，1965)，临界孔隙度的存在使得自洽模型的应用受到一定的限制，因为对于固体相和流体相(孔隙相)组成的混合岩石，当流体相的体积百分比大于60%时，自洽模型计算的剪切模量趋向于0(Berryman，1980)，意味着当流体相的体积百分比大于60%时，固体相失去了连续性.Agnibha和Michael(2009)认为在无限大背景条件下由于自洽理论中固体相和流体相具有对称性，因此当流体相的体积百分比在40%到60%之间时，固体相和流体相组成的有效介质是双联通的.为了解决自洽模型受限的问题，通常的做法是将自洽模型与微分等效介质模型相结合(Hornby et al., 1994；Agnibha and Michael， 2009).本文在自洽模型和微分等效介质模型的基础上，引入Berryman三维孔隙形态，模拟不同孔隙形状对于自洽模型临界孔隙度与岩石速度的影响.在泥页岩中有机质固体充填的思想基础上，利用Brown-Korringa方程(Brown and Korringa， 1975)将有机质作为固体充填物，模拟不同泥质含量条件下，不同孔隙形态中填充的有机质对岩石弹性性质的影响.最后，针对中国南方某地区的一口页岩气井，利用建立的新的岩石物理模型进行速度预测，预测结果与测井曲线测量的结果吻合的很好，验证了该模型的适用性.

富有机质SCA_DEM方程理论是基于Budiansky-Hill(1965)自洽理论(简称SCA)、Zimmerman(1991)微分等效介质理论(简称DEM)、Gassmann方程(1951)、Brown-Korringa方程(1975)、Wu(1966)任意孔隙纵横比、及Berryman(1995)三维特殊孔隙的划分构建成的岩石物理模型.

Budiansky(1965)与Hill(1965)提出了自洽模型(Self-consistent approximation，SCA)，其基本建模思想如下：将要求解的多相介质放置于无限大的背景介质中，该背景介质的弹性参数是任意可调的.通过调整背景介质的弹性参数，使得可调节背景介质弹性参数与多相介质的弹性参数相匹配，当有一平面波入射时，多相介质不再引起散射，此时背景介质的弹性模量与多相介质的有效弹性模量相等.该方法既考虑到孔隙形状的影响，又能够适用于孔隙度较大的岩石.这种方法仍然是计算的包裹物的变形，但是该方法中不再选用多相材料中的一相作为背景介质，而是用要求解的有效介质作为背景介质，通过不断改变基质来考虑包裹物之间的相互作用.因为该方法考虑了包裹物的相互作用，所以能适用于孔隙度较大的岩石(原宏壮，2007).计算公式如下(Wu，1966)：

(1)Wu的任意孔隙纵横比取值方式

设定孔隙纵横比a的初值，当结果与实测结果不匹配时需要不断地调整a的值.孔隙纵横比a是速度预测的重要参数，一般a在0～1之间，越接近于0表示孔隙越扁(裂缝)，接近于1表示圆形孔(孔洞).

(2)Berryman的三维特殊孔隙取值方式

Berryman提出将实际岩石孔隙等效模拟成球形孔隙，针形孔隙，硬币状孔隙，裂缝型孔隙.各种孔隙几何尺寸及三维孔隙分别在砂质和灰质孔隙，泥质孔隙中所占的比例系数如表 1所示.本研究岩石物理模型只选取四种特殊三维孔隙形状中的三个，即球形(表征溶蚀洞)、针形(表征溶蚀孔)、硬币状裂 缝(表征缝).由于DEM理论假设加入的孔隙是有 次序的，因此通过相关资料可以确定孔隙类型主次的情况下，依照这个次序分别加入孔隙，如果无法确定孔隙类型主次，则可以先根据划分的孔隙度大小 排列一个主次顺序，然后按照这个次序迭代求算加入所有孔隙后的干岩石的弹性模量.此时公式(1)和(2)中的P和Q随着加入的孔隙类型而作相应的变化.

表 1(Table 1)

表 1 四种三维孔隙的几何尺寸常量P、Q及其在砂质与灰质、泥质孔隙中的比例系数(Mavko et al., 1998)Table 1 Constants P，Q of Four three-dimensional pore geometry and scale factors of s and y，lime，clay pore

表 1 四种三维孔隙的几何尺寸常量P、Q及其在砂质与灰质、泥质孔隙中的比例系数(Mavko et al., 1998)Table 1 Constants P，Q of Four three-dimensional pore geometry and scale factors of s and y，lime，clay pore

2.2 微分等效介质(DEM)理论

微分等效介质(Differential Effective Medium，DEM)岩石物理模型通过往固体相中逐渐加入填入物来模拟双相混合物(Norris et al., 1985；Zimmerman，1991)，固体矿物是相1，之后逐步加入相2的材料.此过程一直进行到需要的各成分含量达到为止.DEM理论并不是对称的对待每个组成成分，被当成固体矿物或主相的成分可以有不同的选择，且最终的等效模量会依赖于达到最终混合物所采用的路径.用相1作为主相并逐渐加入材料2，与以相2作为主相并逐渐加入材料1，会导致不同的等效属性.公式可以表示为

岩石物理分析中的一个重要问题就是从一种流体饱和的岩石地震速度预测另一种流体饱和的岩石地震速度，即用岩石骨架速度预测饱和岩石速度，反之亦然，这就是流体替换，而流体替换的基础就是Gassmann方程.

Gassmann(1951)提出了饱和流体岩石的弹性模量公式：

由此可得，纵波速度与岩石固体颗粒的弹性模量，流体体积模量，干燥岩石体积模量及孔隙度、岩石密度的关系：

Gassmann方程的基本假设是：(1)岩石(基质和骨架)宏观上是均匀各向同性的；(2)所有的孔隙都是连通的；(3)孔隙中充满着流体；(4)研究中的岩石—流体系统是封闭的(不排液)；(5)当波在岩石中 传播时，流体和骨架之间的相对运动可以忽略；(6)孔隙流体不对固体骨架产生软化或硬化作用.Gassmann 方程没有考虑孔隙形状的变化，适用于低频情况，当频率足够低，使得孔隙流体有足够的时间流动并没有波动诱发产生孔隙压力梯度时才成立，这也就说明了该方程非常适用于地震资料频带(<100 Hz).

Gassmann方程中干岩样的体积模量，一般是未知的，需要实验室测量得到，也有很多人提出了一些简单的计算公式，如Biot(1956)提出了K_d与K_m之间的关系：

Brown-Korringa(1975)利用岩体弹性张量表示流体替代的理论，也称为各向异性的流体替代.经过简单的改造则可以用来进行固体替代计算含干酪根岩石的体积模量和剪切模量，即，

3 三维 SCA_DEM岩石物理模型建模流程

本研究在自洽模型和微分等效介质模型的基础上，引入Berryman三维孔隙形态，模拟不同孔隙形状对于自洽模型临界孔隙度与岩石速度的影响；在泥页岩中有机质固体充填的思想基础上，利用Brown-Korringa方程将有机质作为固体充填物，图 1为富有机质泥页岩三维SCA_DEM岩石物理建模流程图，其主要流程如下所述：

图 1

Fig. 1

图 1 富有机质泥页岩三维SCA_DEM模型建模流程图Fig. 1 Schematic view of 3D SCA_DEM rock physics model for organic-rich shale

(1)给定临界孔隙度，利用SCA理论计算临界孔隙度时干岩石的体积模量和剪切模量；临界孔隙度是指岩石中矿物由承载域转变为悬浮域的孔隙度值，当孔隙度低于这个值时，矿物承载，而当孔隙度大于这个值时，矿物散开并且变成承载液体中的悬浮物.临界孔隙度的值取决于沉积过程中颗粒的分选和磨圆，国外一些学者通过测试得到了一些岩石的临界孔隙度(Mavko et al., 1998)，在计算中，大都是取这些值.

(2)利用DEM理论逐步调整孔隙度至真实孔隙度时的干岩石的体积模量和剪切模量；

(3)利用Gassmann方程计算饱和流体岩石的体积模量和剪切模量；

(4)利用Brown-Korringa方程进行固体替代计算含干酪根岩石的体积模量和剪切模量；

(5)最后再利用纵横波速度与弹性模量(体积模量K和剪切模量μ)、密度(ρ)之间的关系式求得最终的饱和流体岩石的纵波速度V_P和横波速度V_S.

4 基于三维 SCA_DEM岩石物理模型正演分析

图 2和图 3显示了利用SCA_DEM模型模拟不同孔隙形态条件下水饱和砂岩的纵波速度与横波速度随孔隙度的变化趋势.图 2显示了利用二维SCA_DEM模型得到的纵横波速度随孔隙度的变化趋势，其利用Wu的任意二维孔隙纵横比模拟孔隙形态的变化.从图中我们可以看到，随着孔隙度增大，纵横波速度降低；随着孔隙纵横比减小，纵横波速度降低，特别是孔隙纵横比在0.1~0.2之间时，速度的变化非常明显.图 3显示了利用三维SCA_DEM模型得到的纵横波速度随孔隙度的变化趋势，其利用Berryman提出的三维特殊孔隙形态来模拟孔隙形状的变化(球形孔隙、针形孔隙、硬币状裂缝孔隙)，从图中可以看出，随着孔隙度增加，球形孔隙条件下的纵横波速度下降的程度要比针形孔隙及硬币状裂缝孔隙要低.同时与二维情况下类似，硬币状裂缝的纵横比落在0.1~0.2之间时，可以看到纵横 波速度有比较明显的不可忽视的变化，表明不同孔 隙形状对于岩石速度的影响是非常重要的.

图 2

Fig. 2

图 2 基于Wu的任意纵横比的2D SCA_DEM模型的水饱和砂岩的纵横波速度随孔隙度的变化Fig. 2 Forward modeling of 2D SCA_DEM with Wu′s arbitrary aspect ratio

图 3

Fig. 3

图 3 基于Berryman的特殊孔隙形态的三维SCA_DEM模型的水饱和砂岩的纵横波速度随孔隙度的变化Fig. 3 Forward modeling of 3D SCA_DEM rock physics model with Berryman′s 3D special pore geometry

对于固体相和流体相(孔隙相)组成的混合岩石，当流体相的体积百分比大于60%时，自洽模型计算的剪切模量趋向于0(Berryman，1980)，意味着当流体相的体积百分比大于60%时，固体相失去了连续性.Agnibha和Michael(2009)提出由于自洽理论中在无限大背景条件下，固体相和流体相具有对称性，因此当流体相的体积百分比在40%到60%之间时，固体相和流体相组成的有效介质是双联通的.基于Berryman的三维特殊孔隙形态，我们利用三维SCA_DEM模型模拟不同孔隙形状对应的临界孔隙度.图 4a显示对于球形孔隙，随着临界孔隙度的增大，SCA_DEM模型计算得到的纵横波速度降低，特别是当临界孔隙度达到0.5时，纵波与横波都发生明显的扭曲，即其计算的剪切模量趋向于0.对于针形孔隙，其情况与球形孔隙类似，但其变化梯度稍有不同(如图 4b所示).然而对于硬币状裂缝而言，特别是其裂缝纵横比在0.1~0.2之间时，情况变得相当的复杂.图 4c表明裂缝纵横比等于0.1时，当临界孔隙度达到0.35时，纵波速度与横波速度都开始扭曲变形.而当裂缝纵横比达到0.2时，纵波速度与横波速度的变化趋势几乎与针形孔隙变化趋势一致.

图 4

Fig. 4

图 4 不同孔隙类型条件下的临界孔隙度分析(A)球形孔隙；(B)针形孔隙；(C)硬币状裂缝孔隙(a=0.1)；(D)硬币状裂缝孔隙(a=0.2).Fig. 4 Different critical porosity for different pore geometry(A)Spherical pores;(B)Needle-shaped pores;(C)Penny-shaped crack(a=0.1);(D)Crack(a=0.2).

5 基于三维 SCA_DEM模型的有机质模拟分析及速度预测

(1)有机质正演分析

基于正演分析的结果，利用建立的三维SCA_DEM模型模拟有机质对岩石弹性性质的影响.首先在固体背景条件下建立不同孔隙形状充填物的干岩石，然后利用Gassmann方程进行流体替代，最后利用Brown-Korringa方程进行固体替代将有机质充填到孔隙中.设计富有机质的石英-黏土组成的泥页岩，固定其孔隙度为15%，且为100%水饱和(泥质体积模量K_clay=17 GPa，泥质剪切模量μ_clay=8 GPa，石英体积模量K_quartz=36 GPa，石英剪切模量μ_quartz=44 GPa).图 5为球形孔隙条件下，有机质对于岩石弹性性质的影响的交会图.从图中可以看 出，当黏土含量较低的情况下，岩石显示低纵横波速度 比(V_P/V_S)，高杨氏模量(E)以及低泊松比(ν)的特征；黏土含量较高的情况下，有机质对于纵横波速度比，杨氏模量，泊松比的影响降低.图 6为硬币状裂缝孔隙条件下(孔隙纵横比=0.1)，有机质对于岩石弹性性质的影响的交会图，从图中可以看到，当黏土含量较低的情况下，与球形孔的情况类似，但影响幅度略显不同；当黏土含量很高的条件下，岩石显示出高纵横波速度比，低杨氏模量，高泊松比的特征.并且对于纵横波速度比与泊松比，在有机质含量达到一定程度时，黏土含量对于岩石弹性性质的影响会发生突变.通过以上的分析表明，有机质对岩石弹性性质会产生一定的影响，但这种影响敏感程度是黏土含量依赖与孔隙类型依赖的.有机质对于低黏土含量时的岩石弹性特征具有明显的影响，并且这种影响是孔隙类型依赖的；对于高黏土含量时，针对球形孔隙，固体有机质充填物一定程度上被黏土的影响所掩盖；但是对于硬币状裂缝孔隙，其影响是不容忽视的.

图 5

Fig. 5

图 5 有机质对石英-黏土组成的岩石弹性性质的影响孔隙类型为球形孔，φ=15%，水饱和度=100%，Kclay=17 GPa，μclay=8 GPa，Kquartz=36 GPa，μquartz=44 GPa.(a)总碳含量(TOC)与纵横波速度比(VP/VS)；(b)总碳含量(TOC)与杨氏模量(E)；(c)总碳含量(TOC)与泊松比(ν)；(d)纵波阻抗(PI)与纵横波速度比(VP/VS)；(e)杨氏模量(E)与泊松比(ν)；(f)第一拉梅常数(λρ)与第二拉梅常数(μρ)Fig. 5 Model prediction of TOC effect on a shale with quartz-clay mineralogy and spherical pore geometry(Porosity=15% and Water saturation=100%)

图 6

Fig. 6

图 6 有机质对石英-黏土组成的岩石弹性性质的影响孔隙类型为硬币状裂缝(纵横比a=0.1); φ=15%，水饱和度=100%;Kclay=17 GPa，μclay=8 GPa，Kquartz=36 GPa，μquartz=44 GPa.(a)总碳含量(TOC)与纵横波速度比(VP/VS)；(b)总碳含量(TOC)与杨氏模量(E)；(c)总碳含量(TOC)与泊松比(ν)；(d)纵波阻抗(PI)与纵横波速度比(VP/VS)；(e)杨氏模量(E)与泊松比(ν)；(f)第一拉梅常数(λρ)与第二拉梅常数(μρ)Fig. 6 Model prediction of TOC effect on a shale with quartz-clay mineralogy and crack(aspect ratio=0.1)geometry，Porosity=15%， and Water saturation=100%

(2)某页岩气井速度预测结果分析

图 7是利用三维SCA_DEM模型对中国南方的一口页岩气井进行纵横波速度预测成果与实际测井曲线测量结果的对比，讨论该模型在实际页岩气储层中的预测效果.从图中可以看出，预测的纵横波时差三维SCA_DEM与实测的纵横波时差AC、DTS达到了很好的吻合，这就保证了本文提出的三维SCA_DEM模型的适用性及可靠性.因此，这样建立在精细测井分析和岩石物理模型预测技术上得到的纵横波速度数据就为在无横波测井区同样实现叠前反演、时移地震等需要纵、横波信息联合的储层预测或监测手段带来了可能性.

图 7

Fig. 7

图 7 基于三维SCA_DEM模型的纵横波速度预测结果与实测纵横波对比结果(以中国南方某页岩气井为例，左侧为实际测井曲线，右侧为岩性解释剖面及纵横波预测结果，1 ft=0.3048 m)Fig. 7 Comparison between P-wave and S-wave Slowness predicted by 3D SCA_DEM model and those measured by well logging

6 结论

通过三维SCA_DEM岩石物理建模及其正演以及实际应用效果分析，我们可以得到如下结论：

(1)本文提出的三维SCA_DEM岩石物理模型对于富有机质泥页岩速度预测是有效的、可靠的；

(2)SCA理论中虽然将未知的复合岩石作为背景，但其与DEM结合时，固体相和流体相仍然是不 对称的，因此临界孔隙度不一定必须落在0.4~0.6 之间.临界孔隙度对于纵横波速度预测的结果有明显的影响.不同的孔隙类型(球形、针形、硬币状裂缝)对应着不同的临界孔隙度，因此在岩石物理速度预测及岩石物理模型分析过程中，临界孔隙度必须根据岩性、孔隙类型具体选择；

(3)有机质对于岩石弹性性质具有一定的影响，但其敏感度取决于岩石中的黏土含量及孔隙类型.