N2/CO2 注入压力对含瓦斯煤岩中甲烷解吸的影响

李树刚 白 杨 林海飞 严 敏 龙 航 郭豆豆

1. 西安科技大学安全科学与工程学院 2.西部矿井瓦斯智能抽采工程研究中心

0 引言

我国煤层气（俗称瓦斯，主要成分为CH4）资源丰富，但是许多矿区煤层渗透率低，使得采用常规钻孔预抽煤层瓦斯的效果不理想。为了提高煤层瓦斯预抽率（钻孔抽放瓦斯量与钻孔控制范围内煤层瓦斯储量的比值）、缩短预抽时间，应对煤层采取强化增透措施。目前，注气（以N2或CO2为主要成分）是促进瓦斯解吸的主要方法之一[1-2]，通过边注气边抽采来提高低渗透煤层的抽采效率。

为了深入揭示气驱促进瓦斯解吸的机理，国内外学者开展了煤层气吸附/解吸的相关研究[3]。目前，针对注气促进瓦斯解吸的大部分研究都通过室内实验及现场试验来进行[4-6]，而煤是由无机物和有机物构成的，内表面结构特征复杂，因而应当从分子水平来解释注气驱替机制[7-8]。采用分子模拟方法，能够从分子水平来研究吸附质与吸附剂之间的微观作用机理[9]。

Matranga 等[10]运用巨正则系综蒙特卡洛（GCMC）方法模拟CH4、CO2、N2等气体在活性炭上的吸附。Mosher 等[11]、李树刚等[12-13]对煤孔隙结构中CH4、CO2的竞争吸附进行了分子模拟研究。相建华等[14]、曹达鹏等[15]运用蒙特卡洛方法模拟多组分气体在煤分子及活性炭孔隙中的吸附。周婧等[16]采用分子动力学模拟方法，研究含水情况下煤层气在煤表面的吸附/解吸性质。Hu 等[17]、林海飞等[18]、You 等[19]通过分子模拟方法，研究了煤储层中CO2、CH4的解吸、扩散。Sun等[20]、王宝俊等[21]通过改进的MS-BPD 孔隙扩散模型，对注气过程中多组分气体吸附、扩散进行研究，分析了注入气体组分对煤层气预抽率的影响。上述研究主要通过室内实验、现场试验及数值模拟手段，研究活性炭材料及煤结构中不同气体的吸附/解吸机理，而实际煤体本身含有吸附态瓦斯，注入气体类型、注入压力对煤层瓦斯的解吸影响很大，因而目前对于注气过程中注入气体分子与含瓦斯煤间的相互作用机制尚不明确。

为此，笔者基于含瓦斯煤吸附构型，采用GCMC 方法，研究不同注入压力下，注入N2/CO2后含瓦斯煤结构中CH4的解吸行为，以期明确注入N2/CO2促进煤体中CH4解吸的作用机理。

1 模拟方案

1.1 煤结构模型构建

1.1.1 煤分子确定

采用Materials Studio 软件构建煤体三维模型，并对各种无定型、高分子及晶体材料的物理化学性质及分子间相互作用进行研究，建立Wiser 煤结构模型（C184H155N3O20S3）[22-23]，该模型考虑了C、O、H、S 及N 原子形态及煤中芳香核大小，能够较全面地反映模拟煤体三维结构特征。该模型中C、O、H、S、N 原子的质量百分数分别为78.3%、11.3%、5.5%、3.4%、1.5%，与实测数据接近[24]。

采用COMPASS 力场[25]来描述吸附剂和吸附质间的相互作用，体系总能量计算式为：

式中Etotal表示体系总能量，kcal/mol（1 cal=4.19 J，下同）；Evalence表示键键相互作用能，kcal/mol；Enon-bond表示非键相互作用能，kcal/mol。

Evalence计算式为：

式中Eb表示键长伸缩振动能，kcal/mol；Eθ表示键角弯曲振动能，kcal/mol；Eφ表示二面角扭转振动能，kcal/mol；Eω表示面外角的偏移振动能，kcal/mol；Ec表示键间相互作用交叉偶合能，kcal/mol。

Enon-bond计算式为：

式中Ev表示范德瓦尔斯能，kcal/mol；Ee表示库伦相互作用能，kcal/mol；Eh表示氢键相互作用能，kcal/mol。

1.1.2 模型结构优化及合理性验证

首先，构建含2 个煤分子的晶胞，然后，采用Forcite 模块的Geometry Optimization 模块进行结构优化，得到能量最低的煤分子模型（初始构型），再根据体系能量变化，运用Metropoils 运算规则来形成新构型。在三维方向上，通过周期性边界条件进行限制，晶胞为棱长为21.7706 Å（1 Å=0.1 nm，下同）的立方体。

在Geometry Optimization 模块中，选择优化方法为Medium 法；力场选择COMPASS 力场；范德华力、氢键作用力的计算采用Atom based 法，库伦力的计算采用Ewald 方法。分子力学优化后得到局部能量最小化的Wiser 结构模型，在此基础上，基于退火动力学模拟，在整个势能面上对Wiser 结构模型进行全局优化，得到一个全局能量最小化构型；然后，将晶胞扩展，构建1×2×2 的超晶胞，并加入5 nm真空层。模型构建过程如图1 所示。当体系能量达到最低时，计算得到模型密度为1.23 g/cm3。该数值在实际烟煤密度（一般介于1.2 ～1.4 g/cm3）范围内[26]，说明优化后的煤分子空间结构合理。

图1 煤结构模型构建过程示意图

1.2 含瓦斯煤吸附构型筛选

为了确定最低能量含瓦斯煤吸附构型，需要计算超晶胞模型吸附特性参数。Sorption 吸附模块中，设置任务项为Adsorption isotherm，注入压力介于0.01 ～1.00 MPa，温度为303 K，模拟加载至平衡态步骤数为1×106步，总步骤数为2×106步；构型计算选择Metropolis 法，电荷计算采用电荷平衡（QEq）法，力场及范德华力、氢键作用力、库伦力计算方法的设置与前述模型结构优化中的设置一致。Sorption 吸附模块通过逸度计算吸附量，基于PR 状态方程进行逸度和压力的转化[27]。通过3 次模拟计算，得到不同吸附构型的吸附特性参数（表1），由于构型A 能量最低，将其作为含瓦斯煤吸附构型。

表1 不同吸附构型的吸附特性参数统计表

1.3 注气模拟条件

CO2、N2是目前注气促进瓦斯解吸的主要气体[28-29]，注入压力对瓦斯解吸影响较大。在设置注入压力及温度后，确定含瓦斯煤模型中注入CO2、N2的数量（表2）。

Sorption 吸附模块中，设置任务项为Fixed Pressure，设置注入压力介于1 ～5 MPa，温度为303 K，模拟加载至平衡态步骤数为1×106步，总步骤数为2×106步；构型计算方法、电荷、能量计算方法及力场与1.2 部分相应设置一致。

注气促解模拟的初始构型即是在1.2 部分确定的含瓦斯煤吸附构型基础上（图2），分别注入不同数量的CO2、N2，然后进行下一步模拟。在注气促解模拟中，不同注入压力下CO2、N2的数量与表2中CO2、N2的平均吸附量一致。在Forcite 模块对体系进行动力学计算，采用NVT 系统对体系进行优化。Forcite 模块中，设置任务项为Dynamics，力场为COMPASS 力场，注入压力为1 MPa、2 MPa、3 MPa、4 MPa、5 MPa，温度为303 K，温度控制方法为Nose，压力控制方法为Berendsen，步数为10 000，时间步长为1 fs（1 fs=1×10－15s，下同），模拟时间为40 ps（1 ps=1×10－12s，下同）。

表2 不同注入压力下注入气体吸附特性参数统计表

图2 含瓦斯煤体注入N2/CO2 后初始构型图

2 结果与讨论

2.1 不同注入压力下含瓦斯煤体解吸构型能量对比与分析

温度为303 K、注入压力介于1 ～5 MPa，注入N2、CO2的含瓦斯煤体解吸构型如图3-a、b 所示。

如图3 所示，N2和CO2两种气体主要通过占据吸附位来驱替CH4，并且随着体系压力升高，真空层中游离气体分子逐渐增多。较之注入N2的构型，注入CO2的构型中游离的CH4分子明显更多，并且后者是前者的2 倍，可见注入CO2促进CH4解吸的效果优于N2。这是由于吸附CO2后煤表面自由能变化值普遍大于吸附CH4后的煤表面自由能变化值，并且CO2在煤基质表面存在多分子层吸附，而CH4和N2为单分子层吸附。对于孔径分布确定的煤分子结构，吸附质分子动力学直径（CH4为0.380 nm，N2为0.364 nm，CO2为0.330 nm）决定了其能够进入的煤体有效孔径，大于气体分子直径的煤体孔径即成为有效吸附孔径，前述3 种气体中CO2分子直径最小，由于分子直径与被吸附数量成反比，也使得CO2吸附性能优于CH4和N2。

图3 不同注入压力下含瓦斯煤体解吸构型图

图4 注入N2/CO2 初始构型、解吸构型Ev、Ee 统计图

图5 不同注入压力下CH4—N2/CH4—CO2 体系Etotal 统计图

如图4 所示，CH4—N2和CH4—CO2构型中体系Ev均与注入压力呈线性负相关，体系Ee与注入压力呈线性正相关。注入N2后体系Ev和Ee均大于初始构型，而注入CO2后体系Ev减小，而Ee增加，说明CH4—CO2体系中分子间作用力更强。如图5 所示，CH4—N2和CH4—CO2体系初始构型和解吸构型的体系总能量均随着注入压力增加而降低，并且呈线性负相关。CH4—N2和CH4—CO2构型在解吸后，较之初始构型，体系总能量均升高，可见解吸为吸热反应。与初始构型相比，在注入压力介于1 ～5 MPa条件下，注入N2，解吸构型的体系总能量分别升高73.72%、68.63%、63.76%、66.15%、63.86%，而注入CO2，解吸构型的体系总能量分别升高30.68%、18.23%、17.81%、15.20%、4.05%，后者体系总能量增加率低于前者，原因在于吸附为放热反应，并且CO2吸附能力远大于N2。CH4—CO2体系总能量整体低于CH4—N2体系，进而从微观角度说明CH4—CO2体系更稳定。

2.2 不同注入压力下含瓦斯煤体中气体分子扩散规律

2.2.1 均方根位移的计算

笔者利用前述Forcite 模块对体系进行动力学计算，优化后利用Focite 模块的Analysis 任务项进行计算，即可得到不同体系在不同注入压力下的均方根位移，该参数可以表征气体分子的运动规律。均方根位移计算式[21]为：

式中MSD 表示均方根位移，Å；ri(t)、ri(0)分别表示第i 个粒子在t 时刻与初始时刻的位置向量，无量纲；t 表示模拟时间，ps；N 表示吸附质分子数目，Nt表示统计平均的分子动力学步数；t0表示初始时刻。

如 表3 所 示，随 注 入 压 力 增 大，CH4—N2、CH4—CO2体系中CH4、CO2、N2的均方根位移均增大。在相同的注入压力条件下，不同体系中CH4、CO2、N2的均方根位移逐渐减小。CH4—N2、CH4—CO2体系中各气体均方根位移均随注入压力升高而增加。均方根位移表示t 时刻粒子所在位置与初始位置的距离，气体均方根位移随压力增大而增大，说明体系中分子活跃程度随压力增大而增强。在相同注入压力条件下，较之CH4—N2体系，CH4—CO2体系中CH4均方根位移更大，可见CH4—CO2体系中CH4分子活跃程度更高，CO2促进CH4解吸的效果更好。

表3 不同注入压力下CH4—N2、CH4—CO2 体系的均方根位移统计表

2.2.2 不同气体扩散系数的计算

扩散系数采用均方根位移曲线和Einstein 法计算得到，Einstein 法计算式[21]为：

式中D 表示扩散系数，m2/s。

将均方根位移曲线进行线性回归，得到斜率（k），如式（6）所示，则扩散系数（D）计算式可以简化为式（7）。

如图6 所示，CH4—N2/CH4—CO2体系中气体均方根位移均随着模拟时间增加而变大，并且在相同时间下，体系中气体均方根位移随压力升高也逐渐增加。N2/CO2注入煤体后，由于注入压力高于煤层地层压力，在压力梯度影响下吸附态CH4被注入气体驱替，并且和注入气体一起被采出。另外，N2/CO2的注入对真空层的CH4有稀释作用，使煤表面分子层CH4浓度大于真空层，在浓度差作用下煤表面分子层的CH4向真空层扩散。与之相反，由于真空层的游离态N2/CO2浓度大于煤表面分子层中N2/CO2的浓度，部分游离态N2/CO2将从游离态转化为吸附态。因此，在煤体中实际发生的是N2/CO2以扩散方式来置换CH4[30]。

图6 不同注入压力条件下CH4—N2/CH4—CO2 体系中气体均方根位移变化曲线图

如图7 所示，CH4扩散系数均随着N2/CO2注入压力增大而增大，呈线性正相关，并且N2、CO2、CH4扩散系数依次增大；注入压力介于1 ～5 MPa，CH4扩散系数平均为0.100 nm2/ps（CH4—N2体系）、0.119 nm2/ps（CH4—CO2体系），N2扩散系数平均为0.021 nm2/ps，CO2扩散系数平均为0.044 nm2/ps；较之CH4—N2体系，CH4—CO2体系中CH4扩散系数更大，可见CO2的促解吸效果优于N2。

2.3 注气后解吸构型的相对浓度及速度分布

2.3.1 注气后解吸构型的相对浓度分布

为了从微观角度进一步研究CH4、CO2及N2在煤结构模型中的分布，对注入气体后的解吸构型进行了相对浓度分析。如图8 所示，含瓦斯煤结构中注入CO2/N2后，模型真空层中CH4相对浓度均升高。并且，随着注入压力增加，真空层中CO2、N2、CH4的相对浓度均升高，总体上CH4—CO2体系中CH4相对浓度最高，CH4—N2体系中CH4、CH4—CO2体系中CO2相对浓度居中，CH4—N2体系中N2相对浓度最低，这说明体系中一部分CH4扩散到真空层，而CO2和N2仍主要集中于煤分子层，进而从微观角度证明CO2和N2主要是通过占据吸附位来驱替CH4。针对真空层而言，相同注入压力下CH4—CO2体系中CH4平均相对浓度比CH4—N2体系中CH4平均相对浓度增加了1 ～2 倍，可以看出，较之N2，注入CO2促进CH4解吸的效果更优，这与马志宏等[31]开展的实验结果也一致。若注入N2，在注入压力介于2 ～3 MPa 条件下促进CH4解吸的效果更好；若注入CO2，注入压力介于1 ～4 MPa 条件下的促进CH4解吸效果更好。同时，由于所建立的煤结构模型为对称模型，0 ～25 Å、75 ～100 Å 区域为煤分子层，25 ～75 Å 区域为真空层，气体相对浓度也呈对称分布（图8）。比较注气后CH4的扩散系数和相对浓度，注入N2促进CH4解吸的合理注入压力介于2 ～4 MPa，注入CO2促进CH4解吸的合理注入压力介于3 ～4 MPa。

图7 气体扩散系数与注入压力关系曲线图

图8 不同注入压力条件下CH4—N2/CH4—CO2 体系相对浓度分布图

2.3.2 注气后解吸构型的速度分布

如图9 所示，含瓦斯煤结构中注入CO2/N2后，模型真空层中CH4分子运动速度均上升，证明一部分CH4从分子层运动到真空层；随着注入压力增加，真空层CO2、CH4运动速度降低，N2运动速度上升，整体上，真空层CH4、CO2、N2的速度依次降低；在相同注入压力下，较之CH4—CO2体系，CH4—N2体系中CH4在真空层的平均速度偏高。

图9 不同注入压力条件下CH4—N2/CH4—CO2 体系速度分布图

3 结论

1）较之注入N2，含瓦斯煤体在注入CO2后，体系中游离CH4分子更多并且CH4—CO2体系总能量整体低于CH4—N2体系，前者更稳定。

2）在相同注入压力条件下，较之CH4—N2体系，CH4—CO2体系中CH4均方根位移更大，CH4分子活跃程度更高，并且CH4—CO2体系中CH4扩散系数更大，CO2促进CH4解吸的效果更好。

3）含瓦斯煤结构中注入CO2/N2后，随着注入压力增加，真空层中CO2、N2、CH4的相对浓度均升高，总体上CH4—CO2体系中CH4相对浓度最高，CH4—N2体系中CH4、CH4—CO2体系中CO2相对浓度居中，CH4—N2体系中N2相对浓度最低；注入N2促进CH4解吸的合理注入压力介于2 ～4 MPa，注入CO2促进CH4解吸的合理注入压力介于3 ～4 MPa。

4）所取得的研究成果可以为煤层注气增产技术研究提供理论支撑。由于此次研究仅对烟煤模型进行了分子模拟，后续将拓展研究范围，进一步完善研究成果。