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

CO2置换开采CH4水合物的深海地层多场耦合连续介质数值方法研究

蒋明镜, 李子煜, 李承超, 姜朋明

蒋明镜, 李子煜, 李承超, 姜朋明. CO2置换开采CH4水合物的深海地层多场耦合连续介质数值方法研究[J]. 岩土工程学报, 2025, 47(7): 1354-1364. DOI: 10.11779/CJGE20240762
引用本文: 蒋明镜, 李子煜, 李承超, 姜朋明. CO2置换开采CH4水合物的深海地层多场耦合连续介质数值方法研究[J]. 岩土工程学报, 2025, 47(7): 1354-1364. DOI: 10.11779/CJGE20240762
JIANG Mingjing, LI Ziyu, LI Chengchao, JIANG Pengming. Multi-field coupling continuum numerical method for exploiting CH4 by CO2 replacement in deep-sea formation[J]. Chinese Journal of Geotechnical Engineering, 2025, 47(7): 1354-1364. DOI: 10.11779/CJGE20240762
Citation: JIANG Mingjing, LI Ziyu, LI Chengchao, JIANG Pengming. Multi-field coupling continuum numerical method for exploiting CH4 by CO2 replacement in deep-sea formation[J]. Chinese Journal of Geotechnical Engineering, 2025, 47(7): 1354-1364. DOI: 10.11779/CJGE20240762

CO2置换开采CH4水合物的深海地层多场耦合连续介质数值方法研究  English Version

基金项目: 

国家自然科学基金重点项目 52331010

国家自然科学基金重大项目课题 51890911

国家自然科学基金青年项目 52309139

江苏省高校自然科学面上基金项目 22KJB570004

江苏省研究生科研创新计划项目 KYCX23_3330

详细信息
    作者简介:

    蒋明镜(1965—),男,教授,博士生导师,主要从事天然结构性黏土、砂土、非饱和土、太空土和深海能源土宏观微观试验、本构模型和数值分析方面的研究。E-mail: mingjing.jiang@usts.edu.cn

    通讯作者:

    姜朋明, E-mail: pmjiang@usts.edu.cn

  • 中图分类号: TU432

Multi-field coupling continuum numerical method for exploiting CH4 by CO2 replacement in deep-sea formation

Funds: 

Key Program of National Natural Science Foundation of China 52331010

Major Program of National Natural Science Foundation of China 51890911

Youth Program of National Natural Science Foundation of China 52309139

the Natural Science Foundation of the Higher Education Institutions of Jiangsu Province 22KJB570004

the Jiangsu Provincial Graduate Research and Practice Innovation Program KYCX23_3330

  • 摘要: 通过诱导水合物分解的开采方式将致使含天然气水合物沉积物(后称能源土)力学性质劣化,进而可能引发一系列地质灾害。CO2置换法能有效改善地层强度,但目前鲜有研究考虑CO2置换开采的分析方法。基于此,针对CO2置换法开展数值方法研究和应用。通过拓展TOUGH+HYDRATE数学模型、嵌入Chen-Guo模型计算CH4-CO2多元水合物相平衡条件,建立CO2置换的数值模拟器T+MixH V1.0,并与室内试验和已有数值计算结果对比,验证模拟器的可靠性。随后与FLAC3D结合,建立了CO2置换开采工况下的温-压-力-化多场耦合数值分析方法。最后对比分析了南海储层降压、置换开采CH4水合物两种情景以及经历不同开采时间后的能源土地基载荷板试验,结果表明:相比于降压开采方式,置换开采提高了产气量并减缓了地层沉降,其原因在于降压造成的能量损失使CH4水合物分解范围有限,而置换开采在合适的温压条件下自发进行且生成的CO2水合物对地层有一定支撑作用。此外,能源土地基承载特性主要受控于两种因素:孔压变化引起的有效应力变化,水合物分解/生成引起的胶结强度变化。
    Abstract: The extraction method that induces hydrate dissociation may deteriorate the mechanical properties of methane hydrate bearing sediments (MHBS), potentially leading to a series of significant disasters. The carbon dioxide (CO2) replacement method can effectively improve the sediment strength, but few methods have considered its mechanical response. So, a numerical method for the CO2 replacement method as well as its application is studied. Based on TOUGH+HYDRATE, a numerical simulator, T+MixH V1.0 is developed to simulate the CO2 replacement process by incorporating the Chen-Guo hydrate model to calculate the multi-component hydrate phase equilibrium conditions. The reliability of the simulator is validated by comparing with other experimental and numerical results. Then, a thermo-hydro-mechanic-chemical multi-field coupling numerical method for the CO2 replacement is established by coupling with FLAC3D. Finally, the gas hydrate exploitation is simulated through depressurization and CO2 replacement methods in the South China Sea, as well as the plate loading tests on MHBS ground during gas production. The results indicate that the CO2 replacement increases gas production and mitigates subsidence of the ground, which comes from that depressurization causes energy loss, limiting the decomposition of CH4 hydrate, whereas the CO2 replacement method occurs spontaneously, and the generated CO2 hydrate provides support in the ground. Besides, the bearing characteristics of MHBS ground are mainly affected by two factors: the effective stress variation caused by pore pressure changes, and the bond strength variation induced by gas hydrate decomposition/generation.
  • 页岩油是一种在中国存储量较丰富且仅次于石油的油气资源[1]。为了满足中国日益增长的能源需求,大力发展页岩产业和页岩开采技术是十分必要的。射流破岩技术能够大面积且有导向地碎裂储层岩体,方便页岩气和煤层气等油气资源的解吸与渗流,已被广泛应用于页岩的开采领域[2-3]

    目前,国内外学者已对磨料水射流破岩问题展开了大量的研究[4-6],但由于射流和岩体具有不同的物理属性和力学性能,尤其是岩体中还存在着大量的孔隙[7],再加上磨料颗粒对岩体的冲蚀作用,使得岩体的变形破损机制变得极其复杂[8]。因此,开展磨料水射流冲击破岩的研究,并对含有孔隙结构岩体的破损过程进行分析,能够为射流破岩开采技术的优化和完善提供一定的理论依据。

    随着计算机技术的发展,数值模拟方法已经成为研究射流破岩问题的有效手段之一[9-10]。黄中伟等[11]基于有限元法模拟了液氮磨料射流冲击岩体的过程,发现冲击是岩体破损的主要因素。庄欠伟等[12]应用有限元方法模拟研究了磨料水射流冲击钢筋和混凝土的过程。穆朝民等[13]使用网格法构建了岩体,并分析了在受到磨料水射流冲击后岩体的损伤演化过程。由于网格类的数值方法在模拟大变形等问题时存在着网格畸变的情况,因此,近年来无网格方法得到了大量的关注并广泛应用于射流破岩问题的研究。

    光滑粒子流体动力学(smoothed particle hydro- dynamics,SPH)方法作为一种典型的无网格粒子方法,在处理大变形、应力局部化等问题中具有一定的优势,被广泛采纳并应用于冲击大变形和流固耦合等领域。刘勇等[14]运用SPH方法模拟分析了磨料空气射流破岩的冲蚀过程。李世杰等[15]采用SPH方法模研究了不同直径的水射流对土壤冲蚀坑的影响。赵健等[16]将SPH方法和有限元方法进行耦合,模拟分析了磨料水射流冲击下岩体的损伤演化机制。虽然目前有许多学者已采用SPH方法模拟研究了不同的射流破岩问题,但针对射流冲击下含孔隙岩体的破损过程以及孔隙结构对破岩效果影响的研究相对较少。鉴于实际岩体中含有孔隙结构,本文在岩体模型的内部引入了孔隙,模拟研究了磨料水射流冲击含孔隙岩体的破损过程,并分析了多种孔隙特性对射流破岩效果的影响。

    本文依据SPH方法的基本理论,分别引入描述射流和岩体力学特性的NULL型和H-J-C型本构模型,建立磨料水冲击含孔隙岩体的数值模型。通过模拟磨料水射流破岩问题,验证了方法和模型的有效性。在此基础上,加入随机混合的孔隙,构建出磨料水射流冲击含孔隙结构岩体的动力学模型,对射流破岩的动态过程进行模拟研究,分析了岩体的孔隙特性(孔隙尺寸、孔隙度和孔隙填充物)以及射流参数(射流速度和磨料浓度)对岩体破损效果的影响,为深入理解射流开采中岩体的破损机理和提高开采效率提供了一定的理论依据。

    SPH方法是一种通过利用离散点来构造近似函数的纯拉格朗日方法。对于磨料水射流冲击含孔隙岩体的问题而言,可将含磨料颗粒的水射流和具有孔隙结构的岩体离散为一系列的粒子,通过差值求解非稳态N-S控制方程即可得到粒子的物理属性。因此,宏观变量可通过计算域内离散点的积分插值得到,质点的近似函数为[17]

    f(x)=Ωf(x)W(xx,h)dx (1)

    式中f(x)为三维坐标向量x的函数;Ωx的支持域;xx为粒子间距离;h为SPH粒子光滑长度,光滑长度随时间和空间变化,如图 1所示;W(xx,h)为核函数,用辅助函数θ(xx)定义:

    W(xx,h)=1h(xx)dθ(xx) (2)
    图  1  粒子近似法
    Figure  1.  Particle approximation

    式中,d为空间维数。用粒子近似法将连续形式积分方程换算为离散形式的方程,即

    f(x)=ni=1miρi1h(xx)dθ(xx) (3)

    式中,ρi为粒子i的密度,mi为粒子i的质量。

    SPH方法下的N-S方程中质量守恒方程为

    ρi=nj=1mjWijNj=1mjρjWij (4)

    动量守恒方程为

    DvαiDt=Nj=1mj(σαβiρ2i+σαβjρ2j)Wijxβi (5)

    能量守恒方程为

    DeiDt=12Nj=1mj(Pi+Pjρiρj)vβijWijxβi+μi2ρiεαβiεαβj (6)

    式中,v为速度,e为热能,σ为应力张量,P为压力,m为质量,ρ为密度,μ为动力黏度系数,ε为应变率,N为质点总数,x为坐标轴符号,αβ代表不同坐标轴,ij代表不同质点

    为了保证计算中的精度和稳定性,提高计算效率,研究采用了应用较广的三次B样条核函数[17]

    研究采用了LS-DYNA中的NULL本构模型描述水和磨料,塑性材料是一种在外力的作用下,能够产生显著变形而不被破坏的材料,由于在射流破岩问题中,磨料和水较符合塑性材料的特性[2],因此本文在仿真模拟中将两者都视为完全塑性材料并赋予Gruneisen状态方程[18],即

    p=ρ0C2μ[1+(1γ02)μa2μ2][1(S11)μS2μ2μ+1S3μ3(μ+1)2]+(γ0+aμ)Ea (7)

    式中C为冲击波速度与质点速度变化曲线截距;S1S2S3为冲击波与质点速度变化曲线的斜率系数;γ0为Gruneisen常数;aγ0μ=ρ/ρ01的一阶体积修正量。水和磨料颗粒的本构模型参数[16]表 1所示。为了实现磨料粒子在流场中的真实分布情况,从射流粒子中进行随机抽样,并赋予磨料的属性,实现射流中磨料颗粒的随机混合。

    表  1  水和磨料本构模型的参数
    Table  1.  Parameters of constitutive model for water and abrasive
    名称 ρ0 /(g·cm-3) C/(m·s-1) S1 S2 S3
    1.0 1480 2.56 -1.986 0.228
    磨料 3.5 4569 1.49 0 0.460
    下载: 导出CSV 
    | 显示表格

    在射流的冲击作用下,不同性质的岩体会呈现出不同程度的破损。研究选用页岩材料,并使用H-J-C本构模型描述,页岩本构模型的参数[19]表 2所示。H-J-C本构模型综合考虑了岩石材料损伤、应变率和静水压力对屈服力的影响,是一种能够在岩体大变形、高应变率和高压等情况下应用的模型,其等效强度通过应力、应变速率和损伤状态来描述,屈服面的方程为[20]

    σ=[A(1D)+BpN](1Clnˉε) (8)
    表  2  页岩本构模型的参数
    Table  2.  Parameters of constitutive model for shale
    岩石 ρ/(g·cm-3) G/GPa fc/GPa T/GPa μcrush
    页岩 2.44 13.6 0.0048 0.004 0.001
    下载: 导出CSV 
    | 显示表格

    式中σ为岩体冲击载荷作用下等效应力与静态屈服强度之比;P为岩体实际压力与静态屈服强度之比;¯ε为冲击载荷作用下岩体应变率与静态应变率之比;ABNc为材料的强度参数。

    损伤因子D通过等效塑性应变和塑性体积应变累计得到

    D=Δεp+ΔμpD1(p+T)D2 (9)

    式中ΔεpΔμp分别为一个计算循环内的等效塑性应变和塑性体积应变;D1D2为损伤常数;T=T/feT为材料的最大拉伸强度,fe为材料静态抗压强度。

    为了构建与天然岩体更接近并且能够进行数值计算的孔隙结构,研究中参考最大球算法[21]建立孔隙的过程和随机堆积球模型[22-24],将孔隙离散为球体,允许孔隙相互重叠,重叠程度与孔隙自身相比较小,如图 2所示。

    图  2  球体重叠程度
    Figure  2.  Sphere overlap

    在使用SPH方法将岩体与孔隙的离散后,空间中的任意一个代表孔隙中心位置的粒子,与最近岩体粒子的欧氏距离,为球心与岩体表面相切的最大球的半径,计算公式如下:

    R=(x1x2)2+(y1y2)2+(z1z2)2 (10)

    式中,(x1,y1,z1)为孔隙中心点的坐标,(x2,y2,z2)为岩体点的坐标。

    根据岩体与孔隙的质量比α能够确定出页岩的孔隙度β,之后通过β即可确定出孔隙的个数:

    α=mamw+ma=VaρaVwρ0+Vaρa (11)
    β=VwVa=αρ01αρ0 (12)
    na=3dw3β4π R3 (13)

    式中mamw分别为岩体与孔隙的质量;VaVw分别为岩体与孔隙的体积;ρaρ0分别为岩体与孔隙的密度;na为页岩模型中孔隙的个数;dw为岩体的边长,R为孔隙最大内切球的半径。

    页岩等沉积岩的形成过程十分复杂,因此在建立孔隙结构需要有以下特征[25]:①随机性。在岩石的沉积过程中,孔隙在岩体中随机生成,且在任何点都有均等形成的机会。②独立性。孔隙生成在岩体中达到平衡且静止的状态,不会受到其他孔隙的影响。③长周期性。在岩体沉积、形成的过程中,孔隙生成的方向为重力势能梯度最大的方向,选取的孔隙的位置需要符合此规律。研究中采用了蒙特卡洛法随机生成孔隙,对于空间内任意一个SPH粒子,具有相同的概率被采样,能够使需要生成的na个孔隙均匀的分布在岩体内部,被选取的单个孔隙位置表示为Pi(x,y,z),其描述如下:

    Pi(x,y,z)={x,y,z=Random(L)},ina (14)

    式中,Random(L)是能够满足孔隙分布的随机函数,xyz为坐标信息。

    在随机得到的Pi(x,y,z)位置上生成一个SPH粒子之后,沿着XYZ三个方向继续增加粒子,最终生成单个孔隙,演化过程如图 3所示。

    图  3  孔隙生成过程
    Figure  3.  Formation process of pore

    本文离散球体的规则与SILIN等[21]和文献[2223]一致。为了比较不同尺寸的孔隙之间的差异,建立球体大小相等的孔隙模型,如图 4所示。孔隙的本构模型及状态方程与水相同

    图  4  球体的离散与分布
    Figure  4.  Dispersion and distribution of sphere

    经过对比,发现在相同参数下的构建的岩体仿真模型与真实岩体的孔隙分布[26]基本相同,并且通过三轴压缩试验,得到了在载荷作用下孔隙岩体模型的应力应变曲线,如图 5所示。经计算得知,该模型的弹性模量与文献[27]中真实孔隙岩体相差1%,其他力学参数基本一致,说明了构建的孔隙岩体模型能够非常接近真实岩体孔隙结构。

    图  5  孔隙岩体模型的应力应变曲线
    Figure  5.  Stress-strain curve of model for porous rock

    针对所研究的磨料水射流冲击含孔隙岩体的问题,进行如下的简化假设:①磨料水射流破岩过程中忽略气相环境的影响,只包含水射流、磨料颗粒、岩体和孔隙4种物质;②磨料颗粒与水粒子为相同直径的球体。

    为了保证计算效率,采用了1/4的对称模型,对称面设置为SPH_SYMMETRY_PLANE对称约束,并在对称面的基础上增加SPH_NON_REFLECTING约束来模拟岩石无限大的情况。对称后岩体为50 mm×50 mm×25 mm的立方体,孔隙的直径介于0.3~2 mm,射流喷距为1 mm,模型共由424395个光滑粒子构成,其中2394个为射流粒子,包含了水与磨料。磨料水射流冲击孔隙岩体的几何模型如图 6所示。

    图  6  磨料水射流冲击孔隙岩体模型
    Figure  6.  Model for abrasive water jet impacting pore rock

    研究中选取了影响岩体破损效果的以下因素进行分析:孔隙特性(孔隙尺寸ϕ、孔隙度β和孔隙填充物)以及射流参数(射流速度v和磨料浓度δ)。

    为了验证所建模型的有效性,对李井慧[28]所研究的磨料水射流冲击大理石过程进行了模拟,各项参数均与文献[28]保持一致,图 7为不同射流压力下破损坑截面演化图,由图 7可以看出,随着射流压力的增加,大理石的侵彻深度在逐渐的增大。

    图  7  不同射流压力下破损坑演化图
    Figure  7.  Evolution of damaged pits under different jet pressures

    图 8为试验与模拟结果对比图。由图 8可知,由于真实试验的冲击时间大于仿真模拟,且试验中的磨料水射流能够持续冲击岩体,导致试验切割深度整体大于模拟结果,但岩体的切割深度随射流压力的变化趋势基本相同,均随着射流压力的增加而增大,验证了本文所构建的磨料水射流模型的有效性。

    图  8  试验与模拟结果对比
    Figure  8.  Comparison between experimental and simulation results

    本文还对林晓东等[29]所研究的磨料水射流冲击灰岩问题进行了模拟,对岩体破损形状进行分析和验证。射流直径为4 mm,射流高度为75 mm,磨料浓度为30%,均与文献[29]保持一致。射流冲击灰岩破损坑截面演化过程如图 9所示,当射流撞击岩体的瞬间,受到冲击的区域出现了破损坑,随着冲击时间的增加,破损坑逐渐扩大,最终接近子弹型。由于林晓东等[29]采用网格法描述岩体,受到冲蚀作用的岩石网格单元破坏失效后删除,而本文采用SPH方法对岩体进行了离散,使得破损坑的形状略有差异,但冲孔截面形态基本相同,与试验的冲孔扫描结果[29]吻合,说明使用SPH方法替代网格法构建岩体是可行的。

    图  9  灰岩破损坑对比图
    Figure  9.  Comparison of damaged limestone pits

    图 10给出了射流速度模拟结果对比,由图 10可知,当射流速度小于200 m/s时,冲击速度对破损坑半径的影响较小;当射流速度超过200 m/s时,破损坑半径随着射流速度的增加急剧增大,不同模拟方法所得的岩体冲孔半径随冲击速度的变化趋势较为相似,进一步验证了本研究基于SPH方法建立的数值模型的有效性。后续关于磨料水射流冲击含孔隙岩体的研究以此模型为基础,为了更好的观察破损坑形态的变化,射流高度调整为25 mm。

    图  10  射流速度模拟结果对比
    Figure  10.  Comparison of simulated results of jet velocity

    为了全面分析岩体的破损情况,对不同工况下的破岩过程进行了模拟,并提取了岩体破损图像的特征,研究表明破损坑深和破损面积的变化规律基本一致,因此本文选取了表面坑径(r/R)和破损面积(s/S)为岩体破损指标,并进行了无量纲处理,如图 11所示,其中rs分别为岩体破损坑的表面坑径和纵截面积,RS分别是岩体的宽度和纵截面积。

    图  11  破岩指标参数示意图
    Figure  11.  Diagram of rock breaking index parameters

    选取孔隙度β=10%,孔隙尺寸为1 mm3,孔隙填充物为水,射流速度v=500 m/s,磨料浓度δ=5%的数值模型进行Von-Mises应力分析,图 12为射流冲击过程的Von-Mises应力演化图。

    图  12  Von-Mises应力演化图
    Figure  12.  Evolution diagram of Von Mises stress

    图 12可知,射流撞击到岩体的瞬间,受到冲击的岩体区域出现了较大的应力集中现象,当Von-Mises应力超过岩体的强度极限时,会导致岩体与射流的接触区域破损。由于岩体部分采用了脆性较高的页岩,当t=0.03 ms时,应力波开始呈现出裂纹状;距离射流较近的孔隙边缘出现了应力集中现象,并且孔隙的形状发生了改变,因孔隙填充物为水,导致孔隙内部应力较小。随着冲击时间的增加,冲击区域两侧的孔隙持续受到岩体挤压。当t≥0.06 ms时,孔隙被完全压碎,导致无法看出孔隙具体的形状和位置,岩体的破损程度不断加剧,应力波在岩体内部向四周传播并且以裂纹状继续扩散。

    孔隙的尺寸会影响射流破岩的效率和岩体的破损程度。研究中取孔隙尺寸ϕ=1~6(以ϕ=3为例说明,孔隙直径离散为3个粒子),如图 13所示。孔隙度β=15%,孔隙填充物为水,磨料体积浓度δ=5%,射流速度v=500 m/s。

    图  13  孔隙尺寸示意图
    Figure  13.  Diagram of pore size

    图 14给出了不同孔隙尺寸下t=0.1 ms时岩体破损区域与破损指标随冲击时间的变化曲线。由图 14(a)可看出,在相同的孔隙度下,当孔隙尺寸较小时(ϕ≤4),岩体内孔隙数量较多,能够均匀的分布在岩体中,随着孔隙体积的增大,岩体破损加剧;当孔隙尺寸ϕ≥5时,岩体的破损坑的形状变化较大。由图 14(b)可以看出,随着射流冲击时间的增加,岩体破损坑的表面坑径和破损面积都逐渐增大。当孔隙较小时(ϕ≤4),随着孔隙尺寸的增加,破损面积不断增大,但表面坑径之间的差异较小;当孔隙尺寸ϕ≥5时,表面坑径之间出现了较大的差异。

    图  14  不同孔隙尺寸
    Figure  14.  Different pore sizes

    为了探究岩体的破损是否受到孔隙分布位置的影响,研究中对孔隙尺寸ϕ≤4的4种情况进行了多次模拟并进行了对比,发现在相同孔隙尺寸下,岩体的破损情况受到孔隙分布位置的影响较小。在孔隙尺寸较大的两种情况下(ϕ≥5),分别增加3种不同的孔隙随机分布。图 15为不同的孔隙分布位置下t=0.1 ms时岩体破损区域与破损指标随冲击时间的变化曲线。结合图 1415可以得知,当孔隙较大并且在相同的孔隙度下,不同孔隙分布位置的破损坑形状有显著的差别。当孔隙尺寸ϕ≥5时,不同的孔隙分布位置下破损坑的表面坑径以及破损面积差异较大,说明了在磨料水射流冲击的过程中含大孔隙岩体的破损受到孔隙随机分布位置的影响,因此后续的研究中孔隙尺寸选取为1 mm3ϕ=3)。

    图  15  不同孔隙分布位置
    Figure  15.  Distribution location of different pores

    孔隙度是影响岩体破损的重要原因。研究中取孔隙度β=0~25%,孔隙尺寸为1 mm3ϕ=3),孔隙填充物为水,射流速度v=500 m/s,磨料体积浓度δ=5%。图 16给出了不同孔隙度下t=0.1 ms时岩体破损区域与破损指标随冲击时间的变化曲线。从图 16可以看出,随着岩体内部孔隙的增多,在磨料水射流的冲击作用下岩体破损坑的表面坑径和破损面积均增大,尤其是破损面积急剧增加,说明孔隙率越高的岩体越容易破碎。当岩体内不含孔隙(β=0%)时,岩体的破损面积略大于含孔隙较少的岩体(0%<β≤10%),这也说明少量的孔隙及其填充物可能会对磨料水射流产生干扰并起到一定的缓冲作用,降低了射流的冲击破岩能力。

    图  16  不同孔隙度
    Figure  16.  Different porosities

    不同区域的页岩孔隙中蕴藏着不同的填充物质。研究选取了水、轻油、重油和空气4种孔隙填充物,由于空气的密度与液体相比非常小,可近似忽略,因此将孔隙粒子删除视为空气,填充物的参数[30]表 3所示,孔隙尺寸为1 mm3ϕ=3),孔隙度β=15%,孔隙的分布位置完全相同;射流速度v=500 m/s,磨料体积浓度δ=5%。

    表  3  孔隙填充物的参数
    Table  3.  Parameters of pore fillers
    材料 ρ /(g·cm-3) 动力黏度系数/(pa·s)
    1.00 1.009×10-3
    重油 1.10 2.5×10-3
    轻油 0.85 2.5×10-3
    下载: 导出CSV 
    | 显示表格

    图 17给出了不同填充物下t=0.1 ms时岩体破损区域与破损指标随冲击时间的变化曲线。由图 17可以看出,当填充物为液体时,岩体破损坑的形状较为相似,表面坑径以及破损面积差异较小;当孔隙填充物为空气时,岩体破损坑的面积最大,但表面坑径最小,这也说明与其他填充物相比,含气体的岩体更容易沿着磨料水射流的冲击方向发生破损。

    图  17  不同孔隙填充物
    Figure  17.  Different pore fillers

    射流速度是决定开采破损效率的重要因素。射流速度取v=100~500 m/s,磨料浓度δ=5%,孔隙尺寸为1 mm3ϕ=3),孔隙度β=15%,孔隙填充物为水。图 18给出了不同射流速度下t=0.1 ms时岩体的破损区域与破损指标随冲击时间的变化曲线。由图 18可以看出,岩体破损坑的表面坑径和破损面积随着射流速度的增大而增加。当射流速度v>200 m/s时,表面坑径和破损面积急剧增加,说明磨料水射流速度越高,破损孔隙岩体能力越强。

    图  18  不同射流速度
    Figure  18.  Different jet velocities

    磨料浓度也是影响射流破岩效率和破损精度的重要因素之一。根据实际磨料水射流破岩过程中的磨料颗粒体积占比,磨料浓度δ=0~25%,射流速度v=500 m/s,孔隙尺寸为1 mm3ϕ=3),孔隙度β=15%,孔隙填充物为水。图 19给出了不同磨料浓度下t=0.1 ms时岩体的破损区域与破损指标随冲击时间的变化曲线。

    图  19  不同磨料浓度
    Figure  19.  Different abrasive concentrations

    图 19可知,随着磨料浓度的增加,岩体破损坑的纵截面形状发生了明显的改变。当水射流中添加磨料颗粒后,岩体破损坑的表面坑径急剧减小,磨料浓度越大,岩体破损坑的表面坑径和破损面积越小;当磨料浓度增大到10%以上时,各浓度之间的表面坑径以及破损面积差异较小。

    为了深入分析磨料浓度对岩体破损效果的影响,给出了底部坑径(l/R)指标,如图 20所示。

    图  20  底部坑径示意图
    Figure  20.  Diagram of diameter of bottom pit

    图 21为不同磨料浓度下底部坑径随冲击时间的变化曲线,由图 21可知,当水射流中混有磨料颗粒时,岩体破损坑的底部坑径急剧增加,并且随着磨料浓度的增加而增大。结合图 19(a)可知,磨料颗粒主要集中在破损坑四周的表面上并产生冲蚀作用,使破损坑的表面坑径和底部坑径相差较小,说明添加磨料后能够有效的提高射流破岩的精准性。

    图  21  不同磨料浓度下底部坑径随冲击时间的变化曲线
    Figure  21.  Variation curves of diameter of bottom pit with impact time under different abrasive volume concentrations

    (1)模拟分析了射流冲击岩体过程的破损形态和指标,通过与已有研究结果进行对比,验证了本研究所建模型的有效性。

    (2)文中采用PYTHON语言编写离散后的球体孔隙,并实现随机分布,采用该方法可以对数值模拟中的孔隙特性实现精准的控制。

    (3)模拟研究了孔隙特性和射流参数对岩体破损的影响。孔隙尺寸和孔隙分布情况会显著影响岩体的破损形状;当孔隙度介于0~25%时,相同的孔隙尺寸下,岩体孔隙度越大,岩体破损坑的表面坑径和破损面积越大;当孔隙填充物之间的密度以及动力黏度系数相对接近时,填充物对岩体破损的影响相对较小,其中当填充物为气体时,岩体则更容易沿着射流冲击方向发生破损变形。射流的初始速度和磨料浓度也是影响岩体破损情况的主要因素,射流初速度越大,射流的冲击能力越强,岩体破损坑的表面坑径和破损面积也随之增大;对于多孔岩体而言,冲击过程中的磨料颗粒主要集中在破损坑四周的壁面上,随着磨料浓度的增大,磨料颗粒对破损坑壁面的冲蚀效果越强,坑底部的直径明显增大。

    利益冲突声明/Conflict of Interests:所有作者声明不存在利益冲突。All authors disclose no relevant conflict of interest.
    作者贡献/Authors' Contributions:蒋明镜、李子煜参与数值方法框架建立;蒋明镜、李子煜、李承超、姜朋明参与模拟工况设计;李子煜完成模拟方法数值化和工况模拟;蒋明镜、李子煜、李承超、姜朋明参与论文写作和修改、所有作者均阅读并同意最终稿件的提交。The numerical method framework was disigned by Jiang Mingjing and Li Ziyu. The simulation experiment was designed by Jiang Mingjing, Li Ziyu, Li Chengchao and Jiang Pengming. The numerical method and experimental simulation was completed by Li Ziyu. The manuscipt was drafted and revised by Jiang Mingjing, Li Ziyu, Li Chengchao and Jiang Pengming. All the authors have read the last version of paper and consentedfor submission.
  • 图  1   相平衡条件插值计算结果与实测值对比

    Figure  1.   Comparison between interpolation calculation results and measured values of phase equilibrium line

    图  2   CH4回收率的计算值与实测值

    Figure  2.   Predicted and measured recovery ratios of CH4

    图  3   数值模型及本文与White等[13]数值计算结果对比

    Figure  3.   Model and comparison of predicted results

    图  4   CO2置换开采THMC多场耦合流程图

    Figure  4.   Flow chart of CO2 replacement method

    图  5   本文与室内三轴压缩试验结果对比

    Figure  5.   Comparison of triaxial compression results between experiment and this study

    图  6   本文与室内三轴压缩试验结果对比

    Figure  6.   Comparison of triaxial compression results between experiment and this study

    图  7   数值模型及其边界条件

    Figure  7.   Model and boundary conditions for simulation

    图  8   水合物相中CO2水合物质量分数分布演化

    Figure  8.   Evolution of spatial distribution of xch, GH

    图  9   不同生产方案下的产气图

    Figure  9.   Gas production under different production schemes

    图  10   置换过程中地层孔隙率、固有渗透率分布图

    Figure  10.   Distribution of porosity and inherent permeability during replacement process

    图  11   水合物开采过程中能源土地基p-s曲线

    Figure  11.   p-s curves during hydrate exploitation

    图  12   水合物开采过程中能源土地基有效应力分布图

    Figure  12.   Effective stresses during hydrate exploitation

    图  13   水合物开采过程中胶结强度分布图

    Figure  13.   Distribution of bond strength during hydrate exploitation

    图  14   降压开采载荷板试验过程中胶结强度分布图

    Figure  14.   Distribution of bond strength in plate loading tests during depressurization exploitation

    图  15   置换开采载荷板试验过程中胶结强度分布图

    Figure  15.   Distribution of bond strength in exploitation plate loading tests during CO2 replacement

    表  1   不同相态的主变量[6]

    Table  1   Primary variables for different phase states

    相态 主变量1 主变量2 主变量3 主变量4
    G PG χmG χcG T
    A P ωmA ωcA T
    A+G PG SA χcG T
    I+G PG SI χcG T
    A+H P SA χcA T
    I+H P SI χcH T
    H+G PG SG χcG T
    A+I P SA ωcA ωmA
    A+H+G SG SA χcG T
    A+I+G PG SA χcG SG
    A+I+H P SA χcA SI
    I+G+H SG SI χcG T
    I+H+A+G SG SA χcG SI
    下载: 导出CSV

    表  2   储层和模型相关参数

    Table  2   Related parameters of reservoir and model

    参数 数值
    储层固有渗透率k 1.0×10-12 m2
    最大液相饱和度SmxA 1.0
    毛细进气压力P0 1000 Pa
    孔隙分布指数λ 0.77
    残余水饱和度SirA 0.1
    下载: 导出CSV

    表  3   储层模型相关参数

    Table  3   Parameters used in reservoir model

    参数 数值
    储层固有渗透率k1 7.5×10-12 m2
    上覆层、下伏层固有渗透率k2 7.5×10-14 m2
    孔隙率ϕ 0.4
    液相衰减指数n 3.572
    气相衰减指数nG 2.0
    残余水饱和度SirA 0.3
    残余气饱和度SirG 0.05
    孔隙分布指数λ 0.45
    毛细进气压力P0 105 Pa
    下载: 导出CSV

    表  4   深海能源土本构模型参数

    Table  4   Constitutive model parameters used in analyses

    净砂参数 变量 数值 新增模型参数 变量 数值
    弹性参数 G0 125 状态参量参数 mb 0.04
    ν 0.05 x 0
    临界状态参数 Mcc 1.2 临界状态参数 ae 0.3
    c 0.712 初始胶结强度参数 ap 880
    eref 0.934 bp 1.3
    λ 0.019 cp -0.56
    ξ 0.7 dp 0.12
    屈服面参数 m 0.03 胶结破坏参数 hp 5.0
    塑性模量
    参数
    h0 10
    ch 0.85
    nb 1.0
    剪胀参数 A0 0.5
    nd 3.5
    下载: 导出CSV
  • [1] 周守为, 李清平, 朱军龙, 等. 中国南海天然气水合物开发面临的挑战与思考[J]. 天然气工业, 2023, 43(11): 152-163.

    ZHOU Shouwei, LI Qingping, ZHU Junlong, et al. Challenges and considerations for the development of natural gas hydrates in South China Sea[J]. Natural Gas Industry, 2023, 43(11): 152-163. (in Chinese)

    [2]

    OHGAKI K, TAKANO K, SANGAWA H, et al. Methane exploitation by carbon dioxide from gas hydrates: Phase equilibria for CO2-CH4 mixed hydrate system[J]. Journal of Chemical Engineering of Japan, 1996, 29(3): 478-483. doi: 10.1252/jcej.29.478

    [3]

    SCHODERBEK D, FARRELL H, HOWARD J, et al. ConocoPhillips Gas Hydrate Production Test Final Technical Report[R]. Houston: United States Department of Energy, 2013.

    [4]

    LI Q P, LI S X, DING S Y, et al. Numerical simulation of gas production and reservoir stability during CO2 exchange in natural gas hydrate reservoir[J]. Energies, 2022, 15(23): 8968. doi: 10.3390/en15238968

    [5] 张逸群, 杜红星, 王海柱, 等. 双井周期注CO2联合降压法开采天然气水合物分析[J]. 天然气工业, 2024, 44(3): 199-213.

    ZHANG Yiqun, DU Hongxing, WANG Haizhu, et al. Dual-well cyclic CO2 injection assisted gas hydrate depressurization production[J]. Natural Gas Industry, 2024, 44(3): 199-213. (in Chinese)

    [6]

    KAN J Y, SUN Y F, DONG B C, et al. Numerical simulation of gas production from permafrost hydrate deposits enhanced with CO2/N2 injection[J]. Energy, 2021, 221: 119919. doi: 10.1016/j.energy.2021.119919

    [7]

    TIAN H L, YU Z M, XU T F, et al. Evaluating the recovery potential of CH4 by injecting CO2 mixture into marine hydrate-bearing reservoirs with a new multi-gas hydrate simulator[J]. Journal of Cleaner Production, 2022, 361: 132270. doi: 10.1016/j.jclepro.2022.132270

    [8]

    CHEN G J, GUO T M. A new approach to gas hydrate modelling[J]. Chemical Engineering Journal, 1998, 71(2): 145-151. doi: 10.1016/S1385-8947(98)00126-0

    [9]

    GUO Y, LI S X, HUANG X, et al. Numerical simulation and optimization of well parameters for depressurization-assisted CO2 replacement using a PSO algorithm[J]. Energy & Fuels, 2024, 38(11): 9722-9733.

    [10]

    ZHANG S L, MA Y R, XU Z H, et al. Numerical simulation study of natural gas hydrate extraction by depressurization combined with CO2 replacement[J]. Energy, 2024, 303: 131998. doi: 10.1016/j.energy.2024.131998

    [11] 万春燕, 张贺恩, 李磊, 等. 海洋天然气水合物降压开采装备现状与技术探讨[J]. 石油机械, 2024, 52(10): 83-90.

    WAN Chunyan, ZHANG Heen, LI Lei, et al. Current status and techniques of equipment for depressurization exploitation of marine gas hydrate[J]. China Petroleum Machinery, 2024, 52(10): 83-90. (in Chinese)

    [12]

    MORIDIS G J. User's Manual for the Hydrate v1.5 Option of TOUGH+ v1.5: A Code for the Simulation of System Behavior in Hydrate-Bearing Geologic Media[P]. 2014.

    [13]

    WHITE M D, MCGRAIL B P. Numerical simulation of methane hydrate production from geologic formations via carbon dioxide injection[C]// Offshore Technology Conference. Texas, 2008.

    [14]

    PENG D Y, ROBINSON D B. A new two-constant equation of state[J]. Industrial & Engineering Chemistry Fundamentals, 1976, 15(1): 59-64.

    [15]

    ADISASMITO S, FRANK R J III, SLOAN E D Jr. Hydrates of carbon dioxide and methane mixtures[J]. Journal of Chemical & Engineering Data, 1991, 36(1): 68-71.

    [16]

    BELANDRIA V, ESLAMIMANESH A, MOHAMMADI A H, et al. Compositional analysis and hydrate dissociation conditions measurements for carbon dioxide + methane + water system[J]. Industrial & Engineering Chemistry Research, 2011, 50(9): 5783-5794.

    [17]

    ZANG X Y, LIANG D Q. Phase equilibrium data for semiclathrate hydrate of synthesized binary CO2/CH4 gas mixture in tetra-n-butylammonium bromide aqueous solution[J]. Journal of Chemical & Engineering Data, 2017, 62(2): 851-856.

    [18] 王曦. CO2+N2混合气置换开采天然气水合物实验研究及过程模拟[D]. 广州: 华南理工大学, 2017.

    WANG Xi. Experimental Study and Process Simulation of CO2+N2 Gas Mixture Displacement Exploitation of Natural Gas Hydrate[D]. Guangzhou: South China University of Technology, 2017. (in Chinese)

    [19]

    ZHANG A, JIANG M J, WANG D, et al. A bounding surface plasticity model for methane hydrate-bearing sediments in deep seabed[J]. Computers and Geotechnics, 2023, 163: 105720.

    [20] 车竟畅. 深海能源土边坡开采效率与边坡稳定性多场耦合数值模拟分析[D]. 天津: 天津大学, 2024.

    Che Jingchang. Multi-Field Numerical Analysis on Gas Production and Slope Stability of Methane Hydrate Bearing Sediments[D]. Tianjin: Tianjin University, 2024. (in Chinese)

    [21]

    HYODO M, NAKATA Y, YOSHIMOTO N, et al. Triaxial compressive strength of methane hydrate[C]//ISOPE International Ocean and Polar Engineering Conference. ISOPE, 2002.

    [22] 蒋明镜, 董硕, 韩亮, 等. 结构性深海固碳黏土的离散元模拟方法研究[J/OL]. 岩土工程学报. DOI: 10.1179/CJGE20240141.

    JIANG Mingjing, DONG Shuo, HAN Liang, et al. Discrete element method for structural carbon dioxide hydrate bearing clay[J/OL]. Chinese Journal of Geotechnical Engineering. DOI: 10.1179/CJGE20240141.(in Chinese)

    [23]

    HYODO M, LI Y, YONEDA J, et al. A comparative analysis of the mechanical behavior of carbon dioxide and methane hydrate-bearing sediments[J]. American Mineralogist, 2014, 99(1): 178-183.

    [24]

    YAN C L, CHEN Y, TIAN W Q, et al. Effects of methane-carbon dioxide replacement on the mechanical properties of natural gas hydrate reservoirs[J]. Journal of Cleaner Production, 2022, 354: 131703.

    [25]

    WANG Y, FENG J C, LI X S, et al. Evaluation of gas production from marine hydrate deposits at the GMGS2-site 8, Pearl River mouth basin, South China Sea[J]. Energies, 2016, 9(3): 222.

  • 期刊类型引用(8)

    1. 王泽志,王杰,马小刚,李帆,范新亚,陈燕. 用磁性微磨料射流技术光整加工交叉深孔内壁. 金刚石与磨料磨具工程. 2025(02): 245-255 . 百度学术
    2. 王高辉,孔维伟,卢文波,潘鑫豪,舒奕展. 高拱坝坝后侵彻爆炸毁伤效应分析. 武汉大学学报(工学版). 2024(07): 853-862 . 百度学术
    3. 王义江,郁东旭,孙立鹏,朱启银,王建州. 激光照射岩石热裂特性与裂隙分布研究. 岩土工程学报. 2024(09): 1809-1819 . 本站查看
    4. 韩明兴,余锴,段宏兵,熊利荣,徐琨,李淼,刘琦. 基于SPH-FEM的油菜茎秆磨粒气体射流切割仿真与试验. 农业工程学报. 2024(17): 82-92 . 百度学术
    5. 李新旺,王汉青,程立朝,张学栋,温学君,申鹏飞. 高压水射流冲击下煤体初始失效形式分析. 中国科技论文. 2023(01): 64-70+90 . 百度学术
    6. 胡阳,曹安业,秦续峰,薛成春,郭文豪,刘耀琪,彭雨杰. 深部冲击地压煤层水射流割缝卸压参数优化研究. 岩土工程学报. 2023(07): 1509-1516 . 本站查看
    7. 林琳,赵志磊,蒋东岑,张云朋,尤晖. 尖边喷嘴磨料水射流特性分析及对去除函数轮廓的影响. 液压与气动. 2023(09): 116-129 . 百度学术
    8. 王璐,徐绯,杨扬. 完全拉格朗日SPH在冲击问题中的改进和应用. 力学学报. 2022(12): 3297-3309 . 百度学术

    其他类型引用(9)

图(15)  /  表(4)
计量
  • 文章访问数:  207
  • HTML全文浏览量:  20
  • PDF下载量:  64
  • 被引次数: 17
出版历程
  • 收稿日期:  2024-07-30
  • 修回日期:  2025-02-25
  • 录用日期:  2025-03-14
  • 网络出版日期:  2025-03-12
  • 发布日期:  2025-03-15
  • 刊出日期:  2025-06-30

目录

/

返回文章
返回