2. 西安大学地质工程与测绘学院, 西安 710054;
3. 中国科学院地质与地球物理研究所, 北京 100029
2. School of Geological Engineering and Geomatics, Chang'an University, Xi'an 710064, China;
3. Institute of Geology And Geophysics, Chinese Adademy of Sciences, Beijing 100029, China
岩石是矿物颗粒的集合体,矿物晶体颗粒的异质性与排列方式导致岩石在不同方向上具有不同模量,而非理想弹性体.这种矿物颗粒的复杂性常造成理论模型难以应用到实际中.此外,岩石弹性各向异性对岩石压裂拉张破坏过程(Lan et al., 2010, 2019)有着非常重要的影响(胡起等, 2014;刘财等, 2018).
目前,研究已经普遍认识到矿物和裂缝的含量、分布与弹性各向异性之间有着密切联系.在岩芯尺度上,一些具有层理分布的岩石(如:板岩、页岩)表现出明显层状各向异性的特点.在矿物晶体尺度上,某些矿物(云母、石英)(Shea and Kronenberg, 1993; Voltolini et al., 2008)、裂缝(Kwon et al., 2004)、甚至是黏土(Hornby et al., 1994; Sayers, 1995, 2013;Sayers and den Boer, 2016, 2018)都会导致岩石弹性的各向异性.各向异性影响着岩石的应力应变特征(Liu et al., 2019; Liu and Yang, 2018).但是目前,从微观尺度研究宏观各向异性的研究较少.
随着数字岩心技术的出现,岩石的矿物分布和内部结构得到了充分细致的描述(Madonna et al., 2012; Zhu and Shan, 2016).CT技术能够无损探测孔裂缝特征(Spane et al., 1994;Coker et al., 1996),电子背散射(EBSD)能够精确识别矿物晶体的取向(Valcke et al., 2006; Bascou et al., 2001; Burlini and Kunze, 2000; Lloyd and Kendall, 2005; Pera et al., 2003),这些技术为研究微观矿物和岩石弹性参数关系提供了新的技术途径.已有研究揭示(Zhong et al., 2014)晶格取向对弹性各向异性特征影响较形态择优取向更明显,但尚没有建立能够明确刻画不同矿物分布特征对于各向异性的影响方法.
针对这一问题,本研究以龙马溪组页岩为例,利用能谱技术识别矿物分布特征和室内声波技术分析岩石弹性各向异性,并采用空间分析方法建立两者之间的联系,试图从实验室尺度研究岩石矿物空间分布和弹性各向异性之间的联系.
1 方法 1.1 宏观弹性各向异性表述方法根据广义胡克定律(Mah and Schmitt, 2003; Markham, 1957),弹性材料的刚性矩阵可以表示为
(1) |
其中,Ci, j是弹性张量,σx,σy,σz为三个方向的正应力,τxy,τyz,τxz为剪应力.εx,εy,εz为正应变,γxy,γyz,γxz为三个方向的剪应变.根据弹性波理论,一维波动方程的通解公式为
(2) |
其中,Ui为质点的位移量,Ai为振动强度,λ是波长,nr是方向矢量,xr是位置矢量,v是波速,t是时间.联立方程(1)和(2),可得到克里斯托费尔方程,即:
(3) |
其中,
(4) |
联立方程(3)和(4),可以利用最小二乘法求得各个方向的弹性张量(Cij)(Markham, 1957).一般地,从不同角度观察,可得到不同的弹性张量矩阵(Cij),也可以通过旋转变换,获得最简矩阵.Thomsen曾提出了3个指数(Thomsen, 1986)来表征横向各向同性模型(VTI)的各向异性特征,分别为
(5) |
(6) |
(7) |
其中,ε和γ分别为纵波和横波的各向异性指数,σ是慢度的各向异性指数.这三个量是目前普遍采用描述宏观弹性各向异性特征的主要参量.但是,当岩石模型不满足横向各向同性模型(VTI)时,该三个指数可能不能够准确反映岩石的各项异性特征.因此,研究参照这三个指数建立新的波速差异性系数An:
(8) |
其中,ρ为密度,v是任意一方向的波速和v′为与v方向垂直的波速.
1.2 微观矿物异质性表述方法在微观方面,目前尚未有一个相应的指数来表征矿物分布的特征.为此,本文提出一个矿物微观异质性参数的概念.首先定义某种矿物在观测面上某个排列方向上的变异系数:
(9) |
其中,σi为某种矿物在i方向线含量的标准差,μi为某种矿物在i方向线含量的均值.在任一测线s上,某种矿物的线含量C(s)为
(10) |
其中,f(x)为该矿物在x位置的含量,L为测线总长度.其均值和标准差分别为:
(11) |
(12) |
公式中n为该矿物在s射线方向上的射线数量.参照宏观各向异性定义,利用任意的两个正交方向的变异系数来建立微观异质性描述参量,即,矿物微观矿物异质性系数定义为:
(13) |
其中,di为某个排列方向上的变异系数,d′i是排列方向的垂直变异系数.下面四种特殊情况可说明该矿物微观异质性系数的特点.
图 1显示了4种不同的微观分布特例,图 1a和图 1d是没有异质性的,而图 1b和图 1c是存在明显异质性的.计算结果得知,图 1a和图 1d在各个方向上的矿物A和矿物B的异质性指数均为0.而图 1B和图 1c在不同方向上的异质性系数是不同的,图 1b矿物A和矿物B在竖直方向和水平方向的差异最大,最大异质性系数为0.5,而在图 1c上,矿物A和矿物B在45°上的差异为最大的,最大异质性系数也为0.5.由此可见,微观异质性系数可以较好描述岩石矿物的微观异质性特征.
岩样采集地点是重庆涪陵地区(图 2),采样区域在山体内部开挖几十米深处.从现场测量看,露头区域的岩面非常平整,层理性强,节理发育,浓黑色,并有黄铁矿层,岩石层面产状:330°∠35°,样品的密度为2.6 g·cm-3.采集的样品中根据地质图标记出相关的主应力方向,其中x轴为最大主应力轴(σ1方向),y轴为次主应力轴(σ2方向),z轴为最小主应力轴(σ3方向),其中X面为垂直σ1方向的面,Y面为垂直σ2方向的面,Z面为垂直σ3方向的面.在同一岩块上,在X面、Y面和Z面表面分别挖取一块样品,制成3块3 mm×3 mm×1 mm的薄块,其中最大表面分别对应于X面、Y面和Z面,并进行亚离子剖光,用于观测样品在微观下的特征;同时在相同岩块中以主应力方向作为边线取出一个30 mm×30 mm×30 mm的立方体,并把8个角切去部分做成声波测试样,用于声学各向异性测量.
岩石微观测量是在中国科学院地质与地球物理研究所的矿物分析电镜上(AmicSCAN)完成.该设备根据x光能谱判断矿物晶体的分布.在该设备上,对X面,Y面和Z面分别获取样品进行了矿物晶体分布的测量.统计结果显示:岩样主要矿物为石英和伊利石,其他矿物含量较低(表 1).从不同观察面看,石英、伊利石分布不均匀,并且石英多被伊利石包裹;在X面和Y面上,石英和伊利石呈现出隐隐的条带状分布,层状堆积面应该是平行于Z面(图 3和图 4);而在Z面上,矿物的定向分布不明显(图 5).这一点符合现场的层理特征(兰恒星等, 2017, 2018).
此外,孔裂隙分布与石英、伊利石有着密切联系(图 6).石英附近的孔裂隙主要分布在晶体间,穿晶裂缝较少.伊利石附近孔裂隙既存在晶体内部也存在晶体间,且微裂缝常延伸到石英矿物之前.不同观察面的孔裂隙是不同的.在X面上,孔裂隙含量较多,指向y轴方向为主,孔裂隙的面含量最大(6.04%);Y面上次之,指向x轴方向为主,孔裂隙的面含量为中间值(5.63%);Z面上最少,孔裂隙没有明显方向性,裂隙的面含量(3.00%)为最小值.
岩石声学测试是在中国科学院地质与地球物理研究所岩石物性实验室岩石弹性参数测试仪上完成(图 7).脉冲发生接收器为泛美公司的5077PR产品.示波器为Tektronix的TDS210,数据采集速率为50 MHz,时间分辨率为0.02 μs,声波换能器的频率为1 MHz,总体测量精度优于0.5%.实验中,超声换能器分别置于样品的两个平行端面,纵波测量使用一对纵波换能器,横波测量使用一对剪切横波换能器.实验开始前,先用游标卡尺测量2个测量面之间的距离,再将换能器直接对接,确定系统基时,然后再对样品进行声学测量.计算波速公式为
(14) |
其中, l为岩样测量面之间的距离,t为声波走时,t0为对接探头时确定的系统基时.测量时,样品沿最大主应力方向标记为x方向,次主应力的方向为y方向,最小主应力方向为z方向,三个主应力方向相互垂直.基于该定义,不同测点上声波传播方向、振动方向和获得的纵横波速度如表 2所示.
根据表 2测量的结果,研究假定样品满足横向各向同性模型(VTI)(Zhubayev et al., 2016; 汪虎等, 2017; 张冰等, 2018).该样品的刚性矩阵可以表示为
(15) |
利用式(1)—(4)联解,可得该样品的刚性矩阵为
(16) |
其中C矩阵的单位为10 GPa.从获得的刚性矩阵看出,z轴方向的弹性参数(C33)为x轴(C11)弹性参数和y轴弹性参数(C22)的0.75倍,说明岩石存在明显层理特征.
3 讨论 3.1 宏观波速与微观变异系数之间的关系为了研究样品声波和矿物、孔裂缝之间的关系,分别在不同观察面(X,Y和Z)上,按照不同的夹角统计了主要矿物(石英和伊利石)和孔裂隙的变异系数(图 8).研究结果显示主要矿物变异系数随角度的变化在不同面上呈现出明显的异质性.在X面上(图 8a和8d),石英、伊利石和孔裂隙的变异系数随着角度的增大而降低,特别是孔裂隙的降幅最大,在45°附近急降到最小值附近,随后随角度增大维持在低值.在Y面上(图 8b和8e),同样观察到矿物的变异系数随着角度明显变化.石英变异系数在角度小于45°时呈现高值、而在大于45°后为相对低值;伊利石的变异系数先是随角度增大而逐步降低,在45°附近达到低值,随后随角度增大保持低值.孔裂隙变异系数随角度呈V型变化,在45°附近达到最低值.在Z面上(图 8c和8f),石英、伊利石和孔裂隙变异系数随角度没有明显变化.相应地,当垂直层理时(X面和Y面的角度为90°时),纵波波速和横波波速最低,而当平行层理时(X面和Y面的角度为0°时),纵波波速和横波波速最高.Z面为平行层理面,所观察到纵波波速和横波波速基本不变.波速变化特征与先前的研究结果一致(伍宇明等, 2016).
岩石纵波速度与矿物和裂缝之间的变异系数的相关性统计列于表 3.岩石纵波速度与石英变异系数的相关系数在X、Y和Z面上都达到了0.7以上.波速与裂隙、伊利石变异系数的相关系数在X和Z面也很高,但在Y面相对低.孔裂隙和波速之间相关系数在X面上高达0.94,而在Y面上为0.33,与X面相关系数相差接近3倍,这说明孔裂隙对于波速的影响较为复杂.上述数据说明,宏观波速与微观变异系数是有明显相关性的,主要矿物的分布特征是影响波速随方向变化的主要因素.
岩石波速的各向异性可以使用波速差异性系数(公式8)来表征,而岩石的微观分布方向性可以用微观异质性系数来刻画.图 9显示波速差异性系数与微观异质性系数的关系.在X面和Y面上,石英、伊利石和孔裂隙的微观异质性系数呈现明显随角度增大而下降的趋势,孔裂隙的微观异质性系数甚至在45°附近出现断崖式下落.而在Z面上,矿物和孔裂隙的异质性系数变化不明显.这说明在X面和Y面上,矿物和孔裂隙具有明显的定向排列特征,而在Z面上属于均质分布.而纵波速度和横波速度的差异性系数在X和Y面都表现出随角度增大而下降的行为,显示出矿物和孔裂隙异质性系数存在明显的正相关特征.在Z面上纵波速度和横波速度的差异性系数基本不变,与Z面上矿物和孔裂隙的异质性系数的表现相一致.这说明,纵波速度和横波速度的差异性系数是依赖岩石的微观分布的.这和宏观上波速的方向性变化是依赖于层理的认知是一致的.纵横波波速一般在平行层理方向上最大,垂直层理方向上最小(伍宇明等,2016),实验室的三轴测量结果也指出了这种行为特征(Chen et al., 2018).这些结果表明宏观各向异性与微观矿物分布具有高度相关性,说明微观矿物分布特征应该是导致弹性各向异性的主要原因.因此,各向异性特征可能是微观矿物分布的一种宏观表现.对于这种异质性和各向异性的关系,传统岩石物理模型(Reuss,1929;Voigt, 1928)不能够充分表达这一特征,因此,需要发展更加完善的微观岩石物理模型来刻画这种关系.
研究建立了微观异质性指标,并和岩石的宏观弹性各向异性特征建立了联系,得出以下三点结论:
(1) 微观背散射图像显示龙马溪组页岩主要由石英和伊利石构成,且这两种矿物在X和Y面上都有很好的定向排列特征.孔裂隙也表现出定向排列性质.相应地,主要矿物和孔裂隙变异系数在X面和Y面上表现出随角度增大而降低的特征,特别是孔裂隙的降幅最大.而在Z面,主要矿物和孔裂隙变异系数变化不明显,与Z面上矿物和孔裂隙无明显的方向性一致.
(2) 在X和Y面上,波速随角度增加而减小,与主要矿物和孔裂隙变异系数变化趋势相同.而在Z面上,波速变化不大,这与主要矿物和孔裂隙变异系数变化不明显的特征一致,这说明,宏观波速与微观变异系数是有明显相关性的.
(3) 波速的差异性系数与主要矿物和孔裂隙的异质性系数在X和Y面都表现出随角度增大而下降的行为,显示出它们之间存在正相关特征.波速的差异性系数在Z面上基本不变,与Z面上矿物和孔裂隙的异质性系数的表现相一致.这说明,宏观波速差异性系数是依赖微观异质性系数,这和宏观上波速的方向性变化依赖于层理的认知是一致的.宏观各向异性性质与矿物的分布异质性具有高度相似性,说明微观矿物和孔裂隙分布特征应该是导致弹性各向异性的主要原因.因此,各向异性特征是微观矿物分布的一种宏观表现.
Bascou J, Barruol G, Vauchez A, et al. 2001. EBSD-measured lattice-preferred orientations and seismic properties of eclogites. Tectonophysics, 342(1-2): 61-80. DOI:10.1016/S0040-1951(01)00156-1 |
Burlini L, Kunze K. 2000. Fabric and seismic properties of Carrara marble mylonite. Physics and Chemistry of the Earth, Part A:Solid Earth and Geodesy, 25(2): 133-139. DOI:10.1016/S1464-1895(00)00022-3 |
Chen J H, Lan H X, Macciotta R, et al. 2018. Anisotropy rather than transverse isotropy in Longmaxi shale and the potential role of tectonic stress. Engineering Geology, 247: 38-47. DOI:10.1016/j.enggeo.2018.10.018 |
Coker D A, Torquato S, Dunsmuir J H. 1996. Morphology and physical properties of Fontainebleau sandstone via a tomographic analysis. Journal of Geophysical Research:Solid Earth, 101(B8): 17497-17506. |
Hornby B E, Schwartz L M, Hudson J A. 1994. Anisotropic effective-medium modeling of the elastic properties of shales. Geophysics, 59(10): 1570-1583. DOI:10.1190/1.1443546 |
Hu Q, Chen X H, Li J Y. 2014. Shear wave velocity prediction for shale gas reservoirs based on anisotropic rock physics model. Geophysical Prospecting for Petroleum (in Chinese), 53(3): 254-261. |
Kwon O, Kronenberg A K, Gangi A F, et al. 2004. Permeability of illite-bearing shale:1. Anisotropy and effects of clay content and loading. Journal of Geophysical Research:Solid Earth, 109(B10): B10205. DOI:10.1029/2004JB003052 |
Lan H X, Martin C D, Hu B. 2010. Effect of heterogeneity of brittle rock on micromechanical extensile behavior during compression loading. Journal of Geophysical Research:Solid Earth, 115(B1): B01202. DOI:10.1029/2009JB006496 |
Lan H X, Wu Y M, Li Q W, et al. 2017. Reconstruction and fractal analysis for three dimensional fracture network of Longmaxi shale. Journal of Engineering Geology (in Chinese), 25(6): 1557-1565. |
Lan H X, Chen J H, Wu Y M. 2018. Spatial characterization of micro-and nanoscale micro-cracks in gas shale before and after triaxial Compression test. Journal of Engineering Geology (in Chinese), 26(1): 24-35. |
Lan H X, Chen J H, Macciotta R. 2019. Universal confined tensile strength of intact rock. Scientific Reports, 9(1): 6170. DOI:10.1038/s41598-019-42698-6 |
Li Q H, Chen M, Jin Y, et al. 2012. Rock mechanical properties and brittleness evaluation of shale gas reservoir. Petroleum Drilling Techniques (in Chinese), 40(4): 17-22. |
Liu C, Fu W, Guo Z Q, et al. 2018. Rock physics inversion for anisotropic shale reservoirs based on Bayesian scheme. Chinese J. Geophys. (in Chinese), 61(6): 2589-2600. DOI:10.6038/cjg2018L0176 |
Liu X, Yang J. 2018. Shear wave velocity in sand:effect of grain shape. Géoechnique, 68(8): 1-22. |
Liu X, Zhang N, Lan H X. 2019. Effects of sand and water contents on the small-strain shear modulus of loess. Engineering Geology, 260: 105202. DOI:10.1016/j.enggeo.2019.105202 |
Lloyd G E, Kendall J M. 2005. Petrofabric-derived seismic properties of a mylonitic quartz simple shear zone:Implications for seismic reflection profiling. Geological Society, London, Special Publications, 240(1): 75-94. DOI:10.1144/GSL.SP.2005.240.01.07 |
Madonna C, Almqvist B S G, Saenger E H. 2012. Digital rock physics:numerical prediction of pressure-dependent ultrasonic velocities using micro-CT imaging. Geophysical Journal International, 189(3): 1475-1482. DOI:10.1111/j.1365-246X.2012.05437.x |
Mah M, Schmitt D R. 2003. Determination of the complete elastic stiffnesses from ultrasonic phase velocity measurements. Journal of Geophysical Research:Solid Earth, 108(B1): ECV 6-1-ECV 6-11. DOI:10.1029/2001JB001586 |
Markham M F. 1957. Measurement of elastic constants by the ultrasonic pulse method. British Journal of Applied Physics, 8(S6): S56. DOI:10.1088/0508-3443/8/S6/312 |
Pera E, Mainprice D, Burlini L. 2003. Anisotropic seismic properties of the upper mantle beneath the Torre Alfina area (Northern Apennines, Central Italy). Tectonophysics, 370(1-4): 11-30. DOI:10.1016/S0040-1951(03)00175-6 |
Reuss A. 1929. Berechnung der fließgrenze von mischkristallen auf grund der plastizitätsbedingung für einkristalle. ZAMM-Journal of Applied Mathematics and Mechanics/Zeitschrift für Angewandte Mathematik und Mechanik, 9(1): 49-58. |
Sayers C M. 1995. Anisotropic velocity analysis. Geophysical Prospecting, 43(4): 541-568. DOI:10.1111/j.1365-2478.1995.tb00266.x |
Sayers C M. 2013. The effect of kerogen on the elastic anisotropy of organic-rich shales. Geophysics, 78(2): D65-D74. |
Sayers C M, den Boer L D. 2016. The elastic anisotropy of clay minerals. Geophysics, 81(5): C193-C203. DOI:10.1190/geo2016-0005.1 |
Sayers C M, den Boer L D. 2018. The elastic properties of clay in shales. Journal of Geophysical Research:Solid Earth, 123(7): 5965-5974. DOI:10.1029/2018JB015600 |
Shea W T Jr, Kronenberg A K. 1993. Strength and anisotropy of foliated rocks with varied mica contents. Journal of Structural Geology, 15(9-10): 1097-1121. DOI:10.1016/0191-8141(93)90158-7 |
Thomsen L. 1986. Weak elastic anisotropy. Geophysics, 51(10): 1954-1966. DOI:10.1190/1.1442051 |
Valcke S L A, Casey M, Lloyd G E, et al. 2006. Lattice preferred orientation and seismic anisotropy in sedimentary rocks. Geophysical Journal International, 166(2): 652-666. DOI:10.1111/j.1365-246X.2006.02987.x |
Voigt W. 1928. Lehrbuch der Kristallphysik (BG Teubner Leipzig und Berlin) 980 S. Reproduced 1966 Spring Fachmedien Wiesbaden GmbH.
|
Wang H, Guo Y T, Wang L, et al. 2017. An experimental study on mechanical anisotropy of shale reservoirs at different depths. Rock and Soil Mechanics (in Chinese), 38(9): 2496-2506. |
Wu Y M, Lan H X, Gao X, et al. 2016. The relationship between the volume density of cracks and acoustic properties of the shale core samples from Fulin. Chinese J. Geophys. (in Chinese), 59(10): 3891-3990. DOI:10.6038/cjg20161032 |
Zhang B, Guo Z Q, Xu C, et al. 2018. Fracture properties and anisotropic parameters inversion of shales based on rock physics model. Journal of Jilin University (Earth Science Edition) (in Chinese), 48(4): 1244-1252. |
Zhong X, Frehner M, Kunze K, et al. 2014. A novel EBSD-based finite-element wave propagation model for investigating seismic anisotropy:Application to Finero Peridotite, Ivrea-Verbano Zone, Northern Italy. Geophysical Research Letters, 41(20): 7105-7114. DOI:10.1002/2014GL060490 |
Zhu W, Shan R. 2016. Digital core based transmitted ultrasonic wave simulation and velocity accuracy analysis. Applied Geophysics, 13(2): 375-381. |
Zhubayev A, Houben M E, Smeulders D M J, et al. 2016. Ultrasonic velocity and attenuation anisotropy of shales, Whitby, United Kingdom. Geophysics, 81(1): D45-D56. |
胡起, 陈小宏, 李景叶. 2014. 基于各向异性岩石物理模型的页岩气储层横波速度预测. 石油物探, 53(3): 254-261. DOI:10.3969/j.issn.1000-1441.2014.03.002 |
兰恒星, 伍宇明, 李全文, 等. 2017. 龙马溪组页岩三维缝网重构及分形分析. 工程地质学报, 25(6): 1557-1565. |
兰恒星, 陈俊辉, 伍宇明. 2018. 三轴压缩试验前后含气页岩微纳尺度裂隙空间分布特征研究. 工程地质学报, 26(1): 24-35. |
李庆辉, 陈勉, 金衍, 等. 2012. 页岩气储层岩石力学特性及脆性评价. 石油钻探技术, 40(4): 17-22. DOI:10.3969/j.issn.1001-0890.2012.04.004 |
刘财, 符伟, 郭智奇, 等. 2018. 基于贝叶斯框架的各向异性页岩储层岩石物理反演技术. 地球物理学报, 61(6): 2589-2600. DOI:10.6038/cjg2018L0176 |
汪虎, 郭印同, 王磊, 等. 2017. 不同深度页岩储层力学各向异性的试验研究. 岩土力学, 38(9): 2496-2506. |
伍宇明, 兰恒星, 高星, 等. 2016. 涪陵地区井下页岩岩芯裂缝体密度与声学性质关系实验研究. 地球物理学报, 59(10): 3891-3900. DOI:10.6038/cjg20161032 |
张冰, 郭智奇, 徐聪, 等. 2018. 基于岩石物理模型的页岩储层裂缝属性及各向异性参数反演. 吉林大学学报(地球科学版), 48(4): 1244-1252. |