刘学锋, 刁庆雷, 孙宝佃, 李国利, 曾文冲

(1.中国石油大学理学院, 山东 青岛 266580; 2.中国石油集团测井有限公司油气评价中心, 陕西 西安 710077; 3.中国石化胜利油田, 山东 东营 257001)

0 引 言

岩石物理数值模拟已经成为岩石物理研究的重要手段之一[1-2]。与岩石物理实验相比,数值模拟有助于在微观尺度上定量研究各种因素对岩石物理属性的影响[3],为改进现有岩石物理模型或建立新的模型提供理论依据。

三维数字岩心在孔隙尺度上刻画了岩石的微观结构,是岩心的三维数字化图像,已经成为岩石物理数值模拟的基础。构建三维数字岩心的常用方法可分为2大类:三维物理成像[4-5]和基于二维图像的重建算法[6]。重建算法又分为过程法[7]和随机法[8-9]。过程法结合岩石的颗粒粒径分布,通过对沉积类岩石形成过程的模拟建立数字岩心。在成岩作用的模拟中,考虑了石英胶结质的生长和黏土物质的填充作用,因此过程法只适用于成岩过程简单的岩石。随机法以岩石二维图像统计特性为约束函数,采用随机算法构建三维数字岩心,其统计特性与二维图像相同。根据选用的随机算法和统计特性,随机法又分为模拟退火法、高斯场法和顺序指示模拟法等,选用的约束函数大都反映二维图像的两点统计特性,不能精确反映孔隙的长程相关性,导致重建三维数字岩心与真实岩心在孔隙连通性方面具有明显差异。Okabe等[10]采用多点地质统计法基于单方向二维切面图像构建重建三维数字岩心,该算法通过旋转单方向二维切面获得三维情况下的条件概率分布函数,未采用硬数据约束重建过程,该方法仅适用于重建各向同性的三维数字岩心。张伟等[11-12]也采用多点地质统计学在岩石二维图像基础上重建三维数字岩心。大量重建实例说明仅采用单方向二维切面图像,考虑1个方向的条件概率分布函数,重建三维数字岩心的孔隙仅在该方向上具有良好的孔隙连通性,而其余2个正交方向与实际岩心差别较大。

为解决随机法重建三维数字岩心孔隙连通性差异大的问题,本文在岩心3个方向二维切面图像基础上,修正数据模板和条件概率函数选取方法,引入硬数据的约束,采用多点地质统计学重建三维数字岩心,并定量比较了重建岩心和X射线CT扫描得到真实岩心的自相关函数、局部孔隙度函数、渗流概率分布函数,评价重建结果的准确性。

1 多点地质统计法

1.1 基本原理

多点地质统计学最初被用在地质统计学中,用来处理油藏尺度下的连续地质实体。多点地质统计方法的基础是用训练图像代替了两点地质统计法中的变差函数。该方法的基本思路是从训练图像中提取特征的图像模式,然后将这些模式还原到最终的模型上。

训练图像原本是一个地质概念,是能够表述实际储层结构、几何形态及其分布的数字化图像。训练图像的代表性决定了模拟结果的准确性[13]。训练图像包含了待模拟区域的各种特征模式。以三维数字岩心建模为例,若将岩石划分为2种组分:孔隙空间和骨架,则训练图像可视为二值化的数字图像(见图1),图1中黑色和白色分别代表岩石骨架和孔隙。训练图像的特征可被在其上方滑动的窗口所捕获,该窗口被称为数据模板[14]。图1(a)为尺寸为3像素×3像素的正方形数据模板,u表示未取样点。数据模板的几何形态和尺寸是任意的,模板尺寸越大,越能反映长程相关性,但计算效率变低。利用数据模板扫描训练图像得到1个数据事件[见图1(c)]。当训练图像中的1个数据事件与数据模板的数据事件相同时,称为1个重复。扫描整个训练图像后即可获得待模拟点的条件概率分布函数,组成搜索树。然后,建立需重构的三维网格,将已知硬数据点部署到网格节点上。硬数据是对客观存在的事物或现象进行测量和观察的结果。最后,定义一条随机路径模拟所有网格节点,对路径上的每一个待模拟节点,根据已建立的条件概率分布函数,利用Monte Carlo方法确定该节点的随机模拟值,直至完成所有节点的模拟。

图1 扫描训练图像过程

1.2 基于三方向二维岩心切面重建三维数字岩心的基本过程

采用多点地质统计学构建三维数字岩心的基本思路与油藏地质建模类似,不同之处在于三维数字岩心的尺寸小,空间分辨率在微米尺度,主要步骤分为准备训练图像、选择数据模板建立搜索树、添加硬数据和随机模拟。本文提出的重建方法是在三方向岩石二维切面图像基础上构建三维数字岩心,修正了数据模板选取、硬数据提取和待模拟点条件概率分布函数计算方法,其基本过程分为5步。

(1) 建立二维训练图像。获取岩心X、Y、Z这3个方向的二维切面灰度图像,经过灰度图像分割获得岩心3个方向的二值化图像。

(2) 选择合适的数据模板,确定条件概率分布函数,建立搜索树。在3个相互垂直面上分别选取相应数据模板,模板尺寸和结构根据二维平面图像结构确定。数据模板尺寸根据二维图像的自相关函数确定。为了在随机模拟过程中充分利用已经完成模拟点和硬数据点的数据,调整模板中待模拟点的位置,X、Y、Z这3个方向数据模板结构示意图如图2所示,图2中显示模板尺寸为4像素×3像素,待模拟点u位置放在了模板边界上,而非传统的模板中心位置。采用相应模板扫描3张训练图像,分别建立3个方向的搜索树。

图2 3个方向图像所用数据模板

(3) 合理选取硬数据。硬数据是测量和观测的结果,所在的网格节点具有明确的属性值,不需要随机模拟产生。在三维数字岩心随机模拟过程中应充分考虑硬数据的作用,如某一XY平面中存在大尺寸孔隙,则与其相邻的平面内相近位置处也必定为孔隙。但如果硬数据数目过大,则待模拟层不能产生预期的变化或者变化很小,所以选用多尺度算法选取硬数据[15]。

(4) 随机模拟。建立起需要重建的三维网络,定义一条访问待模拟节点的随机路径,获取条件概率。分别在3个相互垂直的面上利用以上数据模板扫描得到待模拟点数据事件,从相对应搜索树上检索到各自的条件概率分布函数。利用3个方向在该点处的条件概率分布函数加权平均即可得到待模拟点u(i,j,k)的总条件概率,计算公式为

p(xijk|x(N30(ijk)))=α{p(xijk|x(Ni,12(ijk)))+p(xijk|x(Nj,12(ijk)))+p(xijk|x(Nk,12(ijk)))}

(1)

式中,p(xijk|x(Ni,12(ijk)))表示u在YZ平面上数据模板大小为12时的概率分布函数;α是平均参数。

Monte Carlo模拟从待模拟点处的条件概率分布中提取一个值作为该处的随机模拟值。更新条件数据,将该模拟值加入到条件数据集中。

(5) 沿已定义的随机路径访问下一个待模拟节点,并重复(3)～(4)步骤直至确定所有待模拟节点的数值,从而产生一个随机模拟实现。

2 实例分析

X射线CT是三维无损成像方法,利用不同岩石组分对X射线吸收系数的差异,获得岩石内部的三维图像。在分辨率足够高的前提下,X射线CT建立三维数字岩心与真实岩心孔隙结构相同。为了定量对比真实岩心与重建数字岩心,选取某一枫丹白露砂岩X射线CT获得3个相互垂直二维切面作为训练图像(见图3)。图3中黑色表示岩石骨架,白色表示岩石孔隙,分别用1和0表示,XY、XZ和YZ这3个切面的孔隙度分别为19.92%、20.81%和21.13%,平均孔隙度为20.62%。根据二维图片的自相关函数,确定数据模板尺寸为11像素×11像素。

图3 训练图像

通过多点地质统计学建立三维数字岩心,采用数学形态学方法将孔隙和骨架边界平滑,剔除孤立骨架点后,岩心切面见图4,白色为孔隙空间,黑色为岩石骨架,孔隙度为19.21%。为了对比基于单一方向岩心二维切面的重建结果,采用Alireza Hajizadeh提出的方法在XY切面图像基础上重建三维数字岩心[15],岩心切面如图4所示,孔隙度为18.86%。实例2基于3张图像重建岩心孔隙度为14.46%,基于1张图像重建岩心孔隙度为14.62%,真实岩心孔隙度为13.15%(见图5)。

图4 重建实例1

图5 重建实例2

通过直观对比2种方法重建的三维数字岩心发现,利用本文提出的基于三方向二维切面重建的三维数字岩心孔隙分布比较均匀,尺寸和形状与训练图像类似,3个方向孔隙发育情况基本相同。但仅采用1个方向切面图片重建的三维数字岩心仅在Z方向体现了较为明显的空间相关性,但其余2个方向差别较大,导致三维孔隙尺寸和形状与实际岩心存在较大差异。

3 准确性评价

重建结果表明,由于岩心二维图像包含的信息量较少,重建的三维数字岩心在孔隙结构上与真实岩心并不完全吻合。为了定量评价重建结果的精度,分别计算三维数字岩心的自相关函数、局部孔隙度分布函数、平均渗流概率函数,分析孔隙空间的形状、均质性和连通性等特征。

3.1 自相关函数

岩心孔隙空间的自相关函数表示岩心中相距为r任意2点分布在孔隙空间中的概率,计算公式为

岩心孔隙空间的自相关函数反映了孔隙空间的形状和尺寸,X方向自相关函数计算结果见图6。3条曲线依次为X射线CT扫描得到真实岩心、基于3个方向二维切面图像重建和基于单方向二维图像重建的三维数字岩心。当滞后为0时,3条曲线的起点较为接近,说明3种三维数字岩心的孔隙度较为接近。基于3个方向二维切面重建的三维数字岩心与X射线CT构建的三维数字岩心具有相似的自相关函数,说明该方法重建的三维数字岩心孔隙尺寸与真实岩心类似。而采用1张图像重建三维数字岩心的自相关函数曲线下降速度快,说明该方法构建的三维数字岩心孔隙尺寸相对较小。

3.2 局部孔隙度分布函数

局部孔隙度分布函数反映岩石孔隙空间的均质性,其核心思想:将整个数字岩心划分为多个小尺寸测量单元,计算测量单元的孔隙度,获得岩心孔隙度的分布曲线。计算公式为

μ(φ,L)=1m∑rδ[φ-φ(r,L)]

(3)

式中,m为测量单元的个数;φ(r,L)为孔隙度;δ(x)为狄拉克函数。若测量单元的大小L一定,则数字岩心的均质性与曲线的宽度成反比,即曲线开口越小,均质性越好。当测量单元大小L=50时,局部孔隙度分布函数见图7。图7中3条曲线的峰值均集中在0.2附近,表示2种方法构建的三维数字岩心孔隙度均在20%左右,与X射线CT扫描得到真实三维数字岩心相符。3条曲线开口较小,且宽度相近,说明2种重建的三维数字岩心与真实岩心具有相似的均质性。

图6 自相关函数

图7 局部孔隙度分布函数

3.3 平均渗流分布函数

渗流分布函数反映了孔隙连通性,其基本思想:通过在某一方向上检测测量单元K(r,L)渗透性反映重建数字岩心的孔隙连通性[16]。渗流特征函数具体定义为

∧α(r,L)=1K(r,L)在α方向上有渗透性

0其他

(4)

α值可以取x、y、z、c、3。∧α(r,L)=1说明K(r,L)在α方向上存在逾渗通道,流体可以从一侧渗透到另外一侧。∧3=1说明K(r,L)在X、Y和Z这3个方向上全具有渗透性。∧c=1说明K(r,L)

至少在1个方向上具有渗透性。利用渗流特征可以获得局部渗流概率函数,经过积分后得到平均渗流概率函数pα(L)

(5)

平均渗流概率函数可表示数字岩心孔隙连通性的优劣,孔隙连通程度与曲线的斜率成正比,即曲线上升的越快,连通性就越好。

三维数字岩心的平均渗流概率函数计算结果见图8。该函数为至少在1个方向具有渗透性的渗流概率。3条曲线中,X射线CT扫描得到真实岩心和基于3张二维图像构建的数字岩心上升速度基本类似,而基于1张图片构建的三维数字岩心曲线上升速度明显偏低,说明本文提出算法的重建结果与真实岩心具有相似的孔隙连通性,而基于1张二维图片重建岩心的孔隙连通性较差。

4 结 论

(1) 在三方向二维切面图像的基础上重建的三维数字岩心与X射线CT扫描得到真实三维数字岩心具有相似的自相关函数、局部孔隙度分布函数、平均渗流概率函数,说明了重建的数字岩心与真实岩心具有相似的孔隙尺寸、均质性和孔隙连通性。同时引入3个方向的条件概率分布函数,改善了随机法重建三维数字岩心的孔隙连通性。

(2) 在1张二维图像的基础上,仅考虑单一平面的条件概率分布函数,采用多点地质统计学重建的三维数字岩心,在单一方向上具有较好的孔隙连通性,但其他2个方向的孔隙连通性较差,这是由于仅考虑某一方向条件概率分布函数,重建结果无法反映岩石真实孔隙空间特征。

参考文献:

[1] 孙建孟, 姜黎明, 刘学锋, 等. 数字岩心技术测井应用与展望 [J]. 测井技术, 2012, 36(1): 1-7.

[2] Liu X F, Sun J M, Wanga H T. Reconstruction of 3-D Digital Cores Using a Hybrid Method [J]. Applied Geophysics, 2009, 6(2): 105-112.

[3] 姚军, 赵秀才, 衣艳静, 等. 储层岩石微观结构性质的分析方法 [J]. 中国石油大学学报: 自然科学版, 2007, 31(1): 80-86.

[4] Dunsmuir J H, Ferguson S R, Damico K L, et al. X-ray Microtomography: a New Tool for the Characterization of Porous Media [C]∥SPE22860, 1991.

[5] Coenen J, Tchouparov A E, Jing X. Measurement Parameters and Resolution Aspects of Micro X-ray Tomography for Advanced Core Analysis [C]∥Proceedings of International Symposium of the Society of Core Analysts, 2004.

[6] Hazlett R D. Statistical Characterization and Stochastic Modeling of Pore Networks in Relation to Fluid Flow [J]. Mathematical geology, 1997, 29(6): 801-822.

[7] Bakke S, Oren P E. 3D Pore-scale Modeling of Sandstones and Flow Simulations in the Pore Networks [J]. SPE Journal, 1997, 2(2): 136-149.

[8] Hidajat I, Rastogi A, Singh M. Transport Properties of Porous Media from Thin Section [J]. SPE Journal, 2002, 7(1): 40-48.

[9] 刘学锋, 孙建孟, 王海涛, 等. 顺序指示模拟重建三维数字岩心的准确性评价 [J]. 石油学报, 2009, 30(3): 391-395.

[10] Okabe H, Blunt M J. Pore Space Reconstruction Using Multiple-point Statistics [J]. Journal of Petroleum Science and Engineering, 2005, 46: 121-137.

[11] 张伟, 林承焰, 董春梅. 多点地质统计学在秘鲁D油田地质建模中的应用 [J]. 中国石油大学学报: 自然科学版, 2008, 32(4): 24-28.

[12] 张挺. 基于多点地质统计的多孔介质重构方法及实现 [D]. 合肥: 中国科学技术大学计算机科学与技术学院, 2009.

[13] Steebelle S. Sequential Simulation Drawing Structures from Training Images [D]. California: Department of Geological and Environmental Science, Stanford University, 2000.

[14] Zhang T F, Switzer P, Journel A G. Filter-based Classification of Training Image Patterns for Spatial Simulation [J]. Mathematical Geology, 2006, 38(1): 63-80.

[15] Alireza H, Aliakbar S, Farhad A, et al. A Multiple-point Statistics Algorithm for 3D Pore Space Reconstruction from 2D Images [J]. Advances in Water Resources, 2011, 34: 1256-1267.

[16] 张丽, 孙建孟, 孙志强, 等. 多点地质统计学在三维岩心孔隙分布建模中的应用 [J]. 中国石油大学学报: 自然科学版, 2012, 36(2): 105-109.