第1章绪论
1.1 岩石物理
石油的需求日益增加及石油勘探开发技术的发展促使致密砂岩、页岩及煤层一类的非常规油气储层成为近年来的研究热点。岩石物理实验一直以来都是对储层岩石物理性质进行评价的重要手段，对油气资源的勘探开发具有重要的意义。然而传统的岩石物理实验应用在复杂储层及非常规储层时却出现了新的难题，如对于非均质性较强的一类储层难以取到具有代表性的岩心、难以对低渗透率的岩石进行驱替等。此外，传统的岩石物理实验也难以研究岩石微观组分对岩石宏观物理性质的影响。
近年来，随着计算机技术的发展，数字岩心成为岩石物理实验发展的一个热门方向。数字岩心技术能够利用计算机通过岩心的微观结构成分在计算上重构出一个可以反映岩石真实孔隙空间的岩心模型，重构后的三维岩心模型通过数值模拟计算等方法可以用来研究岩心的物理属性特征。
在计算机上利用重构的三维数字岩心模型研究岩石的物理属性特征具有很多传统岩石物理实验不可比拟的优势。三维数字岩心模型建立后可以反复使用，可以用同一块数字岩心研究不同的岩石物理属性，如电学、声学和放射性属性特征等，能够比传统岩石物理实验具有更好的统一内在性。数字岩心技术不仅可以用来计算传统岩石物理实验难以计算的一些物理量(如三相渗透率)，还能够很好地用来研究岩心微观组分对岩石物理宏观属性特征的影响。
自发展放射性测井技术以来，该技术已被大量用来评价岩性、孔隙度等储层物理性质并取得了良好效果。目前，国内外模拟实验中所构建的地层模型都较简单，与真实的复杂地层条件有较大差别，难以对中子仪器和地层进行基准校验。由于仪器刻度、地层复杂等，放射性测井在非常规储层会有偏差，通过数字岩心技术研究岩石微观组分对放射性属性特征的影响，能够为放射性技术在复杂储层的研究中起到积极的作用。
1.2 数字岩心
数字岩心建模的准确性直接决定了数值模拟的结果，只有在能够反映真实岩心状况的数字岩心上进行数值模拟，模拟结果才是可靠的。数字岩心建模的依据是根据各组分对岩石物理属性的影响来决定是否在数字岩心中将其表现出来。早期的数字岩心模型大多分为纯岩石骨架和孔隙两部分，但对于页岩这类非常规储层岩石，若要通过数字岩心计算页岩的电阻率，由于页岩骨架成分复杂，骨架成分中包含有对电阻率起影响作用的矿物组分如黄铁矿等，因此需要在对页岩这类岩石进行数字岩心建模时将骨架中的微观结构表现出来。
数字岩心建模的一个重要问题是建立多大尺寸的数字岩心，也就是建立的数字岩心能否反映真实岩心的结构，数值结果只能确保数值模型在相对较大的样本中具有代表性后才能与实验结果进行比较。Arns等[1]、Uribe等[2]和 Faisal等[3]强调了 RVE(代表体积元)尺寸效应在其研究中的重要性。 Shaina等[4]探讨了利用聚焦离子束扫描电子显微镜( FIB-SEM)建立的页岩能否找到代表体积元，模拟结果表明：当 FIB-SEM的图像分辨率低于 5000μm3时，不适合用建立的页岩数字岩心进行渗透率的模拟。Faisal等[3]在不同分辨率和不同孔隙尺寸下进行了岩石弹性模量的测量，研究了选取代表体积元的重要性。数字岩心的尺寸与数字岩心的分辨率是呈负相关的，并不是说建立的尺寸越大越能反映岩心的宏观物理属性，尺寸大了可能引起数字岩心失真，尺寸偏小则在一定程度上不能反映出岩石的各向异性等。
目前，数字岩心建模的方法主要有物理实验方法和数值重建方法两类 [4]。物理实验方法指借助高精度的仪器对岩心进行扫描，将扫描得到的二维图片进行三维重建。数值重建方法是在少量的岩心二维图像上，通过建模采用重建方法建立数字岩心。物理实验方法常用的仪器有光学显微镜、X射线计算机断层扫描成像 (X-CT)、FIB-SEM、宽离子束扫描电子显微镜（ BIB-SEM）等。数值重建方法包括孔隙网络模型、模拟退火算法、过程法、随机法与过程法相结合、顺序指示模拟技术、多点地质统计学方法、马尔可夫链-蒙特卡洛法等[5]。
1.3 岩石物理数值模拟
王克文等[6]、Berg等[7]和 Huo等[8]分别基于数字岩心利用格子 Boltzman方法对岩石内流体传输特性进行研究，同时求取岩石的相对渗透率和绝对渗透率。此外，Sun等[9]针对格子 Boltzman方法计算速度慢的问题，采用改进后并行格子 Boltzman方法求取岩石的渗透率。格子 Boltzman方法用于计算岩石渗透率的步骤是先对分割后的图像网格利用欧拉方程求取连通性指数。连通性指数取决于所求岩石的尺寸，用连通密度表示单位体积内的连通性，再识别岩石的连通性后用流体动力学的方法计算压力场和速度场，然后在有限体积内利用斯托克斯方程求取岩石的渗透率。测量的速度场用于计算流体通过孔隙的总通量；然后利用压力梯度和通量的值基于达西定律计算渗透率[6-8]。格子 Boltzman方法速度分量越多，模拟的结果越准确，但相应的运算量增大，运行速度降低。并行格子 Boltzman方法可用于解决运行速度慢的问题，并行格子 Boltzman方法首先将整个流域分为 N个子区域，计算任务分配给 N个核心，并且将每个子任务的初始流场信息(密度、速度等)和几何信息(网格坐标等)初始化，然后从主程序调度子程序的每个核心开始子任务 [8]。在每次迭代过程中，所有子程序都是独立执行的，同时执行子任务间的数据通信，利用主程序来判断收敛规则和迭代是否完成，通过处理后可以计算结果。
聂昕等[10]在 2016年利用有限元方法基于数字岩心对页岩电学特性模拟进行了研究分析，Wiegmann和 Zemitis[11]在 2006年使用显式跳跃调和平均法研究了岩石的导热性特征，Giorvana等[12]在 2014年和孔强夫等 [13]在 2016年都使用随机游走法研究了岩石的电性特征，结果表明使用随机游走法与有限元方法结果相近。用有限元方法计算岩石的电阻率是将数字岩心的每个像素视为一个元素来处理(离散化)，每个元素映射到相对应的岩石组分的立方体元素上，且满足相应的边界条件，然后施加外加电场计算电压分布 [10， 11]。Knackstedt等[14]在 2007年用有限差分法计算岩石的电阻率，将岩心离散化后用渐进松弛法来处理相邻网格的高对比度电阻率值。显式跳跃调和平均法主要解决了在 X、Y、Z三个不同方向上给定的电位差异的情况[11]。随机游走法求解电阻率是根据岩石内不同流体状态时地层因素不同来计算不同含水饱和度下的电阻增大系数[13]。
基于数字岩心进行声学特性模拟主要是模拟计算岩石的弹性模量。由弹性动力学基础可知，弹性模量可用来计算声波速度。 Saenger[15]在 2008年、姜黎明 [16]在 2012年利用有限元方法对岩石的声学特性进行了模拟。此外， Andra等[17]在 2013年分别利用基于傅里叶的 Lippmann-Schwinger和动态脉冲求取了岩石的弹性模量，Press等[18]在 2007年利用广义麦克斯韦模型和位移应力旋转交错网格等方法计算岩石的弹性模量。用有限元方法计算岩石的弹性模量是将数字岩心的每个像素视为一个元素来直接处理，并假设每个元素的位移为其节点坐标的线性函数，然后根据计算不同方向产生的应力求取弹性模量[15]。这种方法适用于任何数量的不同物质成分。该算法能够计算三维微结构的均匀应变产生的局部应力。
Wang等[19]在 2016年基于 X-CT扫描的数字岩心成功进行了中子传输特性模拟，利用蒙特卡洛方法描述中子在数字岩心中的传输过程，结果表明在探测器上接收到的中子计数与岩心孔隙度的线性关系略微受到各向异性的影响，同时证明了该方法可以用于孔隙结构各向异性的分析。此外，他们还发现利用热中子进行传播特性研究时，阵列探测器的分辨率并不是越高越好。蒙特卡洛方法是一种随机或统计学实验方法，该仿真方法可以根据概率现实地描述物理实验理论。中子的状态可由位置、能量及运动方向来描述，由于相邻碰撞间的时间是极短的，因此可以认为中子每两次碰撞间隔间的状态是不变的，这也是蒙特卡洛方法用于模拟中子传输特性的关键基础。这样中子的状态可由一系列碰撞点的状态来表示，即 S0， S1，. ， SM.1， SM。蒙特卡洛方法模拟中子传输过程是一个由已知确定未知的过程。模拟过程主要分为两步，一是由源分布确定 S0，二是由 SM.1确定 SM。
Arns等[20-22]和 Guo等[23， 24]利用随机行走法基于数字岩心对岩石的核磁共振特性进行了模拟研究。利用随机行走法进行核磁共振特性模拟主要有两步，第一步是获取随机行走粒子的生命曲线，第二步是通过多指数的反演得到 T2谱。模拟过程是将大量粒子随机放在岩心孔隙中，如果粒子遇到骨架表面，将按照一定的概率消失掉 [20-24]，如果粒子遇到骨架表面后继续游走，将会被骨架表面反弹并昀终回到原来的位置，同时将时间更新。反复重复模拟过程获得随机行走粒子的生命时间曲线，对随机行走粒子的生命时间曲线进行多指数的反演就可获得 T2谱。
参考文献
第2章数字岩心建模方法
三维数字岩心重建方法主要有数值重建方法和物理实验方法 [1]。数值重建方法是以少量岩心薄片的二维图像为基础，通过图像处理得到统计数据，采用某些数学方法重建三维数字岩心。物理实验方法建模是利用 X射线计算机断层扫描成像、聚焦离子束扫描电子显微镜、高倍光学显微镜、扫描电子显微镜等仪器，采集来自地下不同位置的岩心，通过重建算法得到二维图像从而得到岩心图像。
2.1 基于 CT数字岩心
用基于 X射线的 CT扫描技术建立数字岩心主要包括三个步骤：一是对岩样进行预处理，之后用 CT扫描获得图像数据；二是选择一种图像重建方法根据获得的图像数据重建数字岩心灰度图像；三是用二值分割方法分割灰度图像中的骨架和孔隙空间，从而建立数字岩心。
CT为计算机断层扫描(computed tomography)的简称，其全称为 X射线计算机断层扫描成像技术。CT技术起源于 1917年，奥地利数学家 Radon证明了可以通过无穷多个投影确定二维或三维物体的密度分布。在此基础上，Cormack证明了通过 X射线投影可以重构图像。这些理论都为 CT技术的发展奠定了基础。
20世纪 80年代研制出世界上第一台 CT仪器，主要应用于医学领域，之后 Dunsmui和 Coene等对 CT技术加以改进并引入石油开发领域中[2， 3]。目前，应用于石油领域的 CT仪器有两种：一种是台式 CT机，另一种是同步加速 CT机，二者的区别在于拥有不同的 X射线发生器。
CT设备发射 X射线穿过实验岩样，根据穿透的射线的衰减程度确定被探测物体的密度分布，从而获得岩心内部的结构信息。根据投影数据重构图像是 CT技术的关键，该方法的实现是用获得的投影数据，通过一定的计算方法确定衰减系数与样本中各组分空间位置的对应关系并构建图像，从而重构物体内部结构信息。
X射线衰减的原因是：当 X射线穿过岩心样品时，会与岩心内的原子发生物理反应而衰变，不同物质对 X射线的衰减系数不同。据此可以通过衰减系数的差异判断样品的内部结构组成。当 X射线穿过某一样品时，它所穿过的路径上的所有物质对衰减都有贡献，都反映在 X射线强度测量结果中，即
(2.1)
