通过数字岩心计算岩石中的流体性质

1、通过数字岩心计算岩石中的流体性质P.E. ?ren, S. Bakke, and H.G. Ruesl?tten 张杰译 翻译 杨晓宁 校对摘要在北海原油储层中， 我们运用地质学重建技术来形成复杂砂岩相的虚拟岩石。 计算机生 成岩石的有效性质(如流体渗透性、 电阻率和弹性模量) ，与这些由实际岩石的微地形图像 决定的岩石性质， 二者可以进行很好的对比。 利用网络建模技术， 我们为计算机生成岩石计 算了注水相对渗透率，并将其和实验数据进行了对比。实验岩心样品与测出的变化范围 0.3-0.7 的 Amott 指标在湿度上显示了很大的不同。这种变化的物理成因还有待研究，然而 我们的结果明显证明，对于

2、固定湿度孔隙水平分布条件，剩油饱和度和 Amott 指标随着初 始水饱和度的升高而升高。 这主要是由混合潮湿孔中残油膜的稳定性引起的。 当前的研究证 实：对于储集岩来说， 结合计算机生成岩石和数字计算来得出重要岩石和流体性质是可能的， 并且是可行的。简介由于计算机技术的飞速发展， 宏观岩石和传输性质现在可以直接来自于岩石微观结构的三维(3D)图像。可被用作数字计算的有效性质包括弹性模量(Arns et al., 2002),电阻率(？ren and Bakke, 2002),渗透率(Jin et al., 2004, Arns et al., 2004)， 核磁共振弛豫( NMR relaxa

3、tion) (?ren et al., 2002, Knackstedtet al.,2004)， 力学性质(Jin et al., 2003)和2-3相流体的构成关系(?ren et al., 1998, Lerdahl, et al., 2000, Valvatne and Blunt, 2004, Piriand Blunt, 2005)。因此,储集岩孔隙结构详细 3D表示的获得对石油工业来说是非常重要的。在过去的几十年里，储集岩微结构的3D图像的获得已经取得了一些进展。一种普遍使用的方法是随机重建 (Adler et al., 1990, Hazlett, 1997, Yeoung a

4、nd Torquato, 1998)。这种技术依赖于先对实际岩石2D微观图像进行测量，然后进行 3D模拟，最后匹配其统计特性。现今对微层析成像随机重建的定量对比已经显 示：尤其对于低孔隙度的样品而言， 随机性模型趋向于忽略孔隙空间的连通性及 渗透率 （Biswal et al., 1999, Manswartet al., 2002, ?ren and Bakke,2003）。近年来，为了达到将储集岩的微构造进行数字化表示， 基于地质学重建技术 的程序已经建立 （Bryant et al., 1993, Bakke and ?ren, 1997, ?ren and Bakke, 2002, J

5、in etal., 2003）。简而言之，这些技术是基于岩石形成过程的地质学直接模拟，即沉积 作用、压实作用和成岩叠加作用。 对于重建来说， 必要的输入参数是从二维薄片 的背散射（BSE）图像中提取的（?ren and Bakke, 2002）。在目前的工作中，我们运用了基于地质学的重建技术，生成了北海油储复杂 砂岩相的虚拟岩石。 我们计算了原料和计算机生成岩石的传输相关有效性质， 并 将它们与来自实际岩石标本进行X射线微形态（微CT）图像的原料和性质进行了 对比。网络建模技术被用来计算重建岩石的注水相对渗透性。 预知的相对渗透性 可以与实验数据进行对比。样品描述现有的储集岩样品为分选不好的砂

6、岩，平均粒径 d在140uni到900 之间， 孔隙率U在0.23到0.27之间，渗透率在1-10达西之间。沉积环境为冲刷河道。因此， 样品是非均质性的， 并且经历了复杂的成岩改造， 即自生粘土矿物的形成 （5-8%） 和斑块状碳酸盐岩胶结物 （见图 1）。对岩石类型的典型岩心分析信息是可用的， 包括Amott湿度测量、对三块岩心柱塞离心测量油的相对渗透率及对两块合成岩 心定态相对渗透率的测量。Amott湿度指数Iwo，对不同岩石类型来说变化范围在 0.3到0.7之间，平均值为 0.5。来自图1BSE （背散射）图像的孔隙度是0.236。薄片下256张这种图像的镶嵌 图被用来提取必要的输入参数

7、， 以期对样品进行地质重建。 通过引入粘土成分和 目标孔隙度上微小的变化， 已经在总共 10种岩石类型上实现。 关于重建法则和应 用程序详细的介绍在其它地方给出 （Bakke and ?ren, 1997, ?ren and Bakke, 2002, ?ren andBakke, 2003）。从一个典型的岩石样品（提取的样品直径为 10mm）中获 得了高分辨率的微CT图像。微层析成像是从澳大利亚国立大学获得的（Arns et al.,2004, Kn ackstedtet al., 2004）。微层析成像(microtomographic)的截面(记为MCT)和重建的标本(记为 PBM)在图2

8、中进行了对比。两块样品都是基于5123大小和5.24叩的分辨率的三 维像素。接下来，我们计算并比较了这些样品的有效材料和传输性质。 岩心柱塞 的测量是从其有效处给出的。13图1.不同类型的北海储层砂岩薄片的背散射 BSE图像。黑色-孔隙；深灰色-粘土;灰色-石英；浅灰色-长石；浅灰白色-碳酸盐胶结物。像素分辨率为3.31卩。传输性质本处连续方程的平均计算将显微结构的影响和宏观尺度上适用的有效物理 传输属性联系了起来。对于均质媒介中的传输属性来说，与con servable量有关的一般通量(如电流、压力或动量)与一般梯度(电场、张力、压力)和比例常数 Ke是线相关关系。对于弹力和电传导来说，这可

9、以写作F(x) = Ke (x) ?G(x),其中 F服从微分方程？ i F(x) = 0。在非均质媒介中，我们希望类似的线相关关系包含 在平均场中F(x)= Ke G (x) and及(1)Ke成为有效参v(x) = Ke g(x)其中 Ke = -k/ 卩分别对于弹力/传导率和流体渗透率来说。方括号表示平均化以使 数，v代表流体速率，g为实用压力梯度，k为绝对渗透性， 伪流体粘滞性图2. MCT微层析图(左边)和PBM样品重建(右边)，非均质的北海储层砂岩, 黑色代表孔隙，深灰色代表粘土，灰色带表石英，亮灰色代表长石，白色代表碳 酸盐岩胶结物。图像侧边长2.68mm，三位图像分辨率为5.2

10、4卩。渗透率不可压缩的牛顿流体的低雷诺系数流体服从恒稳态的斯托克斯方程2战 v = ? p ,(2)? v = 0(3)服从固壁上的边界条件v=0。v和p分别代表速率和压力。A D3Q19晶格波尔兹曼 算法(Jin et al., 2004)被直接用于解决数字化图像上的斯托克斯方程。定向绝对渗 透率k是由i轴常压梯度的应用决定的(i = x, y, z)。宏观通量是通过计算局部流体 速率体积平均值获得的，ki是由达西法则决定的。我们将平均渗透率k定义为定向 渗透率的算法平均值。形成因素(formation factor)关于恒稳态的传导率问题，局部控制方程就变成了拉普拉斯方程? J=0(4)J

11、= cw ?(5)服从固壁上的边界条件？n = 0。J为电流，c为孔隙中流体的电导率， 为电压，n为单位矢量正常到固壁。拉普拉斯方程的数解是通过有限差方法获得 的(?ren and Bakke, 2002)。定向的形成因素Fi被定义为有效电导率Fi = w/ i的反 面。我们定义平均形成因素F为定向相关形成因素的调和均值。弹性模量决定非均质体介质弹性习性的局部方程为弹性静力学的基础方程? t = 0(6)T= &其中 & = 1/2? d+(? d)T(7)其中，和分别表示压力和张力张量，d为位移场，C为刚性张量。上述方程 式被通过有限元方法，用线性弹性波方程的能量表示解决了。 周期性边界条件

12、被 用于数字化图像的外观。有效体积和剪切模量为假定均质的线弹性性质而计算得 出的。本构关系 本构关系，如毛细管压力和相对渗透率曲线，是由模拟计算机生成岩石的孔 隙网络表征中的两相位移(如初次排水，注水，二次排水)而决定的。孔隙网络 的拓扑学和连通性是由提取孔隙中的残余物决定的。 这是通过颗粒的最终扩大完 成的(Bakke and ?ren,1997)。Voronoi多面体的顶点定义了孔喉，同时，许多多面 体的边界确定了喉。以数学纲要作为基础，我们用标准图像分析技术直接测量了 全部孔喉和喉的大小和体积 (?ren and Bakke, 2003)。由于提取的孔隙网络与重建 的孔隙是一对一的对应关

13、系，我们没有引进任何配置或调整参数来匹配诸如孔隙 度和渗透率等宏观参数。在所有多相流体模拟中假定毛细管压力可以左右孔径大小。 模拟毛细管控制 位移的基础是流体在孔隙中的恰当分布。 对于两相流体来说，均衡流体分布受湿 度和毛细管压力的控制，在两相之间对于任何外加的压力差异，杨氏-拉普拉斯方程都可以应用。包括湿度影响在内的所有精确细节清晰而广泛的讨论都已经包 含在模拟中，并且之前已经介绍过了 (?ren et al., 1998, Patzek, 2001, ?ren and Bakke, 2003)。当缺乏粘滞性影响时，流体的移动经过了一系列的均衡流体形成过程。 在每 次形成过程中，我们计算了饱

14、和度，毛细管压力和相对渗透率。在所有孔和喉中， 通过大量拉平相饱和度，宏观相饱和态得以获得。简单地说，毛细管压力Pc就是与末次均衡构型变化有关的入口毛细管压力。 为了计算相对渗透率，网络中每个相态流体必须指明。对于层流来说，介于两个连接孔I和J之间的流体流动速率i由下面方程给出(8)其中是指Iij孔隙中心间距。假定有效电导率gi,ij为喉的电导率及连接的两个半孔的调和平均值。在不同几何学形状和不同流体结构的孔中，水压电导率的表达式来自斯托克斯方程的解（?ren et al., 1998, Patzek and Kristensen, 2001。）在每个孔中，我们调用了质量守恒定律另加=丿（9）

15、其中J通过所有喉连接到孔I。对于可以用一般技术解决的孔压力来说， 方程式（8） 和（9）产生了一组线性方程。最初，当网络充满水达到饱和时，绝对渗透率 k即被估计。通过计算的压力 场我们可以计算总流量的流速，从而可以运用达西法则计算绝对渗透率。 类似地， 相对渗透率可以进行计算得出。假定所有界面都进行了适当的界定， 每个相态中 的压力都是分别计算的。通过入口的流速总和决定了相的宏观流速， 即通过达西 法则计算所得的相的渗透率ki。相的相对渗透率由公式kri = ki/k给出。对于饱和 度每个0.025的变化，给出的位移过程的构成关系通过计算 kr和Pc建立的。结果和讨论对于微CT和计算机生成的孔

16、结构来说，两点相关函数C2的一致性非常好（见 图3）。提出的函数是方向函数的平均值。对于 MCT和PBM样品来说，C2的初始 斜率允许特定的表面面积分别具有 0.0205 和0.0192卩用1的估计值。C2 = 0的 特征长度L=0.72mm在两块样品中都发现了。 C2的末端表明样品中存在更长的距 离结构，这是由于粒径大小变化较大，粘土的存在以及碳酸盐胶结物。100500L i micro ns)Porosity3D局部孔隙度分布图3. MCT和PBM样品的平均两点相关函数（左边）和KL=0.39mm)(右边)表1.计算出的微CT和重建的样品的传输性质SamplePorosity西%FA(mD

17、)&(mD)底(mD)A(mD)MCT0.24012.716.612 213.65853432263345502PBM0.24310.810911.010.96515613260156221两种样品的3D局部孔隙度分布， 讥册，（Hilfer, 1991）见图3 （右边）。卩（L）, 指的是发现局部孔隙度 旳线尺寸L立方体单元中的经验概率。图3中显示的空隙 度分布是用L = 0.39 mm计算得出的。分布的宽度（孔隙度的差异）和峰值（可 能性最大的孔隙度）二者很相似。重建的样品起点值比较高表明更大的基质区域 出现得更频繁了。这与特征长度L* （Biswal et al.,1999）的计算是一致

18、的。计算得 出PBM样品的L*=0.84mm, MCT样品的L*=0.79mm。L*给出了可以适合基质空 间及可以作为最大颗粒大小测量的最大立方体的边长。对于全部的样品大小（即5123的三维像素，5.24 口的分辨率），表1列出了 计算过的方位、一般形成因素和绝对渗透率。PBM样品的传输性质完全是等方性的，同时，MCT样品在y轴上具有比较高的形成因素和较低的渗透率的一向等 方性（an-isotropic）。与MCT样品相比，重建的孔洞结构高估了一般绝对渗透率 13% （5503 mD vs. 6221 mD），低估了一般形成因素 20% （13.6 vs. 10.9）由于两种 样品的孔隙度是相

19、似的，这意味着 MCT样品的弯曲比PBM样品要大。实验方法 也证实了较大柱塞的绝对渗透率是4370 mD。重建的样品的孔隙网络表现出的渗 透率和形成因素计算值分别为6347 mD和 16.5。2500020000LDOOCi5000o.ibsoPwKitya图4. MCT和PBM样品随孔隙度变化而计算的电导率（左图）和绝对渗透率（右 图）。标出了对大小为1283和2563三维像素的副样计算。 MCT2M hlCTIZE rm 354 a mii2s.-1口口Q 040M20-rd它 mqSIQ00.10.2 nD.3D4D0.1QJ2 .O.J0.4PorosiKPok）聶（y图5. MCT和

20、PBM样品随孔隙度变化而计算的体积模量（左图）和剪切模量（右 图）。标出了对大小为1283和2563三维像素的副样计算。图4列出了关于非重叠副样计算的传输性质。立方体副样的侧边长为 0.67mm（1283三维像素）和1.34mm（2563三维像素），从而符合样品确定的相关长度 的范围（L = 0.72 mm for C2 = 0）。对于微CT和计算机生成孔隙结构来说，随着孔 隙度的变化，绝对渗透率的变化曲线是相似的。 随着孔隙度变化，计算的电导率 曲线显示出重建的孔隙结构比微层析图像具有略高的电导率，或预测了较低的形成因素。一个可能的解释可能在PBM样品的小孔隙的预测计算中可以找到。这 些孔对

21、电导率具有很大的影响，但对于绝对渗透率的影响却微乎其微。图5列出了体积和剪切模量的计算。结果再次显示了大小1283和 2563三维像素的副样，并指出了每个样品中包含的可变性。曲线走势表明空隙度接近线性， 并且孔隙度高的具有较大的分散性。MCT和PBM数据显示斜率明显不同。这导致了小的，但对重建孔隙结构弹性模量系统的低估，或反之亦然的对层析图像的高估。弹性模量强烈依赖于颗粒与颗粒接触的表现。不论MCT样品中图像处理（即阈值）的并非决定性的影响，但有迹象支持这种解释图6.模拟的和离心测量的油相对渗透率比较。模拟结果的初始水的饱和度为Swi也0.15 （左图）和Swi也0.29 （右图）。关于重建样

22、品10个实现的注水油和水的相对渗透率进行了计算。首先，通过模拟强水湿度条件下主要排水建立了初始水饱和度Swi。然后，在建议的湿度状态下（Iwo=0.5）,注水进行了模拟。油浸入带孔隙中40%变为油湿的目标Iwo达 到了。油湿孔和喉的前进接触角（a随机分布于iio $160之间，同时，对于水 湿孔和喉来说，（随机分布于20 （60之间。对于不同实现方法来说，计算的 Iwo在0.4和0.6之间变化。在图6 （左）中，将计算的油相对渗透率和离心测量的 油相对渗透率进行了比较。实验数据来自三个不同的柱塞。计算结果倾向于高估 测量的油相对渗透率，尤其是低油饱和度，模拟的剩余油饱和度Sorw比用实验方法得

23、出的小（0.15比0.25）。实验柱塞的Swi值介于0.26和0.33之间。这比模拟值Swi也0.15明显高，也比野 外观察所得的值高。这种偏差的原因很可能是由于岩心样品固结性较差的事实所 导致的实验困难导致的，从而很容易地在高转速时破坏了，而高转速是对达到低 的初始水饱和度所必须的。然而，众所周知，初始水饱和度对后来的注水具有非 常重要的影响，尤其对非水湿样品而言（Jadhunandan and Morrow, 1995。为了 允许与实验数据进行直接的比较，我们用类似于实验值的Swi值重做了所有的模拟。这是通过终止Swi也0.29时初次排水模拟来实现的。湿度的孔隙级别分布（即 解除角度和部分

24、油湿孔）与之前是一样的。结果见图6（右图）。测量的油相对渗透率值和计算值的一致性显著提高了。 湿度在不同等级孔隙中的分布与之前一样，计算的 Swi值增加了，与测量值相似。Swi值的增加主要由混合湿孔中残油膜的稳定性引起的。存在于孔隙倒塌形成的 角落和裂缝中油膜上的负毛细管压力与最大毛细管压力是直接成比例的。因此，Swi值在初次排水时达到了 (Blunt, 1997, ?ren et al., 1998)。Swi值更高(即较小的最 大毛细管压力)，油膜更薄，传导性更差。这依次减少了油饱和度低时油的相对 渗透率，大部分的油被迫通过膜流动。图7将模拟的水、油的相对渗透率和测量地稳定态地相对渗透率进行

25、了比较。 稳定态实验是在两个合成孔中做的，每个合成孔由四个合在一起的岩心柱塞组 成。即使模拟结果走向略高估了水的相对渗透率，图7表明几乎所有的测量数据都落在模拟数据的分布曲线图中。测量的相对渗透率和预测值之间的这种一致性 是非常可喜的，尤其考虑到砂岩岩相复杂的特性和湿度的不确定性。虽然我们需 要对更大批的非均质岩石样品进行调查，但是这些结果表明计算机生成岩石结合 了数值计算，可以成为一种可提前获得重要储集岩属性的可行的方法。图7.预测的和稳态确定的注水相对渗透率二者比较。实验的初始水饱和度及模拟值的Swi也0.19是相似的结论对于北海储集岩中非均质的河床砂岩，我们研究了其材料直接孔径模拟的预

26、测的可能性，并传输了相关有效属性。研究表明，微 CT图像和基于地质而重建 的样品二者的有效属性具有非常好的一致性。对于计算机生成岩石来说，计算的 注水相对渗透率非常符合测量的数据。我们的模拟表明，对于混合的湿样品来说， 油的相对渗透率和剩余油饱和度非常依赖于位移开始时初始水饱和度。这是由于存在于混合湿孔中的油膜的稳定性直接依赖于最大毛细管压力，从而Swi值在初次排水时达到了我们的结果证实了结合计算机生成岩石和数值计算来生成岩石的潜能和可 行性，还证明了储集岩的流动性，来增加实验获得SCAL 数据的分析和解释。对于比实验方法更容易数字构建的属性和 crosSS性来说，这可以导致典型曲线的发 展。

27、这呈现出填补当前存在于精细地质模型和缺乏相关储集岩属性二者之间技术 缺口可喜的可能性。注：此文是为 2006 年 9 月 12-16 日在挪威特隆赫姆举行的国际岩心分析员协会会议而准备的报告。参考文献1. Adler , P.M., Jacquin, C.G., and Quiblier, J.A.,“Flow in simulated pIonrto. uJ.sMmueltdipiahase Flo”w ,(1990), 16, 691-712.2. Arns, C.H., Knackstedt, M.A., Pinczewski, V., and Garboczi, E.J.,line“a

28、rCeolamsptiuctaptriopneortfiesfrom microtomographic images: Methodlogy and agreement between theory and experiments Geo”ph,ysics, (2002), 67, 1396-1405.3. Arns, C.H., Knackstedt, M.A., Pinczewski, V., and Martys, N.S.,permeame“tryVoirntumalicrotomographicimages ”J,. Petroleum Science and Engineering

29、 , (2004), 45, 41-46.4. Bakke, S., and ?ren, P.E., -D Pore“-S3cale Modeling of Sandstones and flow Simulations in the Pore Networks ”S,PE Journal, (1997), 2, 136-149.5. Biswal, B., Manswart, C., Hilfer, R., Bakke, S., and ?ren, P.E.,analysis“oQf euxapnetirtiamtievnetal andsynthetic microstructures o

30、f sedimentary rocks Physica A,”(1,999), 273, 452-475.6. Blunt, M.,“Pore level modeling of the effects of wettabiliStyPE Journal , (1997”), 2, 494-510.7. Bryant, S.L., Cade, C.A., and Melor, D.W.,ility”prePdeircmtioenabfrom geological models ”A, APG Bulletin ,(1993), 77, 1338-1350.8. Hazlett, R.D.,“

31、Statistical characterization and stochastic modeling of pore netiwn orerklastion to fluid flowMath. Geol ., (1997), 29, 801-822.9. Hilfer, R.,“ Geoamndetdriicelectric characterization of porous mediaPhysical Re”vie, w B, (1991), 60-75.10. Jadhunandan, P.P. and Morrow, N.R.,“ Effect of wettability on

32、 waterflood rfeocrocvruedrye-oil/brine/rocksystems ”S,PE Reservoir Engineering , (1995), 10, 40-46.11. Jin, G., Patzek, T.W., and Silin, D.B.,“ Direct Prediction of the AbsoluteoPfeUrmnceoanbsiloitlyidated andConsolidated Reservoir Rock SPE”p,aper 90084, (2004), presented at the SPE ATCE, Houston, T

33、X, September 2004.12. Jin, G., Patzek, T.W., and Sili n, D.B “ Physics based Reconstruction of SedimentaRroycks”S, PE paper 85387, (2003), Long Beach CA., 19-24 May, 2003.13. Knackstedt, M.A., Arns, C.H., Limaye, A., Sakellariou, A, Senden, t.J., Sheppard, A.P., Sok, R.M., andPinczewski, W.V.,” Digi

34、tal CLoarbeoratory: Properties of reservoir core derived from 3D images SPE pape,r87009, (2204), Kuala Lumpur, Malyasia, 2004.14. Lerdahl, T.R., ?ren,P.E., and Bakke, S.,”A predictive network m-opdhealsfeorflothwreine porous mediaSPE paper 59311, (2000), SPE/DOE Symposium on EOR, Tulsa, OK, April 3-

35、5, 2000.15. Manswart, C., Aaltosalmi, U., Koponen, A., Hilfer, R., and Timonen, J.,-Boltzmann an“d Lfinaitteice difference simulations for the permeability for threedimensional porous media ”Ph,ysical Review E, (2002), 66, 016702.16. Patzek, T.W.,“ Verification of a complete network simulator for drainage iamnbdibitionS”P,E Journal,(2001), 6, 144-156.17. Patzek, T.W., and Kristensen, J.D.,“ Shape factor and hydraulic conductance incanpoinllacirriceusl:aIrI.Two- phase creeping flow J. ”Co,lloid & Interface Sci ., (2001), 236, 305-317.18. Piri, M., and Blunt, M.,-dim“eTnhsiroeneal mixed-wet r