﻿ 硬币型裂缝介质的频散与衰减
 石油地球物理勘探  2019, Vol. 54 Issue (4): 836-843  DOI: 10.13810/j.cnki.issn.1000-7210.2019.04.014 0
 文章快速检索 高级检索

### 引用本文

ZHANG Fanchang, SANG Kaiheng, LU Yawei. Frequency dispersion and attenuation of porous rocks with penny shaped fractures. Oil Geophysical Prospecting, 2019, 54(4): 836-843. DOI: 10.13810/j.cnki.issn.1000-7210.2019.04.014.

### 文章历史

Frequency dispersion and attenuation of porous rocks with penny shaped fractures
ZHANG Fanchang , SANG Kaiheng , LU Yawei
School of Geosciences, China University of Petroleum(East China), Qingdao, Shandong 266580, China
Abstract: Conventional integral methods usually have the singularity problem in solving frequency-dependent elastic modulus of penny shaped fractures, and can only obtain a single normal frequency-dependent modulus.In order to quickly and accurately perform numerical simulation of the penny shaped fracture, based on the Gauss-Lobatto discrete integral, a new approach is proposed.First the second Fredholm integral equation is transformed into a discrete integral in a finite interval.Then, the solving the zero-limit of far field multiple scatter equation and the discrete integral procedure are unified with the high-order approximation.To obtain the frequency-dependent elastic tensor, analytic elastic moduli equations at the high and low frequency limits are respectively obtained with the anisotropic Gassmann equation, the linear slip theory, and the Hudson crack model.Based on the analysis of phase velocity and viscoelastic reflection coefficient, the following understanding are obtained:A.The greater the fracture denty is, the graeter the attenuated peak amplitude is; B.The attenuated peak amplitude has lower frequency when the fracture size or the fluid viscosity becomes greater; C.The incident velocity and attenuation in the vertical direction have greater change than that in the horizontal direction; D.The reflection coefficient of gas-bearing media is significantly larger than that of oil-bearing or water-bearing media.
Keywords: penny shaped fracture    discrete integral    anisotropic Gassmann equation    linear slip theory    Hudson crack model
0 引言

1 硬币型裂缝孔隙弹性理论

 $B(x)+\frac{1}{\pi} \int_{0}^{\infty} R(x, y) T(y) B(y) \mathrm{d} y=-p_{0} S(x)$ (1)

 $\begin{array}{l} T(y) = M\left\{ {\frac{{{{\left( {2\alpha g{y^2} - k_2^2} \right)}^2}}}{{2gyHk_2^2(g - 1)\sqrt {{y^2} - k_2^2} }} - } \right.\\ \;\;\;\;\;\;\;\;\left. {\frac{{2yag\sqrt {{y^2} - k_2^2} \left[ {(\alpha g - 2)k_2^2 + 2\alpha g{y^2}} \right]}}{{2gyHk_2^2(g - 1)\sqrt {{y^2} - k_2^2} }}} \right\} \end{array}$ (2)

 ${p_{33}} = \frac{{\rho {\omega ^2}}}{{k_1^2{{\left[ {1 - \frac{{2{\rm{i}}\varepsilon (H - \alpha M)}}{{{a^3}{k_1}\mu H(1 - g)}}\mathop {\lim }\limits_{x \to 0} \frac{{B(x)}}{x}} \right]}^2}}}$ (3)

2 数值求解方法

 $\theta (\zeta ) + \frac{1}{{\rm{ \mathsf{ π} }}}\int_0^a \mathit{\Gamma } (\zeta ,\xi )\theta (\xi ){\rm{d}}\xi = - {p_0}\zeta$ (4)

 $\mathit{\Gamma }(\zeta ,\xi ) = {\rm{ \mathsf{ π} }}\sqrt {\zeta \xi } \int_0^\infty y T(y){J_{\frac{1}{2}}}(\zeta y){J_{\frac{1}{2}}}(\xi y){\rm{d}}y$ (5)

 $B(x)=\frac{2}{\pi} \int_{0}^{a} \theta(\zeta) \sin (\xi x) \mathrm{d} \xi$ (6)

 ${J_{\frac{1}{2}}}(x) = \sqrt {\frac{2}{{{\rm{ \mathsf{ π} }}x}}} \sin x$

 $\mathit{\Gamma }(\zeta ,\xi ) = 2\int_0^\infty T (y)\sin (\zeta y)\sin (\xi y){\rm{d}}y$ (7)

 $\int_{0}^{a} f(x) \mathrm{d} x \approx \sum\limits_{j=1}^{n} A_{f} f\left(x_{j}\right)$ (8)

 $\theta (x) + \frac{1}{\pi }\sum\limits_{j = 1}^n {{A_j}} \mathit{\Gamma }\left( {x,{x_j}} \right)\theta \left( {{x_j}} \right) = F(x)$ (9)

x=x1, …, xn, 并记θi=θ(xi), Γij=Γ(xi, xj), Fi=F(xi), 则式(9)可写为

 ${\theta _i} + \frac{1}{\pi }\sum\limits_{j = 1}^n {{A_j}} {\mathit{\Gamma }_{ij}}{\theta _j} = {F_i}$ (10)

 $\left[ {\begin{array}{*{20}{c}} {{\theta _1}}\\ {{\theta _2}}\\ \vdots \\ {{\theta _n}} \end{array}} \right] + \frac{1}{{\rm{ \mathsf{ π} }}}\left[ {\begin{array}{*{20}{c}} {{A_1}{\mathit{\Gamma }_{11}}}&{{A_2}{\mathit{\Gamma }_{12}}}& \cdots &{{A_n}{\mathit{\Gamma }_{1n}}}\\ {{A_1}{\mathit{\Gamma }_{21}}}&{{A_2}{\mathit{\Gamma }_{22}}}& \cdots &{{A_n}{\mathit{\Gamma }_{2n}}}\\ \vdots & \vdots &{}& \vdots \\ {{A_1}{\mathit{\Gamma }_{n1}}}&{{A_2}{\mathit{\Gamma }_{n2}}}& \cdots &{{A_n}{\mathit{\Gamma }_{nn}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{\theta _1}}\\ {{\theta _2}}\\ \vdots \\ {{\theta _n}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{F_1}}\\ {{F_2}}\\ \vdots \\ {{F_n}} \end{array}} \right]$ (11)

 图 1 当频率为10(a)、100(b)、10000Hz(c)时积分函数随lgy的变化曲线

 $\lim\limits _{x \rightarrow 0} \frac{B(x)}{x}=\lim\limits_{x \rightarrow 0} \frac{\int_{0}^{a} \theta(\xi) \sin (\xi x) \mathrm{d} \xi}{x}=\int_{0}^{a} \xi \theta(\xi) \mathrm{d} \xi$ (12)

 $\int_{0}^{a} x \theta(x) \mathrm{d} x \approx \sum\limits_{j=1}^{n} A_{j} x_{j} \theta\left(x_{j}\right)$ (13)

3 硬币型裂缝介质频变弹性张量

Galvin等[21]只给出裂缝法向的弹性模量，其他方向的弹性模量可以通过Krzikalla等[25]推导的频变弹性模量关系得到

 $p_{i j}(\omega)=C_{i j}^{\operatorname{high}}-\frac{\left[p_{33}(\omega)-C_{33}^{\mathrm{high}}\right]\left(C_{i j}^{\mathrm{high}}-C_{i j}^{\mathrm{low}}\right)}{C_{33}^{\mathrm{low}}-C_{33}^{\mathrm{high}}}$ (14)

3.1 低频极限

 $\Delta_{\mathrm{N}}=\frac{4 \varepsilon}{3 g(1-g)}$ (15)
 $\Delta_{\mathrm{T}}=\frac{16 \varepsilon}{3(3-2 g)}$ (16)

 $\left\{ \begin{array}{l} C_{11}^{{\rm{low}}} = \left( {\lambda + 2\mu } \right)\left[ {1 - \frac{{{\lambda ^2}\Delta_{\rm{N}}}}{{{{\left( {\lambda + 2\mu } \right)}^2}}}} \right] + \\ \;\;\;\;\;\;\;\;\;\frac{{{{\left( {{K_{\rm{g}}} - K + \frac{{\lambda {\Delta _{\rm{N}}}K}}{{\lambda + 2\mu }}} \right)}^2}}}{{\frac{{{K_{\rm{g}}}\phi \left( {{K_{\rm{g}}} - {K_{\rm{f}}}} \right)}}{{{K_{\rm{f}}}}} + {K_{\rm{g}}} - K\left( {1 - \frac{{{\Delta _{\rm{N}}}K}}{{\lambda + 2\mu }}} \right)}}\\ C_{13}^{{\rm{low}}} = \lambda \left( {1 - {\Delta _{\rm{N}}}} \right) + \\ \;\;\;\;\;\;\;\;\;\frac{{\left( {{K_{\rm{g}}} - K + {\Delta _{\rm{N}}}K} \right)\left( {{K_{\rm{g}}} - K + \frac{{\lambda {\Delta _{\rm{N}}}K}}{{\lambda + 2\mu }}} \right)}}{{\frac{{{K_{\rm{g}}}\phi \left( {{K_{\rm{g}}} - {K_{\rm{f}}}} \right)}}{{{K_{\rm{f}}}}} + \left( {{K_{\rm{g}}} - K + \frac{{{\Delta _{\rm{N}}}K}}{{\lambda + 2\mu }}} \right)}}\\ C_{33}^{{\rm{low}}} = H\left[ {1 - \frac{{4\varepsilon {{\left( {H - \alpha M} \right)}^2}}}{{3\mu H\left( {1 - g} \right)}}} \right]\\ C_{44}^{{\rm{low}}} = \mu \left( {1 - {\Delta _{\rm{T}}}} \right)\\ C_{66}^{{\rm{low}}} = \mu \end{array} \right.$ (17)

3.2 高频极限

 $\Delta_{\mathrm{N}}=\frac{4 \varepsilon}{3 g(1-g)\left[1+\frac{1}{\pi g(1-g)}\left(\frac{K_{\mathrm{f}}+\frac{4}{3} \mu_{\mathrm{f}}}{\mu}\right)\left(\frac{a}{c}\right)\right]}$ (18)

 $\left\{ \begin{array}{l} C_{11}^{{\rm{high}}} = H - \frac{{{{\left( {H - 2\mu } \right)}^2}}}{H}{\Delta _{\rm{N}}}\\ C_{33}^{{\rm{high}}} = H\left( {1 - {\Delta _{\rm{N}}}} \right)\\ C_{13}^{{\rm{high}}} = C_{33}^{{\rm{high}}}\left( {1 - 2\frac{\mu }{H}} \right)\\ C_{44}^{{\rm{high}}} = \mu \left( {1 - {\Delta _{\rm{T}}}} \right)\\ C_{66}^{{\rm{high}}} = \mu \end{array} \right.$ (19)
4 黏弹性VTI介质的相速度和反射透射系数

VTI介质中复相速度VP的频散表达式为[31]

 $\begin{array}{l} {V_{\rm{P}}}\left( {\theta ,\omega } \right) = \frac{1}{{\sqrt {2\rho } }} \times \\ \sqrt {\left( {{p_{11}} + {p_{44}}} \right){{\sin }^2}\theta + \left( {{p_{33}} + {p_{44}}} \right){{\cos }^2}\theta + D} \end{array}$ (20)

 $\begin{array}{l} D = \left\{ {{{\left[ {\left( {{p_{11}} - {p_{44}}} \right){{\sin }^2}\theta - \left( {{p_{33}} - {p_{44}}} \right){{\cos }^2}\theta } \right]}^2}} \right.\\ \;\;\;\;\;\; + {\left( {{p_{13}} + {p_{44}}} \right)^2}\sin 2\theta {\} ^{\frac{1}{2}}} \end{array}$

 $\left[ {\begin{array}{*{20}{c}} {{\beta _{{\rm{P1}}}}}&{{\beta _{{\rm{S1}}}}}&{ - {\beta _{{\rm{P2}}}}}&{ - {\beta _{{\rm{S2}}}}}\\ {{\gamma _{{\rm{P1}}}}}&{{\gamma _{{\rm{S1}}}}}&{{\gamma _{{\rm{P2}}}}}&{{\gamma _{{\rm{S2}}}}}\\ {{Z_{{\rm{P1}}}}}&{{Z_{{\rm{S1}}}}}&{ - {Z_{{\rm{P2}}}}}&{ - {Z_{{\rm{S2}}}}}\\ {{W_{{\rm{P1}}}}}&{{W_{{\rm{S1}}}}}&{{W_{{\rm{P2}}}}}&{{W_{{\rm{S2}}}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{R_{{\rm{PP}}}}}\\ {{R_{{\rm{PS}}}}}\\ {{T_{{\rm{PP}}}}}\\ {{T_{{\rm{PS}}}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} { - {\beta _{{\rm{P1}}}}}\\ {{\gamma _{{\rm{P1}}}}}\\ { - {Z_{{\rm{P1}}}}}\\ {{W_{{\rm{P1}}}}} \end{array}} \right]$ (21)

 $\beta = {\left[ {\frac{{{p_{44}}s_x^2 + {p_{33}}s_z^2 - \rho }}{{{p_{11}}s_x^2 + {p_{33}}s_z^2 + {p_{44}}\left( {s_x^2 + s_z^2} \right) - 2\rho }}} \right]^{\frac{1}{2}}}$ (22)
 $\gamma = \pm {\left[ {\frac{{{p_{11}}s_x^2 + {p_{44}}s_z^2 - \rho }}{{{p_{11}}s_x^2 + {p_{33}}s_z^2 + {p_{44}}\left( {s_x^2 + s_z^2} \right) - 2\rho }}} \right]^{\frac{1}{2}}}$ (23)
 $W=p_{55}\left(\gamma_{S_{x}}+\beta s_{z}\right)$ (24)
 $Z=\beta p_{13} s_{x}+\gamma p_{33} s_{z}$ (25)
 $s_{z}=\frac{\left(K_{1} \mp \sqrt{K_{1}^{2}-4 K_{2} K_{3}}\right)^{\frac{1}{2}}}{\sqrt{2}}$ (26)

 $K_{1}=\frac{\rho}{p_{44}}+\frac{\rho}{p_{33}}+\left[\frac{p_{13}}{p_{33} p_{44}}\left(p_{13}+p_{55}\right)-\frac{p_{11}}{p_{44}}\right] s_{x}^{2}$ (27)
 $K_{2}=\frac{p_{11} s_{x}^{2}-\rho}{p_{33}}$ (28)
 $K_{3}=s_{x}^{2}-\frac{\rho}{p_{55}}$ (29)
5 数值模拟结果分析

 图 2 不同裂缝密度(上)、裂缝长度(下)的硬币型裂缝的频散(左)与衰减(右)

 图 3 不同频率的含水介质速度(左)和衰减(右)随入射角的变化

 图 4 不同入射角的含水介质频散(左)和衰减(右)

 图 5 含气介质的速度(左)与能量(右)随入射角和频率的变化

 图 6 含油介质的速度(左)与能量(右)随入射角和频率的变化

 图 7 含水(上)、含油(中)、含气(下)介质频变反射系数(左)与流度(右) 裂缝介质上覆高速、高密度各向同性介质层，密度为2700kg/m3，纵、横波速度分别为3000、1600m/s
6 结论

(1) 通过离散积分和高阶近似方法，将第二类Fredholm积分方程与多重散射方程的远场极限值求解相统一，为数值模拟硬币型裂缝模型的法向频变弹性模量提供快速、准确的方法；结合线性滑移理论和Hudson裂缝模型，得到了频变弹性模量的解析式。

(2) 分析VTI介质相速度与黏弹性反射、透射系数方程可知，裂缝长度决定弛豫向非弛豫过渡区的频带，裂缝密度影响衰减峰值；当裂缝介质中饱含不同流体时，流体的黏度成为决定弛豫向非弛豫过渡的又一因素，饱气时反射系数幅值大于饱油、饱水时，流度的变化趋势与衰减的变化趋势相似。

 [1] 吴国忱, 赵小龙, 罗辑, 等. 基于扰动弹性阻抗的裂缝参数反演方法[J]. 石油地球物理勘探, 2017, 52(2): 340-349. WU Guochen, ZHAO Xiaolong, LUO Ji, et al. Fracture parameter inversion based on perturbation elastic impedence[J]. Oil Geophysical Prospecting, 2017, 52(2): 340-349. [2] 龚明平, 张军华, 王延光, 等. 分方位地震勘探研究现状及进展[J]. 石油地球物理勘探, 2018, 53(3): 642-658. GONG Mingping, ZHANG Junhua, WANG Yanguang, et al. Current situation and recent progress in different azimuths seismic exploration[J]. Oil Geophysical Prospecting, 2018, 53(3): 642-658. [3] 薛姣, 顾汉明, 蔡成国. 准静态孔缝介质广义裂隙弱度研究[J]. 石油地球物理勘探, 2015, 50(6): 1146-1153. XUE Jiao, GU Hanming, CAI Chengguo. General fracture weakness for quasi-static porous fractured media[J]. Oil Geophysical Prospecting, 2015, 50(6): 1146-1153. [4] Rüger A. P wave reflection coefficients for transversely isotropic models with vertical and horizontal axis of symmetry[J]. Geophysics, 1997, 62(3): 713-722. DOI:10.1190/1.1444181 [5] Downton J E, Roure B. Interpreting azimuthal Fourier coefficients for anisotropic and fracture parameters[J]. Interpretation, 2015, 3(3): ST9-ST27. [6] Schoenberg M. Elastic wave behavior across linear slip interfaces[J]. Journal of the Acoustical Society of America, 1980, 68(5): 1516-1521. DOI:10.1121/1.385077 [7] Hudson J A. Overall properties of cracked solid[J]. Mathematical Proceedings of the Cambridge Philosophical Society, 1980, 88(2): 371-384. DOI:10.1017/S0305004100057674 [8] Thomsen L. Weak elastic anisotropy[J]. Geophysics, 1986, 51(10): 1954-1966. DOI:10.1190/1.1442051 [9] Mavko G, Mukerji T, Dvorkin J. The Rock Physics Handbook[M]. New York: Cambridge University Press, 2009. [10] Far M E, Sayers C M, Thomsen L, et al. Seismic cha-racterization of naturally fractured reservoirs using amplitude versus offset and azimuth analysis[J]. Geophysical Prospecting, 2013, 61(2): 427-447. DOI:10.1111/1365-2478.12011 [11] Far M E, Thomsen L, Sayers C M. Seismic characteri-zation of reservoirs with asymmetric fractures[J]. Geophysics, 2013, 78(2): N1-N10. [12] Grechka V, Kachanov M. Seismic characterization of multiple fracture sets: Does orthotropy suffice?[J]. Geophysics, 2006, 71(3): D93-D105. DOI:10.1190/1.2196872 [13] Grechka V, Kachanov M. Effective elasticity of rocks with closely spaced and intersecting cracks[J]. Geophysics, 2006, 71(3): D85-D91. DOI:10.1190/1.2197489 [14] Chapman M. Frequency dependent anisotropy due to mesoscale fractures in the presence of equant porosity[J]. Geophysical Prospecting, 2010, 51(5): 369-379. [15] Chapman M, Liu E, Li X. The influence of fluid sensitive dispersion and attenuation on AVO analysis[J]. Geophysical Journal of the Royal Astronomical Society, 2010, 167(1): 89-105. [16] 吴建鲁, 吴国忱. 一维频率域中观尺度虚岩石物理方法[J]. 石油地球物理勘探, 2018, 53(1): 105-112. WU Jianlu, WU Guochen. A method of 1D mesoscopic digital rock physics in the frequency domain[J]. Oil Geophysical Prospecting, 2018, 53(1): 105-112. [17] Brajanovski M, Gurevich B, Schoenberg M. A model for P wave attenuation and dispersion in a porous medium permeated by aligned fractures[J]. Geophysical Journal of the Royal Astronomical Society, 2005, 163(1): 372-384. DOI:10.1111/j.1365-246X.2005.02722.x [18] Brajanovski M, Müller T M, Gurevich B. Characteristic frequencies of seismic attenuation due to wave-induced fluid flow in fractured porous media[J]. Geophysical Journal of the Royal Astronomical Society, 2006, 166(2): 574-578. DOI:10.1111/j.1365-246X.2006.03068.x [19] Galvin R J, Gurevich B. Interaction of an elastic wave with a circular crack in a fluid-saturated porous medium[J]. Applied Physics Letters, 2006, 88(6). DOI:10.1063/1.2165178 [20] Galvin R J, Gurevich B. Effective properties of a poroelastic medium containing a distribution of aligned cracks[J]. Journal of Geophysical Research Solid Earth, 2009, 114(7). DOI:10.1029/2008JB006032 [21] Galvin R J, Gurevich B. Scattering of a longitudinal wave by a circular crack in a fluid-saturated porous medium[J]. International Journal of Solids & Structures, 2007, 44(22-23): 7389-7398. [22] Johnson D L. Theory of frequency dependent acoustics in patchy-saturated porous media[J]. Journal of the Acoustical Society of America, 2001, 110(2): 682-694. DOI:10.1121/1.1381021 [23] Gurevich B, Brajanovski M, Galvin R J, et al. P-wave dispersion and attenuation in fractured and porous reservoirs poroelasticity approach[J]. Geophysical Prospecting, 2010, 57(2): 225-237. [24] Guo J, Germán R J, Barbosa N D, et al. Seismic dispersion and attenuation in saturated porous rocks with aligned fractures of finite thickness: Theory and numerical simulations-Part 1: P-wave perpendicular to the fracture plane[J]. Geophysics, 2018, 83(1): 1-44. [25] Krzikalla F, Müller T M. Anisotropic P-SV-wave dispersion and attenuation due to inter-layer flow in thinly layered porous rocks[J]. Geophysics, 2011, 76(3): WA135-WA145. DOI:10.1190/1.3555077 [26] Galvin R J, Gurevich B. Frequency-dependent aniso-tropy of porous rocks with aligned fractures[J]. Geophysical Prospecting, 2015, 63(1): 141-150. DOI:10.1111/1365-2478.12177 [27] Yang X, Cao S, Guo Q, et al. Frequency-dependent amplitude versus offset variations in porous rocks with aligned fractures[J]. Pure & Applied Geophy-sics, 2017, 174(3): 1-17. [28] Kong L, Gurevich B, Zhang Y, et al. Effect of fracture fill on frequency-dependent anisotropy of fractured porous rocks[J]. Geophysical Prospecting, 2017, 65(6): 1649-1661. DOI:10.1111/1365-2478.12505 [29] Gassman F. Abiat the elasticity of porous media[J]. Zurich Quarterly Natare Research, 1951, 96(1): 1-23. [30] Bakulin A, Tsvankin I, Grechka V. Estimation of fracture parameters from reflection seismic data-Part Ⅰ: HTI model due to a single fracture set[J]. Geophy-sics, 2000, 65(6): 1788-1802. [31] Rüger A.Reflection Coefficients and Azimuthal AVO Analysis in Anisotropic Media[D]. Colorado School of Mines, 1996. [32] Carcione J M. Reflection and transmission of qP-qS plane waves at a plane boundary between viscoelastic transversely isotropic media[J]. Geophysical Journal International, 1997, 129(3): 669-680. DOI:10.1111/j.1365-246X.1997.tb04502.x