CN107621350B - 一种模拟超临界co2驱替天然气流动的方法 - Google Patents

Abstract

本发明提出了一种模拟超临界CO₂驱替天然气流动的方法，所述方法包括如下步骤：根据气藏地质特征，建立三维地质模型；通过改进状态方程建立流体模型确定气藏流体的基础参数；确定CO₂的扩散系数；建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型；应用模型计算不同CO₂注入量和天然气采出量的对应关系。本发明提出了新的数学模型及处理方法模拟超临界CO₂驱替天然气流动特征，与常规数值模型相比能更好的评价CO₂驱气流动特点及驱替效果。

Description

一种模拟超临界CO2驱替天然气流动的方法

技术领域

本发明属于气藏开发领域，具体地，涉及一种模拟超临界CO₂驱替天然气流动的方法。

背景技术

由于人类对石化燃料(煤、石油、天然气)的过度依赖，工业和人类生活过程中产生的温室气体排放量日益增加，由此导致的空气污染和温室效应正在严重地威胁着人类赖以生存的环境。温室气体二氧化碳(CO₂)使全球气候变暖，对人类的生存和社会经济的发展构成了严重的威胁。CO₂有效利用已成为了全球广泛关注的课题。通过不断的发展及完善，在油藏中注入CO₂成为较为成熟的技术。

但在气藏中注入CO₂相态和驱替机理与油藏完全不同，注CO₂提高气藏采收率成为一项全新的课题。国外在上世纪九十年代中后期开始CO₂EGR研究，国内近两年初步开展室内理论研究。但目前国内外仅开展实验及机理模拟，仍处于探索阶段，未进入现场应用，此项研究仍处于起步阶段。

在CO₂与天然气驱替模拟方面，国内外均应用数值模拟软件开展了CO₂提高采收率机理研究。Oldenburg等人应用TOUGH2软件开展等温多相多组分数值模拟，模拟结果表明底部注气有利于减少CO₂向上锥进和延缓CO₂与CH₄混合，提高平面及垂向驱替效率。TorstenClemens以实例气田为例模拟CO₂提高采收率驱替效率。模拟结果表明，气藏开采到废弃阶段注入CO₂能提高采收率 9.4％。A.AL-Hashami应用机理模型模拟底水气藏CO₂驱替效率，通过数模研究表明CO₂有助于提高气藏采收率，能增产8-11％。目前，普遍采用商业数值模拟软件分析CO₂提高采收率驱替封存效率，模型考虑因素较为简单。未考虑CO₂扩散对驱替效果的影响。CO₂驱气机理不明确，驱替机理论证不系统，不充分。

发明内容

为了解决现有技术的不足，本发明探索提出体积改进状态方程方法，解决低温时高含CO₂气体实验数据拟合难的问题；并且考虑了不同压力条件下CO₂扩散对驱替天然气的影响，建立了考虑扩散的超临界CO₂驱天然气多组分数值模拟模型，模拟超临界CO₂驱替天然气驱替机理及驱替效果。

根据本发明的一个方面，提供一种模拟超临界CO₂驱替天然气流动的方法，所述方法包括如下步骤：

根据气藏地质特征，建立三维地质模型；

通过改进状态方程建立流体模型确定气藏流体的基础参数；

确定CO₂的扩散系数；

建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型；

应用模型计算不同CO₂注入量和天然气采出量的对应关系。

进一步地，确定气藏流体的基础参数，包括：实验测试不同CO₂天然气混合物偏差因子、密度及粘度。

进一步地，改进状态方程，获得描述整个压力温度范围内临界性质的热力学参数。

进一步地，引入体积修正项对状态方程进行修正，对临界点体积参数进行改进，改进项由两部分组成：①临界点改进；②温度改进：

c(T)＝c_cf(T_r) (1)

其中，c(T)体积修正项，m³/mol；c_c临界点改进体积修正项，m³/mol；f(T_r) 温度改进系数，小数。

进一步地，引入临界点改进项：

其中，临界压缩因子Zc为定值，R气体常数，8.314J/mol·K；T_c临界温度， K；P_c临界压力，MPa；v_EOS状态方程计算摩尔体积，m³/mol；v_Exp实验计算摩尔体积，m³/mol；c_c体积修正项，m³/mol。

进一步地，引入温度改进项进行修正，温度改进通过相平衡实验确定Z_c与T 的对应关系，通过实验拟合得方程式(3)-(5)：

f(T_r)＝β+(1-β)exp(0.5β+8.3) (3)

β＝-2.8431exp([-64.2184(0.3074-Z_c)]+0.1735) (4)

β＝T-T_r (5)

其中，T气体所处的温度，K；T_r＝T/T_c对比温度，K；β温度校正项，K。

进一步地，将临界点改进项c(T)＝c_cf(T_r)带入PR状态方程：

状态方程中，ab由临界参数和偏心因子计算。式中，

a＝α(T)·a_c；

m＝0.37464+1.54226ω-0.26992ω²；

T_r＝T/T_c；

b＝0.77796·R·T_c/P_c

其中，p压力，MPa；v摩尔体积，m³/mol；ω偏心因子，无因次。

进一步地，确定CO₂的扩散系数包括：实验测量不同压力条件下CO₂在天然气中的扩散浓度变化规律，应用插值计算方法描述多孔介质中CO₂的扩散系数。

进一步地，建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型包括：应用改进状态方程劈分计算各组分热力学参数，质量守恒方程中考虑扩散项对流动的影响。

进一步地，考虑扩散的超临界CO₂驱天然气多组分数值模拟模型包括，如下：

其中，k渗透率，D；k_rg气相渗透率，小数；μ_g气相粘度mPa·s；ρ_g气相密度kg/m³；y_ig组分i在气相中的摩尔分数，小数；气相的势；φ孔隙度，小数；S_g气相饱和度，小数；q_ik组分i在k层中的产量或流量，m²/d；D_ig为组分i在气相中的扩散系数，m²/d；n_f为井的射孔层段；

扩散张量：

式中，

其中，ρ_i组分i密度kg/m³；S_j j相饱和度，小数；D_ig为组分i在气相中的扩散系数，m²/d；为组分i在j相中的分子扩散系数，m²/d；u_j为速度相量，m²/d；u_jx为x方向速度相量，m²/d；u_jy为y方向速度相量，m²/d；u_jz为 z方向速度相量，m²/d；τ_j迂曲率，小数；j相的沿程弥散系数，小数；α_jt j相的横向弥散系数，小数。

进一步地，应用模型计算不同CO₂注入量和天然气采出量的对应关系包括，根据长岩心驱替实验的物性参数，建立实验尺度数值模型，导入束缚水饱和度、渗透率、孔隙度参数，计算不同CO₂注入量和天然气采出量的对应关系。

本发明解决了超临界CO₂天然气混合物近临界区高压物性拟合难的问题。

常规气体偏差系数计算模型主要有：Dranchuk–Purvis-Robinson(DPR)法、 Hall&Yarborough(HY)法、Sarem方法、Dranchuk-Abu-Kassem(DAK)法、李相方(LXF)法等。采用经验公式对比实验数据，发现计算误差较大；而应用状态方程，高温区域热力学参数计算误差小，但在低温(近临界温度)区且处于近临界区(近临界压力)时，热力学参数误差较大。

应用改进的状态方程，计算偏差因子与实验值误差小于3％、计算不同CO₂含量天然气混合物高压物性参数(密度、粘度与实验值)误差小于3％。

本发明考虑不同压力条件下CO₂扩散对驱替天然气的影响，描述超临界CO₂天然气驱替前沿及驱替效率。

通过长岩心驱替实验，获得CO₂注入量与出口端CO₂摩尔含量及天然气采收率对应关系数据。根据实验数据，建立长岩心驱替模型，在相态拟合的基础上，分别应用常规组分模型和本文考虑吸附的组分模型计算在长岩心中CO₂驱替天然气，使用本发明方法建立的组分模型计算的采收率与实验值更接近。

附图说明

通过结合附图对本公开示例性实施方式进行更详细的描述，本公开的上述以及其它目的、特征和优势将变得更加明显，其中，在本公开示例性实施方式中，相同的参考标号通常代表相同部件。

图1为根据本发明实施例的方法流程图。

图2显示了Z_c与β实验测试拟合图。

图3显示了CO₂在甲烷中的扩散系数。

图4显示了状态方程热力学参数计算结果。

图5显示了不同CO₂含量天然气混合物密度拟合图。

图6显示了不同CO₂含量天然气混合物粘度拟合图。

图7显示了注入CO₂孔隙体积与出口端CO₂摩尔含量计算数据。

图8显示了注入CO₂孔隙体积与采收率计算数据。

具体实施方式

下面将参照附图更详细地描述本公开的优选实施方式。虽然附图中显示了本公开的优选实施方式，然而应该理解，可以以各种形式实现本公开而不应被这里阐述的实施方式所限制。相反，提供这些实施方式是为了使本公开更加透彻和完整，并且能够将本公开的范围完整地传达给本领域的技术人员。

本发明提供了一种模拟超临界CO₂驱替天然气流动的方法，属于注气提高采收率模拟方法领域。本方法以组分模型为基础，基于超临界CO₂气体的相行为特殊，考虑近临界及温度的影响，改进状态方程，探索描述超临界CO₂注气过程中与天然气混合的组分表征方法。在流动表征方面，在流动项表征时考虑不同压力CO₂扩散对驱替天然气的影响，建立了考虑扩散的超临界CO₂驱天然气多组分数值模拟模型，实现了超临界CO₂驱气前缘的模拟，本发明提出了新的数学模型及处理方法模拟超临界CO₂驱替天然气流动特征，与常规数值模型相比能更好的评价CO₂驱气流动特点及驱替效果。

目前，超临界CO₂驱天然气多组分数值模拟模型存在以下问题：

在组分表征方面，主要为经验公式法及状态方程法。经验公式法参变量少，易于计算，但是计算误差较大；而状态方程以热力学理论为基础，描述分子间引力斥力作用，模型通用性很好，计算热力学参数为组分模拟器中各组分劈分的基础，但是，在近临界区高压物性计算精度差。CO₂气体的相行为特殊，在近临界区严重偏移理想气体，CO₂-天然气混合物也具有相似的行为。因此，本发明探索提出描述超临界CO₂注气过程中与天然气混合的组分表征方法。

在流动表征方面，超临界CO₂驱替天然气用达西定律表征，但是，实验表明CO₂扩散对驱替前沿存在影响，流动项表征时需考虑扩散对CO₂驱替的影响。因此，需要在实验的基础上，考虑不同压力条件下CO₂扩散对驱替天然气的影响。

针对以上存在的问题，本发明提出了新的数学模型及处理方法模拟超临界 CO₂驱替天然气流动特征，比常规数值模型相比能更好的评价CO₂驱气流动特点及驱替效果。

参照图1，本发明提出了一种模拟超临界CO₂驱替天然气流动的方法，所述方法包括如下过程：

根据气藏地质特征，建立三维地质模型；

通过改进状态方程建立流体模型确定气藏流体的基础参数；

确定CO₂的扩散系数；

建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型；

应用模型计算不同CO₂注入量和天然气采出量的对应关系。

进一步地，确定气藏流体的基础参数，包括：实验测试不同CO₂天然气混合物偏差因子、密度及粘度；通过引入体积修正项改进状态方程，获得描述整个压力温度范围内临界性质的热力学参数。

优选地，确定扩散系数包括：实验测量不同压力条件下CO₂在天然气中的扩散浓度变化规律，应用插值计算方法描述多孔介质中CO₂的扩散系数。

优选地，建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型，包括：应用改进状态方程劈分计算各组分热力学参数，质量守恒方程中考虑扩散项对流动的影响。

优选地，根据长岩心驱替实验的物性参数，建立实验尺度数值模型，导入束缚水饱和度、渗透率、孔隙度等基本参数，计算不同CO₂注入量和天然气采出量的对应关系，验证模型的实用性。

优选地，根据气藏地质特征，建立三维地质模型。导入气田生产井产气量、注入CO₂量及压力数据，进行流体流动预测分析，计算不同方案生产井CO₂含量及天然气采出程度。

优选地，超临界CO₂-天然气状态方程改进方法，如下：

CO₂相行为特殊，当温度大于31.04C，压力大于7.38MPa，CO₂具有超临界特性，超临界CO₂密度与液体接近，大于天然气密度。在近临界区CO₂物性严重偏移理想气体，CO₂-天然气混合物也具有相似的行为。应用常规状态方程， CO₂-天然气混合物高压物性在近临界区难以拟合，另外，当温度较低，在近临界温度附近时，计算误差较大。

因此引入体积修正项对状态方程进行修正，反映CO₂-天然气混合物相变特征。改进项首先考虑立方体状态方程自身存在的缺陷，立方体状态方程临界压缩因子Zc为定值，最为接近真实气体的PR状态方程Z_c为0.307，而CO₂Z_c 0.274 真实气体Z_c0.23-0.29，存在差异，因此需要对临界点体积参数进行改进。改进项由两部分组成：①临界点改进；②温度改进：

c(T)＝c_cf(T_r) (1)

其中，c(T)体积修正项，m³/mol；c_c临界点改进体积修正项，m³/mol；f(T_r) 温度改进系数，小数。

引入临界点改进项：

其中，临界压缩因子Zc为定值，R气体常数，8.314J/mol·K；T_c临界温度， K；P_c临界压力，MPa；v_EOS状态方程计算摩尔体积，m³/mol；v_Exp实验计算摩尔体积，m³/mol；c_c体积修正项，m³/mol。

引入温度改进项进行修正，温度改进通过相平衡实验确定Z_c与T的对应关系，通过实验拟合得方程式(3)-(5)：

f(T_r)＝β+(1-β)exp(0.5β+8.3) (3)

β＝-2.8431exp([-64.2184(0.3074-Z_c)]+0.1735) (4)

β＝T-T_r (5)

其中，T气体所处的温度，K；T_r＝T/T_c对比温度，K；β温度校正项，K。

进一步地，将临界点改进项c(T)＝c_cf(T_r)带入PR状态方程：

状态方程中，ab由临界参数和偏心因子计算。式中，

a＝α(T)·a_c；

m＝0.37464+1.54226ω-0.26992ω²；

T_r＝T/T_c；

b＝0.77796·R·T_c/P_c

其中，p压力，MPa；v摩尔体积，m³/mol；ω偏心因子，无因次。

应用改进的PR状态方程(6)劈分各组分物性参数，获得超临界CO₂驱替天然气过程中整个压力温度范围内的热力学参数。

进一步地，考虑扩散的超临界CO₂驱天然气多组分数值模拟模型包括，如下：

其中，k渗透率，D；k_rg气相渗透率，小数；μ_g气相粘度mPa·s；ρ_g气相密度kg/m³；y_ig组分i在气相中的摩尔分数，小数；气相的势，MPa；φ孔隙度，小数；S_g气相饱和度，小数；q_ik组分i在k层中的产量或流量，m²/d；D_ig为组分i在气相中的扩散系数，m²/d；n_f为井的射孔层段。

扩散张量：

式中，

其中，ρ_i组分i密度kg/m³；S_j j相饱和度，小数；D_ig为组分i在气相中的扩散系数，m²/d；为组分i在j相中的分子扩散系数，m²/d；u_j为速度相量，m²/d；u_jx为x方向速度相量，m²/d；u_jy为y方向速度相量，m²/d；u_jz为 z方向速度相量，m²/d；τ_j迂曲率，小数；j相的沿程弥散系数，小数；α_jt j相的横向弥散系数，小数。

优选地，CO₂扩散系数的求取，如下：

应用高压PVT仪器测量容器中不同压力条件下不同部位超临界CO₂与天然气混合浓度变化规律，计算CO₂分子扩散系数。扩散系数与压力有关，绘制扩散系数与压力的关系图，见图3。

与压力有关，通过插值确定不同压力时CO₂在甲烷中的分子扩散系数。代入考虑扩散的超临界CO₂驱天然气多组分数值模拟模型计算。

为便于理解本发明实施例的方案及其效果，以下给出一个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。

下面参照图4至图8，描述本发明的一个应用例，将本发明应用于CO₂驱替长岩心实验模拟。

首先，根据气藏地质特征，建立三维地质模型，开展网格划分并粗化处理。应用CMG软件，以长岩心实验为例，建立实验尺度长岩心地质模型，网格确定为横向上20个网格，纵向上1层，网格数量为20×1×1＝20个。横向网格长度为5cm，纵向上为2.5cm。平均孔隙度9.9％，渗透率0.652mD，束缚水饱和度为37.3％。横向上第1个网格处设为注入端，第20个网格处设为产出端。每个网格原始地层压力设为25MPa。

接下来，确定流体参数。实验测试不同CO₂天然气混合物偏差因子、密度及粘度；应用改进状态方程，计算不同CO₂天然气混合物偏差因子、密度及粘度，验证模型的准确性。应用改进状态方程，劈分各组分，获得描述整个压力温度范围内临界性质的热力学参数。将获得的各组分热力学参数录入到CMG软件中。

接下来，确定CO₂扩散系数。实验测量CO₂在天然气中的不同压力条件下不同部位CO₂浓度变化规律。通过实验测试，在扩散作用下，随着压力的不断降低，CO₂与天然气的混合时间越来越短，25MPa时CO₂和天然气混合的时间为210小时左右，15MPa时两者混合的时间为120小时左右，而5MPa下为25 小时左右。进一步定量评价二者的扩散速度，计算多孔介质条件下CO₂在甲烷中的扩散系数。将菲克定律推广到测定气体通过岩石的扩散系数。假定岩心长度为Z(而不是孔径长度)，面积为S(而不是有效面积)，则求得的扩散系数D 应是有效扩散系数。公式为

其中

式中：C₀‘、C₀‘’——定压取样装置气体反应区上部取样点CO₂和下部取样点CO₂的初始浓度；

C_t‘、C_t‘’——t时间定压取样装置气体反应区上部取样点CO₂和下部取样点CO₂的浓度；

V‘、V‘’——CO₂与天然气的初始体积；

S——定压取样装置内部截面积；

Z——定压取样装置上部和下部取样点的距离。

通过计算获得定压取样装置内25MPa时CO₂在天然气中的扩散系数为 2.74×10^{- 4}cm²/s，15MPa时的扩散系数为6.55×10^-4cm²/s，5MPa时的扩散系数为 1.57×10^-3cm²/s。

接下来，建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型。与压力有关，通过插值确定不同压力时CO₂在甲烷中的分子扩散系数，应用于考虑扩散的超临界CO₂驱天然气多组分数值模拟模型(公式6)。在模型应用中，对不同压力步录入扩散系数。组分模型求解方式采用自适应隐式方法。

接下来，模拟多级降压开采，第20个网格处设为产出端，定压力求产。降压间隔均为3MPa，最终压力降至8MPa。第1个网格设为注入端，注入端压力为11MPa，出口端回压为8MPa，驱替压差3MPa，定压力求产。

模型结果能显示不同网格压力，CO₂分布特征，密度粘度的变化规律，及不同CO₂注入量天然气的产出量。通过数模计算，当出口端CO₂为10％时，比衰竭方案提高采收率9.46％。数值模型论证了超临界CO₂提高气藏采收率的可行性。

表1长岩心CO₂驱替实验结果与数值模型计算结果对比表

应用上述方程，计算实验数据，与实验值拟合误差小于5％注入倍数在0.5 左右时出口端气体中CO₂摩尔含量超过10％。

以上已经描述了本公开的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。本文中所用术语的选择，旨在最好地解释各实施例的原理、实际应用或对市场中的技术改进，或者使本技术领域的其它普通技术人员能理解本文披露的各实施例。

Claims (10)

1.一种模拟超临界CO₂驱替天然气流动的方法，其特征在于，所述方法包括如下步骤：

根据气藏地质特征，建立三维地质模型；

通过改进状态方程建立流体模型确定气藏流体的基础参数；

确定CO₂的扩散系数；

建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型；

应用模型计算不同CO₂注入量和天然气采出量的对应关系。

2.根据权利要求1所述的方法，其特征在于，确定气藏流体的基础参数，包括：实验测试不同CO₂天然气混合物偏差因子、密度及粘度；

通过改进状态方程获得描述整个压力温度范围内临界性质的热力学参数。

3.根据权利要求1所述的方法，其特征在于，确定CO₂的扩散系数包括：实验测量不同压力条件下CO₂在天然气中的扩散浓度变化规律，应用插值计算方法描述多孔介质中CO₂的扩散系数。

4.根据权利要求1所述的方法，其特征在于，建立考虑扩散的超临界CO₂驱天然气多组分数值模拟模型包括：应用改进状态方程劈分计算各组分热力学参数，质量守恒方程中考虑扩散项对流动的影响。

5.根据权利要求1所述的方法，其特征在于，应用模型计算不同CO₂注入量和天然气采出量的对应关系包括，根据长岩心驱替实验的物性参数，建立实验尺度数值模型，导入束缚水饱和度、渗透率、孔隙度参数，计算不同CO₂注入量和天然气采出量的对应关系。

6.根据权利要求2所述的方法，其特征在于，引入体积修正项对状态方程进行修正，对临界点体积参数进行改进：

c(T)＝c_cf(T_r) (1)

其中，c(T)体积修正项，m³/mol；c_c临界点改进体积修正项，m³/mol；f(T_r)温度改进系数。

7.根据权利要求6所述的方法，其特征在于，引入临界点改进项：

其中，临界压缩因子Zc为定值，R气体常数；T_c临界温度，K；P_c临界压力，MPa；v_EOS状态方程计算摩尔体积，m³/mol；v_Exp实验计算摩尔体积，m³/mol；c_c体积修正项，m³/mol。

8.根据权利要求7所述的方法，其特征在于，引入温度改进系数进行修正，温度改进系数通过相平衡实验确定，进而建立Z_c与T的对应关系，通过实验拟合得方程式(3)-(5)：

f(T_r)＝β+(1-β)exp(0.5β+8.3) (3)

β＝-2.8431exp([-64.2184(0.3074-Z_c)]+0.1735) (4)

β＝T-T_r (5)

其中，T气体所处的温度，K；T_r＝T/T_c对比温度，K；β温度校正项，K。

9.根据权利要求7所述的方法，其特征在于，将临界点改进项和温度改进系数带入PR状态方程：

状态方程中，ab由临界参数和偏心因子计算，式中，

a＝α(T)·a_c；

α(T)＝[1+m(1-T_r ^0.5)]²；

m＝0.37464+1.54226ω-0.26992ω²；

T_r＝T/T_c；

a_c＝0.457235·R²·T_c ²/P_c；

b＝0.77796·R·T_c/P_c

其中，p压力，MPa；v摩尔体积，m³/mol；ω偏心因子，无因次。

10.根据权利要求8所述的方法，其特征在于，考虑扩散的超临界CO₂驱天然气多组分数值模拟模型包括，如下：

其中，k渗透率，D；k_rg气相渗透率；μ_g气相粘度mPa·s；ρ_g气相密度kg/m³；y_ig组分i在气相中的摩尔分数；气相的势，MPa；φ孔隙度；S_g气相饱和度；q_ik组分i在k层中的产量或流量，m²/d；D_ig为组分i在气相中的扩散系数，m²/d；n_f为井的射孔层段；

扩散张量：

式中，

其中，ρ_i组分i密度kg/m³；S_j j相饱和度；D_ig为组分i在气相中的扩散系数，m²/d；为组分i在j相中的分子扩散系数，m²/d；u_j为速度相量，m²/d；u_jx为x方向速度相量，m²/d；u_jy为y方向速度相量，m²/d；u_jz为z方向速度相量，m²/d；τ_j迂曲率；α_jlj相的沿程弥散系数；α_jt j相的横向弥散系数。