作者：姚军,赵秀才
出版社：石油工业出版社
格式: AZW3, DOCX, EPUB, MOBI, PDF, TXT
数字岩心及孔隙级渗流模拟理论试读：
前言
提高原油采收率是石油工业面临的重要课题之一，从微观尺度研究提高原油采收率机理是发展大幅度提高采收率方法的重要途径。传统的多孔介质渗流理论以介质连续性假设和达西方程为基础，只能描述渗流的宏观规律（即多个孔隙和喉道共同流动时的统计特性），此理论体系为宏观渗流理论。如何在孔隙尺度上研究多孔介质的流动问题对于精细刻画微观流动机制、揭示大幅度提高采收率新的认知等具有重要意义，这种孔隙尺度（微米尺度）上的渗流理论为微观渗流理论。微观渗流理论和宏观渗流理论是两个平行的理论体系，共同构成了多孔介质的渗流理论体系。
利用孔隙级模拟方法研究流体在多孔介质中的流动是揭示其流动机制、获得新认知的重要方法，是目前国际学术界研究多孔介质流动模拟领域的热点。该方法不仅适用于油气田开发领域，还可以应用于环境工程、医学工程等相关领域，凡是涉及多孔介质中流体流动问题，均可应用此方法。
本书是作者近10年来在该领域系统研究的成果总结。全书共分七章。第一章简要介绍了微观渗流理论研究的重要意义、构建数字岩心及孔隙网络模型的必要性、重要性以及研究进展；第二章介绍了建立数字岩心和孔隙网络模型所需要的实验资料；第三章系统介绍了岩心微观图像的数学处理方法，详细介绍了三种可用来高质量分割灰度图像的二值化处理方法；第四章详细介绍了3种数字岩心的建模理论与方法，包括：CT扫描法、模拟退火法和过程模拟法；第五章介绍了用来分析数字岩心微观结构性质的理论与方法，分析了数字岩心的各项微观性质；第六章系统论述了以数字岩心为基础建立孔隙网络模型的理论与方法；第七章介绍了利用三维孔隙网络模型进行油水两相流动模拟的原理与方法，详细阐述了考虑渗流边界层的孔隙级渗流理论。
本书第二章、第五章、第六章、第七章由姚军编写，第一章、第三章、第四章由赵秀才编写，姚军负责全书的统稿工作。
在本书出版之际，向资助本项研究的国家自然科学基金委以及中国石油化工股份有限公司科技事业部表示衷心的感谢。
由于作者水平有限，书中错漏难免，恳请各位专家批评指正。作者2010．7第一章绪 论
当今，石油在世界各国的能源体系中均占据着至关重要的地位，同时作为重要的战略物资，其对国家安全的影响也十分重大。但我国现今石油开发、供给系统的实际情况却不容乐观：东部各大老油田相继进入高含水的后续开发阶段，常规开采技术已无法获得显著效果；西部油田地质条件十分复杂，现有理论在应用过程中产生许多问题，导致无法形成稳定的保障体系；陆地、海上发现新油田越来越困难，同时海上油田的开发风险相当大。可见，在我国经济发展及国防建设对石油资源的需求日益迫切、依赖性逐日增大的情况下而石油开发保障体系的基础并不牢固，因此形势十分严峻。
增大石油产量、保障石油供给的重要途径主要有两个：（1）借助勘探技术、手段不断探明新油田，增大开发储量，为石油开发提供坚强的后盾资源；（2）完善开发理论、改进开发技术、提高原油采收率，在可采储量一定的情况下尽最大努力将石油产量不断提高。如前所述，随着勘探技术的不断发展、勘探密度的逐日加大，国内大多数地质储量均已探明。因此，在我国石油资源储量一定的客观事实下，通过勘探技术增大石油开发后备资源的困难日益加大，所以，以第一种途径来稳定石油保障体系有很大风险。相比而言，第二种途径的可发展空间相当大。迄今，我国各大油田的平均采收率约为32%，仍有约68%的可采储量存于地下，可见，对石油可采储量的动用程度相对偏低，采收率仍有很大的提高空间。目前，我国石油的剩余可采储量84约33×10t，将油田采收率每提高一个百分点即可增油3300×10t，相当于在现有的开采水平下新添了一个中型油田，因此，效益十分可观，这将对保持我国石油供给保障体系的稳定继而推动国民经济的发展和国防体系的建设起到非常巨大的作用。因此，提高原油采收率迫在眉睫！第一节微观渗流理论研究的重要意义
迄今，采收率机理研究都是基于宏观渗流理论体系，宏观渗流理论均以达西方程为基础。由于达西方程是以宏观实验定律为基础的，因此在石油开发领域中应用达西方程建立的渗流方程也只能对油藏流体的运动行为进行宏观表征而无法给出孔隙级别上的精细刻画。众所周知，储层的许多宏观性质，如渗透率、毛细管压力、相对渗透率等，均取决于它的微观结构和组成它的固体及其孔隙空间中流体的物理性质，即微观结构、性质是根本，宏观性质、现象是表象。因此，仅仅通过对宏观性质的观测分析展开研究将永远停留在事物表面而无法对其内部的作用本质有所认识，由此所给出的研究结果也只能对提高采收率的技术手段提供参考，无法从根本上解释问题的产生原因以及其内部的作用机制，也就无法从根本上给出大幅度提高原油采收率的研究方向和研究思路。
由以上分析不难看出，要使原油采收率得以大幅度提高，其理论研究和技术开发决不能仅仅停留在宏观层次上，必须深入石油的储集、运移空间——多孔介质内部，从微观层面上开展研究。通过对多孔介质内部孔隙空间的发育规模、空间分布对流体渗流的影响、流体在其中的分布规律、相互作用机理等一些决定流体在多孔介质中流动的宏观现象的本质问题展开研究，才能够从根本上认识微观与宏观的联系，也只有在此基础上才能真正找到研究提高原油采收率技术的正确方向和应该采取的技术手段，为现场生产提供可靠有力的技术指导。第二节构建数字岩心及孔隙网络模型的必要性和重要性
目前，尽管多孔介质中流体渗流的一些微观机理已经可以通过实验来定性研究，但缺少必要手段和科学合理的方法进行定量描述，即所有的渗流理论仍停留在宏观尺度上，很多微观机理均无法考虑，如润湿性、边界层等。国内外已有学者通过各种方法建立了可以准确刻画岩石孔隙空间分布的数字岩心，它可认为将真实岩心的空间离散，即将真实岩心中的岩石骨架和孔隙均用体素（空间中的小立方体，其边长取决于成像设备的分辨率）表示。由于数字岩心能够给出岩心内部空间分布的精细描述，因而它为准确模拟流体在孔隙空间中的流动、分布等一系列渗流问题奠定了重要基础。
迄今，以数字岩心为基础开展流动模拟时均需采用格子Boltzmann方法，该方法的计算量对模拟介质的规模（即格子的数量）非常敏感，模拟时间随介质规模的增大而急剧增加。由于数字岩心中的每个体素的尺寸即采用格子Boltzmann方法开展流动模拟时对应的格子非常小，因此，若保证模拟介质（岩心孔隙空间）具有一定规模以便获取有代表性的流动规律，则流动模拟的运算量势必异常庞大，这将十分不利于复杂问题的研究。
因此，以数字岩心为基础，构建一种可以准确反映岩心孔隙空间的分布特征，同时又能够反映流体在岩心中渗流特征的孔隙网络模型是十分必要的。以这种孔隙网络模型为基础，不但可以方便地研究流体在其中的分布规律，计算流体的传导性能，还可以充分克服数字岩心计算时间长的缺点，因此它将是后续微观渗流规律研究和提高采收率技术研究的基础平台。所以，以数字岩心为基础快速准确地构建岩心孔隙网络模型具有重要的意义。
综上所述，数字岩心及孔隙网络模型是开展微观渗流规律研究所必需的平台。理论上，借助孔隙网络模型能够充分认识流体在多孔介质中的流动规律、作用机理等微观本质。实际中，以孔隙网络模型为依托，可以在深层次上展开各种研究：水驱过程中油水两相的流动规律和分布特征；聚合物以及各种化学驱过程中，在各种驱替剂的作用下原油在多孔介质中各种物理性质的转变以及由此引起的渗流规律的变化等。因此，岩心孔隙网络模型的开发具有重要的理论和实际意义。第三节数字岩心构建理论与方法研究进展
迄今，数字岩心的构建方法可分为两大类：物理实验方法和数值重建方法。物理实验方法均借助高倍光学显微镜、扫描电镜或CT成像仪等来获取岩心的平面图像，之后对平面图像进行三维重建即可得到数字岩心；数值重建方法则借助岩心平面图像等少量资料，通过图像分析提取建模信息之后采用若干数学方法来建立数字岩心。（1）物理实验方法建立数字岩心。
用以建立数字岩心的物理实验方法主要有序列成像法、聚焦扫描法和CT扫描法。序列成像法建立数字岩心的基本过程如下：将准备好的岩样抛光得到平整的岩样面，用高倍显微镜拍摄岩样抛光面获取岩心微观结构图像；平行于抛光面切除一层岩样薄片，将切割后的岩样面做抛光处理并用显微镜取像；重复切割、抛光、取像的实验过程直至获得一定数量的岩心图像为止；最后，将得到的实验图像分离、组合即可得到数字岩心。Lymberopoulos和Payatakes （1992）、Vogel和Roth（2001）、Tomutsa和Radmilovic（2003）、Tomutsa和Silin（2004）、Tomutsa等人（2007）曾用该方法建立了数字岩心并分析了岩心的孔隙结构（图1－1）。尽管序列成像法可以获得高分辨率的图像（电子显微镜可达到纳米级）；但是，由于岩心切割和抛光处理不仅需要大量时间而且容易破坏岩心孔隙结构，因此，该方法实用性较差。图1－1 Tomutsa等人采用序列成像法建立的石灰岩数字岩心（黑色为孔隙，灰白色为岩石骨架）
采用聚焦扫描微观成像技术获取数字岩心的基本原理如下：向岩心孔隙中灌注染色环氧树脂，借助聚焦扫描仪器对岩心切片逐点扫描并记录环氧树脂对光的反射强度，将强度信号转换为灰度图像从而得到岩心图像。Fredrich等人（1995）曾采用该方法建立了数字岩心（图1－2）。虽然岩心图像具有较高的分辨率，但是，该方法只能对岩心薄片进行成像处理，因此所得到的数字岩心呈薄片状，规模很小。图1－2 Fredrich等人采用聚焦扫描法建立的数字岩心（彩色表示孔隙，岩石骨架透明）
20世纪80年代初，Jim Elliott（1982）研制了世界上第一台CT机，之后CT技术首先在医学领域被广泛应用。Dunsmuir等人（1991）将CT技术加以改进并应用到石油开发领域，Coenen等人（2004）应用CT机对岩心进行扫描得到了分辨率小于1μm的岩心三维图像，CT扫描法建立数字岩心的基本原理将在第二章中系统介绍。目前，应用在石油科技领域的CT机可分为两类：台式CT机和同步加速CT机。尽管台式CT机能够获得分辨率达5μm或更低的图像，但文献中的大多数高质量图像都是通过同步加速CT机扫描得到的，如Conean（2004）获得分辨率为0．7μm的图像。澳大利亚国立大学的Knackstedt和Arns等人（2004）自行研制了一套CT扫描设备。该设备可以对直径达5cm的样品进行扫描并能得到分辨率小于2μm的图像。他们采用1283个节点的并行机来生成CT图像，尺寸为2048体素的三维图像约4小时即可生成（图1－3）。对于弱胶结疏松砂岩，应用台式CT机就可以清楚分辨孔隙结构，但是，对于含有大量微孔结构的致密砂岩或碳酸岩盐，只能采用同步加速CT机来研究孔隙结构。只是，同步加速CT机造价高，试验成本大，不便于大规模应用。图1－3 澳大利亚国立大学XCT实验室采用CT扫描法建立的数字岩心（灰色为孔隙，岩石骨架透明）（2）数值重建方法建立数字岩心。
数值重建法通常以岩心切片图像为基础，借助各种不同的统计方法或模拟岩石的形成过程来建立数字岩心。迄今已发展了多种重建方法，比较典型的有高斯模拟法、模拟退火法、过程模拟法、多点统计法和马尔科夫随机重建法。
高斯模拟法由Joshi于1974年提出，该方法以分析岩石薄片所得的统计资料为基础。首先随机产生一个由相互独立的高斯变量组成的数据集（称为高斯场），集合中所有变量组成的总体满足标准正态分布；之后对高斯场作线性变换，使最初独立的变量具有相关性，该过程中使用孔隙度和两点相关函数作为约束条件。最后通过非线性变换将高斯场转化为数字岩心。受计算工具的限制，Joshi实际只建立起二维岩心。此后，Quiblier（1984）将Joshi提出的算法由二维重建拓展到三维重建，借助其改进后的算法，Quiblier建立了第一个真正的数字岩心。Adler等人（1990）在Quiblier方法的基础上引入了周期性边界条件，使得建模效果有所改善，应用所提出的方法，他们建立了Fontainebleau砂岩的三维孔隙空间。1992年，他们又引入了傅里叶变换的方法使建模速度有所提高。Ioannidis等人在该领域作了更为细致的研究并于1995年指出：用快速傅里叶变换法取代普通傅里叶变换可以大大加快建模速度。
但是仅用孔隙度和两点相关函数无法充分反映孔隙空间的结构特征［Quiblier（1984）；Adler等人（1992）；Roberts（1997）；Adler和Thovert（1998）；Levitz（1998）；Liang等人（1998）；Yeong和Torquato（1998）；Ioannidis等人（1999）；Roberts和Torquato（1999）；Bekri等人（2000）；Ioannidis和Chatzis（2000）］。为更好地描述孔隙空间特征，Hilfer（1991）引入了孔隙度分布、渗流概率分布函数，Torquato和Lu（1993）、Levitz（1998）、Roberts和Torquato（1999）采用孔隙尺寸分布函数来反映孔隙结构特征。这些函数的使用提高了数字岩心的质量，但是孔隙空间的连通性依然较差。
1997年，Hazlett提出了另外一种随机法——模拟退火法，该算法中的系统能量对应着目标函数的取值，目标函数定义为重建介质的性质和待模拟介质的统计性质的差值平方和。通过对系统的不断更新，使之稳定性逐渐增强，从而最终得到数字岩心。模拟退火法较高斯模拟法的优势在于，在建立数字岩心时它可以将反映岩石的更多的信息考虑进来，从而使得所建立的模型与真实多孔介质更接近。1998年，Yeong和Torquato对模拟退火算法建立数字岩心的能力进行了研究，结果表明该方法具有很强的数字岩心建模能力，它可以将参考函数所包含的信息有效地转化到输出的结果系统中。此外，他们还分析了不同建模函数对建模结果的影响。Hidajat等人于2001年综合运用上述高斯模拟法和模拟退火法，他们以高斯模拟法的输出作为模拟退火法的输入，这样由于高斯模拟法计算所得的空间介质要远远优于随机介质，因此使模拟退火法的计算时间大大缩减，同时由于模拟退火法可以将更多的统计性质考虑进来，所以建立的模型更能反映真实孔隙空间。
与上述随机建模方法不同，Bryant等人（1992）提出了通过模拟岩石的地质成岩过程（包括沉积、压实和成岩）来建立数字岩心的方法。研究发现，只有考虑孔隙的空间相关性（即孔隙在空间中的分布不应是任意的，而是具有一定的相互关系）才能对岩石的传导性质做出正确预测。然而，Bryant等人所建立的模型存在巨大局限性，即：只有当所研究岩石中的颗粒尺寸与他们建模采用的圆球尺寸相等时，应用该模型才能做出很好的预测。此后，Bakke、ren等（1997、1998）在该领域做了更深入的研究并给出了一种能更加逼真的模拟真实岩心形成过程的建模方法，该方法建模时不仅考虑了岩石的颗粒粒径分布，还将其他一些通过薄片分析得来的岩石物理性质结合进来。在成岩作用的模拟中，Bakke和ren只考虑了石英胶结质的生长和黏土物质的填充作用。应用这种基于过程模拟的方法，他们建立了Fontainebleau砂岩的数字岩心（图1－4）。通过与该砂岩切片的数字化图像定量比较发现，该模型可以很好地重现真实岩石的几何性质和传导性质。此后，他们又对算法作了进一步改进（2003）：①作为模拟单元的岩石颗粒不仅有球形还有橄榄球形；②模拟黏土物质对孔隙的充填作用时，不仅考虑了长石还考虑了绿泥石和伊利石。此外，Pillotti和Coehlo等还研究了非球形成岩颗粒沉积的模拟方法，从而进一步拓展了采用过程模拟方法建立数字岩心的实用性。过程模拟方法能够建立各向异性的数字岩心且数字岩心具有很好的孔隙连通性；但是该方法无法模拟复杂孔隙系统的成岩过程，因此无法建立这类岩石（如具有微观孔隙和复杂成岩过程的碳酸盐岩）的数字岩心。图1－4 Fontainebleau砂岩CT图像及过程模拟法建立的数字岩心
Okabe和Blunt（2004，2005）提出了用以建立数字岩心的新方法：多点统计法（图1－5）。他们使用［9，9］模板统计并存储岩心切片图像中的孔隙空间结构特征，再把统计得到的信息充分反映到所建立的数字岩心中。他们建立的数字岩心具有良好的孔隙连通性，但是建模速度很慢。之后，Wu等人（2004，2006）开发出一种建立数字岩心的有效方法。该方法以马尔可夫随机滤网统计模型为基础，借助2点及5点邻域模板对孔隙与岩石骨架交接面的特征进行统计并将统计信息映射到所建立的数字岩心中。该方法建立的数字岩心具有良好的孔隙连通性，此外，它的建模速度很快。图1－5 Berea砂岩CT图像及Okabe采用多点统计法建立的该砂岩数字岩心（图中黑色为孔隙，岩石骨架透明）
与过程模拟法相比，上述其他方法可以对任意类型的岩石建模，不足之处是：他们建立数字岩心时都假设岩心是各向同性的。第四节孔隙网络模型构建理论与方法研究进展
迄今，研究人员建立的孔隙网络模型种类繁多，但就所建模型的拓扑性质而言，以往所有网络模型可分为两大类，即：规则拓扑孔隙网络模型和真实拓扑孔隙网络模型。（1）规则拓扑孔隙网络模型。
所谓规则拓扑孔隙网络模型是指组成孔隙网络模型的基本单元（孔隙、喉道）在平面（二维网络模型）上或空间（三维网络模型）中具有十分规整的排布形式。其中，最为典型也是最常用的就是表征孔隙、喉道的基本单元分别沿直角坐标系x、y、z三个主轴方向排布而构成的网络模型，该类模型如图1－6（a）所示。Fatt在1956年最初引入的孔隙网络模型就是这种具有规则拓扑结构的二维网络模型，自那时起直至20世纪90年代初，Bryant等人通过过程模拟法构建起拓扑结构与真实岩心孔隙空间更为逼近的网络模型［见图1－6（b）］的这段时间内，几乎所有网络模型都属于这里所提出的规则拓扑孔隙网络模型。尽管规则拓扑孔隙网络模型中各单元在空间中排布十分规整，但是由于孔隙、喉道单元的表征方法不同以及它们尺寸大小的赋值方法各异而使得该类模型仍然具有十分丰富的表现形式。图1－6 孔隙网络示意图
①规则拓扑孔隙网络模型中孔隙、喉道的几何表征方法。
最初，Fatt建立网络模型时引入了圆柱形毛细管并用其表征喉道，孔隙仅作为喉道的连结点。尽管这种圆柱形毛管网络适合于毛细管压力的计算，也可用于对润湿滞后等现象的定性解释，但是，由于每个毛细管中仅允许单相流体存在，故对于研究两相、多相流动以及润湿性等问题十分不利。因此，以后的学者们提出了其他形状的孔隙、喉道单元。发展较成熟的包括：球形孔隙［梅文荣、邓传光（1994）］、星形孔道［Mason、Morrow（1986）］、变尺寸星形孔道［Man、Jing（1999）］及截面为不规则三角形［Piri、Blunt（2002）］、正方形［Lopez等（2003）］的孔隙、喉道。其实，无论星形、三角形或正方形孔隙、喉道，只要其存在边角结构就能够为非湿相流过时润湿相的残留提供空间和条件，所以，尽管这些几何表征方法与真实岩心孔隙空间相差较大，但它们至少在解释两相、多相流及润湿性影响等问题时具有很好的应用效果。
②规则拓扑孔隙网络模型中孔隙、喉道尺寸的赋值方法。
网络模型中孔隙、喉道尺寸及喉道长度的赋值方法也是研究人员较为关注的问题。最初，Fatt采用完全随机的方式为喉道半径赋值，之后很多学者研究发现：多孔介质中的孔隙、喉道尺寸、喉道长度等均大致满足一定的分布规律，因而可用某种分布函数来表征。迄今，应用较多的主要有以下几种：对数正态分布［胡雪涛、李允（2000）］、Haring－Greenkorn概率分布［Sorbie等（1989）］、Rayleigh分布［Sorbie等（1989）］、威布尔概率分布、截断式威布尔分布［候健等（2005）］等。图1－7 Bethe网络示意图（Z＝5）
③Bethe网络。
Bethe网络（Mohammadi等，1990）是由交点和连线构成、具有无穷分支结构、没有闭环的网络，其特征参数是配位数Z，模型如图1－7所示。该网络空间中用连线来代表孔隙，连线具有一定的半径、体积和流动阻力；连线之间的交点视为没有体积、阻力，只是起到连接的作用。用Bethe网络表征多孔介质的一大优势在于：可以根据研究需要轻易调节网络的拓扑结构，此外，网络的传导性也很容易求出。Bethe网络的几何结构特征由连线的尺寸、体积分布来描述。尽管Bethe网络从表面上看具有较复杂的结构，但是该网络的拓扑结构仅由一个特征参数即配位数Z来控制，因此，这里将它划归到规则拓扑孔隙网络模型一类。（2）真实拓扑孔隙网络模型。
这类孔隙网络模型是在数字岩心的基础上建立的，它们具有与数字岩心孔隙空间等价的拓扑结构。数字岩心若采用物理实验方法建立，则与之对应的孔隙网络模型具有真实岩心孔隙空间的拓扑结构；若采用数值方法建立，则在其基础上建立的孔隙网络模型的拓扑结构与真实岩心的拓扑结构仍有差异，但由于模型的拓扑性质较规则网络模型有较大改进，故此将它们统称为真实拓扑孔隙网络模型。建立这类模型的方法有多种，包括：多向扫描法、居中轴线法、Voronoi多面体法、最大球体法。
①多向扫描法。
Zhao等人（1994）提出了多向扫描方法，通过对孔隙空间进行多方向切片扫描来搜索孔隙、喉道。来自不同方向的扫描切片相交处标志孔隙空间的局部最小位置，将它们定义为喉道。但是该方法很难准确探测孔隙。此后，Baldwin等人（1996）及Liang等人（2000）曾在缩减算法定义的孔隙空间居中轴线的基础上利用这种方法来测量孔隙空间的水利半径。
②孔隙空间居中轴线法。
岩心孔隙空间类似于贯穿于岩心内部的相互连通的空心管道（显然，这些管道的截面形状极不规则），这些空心管道的中轴线相互连接即构成了这里所谓的孔隙空间居中轴线。尽管孔隙空间居中轴线无法表征孔隙空间的特征，但它们准确保留了孔隙空间的拓扑结构，即孔隙空间的结构特征。以居中轴线为基础，通过对整个孔隙空间进行合理分割和简化便得到孔隙网络模型，如Lindquist等人（1996）、Lindquist和Venkatarangan（1999）、Sheppard等人（2005）、Prodanovic等人（2006）等将中轴线的节点定义为孔隙，将中轴线上孔隙空间的局部最小区域作为喉道来建立孔隙网络模型（图1－8、图1－9）。
孔隙空间居中轴线可以通过缩减算法［Thinning algorithm，Baldwin等人（1996）；Liang等人（2000）］或烈火模拟算法［Burning algorithm，Lindquist等人（1996）］运算得到。由于孔隙空间结构复杂，所建立的居中轴线中含有一些因边界噪音而引起的冗余枝节或在单个孔隙中存在多个居中轴线节点的错误结构，因此，在建立孔隙网络模型之前需要将冗余枝节删除［Venkatarangan（2000）；Shin等人（2005）］并将不合理的居中轴线节点合并［Sheppard等人（2005）；Shin和Lindquist（2005）］以提高模型的准确性。Knackstedt科研组对孔隙网络建模方法也作了研究，在Lindquist模型中已有的孔隙空间局中轴线的基础上，他们开发了孔隙、喉道分割及孔隙的合并等方法。由于他们仍在尝试开发更合理的建模技术，因此在公开发表的文献中并未详细介绍已有的建模方法。然而，他们指出：对孔隙进行适当的合并将提高模型对流动性质预测的准确性［Sheppard等人（2005）］。此外，Al－Raoush等人同样以孔隙空间居中轴线为基础建立了孔隙网络模型。他们提出了识别孔隙的新方法，并以此来分割孔隙空间建立网络模型。在他们模型中对孔隙的定义如下：以三条或更多条孔隙居中轴线的连接处作为孔隙中心；在孔隙中心放置圆球，逐步增大圆球半径直至球表面接触岩石骨架，此时的圆球半径即为孔隙半径。用所提出的方法，他们建立了玻璃珠模型及砂岩的孔隙网络模型。分析模型本身发现，所建立的孔隙网络模型具有良好的孔隙喉道的空间分布、尺寸大小及连通性；由于他们并未将结果与其他方法的建模结果或实验资料对比，因此，模型能否对单相、多相流动性质进行预测不得而知。尽管如此，他们研究发现，从岩心微观图像出发建立孔隙网络模型是一种有效的建模途径，但图像分辨率对建模结果影响很大，因此在这方面需要进一步研究。图1－8 Lindquist等人提取的真实岩心孔隙居中轴线图图1－9 Sheppard等人建立的孔隙网络模型
③Voronoi多面体法。
采用过程模拟方法建立的数字岩心，岩心内每个成岩颗粒的位置都是已知的，增大每个颗粒的半径直至整个孔隙空间被完全填充，在此过程中记录颗粒的交界点，所有交界点连接起来便得到Voronoi多面体。多面体中的顶点对应孔隙，顶点间的连线对应喉道，对孔隙、喉道的相关参数统计便得到孔隙网络模型［Bryant和Blunt（1992）；Bryant等人（1993）；Bryant和Raikes （1995）；ren等人（1998）；ren和Bakke（2002、2003）］，如图1－10所示。Delurue（2002）提出了一种可以在任意类型孔隙空间中建立Voronoi多面体继而生成孔隙网络的方法。Okabe 和Blunt（2004）用这种方法建立了Berea砂岩的网络模型，但是该模型的拓扑结构很差，因此，Voronoi多面体法只能应用到以过程模拟方法建立的数字岩心而不适用于普通数字岩心的网络建模。
④最大球体法。
Silin等人（2003）采用最大球体法对岩心孔隙空间进行表征，他们建立了等径球模型和Fontainebleau砂岩的孔隙网络模型并计算了压汞毛管力曲线，结果表明他们建立的模型能够给出合理结果（图1－11）。Silin模型的构建思路如下：对于孔隙空间中的任意一点，以它为中心放置一个微半径圆球，之后增大球半径直至球表面接触岩石壁面，这样便得到了以该孔隙点为中心且能够放置在孔隙空间中的最大球体。于是，岩心孔隙空间中充满了一系列相互交叠及包含的球体。将所有被包含的小尺寸球体删除并将剩余的相互交叠的球体分为主、仆两类：对于任意两个相交的球体，大半径球体为主球体，小半径球体为仆球体；由于一个球体可能同时与多个球体相交，因此，主球体可以是与之相交的更大球体的仆球体，同样，仆球体也可以是与之相交的更小球体的主球体。最后，用所有局部最大的主球体来表示孔隙，连接相邻孔隙的所有球体表征喉道。Silin孔隙网络模型对孔隙定义合理，但是其对喉道的定义十分模糊。尽管Silin并没有建立起孔隙网络模型，但这种最大球体法为以后的许多学者提供了重要的研究思路。图1－10 en等人建立的Bentheimer砂岩的孔隙居中轴线图及孔隙网络模型图1－11 Fontainebleau砂岩孔隙空间图示及Silin在其基础上建立的简易模型
Al－Kharusi和Blunt（2007）进一步发展了这种方法并用它建立了砂岩及碳酸盐岩的孔隙网络模型，他们对球体间的相互关系作了更全面的定义。Silin仅定义了主球体和仆球体两类，分别表示两个相交或相切球体中半径较大和较小的球体。Al－Kharusi和Blunt在此基础上引入了等径球体这一新类型，用来表示相同半径的两个相交或相切球体。显然，他们对球体的定义更全面、合理。但是Al－Kharusi和Blunt所采用的这种建模方法需要极大的内存，因此，这种方法仅能针对小尺寸数字岩心建立网络模型，此外，他们建立的网络模型中孔隙的配位数偏高。Hu Dong（2007）进一步完善了Al－Kharusi和Blunt的建模方法：改进了生成最大球体的方法，极大提高了建模速度；在描述各个球体的关系时，不再采用主、仆和等径球体的分类方法，而是建立一种树状结构，使孔隙喉道间的相互连通关系更加清晰（图1－12）。图1－12 砂岩CT图像及Hu Dong应用最大球体法建立的孔隙网络模型第五节本书内容
自2003年起，中国石油大学（华东）以姚军教授为首的课题组在中国石油化工股份有限公司科技部（孔隙尺度微观渗流机理与应用研究，项目编号P05056）和国家自然基金委（基于数字岩心的孔隙网络建模方法研究，项目批准号：90610015）的资助下在数字岩心、孔隙网络建模及流动模拟领域开展了大量研究并取得了卓有成效的研究成果。迄今，课题组已经建立了一套较完整的构建数字岩心和孔隙网络模型的系统的理论与方法体系，以所研究的理论方法为技术手段可以开展以下各项工作：（1）借助计算机图像处理技术对岩心微观结构图像进行处理，分析岩心的微观结构特征，获取反映孔隙空间分布、连通性等的重要信息。（2）能够分别采用高斯模拟方法、模拟退火方法、过程模拟方法和CT扫描方法建立数字岩心。（3）对各种数值算法构建的数字岩心，能够采用一系列的评价函数及格子Botlzmann方法等进行建模可靠性评价，从而优选数字岩心，为网络模型的构建奠定基础。（4）以数字岩心为基础，采用一系列理论方法，包括孔隙空间的优化处理、孔隙空间居中轴线系统的构建及优化、孔隙喉道的自动识别及合理分割，孔隙喉道相应参数的统计分析等，形成了构建孔隙网络模型的实用方法。
此外，课题组还对以数字岩心和孔隙网络模型为基础的孔隙级流体流动模拟技术开展了系统研究并形成了相应的理论方法。
本书将全面系统介绍数字岩心和孔隙级网络模型的构建理论与方法及其简单的应用。参考文献
［1］Lymberopoulos，D．P．，Payatakes，A．C．．Derivation of topological，geometrical，and correlational properties of porous media from pore－chart analysis of serial section data［J］．Journal of Colloid and Interface Science，1992，150（1）：61－80．
［2］Vogel H．J．，Roth K．Quantitative morphology and network representation of soil pore structure ［J］．Advances in Water Resources，2001，24（3－4）：233－242．
［3］Tomutsa L．，Radmilovic V．Focused ion beam assisted three－dimensional rock imaging at submicron－scale［C］．Proceedings of International Symposium of the Society of Core Analysts，2003，Pau，France．
［4］Tomutsa L．，Silin D．Nanoscale Pore Imaging and Pore Scale Fluid Flow Modeling in Chalk．Lawrence Berkeley National Laboratory：Paper LBNL－56266，2004，http：／／repositories．cdlib．org／lbnl／LBNL－56266．
［5］Tomutsa L．，Silin D．，Radmilovic V．Analysis of chalk petrophysical properties by means of submicron－scale pore imaging and modeling［J］．SPE Reservoir Evaluation＆Engineering，2007，10：285－293．
［6］Fredrich J．T．，Menendez B．，Wong T．F．Imaging the pore structure of geomaterials［J］．Science，1995，268：276－279．
［7］Elliott J．C．，Dover S．D．X－ray micro－tomography［J］．Journal of Microscopy，1982，126：211－213．
［8］Dunsmuir J．H．，Ferguson S．R．，D′Amico K．L．et al．X－ray micro－tomography：a new tool for the characterization of porous media［C］．Paper SPE 22860，Proceedings of 66th Annual Technical Conference and Exhibition of the Society of Petroleum Engineers，1991，Dallas，TX．
［9］Coenen J．，Tchouparova E．，Jing，X．Measurement parameters and resolution aspects of micro Xray tomography for advanced core analysis［C］．Proceedings of International Symposium of the Society of Core Analysts，2004，AbuDhabi，UAE．
［10］Arns C．H．，Knackstedt M．A．，Pinczewski W．V．et al．Virtual permeametry on microtomographic images［J］．Journal of Petroleum Science and Engineering，2004a，45：41－46．
［11］Joshi M．A class of stochastic models for porous media［D］．University of Kansas，Kansas：1974．
［12］Quiblier J．A．A new three－dimensional modeling technique for studying porous media［J］．Journal of Colloid and Interface Science，1984，98（1）：84－102．
［13］Adler P．M，Jacquin C．G．，Quiblier J．A．Flow in simulated porous media［J］．International Journal of Multiphase Flow，1990，16（4）：691－712．
［14］Ioannidis M，Kwiecien M，Chatzis I．Computer Generation and Application of 3－D Model Porous Media：from Pore－Level Geostatistics to the Estimation of Formation Factor［C］．SPE 30201，Proceedings of SPE Computer Conf，Houston，1995，185－194．
［15］Adler P．M．，Jacquin C．G．，Thovert J．F．The formation factor of reconstructed porous－media ［J］．Water Resources Research，1992，28：1571－1576．
［16］Roberts A．P．Statistical reconstruction of three－dimensional porous media from two－dimensional images［J］．Physical Review E，1997，56：3203－3212．
［17］Adler P．M．，Thovert J．F．Real porous media：Local geometry and macroscopic properties［J］．Applied Mechanics Review，1998，51：537－585．
［18］Levitz P．Off－lattice reconstruction of porous media：critical evaluation，geometrical confinement and molecular transport［J］．Advances in Colloid and Interface Science，1998，77：71－106．
［19］Liang Z．R．，Fernandes C．P．，Magnani F．S．et al．A reconstruction technique for three－dimensional porous media using image analysis and Fourier transforms［J］．Journal of Petroleum Science and Engineering，1998，21：273－283．
［20］Yeong C．L．Y，Torquato S．Reconstructing Random Media［J］．Phys．Rev．E，1998，57（1）：495－505．
［21］Yeong C．L．Y，Torquato S．Reconstructing Random Media．II．Three－Dimensional Media from Two－Dimensional Cuts［J］．Phys．Rev．E，1998，58（1）：224－233．
［22］Ioannidis M．A．，Chatzis I．，Kwiecien M．J．Computer enhanced core analysis for petrophysical properties［J］．Journal of Canadian Petroleum Technology，1999，38：18－24．
［23］Roberts A．P．，Torquato S．Chord－distribution functions of three－dimensional random media：approximate first－passage times of Gaussian processes［J］．Physical Review E，1999，59：4953－4963．
［24］Bekri S．，Xu K．，Yousefian F．et al．Pore geometry and transport properties in North Sea chalk ［J］．Journal of Petroleum Science and Engineering，2000，25：107－134．
［25］Ioannidis M．A．，Chatzis I．On the geometry and topology of 3D stochastic porous media［J］．Journal of Colloid and Interface Science，2000，229：323－334．
［26］Hilfer R．Geometric and dielectric characterization of porous－media［J］．Physical Review B，1991，44：60－75．
［27］Torquato S．，Lu，B．Chord－length distribution function for 2－phase random media［J］．Physical Review E，1993，47：2950－2953．
［28］Hazlett R．D．Statistical Characterization and Stochastic Modeling of Pore Networks in Relation to Fluid Flow［J］．Mathematical Geology，1997，29（6）：801－822．
［29］Hidajat I．，Rastogi，A．，Singh，M．et al．Transport Properties of Porous Media from thin sections ［J］．SPE Journal，March 2002，40－48．
［30］Bryant S．，Blunt M．Prediction of relative permeability in simple porous media［J］．Physical Review A，1992，46（4）：2004－2012．
［31］Bakke S．，ren P．－E．3－D pore－scale modelling of sandstones and flow simulations in the pore networks［J］．SPE Journal，1997，2：136－149．
［32］ren P．－E．，Bakke S．，Arntzen，O．J．Extending Predictive Capabilities to Network Models ［J］．SPE Journal，1998，3：324－336．
［33］ren P．E，Bakke S．Process based reconstruction of sandstones and predictions of transport properties［J］．Transport in Porous Media，2002，46（2－3）：311－343．
［34］ren P．E，Bakke S．Reconstruction of Berea sandstone and pore－scale modelling of wettability effects［J］．Journal of Petroleum Science and Engineering，2003，39（3－4）：177－199．
［35］Pillotti，M．Reconstruction of clastic porous media［J］．Transport in Porous Media，2000，41（3）：359－364．
［36］Coehlo D，Thovert J．F，Adler PM．Geometrical and transport properties of random packings of spheres and aspherical particles［J］．Physical Review E，1997，55：1959－1978．
［37］Okabe H．，Blunt M．J．Prediction of permeability for porous media reconstructed using multiplepoint statistics［J］．Physical Review E，2004，70：066135．
［38］Okabe H．，Blunt M．J．Pore space reconstruction using multiple－point statistics［J］．Journal of Petroleum Science and Engineering，2005，46：121－137．
［39］Wu K．，Nunan N．，Crawford J．W．et al．An efficient Markov Chain model for the simulation of heterogeneous soil structure［J］．Soil Science Society of America，2004．68：346－351．
［40］Wu K．，Dijke M．I．J．V．，Couples G．D．et al．3D stochastic modelling of heterogeneous porous media－applications to reservoir rocks［J］．Transport in Porous Media，2006，65（3）：443－467．
［41］Fatt I．The network model of porous media I．Capillary pressure characteristics［J］．Trans AIME，1956，207：144－159．
［42］Fatt I．The network model of porous media II．Dynamic Properties of a Single Size Tube Network ［J］．Trans AIME，1956，207：160－163．
［43］Fatt I．The network model of porous media III．Dynamic Properties of networks with tube radius distribution［J］．Trans AIME，1956，．207：164－181．
［44］梅文荣，邓传光．颗粒运移的计算机模拟［J］．西南石油学院学报，1994；16（1）：82－87．
［45］Mason G．，Morrow N．R．Meniscus displacement curvatures of a perfectly wetting liquid in capillary pore throats formed by spheres［J］．J．Colloid Interface Sci．1986，109（1）：46－56．
［46］Man H．N，Jing XD．Network modelling of wettability and pore geometry effects on electrical resistivity and capillary pressure［J］．J．Petroleum Science and Engineering，1999，24：255－267．
［47］Piri M．，Blunt M．Pore－scale modeling of three－phase flow in mixed－wet systems［C］．Proceedings of the SPE Annual Technical Conference and Exhibition，San Antonio，TX，September，2002．
［48］Lopez X．，Valvatne P．H．，Blunt M．J．Predictive network modeling of single－phase non－newtonian flow in porous media［J］．J．Colloid and Interface Science，2003，264（1）：256－265．
［49］胡雪涛，李允．随机网络模拟研究微观剩余油分布［J］．石油学报，2000．7，21（4）：46－51．
［50］Sorbie K．S．，Clifford P．J．，Jones E．R．W．The rheology of pseudoplastic fluids in porous media using network modeling［J］．J．Colloid and Interface Science．1989，130（1）：508－534．
［51］候健，李振泉，关继腾，王克文，陈月明．基于三维网络模型的水驱油微观渗流机理研究［J］．力学学报，2005，37（6）：783－787.
［52］Mohammadi S，Sorbie KS，Danesh A．et al．Pore－Level Modelling of Gas－Condensate Flow Through Horizontal Porous Media［J］．SPE 20479 1990：149－158．
［53］Zhao H．Q．，Macdonald I．F．，Kwiecien M．J．Multi－orientation scanning：a necessity in the identification of pore necks in porous media by 3－D computer reconstruction from serial section data［J］．Journal of Colloid and Interface Science，1994．162（2）：390－401．
［54］Baldwin C．A．，Sederman A．J．，Mantle M．D．et al．Determination and characterization of the structure of a pore space from 3D volume images［J］．Journal of Colloid and Interface Science，1996，181：79－92．
试读结束[说明：试读内容隐藏了图片]