考虑气体滑脱及Knudsen扩散的低渗致密砂岩微观流动模拟

曹廷宽,刘成川,卜 淘,段永刚,方全堂，孙连浦

[1.中国石化 西南油气分公司 博士后科研工作站，四川 成都 610041; 2.中国石化 西南油气分公司 勘探开发研究院，四川 成都 610041; 3.西南石油大学 石油与天然气工程学院，四川 成都 610500; 4.华夏吉泰(北京)科技有限公司，北京100020]

低渗致密砂岩气藏分布广泛、资源量丰富，是天然气勘探开发的重要领域[1-2]。与常规气藏相比，该类气藏在储层特征、孔隙结构及产出规律等方面都均有其特殊性[3-4]。研究表明，低渗致密砂岩矿物多样、孔隙结构复杂，微观非均质性极强[5-8]。受孔隙结构控制，天然气的流动具有很强的非线性，滑脱、扩散等非线性输运机制对渗流影响明显，难以利用Darcy定律等现有理论及实验方法加以描述[9-11]。鉴于低渗致密砂岩孔隙结构的复杂性及运移机制的多样性，结合数字岩心技术与微观流动模拟方法来综合评价天然气在此类岩石中的渗流机理。以川西地区沙溪庙组储层岩石为例，利用CT扫描获取不同物性岩样的灰度图像，通过图像处理、孔喉识别、参数统计建立与实际岩样孔隙空间高度一致的真实孔隙模型；基于真实孔隙模型，调整喉道大小，形成包含特定喉道分布的三维孔隙模型。以孔隙模型包含的孔喉单元为对象，建立耦合气体滑脱、Knudsen扩散及Darcy渗流流动控制方程，并通过对流动方程进行离散求解，形成微观流动模拟方法。在此基础上，开展模拟研究，全面分析不同孔隙压力下各输运机制对天然气运移的影响，进而评价该类气藏开采过程中地层压力降低对渗流规律的影响。本文研究成果有助于深化对低渗致密砂岩气藏渗流规律认识，从而为该类气藏的高效开发提供理论依据。

1 三维孔隙模型建模

1.1 实验原理

CT扫描(X-Ray Computerized Tomography，简称CT)是计算机和X-射线技术相结合的产物，该仪器通过让X-射线从各个方向穿过目标样品，利用探测器检测沿不同路径穿透被测样品的X-射线穿透量，以获得样品各个断面的投影数据；再利用相应重构算法对投影数据体进行处理，建立样品内部断层图像[12]。由于CT扫描可准确获取岩石孔隙空间三维分布信息，已经成为研究岩石微观结构及渗流机理的重要技术手段[13-14]。本文采用XCT-400扫描系统来对目标岩样进行扫描，该实验仪器结构示意图如图1所示。

图1 XCT-400扫描系统结构示意图Fig.1 Diagram showing the structure of XCT-400 scanning system

1.2 CT扫描结果

利用上述仪器，对取自川西地区沙溪庙组储层的3块砂岩样品(编号S1-S3)进行实验分析。岩样S1,S2和S3的实测孔隙度分别为13.78%,12.54%和7.55%，相应的气测渗透率为2.34×10-3,1.84×10-3和0.07×10-3μm2，属于低渗致密砂岩范畴。根据所选岩样的物理性质，结合实验仪器的性能，将岩样处理为直径5 mm、高15 mm的规则圆柱体进行CT扫描。通过对岩样进行扫描，获得各岩样的三维灰度图像，其中单张切片图像分辨率为2.15 μm，大小均为980×1005像素(图2)。

从岩样S1,S2和S3的三维图像中分别截取500×500×500像素的体积进行研究。通过图像滤波过滤噪点，并利用图像分割方法分割孔隙相和岩石固，建立了体积为1.08 mm×1.08 mm×1.08 mm的岩样三维数字岩心。

1.3 三维孔隙模型构建

利用CT扫描建立岩样三维数字岩心，便已获得了岩石的孔隙分布及相关结构信息。但直接利用数字岩心进行模拟计算效率很低，因此需对数字岩心中的真实孔隙空间进行抽提简化，建立结构更为简单、便于后续流动模拟的三维孔隙模型[15]。针对川西地区的低渗砂岩，笔者曾应用“最大球”法建立了可靠的三维孔隙模型[16]，建模主要步骤如下：①在数字岩心中搜索孔隙边界，建立由孔隙像素点组成的“最大球”；②剔除完全被包含的“最大球”，并根据球间的连接关系建立由“最大球”构成的孔喉链；③确定孔隙-喉道边界，将孔喉链分割为“孔隙节点”和“喉道单元”，计算孔喉空间坐标、半径、连通性等几何特征参数，最终形成三维孔隙模型。在这类孔隙模型中，孔隙空间由“孔隙节点”表示，代表流体储集的主要场所；“喉道单元”则表征连接孔隙间的狭小部分，是流体在孔隙间运移的通道。

图2 低渗致密砂岩CT切片图像Fig.2 CT slices of low-permeability tight sandstonea.S1岩样CT图像;b.S2岩样CT图像;c.S3岩样CT图像

孔隙模型编号L1L2L3孔隙数83311088911366喉道数131951593911550非连通孔喉比/%6．338．4129．71孔隙半径/μm分布0．89～27．880．89～25．840．89～25．94平均值5．655．104．13喉道半径/μm分布0．21～22．140．22～17．450．21～18．58平均值2．622．391．97孔喉比分布0．11～435．460．12～340．370．15～349．16平均值11．5310．739．46配位数分布0～320～240～14平均值3．122．891．98孔隙度/%实测13．7812．547．55计算13．5011．580．67误差/%2．037．669．11渗透率/(10-3μm2)实测2．341．840．07计算2．151．62—误差/%8．1111．96—

根据上述流程，建立了3块岩样的真实孔隙模型(编号L1—L3)，并预测了模型孔隙度和渗透率，结果如表1所示。从结果来看，模型L1和L2预测的孔隙度、渗透率与岩样实测结果较为吻合，说明模型可靠性高。而更为致密砂岩的孔隙模型L3，受CT图像分辨率的限制，仅能获取孔隙分布，难以识别连通孔喉信息，造成孔隙不连通，仅能预测其孔隙度。为提高致密砂岩的建模精度，可通过提高图像分辨率来实现，例如利用Nano-CT和SEM来建立岩样数字图像。但图像分辨率与观测视场成反比，分辨率越高，所获得的图像尺寸越小，因此，具有纳米尺度分辨率的图像往往难以满足真实孔隙模型构建的需要[17-18]。此外，若采用更高分辨率的岩样数字图像进行三维孔隙建模，将大幅增加运算量及模型构建时间，并增大模拟难度。

大量研究表明，低渗致密砂岩中微米尺度孔隙是流体储集的主要场所，而流体运移主要受狭小喉道(纳米-微米尺度)的控制[19-20]。为定量分析天然气在致密砂岩孔隙空间内的运移规律，以低渗砂岩三维真实孔隙模型L2为基础，保持孔隙分布及拓扑结构不变，等比例调模型喉道大小，建立平均喉道半径分别为1.50,1.00,0.50,0.10 μm的孔隙模型(编号D1—D4)。不同孔隙模型中喉道分布如图3。由于衍生的孔隙模型D1—D4与L2的孔隙空间具有相同的拓扑结构，仅以L1和L2为例进行三维显示(图4a，b)。三维视图中，小球为“孔隙节点”，与数字岩心中孔隙具有相同的空间分布，圆柱则表示“喉道单元”。

图3 三维孔隙模型中喉道分布曲线Fig.3 Throat distribution curves of the 3Dpore network models

图4 孔隙模型三维视图Fig.4 3D view of the pore network modelsa.孔隙模型L1三维视图;b.孔隙模型L2三维视图

所建立三维孔隙模型中，模型L3内孔喉单元连通性差，无法开展后续模拟计算，而L1与L2的孔隙结构及孔喉分布较为接近，因此仅利用L2以及基于L2所构建的孔隙模型D1—D4进行后续微观流动模拟。用于开展模拟的5个孔隙模型具有相同的孔隙拓扑结构，可排除孔隙分布及连通关系等差异对模拟结果的影响。

2 流动控制方程

低渗致密砂岩具有复杂的孔隙结构和狭小的孔喉，流体在该类岩石孔喉中运移不再满足Darcy定律，而是受渗流、滑脱及扩散等输运机制综合作用[11]。以单组份甲烷为模拟流体，考虑甲烷在纳米-微米尺度孔喉的多种运移机制，建立耦合Darcy渗流、气体滑脱及Knudsen扩散的流动控制方程，用以全面描述天然气在纳米-微米尺度孔喉中的流动。

2.1 Darcy渗流

多孔介质中，连续流体在压力差驱动下发生的流动称为Darcy渗流(也称粘性流)，是油气产出运移的主要机制。根据Hagen-Poiseuille定律[21]，当仅考虑粘性力作用时，甲烷流经孔喉单元的质量通量与流动压差间存在如下关系：

(1)

式中：Jadv为甲烷流经孔喉的Darcy渗流质量通量，kg/(m2·s)；r为孔喉半径，m；M为气体摩尔质量，kg/mol；p为孔喉平均压力，Pa；dp/dl为孔喉两端的压力梯度，Pa/m；μ为气体粘度，Pa·s；Z为气体偏差因子，无因次；R为理想气体常数，8.314472 m3·Pa/(K·mol)；T为模拟温度，K。

2.2 气体滑脱

狭窄孔隙空间中，壁面处的气体分子可能处于运动状态，并且由于动量交换，相邻层的气体分子可连同壁面处的气体分子一起沿管壁流动，即气体的流动存在滑脱效应[22]。当气体流动不再满足非滑脱边界条件时，需要利用滑脱修正系数对公式(1)进行校正：

Jslip=FJadv

(2)

气体滑脱效应校正系数F受孔隙壁面及流体物性的综合影响，Brown等人[23]建立了该系数与孔喉尺寸、壁面物理性质及流体性质间的相关性：

(3)

式中：Jslip为甲烷流经孔喉的滑脱流动质量通量，kg/(m2·s)；F为气体滑脱效应校正系数，无因次；α为孔隙壁面处的切向动量供给系数，无量纲，α与孔隙壁面形态、气体性质及温度和压力有关，通常取(0，1][22，24]。

2.3 Knudsen扩散

对于单相气体，若气体分子的平均自由程大于孔喉的平均半径，则不能忽略分子与壁面碰撞作用对气体运移的影响[25]。这种情况下气体流经孔喉时需要考虑Knudsen扩散的影响，其质量通量描述如下式：

(4)

式中：JD为甲烷流经孔喉的Knudsen扩散质量通量，kg/(m2·s)；ρavg为孔喉内甲烷平均密度，kg/m3；DK为Knudsen扩散系数，m2/s；cg为气体的压缩系数，Pa-1。

Javadpour[11]建立了Knudsen扩散系数与孔喉几何尺寸、气体物性间的关系如下：

(5)

2.4 多机制流动控制方程

在纳米-微米尺度孔喉中，甲烷的流动总通量可以表示为Darcy渗流通量、滑脱流动通量和Knudsen扩散通量的线性叠加[11,26]：

Jt=Jadv+Jslip+JD

(6)

综合公式(1)，(2)，(4)和(6)，可得到耦合Darcy渗流、气体滑脱及Knudsen扩散影响的流动方程(7)。

(7)

式中：Jt为甲烷流经孔喉的总质量通量，kg/(m2·s)。

由公式(7)可知，当考虑不同输运机制及流体物性变化时，天然气在储层孔喉内的运移将受到温度压力、流体物性及孔隙尺寸等因素的综合影响。下面结合三维孔隙模型，进一步探讨影响天然气流动运移的主要因素及变化规律。

3 模拟结果

3.1 模拟流程

以三维孔隙模型为流动场所，利用流动控制方程(7)来描述流体在孔喉单元中的流动运移，通过模拟计算可综合评价温度压力、流体物性、孔喉参数及输运机制等对流体运移的影响。考虑到气藏开发过程中，地层压力下变化幅度大，而地层温度大致恒定，因此主要通过设定不同的孔隙压力来模拟不同开发阶段中天然气的产出运移机制。具体模拟流程：①在三维孔隙模型中饱和模拟流体，设定孔隙压力及流动压差，计算流体高压物性参数；②基于孔喉单元，根据公式(7)建立以孔隙节点压力为变量的流动方程组；③优选数值计算方法迭代求解方程组，得到孔隙压力分布；④根据孔隙压力，计算设定压差下流经孔隙模型的流体流量或通量；⑤改变模型参数或模拟条件，进行敏感性分析。

3.2 结果分析

3.2.1 孔隙压力对总质量通量的影响

以建立的低渗致密砂岩三维孔隙模型为平台，设定模拟温度为80 ℃，模型两端压力梯度恒定为0.5 MPa/m，预测不同孔隙压力下甲烷流经模型的总通量，结果见图5。模拟结果表明，甲烷的总质量通量与孔隙模型的平均喉道半径(图中用ravg表示)成正比，喉道越大，气体流动能力越强，流动通量也越大。而对于不同孔隙模型，孔隙压力变化对气体总通量的影响规律相同：当孔隙压力高于10 MPa时，压力降低对气体总通量影响较小，甲烷的流动基本呈线性流动；而在孔隙压力低于10 MPa范围，甲烷的总质量通量随孔隙压力的降低而减小，即该孔隙压力范围内甲烷的流动具有较强非线性。下面进一步分析不同孔隙模型中孔隙压力与各气体输运机制的关系。

3.2.2 Darcy渗流对流动的影响

计算不同孔隙压力下Darcy渗流通量占总通量的百分比，结果如图6。相同压差下，不同孔隙模型中Darcy渗流通量与平均喉道半径成正比，喉道越大，Darcy渗流对天然气运移的贡献越大。孔隙压力也与Darcy渗流通量成正相关性，孔隙压力越高，Darcy渗流通量越大。对于平均喉道半径大于0.1 μm的岩样孔隙模型(L2,D1—D3)，当孔隙压力高于10 MPa时，Darcy渗流通量占总通量的90%以上，说明Darcy渗流是主要的输运机制。而天然气流经更为致密的岩石时，Darcy渗流对流动的贡献随着平均喉道的缩小、孔隙压力的降低而减小，在这类储层岩石中应当考虑非Darcy运移机制的影响。

图5 不同条件下天然气流动总质量通量Fig.5 Total mass flux of gas flow at different conditions

图6 Darcy渗流对天然气流动的影响Fig.6 Influences of Darcy flow on gas flow

3.2.3 气体滑脱对流动的影响

不同压力下气体滑脱效应对流动的影响如图7所示。从模拟结果来看，孔隙压力越低，喉道越小，气体滑脱效应越显著，滑脱流动通量占总通量比例相应更大。例如，孔隙压力低至1 MPa时，模型D4中气体滑脱对流动的贡献为18.42%；而相同压力下，模型D3中气体滑脱通量只占总通量的5.56%。综合来看，当平均喉道半径小于0.1 μm且孔隙压力低于10 MPa时，气体滑脱是一种重要的输运机制，对流动的贡献不可忽略。

3.2.4 Knudsen扩散对流动的影响

同样地，根据不同模型预测Knudsen扩散通量/总通量与孔隙压力的关系，结果如图8。对比图7和8，喉道分布、孔隙压力对Knudsen扩散和气体滑脱的影响规律类似，也是主要影响低压(孔隙压力低于10 MPa)、致密砂岩(平均喉道半径小于0.1 μm)中气体的流体运移，但在相同孔隙压力下，Knudsen扩散对流动的贡献要大于气体滑脱效应。

图7 气体滑脱效应对天然气流动的影响Fig.7 Influences of gas slippage on gas flow

图8 Knudsen扩散对天然气流动的影响Fig.8 Influences of Knudsen diffusion on gas flow

4 结论

1) 针对川西低渗致密砂岩的孔隙结构及渗流特征，基于数字岩心理论构建了包含不同尺度孔喉的三维孔隙模型，通过在模型中耦合Darcy渗流、气体滑脱和Knudsen扩散3种气体输运机制，形成了低渗致密砂岩复杂渗流规律微观流动模拟方法。

2) 低渗致密砂岩中喉道分布及孔隙压力是影响天然气传输运移的主要因素，孔喉越小、压力越低，天然气流动的非线性越显著。在平均喉道半径大于0.1 μm的砂岩中，压力高于10 MPa，天然气基本满足Darcy渗流；而岩性更为致密(平均喉道半径小于0.1 μm)时，当孔隙压力低于10 MPa，气体滑脱及Knudsen扩散均是天然气运移的重要机制。

[1] 贾承造,郑民,张永峰.中国非常规油气资源与勘探开发前景[J].石油勘探与开发,2012,39(2):129-136.

Jia Chengzao,Zheng Min,Zhang Yongfeng.Unconventional hydrocarbon resources in China and the prospect of exploration and development[J].Petroleum Exploration and Development,2012,39(2):129-136.

[2] 戴金星,倪云燕,吴小奇.中国致密砂岩气及在勘探开发上的重要意义[J].石油勘探与开发,2012,39(3):257-264.

Dai Jinxing,Ni Yunyan,Wu Xiaoqi.Tight gas in China and its significance in exploration and exploitation[J].Petroleum Exploration and Development,2012,39(3):257-264.

[3] 陈恭洋,晋静,罗迎春,等.沉积相控制下致密砂岩气藏分布特征——以四川盆地川西坳陷中段蓬莱镇组为例[J].石油与天然气地质,2017,38(8):467-477.

Luo Gongyang,Jin Jing,Luo Yingchun,et al.Characteristics of tight sandstone gas reservoir distribution under control of sedimentary facies:A case study from the Penglaizhen Formation in the central western Sichuan depression[J].Oil & Gas Geology,2017,38(8):467-477.

[4] 李杪,侯云东,罗静兰,等.致密砂岩储层埋藏-成岩-油气充注演化过程与孔隙演化定量分析——以鄂尔多斯盆地东部上古生界盒8段天然气储层为例[J].石油与天然气地质,2016,37(6):882-892.

Li Miao,Hou Yundong,Luo Jinglan,et al.Burial,diagenesis,hydrocarbon charging evolution process and quantitative analysis of porosity evolution:a case study from He 8 tight sand gas reservoir of the Upper Paleozoic in eastern Ordos Basin[J].Oil & Gas Geology,2016,37(6):882-892.

[5] 任晓霞,李爱芬,王永政,等.致密砂岩储层孔隙结构及其对渗流的影响——以鄂尔多斯盆地马岭油田长8储层为例[J].石油与天然气地质,2015,36(5):774-779.

Ren Xiaoxia,Li Aifen,Wang Yongzheng,et al.Pore structure of tight sand reservoir and its influence on percolation——taking the Chang 8 reservoir in Maling Oilfield in Ordos Basin as an example[J].Oil & Gas Geology,2015,36(5):774-779.

[6] 刘广峰,白耀星,王文举,等.致密砂岩储层微观孔喉结构及其对渗流特征的影响——以鄂尔多斯盆地周长地区长8储层为例[J].科学技术与工程,2017,17(5):29-34.

Liu Guangfeng,Bai Yaoxing,Wang Wenju,et al.Pore structure of tight sand reservoir and its influence on percolation:a case of Chang 8 Formation in Zhouchang area,Ordos Basin[J] Science Technology and Engineering,2017,17(5):29-34.

[7] Nelson P H.Pore-throat sizes in sandstones,tight sandstones,and shales[J].AAPG Bulletin,2009,93(3):329-340.

[8] 李易霖,张云峰,尹淑丽,等.致密砂岩储集空间多尺度表征——以松辽盆地齐家地区高台子油层为例[J].石油与天然气地质,2016,37(6):915-922.

Li Yilin,Zhang Yunfeng,Yin Shuli,et al.Characterization of the pore space in tight sandstone reservoirs from macroscopic and microscopic perspectives:a case study of Gaotaizi reservoir in Qijia area of the Songliao Basin[J].Oil & Gas Geology,2016,37(6):915-922.

[9] Hadjiconstantinou N G.The limits of Navier-Stokes theory and kinetic extensions for describing small-scale gaseous hydrodynamics[J].Physics of Fluids,2006,18(11):111301.

[10] 杨建,康毅力,李前贵,等.致密砂岩气藏微观结构及渗流特征[J].力学进展,2008,38(2):229-236.

Yang Jian,Kang Yili,Li Qiangui,et al.Characters of micro-structure and percolation in tight sandstone gas reservoirs[J].Advances in Mechanics,2008,38(2):229-236.

[11] Javadpour F.Nanopores and apparent permeability of gas flow in mud rock (shale and siltstone)[J].Journal of Canadian Petroleum Technology,2009,48(8):16-21.

[12] 苏娜,段永刚,于春生.微CT扫描重建低渗气藏微观孔隙结构——以新场气田上沙溪庙组储层为例[J].石油与天然气地质,2011,32(5):792-796.

Su Na,Duan Yonggang,Yu Chunsheng.Reconstruction of microscopic pore structure in low permeability gas reservoirs by micro-CT scanning:an example from upper Shaximiao Formation in Xinchang Field[J].Oil & Gas Geology,2011,32(5):792-796.

[13] Blunt M J,Bijeljic B,Dong H,et al.Pore-scale imaging and modeling[J].Advances in Water Resources,2013,51(1):197-216.

[14] 白斌,朱如凯,吴松涛,等.利用多尺度CT成像表征致密砂岩微观孔喉结构[J].石油勘探与开发,2013,40(3):329-333.

Bai Bin,Zhu Rukai,Wu Songtao,et al.Multi-scale method of nano(micro)-CT study on microscopic pore structure of tight sandstone of Yanchang Formation,Ordos Basin[J].Petroleum Exploration and Development,2013,40(3):329-333.

[15] Dong H,Blunt M J.Pore-network extraction from micro-computerized-tomography images[J].Physical Review E,2009,80(3):036307.

[16] Cao T K,Duan Y G,Yu B M,et al.Pore-scale simulation of gas-water flow in low permeability gas reservoirs[J].Journal of Central South University,2014,21(7):2793-2800.

[17] 屈乐,孙卫,杜环虹,等.基于CT扫描的三维数字岩心孔隙结构表征方法及应用——以莫北油田116井区三工河组为例[J].现代地质,2014,28(1):190-196.

Qu Le,Sun Wei,Du Huanhong,et al.Characterization technique of pore structure by 3D digital core based on CT scanning and its application:an axample from Sangonghe Formation of 116 Well Field in Mobei Oilfield[J].Geoscience,2014,28(1):190-196.

[18] 邹才能,朱如凯,白斌.中国油气储层中纳米孔首次发现及其科学价值[J].岩石学报,2011,27(6):1857-1864.

Zou Caineng,Zhu Rukai,Bai Bing,et al.First discovery of nano-pore throat in oil and gas reservoir in China and its scientific value[J].Acta Petrologica Sinica,2011,27(6):1857-1864.

[19] Mousavi M A,Bryant S L.Connectivity of pore space as a control on two-phase flow properties of tight-gas sandstones[J].Transport in Porous Media,2012,94(2):537-554.

[20] Gao H,Li H A.Pore structure characterization,permeability evaluation and enhanced gas recovery techniques of tight gas sandstones[J].Journal of Natural Gas Science and Engineering,2016,28:536-547.

[21] Sutera S P,Skalak R.The history of Poiseuille's law[J].Annual Review of Fluid Mechanics,1993,25(1):1-20.

[22] Civan F.Effective correlation of apparent gas permeability in tight porous media[J].Transport in porous media,2010,82(2):375-384.

[23] Brown G P,DiNardo A,Cheng G K,et al.The flow of gases in pipes at low pressures[J].Journal of Applied Physics,1946,17(10):802-813.

[24] MillikanR A.Coefficients of slip in gases and the law of reflection of molecules from the surfaces of solids and liquids[J].Physical Review,1923,21(3):217-38.

[25] Javadpour F,Fisher D,Unsworth M.Nanoscale gas flow in shale gas sediments[J].Journal of Canadian Petroleum Technolog,2007,46(10):55-61.

[26] Xiong X,Devegowda D.A fully-coupled free and adsorptive phase transport model for shale gas reservoirs including non-Darcy flow effects[J].SPE 159758,2012.