Multi-field coupling continuum numerical method for exploiting CH4 by CO2 replacement in deep-sea formation
-
摘要: 通过诱导水合物分解的开采方式将致使含天然气水合物沉积物（后称能源土）力学性质劣化，进而可能引发一系列地质灾害。CO2置换法能有效改善地层强度，但目前鲜有研究考虑CO2置换开采的分析方法。基于此，针对CO2置换法开展数值方法研究和应用。通过拓展TOUGH+HYDRATE数学模型、嵌入Chen-Guo模型计算CH4-CO2多元水合物相平衡条件，建立CO2置换的数值模拟器T+MixH V1.0，并与室内试验和已有数值计算结果对比，验证模拟器的可靠性。随后与FLAC3D结合，建立了CO2置换开采工况下的温-压-力-化多场耦合数值分析方法。最后对比分析了南海储层降压、置换开采CH4水合物两种情景以及经历不同开采时间后的能源土地基载荷板试验，结果表明：相比于降压开采方式，置换开采提高了产气量并减缓了地层沉降，其原因在于降压造成的能量损失使CH4水合物分解范围有限，而置换开采在合适的温压条件下自发进行且生成的CO2水合物对地层有一定支撑作用。此外，能源土地基承载特性主要受控于两种因素：孔压变化引起的有效应力变化，水合物分解/生成引起的胶结强度变化。
-
关键词:
- 深海能源土 /
- 二氧化碳置换法 /
- Chen-Guo水合物模型 /
- 多场耦合模拟 /
- 地基承载力
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. -
0. 引言
近年来伴随城市地下空间与地下轨道交通的建设发展，矩形及类矩形顶管在软土地铁车站等大型地下结构建设的应用越来越多，逐渐向大断面、深埋方向发展[1-2]。对于软土深埋超大断面矩形及类矩形顶管，准确评估其结构响应对工程稳定性及安全性尤为重要，其关键在于明确结构水土荷载分布特征及土-结构相互作用关系。
因此，很多学者通过理论分析、现场监测和数值模拟等手段对顶管外荷载分布和幅值展开了研究。在土压力计算理论方面，现行顶管设计规范[3-5]建议根据覆土深度和土层性质的不同采用土柱理论或太沙基理论计算结构竖向土压力；杨仙等[6]考虑土拱效应和土体黏聚力产生的挟持作用，提出了改进的深埋圆形顶管竖向土压力计算模型；贾连辉[2]通过研究矩形顶管施工过程，指出矩形顶管上覆土较难形成稳定土拱；林越翔等[7]提出了基于普式压力拱理论的类矩形顶管土压力的计算方法。综上，目前已有较多针对矩形及类矩形顶管土压力计算的研究，但多针对小断面顶管，且模型的合理性多通过与施工顶推力实测对比而进行间接验证，缺乏与实测水土荷载的直接对比分析。
在结构响应特性方面，雷晗等[8]研究了顶管施工过程中管-土接触压力的变化及对结构内力的影响；Zhen等[9]分析了轴压作用下环向和正交加劲钢顶管的屈曲变形特性；柳献等[10]推导了顶管隧道纵向力随顶进里程的变化关系，并通过对比现场实测进行了模型验证；刘红波等[11]对大断面矩形顶管工程结构应变进行监测，研究了大断面顶管结构的施工力学特性；潘伟强等[12]通过现场实测，研究了大断面顶管结构在不同工况下的结构受力性能。目前，针对顶管结构受力特性的研究主要集中在圆形顶管上。虽然已有部分学者开始探索大断面矩形顶管的前沿研究，但水土荷载模式及土-结构相互作用对顶管受力影响的规律仍有待进一步深入研究，尤其是在实测数据与理论分析的对比方面。
本文的研究目标是揭示软土深埋大断面类矩形顶管的水土压力分布特征，分析水土荷载模式及土-结构相互作用对类矩形顶管受力响应的影响规律，并提出针对顶管结构的设计建议。为此，依托上海轨道交通14号线静安寺站顶管工程，开展了水土压力及结构响应的实测研究，结合不同的土压力理论和地层抗力模型，进行了顶管结构的受力分析，研究了大断面类矩形顶管的水土压力与结构受力响应的空间分布，并对设计理论的适用性进行了评估。
1. 工程概况
1.1 工程简介
上海市轨道交通14号线静安寺站主体结构埋深约24 m，B区采用暗挖顶管法施工，包括两根站台顶管（如图 1所示）和一根站厅顶管。顶管段按照下行线站台顶管、上行线站台顶管和站厅层顶管的顺序依次顶进施工。本文基于下行线站台顶管的实测数据开展研究。下行线顶管段每环宽度为2 m，共41环，纵向总长度为82 m，管顶覆土深度为15.17~15.37 m。顶管管片采用厚度为400 mm的类矩形钢结构，断面尺寸为8.7 m×9.9 m（见图 1）
![]()
图 1 上海轨交14号线静安寺车站顶管工程剖面图Figure 1. Profile of Jing'an Temple Station Project of Shanghai Metro Line 14顶管工程所涉及的土层为典型上海软土，顶管结构主体穿越的主要土层包括④层淤泥质黏土、⑤1-1层黏土和⑤1-2层粉质黏土，水位埋深为0.8~2.5 m。场地主要土层分层及物理力学参数详见图 1和表 1。
表 1 土层物理力学参数Table 1. Physical and mechanical parameters of soil layers地层名称 h/m γ/
(kN·m-3)c′/
kPaφ′/
(°)k/
(kPa·m-1)①填土 2.0 18.8 10 10.0 — ②1黏土 4.0 18.3 3 31.8 1.67×105 ③淤泥质粉质黏土 9.0 17.3 2 29.2 1.31×105 ④淤泥质黏土 17.5 16.7 2 24.5 1.27×105 ⑤1-1黏土 21.8 17.8 3 29.9 1.73×105 ⑤1-2粉质黏土 30.3 18.0 5 31.9 1.72×105 注：h为土层底部埋深；γ为总重度；c′为有效黏聚力；φ′为有效内摩擦角；k为侧向基床反力系数，按照林华国等[13]选取。 1.2 监测方案
为了研究顶管管节水土荷载分布与结构响应，选择下行线顶管的第8环和第21环作为测试管节（见图 2）。其中，第8环位具有最大埋深，而第21环则是距离临近高架桥群桩基础最近的管节。在下行线顶管施工期间，上行线尚未开始施工。因此，选择下行线顶管开展实测研究，有助于揭示在没有其他施工干扰下，单一大断面类矩形顶管的水土荷载及结构响应特征。
顶管结构水土荷载及应变监测点的典型布设位置如图 3所示。其中，水土压力计共布设在8个测试截面（见图 3（a）中W、T测点），分别分布于管节的顶部、底部、两腰、肩部和趾部。每个测试截面设置1个土压力计（见图 3（a）中T N-1至T N-8）和1个孔隙水压力计（见图 3（a）中W N-1至W N-8）, 其中N为管节环号（8与21）。
![]()
图 3 顶管水土荷载与结构应变测点布置图Figure 3. Configuration of water and earth pressures and strain gauges for pipe jacking structures考虑避免施工等干扰因素，结构应变监测断面选择在管节的顶部、底部、两腰以及右上肩部和左下趾部，共计6处（见图 3（a）中GN测点位置）。图 3（b）展示了应变监测截面的布置图（以顶管腰部的GN-6断面为例）。每个监测截面布置6组钢板应变计，分别位于背板和肋板的内侧，编号为GN-n-1至GN-n-6，其中n（1~6）为断面编号。
1.3 现场测试结果分析
（1）水土压力
图 4展示了顶管第8环和第21环外侧水土压力的典型时程曲线。从中可见，随着顶管从始发井出发、经过洞口加固区、进入未加固原状地层，叠加地层差异、注浆扰动、水土压力重分布等因素的影响，管节的水土荷载随时间发生了显著变化。
![]()
图 4 实测顶管水土压力的典型时程曲线Figure 4. Typical time-history curves of measured water and earth pressures on pipe jacking structures如图 4（a）所示，第8环顶部总压力与水压力的时程变化规律总体保持一致。当管节位于洞口加固区内时，顶部水土总压力维持在100~200 kPa。随着管节进入非加固土层，顶部水土总压力迅速增大至约300~320 kPa，并由于注浆和地层变化等因素出现小幅波动。顶管贯通后，水土总压力开始下降并逐渐趋于稳定，最终稳定在约260 kPa。第8环右腰部的总压力变化趋势（见图 4（b））与顶部水土荷载变化相似，即当管节进入未加固土层后，地层压力达到最大值，而在顶进完成后荷载下降并逐渐稳定。
第21环的水土荷载监测结果（见图 4（c），（d））与第8环基本一致，都在进入非加固区后荷载迅速增大。随着顶管的顶进，水土总压力逐渐减小并趋于稳定。另一方面，第21环顶部的水压力显著低于第8环测试结果。鉴于第21环监测数据显示，水压力并未呈现与总压力同步的局部峰值（即由注浆等施工因素引起的压力突增），本文推测孔压计可能发生了阻塞或密封失效，从而低估了实际的水压力。需要指出的是，后续对顶管结构荷载分布的分析主要基于现场试验中土压力计的测试结果（即水土总压力），相比于水压力计，土压力计的测试结果相对更为稳定。
图 5展示了两个测试环水土总压力值沿结构横断面的分布情况。鉴于顶推过程中可能由于“背土”效应以及注浆等因素引起的荷载短暂性增大，图 5中选取了图 4中最终稳定的水土荷载值。从图 5可见，两环的水土荷载在空间分布上基本一致，随着深度的增加，压力值逐渐增大。在具体数值方面，两环在顶部和左肩的实测压力值基本相同，在左腰、底部以及两侧趾部，第21环的压力值略高于第8环，而在右腰和右肩部位置，第21环的压力值显著偏低，两环的实测土压力差值可达到100 kPa。结合图 1中车站剖面图进行分析，上述显著差异可能是由于21环管节与右侧高架桥群桩基础相距较近。群桩基础对地层运动的约束作用，使得该测点的压力值显著低于第8环。
![]()
图 5 顶管结构水土压力实测值沿横面分布图Figure 5. Distribution of measured water and earth pressures along transverse section of pipe-jacking structures（2）管节内力
图 6展示了顶管管节应变时程关系的典型监测结果，各测点的位置详见图 3。各测点的管节应变与水土压力实测值的变化趋势一致。在顶进初期，应变值较小；当管节进入非加固地层后，应变值逐渐增加并达到峰值，随后逐渐降低并趋于稳定。与管节水土压力时程关系（见图 4）相比，管节应变的峰值具有一定滞后性。具体表现为，当管节进入非加固区后，水土压力会迅速达到峰值，而管节应变则在顶进约12~15 m（相当于6~8个管节）后达到峰值。这一现象表明，各管节在纵向上存在相互作用，管节的应变反映了顶管在一定长度范围内的整体结构响应。
![]()
图 6 实测顶管结构应变的典型时程曲线Figure 6. Typical time-history curves of measured structural strain of pipe-jacking structures由图 6可知，不同测点的应变稳定数值与测点到截面形心轴的距离（见图 3）呈正相关关系，即符合平截面假定。在此基础上，本文依据管节应变计算了结构弯矩，结构如图 7所示（以结构内侧受拉为正）。
![]()
图 7 实测管片弯矩的横断面分布图Figure 7. Distribution of measured bending moments along transverse section从图 7可见，两监测管节的横向弯矩的空间分布基本一致。在管节顶部和底部位置，弯矩较大，并表现为内侧受拉；在两腰部位置，弯矩较小，且为外侧受拉。这种弯矩分布与圆形顶管显著不同。基于弹性理论，上述弯矩分布表明，类矩形顶管的变形模式为竖向内凹、横向整体外凸，并且腰部向内收缩。这一变形模式与刘红波等[11]、邓婷等[14]数值模拟分析结果一致。
图 7中顶管顶部的弯矩值约为600 kN·m，这一数值明显小于一般预期（即考虑顶管顶部较大的水土压力（见图 5））。上述较低的弯矩可能由以下两方面原因导致。首先，结构顶部呈拱形（见图 3），有利于降低结构内力；其次，结构侧向外凸变形会挤压土体并引起地层水平抗力，从而为结构提供额外侧向支撑，导致顶部弯矩的下降。后文将对第二方面原因进行详细分析。
从图 7中可见，结构顶部和趾部的弯矩在两个测试环之间存在较大的差距。其中，两个测试环在趾部的弯矩分别为224，446 kN·m，彼此相差近100%。吴列成等[15]给出了静安寺顶管项目中结构外侧的注浆量。数据显示，第21环的注浆量显著大于第8环。由此推测，第21环周围可能存在土-结构脱开的情况，导致较低的土体抗力和较高的结构弯矩，从而造成上述弯矩差异。
2. 顶管结构受力计算分析
2.1 水土荷载模型
图 8展示了类矩形顶管结构受力分析中考虑的外荷载，包括：顶部竖向水土压力、侧向水土压力、底部基地反力和地层水平抗力。其中，土压力计算采取水土分算方法，即水压力按静水压力计算，而土压力则基于有效强度参数确定。
![]()
图 8 类矩形顶管结构荷载分布图Figure 8. Load types and distribution for quasi-rectangular pipe-jacking structures（1）顶部竖向土压力（图 8中p1与p2）
针对依托工程所处的上海软土地基，常用的竖向土压力模型有土柱理论和太沙基理论。这两种理论也是现行顶管结构设计规范[3-4]中推荐的竖向土压力计算方法。
土柱理论将结构正上方至地表的土体重度的累积作为竖向土压力，计算表达式如下：
p1=n∑i=1γihi。 (1) 式中：γi为土体有效重度；hi为土体厚度；i为土层数。
太沙基理论考虑了土拱效应，认为上覆土自重的部分通过剪应力传递给两侧土体，计算表达式如下：
p2=B1γ−2c2Katanφ[1−exp(−2KatanφHB1)]。 (2) 式中：H为管顶的埋深；c为土体有效黏聚力，φ为土体有效内摩擦角；Ka为主动土压力系数；B1为塌落体宽度，按下式计算：
B1=B0+2H0tan(45∘−φ2)。 (3) 式中：B0为顶管结构断面宽度；H0为断面高度。
（2）侧向土压力（图 8中p3）
顶管顶进过程中，可能会出现土体超挖，同时管土之间摩擦会导致侧方土体扰动并向顶管方向发生位移。鉴于此，并参考现有的顶管结构设计规范[3-4]建议，本文采用朗肯主动土压力理论来计算类矩形顶管的侧向土压力。
（3）底部地基反力（图 8中p4）
底部地基反力需要平衡顶管顶部竖向水土压力与管节自重。为了简化计算，本文假设结构底部的地基反力为均匀分布，计算表达式如下：
p4=p1+G或p4=p2+G (4) 式中：G为管片自重在结构水平投影上的平均值。
（4）地层水平抗力（图 8中p5）
上文中实测结构响应与既有数值分析[11, 14]均表明，矩形及类矩形顶管可能会出现明显的侧向外凸水平变形。本文在管片两侧布设了水平弹簧，从而模拟管片与地基间的相互作用及由此产生的地层水平抗力。图 8展示了水平弹簧位置，即从两侧直线段直至与之相接的圆弧的1/2处。由于针对矩形顶管的相关研究较少，上述地层抗力的施加位置参考了盾构隧道结构受力分析方法（即地层抗力弹簧分布于水平中轴线以上与以下45°范围）[16]，并考虑数值分析[11]中顶管结构水平变形较大的位置。抗力计算表达式如下：
p5=kδ。 (5) 式中：k为侧向基床反力系数，取值见表 1；δ为弹簧压缩值（即管片水平位移值），为未知响应，根据水土荷载作用下顶管结构与水平弹簧反力的静力平衡进行求解。需要特别指出的是，大断面类矩形顶管的水平变形分布模式及其与土层相互作用的研究仍较为有限。为此，出于简化模型的考虑，本文假设侧向基床反力系数k沿竖向为均匀分布。在不引入过多复杂因素与不确定性的前提下，上述简化有利于从定性角度研究结构的侧向变形及由此导致的地层水平抗力对深埋大断面类矩形顶管内力的影响。
2.2 水土总压力实测与理论对比
表 2列出了后续对比分析中所考虑的3种荷载组合模式。其中，荷载模式1和荷载模式2为现有顶管结构设计规范[3-4]中采用的荷载组合，而模式3则进一步考虑了可能存在的地层水平抗力。各荷载的具体计算方式如式（1）~（5）所示。
表 2 荷载组合模式表Table 2. Load combinations used in structural analysis名称 荷载组合 荷载模式1* 全覆盖土层竖向土压力、底部基地反力、
侧向主动土压力、结构自重、静水压力荷载模式2* 太沙基竖向土压力、底部基地反力、
侧向主动土压力、结构自重、静水压力荷载模式3 全覆盖土层竖向土压力、底部基地反力、
侧向主动土压力、结构自重、
静水压力、地层水平抗力注：*《顶管工程设计标准（DG/TJ 08—2268—2019）》推荐。 图 9和表 3对比了基于3种荷载模式计算所得的水土压力理论值与实测值。为便于与实测压力值进行对比，本文将荷载计算结果转化为各点的法向压力。由图 9和表 3可知，基于荷载模式2（太沙基理论）得到的顶部竖向水土压力明显低于实测数值（约为实测压力的83%），其他位置的压力值也存在不同程度的低估。相比之下，基于荷载模式1和3（均采用土柱理论）计算得到的顶部竖向压力与实测值更为一致，这表明顶管结构上方可能没有形成明显的土拱效应。
![]()
图 9 水土压力理论值与实测值对比Figure 9. Comparison between theoretical and measured values of water and earth pressures表 3 水土压力理论值与实测值Table 3. Theoretical and measured values of water and earth pressures测点 实测均值/
kPa误差百分比/% 模式1 模式2 模式3 顶部 265 -0.8 -14.7 -0.8 底部 292 0.3 -13.7 0.3 左肩 207 19.3 6.3 36.2 右肩 253 -5.2 -13.0 5.9 左腰 266 -8.6 -12.0 12.8 右腰 323 -24.8 -27.6 -7.1 左趾 301 -4.0 -12.3 2.0 右趾 307 -9.1 -14.0 0.2 注：误差百分比=（理论值-实测值）/实测值×100%。 从图 9和表 3中可以看出，在管片两腰位置，基于荷载模式1和2（均忽略了地层水平抗力）计算得到的总压力明显小于实测值（约25%）。相比之下，基于荷载模式3计算得到的荷载与实测更为接近（误差约为7%）。在管片肩部与趾部区域，基于荷载模式3得到的计算结果同样更接近实测数据。上述对比表明，考虑地层水平抗力有助于更准确地预测顶管结构的侧向荷载分布与大小。由于缺乏顶管水平位移实测数据，无法明确地层水平抗力的区域及大小，后续2.4节将通过数值模拟对此进行进一步分析。
最后，图 9展示了荷载模式3计算所得的荷载分布在顶管肩部和趾部存在锯齿状跳跃。该不光滑分布可能是由于采用了简化的侧向抗力模型（即均匀分布的弹簧系数），从而无法反映实际顶管结构侧向变形的趋势，即腰部最大而在肩部与趾部逐渐衰减至零。需要指出，该跳跃突变仅占总压力的5%，因此预期不会对荷载分布和结构响应产生显著影响。
2.3 结构内力实测与理论对比分析
在上一节水土压力对比的基础上，本节将进一步对比分析顶管横断面弯矩的实测与计算结果。鉴于上一节的分析结果，太沙基理论显著低估了顶管顶部和底部的水土压力，本节将重点分析基于荷载模式1和模式3（均采用土柱理论）所得的结构内力及其与实测的对比。其中，荷载模式3和荷载模式1的区别是前者进一步考虑了结构侧向变形引起的地层水平抗力（如图 8所示）。
如图 10和表 4所示，通过荷载模式1和3计算所得的弯矩分布趋势均与实测结果较为一致。具体表现为，顶管的顶部和底部位置有较大的正弯矩，肩部和趾部有较大的负弯矩，而腰部弯矩较少，整体呈现空间“蝴蝶状”分布。从定量角度，基于荷载模式3计算得到的结构内力与实测值更吻合。与之相比，荷载模式1显著高估了结构内力，特别是在顶管两腰位置，平均误差达到247%。因此，针对依托工程中大断面类矩形顶管，考虑结构侧向水平变形及由此诱发的地层抗力对准确计算结构内力具有关键影响。
![]()
图 10 弯矩理论值和实测值对比Figure 10. Comparison between theoretical and measured values of bending moments表 4 弯矩理论值与实测值Table 4. Theoretical and measured values of bending moment测点 实测均值/
(kN·m)误差百分比/% 模式1 模式3 顶部 588 36.9 5.4 底部 679 15.5 -11.5 右肩 -419 1.2 1.7 右腰 -66 284.8 -106.1 左腰 -82 209.8 -104.9 左趾 -335 66.6 24.5 需要注意的是，不同荷载模式计算得到的结构顶部弯矩均大于底部，而实测结果中顶部弯矩略小于底部的数值（底部弯矩约为顶部的1.1~1.2倍）。上述差异可能源于本文的结构受力分析将结构顶部土压力和底部反力均简化为均布模式，因此有必要在后续研究中对此进行进一步探讨。
2.4 结构变形数值模拟
鉴于现场实测中缺少顶管结构变形数据，为了验证是否可能存在上述地层水平抗力，本节开展了数值模拟。针对依托工程项目，使用Plaxis软件，建立了含有顶管结构和地层的二维对称有限元模型。具体的土层划分、结构尺寸以及单元网格如图 11所示。土体采用HS硬化土模型，土层分布自上而下依次为上海②、③、④和⑤层软土，各土层参数根据表 1和王卫东等[17]的研究选取。顶管管片采用板单元模拟，并参考Chen等[18]建议的方法，通过顶管结构边界收缩模拟由于土体超挖、结构变形等引起的土压力重分布。
图 12展示了有限元数值模拟给出的顶管结构的水平变形响应。结果表明，顶管肩部和趾部之间的区域存在较为显著的水平外凸变形，与上文的荷载结构法分析结果具有较好的一致性。这也解释了考虑地层水平抗力后，计算分析能够更合理地反映现场实测获得的顶管结构在腰部、顶部及底部的弯矩分布特征。
![]()
图 12 顶管结构相对水平变形Figure 12. Relative horizontal deformations of pipe-jacking structures3. 结论
本文依托上海市轨道交通14号线静安寺站顶管工程，结合现场实测和理论分析，研究了软土大断面类矩形顶管的水土荷载与结构受力的空间分布特征，并对比分析了不同土压力理论与土-结构相互作用模式对顶管结构内力响应的影响。主要得出以下3点结论。
（1）结构顶部竖向土压力的实测值接近土柱理论的计算值，且基于土柱土压力计算的结构弯矩与实测值较为吻合，表明本工程中顶管顶部的土拱效应可能较弱。
（2）类矩形顶管的弯矩呈现出“蝴蝶状”分布，表明结构变形模式为竖向内凹及横向外凸，其中横向外凸可进一步诱发地层水平抗力。
（3）地层水平抗力对大断面类矩形顶管结构的内力可产生显著影响，而现有设计规范往往忽略这一因素，可能导致结构弯矩被严重高估（在依托工程中，顶管腰部弯矩的计算误差可达2倍）。
最后有必要指出，尽管本文实测数据与理论分析表明软土大断面类矩形顶管上方的土拱效应可能不显著，但考虑到顶管水土荷载受多种因素及施工状态影响的复杂性，上述初步结论仍需通过进一步的实测和理论分析加以检验。
利益冲突声明/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
![]()
图 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
![]()
图 10 置换过程中地层孔隙率、固有渗透率分布图
Figure 10. Distribution of porosity and inherent permeability during replacement process
![]()
图 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
净砂参数 变量 数值 新增模型参数 变量 数值 弹性参数 125 状态参量参数 0.04 0.05 x 0 临界状态参数 1.2 临界状态参数 0.3 c 0.712 初始胶结强度参数 880 0.934 1.3 0.019 -0.56 0.7 0.12 屈服面参数 m 0.03 胶结破坏参数 5.0 塑性模量
参数10 0.85 1.0 剪胀参数 0.5 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.
-
其他相关附件
计量
- 文章访问数: 207
- HTML全文浏览量: 20
- PDF下载量: 64
