由于我国“富煤、贫油、少气”的能源资源禀赋，煤层气的开发利用不仅可以促进煤炭资源的安全开采，而且对于增强我国天然气自主保障能力大有裨益^[1-3]。煤对CO₂的亲和力比CH₄更强，在正常的储层压力下，煤对CO₂的吸附能力是CH₄的2~10倍^[4-5]。通常情况下，负压抽采瓦斯的资源回收率只能达到30%~70%，在煤层中注入CO₂以驱替CH₄（CO₂-ECBM，CO₂ enhanced coal bed methane recovery）可以显著提高煤层气资源的采收率，是一种有效的煤层气增产技术，该技术集温室气体地质封存和清洁能源开发于一体，应用前景广阔^[6-8]。

CO₂-ECBM涉及煤体变形、CO₂/CH₄竞争吸附、流体运移、温度变化等复杂多物理场耦合过程^[9-11]。因此，建立准确的预测模型以降低投资风险已成为CO₂-ECBM相关领域专家学者研究的热点。WANG等^[12]基于非线性弹性变形力学和气体渗流理论对CO₂在煤层中的地质封存进行了数值模拟研究，结果表明，CO₂封存在煤层中促进了煤层气向产气井快速运移，且产气井生产周期较长，煤层气采收率可达到80%以上；HOU等^[13]利用COMSOL Multiphysics软件模拟了CO₂注入温度和压力对CH₄采收率的影响，得出CH₄采收率与CO₂注入温度成反比，而与注入压力成正比；FANG等^[14-15]在CO₂-ECBM建模方面做了大量的工作，通过分析煤基质收缩、膨胀和有效应力的变化，推导了煤层渗透率的演化模型，建立了气体扩散-吸附-渗流-换热场耦合方程，对CO₂封存量、CH₄产量以及CO₂注入的有效影响半径进行了讨论和分析；LIU等^[16]建立了CO₂-ECBM的全耦合THM数学模型，详细讨论了不同CO₂注入压力下工程井组的CO₂注入和CH₄产出特性，并得出CO₂注入引起的渗透率降低主要影响CO₂注入速率和储气容量，对CH₄产出影响较弱。

煤层渗透率是影响CO₂注入、煤层气产出的关键因素，前人的研究成果多集中在分析CO₂注入压力、温度等因素对CH₄产出效果和渗透率的影响，缺少对CO₂-ECBM过程中煤层渗透率时空演化规律的深入探讨。为此，建立了CO₂-ECBM的流-固-热多物理场全耦合方程，通过数值模拟方法研究CO₂-ECBM过程中煤层渗透率在时间和空间上的动态演化规律；研究结果对于及时调整CO₂-ECBM工程计划、指导煤层气开发具有重要的理论意义。

1.   CO₂-ECBM多物理场耦合方程

由于煤层瓦斯的主要成分为CH₄，在建立CO₂-ECBM多物理场耦合方程时，用CH₄代表煤层瓦斯。

1）气体运移控制方程。煤体对CO₂和CH₄的吸附/解吸过程用二元气体Langmuir等温吸附方程描述，气体在浓度梯度、压力梯度作用下扩散渗流，气体运移控制方程为：

式中：t为时间，s；l为气体组分；M_l为气体摩尔质量，kg/mol；φ为煤层孔隙率；C_l为气体浓度，mol/m³；ρ_s为煤体密度，kg/m³；$\rho_{{\mathrm{ga}}_l} $为气体在标况下的密度，kg/m³；V_Ll为Langmuir体积常数，m³/kg；R为普适气体常数，J/（mol·K）；T为温度，K；P_Ll为Langmuir压力常数，Pa；k为煤层绝对渗透率，m²；μ_l为气体动力黏度，Pa·s；D_l为气体扩散系数，m²/s；Q_l为气体质量源项，kg/（m³·s）。

2）煤体变形控制方程。基于线性热弹性假设建立煤体本构关系，即煤体的总应变为热应变、气体吸附/解吸应变、地应力诱导的应变之和。煤体应变控制方程为：

式中：E为煤体杨氏模量，Pa；u_i,jj、u_j,ji分别为煤体沿i和j方向的位移分量，m；ν为煤体泊松比；K_s为煤体体积模量，Pa；α_sg为煤体吸附应变系数，kg/m³；ΔV_ci为单位质量煤体所吸附气体的体积变化量，m³/kg¹；α_T为煤体热应变系数，K⁻¹；α为煤体的Biot系数；C_1,i、C_2,i分别为CO₂和CH₄浓度，mol/m³；F_i为沿i方向的体积力分量，Pa。

3）热量传输控制方程。在CO₂-ECBM过程中，煤体与流体共存于同一体积空间内，在不考虑热能与机械能转换的条件下，代表性体积单元上的热守恒方程可以表示为：

式中：（ρC_p）_eff为有效比热容，J/（m³·K）；η_eff为有效对流传热系数，W/（m²·K）；λ_eff为有效导热系数，W/（m·K）；q_l为气体组l的等量吸附热，J/mol；α_gl为煤体的吸附应变系数，kg/m³；V_l为单位质量煤体对气体组分l的吸附体积，m³/kg。

4）耦合项。CO₂-ECBM过程中煤层渗透率的演化受气体竞争吸附、孔隙压力、热量传输、煤体形变等复杂多物理场的影响，煤层渗透率演化的控制方程可表示为：

式中：k₀为煤层初始渗透率，m²；φ₀为煤层初始孔隙率；ε_e为煤体的体积应变；Δε_p为孔隙压力引起的煤基质应变；Δε_s为气体吸附引起的煤基质应变；Δε_T为温度变化引起的煤基质应变。

2.   CO₂-ECBM数值模拟模型与参数

基于COMSOL Multiphysics有限元软件，通过求解偏微分方程组实现对CO₂-ECBM真实物理过程的模拟。考虑数值计算的可行性和有效性，将CO₂-ECBM简化为二维过程，采用五点布井法，建立地质模型。CO₂-ECBM数值模拟地质模型如图1所示。

图  1  CO₂-ECBM数值模拟地质模型

Figure  1.  Geological model for CO₂-ECBM numerical simulation

下载: 全尺寸图片 幻灯片

地质模型尺寸为200 m×200 m，CO₂注入井位于模型中心，瓦斯抽采井与CO₂注入井相距100 m，钻孔直径均为0.1 m。为了方便观测煤层渗透率的动态演化规律，设置了观测线AB，以及观测点C（145,145）。地质模型除钻孔边界以外，均为零通量辊支承边界。

根据某矿的煤层瓦斯地质数据设置数值模拟物理参数，煤层初始瓦斯压力1.60 MPa，初始温度300.00 K，初始渗透率5.64×10⁻¹⁶ m²，抽采负压为13 kPa。其他关键参数为：煤体杨氏模量2.80×10³ MPa，煤体体积模量8.47×10² MPa，煤体密度1 470 kg/m³，二氧化碳动力黏度2.22×10⁻⁵ Pa·s，二氧化碳Langmuir压力常数1.38 MPa，二氧化碳Langmuir体积常数2.56×10⁻² m³/kg，二氧化碳比热容6.51×10² J/（kg·K），煤体比热容1.35×10² J/（kg·K），二氧化碳导热系数0.15×10⁻¹ W/（m·K），二氧化碳等量吸附热35.00 J/mol，煤体泊松比0.34，煤层初始孔隙率3.50%，煤体的热膨胀系数2.40×10⁻⁵ K⁻¹，甲烷动力黏度1.84×10⁻⁵ Pa·s，甲烷Langmuir压力常数2.07 MPa，甲烷Langmuir体积常数0.48×10⁻¹ m³/kg，甲烷比热容1.62×10³ J/（kg·K），煤体导热系数0.19 W/（m·K），甲烷导热系数0.31×10⁻¹ W/（m·K），甲烷等量吸附热33.40 J/mol。

在CO₂注入压力和注入温度分别为6 MPa和300 K的条件下，研究煤层渗透率在CO₂-ECBM过程中的演化规律。

3.   CO₂-ECBM过程中的煤层渗透率演化

煤体对CO₂分子的吸附性高于CH₄，随着CO₂-ECBM工程的进行，CO₂持续由注入井向瓦斯抽采井运移，在竞争吸附的作用下，CO₂分子占据了更多的吸附位，导致煤层中的吸附态CH₄分子被逐渐驱替，并解吸成为游离态，在压力梯度的作用下，不断向抽采井运移。

由煤层中CH₄浓度分布（图2）和CO₂浓度分布（图3）可知：驱替时间为500 d时，煤层中CH₄和CO₂的平均浓度分别为478.22、39.75 mol/m³；1000 d时，煤层中CH₄和CO₂的平均浓度分别为386.69、70.14 mol/m³；1500 d时，煤层中CH₄和CO₂的平均浓度分别为325.19、96.57 mol/m³；2000 d时，煤层中CH₄和CO₂的平均浓度分别为280.27、120.51 mol/m³。

图  2  煤层中CH₄浓度分布

Figure  2.  CH₄ concentration distribution in coal seam

下载: 全尺寸图片 幻灯片

图  3  煤层中CO₂浓度分布

Figure  3.  CO₂ concentration distribution in coal seam

下载: 全尺寸图片 幻灯片

由CO₂影响半径随时间变化图（图4）可知：CO₂影响半径随注气时间的增加而逐渐扩大；前600 d影响半径增长速率较大，CO₂的平均运移速率为0.053 m/d；之后，影响半径增长速率逐渐变缓，平均运移速率降至0.023 m/d；注气2000 d时，CO₂影响半径达81.10 m，此时CO₂仍没有突破煤层进入抽采井。

图  4  CO₂影响半径随时间变化图

Figure  4.  Change of CO₂ influence radius over time

下载: 全尺寸图片 幻灯片

煤层渗透率与孔隙率关系密切，即孔隙率越高，渗透率越高。煤基质的膨胀和收缩对孔隙率的变化起着决定性的作用。由于煤基质对CO₂的吸附性强于CH₄，在CO₂驱替CH₄的过程中，煤基质吸附CO₂诱导的吸附膨胀及热膨胀变形占主导地位，使得煤层孔隙率减小，因此CO₂的注入会降低注入井周围的煤层渗透率。

煤层中渗透率随时间的分布云图如图5所示，煤层中渗透率沿观测线AB的分布如图6所示，煤层中渗透率、CO₂浓度、CH₄浓度在观测点C处随时间的变化如图7所示。

图  5  渗透率分布云图

Figure  5.  Cloud diagrams of permeability distribution

下载: 全尺寸图片 幻灯片

图  6  渗透率沿观测线AB的分布

Figure  6.  Permeability distribution along the observed line AB

下载: 全尺寸图片 幻灯片

图  7  观测点渗透率/CO₂浓度/CH₄浓度变化

Figure  7.  Changes of permeability/CO₂ concentration/CH₄ concentration at observation points

下载: 全尺寸图片 幻灯片

由图5可知：500 d时，CO₂影响半径内、外的平均煤层渗透率分别为3.41×10⁻¹⁶、6.11×10⁻¹⁶ m²；1 000 d时，CO₂影响半径内、外的平均煤层渗透率分别为3.38×10⁻¹⁶、6.33×10⁻¹⁶ m²；1 500 d时，CO₂影响半径内、外的平均煤层渗透率分别为3.33×10⁻¹⁶、6.51×10⁻¹⁶ m²；2 000 d时，CO₂影响半径内、外的平均煤层渗透率分别为3.41×10⁻¹⁶、6.67×10⁻¹⁶ m²。可见，注入井周围的煤层渗透率低于初始渗透率，随CO₂影响半径的增加，低渗透率区域的范围也增大，但是CO₂影响半径内的平均煤层渗透率随时间的变化较小；在CO₂影响半径之外，在抽采负压的作用下，煤层中CH₄解吸诱导的煤基质收缩变形占优势，煤层渗透率随时间逐渐增加。煤层中CO₂的注入和CH₄的采出，对渗透率具有显著的影响。

由图6可知：受CO₂注入的影响，注入井周围的渗透率明显低于初始渗透率，煤层渗透率沿观测线A→B方向先快速增加，而后增长速率变缓，在接近CO₂影响范围边缘的区域内，煤层渗透率增长速率再次增加；注气500、1 000、1 500、2 000 d，CO₂的影响半径分别为43.85、60.70、72.43、81.00 m；在CO₂影响范围之内，煤层渗透率的增长速率随时间的增加而减小；在CO₂影响范围之外，煤层渗透率受CH₄解吸的影响，开始以非常缓慢的速率增加，在靠近抽采井的范围内，才发生显著的增长，总体而言，距离注入井越远，煤层渗透率越大。

由图7可知：自CO₂-ECBM工程之始，观测点C处的CH₄浓度一直呈单调下降趋势，CH₄的初始浓度为641.49 mol/m³，1 130 d时CH₄浓度为373.08 mol/m³；2 000 d时CH₄浓度为282.62 mol/m³，其下降速率逐渐减小；在CO₂未运移到观测点C处时，该点处的渗透率呈单调上升趋势，这主要是因为CH₄解吸导致的煤基质收缩变形，使得煤层孔隙率增加，进而使渗透率也随之增加，1 130 d时煤层渗透率达到最大值，为6.12×10⁻¹⁶m²，较初始渗透率增加了8.48%；1 130 d之后，CO₂运移至观测点C处，虽然煤层中的CH₄仍不断解吸，但是煤基质对CO₂的大量吸附，造成煤基质的吸附膨胀变形处于优势地位，煤基质占据了更多的孔隙空间，煤层渗透率急剧下降，2 000 d时煤层渗透率为4.78×10⁻¹⁶ m²，较初始渗透率下降了15.20%，较渗透率最大值下降了21.84%。

4.   结　语

1）建立了CO₂-ECBM流-固-热多物理场耦合数值模型，采用五点布井法构建地质模型，对CO₂-ECBM过程进行了数值仿真。驱替时间为2 000 d时，煤层中CO₂平均浓度增加至120.51 mol/m³、CH₄浓度则由初始的641.49 mol/m³下降至280.27 mol/m³，降幅达56.31%。

2）随着CO₂的持续注入，其影响半径不断增加，且增长速率逐渐变缓，前600 d，平均增长速率为0.053 m/d，之后逐渐变缓，平均增长速率降至0.023 m/d；CO₂影响范围内的煤层渗透率低于初始渗透率，低渗透率区域随CO₂影响半径的增加而增加；CO₂影响范围之外的区域内，煤层渗透率随CH₄解吸逐渐增加。

3）沿CO₂注入井至瓦斯抽采井，煤层渗透率呈上升趋势，在CO₂影响范围之内，煤层渗透率的增长速率随时间的增加而减小。煤层在不受CO₂影响时，CH₄解吸诱导的煤基质收缩变形占优势地位，煤层孔隙率增加，煤层渗透率随时间单调递增；煤层受CO₂影响之后，煤基质对CO₂吸附诱导的煤基质膨胀占主导地位，煤层孔隙率降低，煤层渗透率急剧下降。