张鹏伟等：基于微观机理的页岩气运移分析 ·139· a 影响区域进行了合理概化，认为0.5m范围内岩体 已经极度破碎，此区域可以作为裂隙处理，即裂隙的 竖直井 尺寸（影响区域）为0.5m×5m×5m.图2(c)中有 主裂隙 限元网格采用四面体单元，水平井附近单元进行加 次级裂隙王 密处理，确保计算收敛性 水平井 裂隙间距 孔眼 3数值算例 本文数值模型是基于COMSOL multiphysics有 (b) (c) 限元多场耦合软件进行二次开发实现.数值算例首 15m 15m 先对模型的有效性与解析解进行对比验证，其次对 影响气体运移的各物理机制进行深入分析. 10m 5m 3.1数值模型验证 5 m 为了验证本文数值模型模拟气体在低渗透多孔 15m 介质中流动的有效性，将推导的多物理场耦合模型 影响区城 简化为只考虑Klinkenberg滑移效应的一维稳态控 0.5m 5 m 制方程，将此控制方程数值解与解析解进行对比验 10m 证.式(13)只考虑滑移项的一维稳态连续方程为： 水平井 50 品(4“))=0 (16) 5m 式(16)与Wu等2)考虑气体Klinkenberg滑移一维 图2。页岩气开采三维物理模型与有限元网格.()页岩气开采 示意图：(b)代表性单元体物理模型：（©）有限元网格 稳态流控制方程一致，其一维模型的边界条件为：一 Fig.2 Three dimensional physical model and finite element mesh of 端以恒定速率(v)注入气体，另一端为恒定边界压 shale gas exploitation model:(a)sketch map of shale gas exploita- 力(P). tion;(b)physical model of representative element volume (REV); 为了验证该模型在低渗透储层中气体流动的有 (c）finite element mesh 效性，渗透率取值参考页岩渗透率范围，从5×10-20 过程中会在垂直主裂隙方向产生不规则的次级裂 m到5×l0-8m2,孔隙度取值参考Marcellus页岩 隙，增加裂隙附近岩体的破碎程度，模型对裂隙横向 数据.具体计算参数见表1. 表1一维模型计算参数 Table 1 Parameters in one dimensional model 气体注人速率， 边界气压， 气体黏滞系数， 渗透率， 孔隙度， linkenberg滑移 压缩因子， a/(kg's-1) P/Pa /(Pa-s) k/m2 p/% 因子，b/Pa B/(kg-Pa-1.m-3) 1×10-6 1×105 1.84×10-5 5×10-18~5×10-20 8.67 7.6×105 1.18×10-5 图3为数值计算结果与解析解的对比，结果显 25 k=5x100m2 一解析解（文献26） 示出本数值模型与解析解具有很好的一致性，验证 ◆数值解 了本数值模型的合理性 20 k=1×10-19m2 3.2三维模型边界条件及计算参数 ◆ 15 三维模型的边界条件有两类 (I)Dirichlet边界条件： k=5x1019m2 p(x,y,z,t)=p(x,y,z,t),t>0,(x,y,z)EI k=1x1018m+ (17) k-5×108m2 ◆ (2)Newman边界条件： ◆ 4 6 10 n(货-0-p）-,.,1>0. 距离m 图3考虑滑移效应的一维稳态气体压力分布 (x,y,z)∈T (18) Fig.3 ID static gas pressure distribution considering gas slippage ef- 式中，t为时间，n为x,y,z方向矢量记号，k为渗透 ffect张鹏伟等: 基于微观机理的页岩气运移分析 图 2 页岩气开采三维物理模型与有限元网格. (a) 页岩气开采 示意图; (b) 代表性单元体物理模型; (c) 有限元网格 Fig. 2 Three dimensional physical model and finite element mesh of shale gas exploitation model: ( a) sketch map of shale gas exploita鄄 tion; (b) physical model of representative element volume (REV); (c) finite element mesh 过程中会在垂直主裂隙方向产生不规则的次级裂 隙,增加裂隙附近岩体的破碎程度,模型对裂隙横向 影响区域进行了合理概化,认为 0郾 5 m 范围内岩体 已经极度破碎,此区域可以作为裂隙处理,即裂隙的 尺寸(影响区域)为 0郾 5 m 伊 5 m 伊 5 m. 图 2(c) 中有 限元网格采用四面体单元,水平井附近单元进行加 密处理,确保计算收敛性. 3 数值算例 本文数值模型是基于 COMSOL multiphysics 有 限元多场耦合软件进行二次开发实现. 数值算例首 先对模型的有效性与解析解进行对比验证,其次对 影响气体运移的各物理机制进行深入分析. 3郾 1 数值模型验证 为了验证本文数值模型模拟气体在低渗透多孔 介质中流动的有效性,将推导的多物理场耦合模型 简化为只考虑 Klinkenberg 滑移效应的一维稳态控 制方程,将此控制方程数值解与解析解进行对比验 证. 式(13)只考虑滑移项的一维稳态连续方程为: 鄣 鄣 ( ( x F·k pM RT 滋 ) 鄣p 鄣 ) x = 0 (16) 式(16)与 Wu 等[26] 考虑气体 Klinkenberg 滑移一维 稳态流控制方程一致,其一维模型的边界条件为:一 端以恒定速率(vm )注入气体,另一端为恒定边界压 力(pb ). 为了验证该模型在低渗透储层中气体流动的有 效性,渗透率取值参考页岩渗透率范围,从 5 伊 10 - 20 m 2到 5 伊 10 - 18 m 2 ,孔隙度取值参考 Marcellus 页岩 数据. 具体计算参数见表 1. 表 1 一维模型计算参数 Table 1 Parameters in one dimensional model 气体注入速率, vm / (kg·s - 1 ) 边界气压, pb / Pa 气体黏滞系数, 滋 / (Pa·s) 渗透率, k / m 2 孔隙度, 渍/ % Klinkenberg 滑移 因子,bk / Pa 压缩因子, 茁 / (kg·Pa - 1·m - 3 ) 1 伊 10 - 6 1 伊 10 5 1郾 84 伊 10 - 5 5 伊 10 - 18 ~ 5 伊 10 - 20 8郾 67 7郾 6 伊 10 5 1郾 18 伊 10 - 5 图 3 为数值计算结果与解析解的对比,结果显 示出本数值模型与解析解具有很好的一致性,验证 了本数值模型的合理性. 3郾 2 三维模型边界条件及计算参数 三维模型的边界条件有两类. (1) Dirichlet 边界条件: p(x,y,z,t) = p(x,y,z,t), t > 0, (x,y,z)沂祝1 (17) (2) Newman 边界条件: n·( k 滋 p - 兹D· 驻 p ) = q(x,y,z,t), t > 0, (x,y,z)沂祝2 (18) 式中,t 为时间,n 为 x, y, z 方向矢量记号,k 为渗透 图 3 考虑滑移效应的一维稳态气体压力分布 Fig. 3 1D static gas pressure distribution considering gas slippage ef鄄 ffect ·139·