en
×

分享给微信好友或者朋友圈

使用微信“扫一扫”功能。
作者简介:

严侠(1989-),男,博士,博士后,研究方向为非常规油气藏流固耦合数值模拟。E-mail:jsyanxia1989@163.com。

通讯作者:

姚军(1964-),男,教授,博士,博士生导师,山东省泰山学者攀登计划专家,研究方向为油气田开发工程。E-mail:RCOGFR@126.com。

中图分类号:TE319

文献标识码:A

文章编号:1673-5005(2021)01-0101-09

DOI:10.3969/j.issn.1673-5005.2021.01.012

目录contents

    摘要

    页岩基质通常由微尺度的有机质和无机质组成,由于这两种介质的力学性质差异较大,且气体在两者中的赋存方式和流动机制均不相同,要准确刻画气体在页岩基质中的流固耦合过程,需针对两种介质建立不同的微尺度模型,但由于计算效率问题,微尺度模型无法直接用于宏观模拟。 基于均匀化理论建立一种页岩基质流固耦合尺度升级方法,将微尺度上有机质和无机质各自特征有效地表征到宏观模拟中。 首先,将页岩基质看成是由有机质和无机质组成的非均质多孔弹性介质,考虑真实气体在这两种介质中不同的赋存方式和流动机制,建立微尺度流固耦合模型;然后,采用均匀化理论进行尺度升级,推导得到页岩基质的宏观等效流固耦合模型,并给出相关等效参数的定义和计算方式;最后,通过数值算例验证方法的正确性,并分析有机质力学性质、体积分数和分布对页岩气藏宏观流固耦合数值模拟的影响。 结果表明:有机质弹性模量越小、无机质连通性越好,累积产气量越高,而有机质体积分数越高,早期的累积产气量越低,后期的累积产气量越高。

    Abstract

    In shale gas reservoirs, gas is mainly stored in shale matrix. Therefore, to establish a hydro-mechanical coupling (HM) model in shale is important for the macroscopic HM simulation of shale gas production. At microscopic scale, shale matrix is usually composed of organic and inorganic matters. The mechanical properties of these two media are quite different, and so that both gas storage type and transport mechanisms are different in the two media. Thus it is necessary to develop different microscale models to accurately describe the HM process of gas in shale matrix. However, the microscale models cannot be straightly applied to macroscale simulation due to huge calculation time and data storage requirement. In this study, an efficient upscaling method based on the homogenization theory was proposed for the HM process in shale matrix, which can accurately represent the microscale characteristics of the organic and inorganic matters in macroscale simulations. Firstly, the shale matrix is assumed as a heterogeneous porous elastic medium composed of organic and inorganic matters, and the microscale HM models were developed according to different storage type and transport mechanism of real gas in these two media. Then, the microscale HM models were homogenized to obtain the equivalent macroscopic HM models for the shale matrix, along with the definition and calculation method of equivalent parameters in the micro and macro models. Lastly, the accuracy and application of the proposed method were verified via several case studies, and the influence of the mechanical properties, the content and distribution of organic matters on the macro HM simulation in shale gas production were also analyzed. The results indicate that the cumulative gas production increases with the decrease of the elastic modulus of the organic matter, and it correlates positively with the connectivity of the inorganic matter. Besides, at early period, gas production decreases with the increase of the total organic content, but it will increase at late period of the process.

  • 页岩气资源分布广、储量大,但其基质渗透率极低,通常需要经过水力压裂才能进行商业化开采[1-5], 压裂后的页岩储层具有多尺度孔缝特征[6-7],且处于复杂地应力场和渗流场的共同作用中,流固耦合效应明显[8-10]。其中页岩基质作为气体的主要储集空间,其流固耦合模型的建立是页岩气藏宏观流固耦合数值模拟的基础。由于页岩基质通常由微尺度的有机质和无机质组成[11-13],而这两种介质的力学性质差异较大[14-15],且气体在两者中的赋存方式和流动机制均不相同[12,16],要准确刻画气体在页岩基质中的流固耦合过程,需针对两种介质建立不同的微尺度模型,但由于计算效率问题,微尺度模型无法直接用于宏观模拟。因此有必要建立一种尺度升级方法,将微尺度上有机质和无机质各自特征有效地表征到页岩气藏宏观流固耦合模拟中。近年来,国内外学者分别采用流量等效法[17-20] 和均匀化理论[21-24],对页岩基质中有机质分布对气体宏观渗流规律的影响进行了研究[25-26],但均未考虑流固耦合作用。对此,笔者基于均匀化理论,建立考虑微尺度上有机质和无机质各自特征的页岩基质流固耦合尺度升级方法,通过数值算例验证尺度升级方法的正确性,并分析有机质力学性质、体积分数和分布对页岩气藏宏观流固耦合数值模拟的影响。

  • 1 微观尺度流固耦合模型建立

  • 扫描电子显微镜( scanning electron microscope, SEM)扫描研究发现,在微观尺度上,有机质分散嵌入在无机质中[17]。假设有机质和无机质均满足连续介质假设[11-13],将页岩基质看成是由有机质和无机质组成的非均质多孔弹性介质[14],其中无机质中仅存在游离气,气体的运移机制主要为黏性流和Knudsen扩散,有机质中同时存在吸附气和游离气, 气体的运移机制主要为黏性流、Knudsen扩散和表面扩散[16]

  • 为了便于书写,采用以下广义的流固耦合模型进行描述[27],表达式为

  • σ=0,σ=σs-αpI,σs=Ceus,
    (1)
  • βpt+δ(1-φ)mgt+αρgust+ρgv=0
    (2)
  • v=-kaμp
    (3)
  • 其中

  • eus=12us+usTβ=α-φKsρg+φMgRT1Z-pZZ2p,ρg=MgpZRT

  • 式中,σσs 分别为总应力张量和有效应力张量;p 为孔隙压力;α 为Biot系数;I 为单位张量;C 为弹性张量;u s 为骨架位移;e(us)为应变张量;β 为综合压缩系数;t 为时间;δ 为克罗内克符号;φ 为孔隙度; mg 为吸附/解吸项;ρg 为气体密度;v 为气体流速;k a 为气体的视渗透率;μ 为气体黏度;K s 为基质骨架体积模量;M g 为气体摩尔质量;Z 为气体压缩因子;R 为通用气体常数;T 为储层温度。

  • 根据Langmuir等温吸附[28],气体吸附/解吸项mg 的表达式[6]

  • mg=ρrρgstdVLp/ZpL+p/Z
    (4)
  • 式中,ρr 为岩石密度;ρgstd 为标准条件下的气体密度;V Lp L 分别为有机质的Langmuir体积和Langmuir压力。

  • 无机质中气体的运移机制主要为黏性流和Knudsen扩散,其视渗透率k ai [12,16]

  • kai=k1+12815π2tan-14.0Kn0.4Kn1+4Kn1-bKn
    (5)
  • 式中,K n 为knudsen数。

  • 有机质中气体的运移机制主要为黏性流、Knudsen扩散和表面扩散,其视渗透率k ak [12,16]

  • kak=k1+12815π2tan-14.0Kn0.4Kn1+4Kn1-bKn+μMgDsCmaxdθdp1-reffrh2
    (6)
  • 其中

  • k=rh2φ8τ,Kn=λreff ,λ=πZRT2MgμpCmax=ρrVLVstd,θ=p/ZpL+p/Z

  • 式中, rh τ 分别为孔隙半径和迂曲度;λ 为气体分子的平均自由程;reff 为孔隙有效半径;b=-1 为滑移系数;D s 为吸附气的表面扩散系数;Cmaxθ 分别为有机质的最大吸附量和气体覆盖率[12,16]

  • 以上页岩基质流固耦合模型的非均质系数在无机质中取值为

  • α=αi,C=Ci,φ=φi,δ=0,ka=kairh=rhi,τ=τi,reff=rhi.

  • 在有机质中取值为

  • α=αk,C=Ck,φ=φk,δ=1,ka=kak,rh=rhkτ=τk,reff=rhk-dmp/ZpL+p/Z

  • 式中,下标i和k分别表示无机质和有机质;d m 为甲烷分子直径。

  • 2 基于均匀化理论的尺度升级

  • 在均匀化理论中,物理量特征函数依赖于两种不同尺度:一种是宏观尺度,用来表征介质特征函数在大尺度上的平缓变化;另一种是微观尺度,用来表征特征函数在小尺度上的快速振荡变化[29]。考虑一个特征长度为 L 的页岩基质,其由大量的周期性微尺度元胞组成,微尺度元胞的特征长度为 l,引入不同特征长度之间的比值 ε=l/L。上述页岩基质流固耦合模型中的变量 σusvp 可由以下渐进展开式近似:

  • φ=φ(0)(x,y,t)+εφ(1)(x,y,t)+ε2φ(2)(x,y,t)+
    (7)
  • 式中,y=x/ε,y 为微观尺度坐标;x为宏观尺度坐标;φ(i)为关于坐标x、y以及时间 t 的周期性变量。

  • 此外,根据链式准则,可以得到:

  • x+ε-1y,2x2+2ε-1yx+ε-2y2
    (8)
  • 把方程(7)代入方程(1)~(3)中,并采用式(8) 进行变换,比较方程中不同阶数 ε 的系数。

  • 通过比较方程(3) 中 ε-1 的系数以及方程(1) 中 ε-2 的系数可得

  • yp(0)=0,yCeyus(0)=0
    (9)
  • 由此可知 p(0)us(0) 仅在宏观尺度x上变化, 与微观尺度y无关,即

  • p(0)=p(0)(x,t),us(0)=us(0)(x,t)
    (10)
  • 通过比较方程(2) 中 ε-1 的系数以及方程(3) 中 ε 0 的系数可得

  • yρgp(0)v(0)=0μp(0)kap(0)v(0)+yp(1)=-xp(0)
    (11)
  • 分别对 v(0)p(1)进行分离变量分析:

  • v(0)=-1μp(0)ω(y)xp(0),p(1)=-χ(y)xp(0)
    (12)
  • 式中,ωχ 为关于 y 的周期性变量。

  • 把式(12)代入式(11),可以得到元胞辅助方程为

  • yωi=01kap(0)ωi+χi=ei
    (13)
  • 式中,ei 为笛卡尔坐标系中的i-方向的单位矢量。

  • 对方程(12)中的 v(0) 在元胞上进行体积平均, 可以得到宏观渗流方程为

  • v(0)=-1μp(0)1|Ω|Ω ωp(0)dΩxP(0)=-kequp(0)μp(0)xp(0)
    (14)
  • 式中, |Ω| 为整个元胞的体积。

  • 〈a〉表示变量a( a指本文中涉及到的变量)的体积平均,k equ为等效视渗透率张量,表达式分别为

  • a=1|Ω|Ω adΩ
    (15)
  • kequijp(0)=1|Ω|Ω ωip(0)jdΩ
    (16)
  • 比较方程(1)中 ε -1 的系数可得

  • yCexus(0)+eyus(1)=0
    (17)
  • 由此可知, us(1)可以写成:

  • us(1)=ξ(y)exus(0)+u¯s(1)(x,t)
    (18)
  • 式中,u¯s(1)(x,t)为与坐标y无关的任意变量。

  • 与坐标y相关的周期性变量 ξpq 满足:

  • yCijlle~kl+eyklξpq=0,e~k=12δkpδlq+δkqδlp
    (19)
  • 比较方程(1)中 ε 0 的系数可得

  • xCexus(0)+eyus(1)-αp(0)I+yCexus(1)+eyus(2)-αp(1)I=0
    (20)
  • 将方程(20)在元胞上进行体积平均,并采用散度定理和周期性边界条件,得到

  • xCexus(0)+eyus(1)-αp(0)I=0
    (21)
  • 根据体积平均定理和周期性边界条件,对方程(21)进行整理,可以得到宏观应力平衡方程为

  • xσ(0)=0,σ(0)=Cequ exus(0)-αequ p(0)I.
    (22)
  • 其中

  • Cequijkl =Cijkl+Cijmneymnξkl,αequ =α

  • 比较方程(2)中 ε 0 的系数可得

  • βp(0)p(0)t+δ(1-φ)mgp(0)t+αρgp(0)exiius(0)+eyiius(1)t+xρgp(0)v(0)+yρgp(0)v(1)=0
    (23)
  • 将方程(23)在元胞上进行体积平均,并采用散度定理、体积平均定理和周期性边界条件,可以得到宏观质量守恒方程为

  • βequ p(0)p(0)t+mgequ p(0)t+αequ ρgp(0)×xus(0)t+xρgp(0)v(0)=0
    (24)
  • 其中

  • βequ =ρgp(0)α-φKs+φMgRT1Zp(0)-p(0)Zp(0)Zp(0)2p(0)mgequ =aTOC 1-φkρrρgstdVLp(0)/Zp(0)pL+p(0)/Zp(0)

  • 式中,aTOC 为元胞内有机质的体积分数。

  • 最后,采用下标hm表示下标(0),将宏观等效连续介质流固耦合模型表示为

  • σhm=0,σhm=Cequeuhm-αequphmI,
    (25)
  • βequphmt+mgequt+αequρguhmt+ρgvhm=0
    (26)
  • vhm=-kequμphm
    (27)
  • 可以看出宏观等效连续介质流固耦合模型(25)~(27)的形式与方程(1)~(3)一致,但是由于同时考虑了有机质和无机质的存在,模型的系数发生了变化。元胞辅助方程(13)和(19)很难求得解析解,本文中采用有限单元法[10, 30]进行数值求解。

  • 3 算例分析

  • 3.1 模型验证

  • 图1 为一理想页岩基质模型,尺寸为2 mm×2 mm,内部均匀分布100 个半径为0.05 mm的圆形有机质。采用平面应变模型[30],力学边界条件如图1(a)所示,上边界为流体压力0.1 MPa的定压边界,其余边界为封闭边界,初始流体压力为20 MPa。模型参数如下:有机质的孔隙度、弹性模量和孔隙半径分别为0.004、5.0 GPa和2.51 nm;无机质的孔隙度、弹性模量和孔隙半径分别为0.04、40.0 GPa和25.1 nm;有机质表面扩散系数为1.0×10-8 m 2/s; Langmuir体积和标准摩尔体积分别为0.016 kg/mol、0.41 nm、13.79 MPa、0.028 63 m 3/kg和0.024 m 3/mol;通用气体常数为8.314 J/(mol·K);初始地层压力和温度分别为20.0 MPa和323.0 K。分别基于细网格和粗网格采用有限单元法对该模型进行模拟,网格划分结果如图1( b) 和( c) 所示。其中,细网格对有机质进行精细刻画,其模拟结果视为参考解,每个粗网格单元的等效参数通过提出的尺度升级方法计算得到,其中通过采用有限单元法求解方程(19)并结合方程(23),可得等效弹性张量为36.314.0014.036.30009.8 GPa ,等效视渗透率随压力的变化如图2 所示。图3 为细网格和粗网格两者模拟1 s后的压力场和y方向位移场对比,可以看出两者结果基本一致。

  • 图1 有机质均匀分布模型示意图及不同网格划分结果

  • Fig.1 Schematic of organic matter uniform distribution model and results of different grids

  • 图2 粗网格单元等效视渗透率随压力的变化

  • Fig.2 Relationship between equivalent apparent permeability and pressure for coarse grids

  • 为了进一步验证提出的尺度升级方法的正确性,分别采用细网格和粗网格对有机质半径分别为0.02、0.04 和0.05 mm(对应的有机质体积分数分别为0.031、0.126 和0.196)的情况进行模拟,对比3 种情况下两种方法计算得到的压力观测点(1.0 mm, 0 mm)处的压力和位移观测点(1.0 mm, 2.0 mm)处的y方向位移。压力和位移的对比结果如图4 所示,可以看出两种方法计算得到的结果基本一致,验证了提出的尺度升级方法的正确性。此外,还可以看出随着有机质体积分数的增大,观测点压力下降速度减小,y方向位移增大,这是因为有机质的视渗透率和弹性模量均小于无机质,有机质体积分数增大导致基质整体的渗透性和弹性刚度降低,同时有机质体积分数增大还能增加气体吸附量。

  • 图3 细网格和粗网格模拟1 s后的压力场和位移场对比

  • Fig.3 Comparisons of pressure and displacement fields for fine grid and coarse grid after 1 s

  • 图4 不同有机质体积分数下两种方法计算得到的压力和位移对比

  • Fig.4 Comparisons of pressure and displacement of two different methods for various aTOC

  • 3.2 实例计算及影响因素

  • 3.2.1 有机质力学性质对气体生产的影响

  • 如图5 为某一100 m×100 m的页岩气藏模型及其内部有机质分布,假设其基质的表征单元体为1 mm×1 mm,有机质体积分数为0.138。气藏模型满足平面应变假设,力学边界条件如图5 所示,流动边界为封闭边界,初始地层压力为20 MPa,生产井位于气藏模型左下角,井底流压为5 MPa,其余模型参数与3.1 节算例一致。

  • 图5 页岩气藏模型及其基质表征单元体

  • Fig.5 Schematic of shale gas reservoir and its REV element

  • 为研究有机质力学性质对气体生产的影响,将有机质弹性模量分别设置为1、10 和40 GPa。采用提出的尺度升级方法计算表征单元体的等效参数, 等效视渗透率如图6 所示。等效弹性张量分别为33.5413.110.13813.1135.010.1560.1380.15610.3444.2118.630.03718.6344.530.0280.0370.02812.7753.8523.08023.0853.8500015.38 GPa 。可以看出,对于图5 所示的基质表征单元体,有机质导致的基质等效视渗透率各向异性较弱,且等效视渗透率张量中的非主对角线值远小于主对角线值。此外,还可以看出有机质弹性模量对基质表征单元体的等效弹性张量有较大影响。

  • 图6 基质表征单元体等效视渗透率随压力的变化

  • Fig.6 Relationship between equivalent apparent permeability and pressure for REV element

  • 根据计算得到的表征单元体等效参数,采用有限单元法对图5 中的页岩气藏进行宏观流固耦合数值模拟,模拟过程中采用以下动态渗透率模型考虑岩石变形对渗透率的影响[31] :

  • k=k01+αφ0Δεv+ΔpKs3
    (28)
  • 式中,下角标0 表示初始状态; Δεv 和 Δp 分别为体积应变和气体压力变化量。

  • 图7 为有机质弹性模量取不同值时计算得到的累积产气量对比,可以看出有机质弹性模量越小累积产气量越高。这是由于有机质弹性模量降低,导致页岩基质的可压性增强,弹性开采过程中,相同压降导致的基质孔隙度减小越多,因此累积产气量越高。此外,从图7 中还可以看出页岩基质的应力敏感程度较弱,即岩石变形导致渗透率降低进而降低累积产气量的效果不明显,这也是有机质弹性模量降低导致累积产气量增高的另一原因。

  • 图7 不同情况下页岩气藏累积产气量对比

  • Fig.7 Comparisons of cumulative production for different cases

  • 3.2.2 有机质体积分数对气体生产的影响

  • 图8 为3 个不同有机质体积分数的基质表征单元体示意图,尺寸均为1 mm×1 mm,其余模型参数与3.1 节中算例一致。采用提出的尺度升级方法计算各个表征单元体的等效参数,其中x方向等效视渗透率如图9 所示。等效弹性张量分别为51.1721.77-0.00521.7751.18-0.012-0.005-0.01214.67945.7519.09-0.0619.0945.51-0.06-0.06-0.0613.2340.13 16.57 0.07016.57 40.78 0.0620.070 0.062 11.78 GP 。可以看出,有机质体积分数越高,表征单元体的等效视渗透率和等效弹性张量越小。

  • 根据计算得到的各个表征单元体的等效参数, 考虑岩石变形对渗透率的影响,对图5 中的页岩气藏进行宏观流固耦合数值模拟。图10 为不同有机质体积分数情况下计算得到的累积产气量对比。可以看出,有机质体积分数越高,早期的累积产气量越低,而后期的累积产气量越高。这是由于有机质体积分数增大虽然增加了储层中气体储量,但也导致了基质等效视渗透率减小,以及基质的应力敏感程度增大。由于早期累积产气量受视渗透率的影响较大,因此有机质体积分数越高早期累积产气量越低,而后期累积产气量主要受气体储量的影响,因此有机质体积分数越高后期累积产气量越高。

  • 图8 不同有机质体积分数的基质表征单元体示意图

  • Fig.8 Schematic of REV elements with different aTOC

  • 图9 不同有机质体积分数的基质表征单元体等效视渗透率随压力的变化

  • Fig.9 Relationships between equivalent apparent permeability and pressure for REV elements with different aTOC

  • 3.2.3 有机质分布对气体生产的影响

  • 图11 为3 个不同有机质分布的基质表征单元体示意图,尺寸均为1 mm×1 mm,有机质体积分数为0.138,模型参数与3.1 节中算例一致。采用提出的尺度升级方法计算各个表征单元体的等效参数, 其中等效弹性张量分别为33.4611.27011.2733.460008.0343.6811.74011.7427.390007.8340.1316.570.07016.5740.780.0620.0700.06211.78GPa,等效视渗透率如图12 所示。可以看出,有机质分布对基质表征单元体的等效视渗透率和等效弹性张量均有较大影响,是导致页岩基质渗流和力学性质各向异性的重要因素。

  • 图10 不同有机质体积分数下页岩气藏累积产气量对比

  • Fig.10 Comparisons of cumulative production for shale gas reservoir with different aTOC

  • 图11 有机质分布不同的基质表征单元体示意图

  • Fig.11 Schematic of REV elements with different organic matter distribution

  • 同样,根据计算得到的3 个表征单元体的等效参数,考虑岩石变形对渗透率的影响,对图5 中的页岩气藏进行宏观流固耦合数值模拟。图13 为不同有机质分布情况下计算得到的累积产气量对比,结合3 种分布情况下的无机质连通情况,可以看出,无机质连通性越差,累积产气量越低。

  • 图12 不同有机质分布的基质表征单元体等效视渗透率随压力的变化

  • Fig.12 Relationships between equivalent apparent permeability and pressure for REV elements with different organic matter distribution

  • 图13 不同有机质分布下页岩气藏累积产气量对比

  • Fig.13 Comparisons of cumulative production for shale gas reservoir with different organic matter distribution

  • 4 结论

  • (1)采用均匀化理论可以对包含微尺度有机质和无机质的页岩基质进行流固耦合尺度升级。

  • (2)有机质弹性模量越小、无机质连通性越好, 累积产气量越高,而有机质体积分数越高,早期的累积产气量越低,后期的累积产气量越高。

  • 参考文献

  • 参考文献

  • 版权所有 中国石油大学学报(自然科学版)编辑部 Copyright©2008 All Rights Reserved
    主管单位:中华人民共和国教育部 主办单位:中国石油大学(华东)
    地址: 青岛市黄岛区长江西路66号中国石油大学期刊社 邮编:266580 电话:0532-86983553 E-mail: journal@upc.edu.cn
    本系统由:北京勤云科技发展有限公司设计