Loading [MathJax]/jax/output/SVG/fonts/TeX/Main/Italic/Main.js
  • 全国中文核心期刊
  • 中国科技核心期刊
  • 美国工程索引(EI)收录期刊
  • Scopus数据库收录期刊

基于岩土介质三维孔隙结构的两相流模型

张鹏伟, 胡黎明, MeegodaJay N, CeliaMichael A

张鹏伟, 胡黎明, MeegodaJay N, CeliaMichael A. 基于岩土介质三维孔隙结构的两相流模型[J]. 岩土工程学报, 2020, 42(1): 37-45. DOI: 10.11779/CJGE202001004
引用本文: 张鹏伟, 胡黎明, MeegodaJay N, CeliaMichael A. 基于岩土介质三维孔隙结构的两相流模型[J]. 岩土工程学报, 2020, 42(1): 37-45. DOI: 10.11779/CJGE202001004
ZHANG Peng-wei, HU Li-ming, Meegoda Jay N, Celia Michael A. Two-phase flow model based on 3D pore structure of geomaterials[J]. Chinese Journal of Geotechnical Engineering, 2020, 42(1): 37-45. DOI: 10.11779/CJGE202001004
Citation: ZHANG Peng-wei, HU Li-ming, Meegoda Jay N, Celia Michael A. Two-phase flow model based on 3D pore structure of geomaterials[J]. Chinese Journal of Geotechnical Engineering, 2020, 42(1): 37-45. DOI: 10.11779/CJGE202001004

基于岩土介质三维孔隙结构的两相流模型  English Version

基金项目: 

清华大学自主科研计划项目 THZ-2016-02

中央高校基本科研业务费专项资金项目 2019RC054

详细信息
    作者简介:

    张鹏伟(1990— ),男,工学博士,主要从事岩土介质孔隙尺度多相流、多组分运移数值计算,以及油气藏、水合物三维、四维地质力学计算等方面的研究工作。E-mail:zpw12@tsinghua.org.cn

    通讯作者:

    胡黎明, E-mail:gehu@tsinghua.edu.cn

  • 中图分类号: TU43

Two-phase flow model based on 3D pore structure of geomaterials

  • 摘要: 地下储层中岩土介质一般具有较低的孔隙连通性,宏观流动模拟一般忽略微观尺度的孔隙连通性,通过渗透率、弯曲度等参数反映储层的整体特性。但岩土介质的多孔性及孔隙间复杂的连通性,使得宏观描述流体在岩土介质中流动不能反映其内在流动特征。孔隙结构模型的建立可以反映岩土介质中孔隙的几何形态及空间连通性,为解释流体在复杂多孔介质中的流动特性提供有效手段。通过考虑岩土介质孔隙尺寸分布、孔隙孔喉空间相关性、孔隙连通性等特征参数,建立了反映不同岩土介质连通性、各向异性特征的等效孔隙网络模型。等效孔隙网络模型通过水力特征参数等效的方式反映岩土介质三维微观孔隙结构,通过渗透率计算验证了模型的有效性。此外,基于建立的孔隙结构模型,开发了孔隙尺度动态两相流计算模型,模型可以反映孔隙内弯液面的动态运动过程,直观反映多孔介质中的优势渗流,可以为不同孔隙尺度岩土介质提供表观渗透率、击穿曲线、相对渗透率曲线等宏观计算参数。将孔隙尺度两相流模型应用于页岩气开采中水力阻滞特性研究,结果表明:页岩基质的残余饱和度约为30%,随着平均配位数的增加,残余饱和度显著降低。
    Abstract: Geomaterials normally have low pore-connectivity in underground reservoir, and the macro-scale flow simulation normally ignores the micro pore connectivity and uses macro parameters such as permeability and tortuosity to reflect the conductivity of underground reservoir. However, due to the complex pore structure and pore connectivity of geomaterials, the macro-scale method cannot reflect the micro flow mechanisms. The pore-structure model provides an effective way to reflect the micro-flow mechanisms for complex porous media since the pore geometry and pore connectivity can be included in the model itself. In this work, an equivalent pore-network model (EPNM) is established considering pore-size distribution, spatial correlation and pore-connectivity. EPNM aims at reflecting 3D pore structure of geomaterials by the equivalent hydraulic parameters, and the effectiveness is verified by permeability tests. Furthermore, a dynamic two-phase flow model is developed based on EPNM, and simulate the dynamic invasion of each phase, reflect the preferential flow in porous media, and it can provide apparent permeability, relative permeability curve, breakthrough curve for macro-scale simulation. Finally, the dynamic two-phase flow model is applied to the wetting phase trap during shale gas exploitation. The results show that the residual saturation in shale matrix is around 30%, and this residual saturation decreases significantly with the increase of the average pore coordination number.
  • 孔隙结构模型的发展为揭示岩土介质微观孔隙结构提供了有效手段,可以反映岩土介质孔隙尺寸分布、连通性,摒弃对弯曲度等宏观概化参数的依赖;此外,还可以解释动态多相流、阻滞、指形效应(非均质)等微观流动特性。

    等效孔隙网络模型是孔隙结构模型的一种,是基于岩土介质的物理性质参数统计结果,建立反映其孔隙分布、空间连通性等的孔隙结构。等效孔隙网络模型最早由Fatt[1]提出,其基于毛细管模型提出二维网络模型,并结合球颗粒堆积模型反映孔隙间的连通性。Nicholson等[2]建立了三维等效孔隙网络模型,并进行准静态两相流计算。之后许多学者对孔隙网络模型的孔隙形态进行研究,Dullien[3]将孔隙间的连通喉道简化为连续的毛细管段。Reeves等[4]建立了恒定配位数等效孔隙网络模型,并进一步考虑孔喉尺寸的线性变化。Acharya等[5]提出BACON键模型反映孔隙间配位键断面面积的变化,模型中配位键的半径随位置变化,满足Beti’s影响线数学表达式。Raoof等[6]为了能够反映孔隙间连通的多向性及无序性,将孔隙的连通方向扩展至13个不同的方向,即每个孔隙初始情况下与26个相邻的孔隙相连,为了反映不同岩土介质的孔隙连通性,基于逾渗理论提出了稀疏算法。Gao等[7]建立了高连通性砂土的等效孔隙网络模型,有效模拟了地下水曝气修复中的气液两相流动过程。王晨晨等[8]针对碳酸岩盐建立了多尺度耦合等效孔隙网络模型。

    总结以上孔隙结构模型研究发现,等效孔隙网络模型可以有效模拟岩土介质孔隙结构,但目前对于不同岩土介质连通性、各向异性特征考虑较少。

    相对渗透率-饱和度关系曲线、基质吸力-饱和度关系曲线是研究岩土介质非饱和渗流特性的基本特征曲线,由于试验获得这些基本水力特征参数耗时较长,尤其对于低渗透岩土介质,如致密砂岩、页岩等,获取非饱和特征参数更加困难。等效孔隙网络模型可以从微观尺度模拟两相流动态运动过程,进而获得相应岩土介质的水力特征参数。

    Fatt等[1]首次运用等效孔隙网络模型来模拟两相流过程中的相对渗透率随饱和度变化关系。随后许多学者对孔隙网络模型模拟多相流问题展开研究,但大多研究工作集中在准静态孔隙网络模型[9-11],模型中忽略黏滞力的作用,两相流过程中只考虑基质吸力。Joekar-Niasar等[11]采用管道模型和球管模型分别模拟气液两相流,结果表明管道模型无法反映吸水-排水滞回特征,此外,通过球管模型获得了基质吸力-饱和度-界面面积三者之间的动态关系。

    动态孔隙网络模型由Koplik等[12]首次提出,由于动态模型可以计算相界面位置,更真实地反映实际多孔介质中非饱和流动特性,成为孔隙尺度多相流研究的热点问题[13]。Gao等[14]建立三维动态孔隙网络模型,模拟污染场地地下水曝气修复过程中气液两相流特性。Huang等[15]建立了页岩基质动态两相流模型,模型可以反映气液界面在孔隙内扩张过程,由于孔隙为立方体孔隙,对于润湿相在孔隙边角的滞留现象可以反映,其气液界面进入孔隙内类似于等界面推进。

    总结以上孔隙尺度两相流研究可知:在两相流计算模型方面的进展可以看出,对于毫米、微米级孔隙的岩土介质如砂土、砂岩等建立了较完善的动态两相流模型,但对于微纳米孔隙介质研究较少;此外,两相流模型中对于动态反映弯液面扩展研究较少。

    为此,本文在等效孔隙网络模型建立方面,通过配位键阈值的判断反映不同岩土介质连通性,并定义各向折减系数体现岩土介质孔隙结构多向性、各向异性,通过相连孔隙计算孔喉尺寸反映孔隙、孔喉空间相关性;另一方面,在两相流计算方面,本文提出弯液面动态扩展模型模拟多相动态入侵过程。

    图1所示为等效孔隙网络模型建立的流程图,等效孔隙网络模型建立需要相应岩土介质的结构信息,包括孔隙度、孔隙尺寸分布、孔喉尺寸分布、孔隙间连通性(平均配位数)。

    图  1  等效孔隙网络模型构建流程图
    Figure  1.  Flow chart of establishing EPNM

    获得岩土介质的孔隙结构信息后,根据孔隙尺寸分布、孔隙度可以计算孔隙网络模型中孔隙间距,即配位键的长度,计算的孔隙间距可确定各孔隙的空间坐标。初始孔隙配位键的设置参考Raoof等[6]提出的方法,各孔隙通过26个配位键与相邻的孔隙连接(图2),节点的拓扑编号原则依次为从下到上、从前到后、从左到右。

    图  2  孔隙连通代表性单元
    Figure  2.  Representative element of pore connectivity

    基于以上步骤可以建立每个孔隙具有26个配位键的全连通孔隙网络模型。该立方体型孔隙网络模型中具有6个边界面,考虑到实际岩土介质中渗流问题,边界面处由于与外界存在压力梯度不会出现平面流动,而是流体直接流至裂隙通道。因此模型需要对进出边界面处孔隙间的配位键进行剔除,选取单元体进行计算时,单元体6个边界面内无平面流动发生。

    为了反映实际岩土介质连通性、各向异性特征,需要降低介质内的孔隙连通性。本研究中对每个配位键进行随机判断,每个配位键设置阈值反映该键存在的概率,每个配位键的阈值(pt)与随机生成的剔除数(nr)比较,小于对应剔除数时,该连通方向的配位键剔除。配位键阈值的计算公式为

    pt=min[Nci/Ni,Ncj/Nj], (1)

    式中,Nci,Nci为由正态分布随机分配给相连孔隙i,j的配位键数,Ni,Ni为孔隙i,j连通的方向数。为了反映孔隙连通性的变化以及满足模型稳定性要求,研究中针对剔除数定义了折减因子α,α取值在0~1之间,其与孔隙的平均配位数、总孔隙数成负相关,即孔隙连通性高的岩土介质,α取值较小,则模型生成的孤立孔隙数减少;此外,α取值需要满足统计规律,即多组等效孔隙网络模型具有稳定性。

    最后,孔隙网络模型建立后,需要对流体流动的孔隙主骨架进行搜索提取,剔除孤立孔隙和孔隙串(图3)。

    图  3  三维孔隙网络模型
    Figure  3.  3D pore-network model

    图4所示为孔隙网络模型总孔隙数、配位键数统计结果,可以看出孔隙数以及配位键总数的变异系数(反映每一组数偏离平均值程度)基本在5%内,表明模型稀疏算法具有较好的稳定性。

    图  4  孔隙网络模型稳定性验证结果
    Figure  4.  Stability analysis of pore-network model

    本节主要选取砂土、砂岩以及页岩三类岩土介质展开讨论,对模型的有效性进行验证。以上三类岩土介质孔隙尺寸由毫米级别逐渐过渡到纳米尺度;另一方面,其孔隙间的连通性也由高连通性逐渐向低连通性转变。

    砂土介质一般具有较强的渗透性,孔隙尺寸一般在几十微米到几毫米之间。Mahmoodlu等[16]对不同颗粒级配的砂土介质进行球颗粒堆积模拟,并抓取孔隙结构。抓取获得的孔隙结构信息包括:孔隙半径分布在20~250 μm,孔喉尺寸在20~200 μm。建立等效孔隙网络模型需要岩土介质孔隙连通性信息,本文通过文献调研总结出了砂土、砂岩到页岩的平均配位数范围的重要参数[8-9, 17],如表1所示。

    表  1  不同岩土介质孔隙连通性
    Table  1.  Pore connectivity of different geomaterials
    土壤类别孔隙平均配位数
    砂土6.5~12
    砂岩3.5~4.5
    页岩≤3.5
    下载: 导出CSV 
    | 显示表格

    砂土等效孔隙网络模型建立基于以上参数信息,砂土孔隙结构具有较好的连通性。文中选取砂土介质作为典型的高连通性各向同性岩土介质,对比Mahmoodlu等[16]的数据,本文分别建立了不同孔隙度、不同平均孔隙尺寸的砂土介质孔隙网络模型,通过对固有渗透率的计算来验证模型有效性,同时可以为孔隙度、平均孔隙尺寸对砂土介质固有渗透率影响分析提供参考。表2计算的固有渗透率与Mahmoodlu等[16]通过抓取孔隙结构计算的结果处在同一数量级(10-10 m2),图5为选取孔隙度0.30时生成的砂土介质孔隙网络模型。

    表  2  不同孔隙度砂土介质固有渗透率
    Table  2.  Intrinsic permeability for sand soil with different porosities
    孔隙度 Φ固有渗透率k/(10-10 m2)
    0.3004.72
    0.3254.75
    0.3504.97
    0.3755.48
    0.4006.20
    下载: 导出CSV 
    | 显示表格
    图  5  砂土介质孔隙网络模型(高连通性各向同性)
    Figure  5.  Pore-network model for sand (high connectivity isotropy)

    许多学者针对Berea和Fontainebleau砂岩做了大量试验,包括CT扫描重构三维孔隙网络模型[10, 18-21],获取了大量基础数据。砂岩的平均配位数可根据表1选取。本文为了对比方便,基于Berea和Fontainebleau砂岩CT扫描后重构获得数据建立等效孔隙网络模型。

    等效孔隙网络模型的基本参数如表3所示[10, 18-21],表3中同时对比了等效孔隙网络模型计算的渗透率以及与试验测得的渗透率。

    表  3  砂岩等效孔隙网络模型基本参数
    Table  3.  Parameters of equivalent pore-network model for sandstone
    砂岩类别模型尺寸孔隙度最大孔隙尺寸/μm最小孔隙尺寸/μm平均孔隙尺寸/μm平均配位数计算值渗透率/(10-15 m2)试验值渗透率/(10-15 m2)
    Berea砂岩15 m×15 m×15 m0.20235.10.01403.0403350
    Fontainebleau砂岩15 m×15 m×15 m0.050.8710.010.132.50.001520.00120
    下载: 导出CSV 
    | 显示表格

    图6所示为依据Berea砂岩试验数据建立的等效孔隙网络模型。图7所示为依据Fontainebleau砂岩试验建立的等效孔隙网络模型。表3所计算的Berea砂岩和Fontainebleau砂岩渗透率与试验结果吻合较好,也验证了该等效孔隙网络模型对于低连通性岩土介质水力特征参数计算的有效性。

    图  6  Berea砂岩等效孔隙网络模型(低连通性各向同性)
    Figure  6.  Equivalent pore-network model for Berea sandstone (low connectivity isotropy)
    图  7  Fontainebleau砂岩等效孔隙网络模型(低连通性各向同性)
    Figure  7.  Equivalent pore-network model for Fontainebleau sandstone (low connectivity isotropy)

    各向异性是岩土介质典型特性,尤其对于低渗透岩土介质,主渗透方向渗透系数与其余方向相差较大。文中选取页岩作为低连通性各向异性岩土介质原型,建立页岩基质的等效孔隙网络模型,通过等效孔隙网络模型反映各向异性特性。Barnett、Macellous、Eagle Ford以及国内的龙马溪组页岩为典型的页岩气产地,许多学者也针对上述页岩做了大量研究[22-25],本文所选取的页岩基质的物理性质数据也源于以上产地的页岩典型参数(表4)。

    表  4  典型页岩基本岩石物理参数
    Table  4.  Physical parameters of typical shale
    基本参数最大孔隙尺寸/nm最小孔隙尺寸/nm平均孔隙尺寸/nm孔隙度平均配位数X方向折减因子Y方向折减因子Z方向折减因子
    取值500503000.0730.300.350.45
    下载: 导出CSV 
    | 显示表格

    页岩储层水平向具有较好的连通性,垂直向相对连通性较弱,在建立页岩基质等效孔隙网络模型时需要考虑到此特性。图8所示为建立的页岩基质孔隙网络模型,X方向为主流动方向,连通性高于Y,Z方向。计算的典型页岩基质渗透率分别为:kx=5.53×10-17 m2,ky=4.38×10-17 m2,kz=8.83×10-18 m2

    图  8  页岩等效孔隙网络模型(低连通性各向异性)
    Figure  8.  Equivalent pore-network model for shale (low connectivity anisotropy)

    本文孔隙尺度动态两相流模型基于自主编制的MATLAB计算程序实现。孔隙尺度两相流动模型可以直观反映孔隙内不同流体的动态分布及流动过程中的优势流动、阻滞等现象。本文研究集中在微纳米孔隙尺度岩土介质中两相流动特性,此类岩土介质一般孔隙与孔喉间的尺寸差异较大,在建立两相流动模型时可以忽略孔喉的体积,孔喉只起到流动通道的作用。但需要考虑孔喉断面面积以及孔隙间喉道的长度。动态流动中考虑气、液界面在孔隙内的动态扩张。

    本研究模型在计算过程中采用双压力模型,即孔隙内同时考虑润湿相与非润湿相,气液界面流动非单一的活塞流,弯液面动态扩展可以得到较好反映(图9所示)。计算的基本方程包括:计算压力项的流量连续方程(式2),不考虑压缩性;求解各相态饱和度的各相质量守恒方程(式6)。

    图  9  弯液面动态入侵过程
    Figure  9.  Dynamic invasion process of meniscus
    nj=1(Qij,w+Qij,nw)=0, (2)
    Qij,α=Kij,αΔpij,α   (α=w,nw), (3)

    式中,Kij为喉道导流系数,润湿相、非润湿相的导流系数表达式主要取决于各自的黏滞系数:

    Kij,w=πr4ij8μwlij, (4)
    Kij,nw=πr4ij8μglij, (5)
    ViΔSi,αΔt=nj=1Qij,α      (α=w,nw), (6)

    式中,润湿相与非润湿相饱和度之和为1,即

    Si,w+Si,nw=1 (7)

    研究中针对孔隙内气液界面的扩展提出了相应的基质吸力-饱和度动态变化关系。两相流动计算中考虑到页岩基质内孔隙、孔喉尺寸差异大,假设非润湿相界面在孔隙内呈现出球形气泡状扩展。图9可以看出在该假设条件下,弯液面与孔隙边壁的接触角为90°,弯液面的半径在孔喉半径(rij)与孔隙半径(ri)之间。式(8)中给出了基质吸力与弯液面半径间的动态关系,通过本文的上述假设,弯液面半径与孔隙内饱和度的动态关系可以建立,建立后以饱和度代替弯液面半径推导出基质吸力与孔隙内饱和度间的动态关系。

    pc=pnwpw=2σcosαr, (8)

    式中,σ为气液界面张力系数,一般气、水界面在20℃时取σ=0.0728 N/m,α为接触角。

    弯液面半径与饱和度间动态关系需要定解条件,本文中初始时假设弯液面半径为孔喉ij半径,最终非润湿相充满孔隙时,弯液面半径为孔隙i半径。

    r=rij,    Sw=1r3ij2r3i , (9)
    r=ri,    Sw=0  (10)

    饱和度Sw与孔隙内液相体积分数相关,液相体积又与孔隙半径呈立方关系,因此可以假设弯液面半径(rc)与饱和度关系为

    rc=aS1/3w+b  (11)

    联合式(8)可得任一孔隙i中两相流动基质吸力与饱和度的关系式

    pc=2σcosαriri(ririj)3r3i0.5r3ijS1/3w   (12)

    获得基质吸力与饱和度关系后,联立式(2),(6),(7)可进行整个孔隙网络模型中非润湿相、润湿相压力以及饱和度分布。计算中压力采用隐式计算、饱和度采用显式计算。获得各孔隙非润湿相压力后,根据式(8)可计算出润湿相压力。获得压力项后,需要更新相应时间步各孔隙内饱和度,根据式(6),将润湿相作为计算项,展开后可得

    ViΔSi,wΔt=nj=1Kij,w[(pi,nwpi,c)(pj,nwpj,c)] (13)

    显式计算对于时间步要求较高,需要满足孔隙间的协调稳定性,本文中任一孔隙可出现吸水、排水两种情况,依据以上两种情况分别求解各孔隙吸水、排水时间步:式(14)为排水孔隙时间步计算,考虑单一时间步孔隙排水完全;式(15)为吸水孔隙时间步计算,考虑单一时间步孔隙内充满水。上述时间步针对网络模型内任一孔隙进行计算,当获得所有孔隙排完水或充满水时间步后,选取所有时间步最小值作为整个孔隙网络模型计算时间步,该方法可以确保每一时间步只有单个孔隙充满水或充满气体。

    Δti=ViSi,wnj=1Qij,nw, (14)
    Δti=Vi(1Si,w)nj=1Qij,w (15)

    本文模型在计算过程中为了提高计算效率,非润湿相在单一时间步可以同时充满多个孔隙。整个孔隙网络模型计算时间步需要满足预设的误差限,误差限的设置可以确保总的计算物理时长。多个孔隙同时填充时,需要判断该孔隙是否满足非润湿相入侵条件,即孔隙两端非润湿相与润湿相压差是否大于孔喉处的基质吸力,若大于则该孔隙可以填充,否则搜索下一孔隙。

    Δti(Δt)min(Δt)min<Δerror (16)

    已有页岩储层数据表明,页岩储层中初始水饱和度约为30%,伴随着水力劈裂过程中大量劈裂液注入,储层中液相饱和度增加,页岩气体的流动会呈现出显著的两相流动特性,尤其在裂隙附近区域。页岩基质非常致密,孔隙连通性低,水力劈裂后基质内吸水深度在0.5 m左右,因此水力劈裂后页岩气在裂隙涵盖区域呈现出非润湿相冲破润湿相阻隔的现象。

    由于页岩中有机质部分呈现疏水性,而无机质部分则因含有黏土矿物表现出亲水性,这使得页岩储层呈现出混合润湿特性。此外,页岩中纳米级孔道以及天然微裂隙在气液两相界面具有较大基质吸力,使得劈裂液大部分滞留在储层中无法回收。现场监测数据表明,劈裂液的回收率只有30%±10%[26]。这也说明页岩储层在水力劈裂后呈现出高液相饱和度特征,本文选取裂隙附近区域作为页岩基质内两相流动特性物理模型,研究页岩气开采中两相流动过程、优势渗流特征以及页岩基质相对渗透率曲线,为实际页岩气开采宏观计算模拟提供基础数据。

    页岩基质模型选取如图10所示,页岩基质紧邻水力劈裂裂隙,左侧为含气储层,右侧为水饱和裂隙,本节所指水饱和即指劈裂液饱和。该物理模型针对的基本问题为气相击穿水饱和页岩基质的水动力特性。页岩基质等效孔隙网络模型基本物理参数与2.3节相同,包括平均孔隙尺寸分布300 nm,平均孔喉尺寸6 nm,孔隙度为7%。模型的几何尺寸在X,Y,Z方向为8 μm×5 μm×5 μm,渗流主轴沿着X方向。

    图  10  页岩基质两相流动物理模型
    Figure  10.  Physical model for two-phase flow in shale matrix

    基质内初始条件:初始孔隙水压力为5 MPa,初始状态为水饱和,即Sw=1。

    模型边界条件:流动方向如图沿X轴,上游边界(左侧)为定压边界,设定压力为20 MPa,上游边界连接含气储层,设定水饱和度为零,即Snw=1。下游边界连接裂隙,为恒定压力边界,压力为5 MPa,该压力与初始页岩基质内等压。下游裂隙为水饱和状态(Sw=1)。

    图11所示为页岩基质不同孔隙连通性对于残余饱和度影响,可以看出随着平均配位数(Nave)的增加,页岩基质内残余饱和度降低。平均配位数增加到5时,残余饱和度下降为16%,较低残余饱和度更有利于页岩气体形成连续通道产出。在计算特定岩土介质相对渗透率曲线时(图12),孔隙尺度动态两相流模型可以摒弃对经验公式的依赖,各相态的相对渗透率为各相态出流量与单相流出流量的比值。为了避免边界条件的影响,选取内部多层孔隙网络进行相对渗透率特性计算,计算结果反映出较好的稳定性,残余饱和度基本维持在30%左右。由于该孔隙网络模型中孔隙为球形,孔喉通道为圆柱形,润湿相沿孔隙、孔喉壁面处的流动文中不予考虑,因此计算的相对渗透率曲线即使在残余水饱和度30%情况下,非润湿相的相对渗透率非常高。

    图  11  孔隙连通性与饱和度变化关系
    Figure  11.  Variation of saturation with pore connectivity
    图  12  相对渗透率关系曲线
    Figure  12.  Relative permeability curves

    本文建立的三维孔隙网络多相渗流计算模型可以直观反映不同相流体在孔隙介质中动态运移及分布情况。图13(a),(b)分别表示饱和度分布以及最终残余饱和度分布,图13(c)为非润湿相压力分布,可以看出岩土介质中由于孔隙尺寸分布、连通性等差异会导致明显的优势渗流通道,图13(a)中基质上部存在明显优势通道,饱和度快速下降,图13(c)中的非润湿相压力也与图13(a)饱和度分布结果一致。此外,本文对比了平均配位数为4时孔隙网络模型饱和度及残余饱和度三维时空分布图(图14),可以看出随着平均配位数的增加,残余饱和度显著减小。基质内的残余饱和度由平均配位数为3时的Sr=31%下降到平均配位数为4时的Sr=16%。

    图  13  两相流空间展布(平均配位数=3)
    Figure  13.  Spatial distribution of two-phase flow (Nave=3)
    图  14  两相流空间展布(平均配位数=4)
    Figure  14.  Spatial distribution of two-phase flow (Nave=4)

    本文建立了复杂岩土介质孔隙结构模型。将不同岩土介质孔隙结构特性概化为孔隙尺寸分布、孔喉尺寸分布、连通性、孔隙度等基本物理参数,构建相应的等效孔隙网络模型,模型通过对各配位键设置阈值来反映其存在概率,体现岩土介质孔隙连通性;定义连通性折减因子反映岩土介质孔隙结构的多向性以及各向异性特性。孔喉尺寸由相连孔隙决定,满足孔隙、孔喉的空间相关性;通过计算不同岩土介质渗透率与试验数据对比验证模型的有效性。

    在孔隙结构模型基础上,发展了考虑弯液面动态扩展的两相流数学模型,该模型可以较好地反映多相流中非润湿相的动态扩张过程,利用模型分析了页岩气开采中水饱和区域的两相流动过程,获得了页岩气两相流动中击穿曲线以及相对渗透率曲线等水力特征参数;此外,对孔隙连通性进行了敏感性分析,计算结果表明:随着平均配位数增加,基质内残余饱和度降低。

  • 图  1   等效孔隙网络模型构建流程图

    Figure  1.   Flow chart of establishing EPNM

    图  2   孔隙连通代表性单元

    Figure  2.   Representative element of pore connectivity

    图  3   三维孔隙网络模型

    Figure  3.   3D pore-network model

    图  4   孔隙网络模型稳定性验证结果

    Figure  4.   Stability analysis of pore-network model

    图  5   砂土介质孔隙网络模型(高连通性各向同性)

    Figure  5.   Pore-network model for sand (high connectivity isotropy)

    图  6   Berea砂岩等效孔隙网络模型(低连通性各向同性)

    Figure  6.   Equivalent pore-network model for Berea sandstone (low connectivity isotropy)

    图  7   Fontainebleau砂岩等效孔隙网络模型(低连通性各向同性)

    Figure  7.   Equivalent pore-network model for Fontainebleau sandstone (low connectivity isotropy)

    图  8   页岩等效孔隙网络模型(低连通性各向异性)

    Figure  8.   Equivalent pore-network model for shale (low connectivity anisotropy)

    图  9   弯液面动态入侵过程

    Figure  9.   Dynamic invasion process of meniscus

    图  10   页岩基质两相流动物理模型

    Figure  10.   Physical model for two-phase flow in shale matrix

    图  11   孔隙连通性与饱和度变化关系

    Figure  11.   Variation of saturation with pore connectivity

    图  12   相对渗透率关系曲线

    Figure  12.   Relative permeability curves

    图  13   两相流空间展布(平均配位数=3)

    Figure  13.   Spatial distribution of two-phase flow (Nave=3)

    图  14   两相流空间展布(平均配位数=4)

    Figure  14.   Spatial distribution of two-phase flow (Nave=4)

    表  1   不同岩土介质孔隙连通性

    Table  1   Pore connectivity of different geomaterials

    土壤类别孔隙平均配位数
    砂土6.5~12
    砂岩3.5~4.5
    页岩≤3.5
    下载: 导出CSV

    表  2   不同孔隙度砂土介质固有渗透率

    Table  2   Intrinsic permeability for sand soil with different porosities

    孔隙度 Φ固有渗透率k/(10-10 m2)
    0.3004.72
    0.3254.75
    0.3504.97
    0.3755.48
    0.4006.20
    下载: 导出CSV

    表  3   砂岩等效孔隙网络模型基本参数

    Table  3   Parameters of equivalent pore-network model for sandstone

    砂岩类别模型尺寸孔隙度最大孔隙尺寸/μm最小孔隙尺寸/μm平均孔隙尺寸/μm平均配位数计算值渗透率/(10-15 m2)试验值渗透率/(10-15 m2)
    Berea砂岩15 m×15 m×15 m0.20235.10.01403.0403350
    Fontainebleau砂岩15 m×15 m×15 m0.050.8710.010.132.50.001520.00120
    下载: 导出CSV

    表  4   典型页岩基本岩石物理参数

    Table  4   Physical parameters of typical shale

    基本参数最大孔隙尺寸/nm最小孔隙尺寸/nm平均孔隙尺寸/nm孔隙度平均配位数X方向折减因子Y方向折减因子Z方向折减因子
    取值500503000.0730.300.350.45
    下载: 导出CSV
  • [1]

    FATT I. The network model of porous media: I capillary pressure characteristics[J]. Trans. AIME, 1956, 207: 144. doi: 10.2118/574-G

    [2]

    NICHOLSON D, PETROPOULOS J H. Capillary models for porous media: Ⅲ two-phase flow in three-dimensional network with Gaussian radius distribution[J]. Journal of Physics D: Applied Physics, 1971, 4(2): 181-189. doi: 10.1088/0022-3727/4/2/302

    [3]

    DULLIEN F A. Porous Media-fluid Transport and Pore Structure[M]. New York: Academic Press, 1979.

    [4]

    REEVES P C, CELIA M A. A function relationship between capillary pressure, saturation, and interfacial area as revealed by a pore-scale network model[J]. Water Resources Research, 1996, 32(8): 2345-2358. doi: 10.1029/96WR01105

    [5]

    ACHARYA R C, SJOERD E A T M, LEIJNSE A. Porosity-permeability properties generated with a new 2-parameter 3D hydraulic pore-network model for consolidated and unconsolidated porous media[J]. Advances in Water Resources, 2004, 27: 707-723. doi: 10.1016/j.advwatres.2004.05.002

    [6]

    RAOOF A, HASSANIZADEH S M. A New method for generating pore-network models of porous media[J]. Transport in Porous Media, 2009, 81(3): 391-407.

    [7]

    GAO S Y, MEEGODA J N, HU L M. Two methods for pore network of porous media[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2012, 36: 1954-1970. doi: 10.1002/nag.1134

    [8] 王晨晨, 姚军, 杨永飞, 等. 基于规则网络的碳酸盐岩多尺度网络模型构建方法研究[J]. 计算力学学报, 2013, 30(2): 0231-0235. https://www.cnki.com.cn/Article/CJFDTOTAL-JSJG201302011.htm

    WANG Chen-chen, YAO Jun, YANG Yong-fei, et al. The construction of carbonate multiscale network model based on regular network[J]. Chinese Journal of Computational Mechanics, 2013, 30(2): 0231-0235. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-JSJG201302011.htm

    [9]

    ARNS J Y, ROBINS V, SHEPPARD A P, et al. Effect of network topology on relative permeability[J]. Transport in Porous Media, 2004, 55: 21-46. doi: 10.1023/B:TIPM.0000007252.68488.43

    [10]

    BLUNT M J, JACKSON M D, PIRI M, et al. Detailed physics, predictive capabilities and macroscopic consequences for pore-network models of multiphase flow[J]. Advances in Water Resources, 2002, 25: 1069-1089. doi: 10.1016/S0309-1708(02)00049-0

    [11]

    JOEKAR-NIASAR V, HASSANIZADEH S M, LEIJNSE A. Insights into the relationships among capillary pressure, saturation, interfacial area and relative permeability using pore-network modeling[J]. Transport in Porous Media, 2008, 74: 201-219. doi: 10.1007/s11242-007-9191-7

    [12]

    KOPLIK J, LASSETER T J. One and two-phase flow in network models of porous media[J]. Chemical Engineering Communications, 1984, 26: 285-295. doi: 10.1080/00986448408940216

    [13]

    JOEKAR-NIASAR V, HASSANIZADEH S M, DAHLE H K. Non-equilibrium effects in capillarity and interfacial area in two-phase flow: dynamic pore-network modelling[J]. Journal of Fluid Mechanics, 2010, 655: 38-71. doi: 10.1017/S0022112010000704

    [14]

    GAO S Y, MEEGODA J N, HU L M. A dynamic two-phase flow model for air sparging[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2013, 37: 1801-1821. doi: 10.1002/nag.2109

    [15]

    HUANG X W, BANDILLA K W, CELIA M A. Multi-physics pore-network modeling of two-phase shale matrix flows[J]. Transport in Porous Media, 2016, 111: 123-141. doi: 10.1007/s11242-015-0584-8

    [16]

    MAHMOODLU M G, RAOOF A, SWEIJEN T, et al. Effects of sand compaction and mixing on pore structure and the unsaturated soil hydraulic properties[J]. Vadose Zone Journal, 2016, 15(8): 1-11. doi: 10.2136/vzj2015.12.0157

    [17]

    HAUGHEY D P, BEVERIDGE G G. Local voidage variation in arandomly packed bed of equal-sized spheres[J]. Chemical Engineering Science, 1966, 21(10): 905-915. doi: 10.1016/0009-2509(66)85084-4

    [18]

    ØREN P E, BAKKE S. Process based reconstruction of sandstones and prediction of transport properties[J]. Transport in Porous Media, 2002, 46: 311-343. doi: 10.1023/A:1015031122338

    [19]

    VALVATNE P H, PIRI M, LOPEZ X, et al. Predictive pore-scale modeling of single phase and multiphase flow[J]. Transport in Porous Media, 2005, 58: 23-41. doi: 10.1007/s11242-004-5468-2

    [20]

    MADONNA C, QUINTAL B, FREHNER M, et al. Synchrotron-based X-ray tomography microscopy for rock physics investigations[J]. Geophysics, 2013, 78(1): 53-64. doi: 10.1190/geo2012-0113.1

    [21]

    SHARQAWY M H. Construction of pore network models for Berea and Fontainebleau sandstones using non-linear programing and optimization techniques[J]. Advances in Water Resources, 2016, 98: 198-210. doi: 10.1016/j.advwatres.2016.10.023

    [22]

    SOEDER D J. Porosity and permeability of Eastern Devonian gas shale[J]. SPE Formation Evaluation, 1988, 3: 116-124. doi: 10.2118/15213-PA

    [23]

    CURTIS M E, SONDERGELD C H, AMBROSE R J, et al. Microstructural investigation of gas shales in two and three dimensions using nanometer-scale resolution imaging[J]. AAPG Bulletin, 2012, 96(4): 665-677. doi: 10.1306/08151110188

    [24]

    CHALMERS G R, BUSTIN R M, POWER I M. Characterization of gas shale pore systems by porosimetry, pycnometry, surface area, and field emission scanning electron microscopy/transmission electron microscopy image analyses: examples from the Barnett, Woodford, Haynesville, Marcellus, and Doig units[J]. AAPG Bulletin, 2012, 96(6): 1099-1108. doi: 10.1306/10171111052

    [25] 郭为, 熊伟, 高树生, 等. 页岩气等温吸附/解吸特征[J]. 中南大学学报(自然科学版), 2013,44(7): 2836-2840. https://www.cnki.com.cn/Article/CJFDTOTAL-ZNGD201307029.htm

    GUO Wei, XIONG Wei, GAO Shu-sheng, et al. Isothermal adsorption/desorption characteristics of shale gas[J]. Journal of Central South University (Science and Technology), 2013, 44(7): 2836-2840. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-ZNGD201307029.htm

    [26]

    ROYCHAUDHURI B, TSOTSIS T, JESSEN K. An experimental investigation of spontaneous imbibition in gas shales[J]. Journal of Petroleum Science and Engineering, 2013, 111: 87-97. doi: 10.1016/j.petrol.2013.10.002

  • 期刊类型引用(12)

    1. 赵明凯,孔德森,滕森,邓美旭. 应力作用下岩石介质两相流束缚流体饱和度分形预测模型研究. 岩土工程学报. 2024(04): 871-879 . 本站查看
    2. 叶美娜,仲海书,张震,李保珠,田向盛,张英. 甘肃早子沟金矿典型岩石微观结构特征及其渗透性研究. 地质与勘探. 2024(06): 1142-1152 . 百度学术
    3. 孔德森,赵明凯,时健,滕森. 基于分形维数特征的岩石介质气-水相对渗透率预测模型研究. 岩土工程学报. 2023(07): 1421-1429 . 本站查看
    4. 王萌,于群丁,肖源杰,华文俊,李文奇. 振动荷载下不同级配的基床粗粒土填料孔隙连通性特征研究. 中南大学学报(自然科学版). 2023(11): 4436-4448 . 百度学术
    5. 周新超,马小晶,廖翔云,齐思维,李宏煜. 磨料水射流冲击孔隙岩体的SPH模拟研究. 岩土工程学报. 2022(04): 731-739 . 本站查看
    6. 朱维耀,李华,邓庆军,马启鹏,刘雅静. 多孔介质细观流动理论研究进展. 工程科学学报. 2022(05): 951-962 . 百度学术
    7. 李博,王晔,邹良超,杨磊. 岩石裂隙内浆液–水两相流可视化试验与驱替规律研究. 岩土工程学报. 2022(09): 1608-1616+2-3 . 本站查看
    8. 张兴昊,林丹彤,胡黎明. 基于等效孔隙网络模型的水动力弥散数值模拟. 清华大学学报(自然科学版). 2022(12): 1906-1914 . 百度学术
    9. 叶曾云,郑芳云,谈勋勋,石旺. 干密度对花岗岩全风化土的渗流特性研究. 建筑施工. 2022(10): 2482-2486+2494 . 百度学术
    10. 王志兵,谈勋勋,王远,孙广,刘金明. 不同压实度花岗岩残积土的渗流模拟. 科学技术与工程. 2021(14): 5913-5921 . 百度学术
    11. 蔡沛辰,阙云,李显. 非饱和花岗岩残积土水-气两相驱替过程数值模拟. 水文地质工程地质. 2021(06): 54-63 . 百度学术
    12. 蔡沛辰,阙云,李显. 虚拟仿真APP开发在土体孔隙细观两相渗流研究中的应用. 水资源与水工程学报. 2021(06): 207-214 . 百度学术

    其他类型引用(10)

图(14)  /  表(4)
计量
  • 文章访问数:  524
  • HTML全文浏览量:  34
  • PDF下载量:  318
  • 被引次数: 22
出版历程
  • 收稿日期:  2019-04-19
  • 网络出版日期:  2022-12-07
  • 刊出日期:  2019-12-31

目录

/

返回文章
返回