基于孔隙网络模拟的致密砂岩气充注与微观气水赋存特征

2022-06-07吴克柳何敏侠李相方贾承造

天然气工业 2022年5期

赵文吴克柳姜林何敏侠李相方贾承造,4

1. 油气资源与探测国家重点实验室·中国石油大学（北京） 2. 中国石油勘探开发研究院3. 中国石油长庆油田公司勘探开发研究院 4. 中国石油天然气集团有限公司

0 引言

非常规油气的成功开发推动了全球油气产量增长，致密砂岩气是非常规天然气勘探开发的重要对象，是我国油气接替上产的重要资源[1-2]。在致密砂岩气的勘探开发过程中，人们发现致密气储层中气水分布十分复杂，这极大制约了地质储量的有效动用[3-4]。复杂气水分布特征的形成实质上是充注过程中动力学机制作用结果。致密砂岩气充注机制与常规气藏存在明显差异[5-6]，深化致密气充注机制研究可以深入致密气成藏机理认识[7-8]，提升对非常规资源的预测和评价能力，对进一步高效勘探开发我国致密砂岩气具有重要意义。

常规油气成藏以圈闭富集保存油气为核心，其成藏动力主要为浮力；非常规油气则以连续性聚集为特征，是一种非浮力聚集，成藏主要动力为生烃膨胀力，主要阻力是毛细管力。但目前对于微纳米孔隙中毛细管动力学机制与流体赋存特征认识尚存在局限性。毛细管力实质是孔隙壁面与孔隙内流体之间的分子间作用力的宏观表达[9]。目前针对毛细管力的研究常用Young-Laplace方程进行描述，由Young-Laplace方程可知，两相流体的界面张力、毛细管半径和岩石润湿性这三个因素共同控制了地层中毛细管压力[10]。成藏过程中温度[11]、压力[12]、岩石润湿性[13]等因素都将对毛细管压力造成影响。深入理解致密气成藏充注过程中毛细管压力作用变化特征，有利于理解致密气形成分布特征，指导大规模高效开采，有助于提高对非常规油气开发生产机理和潜在生产能力的认识。

针对致密砂岩气成藏过程中毛细管力研究，国内外学者从理论研究[14-15]、物理模拟[16-17]与数值模拟[18-21]3个方面开展了大量工作。在理论研究上，人们着重于对毛细管力的自身特性[14]与其在油气充注运移[15]过程中发挥的作用进行描述，缺少针对微观机理的深入研究。而物理模拟主要包括砂箱模拟[16]和岩心渗流实验[17]。通过砂箱模拟，观测到油气充注过程中的流体流动特征，但由于实验条件与地质条件相差较大，实验结果难以反映真实流动过程。岩石渗流实验虽然还原了地质条件，但实验理论分析均建立在间接数据处理上，缺少直观可视化分析。目前针对油气成藏过程的数值模拟主要集中在毛细管束模型下的数值计算[18]，表征单元体尺度下的格尔玻尔兹曼模拟[19]，以及盆地尺度下的盆地模拟[20]，在孔隙网络尺度下的数值模拟研究较少，同时微观上的理论机理难以在数值模拟中体现。

为此，本文从分子间作用力出发，探讨了致密砂岩储层中孔隙壁面—天然气—水之间的相互作用及其在毛细管力上的宏观表达，并基于苏里格气田致密砂岩气藏孔隙结构特征，利用数字岩心技术生成了孔隙网络模型，开展了岩心尺度下的天然气充注数值模拟研究，揭示了含气性增长机制，进行了含气性预测与评价。本研究成果对深入理解致密油气成藏机理，开展气水分布评价，指导致密砂岩气勘探开发具有重要指导和借鉴意义。

1 典型致密砂岩气地质特征

鄂尔多斯盆地苏里格气田属上古生界含气系统，烃源岩主要为二叠系太原组和山西组的一套海陆交互相含煤层系，包括暗色泥岩和煤层；储集层主要为山西组的致密砂岩段。从源—储空间分布关系看，属于源—储垂向叠置的下生上储式配置关系（山2段—太原组为烃源层，盒8段—山1段为储集层），这是苏里格气田主要成藏模式。

在底部烃源岩生烃超压的驱动下，油气沿网状输导通道向山1段、盒8段储集层运移，并在相对发育的储集体中富集成藏。随着下部气源的不断补充，天然气逐渐向顶部的盒8段储集层运移、富集。在充注与运移过程中，天然气在生烃膨胀力推动下，沿网状输导体系向储集层充注运移，聚集成藏。天然气自烃源岩进入输导体系之前，由于岩石水润湿性，因此油气必须克服输导体系内毛细管压力才能进入到输导体系中。与常规油气藏不同的是，常规油气藏中的动力（浮力）是具有方向性的力，浮力始终给天然气一个垂直向上的力以推动气体运移，但是在致密储层中，作为动力的生烃膨胀力不具备方向性，这使得致密储层内的天然气运移过程和机理与常规储层之间存在巨大差异。

研究区山西组包裹体均一温度最高为180 ℃，平均温度为115 ℃，目前地温为83～88 ℃。苏里格地区成藏时期地层压力为65 MPa，而目前地层压力介于23～26 MPa[17]，这种温度与压力的差异将使毛细管压力在成藏过程与目前生产过程存在差异。本文将具体讨论地质条件下温度与压力变化对孔隙壁面—天然气—水之间的相互作用，及其在毛细管力上的宏观表达和变化特征。

2 孔隙网络模拟参数确定

天然气—水—岩石的润湿性和界面张力特征，反映了孔隙中流体—岩石、流体—流体之间的相互作用[21]。本文将建立一种半解析的接触角和界面张力随地层温度压力变化的计算方法，并将其嵌入到孔隙网络数值模拟方法中。

2.1 接触角特征

润湿性反映了孔隙壁面与流体之间的相互作用，一般利用接触角来描述宏观上岩石对流体的润湿性特征。在微观上，孔隙壁面与流体之间的作用主要发生在与孔隙壁面接触的薄液膜之间，一般用分离压来表征薄液膜中由于分子间作用力产生的额外压力[22-23]。Frumkin—Derjaguin方程表述了宏观接触角与微观孔隙壁面和孔隙内流体的作用关系[24]，表达式如下：

式中θ表示接触角，（°）；Π表示分离压，MPa；σ表示界面张力，mN/m；f表示液膜厚度，m。

基于扩散双电层模型基础上提出的疏液溶胶稳定性理论（Derjaguin-Landau-Verwey-Overbeek，DLVO理论）[25]，孔隙内流体与壁面之间的分离压力由三种分子间作用力构成，分别为静电力（Fe）、结构力（Fs）和范德华力（FvdW）。在这三种分子间作用力的共同作用下，孔隙壁面与流体之间作用合力可以表现为一种吸引力，亦可以表现为一种排斥力。

致密砂岩气藏中，由于天然气的成分相对简单，极性物质较少，因此孔隙内流体与壁面之间静电力和结构力较弱，在实际计算中可省略[26]。因此，分离压可以直接用范德华力描述，范德华力的大小为：

式中FvdW表示范德华力，N；HvdW表示哈梅克常数（Hamaker Constant，表示物质之间范德华吸引能的大小[27]），J。

在一定的温度和压力条件下，哈梅克常数可以通过折射率与介电常数得到[24]。宏观接触角与微观孔隙壁面和孔隙内流体的作用关系可以表示为：

Hamaker常数随地层不同而表现出差异。为了使计算更接近实际地层情况，在实验室条件下，利用润湿角测定仪，对苏里格地区山西组实际岩心进行了接触角测量。测量过程中实验温度为295.15 K，实验压力为0.1 MPa，平均接触角为38°，据此确定了计算过程中使用的哈梅克常数值。

2.2 界面张力特征

界面张力这一参数表征了气水界面之间的相互作用，定量描述了气水界面动力学特征。从微观上看，界面张力反映了在一定温度和压力条件下气体分子与水分子之间的相互作用，以及它们在化学性质与分子密度间的差异。前人已经通过大量实验测量[28-30]和分子动力学模拟[31]阐明了不同温度压力条件下流体系统中界面张力特征。本文利用拟合的半经验—半解析方程给出了气水界面张力与气水密度差和对比温度之间的关系[式（4）]。该方程结构简单，方便应用于工程计算与理论分析，并可嵌入到孔隙网络数值模拟程序中，进而实现模拟界面张力随温度、压力等地层条件变化而改变的过程。

式中σgw表示气水界面张力，mN/m；ρg、ρw分别表示气、水密度，g/cm3；Tr表示地层对比温度，K；A、B、C和D为拟合得到的4个常数。

本文以天然气为研究对象，利用167例甲烷/水的界面张力岩样测定案例[32-34]，基于非线性回归方法拟合得到式（4）中4个常数的最优值（图1），即A=82.61、B=6.64、C=－2.25和D=41.85，图 1中对角线为半解析—半经验模型计算结果与已发表的167个案例实验结果的比较线[35-39]。对比结果良好，主要数据偏差不超过5%，证实可靠。其中167例的温度范围介于273～ 373 K，压力范围介于0.1～ 100 MPa。

图1 界面张力结果与前人案例的对比图

天然气密度可以通过气体状态方程计算得到[35]：

式中M表示气体的摩尔质量，kg/mol，甲烷的摩尔质量为16×10－3kg/mol；R表示气体常数，取8.314 J/(K·mol)；Z表示气体偏差系数，根据地层温度和压力计算得到[40]：

式中Z表示偏差因子，无量纲；pcc表示临界压力，MPa；Tc表示临界温度，K；pr表示对比压力，无量纲；Tr表示对比温度，无量纲；p表示地层压力，MPa；T表示地层温度，K。

水的密度随地层温度和压力变化较小，可以通过下式计算得到[41]：

式中ρw,0表示地层压力p0=1.1 MPa和地层温度T0= 0 ℃条件下水的密度，取999.8 kg/m3；需要注意的是，T1表示地层的摄氏温度（T1=T－273.15）；T表示地层的开尔文温度，K；E表示水的弹性模量，MPa；β表示水的体积温度膨胀系数，此处取2.04×10－8℃[42]。

2.3 流体黏度

流体黏滞力作为成藏充注过程中的阻力之一，同样会影响流体分布与最终的含气饱和度。流体黏滞力可以用流体黏度来表征充注过程中黏滞力的作用。本文通过解析计算方法，得到不同地层温度和压力条件下气水黏度。

天然气黏度计算：

式中μg、μg分别表示不同地层温度和压力条件下和理想状况下的天然气黏度，mPa·s；A1、A2和A3表示计算常数，其中A1=7.9，A2=9×10－6和A3=0.2[41]。

水的黏度计算：

式中μwr表示不同地层温度和压力条件下水的黏度，mPa·s；ΔT表示环境温度，℃，由ΔT=T－273 K计算得[11]。

2.4 地层温度和压力的影响

成藏期间地层温度与压力均存在差异，本文以鄂尔多斯盆地苏里格地区地层为例，探讨地层温度与压力对物性参数和界面作用参数的影响。苏里格地区上古生界致密砂岩气成藏期地层压力约65 MPa，目前地层压力约23 MPa；流体包裹体测温成藏期最高地层温度为180 ℃（453.15 K），平均地层温度为110 ℃（388.15 K），目前地层温度约为83 ℃（356.16 K）。

2.4.1 流体黏度变化特征

如图2-a所示，当压力不变的情况下，天然气的黏度随温度的升高整体呈降低趋势。当环境压力p=65 MPa时，在温度升至353 K以前，随着温度升高，天然气黏度下降较快；而当温度超过353 K后，随着温度升高，天然气黏度下降速度逐渐减慢。且压力越高，随着温度的升高，天然气黏度下降的越快。如图2-b所示，在温度不变的情况下，天然气的黏度随着压力的升高而增大，且几乎呈线性关系。温度越低，天然气黏度随压力升高而增大的速度越快。

图2 天然气黏度随温度和压力的变化图

图3-a为不同地层压力条件下，水的黏度随地层温度变化的特征，即随着地层温度的升高，水的黏度持续降低。在温度升高初期，水的黏度迅速降低，当温度超过333 K后，水的黏度依旧随着温度的升高而降低，但降低速度逐渐放缓。压力对水的黏度影响较小，在相同的温度条件下，随着压力的升高，水的黏度几乎不发生改变（图3-b）。因此水的黏度主要受温度的影响，而对压力的变化不敏感，这是由于水分子之间形成的氢键网络，其密度对温度异常敏感，而对压力不敏感造成的。

图3 水的黏度随温度和压力的变化图

2.4.2 气水界面张力变化特征

图4-a为不同压力条件下，气水界面张力随温度变化特征，即在一定压力条件下，气水界面张力随着温度的升高而降低。温度越高，气水界面张力随温度升高而降低的速度越慢。图4-b为不同温度条件下，气水界面张力随压力变化特征，即在一定温度条件下，气水界面张力随压力的升高而降低。压力越高，气水界面张力随温度升高而降低的速度越慢。

图4 气水界面张力随温度和压力变化图

2.4.3 接触角（润湿性）变化特征

通常用接触角来表征岩石润湿性的变化特征，此次以实验测得的鄂尔多斯盆地苏里格地区致密砂岩气水润湿角为基础。图5-a为在一定压力条件下，接触角随温度变化特征，即随着温度升高，接触角逐渐增大，表明润湿性逐渐由强水湿向弱水湿方向变化。图5-b为在一定温度条件下，接触角随压力变化特征，即随着压力升高，接触角逐渐增加，且在80 MPa以前接触角增加较快，当压力超过80 MPa以后，接触角增加的速度逐渐放缓。

图5 接触角（润湿性）随温度和压力变化特征图

3 孔隙网络数值模拟

以鄂尔多斯盆地苏里格地区山西组2段致密砂岩为例，取样深度3 292.3 m，具体物性如表1所示。本文基于最大球法[34-35]生成了致密砂岩孔隙网络模型，进行了致密砂岩储层天然气孔隙尺度上的充注模拟，探讨了天然气充注机制，分析了充注过程中气水分布特征与形成过程，并基于孔隙流体系统中的界面作用和物性特征，讨论了不同地质条件下充注过程中含气性差异。

表1 实际岩心与孔隙网络参数表

3.1 孔隙网络特征

由于微米CT扫描精度有限，为了更好地模拟天然气充注过程，剔除了孔隙网络模型中不连通孔隙[43-44]。图6-a为孔隙网络模型，孔隙网络中圆球表示孔隙，其直径表示孔隙大小；圆柱表示喉道，其长度及直径表示喉道长短与大小。图6-b为鄂尔多斯盆地山西组致密砂岩岩心压汞实验测得的毛细管压力曲线与孔隙网络模型模拟，获得的毛细管压力曲线对比图，两条毛细管压力曲线具有较高的相似度，表明本文采用的孔隙网络模型与实际情况吻合。

图6 修正的孔隙网络模型与孔隙结构特征图

3.2 天然气充注机制

在致密砂岩气充注过程中，主要动力为源—储压差，主要阻力为毛细管压力，图7展示了两种可能的毛细管力变化，图7中气体（红色）均向上充注，图7-a中气体由大孔向小孔充注，其阻力f为：

由式（1）Young-Laplace公式可得：

式中f表示阻力，N；pc1、pc2分别表示小孔和大孔的毛细管压力，MPa；σ表示界面张力，mN/m；θ表示接触角，（°）。

阻力的大小取决于r1与r2之间的差值，若充注压力大于阻力则气体可以充注进入小孔隙中，若充注压力不足以克服毛细管阻力，则小孔隙无法被充注，表现为高含水。图7-b为气体从小孔向大孔运移，由上述分析可知，此时气体易于进入大孔中，气体更容易在大孔中聚集。

图7 两种充注过程中毛细管内的受力分析图

在致密砂岩孔隙网络天然气充注模拟过程中，天然气同样需要克服孔隙内的毛细管压力才能充注进入储层孔隙中[45]。在初始条件下，孔隙网络模型首先饱和水，然后天然气从底部充注进入孔隙网络模型。

图8为充注过程中不同含气饱和度下的气体分布特征。本模拟试验实现了结果的可视化表征，为了更清晰地观察气体的充注过程与分布特征，图8中不显示未被天然气充注的孔隙（即仍饱含水的孔隙），并用红色表示已经被气体充注的孔隙和喉道，同时保留了孔隙网络的几何结构特征（圆球表示孔隙，其直径表示孔隙大小；圆柱表示喉道，其长度与直径表示喉道的长短与宽窄）。

在气体充注初期，含气饱和度为10%时（图8-a），气体首先由连通的较大孔隙充注进入孔隙网络中，并产生充注前缘（图8-a中A区域）；底部的部分孔隙虽然也被气体充注，但是由于与其相连接的喉道尺度较小，毛细管阻力大，因此气体充注进入后无法立即向前推进（图8-a中B区域）。随着充注的进行，更多的气体进入了孔隙网络中，含气饱和度增加，此时图8-a中B区域的部分气体继续向前推进（图8-b中的B区域的蓝色箭头），并与图8-b中C区域的部分气体相连通，此时气体通过较大孔隙与喉道组成的高渗通道向前推进，形成新的充注前缘（图8-b中的D区域），且充注过程中天然气的流型呈现出“指进”的特征。随着充注的继续进行，含气饱和度进一步升高，当含气饱和度达到50 %时，此时较大孔隙基本被完全充注（图8-c）。值得注意的是，并非所有的大孔隙都会首先被充注，这是由于孔隙网络的非均质性特征，与小孔隙或者窄喉道相连接的大孔隙，其只有在充注后期，与其相连接的小孔隙或者窄喉道被充注后，它才能被气体充注（图8-d），甚至始终无法被充注。

图8 天然气充注过程中模拟的气体分布特征图

3.3 充注过程中微观气水分布特征

图9为含气饱和度10%、20%、50%与70%情况下的孔隙网络内气水分布特征，其中红色表示已经被气体充注的孔隙或喉道，蓝色表示未被气体充注仍饱含水的孔隙或喉道。图9完整地展示了充注过程中的气水分布变化特征。当充注开始时，气体首先从底部进入较大的孔隙中（图9-a）；随着充注的进行，气体逐渐占据较大的孔隙（图9-b）与中等大小的孔隙（图9-c），连通的中大孔隙内水首先被驱替出；充注结束后，残余水主要赋存在小孔隙与窄喉道中，而少量的残余水赋存在与小孔隙或窄喉道相连接的大孔隙中（图9-d）。值得注意的是，少量水也残留在与小孔隙或窄喉道相连接的大孔隙中，这部分大孔隙无法被天然气充注。

图9 天然气充注过程中模拟的气水分布特征图

3.4 不同地质条件下的充注特征

根据前文分析，温度与压力等地层条件变化在微观上会影响孔隙内流体分子之间、流体分子与孔隙壁面之间的分子间相互作用力，从而在宏观上影响毛细管压力的大小。本文根据典型致密砂岩气藏成藏特征，选取了三个地层温度和地层压力条件（表2），进行了天然气充注的孔隙网络模拟。表2展示了不同地层温度和压力条件下的孔隙网络模拟结果。从条件1到条件3，地层温度和地层压力逐渐升高，孔隙中流体的界面张力也随之降低，接触角逐渐增大；由式（1）可知，此时天然气充注进入孔隙所需要克服的毛细管压力也随之降低。在相同的充注压力下，气体可以充注进入更小的孔隙中，进而提高充注后的含气饱和度。

表2 不同温度压力条件下的孔隙网络模拟结果表

图10为天然气充注进入孔隙网络过程中需要克服的毛细管压力与含水饱和度之间的关系。天然气充注进入孔隙网络中后，含气饱和度逐渐升高，含水饱和度逐渐降低。气体首先充注进入相互连通的中大孔隙，随着气体的逐渐充注，连通的大孔隙基本都被完全充注，气体将逐渐进入更小的孔隙中，含水饱和度进一步降低，此时充注过程需要克服更大的毛细管压力。

如图10所示，随着含气饱和度的升高，充注过程所需要克服的毛细管压力逐渐增大。图10中条件3（蓝色曲线）相比于条件1（红色曲线），具有更高的地层温度和地层压力，因此在条件3情况下，气体充注进入相同尺度的孔隙或喉道中，所需要克服的毛细管压力要小于条件1情况，直观的表现为条件3的毛细管压力曲线位于条件1的毛细管压力曲线的下方。模拟过程中设置最大充注压力为10 MPa，在充注结束后，条件1情况下最终含气饱和度为62.9%，条件3含气饱和度为67.4%，含气饱和度提高了4.5%。

图10 天然气充注过程中需要克服的毛细管压力与含气饱和度关系图

地质条件的改变将在微观上导致孔隙内流体分子之间，流体分子与孔隙壁面之间的分子间相互作用力的改变，在宏观上表现为孔隙内毛细管压力的变化，这种变化导致了不同地层温度和压力储层中天然气充注下限存在差异，最终造成储层含气性的差异。

4 应用前景

基于孔隙网络的致密气充注数值模拟方法，拓展了国内针对成藏数值模拟手段，揭示了连通孔隙中的自由水的成因机制。与此同时，孔隙网络模拟技术同样可以充分利用数字岩心资源，将CT扫描技术与数值模拟方法相结合，进一步探究成藏过程中的微观现象，揭示微观机理，更好地进行含气性与气水分布的预测。

在致密砂岩气开发过程中往往面临气水分布复杂的难题，孔隙网络模拟技术有助于深化对致密砂岩气藏成藏机理的认识，深入理解复杂气水分布关系的形成过程与微观渗流特征，有助于提高对致密砂岩气藏的生产机理和潜在生产能力的认识。本文中使用的岩心取自苏里格地区某井的储集层段，根据主充注期的成藏压力与温度环境下模拟的含气饱和度为65.2%，测井解释含气饱和度为68%，两者具有较好的吻合程度。

诚然，主成藏期的充注过程只是成藏过程的一个阶段，储层的含气性在后期地质变化过程中均有可能发生改变。但基于孔隙网络的致密砂岩气充注数值模拟方法拓展了研究人员们对数字岩心和CT扫描技术的应用，丰富了对充注过程的模拟手段，有利于深入理解致密油气，以至于页岩油气的成藏机理。