CO2在饱和稠油岩心中的扩散及仿真模拟

杜 林,王仕莉,袁京素,马剑坤,孙来喜

(1.油气藏地质及开发工程国家重点实验室(成都理工大学)，成都 610059;2.中国石油新疆油田公司 勘探开发研究院，新疆 克拉玛依 834000;3.中国石油长庆油田分公司 第八采油厂，西安 710018)

常规油气资源的减少使得稠油油藏的开发越来越受到重视。据统计，世界上稠油的储量约为100×109t[1]。随着常规石油资源的日益枯竭和日益增加的油气消耗量，稠油资源的开采比重也会不断增加。稠油具有黏度高、密度大、开采过程中流动阻力大等特点[2]，因此开采过程中体积扫油效率较低。稠油油藏开发的重点在于降低稠油的黏度，增加稠油的流动性。

气体在稠油中的扩散性能对于稠油冷采，尤其是注气提高采收率项目设计具有重要的意义[3-9]。气体与稠油接触的过程中，气体通过扩散进入稠油，进而溶解于稠油中，降低稠油的黏度，使得稠油易于开采。目前关于气体在稠油中的扩散研究多集中于扩散系数的确定[10-12]和数学模型的建立[13-15]等，几乎没有涉及其扩散过程的仿真模拟。仿真模拟可更直观地观察气体在稠油中扩散过程及扩散结果，以及可通过调整参数进而观察扩散过程的变化，因此对于稠油油藏注气提高采收率具有一定的指导意义。

在进行气体在油相中扩散过程的实验研究中，最常见的方法是压降分析方法(PVT)。气体向油相扩散过程中，扩散会引起气体的压力降低。将扩散系数与压力进行耦合，通过监测压力的变化，可计算出扩散系数的值[3-4,6,11,13-14,16-20]。本文采用压降分析法求取CO2在饱和稠油岩心中的扩散系数。但PVT方法易受到流体对流产生的影响，如气体导入速率过快等引起的气体在液相中流动。为了消除对流的影响，最常用的方法是采用极缓慢的注入速率，以最大程度降低对流的影响[3,21]。

本文针对CO2在饱和稠油的一维岩心中的扩散建立物理模型及数学模型，通过求解扩散数学模型，得到扩散过程中的浓度分布函数，将其与气体状态方程耦合，建立CO2在饱和稠油的一维岩心中扩散的压降模型。在此基础上，设计并进行CO2在饱和稠油的一维岩心中的扩散实验，利用构建的数学模型对实验过程中的压力变化图像进行拟合，求得此扩散过程中的扩散系数。利用求得的扩散系数基于有限元的方法进行仿真，基于仿真结果，对CO2在饱和稠油的岩心中扩散规律进行研究，并对影响CO2在稠油中扩散过程的因素进行分析。

1 扩散模型

1.1 物理模型

为了建立气体在饱和稠油的多孔介质中扩散过程的物理模型，考虑在一个PVT室中，上部为气体、下部为饱和稠油的多孔介质，二者相互接触，则气体在浓度差的作用下，会自发地向稠油中扩散。在扩散过程中，气体会溶解于原油之中，进而使得稠油的黏度下降。据此建立的物理模型如图1所示。

图1 物理模型图Fig.1 Diagram demonstrating the physical model of the CO2 diffusion in porous media saturated by heavy oil

由图1可以看出上部的气体由于浓度高于下部多孔介质中气体的浓度，在此浓度差的作用下，气体扩散进入下部多孔介质中的稠油。其中x为液相中的坐标，hl、hg分别为液相和气相的高度。

1.2 数学模型

基于图1的物理模型，建立气体在饱和稠油的多孔介质中扩散的数学模型。为了增加扩散方程的可解性，必要的假设条件如下：①忽略由于气体溶解于稠油而引起的稠油膨胀；②扩散室中保持恒温恒容状态；③扩散系数为定值，不会随着扩散的进行而改变。

此外，在此过程中，在气体和稠油的界面处，即x=hl处，气体的浓度为瞬时平衡浓度，此边界条件是由Y.P.Zhang[3]于2000年首次提出来，并在其后相关的扩散研究中得到广泛的应用。另一方面，在包含多孔介质的稠油底部截面处，即x=0处，没有CO2通过此截面，则CO2的扩散通量为0。在扩散开始之前，稠油中无气体存在。则根据菲克定律，该过程的数学模型为

(1)

其中：c为气体在稠油中的浓度(mol/L);Deff为CO2在多孔介质中的有效扩散系数(m2/s);t为扩散时间(s);x为位置(m)。

利用分离变量法[22]，解得上述扩散方程的解析解为

(2)

式(2)即为CO2在饱和稠油的一维岩心中扩散的数学模型的解，该式可用于求解在不同时刻不同位置处的浓度。

对式(2)在区间0～hl上进行积分，可得到在不同时间稠油中扩散的CO2的量，其表达式如下

(3)

在PVT室中，CO2向稠油中的扩散使得CO2在上层气相中的量减少，进而导致上层气体的压力减小。显然体系的压力与该扩散过程存在一定的关系，该关系可通过气体的状态方程进行表征。根据气体状态方程，气体的量、温度、压力之间的关系

pV=ZnRT

(4)

其中：p为气体的压力(Pa);V为气体的体积(m3);Z为气体的压缩系数，理想气体Z=1;n为气体的量(mol);R为普适气体恒量，为8.31 Pa·m3/(K·mol);T为热力学温度(K)。

假设气体的压缩因子和温度是常数，则设气体在时间t内因扩散物质的量减少量为Δn，则状态方程可化为

ΔpV=ZΔnRT

(5)

变形为

(6)

由质量守恒可知，气相中CO2的减少量等于CO2扩散进入原油中的量，即m=Δn，因此，由式(3)和式(6)可得

(7)

式(7)即为CO2在饱和稠油的一维岩心中扩散过程的压降模型。从该式可以看出，将实验过程中的压力-时间关系图像进行拟合，可求得CO2在饱和稠油的一维岩心中扩散过程中的扩散系数。

2 扩散实验

2.1 实验样品

实验中使用的岩心为人造岩心，直径为2.526 cm，长度为20.612 cm。实验中采用的CO2的纯度(体积分数)可达到99.9%，在50℃时原油的黏度为2 Pa·s，相对密度为0.957。实验过程中的条件及岩心数据如表1所示。

2.2 实验仪器及步骤

基于CO2在饱和稠油的一维岩心中扩散过程中的压降模型进行CO2扩散实验。实验采用的装置如图2所示。

表1 扩散实验条件Table 1 Experimental parameters of the CO2 diffusion in porous media saturated by heavy oil

图2 CO2在饱和稠油的一维岩心中扩散实验装置图Fig.2 Sketch showing experimental equipment of the CO2 diffusion in porous media saturated by heavy oil

实验装置主要包括CO2气瓶、高压活塞泵、中间容器、恒温箱、扩散PVT室和传感器及计算机数据采集系统。其中CO2气瓶、高压活塞泵和中间容器用于将CO2导入扩散室中，恒温箱用于将实验过程中的温度稳定在50℃，PVT室是扩散实验的核心，传感器及计算机数据采集系统用于采集扩散实验过程中的压力数据。

①实验前准备：将岩心洗净、烘干，备好CO2气样，连接实验装置。

②稠油饱和岩心：将洗净烘干后的岩心饱和稠油。

③安放岩心：将饱和稠油的岩心垂直放置在扩散室中，拧紧扩散室的盖子，安装好实验装置，将实验装置上的阀门关闭。

④升温至设置温度：打开恒温箱，给扩散室加温，并在实验温度下放置4～6 h。

⑤将CO2导入扩散室：先打开阀门1和4，将CO2导入中间容器中，关闭阀门。再打开阀门6和4，利用高压活塞泵以极其缓慢的速率(0.1 mL/min)将CO2导入扩散室中直至压力升高至实验压力，关闭阀门4。当CO2与稠油接触时，CO2在浓度差的作用下向稠油中扩散，并开始采集压力数据。

⑥整理仪器，结束实验：待采集完成压力数据后，缓慢打开阀门3，慢慢释放扩散室中的CO2，并取出岩心，洗净烘干，同时整理好仪器，结束实验。

值得注意的是，实验过程中的第⑤步将CO2导入扩散室的过程中，先将CO2导入中间容器中，再利用高压活塞泵将中间容器中的CO2缓慢导入扩散室中使其与饱和稠油的岩心接触，并缓慢注入加压至实验设置压力。此过程的目的是为了防止CO2从高压气瓶直接导入扩散室时，高压引起的高流速产生对流，进而产生多相流，对扩散过程产生影响。

2.3 实验结果及分析

通过压力传感器及计算机采集系统采集到上述实验过程中的压力曲线如图3所示。

图3 实验过程中的压力与时间关系图Fig.3 The pressure verse time curve of the diffusion experiment

从图3可以看出，CO2在饱和稠油岩心的扩散过程中，压力值逐渐减小，证明CO2向岩心中发生了扩散。此外，在扩散的初期，由于上部CO2气区与岩心中CO2浓度差较大，扩散速率较大，表现为曲线的斜率较大；随着扩散的进行，CO2的浓度差逐渐减小，扩散速率逐渐减小，表现为曲线的斜率逐渐减小，其浓度分布也逐渐趋于平缓。利用式(7)对图3进行拟合，求得CO2在饱和稠油岩心中的扩散系数Deff=1.025×10-10m2/s。将此值与相关文献中[23- 24]的扩散系数进行对比，均处于10-10m2/s的数量级，一定程度上证明了模型的正确性。

3 扩散仿真

上述实验所用岩心为砂岩岩心，其孔隙度(q)为26.32%。将该岩心部分制作成薄片，在显微镜下拍摄岩心薄片的照片；利用图像处理软件对薄片中基质边界进行识别，得到灰度图像；再利用拓扑软件对其进行矢量化，得到岩心薄片的矢量化图片。岩石基质边界识别和矢量化图像如图4所示。

图4 薄片基质边界识别及矢量化图像Fig.4 The identification of boundary matrix and vectorized image

将矢量化图像导入仿真软件COMSOL Multiphysics之中，结合其识别的岩石基质，利用全视域对岩石基质作差集，得到岩石薄片的孔喉分布图像(图5)。

图5 岩石薄片孔喉分布图像Fig.5 The image showing distribution of pore and throat in the porous media

图5的大小是640 μm×320 μm。利用图像分析软件对该视域下的孔喉分布进行分析，得到其面孔率为41.48%，平均孔隙半径为26.35 μm，平均喉道半径为17.21 μm。值得一提的是，为了使得仿真结果收敛，孔喉识别过程中对部分边界作了圆滑处理，并且在进行岩石基质边界识别及勾画过程中存在一定的人为操作误差。

得到岩石孔喉分布图像后，对图5进行边界条件和初始条件的设置，设置仿真温度为50℃，仿真初始压力为7.5 MPa，孔喉中稠油的黏度为2.0 Pa·s。设置图5左侧为饱和稠油的岩石与CO2接触的区域，即边界浓度为瞬时平衡浓度ceq，由气体状态方程得出ceq=2.8 mol/m3；设置最右侧为扩散通量为0的边界，即浓度对位置的一阶偏导数为0，并且设置在孔喉中CO2的初始浓度为0。在设置好上述条件后，对仿真岩石薄片进行网格剖分。为提高计算精度，采用极细化三角形网格剖分，剖分网格的结果如图6所示。

图6 对仿真孔喉进行网格剖分Fig.6 The mesh generation of the porous media for the simulation

从图6中可以看出，网格剖分时在孔喉曲率变化较大处采用了网格加密，此加密方法一方面可以提高计算精度，另一方面可以保证计算过程收敛。在进行参数输入和网格剖分后，基于有限元的方法对该模型进行仿真分析。由于岩心薄片的长度仅有640 μm，CO2扩散前缘会在较短的时间内到达最右端，因此仅设置扩散仿真模拟时间为 1 800 s，并且设置时间步长为1 s，得到在扩散初始时刻(0 s)、1 s、500 s、1 000 s、1 500 s时刻和 1 800 s时刻的扩散形态图像。其中，1 s、500 s、1 000 s、1 500 s时刻的扩散图像中添加了扩散矢量方向箭头(图7)。

图中的颜色分布代表CO2在多孔介质中的浓度，且从红到蓝CO2浓度逐渐降低。图7-A为扩散初始时刻(0 s)的扩散形态，岩心孔隙中的CO2浓度为0，全部为稠油。随着扩散的进行，在扩散的结束时刻(图7-F)，CO2在岩心孔隙中的不同位置处均呈有分布。在左端，岩心中稠油直接与CO2相接触，CO2浓度较高；在岩心的右端，CO2浓度较低。从左端到右端，CO2浓度逐渐降低，这是由于CO2在饱和稠油的岩心中扩散过程较慢且稠油黏度较大，对CO2扩散过程的阻碍也较大，因此CO2从最左端扩散到最右端需要一定的时间。此时，CO2在最右端的浓度为0.157 mol/m3，证明在第1 800 s时刻，CO2的扩散前缘已经到达了岩心最右端。

对仿真扩散过程进行分析，可以得到一些CO2在饱和稠油岩心中的扩散规律：①图7-B～图7-E中的扩散矢量方向表明，在岩石基质表面的扩散矢量方向沿基质表面的切线方向，且指向CO2浓度势降低的方向。在第1 s时刻，扩散的矢量方向垂直于CO2气相边界且只在边界处存在，证明扩散刚发生时扩散方向为垂直于边界指向浓度势降低方向，上述现象是由于扩散是在浓度势差的驱动下产生的。②扩散速度的大小受稠油与CO2的接触面积大小影响(图7-C)。在该区域上半部分，即红线以上区域，岩心中稠油与CO2气相区仅通过一个孔隙通道相连接；而下半部分岩心中，即红线以下区域，稠油与CO2气相区通过5个孔隙通道相连接：上半部分的接触面积远小于下半部分，因此，从图中可以看出，下半部分的CO2浓度明显高于上半部分。③扩散速度的大小受孔径大小影响。如图7-D中区域1和区域2所示，区域1的平均CO2浓度为1.32 mol/m3，区域2的平均CO2浓度为1.18 mol/m3，孔径更小的孔隙通道处，CO2浓度更高，扩散的速度更大。这是由于更小的孔隙通道对CO2的质量分散作用越小，CO2扩散的物质的量更大，因此CO2扩散的速度越快，浓度越高。

图7 不同时刻的扩散形态图像和扩散矢量方向箭头Fig.7 The diffusion situation and directional arrow of diffusion at different times

4 结 论

本文针对CO2在饱和稠油的一维岩心中的扩散，建立了该过程中的物理模型和数学模型，将其与气体状态方程进行耦合，建立了压降模型。在此基础上，设计并进行了CO2在饱和稠油的一维岩心中的扩散实验，利用模型对实验结果进行拟合，求得了该过程的扩散系数，并对CO2在饱和稠油的一维岩心中的扩散进行了仿真分析，得到的主要结论如下：

a.CO2在饱和稠油的一维岩心中的扩散系数为1.025×10-10m2/s，其所处的数量级契合现有研究中扩散系数的数量级。

b.通过对扫描电镜下的岩石薄片进行基质边界识别、矢量化后，可以对其中的扩散情形进行模拟仿真，仿真时参数设置和边界条件设置应与实验相契合。

c.诸多因素会影响CO2在饱和稠油的岩心中的扩散特征，如孔隙直径的大小、稠油与CO2的接触面积等。

d.在扩散过程中，岩心基质表面的扩散方向为基质表面的切向方向，且其水平方向的分量指向CO2浓度势降低的方向。