2017, Vol. 32
2. 中国石油大学地球物理与信息工程学院, 北京 102249
2. College of Geophysics and Information Engineering, China University of Petroleum, Beijing 102249, China
在地震勘探领域中，对地下介质的认识和模拟是开展勘探工作及后续处理解释的基础.随着油气勘探目标从简单的构造性油气藏，深入到复杂的岩性油气藏和裂缝性油气藏等，学者们发现传统的弹性介质理论，渐渐不能满足日益复杂日益困难的勘探要求.尤其是对于油气储层，较高的孔隙度，以及孔隙中水、油、气等流体的地球物理特征与固体岩石差异极大，所以考虑储层是多孔双相系统是极有必要的.对于一个流体饱和岩石系统，其体积模量在一定频率范围内 (10-2~106 Hz) 应该具有2次改变 (Clear, 1978; Sarout, 2012)，分别是岩石体积模量由松弛状态改变为未松弛状态，由饱和岩石状态改变为干岩石状态.这两次改变分别由黏弹机制和饱和岩石孔隙流体流动机制控制.
岩石黏弹机制的研究起始于19世纪40年代，Stokes (1845)建立了Stockes方程来描述波在黏滞性介质内传播时的能量消耗，随后许多学者也投入到黏弹性介质的研究中来.描述岩石黏弹性机制的模型主要有3种，早期对黏弹性介质的研究主要基于Voigt模型和Maxwell模型.随着研究的深入，Zener等 (1956)发现，Voigt模型并不能合理的描述实际情况，并提出了标准线性体 (Standard linear solid) 模型，并指出标准线性提模型更适合描述波在地下介质中衰减的物理机制，更接近真实情况.之后，Carcione (1988)在此基础上发展了基于广义Zener线性体的方法，并根据流变学原理将黏弹性介质中的应力表示为松弛函数与应变的卷积，得到了各向异性黏弹性介质的本构方程，建立了线性黏弹性各向异性介质的波动方程.
多孔岩石系统的研究起始于Gassmann方程，Gassmann (1951)首次建立了流体饱和岩石体积模量和干岩石体积模量的关系.之后，Biot (1956a, b, 1962) 发展了Gassmann方程，提出的Biot理论奠定含饱和流体孔隙介质中波传播的研究基础.随后，大量学者对Boit理论进行了发展和研究 (Plona, 1980; Berryman and Wang, 1995; Pride et al., 2004).随着研究的深入，研究者发现Biot理论不能合理的解释地震波高频散高衰减现象，Dvorkin和Nur等 (Dvorkin and Nur, 1993; Dvorkin et al., 1994) 引入了squirt喷射流动理论，很好地解释了实验中观测到的高频散和强衰减现象，建立了统一的BISQ模型，并研究了关键参数对于弹性波相速度频散和衰减的影响，有利地促进了含流体孔隙介质中波传播理论的发展.
多孔岩石系统具有复杂的孔隙结构，不同尺度上的孔隙结构都会对地震波的衰减和频散特征造成复杂的影响 (邓继新等，2015；吴建鲁等，2015)，单一的孔隙介质模型侧重于流体流动引起的衰减和频散 (吴国枕等，2014；孙卫涛等，2015)，而忽略了骨架黏弹特征带来的影响.研究表明，单一的模型都不能全面的描述在一定频率范围内地震波的衰减和频散特征.所以很多学者进行了大量研究，试图将两种机制结合，建立统一的黏弹孔隙介质模型 (Carcione and Helle, 1999).Nie等 (2012)将Voigt模型引入了BISQ模型，但相对于Voigt模型，Zener模型更加贴近地下介质黏弹响应的真实情况.本文基于Zener模型，分析了体积模量由松弛到未松弛的变化机理，及常Q衰减机制的建立；然后，我们通过将Zener线性体的本构关系引入BISQ理论，建立了统一的Zener-BISQ黏弹孔隙介质模型，给出了该模型的平面波解，并在此基础上，计算了Zener-BISQ模型的相速度衰减和频散曲线，结果表明:该模型很好的刻画了岩石体积模量由松弛到未松弛，由饱和岩石到干岩石的两次变化；最后，通过三个模型的对比试验，分析了Zener黏弹机制在Zener-BISQ模型中的作用方式和控制因素.
1 Zener线性体黏弹介质模型Zener线性体模型的本构方程可以用波尔兹曼叠加原理构建得到，其表达式为
|
(1) |
其中，TI表示应力向量，SJ表示应变向量，ΨIJ为弛豫函数，“*”表示时间的卷积.弛豫函数的表达式为
|
(2) |
其中，τεl和τσl分别是材料的第l个机制的应变弛豫时间和应力弛豫时间，表示为
|
(3) |
fl是中心频率，Ql是中心频率对应的品质因子数，L是松弛机制的总数，每组τεl和τσl表示一种衰减机制.
对公式 (2) 进行傅里叶变换，即可得到复平面波模量的表达式为
|
(4) |
其中，MR为松弛模量，通过令角频率ω→0时，计算得到.而当ω→∞时，可以得到未松弛模量的表达式为
|
(5) |
图 1给出了弛豫函数的示意图，参数为f=25 Hz, Q=12.8, MR，可以看出，未松弛的模量MU需要一定时间才能完全松弛为MR，这就说明，如果频率高到一定程度时，岩石骨架没有足够的时间回复形变，表现出未松弛时的弹性性质.图 2a, b分别给出了Zener线性体模型的相速度和逆品质因子的随频率变化曲线.如图所示，在低频时，相速度表现为岩石松弛时的
|
图 1 Zener模型弛豫函数曲线 Figure 1 Relaxation function of Zener model |
|
图 2 相速度 (a) 和逆品质因子 (b) 曲线 Figure 2 Phase velocity (a) and inverse quality factor (b) curves |
孔隙介质理论的起始是从Gassmann方程开始的，它建立了饱和岩石体积模量Ku与干岩石体积模量Kd之间的联系, 公式为
|
(6) |
其中，Ks为固体颗粒的体积模量，Kf为孔隙中流体的体积模量，φ为孔隙度，α为Biot特征系数，公式为
|
(7) |
许多学者在Gassmann理论的基础上进行了大量研究，最终建立起了同时包括Biot流动机制和Squirt喷射流动机制的BISQ模型，其本构方程可以表示为
|
(8) |
其中，P为孔隙流体压力，求取公式为
|
(9) |
其中
|
(10a) |
|
(10b) |
|
(10c) |
|
(10d) |
BISQ模型的建立过程较为复杂，而且在很多文章中都有详细的论述 (Dvorkin and Nur, 1993; Dvorkin et al., 1994)，在本文中不做重点讨论.
3 Zener-BISQ黏弹孔隙介质模型的本构方程及平面波解Zener线性体模型描述了岩石固体弛豫作用引起的，平面波模量由松弛向未松弛之间的变化，BISQ孔隙模型描述了流体流动引起的，岩石体积模量由干岩石向饱和岩石之间的变化.这两种岩石体积模量的改变都会影响地震波传播和衰减频散特征.我们通过将Zener模型和BISQ模型结合，构建了统一的Zener-BISQ模型.
Zener线性体模型描述的弛豫现象一般适用于固体部分，所以我们使用Zener模型的本构关系代替BISQ模型中的固体骨架部分的本构关系，即得到Zener-BISQ模型的本构方程为
|
(11) |
为了求取其平面波解，我们将方程 (11) 转换到频率域，得
|
(12) |
根据平面波解的假设，x方向的波场可以表示为以下形式：
|
(13a) |
|
(13b) |
|
(13c) |
将上述公式带入本构方程，并求解得到如下形式的相速度和逆品质因子的解析解为
|
(14) |
|
(15) |
其中
|
(16) |
其中
|
(17a) |
|
(17b) |
|
(17c) |
我们使用饱和含水砂岩模型来进行数值实验，其基本参数如表 1所示，而对于Zener黏弹机制，这里使用的参数为Q=20, f=100 Hz.图 3给出了快P波和慢P波的相速度和逆品质因子的频率依赖曲线，从图中可以看出，快P波相速度和逆品质因子曲线具有两个过渡带，分别由两个过渡频率控制，第一个由Zener模型控制，表示岩石体积模型由松弛状态过渡到未松弛状态；第二个过渡频率由BISQ孔隙模型控制，表示，岩石体积模量由饱和岩石状态过渡到干岩石状态.因为Zener线性体对BISQ模型的改造，仅限于固体骨架部分，所以如图所示，BISQ模型中包含的慢P波的相速度和逆品质因子曲线没有受到Zener模型影响，符合本文中Zener-BISQ模型的构建方式.
|
表 1 含水饱和砂岩模型参数 Table 1 Parameters of water saturated sandstone |
|
图 3 Zener-BISQ模型相速度和逆品质因子曲线，上两图为快P波，下两图为慢P波 Figure 3 Phase velocity and inverse quality factor curves of Zener-BISQ model. Top: fast P wave. Bottom: slow P wave |
为了进一步研究黏弹骨架参数弛豫品质因子和弛豫中心频率对模型衰减频散特征的作用方式，我选择了3组不同的参数进行数值实验，分别为Q1：Q1=20, f1=140 Hz；Q2：Q2=36, f2=1460 Hz；而模型Q3为BISQ模型，即Q3=∞.因为黏弹机制并不作用于慢P波，所以只给出了3组模型的快P波相速度及逆品质因子曲线，如图 4所示.从图中可以看出，当Q取值为无穷大时，模型Q3为退化为基本的BISQ模型，逆品质因子曲线只有一个峰值，相速度也只在一个过渡频率区间发生变化.模型Q1和Q2的逆品质因子曲线均具有两个峰值，与前文中模拟保持一致，其过渡频率由中心频率f控制，如图中所示，其过渡频率为通过中心频率计算得到的角频率2πf.其中模型Q1的过渡频率为2πf1=849 Hz，模型Q2的过渡频率为2πf2=9173 Hz.而模型Q1的弛豫品质因子为20小于模型Q2的36，所以模型Q1表现出更大的衰减特征.所以，对于Zener模型，我们可以较为直观和便捷的控制其衰减和频散机制的频率范围和衰减程度.同时，固体骨架弛豫机制的引入，较好的补充了BISQ模型在刻画地震波低频衰减现象的不足.
|
图 4 3个模型相速度和逆品质因子曲线，其中黏弹骨架参数为Q1：Q1=20, f1=140 Hz；Q2：Q2=36, f2=1460 Hz；Q3：Q3=∞ Figure 4 Phase velocity and inverse quality factor curves of three models with viscoelastic skeleton parameters: Q1: Q1=20, f1=140 Hz; Q2:Q2=36, f2=1460 Hz; Q3:Q3=∞ |
储层模型是包含固体骨架、孔隙结构和流体的复杂介质模型，固体骨架的黏弹性性质，流体的物理性质，以及不同尺度的孔隙结构引起的流动机制都会对地震波的衰减和频散产生极大的影响.所以，很难使用一个通用的模型来描述所有影响因素，在对于复杂储层的研究模拟中，要根据实际储层的岩石物理特征，填充流体类型，建模任务需求等具体情况，选择适合的模型.
5.2本文介绍的Zener-BISQ黏弹孔隙介质模型，较好的结合了固体骨架黏弹机制和饱和流体的流动机制，可以有效的模拟地震波传播经过储层时，岩石体积模量变化的两个明显的阶段，由松弛的状态变化为未松弛状态，由饱和岩石过渡到干岩石.在该模型中，这两个阶段分别由Zener特征频率和BISQ特征频率控制，可以通过改变特征频率适应不同类型储层的模拟.数值实验证明，Zener模型较好的补充了BISQ模型在解释低频衰减现象时的不足，而且其灵活性较高，只需两个参数f和Q控制.所以，Zener-BISQ模型可以灵活、准确的模拟地震波在传播经过储层时的衰减和频散特征，成为储层模拟技术新的补充.当然，地下储层的结构复杂，该模型的适用性还需要结合实际资料进一步研究.
致谢 感谢审稿专家提出的修改意见和编辑部的大力支持！| [] | Berryman J G, Wang H F. 1995. The elastic coefficients of double-porosity models for fluid transport in jointed rock[J]. Journal of Geophysical Research, 100(B12): 24611–24627. DOI:10.1029/95JB02161 |
| [] | Biot M A. 1956a. Theory of propagation of elastic waves in a fluid-saturated porous solid. Ⅰ. Low-frequency range[J]. J. Acoust. Soc. Am., 28(2): 168–178. DOI:10.1121/1.1908239 |
| [] | Biot M A. 1956b. Theory of propagation of elastic waves in a fluid-saturated porous solid. Ⅱ. higher frequency range[J]. J. Acoust. Soc. Am., 28(2): 179–191. DOI:10.1121/1.1908241 |
| [] | Biot M A. 1962. Generalized theory of acoustic propagation in porous dissipative media[J]. J. Acoust. Soc. Am., 34(5): 1254–1264. |
| [] | Carcione J M. 1988. Viscoacoustic wave propagation simulation in the earth[J]. Geophysics, 53(6): 769–777. DOI:10.1190/1.1442512 |
| [] | Carcione J M, Helle H B. 1999. Numerical solution of the poroviscoelastic wave equation on a staggered mesh[J]. Journal of Computational Physics, 154(2): 520–527. DOI:10.1006/jcph.1999.6321 |
| [] | Cleary M P. 1978. Elastic and dynamic response regimes of fluid-impregnated solids with diverse microstructures[J]. International Journal of Solids and Structures, 14(10): 795–819. DOI:10.1016/0020-7683(78)90072-0 |
| [] | Deng J X, Zhou H, Wang H, et al. 2015. The influence of pore structure in reservoir sandstone on dispersion properties of elastic waves[J]. Chinese Journal of Geophysics (in Chinese), 58(9): 3389–3400. DOI:10.6038/cjg20150931 |
| [] | Dvorkin J, Nolen-Hoeksema R, Nur A. 1994. The squirt-flow mechanism:Macroscopic description[J]. Geophysics, 59(3): 428–438. DOI:10.1190/1.1443605 |
| [] | Dvorkin J, Nur A. 1993. Dynamic poroelasticity:A unified model with the squirt and the Biot mechanisms[J]. Geophysics, 58(4): 524–533. DOI:10.1190/1.1443435 |
| [] | Gassmann F. 1951. Elastic waves through a packing of spheres[J]. Geophysics, 16(4): 673–685. DOI:10.1190/1.1437718 |
| [] | Nie J X, Ba J, Yang D H, et al. 2012. BISQ model based on a Kelvin-Voigt viscoelastic frame in a partially saturated porous medium[J]. Applied Geophysics, 9(2): 213–222. DOI:10.1007/s11770-012-0332-6 |
| [] | Plona T J. 1980. Observation of a second bulk compressional wave in a porous medium at ultrasonic frequencies[J]. Applied Physics Letters, 36(4): 259–261. DOI:10.1063/1.91445 |
| [] | Pride S R, Berryman J G, Harris J M. 2004. Seismic attenuation due to wave-induced flow[J]. Journal of Geophysical Research, 109(B1): B01201. |
| [] | Sarout J. 2012. Impact of pore space topology on permeability, cut-off frequencies and validity of wave propagation theories[J]. Geophysical Journal International, 189(1): 481–492. DOI:10.1111/gji.2012.189.issue-1 |
| [] | Stokes G G. 1845. On the theories of the internal friction of fluids in motion, and of the equilibrium and motion of elastic solids[J]. Trans. Camb. Philos. Soc., 8: 287–305. |
| [] | Sun W T, Liu J W, Ba J, et al. 2015. Theoretical models of elastic wave dispersion-attenuation in porous medium[J]. Progress in Geophysics (in Chinese), 30(2): 586–600. DOI:10.6038/pg20150215 |
| [] | Wu G C, Wu J L, Zong Z Y. 2014. The attenuation of P wave in a periodic layered porous media containing cracks[J]. Chinese Journal of Geophysics (in Chinese), 57(8): 2666–2677. DOI:10.6038/cjg20140825 |
| [] | Wu J L, Wu G C, Zong Z Y. 2015. Attenuation of P waves in a porous medium containing various cracks[J]. Chinese Journal of Geophysics (in Chinese), 58(4): 1378–1389. DOI:10.6038/cjg20150424 |
| [] | Zener C M, Siegel S. 1949. Elasticity and anelasticity of metals[J]. Journal of Physical Chemistry, 53(9): 1468. DOI:10.1021/j150474a017 |
| [] | 邓继新, 周浩, 王欢, 等. 2015. 基于储层砂岩微观孔隙结构特征的弹性波频散响应分析[J]. 地球物理学报, 58(9): 3389–3400. DOI:10.6038/cjg20150931 |
| [] | 孙卫涛, 刘嘉玮, 巴晶, 等. 2015. 孔隙介质弹性波频散一衰减理论模型[J]. 地球物理学进展, 30(2): 586–600. DOI:10.6038/pg20150215 |
| [] | 吴国忱, 吴建鲁, 宗兆云. 2014. 周期性层状含孔隙、裂隙介质模型纵波衰减特征[J]. 地球物理学报, 57(8): 2666–2677. DOI:10.6038/cjg20140825 |
| [] | 吴建鲁, 吴国忱, 宗兆云. 2015. 含混合裂隙、孔隙介质的纵波衰减规律研究[J]. 地球物理学报, 58(4): 1378–1389. DOI:10.6038/cjg20150424 |