2020, Vol. 63
2. 中国石油大学(华东)深层油气重点实验室, 青岛 266580;
3. 南京工业大学交通运输工程学院, 南京 210009
2. Key Laboratory of Deep Oil and Gas, China University of Petroleum(East China), Qingdao 266580, China;
3. School of Transportation Engineering, Nanjing Tech University, Nanjing 210009, China
龙马溪组富有机质页岩主要由矿物骨架、孔隙结构及有机质(TOC)构成，共同决定了其弹性特征，尤其受TOC含量、分布及有机孔的影响极大.研究TOC相关特征对弹性模量的影响对页岩具有重要指导作用.与常规的等效介质解析方法及有限元数值方法不同，本文采用更加精细的晶格弹簧-随机孔隙耦合模型(LSM-RVM)，研究TOC含量及成熟度对页岩弹性模量的影响.
富有机质页岩经历漫长的地质活动进程，岩石内部组分复杂(Zhao et al., 2016).特别是页岩中的有机质干酪根，其成熟过程涉及大量物理及化学变化，对页岩力学性质产生重要影响.页岩实验研究(Vernik and Nur, 1992；Vernik and Liu, 1997)建立了TOC含量及其成熟度与弹性波传播特性的量化关系，分析了矿物组分、裂缝及节理排布对Bakken页岩弹性各向异性特征的影响.Sondergeld等(2000)、Sondergeld和Rai(2011)基于声学实验分析了TOC含量对Kimmeridge页岩弹性参数的影响，结果表明随着TOC含量的增多，页岩样本呈现出更显著的各向异性特征.邓继新等(2004)基于超声实验研究了干燥和饱油条件下，不同围压及孔隙流体对泥岩和页岩各向异性的影响.一般而言，实验研究是最直观有效的方法，然而，由于研究对象的错综复杂以及对实验条件要求苛刻，因而也具有一定的局限性.
理论分析方法广泛用于富有机质页岩的等效弹性建模.Hornby等(1994)综合利用自洽(SCA)模型与微分等效介质(DEM)模型，计算了Kimmeridge页岩的等效弹性参数.Vernik和Landis(1996)改进了Backus平均方法(Backus, 1962)，进而分析页岩中干酪根网络对弹性各向异性特征影响.Carcione等(2011)利用Backus平均方法计算了页岩中黏土颗粒定向排列对其地震弹性性质及各向异性的影响.Guo等(2013)综合利用各向同性SCA模型和Backus平均方法分析了黏土和干酪根颗粒的分布方向对Barnett页岩脆性指标的影响.Zhao等(2016)基于DEM模型分析了页岩弹性参数随不同成熟度TOC含量的变化趋势.
作为衔接理论分析与实验测量之间的有效方法，数值模拟在科学研究领域发挥重要作用.页岩气储层微观尺度建模需要考虑矿物组分、岩石颗粒、孔隙/裂缝网、TOC分布及有机孔、岩石层理等，裂缝网得益于成像技术的进步，采用扫描电镜(SEM)和高精度X射线CT成像技术，可成像低至微、纳米尺度的细观结构，数字重构包括孔隙流及晶粒形态等岩石中的各种不均匀体.基于高精度数字岩心数据，采用有限差分(FDM)或有限单元(FEM)数值方法可以研究岩石的弹性性质和声响应特征.Arns等(2002)采用FEM数值模拟Fontainebleau砂岩的弹性性质，并与实验数据进行了比较.数值模型引入更多的岩石微结构细节有助于理解富含有机质页岩的岩石物理性质(Dvorkin et al., 2003；Saenger et al., 2004)、电学性质(Liu et al., 2009)、孔隙流特征(Kewen and Ning, 2008)、储层特性(Ceron et al., 2013)等.Zhang等(2016)基于高精度数字岩心数据，分析了龙马溪组页岩孔隙度和干酪根含量对弹性性质的影响.Liu和Fu(2020)基于LSM模型，详细分析了微观结构属性对岩石各向异性特征的作用程度.
高精度的数值建模方法有助于捕捉页岩中细微结构(例如岩石颗粒形状、TOC分布形态、微裂缝和层理等)对样本弹性性质的影响.传统数值方法基于连续介质理论，忽略了页岩内部矿物组分离散分布特点和形态特征.页岩本身结构具有离散性，不同矿物晶粒结构之间的局部应力也呈现出不连续性，需要基于离散介质理论的数值建模方法.基于上述考虑，本文提出了一种离散数值建模方法——晶格弹簧-随机孔隙耦合模型(LSM-RVM).首先，通过高精度CT成像技术，得到含多种矿物组分和精细微观结构的数字岩心.引入Lin等(2011)提出的随机孔隙模型(RVM)，调整相关系数，在数字岩心岩石骨架中生成与实验数据统计特点符合的TOC分布，据此建立相关LSM模型.利用该LSM-RVM耦合模型进行数值模拟，最后讨论TOC含量及成熟度对页岩弹性参数的影响.
1 基本方法 1.1 LSM模型介绍基于原子晶格结构发展起来的晶格弹簧模型(Lattice Spring Model，LSM)将研究对象离散为分散的质量点，质量点之间采用弹簧相连(Ostoja-Starzewski, 2002).晶格模型的起源可以追溯到1941年(Hrennikoff, 1941)，最早的晶格模型只使用拉压弹簧，节点间通过中心力作用，曾广泛用于含缺陷介质的弹性和断裂行为研究(Ray and Chakrabarti, 1985；Sahimi and Goddard, 1986；De Arcangelis et al., 1989；Liu J X et al., 2007；Liu N et al., 2020).
在本文中，晶格弹簧模型两质量点间通过弹性刚度分别为Kn和Ks的拉压弹簧和剪切弹簧发生作用，受力模型如图 1所示.其中两节点间的拉压力Fn，可以写为如下形式：
|
图 1 线弹性弹簧模型(包括正向力和切向力) Fig. 1 Linear elastic spring model including normal and shear forces |
|
(1) |
这里，
|
(2) |
|
(3) |
其中，leq、l、Δl、x1、x2和l12分别是两节点间弹簧平衡距离，两节点距离，法向位移增量，两节点的位置坐标以及两节点间弹簧的法向向量.对于切向力Fs，有
|
(4) |
|
(5) |
这里ls是切向方向向量；Δs是切向位移增量.
在离散晶格模型中，通过如图 2所示的Verlet跳蛙积分算法来计算晶格节点的运动，并由此计算晶格节点的速度，进而更新晶格节点位置坐标等.从图中左侧的三元组
|
图 2 Verlet跳蛙积分算法 Fig. 2 Verlet leapfrog integration |
|
(6) |
|
(7) |
|
(8) |
|
(9) |
其中，mi为晶格节点i的质量，Fi(t)和ai(t)是对应节点在t时刻的合力和加速度，Δt是时间步长，
为了保证上述(7)和(8)式显式积分算法的稳定性，Δt有一个上限值，需满足如下的关系式(Šmilauer et al., 2010)：
|
(10) |
|
(11) |
其中ωmax是系统最大特征频率值.这里，如图 3所示，每个节点受到来自多个连接的作用，节点间作用既包含拉压弹簧，也包含切向弹簧，因而将节点位移刚度矩阵设为Kij，且令切向弹簧刚度Ks=ξKn，ξ为切向弹簧与拉压弹簧的刚度比值.则作用在该节点上的位移刚度矩阵Kij可以表示为(Chareyre and Villard, 2005)
|
图 3 晶格弹簧模型中节点受力情况 Fig. 3 Nodal forces in LSM |
|
(12) |
在笛卡尔坐标系下，i和j∈{x, y, z}.在数值模拟中，(11)式和(12)式共同决定Δtcr的大小.
为了满足准静态条件，需要引入非黏性阻尼，通过数值阻尼压制节点速度改变引起的振荡作用力(Šmilauer et al., 2010)，单位时间步长Δt内的衰减比率可表示为
|
(13) |
其中，ζ为数值阻尼.因此，整个晶格弹簧模型的计算执行流程可总结为图 4，图中的停止条件可设定为给定时间步或给定拉伸(压缩)应变等.
|
图 4 LSM模型数值模拟流程图 Fig. 4 Flowchart of LSM |
随机介质理论方法广泛用于模拟岩石中矿物颗粒的随机分布(Korn, 1993；Ikelle et al., 1993)，例如基于高斯型和指数型自相关函数的随机介质模拟.Lin等(2011)提出了一种基于指数型自相关函数的二维随机孔模型(Random Void Model，RVM)，用于描述介质中的随机孔分布，其表达式如下：
|
(14) |
其中a和b分别是自相关长度；当a≠b时，(14)式为椭圆型自相关函数.图 5给出了生成随机孔模型(RVM)的具体实现流程.
|
图 5 生成二维随机孔模型流程图：FFT代表快速傅里叶变换；IFFT代表快速傅里叶逆变换(Lin et al., 2011)；图中ε为随机扰动参数 Fig. 5 Flowchart of generating the 2D random void model: FFT and IFFT represent fast-Fourier transform and inverse fast-Fourier transform, respectively (Lin et al., 2011); ε is the random variance |
一般而言，孔隙岩石弹性建模需要考虑各种矿物组分的弹性模量、矿物比例、矿物颗粒几何形态等.忽略岩石中微观孔/裂隙结构的影响，采用最小势能和最小余能原理可分别预测岩石等效弹性参数的上限和下限.Voigt (1889)提出的体积弹性模量KVoigt和剪切弹性模量μVoigt，表示如下：
|
(15) |
|
(16) |
式中，N代表该材料具有N种组分，第i种成分的体积模量和剪切模量分别为Ki和μi，所占比例用Vi表示.该方法实质为刚度均匀化等效.Reuss(1929)给出的材料弹性参数估计下边界表达式如下：
|
(17) |
|
(18) |
这种方法实质为柔度均匀化等效.
Hashin和Shtrikman(1962)基于变分原理给出了范围更小的上下界限理论.随后，一些学者(如Hill, 1965；Budiansky, 1965；Mori and Tanaka, 1973)基于Eshelby(1957)提出的含孔介质等效弹性原理又进一步完善了复杂介质等效弹性参数的估计模型.其中，Voigt-Reuss边界理论方法由于形式简单，被广泛地用于地球物理领域.
2 龙马溪组页岩数字岩心 2.1 龙马溪组页岩近年来，随着四川盆地及其周缘页岩气勘探开发的不断深入，盆地内下奥陶系和上志留系的五峰-龙马溪组有机质演化程度高，逐渐成为了国内勘探突破的首选区域(邓继新等, 2015).龙马溪组页岩与美国Barnett硅质页岩生物成因相类似，均是沉积于深水陆棚闭塞环境，海平面上升控制形成优质页岩的缺氧环境.二者都相对富集脆性矿物，但龙马溪组页岩脆性矿物含量变化较大，且矿物含量普遍低于Barnett页岩，黏土矿物含量相对较高，脆性矿物以陆源石英、长石为主(Guo et al., 2013).因此，龙马溪组页岩具有与Barnett页岩相似的岩相学特征，但岩石矿物成分和岩性特征差异较大.
2.2 高精度X射线CT扫描技术与数字岩心本研究采用上海光源第三代同步辐射装置对直径为800 μm、长度为3 mm的龙马溪组页岩样本进行高精度X射线CT(X-ray Micro-CT，XMCT)扫描成像，实验装置及岩心样本如图 6所示.采用有序子集最大期望值(Ordered Subset Expectation Maximization，OSEM)重建算法(Hudson and Larkin, 1994)对数字岩心数据进行处理，提高图像精度.图 7为3.5 s曝光时间条件下，通过1080个投影CT图重构的三维数字岩心三视图.为了更精确地捕获样品内部的矿物组分和微观结构，我们从中提取了一个尺寸为390 μm×390 μm×455 μm，包含600×600×700个像素体的立方块，其中数字图像分辨率为0.65 μm/像素.
|
图 6 高精度X射线CT扫描设备(左)；直径为800 μm的龙马溪岩心样本(右) Fig. 6 X-ray micro-CT scanner at Laboratory (left) and a Longmaxi shale sample with diameter of 800 μm (right) |
|
图 7 三维数字岩心的三视图 Fig. 7 Three views of the 3D digital core |
页岩微观结构复杂，包括如石英、长石、黄铁矿、碳酸钙等脆性矿物，以伊利石为主的黏土类矿物，以及不规则孔隙结构等(邓继新等, 2015).XMCT技术根据不同介质对X射线透射率的差异性来成像页岩内部的微观特征.表 1列出了页岩中常见的几种矿物吸收系数，在重构的图像上表现为像素点的不同灰度值，对于8位二进制图像，对应像素值为0~255.从表 1中可以看出，黄铁矿具有较高的吸收系数，而其他矿物对应的值较低.
|
表 1 龙马溪组页岩主要矿物在20 keV射线光源 Table 1 Absorption coefficients of X-ray with 20 keV energy for common minerals in the Longmaxi shale sample |
脆性矿物与黏土类矿物吸收系数相近，需要高精度的图像分辨处理技术.本文基于大津算法进行多级图像分割，将数字岩心按灰度值划分为三组：包括不规则孔隙结构(0~36)、黄铁矿(131~255)和其他矿物(37~130).图 8是以3D数字岩心的2D切片为例的矿物组分划分结果，其中黑色为背景色.脆性矿物与黏土类矿物构成岩石骨架部分，呈草莓状分布的黄铁矿由多个紧密堆积的晶粒构成.随机分布的孔隙大部分可近似为圆孔状，少部分由于受水平方向主裂缝影响而呈偏平状.基于Reuss等效弹性模量边界理论方法((17)和(18)式所示)，我们计算了不同黏土类矿物含量的龙马溪组页岩样本的弹性参数，结果如表 2所示.
|
图 8 2D数字岩心切片及矿物组分划分：骨架(37~130)；黄铁矿(131~255)；孔隙结构(0~36) Fig. 8 2D digital core slice and mineral classification: matrix (37~130); pyrites (131~255); pores (0~36) |
|
|
表 2 龙马溪组页岩样本主要组分的材料常数 Table 2 Material constants of the main constituents in the Longmaxi shale samples |
邓继新等(2015)在其研究中曾指出页岩储层孔隙类型主要包括由脆性矿物、黏土类矿物颗粒间微孔形成的原生孔隙；页岩中有机质生烃、黏土矿物脱水、伊利石和长石等不稳定矿物溶蚀作用形成的微小孔洞和微裂隙，即次生孔隙等.TOC与黄铁矿主要分布于原生粒间孔隙中，孔隙结构决定了二者的形态.图 9统计了2D数字岩心切片中的孔隙结构、黄铁矿的平均半径，以及纵横比.从该图中可以看出孔隙与黄铁矿的平均半径集中在1.3~3.2 μm之间，纵横比近似为1.0.
|
图 9 2D数字岩心切片中孔隙及黄铁矿平均半径与纵横比统计分布 Fig. 9 Statistical distributions of average radius and aspect ratio for pores and pyrites in the 2D digital core slice |
龙马溪组页岩的TOC主要分布于原生粒间孔隙中，呈随机斑块状分布，无明显的定向性.由于含量较低，一般不作为岩石骨架的组成部分来承担作用力(邓继新等, 2015).因而，本文采用RVM模型在数字岩心中随机生成干酪根，将(14)
式中自相关系数a与b的比值设为1.0，a与b的值在3.0~9.0之间变化，据此重建不同TOC含量的数字岩心样本.通过设定不同阈值系数T进行图像分割，得到与大多数龙马溪组页岩样本TOC含量相符的LSM-RVM模型，具体流程如图 10所示.
|
图 10 生成符合TOC实际分布特征的LSM-RVM流程图 Fig. 10 Flowchart for generating LSM-RVMs with the actual distribution characteristics of TOC |
图 11为根据数字岩心生成的包含不同TOC体积含量的LSM-RVM.从该图中可以看出随着自相关系数a、b的增大，随机生成斑块的尺度增大；RVM中随机扰动系数ε2(见图 5)与斑块尺寸同样呈正相关.
|
图 11 数字岩心的LSM-RVM模型 (a) a=b=3, ε2=0.1, T=0.4；(b) a=b=9, ε2=0.1, T=0.4；(c) a=b=9, ε2=0.5, T=0.4. Fig. 11 LSM-RVM models for digital cores |
假设页岩样本中的每组矿物是各向同性的，因而采用三角晶格排布的LSM模型，该模型的弹簧刚度与其晶格单胞的等效弹性参数存在如下关系：
|
(19) |
|
(20) |
如图 12所示，这里采用沿不同方向的单轴压缩试验估计LSM模型的等效弹性参数，得到的值与理论值比较，用于调整输入参数与验证该模型的有效性.表 3中给出了不同成熟度干酪根的弹性模量(Lucier et al., 2011；Yan and Han, 2013；Qin et al., 2014).
|
图 12 LSM-RVM模型的加载方式 (a)沿竖直轴压缩；(b)沿水平轴压缩. Fig. 12 Loading paths of LSM-RVM model (a) Uniaxial compression along vertical axis; (b) Uniaxial compression along horizontal axis. |
|
|
表 3 不同成熟度有机质的弹性模量(Zhao et al., 2016) Table 3 Elastic moduli of organic matter at different maturity stages (Zhao et al., 2016) |
LSM模型中，拉压弹簧和剪切弹簧刚度比值决定了所表征材料的泊松比，同时，该材料常数也是验证晶格模型的重要指标.图 13显示的是表 3中列出的三种材料的LSM模型在不同应变速率0.01/Δts-1、0.05/Δts-1、0.10/Δts-1加载条件下得到的泊松比随着数值阻尼的变化趋势：材料1和材料3，随着数值阻尼的增大而泊松比减小；材料2的泊松比随着数值阻尼增大，呈现先上升后下降趋势.
|
图 13 三组材料在应变速率分别为0.01/Δts-1，0.05/Δts-1，0.10/Δts-1条件下的泊松比随数值阻尼变化趋势 Fig. 13 Poisson′s ratios of three groups of materials change with increasing numerical damping at the strain ratios of 0.01/Δts-1, 0.05/Δts-1, 0.10/Δts-1 |
这里采用NL2SOL算法估计整体误差，表达式如下：
|
(21) |
式中，
|
(22) |
上标j表示如表 3所示的三种材料中的第j种介质；这里，数值阻尼ζ=0, 0.2, 0.4, 0.6, 0.8共五个值，加载速率
|
|
表 4 不同加载条件下得到结果的误差值 Table 4 Result errors at different loading conditions |
图 14为岩石基体分别包含40%、60%和80%黏土类矿物数字岩心样本的杨氏模量随不同成熟度TOC含量变化图，并与Voigt-Reuss等效边界理论进行比较.从该图中可以看出，沿x轴压缩得到的杨氏模量Exx普遍高于对应的Eyy；过度成熟的岩心样本弹性参数落在Voigt-Reuss上下边界包围的区间内；随着TOC成熟度的减小，弹性参数偏离基于各向同性假设的等效边界理论，反映了各向异性程度的增大；岩石骨架中随着黏土含量增大，岩心样本的刚度呈上升趋势.
|
图 14 比较不同成熟度TOC含量对杨氏模量的影响 (a)未成熟；(b)成熟；(c)过熟. Fig. 14 Comparisons of Young′s moduli along two directions with various TOC content at different maturity stages (a) Immature; (b) Mature; (c) Over-mature. |
图 14显示了基于弹性假设，数字岩心沿不同方向加载得到的杨氏模量值存在微小差异，说明了该数字岩心表现为弱各向异性特点.为了更直观地呈现这种差异性，图 15为包含7.18%成熟TOC，岩石骨架中黏土占比60%的数字岩心分别沿x轴、y轴压缩得到的偏应力场分布.从该图中可以看出，二者有显著差异，局部最大应力出现在垂直压缩轴方向的孔隙尖端、不同矿物介质分界面.为更清楚地表征不同矿物组分(图 16a与图 16b对比)、TOC成熟度(图 16a与图 16c对比)、TOC含量(图 16a与图 16d对比)对结构局部各向异性的影响差异，图 16显示了4种数字岩心样本在外力加载条件下的局部偏应力场分布.
|
图 15 ε=0.008时，含7.18%成熟TOC，骨架黏土占比60%的数字岩心沿x、y轴加载的偏应力场 Fig. 15 The deviatoric stress fields of digital shale cores with 7.18% mature TOC and 60% clay-matrix at εxx=0.008 and εyy=0.008 |
|
图 16 εyy=0.008时，不同成熟度TOC、黏土含量页岩样本的偏应力场 Fig. 16 When εyy=0.008, the deviatoric stress fields of digital shale cores with various TOC and clay contents at different maturity stages |
本文以直径为800 μm钻心取样自龙马溪组页岩储层的岩心样本为例，说明了一种离散数值模拟方法：首先，利用高精度CT成像技术得到龙马溪组页岩的数字岩心；根据介质的光学、物理、化学特性等，采用大津算法进行多级图像划分，进而识别多种
矿物组分，并确定每个矿物组的材料参数；综合实验数据及RVM方法生成TOC分布；建立相关LSM-RVM模型.基于该模型，量化研究等效弹性参数和各向异性特性随TOC含量及成熟度的变化趋势，得到以下结论.
(1) 通过与理论值比较验证了LSM的有效性；在LSM建立模型过程中，为获取准静态加载条件，需引入数值阻尼；数值模拟结果显示，数值阻尼与加载应变速率会影响计算精度，随着应变速率增大，等效泊松比值与理论值误差增大；数值阻尼变化趋势则与目标参数值有关.
(2) 岩心中TOC成熟度、含量以及其他矿物组分共同决定弹性参数变化趋势，随着TOC成熟度增大，样本表现为较弱的各向异性；得到的变化趋势与实验结果相一致，说明LSM-RVM能够生成符合TOC及多种矿物实际分布特征的数值模型，是一种精细数值建模方法.
Arns C H, Knackstedt M A, Pinczewski W V, et al. 2002. Computation of linear elastic properties from microtomographic images:Methodology and agreement between theory and experiment. Geophysics, 67(5): 1396-1405. DOI:10.1190/1.1512785 |
Backus G E. 1962. Long-wave elastic anisotropy produced by horizontal layering. Journal of Geophysical Research, 67(11): 4427-4440. DOI:10.1029/JZ067i011p04427 |
Budiansky B. 1965. On the elastic moduli of some heterogeneous materials. Journal of the Mechanics and Physics of Solids, 13(4): 223-227. DOI:10.1016/0022-5096(65)90011-6 |
Carcione J M, Helle H B, Avseth P. 2011. Source-rock seismic-velocity models:Gassmann versus Backus. Geophysics, 76(5): N37-N45. DOI:10.1190/geo2010-0258.1 |
Ceron M R, Martinez J F, Diaz E, et al. 2013. Digital rock physics for reservoir characterization.//13th International Congress of the Brazilian Geophysical Society & EXPOGEF. Rio de Janeiro, Brazil: Society of Exploration Geophysicists and Brazilian Geophysical Society, 1134-1136.
|
Chareyre B, Villard P. 2005. Dynamic spar elements and discrete element methods in two dimensions for the modeling of soil-inclusion problems. Journal of Engineering Mechanics, 131(7): 689-698. DOI:10.1061/(ASCE)0733-9399(2005)131:7(689) |
De Arcangelis L, Hansen A, Herrmann H J, et al. 1989. Scaling laws in fracture. Physical Review B, 40(1): 877-880. DOI:10.1103/PhysRevB.40.877 |
Deng J X, Shi G, Liu R X, Yu J. 2004. Analysis of the velocity anisotropy and its affection factors in shale and mudstone. Chinese Journal of Geophysics (in Chinese), 47(5): 862-868. |
Deng J X, Wang H, Zhou H, et al. 2015. Microtexture, seismic rock physical properties and modeling of Longmaxi Formation shale. Chinese Journal of Geophysics (in Chinese), 58(6): 2123-2136. DOI:10.6038/cjg20150626 |
Dvorkin J, Walls J, Tutuncu A, et al. 2003. Rock property determination using digital rock physics.//73rd Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1660-1663.
|
Eshelby J D. 1957. The determination of the elastic field of an ellipsoidal inclusion, and related problems. Proceedings of the Royal Society of Series A:Mathematical, Physical and Engineering Sciences, 241(1226): 376-396. |
Guo Z Q, Li X Y, Liu C, et al. 2013. A shale rock physics model for analysis of brittleness index, mineralogy and porosity in the Barnett Shale. Journal of Geophysics and Engineering, 10(2): 025006. DOI:10.1088/1742-2132/10/2/025006 |
Hashin Z, Shtrikman S. 1962. A variational approach to the theory of the effective magnetic permeability of multiphase materials. Journal of Applied Physics, 33(10): 3125-3131. DOI:10.1063/1.1728579 |
Hill R. 1965. A self-consistent mechanics of composite materials. Journal of the Mechanics and Physics of Solids, 13(4): 213-222. DOI:10.1016/0022-5096(65)90010-4 |
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 |
Hrennikoff A. 1941. Solution of problems of elasticity by the framework method. Journal of Applied Mechanics, 8(4): 169-175. |
Hudson H M, Larkin R S. 1994. Accelerated image reconstruction using ordered subsets of projection data. IEEE Transactions on Medical Imaging, 13(4): 601-609. DOI:10.1109/42.363108 |
Ikelle L T, Yung S K, Daube F. 1993. 2-D random media with ellipsoidal autocorrelation functions. Geophysics, 58(9): 1359-1372. DOI:10.1190/1.1443518 |
Kewen W, Ning L. 2008. Numerical simulation of rock pore-throat structure effects on NMR T2 distribution. Applied Geophysics, 5(2): 86-91. |
Korn M. 1993. Seismic waves in random media. Journal of Applied Geophysics, 29(3-4): 247-269. DOI:10.1016/0926-9851(93)90007-L |
Lin L, Chen J, Zhang X, et al. 2011. A novel 2-D random void model and its application in ultrasonically determined void content for composite materials. NDT & E International, 44(3): 254-260. |
Liu J X, Deng S C, Zhang J, et al. 2007. Lattice type of fracture model for concrete. Theoretical and Applied Fracture Mechanics, 48(3): 269-284. DOI:10.1016/j.tafmec.2007.08.008 |
Liu N, Fu L Y, Tang G, et al. 2020. Modified LSM for size-dependent wave propagation:comparison with modified couple stress theory. Acta Mechanica, 231(4): 1285-1304. DOI:10.1007/s00707-019-02580-y |
Liu N, Fu L Y. 2020. Elastic characteristics of digital cores from Longmaxi shale using lattice spring models. Communications in Computational Physics, 28(1): 518-538. DOI:10.4208/cicp.OA-2018-0085 |
Liu X F, Sun J M, Wang H T. 2009. Numerical simulation of rock electrical properties based on digital cores. Applied Geophysics, 6(1): 1-7. |
Lucier A M, Hofmann R, Bryndzia L T. 2011. Evaluation of variable gas saturation on acoustic log data from the Haynesville Shale gas play, NW Louisiana, USA. The Leading Edge, 30(3): 300-311. DOI:10.1190/1.3567261 |
Milauer V, Catalano E, Chareyre B, et al. 2010. Yade Reference Documentation. Yade Documentation.
|
Mori T, Tanaka K. 1973. Average stress in matrix and average elastic energy of materials with misfitting inclusions. Acta Metallurgica, 21(5): 571-574. DOI:10.1016/0001-6160(73)90064-3 |
Ostoja-Starzewski M. 2002. Lattice models in micromechanics. Applied Mechanics Reviews, 55(1): 35-60. DOI:10.1115/1.1432990 |
Qin X, Han D H, Zhao L X. 2014. Rock physics modeling of organic-rich shales with different maturity levels.//84th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 2952-2957.
|
Ray P, Chakrabarti B K. 1985. A microscopic approach to the statistical fracture analysis of disordered brittle solids. Solid State Communications, 53(5): 477-479. DOI:10.1016/0038-1098(85)91061-0 |
Reuss A. 1929. Berechnung der fließgrenze von mischkristallen auf grund der plastizitätsbedingung für einkristalle. ZAMM-Journal of Applied Mathematics and Mechanics, 9(1): 49-58. |
Saenger E H, Krüger O S, Shapiro S A. 2005. Numerical considerations of fluid effects on wave propagation.//Krause E, Jäger W, Resch M eds. High Performance Computing in Science and Engineering. Berlin, Heidelberg: Springer, 385-394.
|
Sahimi M, Goddard J D. 1986. Elastic percolation models for cohesive mechanical failure in heterogeneous systems. Physical Review B, 33(11): 7848-7851. DOI:10.1103/PhysRevB.33.7848 |
Sondergeld C H, Rai C S, Margesson R W. 2000. Ultrasonic measurement of anisotropy on the Kimmeridge Shale.//70th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 1858-1861.
|
Sondergeld C H, Rai C S. 2011. Elastic anisotropy of shales. The Leading Edge, 30(3): 324-331. DOI:10.1190/1.3567264 |
Vernik L, Nur A. 1992. Ultrasonic velocity and anisotropy of hydrocarbon source rocks. Geophysics, 57(5): 727-735. DOI:10.1190/1.1443286 |
Vernik L, Landis C. 1996. Elastic anisotropy of source rocks:Implications for hydrocarbon generation and primary migration. AAPG Bulletin, 80(4): 531-544. |
Vernik L, Liu X Z. 1997. Velocity anisotropy in shales:A petrophysical study. Geophysics, 62(2): 521-532. |
Voigt W. 1889. Ueber die Beziehung zwischen den beiden Elasticitätsconstanten isotroper Körper. Annalen der Physik, 274(12): 573-587. DOI:10.1002/andp.18892741206 |
Yan F Y, Han D H. 2013. Measurement of elastic properties of kerogen.//83rd Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts, 2778-2782.
|
Zhang G Z, Chen J J, Chen H Z, et al. 2015. Prediction for in-situ formation stress of shale based on rock physics equivalent model. Chinese Journal of Geophysics (in Chinese), 58(6): 2112-2122. DOI:10.6038/cjg20150625 |
Zhang W H, Fu L Y, Zhang Y, et al. 2016. Computation of elastic properties of 3D digital cores from the Longmaxi shale. Applied Geophysics, 13(2): 364-374. |
Zhao L X, Qin X, Han D H, et al. 2016. Rock-physics modeling for the elastic properties of organic shale at different maturity stages. Geophysics, 81(5): D527-D541. DOI:10.1190/geo2015-0713.1 |
Zhou H, Fu L Y. 2018. Scattering and intrinsic components of attenuation through the spectral ratio method in ultrasonic laboratory experiment. Chinese Journal of Geophysics (in Chinese), 61(3): 1083-1094. DOI:10.6038/cjg2018L0359 |
邓继新, 史謌, 刘瑞珣, 等. 2004. 泥岩、页岩声速各向异性及其影响因素分析. 地球物理学报, 47(5): 862-868. DOI:10.3321/j.issn:0001-5733.2004.05.018 |
邓继新, 王欢, 周浩, 等. 2015. 龙马溪组页岩微观结构、地震岩石物理特征与建模. 地球物理学报, 58(6): 2123-2136. DOI:10.6038/cjg20150626 |
张广智, 陈娇娇, 陈怀震, 等. 2015. 基于页岩岩石物理等效模型的地应力预测方法研究. 地球物理学报, 58(6): 2112-2122. DOI:10.6038/cjg20150625 |
周浩, 符力耘. 2018. 超声实验中谱比法衰减的散射与本征吸收特性. 地球物理学报, 61(3): 1083-1094. DOI:10.6038/cjg2018L0359 |