基于裂缝网络模型的CO₂地质封存稳定性评价

Stability Assessment of CO₂ Geological Storage Based on a Fracture Network Model

摘要

裂缝作为CO₂地质封存的主要储集空间和渗流通道，直接影响CO₂的封存效率及长期封存的安全性.基于无人机倾斜摄影技术，构建了鄂尔多斯盆地裂缝性砂岩储层离散裂缝网络模型，并基于多物理场耦合仿真软件COMSOL Multiphysics建立了考虑基质‒层理‒裂缝的CO₂-水两相流固耦合数值模型.研究发现，CO₂优先沿高渗透性的层理和裂缝运移，横向层理与低倾角、低连通性的自然裂缝阻滞垂向渗流，降低了CO₂逃逸到盖层的风险；裂缝网络加速了压力传递，诱发显著位移响应，初始位移增长率是基质模型的6.2倍.因此，考虑基质‒层理‒裂缝多重介质系统对CO₂地质封存的稳定性具有重要意义.

Abstract

Fractures serve as the primary storage space and seepage pathways for CO₂ geological storage, directly influencing storage efficiency and long-term containment security. This study employed UAV oblique photogrammetry to construct a discrete fracture network (DFN) model of a fractured sandstone reservoir within the Ordos Basin. Subsequently, a fluid-solid coupling numerical model for CO₂-water two-phase flow was established using the multiphysics simulation software COMSOL Multiphysics, explicitly accounting for the matrix-bedding-fracture system. Key findings reveal that CO₂ preferentially migrates along high-permeability bedding planes and fractures. Horizontal bedding, combined with low-dip, low-connectivity natural fractures, impedes vertical seepage, thereby reducing the risk of CO₂ escape into the caprock. The fracture network accelerates pressure transmission, inducing significant displacement responses. The initial rate of displacement increase was found to be 6.2 times higher than that observed in the matrix model. Consequently, accurately representing the multi-porosity system encompassing the matrix, bedding, and fractures is crucial for assessing the stability of CO₂ geological storage.

Graphical abstract

关键词

CO₂地质封存 / 裂缝网络模型 / 封存潜力 / 环境科学.

Key words

CO₂ geological storage / fracture network model / storage potential / environmental science

随着全球工业化持续推进，化石燃料的消耗量大幅上升，导致大量二氧化碳排放，加剧温室效应，对环境和生态系统造成严重危害.减少二氧化碳排放、实现碳中和成为全球共识.我国在2020年提出了“碳达峰、碳中和”的双碳目标，正式开启了社会经济全面绿色低碳转型的新纪元（Cai et al.，2023）.碳捕集、利用与封存（CCUS）技术是实现CO₂深度减排和可持续发展的重要途径（Rui et al.，2025）.其中，二氧化碳地质封存（CO₂ Geological Storage）是CCUS技术链条中的核心环节：将捕集的CO₂压缩至超临界状态（温度>31.1 ℃、压力>7.38 MPa），兼具气体的流动性与液体的高密度（Jiang et al.，2025），然后通过管道注入深部咸水层、枯竭油气藏等地质构造层中进行长期封存，从而避免其进入大气层加剧温室效应（廖志伟等，2024）.

受沉积作用和构造作用的影响，地下储层广泛发育层理、节理、断层等多尺度裂缝，形成复杂裂缝网络系统，对CO₂地质封存产生重要影响.一方面，裂缝网络提升了储层的封存潜力与CO₂注入效率（Li et al.，2023）.这是由于裂缝网络通过连通孔隙率较低的岩石基质孔隙（仅2%~5%），有效扩大了CO₂的储集空间.此外，裂缝的渗透率（10⁻³~ 10⁻¹ μm²）较岩石基质的渗透率（10⁻⁶~10⁻³ μm²）高2~4个数量级，裂缝贡献了储层绝大部分渗透率，形成了优势渗流通道，从而降低了注入压力并提高了注入效率（杨圣奇等，2022， 2023）.另一方面，裂缝网络系统对CO₂封存安全性与稳定性构成威胁.当裂缝网络与盖层（封存的关键屏障）连通，或存在延伸至浅部含水层/地表的断层时，可能诱发CO₂泄漏（Sun et al.，2024）.例如，美国Cranfield油田曾因裂缝‒断层导通，导致CO₂垂向运移800 m侵入浅层含水层.此外，裂缝的分布（密度、走向、连通性）受构造运动、岩性差异影响显著，导致储层渗透率在空间上剧烈变化，裂缝空间的非均质性会导致CO₂窜流（即沿高渗透裂缝快速突破），使其偏离目标封存区并增加监测难度.裂缝的高渗透性会加速储层压力传递，若注入速率过高，可能引发裂缝进一步扩展或盖层剪切破裂（ΔP>0.7 σ_v时库仑失效），破坏封存结构（Zheng et al.，2025）.因此，鉴于裂缝网络普遍存在于地质环境中且主导烃类运移‒圈闭过程，评估裂隙储层中CO₂封存的过程至关重要.

目前，针对裂缝型储层碳封存的研究已经取得了一定的成果.在实验研究方面，通过室内岩心驱替和现场中试（如Cranfield、Otway项目），对CO₂注入性、盖层密封性及存储机理进行了研究，揭示了裂缝系统对二氧化碳存储的重要作用.例如，研究发现超临界CO₂压裂增产技术能够利用其低粘度高扩散性，进入更小的微裂缝，有效扩展微裂缝网络，激活次生储集空间（基质波及效率提升15%~40%），并将返排率抑制至5%以下（Han et al.，2025）.在数值模拟方面，建立了多种模型来模拟CO₂在裂缝型储层中的流动、扩散以及与岩石的相互作用等过程，为碳封存的方案设计和优化提供了理论支持（袁哲等，2025）.关键迁移机理表明，CO₂（无论是超临界状态、气态还是溶解态）具有流动性，在储层内会发生迁移，其动态分布受控于储层内渗透率级差（裂缝渗透率较基质高2~4个数量级），尤其在低渗透储层中，CO₂会优先驱替裂缝网络中的水，这些裂隙成为流体运移的主要通道和封存场所.此外，层理诱导的渗透率各向异性使得渗透率在平行和垂直层理方向的存在明显差异，影响了储层的非均质性，对CO₂的流动和分布产生影响.

将CO₂储存于深层地质构造层中涉及复杂的物理过程，包括CO₂⁃水两相渗流、CO₂在矿物界面的竞争吸附以及超临界二氧化碳（scCO₂）诱导的岩石基质和裂隙的变形等（Ren et al.，2017）.而在地质构造中裂缝与断层网络广泛发育，这些不连续的构造强化了地质结构的非均质性，加剧了CO₂封存的复杂性，对CO₂的运移和圈闭起着关键作用.Hyman et al. （2020）基于三维离散裂缝网络模型（DFN）开展多相流模拟，探究scCO₂在静水压力条件下的驱替行为，并对比水力特性和结构非均质性对其在盖层初始运移的影响，分析表明裂隙网络结构（连通性、悬挂裂隙）比局部水力特性对scCO₂驱水控制作用大.Bigi et al. （2013）提出了一种模拟裂隙岩石模型并进行尺度提升的类比模型，利用商业软件COMSOL研究了CO₂通过裂隙网络的泄漏情况.Lee and Ni（2015）开发了一种网格生成器，为TOUGH系列模型生成三维离散裂隙网格，用于研究向裂隙性地层中注入CO₂的情况.Oh et al. （2013）通过实验和数值模拟相结合的方式，研究了岩心尺度下向裂隙性岩石样品中注入CO₂的行为.Ren et al.（2017）提出了两相流管网法用于研究裂隙网络对含水层中CO₂运移的影响，发现裂隙可显著影响含水层中 CO₂的分布和演化.Wang et al. （2019）采用综合方法，融合有限元法（FEM）和统一管网方法（UPM），模拟干热岩储层中二氧化碳地质封存和地热开采过程.Wang et al. （2020）运用通用油藏模拟器CMG⁃GEM，评估CO₂通过高渗透裂隙区泄漏的风险.Sun et al. （2024）采用扩展数值流形法来模拟CO₂在含两组裂隙的盖层中的运移，分析了裂缝密度和连通性对CO₂地质封存的影响.Qin et al. （2025）建立页岩储层裂隙中CO₂⁃水‒油流动的孔隙尺度模型，注入二氧化碳后进行注水，由于毛细管效应抑制了二氧化碳窜流，从而提升采收率与封存效率.然而，当前主流的CO₂地质封存模型普遍对储层裂缝系统进行了高度简化，未能充分考虑层理与天然裂缝网络交互效应，导致封存容量等关键参数的预测与实际情况存在显著偏差，严重制约了评估结果的可靠性.

鄂尔多斯盆地作为中国实现“双碳”目标的核心战略区域，其二氧化碳地质封存技术已通过示范工程（10⁵ t/a的CO₂捕集与封存）验证了深部咸水层封存的可行性，并展现出巨大的封存潜力.盆地内的主力封存层以低孔渗砂岩和碳酸盐岩为主，包括奥陶系马家沟组、三叠系延长组、二叠系石盒子组/石千峰组等（付金华等，2020）.然而，盆地内普遍存在的低孔渗、高非均质性储层以及天然裂缝发育等复杂的地质结构，导致CO₂运移与矿物反应存在显著不确定性，制约了技术的规模化应用.本文基于无人机倾斜摄影技术，构建了鄂尔多斯盆地延长组裂缝性砂岩储层高精度地质结构模型，并基于多物理场耦合仿真软件COMSOL Multiphysics，建立了基质‒层理‒裂缝耦合的CO₂⁃水两相流固计算模型.通过精细化模拟CO₂在复杂裂隙网络中的运移过程，揭示了层理结构控制CO₂横向扩散范围、裂缝系统主导垂向运移速率的规律，量化了CO₂扩散对盖层密封效率的削弱效应.本研究为裂缝型储层封存安全评估提供了理论依据与预测工具.

1 地质背景与裂缝表征

1.1　地质背景

鄂尔多斯盆地属于华北克拉通基础上发育的大型多旋回叠合盆地，是中国第二大沉积盆地，现可划分为伊盟隆起、渭北隆起、西缘冲断带、晋西挠褶带、天环坳陷、伊陕斜坡六大构造单元，呈现西深东浅、西陡东缓的不对称向斜结构（图1a）.作为稳定的克拉通盆地，鄂尔多斯盆地构造活动长期较弱，这使得其内延长组地层未受到强烈的构造运动改造，得以保持较好的完整性.因此，潜在的二氧化碳泄漏通道较少，为二氧化碳长期地质封存提供了稳定的地质环境.

三叠系延长组是鄂尔多斯盆地最重要的含油层系，其沉积经历了大型内陆湖盆从形成、鼎盛到消亡的全过程.根据沉积特征，延长组自下而上划分为长10至长1十个储层组（图1b），其砂岩储层在整个鄂尔多斯盆地内分布较为广泛，具有良好的横向连续性.其中，长6、长7、长8段等油层组以细砂岩为主，孔隙度普遍为 8%~15%，渗透率为0.1~ 5 mD，属典型低渗储层（宿晓岑等，2023）.但广泛发育的天然裂缝系统（密度2~8条/m）显著增加了有效孔隙空间（裂缝孔隙度贡献达 0.5%~1.2%），为CO₂提供了有效的赋存空间.同时，长7段页岩连续分布（厚度 40~120 m，突破压力 >8 MPa），且富含伊蒙混层矿物（含量>35 %），遇CO₂⁃水作用后自愈合能力增强，可有效阻隔CO₂垂向迁移，构成垂向封盖层.凭借砂泥岩互层的储盖组合、稳定的构造背景及上述地质特性，延长组成为盆地内理想的CO₂封存靶层.其低渗基质中的裂缝系统通过科学调控可转化为封存增效通道，而矿物固碳潜力进一步保障了长期封存安全性，为规模化碳封存工程提供了不可替代的地质载体.然而，储层强烈的非均质性，尤其是层理与裂缝网络交织形成的复杂结构，显著增加了CO₂在地下空间分布与运移路径的不确定性.这种不确定性直接关系到封存体的密封性、有效容量利用及长期安全性，因此，深入探究CO₂在裂缝性储层中的运移规律与主控因素，并据此优化注入参数以精准调控其分布范围，是实现该区域高效、安全封存亟待解决的核心科学问题.

1.2　裂缝表征

本文研究对象为位于宜山斜坡构造单元上三叠系延长组长7地层剖面.剖面主体岩性为砂岩，厚度为22 m，砂岩上部和下部均为黑色泥岩.基于大疆幻影4实时动态定位（RTK）无人机平台开展裂缝测绘工作（图2）.该平台搭载视场角84°的航摄镜头，配合1英寸CMOS影像传感器提供2 000万有效像素.RTK系统采用双备份高精度全球导航卫星系统（GNSS），实现垂直方向1.5 cm、水平方向1 cm的定位精度.通过多航次重叠扫描获取319张高分辨率影像（5 472×3 648像素），并基于影像数据构建露头数字照片.

通过无人机采集的图片，在Engineering Profile System软件中通过人工识别数字模型中的裂缝，并提取天然裂缝上具有位置信息的特征点.理论上，通过选取裂缝面上任意3个非共线特征点即可表征裂缝平面形态.本文通过MATLAB编写程序计算节理的产状、迹长等参数.计算原理如下.

将裂隙视为平面，其平面方程为：

Z=AX+BY+C，(1)

该平面的法向量为：

\mathit{n}=\left(-A,-B,\mathrm{1}\right)，(2)

假如在某一裂隙上共测得n个 （n≥3）不共线特征点P_i （X_i， Y_i， Z_i ），通过最小二乘法可求解出A、B、C （Liu et al.，2022）：

\left[\begin{array}{l}A\\ B\\ C\end{array}\right]={\left({\left[\begin{array}{ccc}{X}_{\mathrm{1}}& Y_{\mathrm{1}}& \mathrm{1}\\ \begin{array}{l}{X}_{\mathrm{2}}\\ ⋮\end{array}& \begin{array}{l}{Y}_{\mathrm{2}}\\ ⋮\end{array}& \begin{array}{l}\mathrm{1}\\ ⋮\end{array}\\ X_n& Y_n& \mathrm{1}\end{array}\right]}^{\mathrm{T}}\left[\begin{array}{ccc}{X}_{\mathrm{1}}& Y_{\mathrm{1}}& \mathrm{1}\\ \begin{array}{l}{X}_{\mathrm{2}}\\ ⋮\end{array}& \begin{array}{l}{Y}_{\mathrm{2}}\\ ⋮\end{array}& \begin{array}{l}\mathrm{1}\\ ⋮\end{array}\\ X_n& Y_n& \mathrm{1}\end{array}\right]\right)}^{-\mathrm{1}}{\left[\begin{array}{ccc}{X}_{\mathrm{1}}& Y_{\mathrm{1}}& \mathrm{1}\\ \begin{array}{l}{X}_{\mathrm{2}}\\ ⋮\end{array}& \begin{array}{l}{Y}_{\mathrm{2}}\\ ⋮\end{array}& \begin{array}{l}\mathrm{1}\\ ⋮\end{array}\\ X_n& Y_n& \mathrm{1}\end{array}\right]}^{\mathrm{T}}\left[\begin{array}{l}{Z}_{\mathrm{1}}\\ Z_{\mathrm{2}}\\ ⋮\\ Z_n\end{array}\right],(3)

该结构面的倾向β与倾角α（图3）由法向量 n 来计算：

（1）A=0时，

，(4)

（2）A≠0时，

.(5)

裂隙迹线长度被定义为同一节理上任意两个特征点相连所组成的线段中的最大线段长度（图4）（Kong et al.，2020）， 可通过欧氏距离公式计算得到：

d_{ij}=\sqrt[]{{\left(X_i-X_j\right)}^{\mathrm{2}}+{\left(Y_i-Y_j\right)}^{\mathrm{2}}+{\left(Z_i-Z_j\right)}^{\mathrm{2}}}，(6)

式中： i=1，2，3，…，（n-1），j=i+1，i+2，…，n.

从数字模型上一共提取了223个节理和10个层理数据（图5）.节理划分为两组（Sets 1 and 2），其中，Set 1包含75条节理，平均产状为81°/89°；Set 2有148条节理，平均产状为170°/87°（图6a）.使用Kolmogorov⁃Smirnov （KS）检验进行统计分析，发现两组节理的迹长均服从指数分布（图6b）.

本研究基于地表露头数据构建裂缝网络模型.尽管卸荷、风化等后期表生改造作用会使地表岩体发育次生裂隙，导致其密度与连通性等同深部储层存在差异，但本研究在裂缝识别阶段已对此类次生裂隙进行了有效辨识与剔除.由于地表露头本质上是深部地层在地表的延伸，其所保留的裂缝网络主要继承了区域构造应力场历史与岩性组合的控制.因此，经筛选后的地表裂缝在优势方位、发育组数与基本形态等核心结构信息上，与深部储层具有显著的同源性.鉴于此，利用地表数据所建立的模型能够反映储层地质结构特征.本研究旨在揭示“基质‒层理‒裂缝”复杂多重介质系统中CO₂运移、压力传递及地质力学响应的内在机理，所探究的相对关系（如裂缝渗透率>基质渗透率）与变化趋势具有一定的普适性，对于理解深部流体潜在的优势运移方向具有一定的参考价值.

2 物理模型与数学模型

2.1　物理模型

本研究基于无人机勘测数据构建了地质结构模型，并据此建立了储盖层系统的CO₂地质封存二维数值模型.该模型自上而下由泥岩盖层、砂岩储层和低渗透层组成，其中砂岩储层包含岩石基质、层理和天然裂缝等多重介质.储层埋深为2 021 m，储层温度约为75 °C，模型尺寸简化为40 m×42 m，采用有限元网格离散，共包含41 942个域单元和 3 407个边界单元，如图7所示.模型上下边界为无流动边界，左右边界为静水压力边界.储层初始饱和度为1，初始压力为22 MPa；注入井位于模型中心，定流量注入，模拟方案包括25 d注入期和25 d扩散期.由于本研究区域为“冷盆地”特征，储层埋深为2 021 m，因此CO₂为超临界状态（Nordbotten et al.，2005）.为探究裂隙系统对CO₂运移的影响，本研究设计了三种对比模型，分别为岩石基质模型、层理模型、层理和裂隙模型.模型建立基于以下基本假设：（1）储层为咸水‒超临界CO₂两相体系，流体运移遵循达西定律；（2）暂不考虑长期地质封存过程中的溶解封存和矿物封存；（3）系统温度保持恒定；（4）气体吸附服从Langmuir吸附

2.2　数学模型

2.2.1　基岩两相渗流方程

COMSOL Multiphysics是一款基于有限元方法的多物理场耦合仿真软件，其多孔介质多相流模块（Porous Media and Subsurface Flow Module）可有效模拟气‒水两相流动过程.超临界CO₂⁃盐水在储层中的运移满足质量守恒方程和动量守恒方程，构成了一个耦合的非线性方程组.通过配合初始条件与边界条件可实现数值求解.sc⁃CO₂相和盐水相的质量守恒方程为：

\frac{\partial }{\partial t}\left(\phi {\rho }_i{S}_i\right)+\nabla \left({\rho }_i{u}_i\right)=Q_i，(7)

式中：φ为多孔介质的孔隙度，无量纲；S_i 为i相的饱和度，无量纲；ρ_i 为i相的密度，kg/m³；Q_i 为i相的相变项，本文忽略sc⁃CO₂和盐水的相变，因此Q_i=0；u_i 为i相的达西速度，m/s.

达西定律常被应用于研究单相流体在多孔介质中的流动，但不适用于两相流体.孔祥言等（1999）提出了扩展的达西定律方程，通过引入相对渗透率描述相流动能力随饱和度的变化，结合毛管压力刻画相间压力差异，形成了适用于多相流的流动方程.流体流动方程由达西定律表示：

u_i=-\frac{k{k}_{ri}}{{\mu }_i}\nabla \left(P_i-{\rho }_i\mathrm{g}\right)，(8)

式中：k为多孔介质渗透率，m²；k_ri 为i相的相对渗透率，无量纲；μ_i 为i相的动力粘度，Pa∙s；P_i 为i相压力，MPa；g为重力加速度，m/s².

将（8）式代入（7）式，可以得到两相流平衡方程：

\frac{\partial }{\partial t}\left(\phi {\rho }_i{S}_i\right)-\nabla \left[{\rho }_i\frac{k{k}_{ri}}{{\mu }_i}\nabla \left(P_i-{\rho }_i\mathrm{g}\right)\right]=\mathrm{0}.(9)

本文假设储层中只含有水和sc⁃CO₂两相（i=w，n），因此水和sc⁃CO₂共同完全地填充了孔隙空间，两相饱和度之和为1，储层的总饱和度方程为：

S_{\mathrm{w}}+S_{\mathrm{n}}=\mathrm{1}.(10)

毛细压力P_c可以被定义为两相压力之差：

P_{\mathrm{c}}=P_{\mathrm{n}}-P_{\mathrm{w}}.(11)

考虑到咸水层中封存的CO₂部分以吸附态存在于岩石骨架中，因此封存在储层中的CO₂总质量为CO₂自由相质量和CO₂吸附相质量之和（盛丹娜等，2025）.因此，单位总体积的CO₂总质量变为：

m_{\mathrm{n}}=\phi {\rho }_{\mathrm{n}}{S}_{\mathrm{n}}+\left(\mathrm{1}-\phi \right){\rho }_{\mathrm{c}}{\rho }_{\mathrm{n}\mathrm{s}}\frac{V_{\mathrm{L}}{P}_{\mathrm{n}}}{P_{\mathrm{n}}+P_{\mathrm{L}}}，(12)

式中：P_L是Langmuir吸附压力，MPa；V_L是Langmuir吸附体积，m³/kg.

毛细压力P_c与饱和度的关系采用Van⁃Genuchten模型，平衡方程中主变量为盐水相的饱和度S_w和压力P_w.因此，将上述公式联立得到：

盐水相平衡方程：

\frac{\partial }{\partial t}\left(\phi {\rho }_{\mathrm{w}}{S}_{\mathrm{w}}\right)-\nabla \left[{\rho }_{\mathrm{w}}\frac{k{k}_{\mathrm{r}\mathrm{w}}}{{\mu }_{\mathrm{w}}}\nabla \left(P_{\mathrm{w}}-{\rho }_{\mathrm{w}}\mathrm{g}\right)\right]=\mathrm{0}.(13)

sc⁃CO₂相平衡方程：

\begin{array}{l}\frac{\partial }{\partial t}\left(\phi {\rho }_{\mathrm{n}}{S}_{\mathrm{n}}\right)+{\rho }_{\mathrm{n}\mathrm{s}}{\rho }_{\mathrm{c}}{V}_{\mathrm{L}}\left[\left(\mathrm{1}-\phi \right)\frac{P_{\mathrm{L}}}{{\left(P_{\mathrm{n}}+P_{\mathrm{L}}\right)}^{\mathrm{2}}}\frac{\partial P_{\mathrm{n}}}{\partial t}-\frac{P_{\mathrm{n}}}{P_{\mathrm{n}}+P_{\mathrm{L}}}\frac{\partial \phi }{\partial t}\right]\\ -\nabla \left[{\rho }_{\mathrm{n}}\frac{k{k}_{\mathrm{r}\mathrm{n}}}{{\mu }_{\mathrm{n}}}\nabla P_{\mathrm{w}}+\nabla \left({\rho }_{\mathrm{n}}\frac{k{k}_{\mathrm{r}\mathrm{n}}}{{\mu }_{\mathrm{n}}}\cdot \frac{\partial P_{\mathrm{c}}}{\partial S_{\mathrm{n}}}\nabla S_{\mathrm{n}}-{\rho }_{\mathrm{n}}\mathrm{g}\right)\right]=\mathrm{0}.\end{array}(14)

2.2.2　裂隙两相渗流方程

考虑基质‒层理‒裂缝的两相流动，本文利用COMSOL自带的达西定律‒裂缝接口，将物理模型视为离散裂缝网络模型，裂缝看作是渗流区域的内边界，求解裂隙流时可通过对所有裂隙边界增加弱形式的方式表现，其渗透率k_f与开度d_f满足关系：

k_{\mathrm{f}}=\frac{{d_{\mathrm{f}}}^{\mathrm{2}}}{\mathrm{12}}.(15)

其渗流守恒方程转换为：

盐水相平衡方程：

d_{\mathrm{f}}\frac{\partial }{\partial t}\left({\phi }_{\mathrm{f}}{\rho }_{\mathrm{w}}\right)+{\nabla }_{\mathrm{T}}\left(-\frac{{d_{\mathrm{f}}}^{\mathrm{3}}}{\mathrm{12}}{\rho }_{\mathrm{w}}{\nabla }_{\mathrm{T}}{P}_{\mathrm{w}}\right)=d_{\mathrm{f}}{Q}_{\mathrm{w}},(16)

sc⁃CO₂相平衡方程：

d_{\mathrm{f}}\frac{\partial }{\partial t}\left({\phi }_{\mathrm{f}}{\rho }_{\mathrm{n}}\right)+{\nabla }_{\mathrm{T}}\left[-\frac{{d_{\mathrm{f}}}^{\mathrm{3}}}{\mathrm{12}}{\rho }_{\mathrm{n}}\left({\nabla }_{\mathrm{T}}{P}_{\mathrm{w}}+\frac{\partial P_{\mathrm{c}}}{\partial S_{\mathrm{n}}}{\nabla }_{\mathrm{T}}{S}_{\mathrm{n}}\right)\right]=d_{\mathrm{f}}{Q}_{\mathrm{n}}，(17)

式中：k_f为裂缝渗透率，m²；d_f为裂缝开度，m；φ_f为裂缝孔隙率（考虑到裂缝内部无固体骨架，因此φ_f=1），无量纲；层理选择同一接口，只设置渗透率.

2.2.3　储盖层变形控制方程

孔隙率模型：

\phi =\mathrm{1}-\left(\mathrm{1}-{\phi }_{\mathrm{0}}\right)\mathrm{e}\mathrm{x}\mathrm{p}\left[-\left({\epsilon }_{\mathrm{v}}-{\epsilon }_{\mathrm{v}\mathrm{0}}\right)-\frac{\mathrm{1}}{K_{\mathrm{s}}}\left(P-P_{\mathrm{0}}\right)\right]，(18)

式中：φ₀为初始孔隙度，无量纲；ε_v为体积应变，ε_v0为初始体积应变；K_s为岩体骨架的体积模量；P₀为初始地层压力，MPa.其中，固体力学模块中因变量x、y两个方向的位移场分量为u、v，体积应变可表示为：

{\epsilon }_{\mathrm{v}}=\frac{\mathrm{d}u}{\mathrm{d}x}+\frac{\mathrm{d}v}{\mathrm{d}y}.(19)

渗透率模型采用立方定律表示为：

\frac{k}{k_{\mathrm{0}}}={\left(\frac{\varphi }{{\varphi }_{\mathrm{0}}}\right)}^{\mathrm{3}}，(20)

式中：k₀为多孔介质的初始渗透率，m².

在弹性力学中，对于各项同性材料，弹性力学平衡方程为：

\sum _{j=x,y}\frac{\partial {\sigma }_{ij}}{\partial j}+F_i=\mathrm{0}，(21)

式中：σ_ij 为应力张量的分量，Pa；F_i 为体积力在i方向的分量，N/m³.

根据广义的胡克定律，应力与应变的关系为：

{\sigma }_{ij}=\lambda {\epsilon }_{\mathrm{v}}{\delta }_{ij}+\mathrm{2}G{\epsilon }_{ij}.(22)

式中：λ为拉梅常数；δ_ij 为kronecker delta符号（当i=j时，δ_ij =1，当i≠j时，δ_ij =0）；G为剪切模量，GPa.

应变张量与位移的关系可表示为：

{\epsilon }_{ij}=\frac{\mathrm{1}}{\mathrm{2}}\left(u_{i,j}+u_{j,i}\right).(23)

考虑孔隙流体压力的影响，得到有效应力与压力的关系为：

\sigma {\text{'}}_{ij}={\sigma }_{ij}-\alpha \cdot {\delta }_{ij}\cdot \left(\mathrm{\Delta }{P}_{\mathrm{w}}+\mathrm{\Delta }{P}_{\mathrm{n}}\right)，(24)

式中：α为Bito系数.

考虑到吸附会导致介质发生膨胀或收缩变形，该吸附应变同样满足广义胡克定律，通过计算，得到的吸附体积力为：

{F_i}^{\mathrm{a}\mathrm{d}\mathrm{s}}=-\left(\lambda +\mathrm{2}G\right)\frac{\partial }{\partial x_i}\left({\epsilon }_{\mathrm{L}}\frac{P}{P+P_{\mathrm{L}}}\right).(25)

将式（22）~式（25）代入式（21）中，得到x，y方向上的应力场平衡方程：

\begin{array}{l}\lambda \frac{\partial {\epsilon }_{\mathrm{v}}}{\partial x}+G\left(\mathrm{2}\frac{{\partial }^{\mathrm{2}}u}{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}u}{\partial y^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}v}{\partial x\partial y}\right)-\\ \frac{\partial \left[\alpha \cdot \left(\mathrm{\Delta }{P}_{\mathrm{w}}+\mathrm{\Delta }{P}_{\mathrm{n}}\right)\right]}{\partial x}+{F_x}^{\mathrm{a}\mathrm{d}\mathrm{s}}=\mathrm{0}\end{array}，(26)

\begin{array}{l}\lambda \frac{\partial {\epsilon }_{\mathrm{v}}}{\partial y}+G\left(\frac{{\partial }^{\mathrm{2}}v}{\partial x^{\mathrm{2}}}+\mathrm{2}\frac{{\partial }^{\mathrm{2}}v}{\partial y^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}u}{\partial x\partial y}\right)-\\ \frac{\partial \left[\alpha \cdot \left(\mathrm{\Delta }{P}_{\mathrm{w}}+\mathrm{\Delta }{P}_{\mathrm{n}}\right)\right]}{\partial y}+{F_y}^{\mathrm{a}\mathrm{d}\mathrm{s}}=\mathrm{0}\end{array}.(27)

用于模拟的具体参数见表1.

3 模型验证

本研究基于Pavan et al.（2024）的模型进行验证.该模型为长5 000 m、宽50 m的矩形.模型介质被视为各向同性且均质，初始时被润湿盐水相饱和.在模型左侧边界以恒定速率注入二氧化碳，模型的上、下边界均为无流动边界.图8展示了使用相同参数下，利用本文研究方法和Pavan等人的方法得到的CO₂体积分数在上边界变化.在初始阶段，二氧化碳进入储层，将井周围的水置换出去.随后，二氧化碳继续向远处迁移，并在重力作用下不断在储层上部边界积累.二氧化碳饱和度从上边界到下边界逐渐减少，并且从注入边界到二氧化碳羽流前沿也逐渐减少.但由于本文未考虑温度变化对CO₂运移的影响，因此体积分布图存在轻微差异.此外，1 500 d后注入的模拟二氧化碳羽流与Pavan等人预测的Caseii（a）进行了对比，如图9a和图9b所示，两种结果的二氧化碳羽流空间分布相似，储层上边界处二氧化碳迁移的距离分别为601 m和660 m.结果表明，采用多孔介质两相流动模型在模拟地下储层中的二氧化碳迁移过程方面具有足够的可靠性.

4 结果

本文通过数值模拟研究了层理和裂缝等关键参数对CO₂运移及表面压力的影响.基于模拟结果，定量分析了CO₂注入总量及扩散至盖层的体积分数，并通过计算逃逸比，评估了层理和裂缝对CO₂封存效率及长期安全性的影响.

4.1　层理和裂缝对CO₂饱和度的影响

图10a展示了CO₂在岩石基质模型中的气体饱和度分布演化过程.注入初期，CO₂在注入井附近聚集，形成局部高饱和度区域（S_g>0.2）.注入至t= 25 d时，CO₂羽流范围显著扩展.CO₂的运移主要受以下因素控制：（1）下部低渗透层（K<1 mD）对CO₂垂向运移的阻隔作用；（2）上部泥页岩盖层的高毛细管压力导致突破压力低于临界值.在此条件下，CO₂优先沿水平方向扩散，并在浮力作用下向储层顶部富集，形成顶部（S_g=0.3±0.05）显著高于底部（S_g=0.25±0.05）的饱和度分布特征.停止注入后的25 d监测期内，CO₂饱和度呈现系统性降低趋势，其原因可归纳为：（1）注入压力消散引发流体重新分布；（2）浮力驱动下的持续上移导致气相分散；（3）砂岩孔隙毛细力对CO₂的捕集作用形成残余气饱和度.在盖层封堵作用下，CO₂最终在储层顶部形成稳定分布的薄层气带.当以S_g>0.12为羽流边界时，基质模型中CO₂羽流的最大水平扩展距离为10.8 m.

图10b展示了含层理构造的储层中CO₂饱和度（S_g）的时空分布特征.由于注入口靠近高渗透性的层理，在注入压力的驱动下，CO₂优先进入高渗透层理，形成快速运移通道.层理模型中CO₂羽流的最大水平扩展距离增至20.8 m，较基质模型增加43.45%.连续注入25 d后，CO₂在浮力驱动下垂向运移至上部层理，并通过层理高渗带横向扩展，显著增强了储层的水平连通性，形成以高渗层理为核心的椭球状分布模式.

图10c揭示了裂缝‒层理耦合系统中CO₂运移的动力学特征.由于裂缝的渗透率高于基质和层理，CO₂在注入压力梯度驱动下优先沿裂缝网络形成优势渗流通道.裂缝‒层理模型中CO₂羽流的最大水平扩展距离进一步增至25.2 m，分别是基质模型和层理模型的1.74倍和1.21倍.裂缝系统的特殊几何特征导致：高角度裂缝穿透层理界面，显著降低垂向封隔效率，使CO₂突破垂向封隔，快速向上运移；低角度裂缝构成横向输导网络，促进CO₂横向扩散.连续注入25 d后，裂缝成为CO₂运移的主干通道，大量CO₂通过裂缝‒层理网络运移至远端，但裂缝连通性差导致CO₂在裂缝节点处聚集，形成“骨架‒枝叶”复合羽流结构.

图11展示了关井25 d后，三种地质模型（基质模型、层理模型、基质‒层理‒裂缝模型）的CO₂饱和度随剖面高度变化.由于剖面是以注入口为中心的垂向对称剖面，为了方便观察，将剖面高度 （-2 042 m<h<-2 000 m）简化为（-21 m<h<21 m）.25 d扩散期后，在CO₂注入口上方（0 m<h<11 m），CO₂饱和度逐渐下降，但比注入口下方（h<0 m）的CO₂饱和度高，这是由于注入的超临界CO₂密度显著低于储层盐水密度，从而产生浮力，驱动CO₂自发向上运移.随着剖面高度的增加（h>11 m，即储盖层界面），CO₂饱和度骤降，这表明盖层有效封堵了CO₂的扩散.在基质模型中，储层均质中的毛细管力有效限制了CO₂垂向扩散.相比之下，含层理的模型中层理的高渗透性成为CO₂的主要运移通道，CO₂快速扩散，注入口的CO₂垂向饱和度较基质明显降低.而层理‒裂缝模型峰值的CO₂饱和度最低，但在注入口上方（h>0 m）随着剖面高度的增加，CO₂饱和度比含层理模型的饱和度高.这表明裂缝成为主要运移通道，构成的裂缝网络促进了CO₂的垂向扩散，从而增加了CO₂的垂向饱和度.

4.2　层理和裂缝对表面压力的影响

如图12a所示为岩石基质模型中注入CO₂后的水压力分布特征.在注入初期，注入井周边迅速形成高压区，压力值快速抬升，压力等值线呈现同心圆形态，表明在均质介质中以径向流动为主.当持续注入至 25 d时，高压区扩展至井周20 m 范围.受盖层与低渗透层的阻隔作用，垂向压力传递受到抑制，而左右两侧的静水压力边界则形成压力 “泄流区”，导致侧向压力发生耗散，压力等值线随之转变为椭圆扁平状.在停止注入后的25 d，系统压力出现骤降；同时，由于重力分异作用，CO₂羽流因密度较低在盖层下方聚集，咸水则向中下部区域运移，但残余CO₂气泡产生的毛细阻力对压力恢复过程形成阻碍.

如图12b所示为含层理多孔基质中CO₂注入过程的水压力分布特征.注入初期，压力沿高渗透层理快速延伸，层理作为压力传导的优势通道加速了压力扩散进程.受盖层和低渗透层的阻挡，垂向压力传递受到强烈抑制，形成有效的压力屏障，导致垂向压力梯度显著低于水平方向.持续注入25 d后，受重力分异作用影响，CO₂逐渐在盖层下方聚积，导致储层顶部形成局部高压区，压力等值线呈现明显的上凸特征.水平方向上高渗层理促进了压力耗散，使得压力场呈带状分布.停注25 d后，层理方向的CO₂持续扩散，导致近井地带压力快速降低，但盖层下方的CO₂积聚区压力衰减较慢.受浮力驱动，残余CO₂沿高渗透层理向上运移，并在盖层下方重新聚集，形成长期存在的超压带，显著延缓了储层上部的压力恢复.以上分析表明，高渗透层理显著改变压力传播路径.

如图12c所示为含层理‒裂缝的多孔基质中CO₂注入过程的水压力分布特征.注入初期，裂缝作为CO₂运移的主导通道，使CO₂在网络结构内迅速形成超高压带.其中垂向裂缝承担压力向上传递的功能，横向裂缝则与层理共同构成横向压力传导路径.持续注入25 d后，裂缝与层理交织形成的贯通网络使压力分布呈现 “主干‒分支”特征.同时，受重力分异作用驱动，CO₂持续向上运移，导致盖层底部逐渐形成高压区.停注25 d后，压力沿裂缝‒层理网络持续耗散，与仅含层理的情形相比，盖层底部的CO₂积聚量显著减少，这是由裂缝进一步增强了压力分散效应所致.

图13展示了关井25 d后，三种地质模型（基质模型、层理模型、基质‒层理‒裂缝模型）的CO₂压力随剖面高度变化.图中显示储层（-11 m< h<11 m）中的CO₂压力最高，而盖层（11 m<h<21 m）中的CO₂压力明显较在储层中的压力低，这表明盖层有效阻止了CO₂的扩散.图中显示在注入点为压力峰值，向上下地层压力快速衰减，但注入口上方压力低于注入口下方压力，这表明CO₂沿着孔隙、层理、裂缝向上运移，扩散范围变广，压力被分散.而层理模型和层理‒裂缝模型的CO₂压力要高于均质模型，表面层理和裂缝成为CO₂的主要运移通道，CO₂分布范围更广，整体压力更高.而裂缝渗透性最强，CO₂通过裂缝快速扩散至整个储层，甚至局部贯穿储层上下边界，导致压力最高.

注入压力直接影响CO₂相态在孔隙中的流动能力，如图14所示不同地质条件下注入口压力随时间的变化.在注入初期，注入口压力达峰值3.6×10⁷ Pa.随着注入时间的增加，压力逐渐扩散，但基质‒层理模型和基质‒层理‒裂缝模型均高于基质模型，表明层理与裂缝的高渗透性加速了压力耗散，层理‒裂缝下注入口压力较基质模型降低了约40%，表明层理和裂缝显著缓解了局部压力积聚.停止注入后，注入口压力骤降，并低于初始压力，但随着CO₂的扩散，表面压力趋于稳定，注入口压力缓慢回升，最后到达初始压力.

4.3　层理和裂缝对位移的影响

如图15所示为不同条件下的剖面位移. 在基质模型中，CO₂注入初期，基质作为相对均一的介质，CO₂在其中缓慢扩散，对周围岩石产生较小的压力，从而引起剖面位移的缓慢增加.由于基质中没有明显的薄弱结构，二氧化碳的运移相对较为均匀，所以位移变化较为平缓.在注入中期，位移继续上升，但斜率逐渐变小，呈现出趋近于稳定增长的趋势.表明随着CO₂的不断注入，基质中的孔隙逐渐被CO₂填充，CO₂的扩散速度减慢，对岩石的压力增加也逐渐变缓，导致剖面位移的增长速度降低.停止注入CO₂后，已有的CO₂在孔隙中缓慢扩散，压力变化较小，位移逐渐趋近于一个稳定值，在基质‒层理模型中，位移起始值略高于基质模型，且上升速度稍快.这是由于层理作为岩石中的薄弱面，在注入期为CO₂的扩散提供了相对通畅的通道，CO₂更容易沿着层理面扩散，相比之下，单纯的基质模型没有这种层理通道的优势，位移变化速率稍慢，在相同时间内，CO₂在基质‒层理模型中能够更快地传播并对周围岩石产生压力，导致剖面位移增加较快.随着注入的稳定，位移持续上升，增长速度比基质模型快，但同样斜率逐渐变小.层理面虽然有助于CO₂扩散，但随着注入的进行，层理面附近的压力逐渐平衡，CO₂的扩散也受到一定限制，位移增长速度放缓.停止注入CO₂后，位移增长明显变缓，最终趋于稳定.层理面的存在使得CO₂在扩散阶段仍然能够在一定范围内继续运移，但随着时间推移，压力逐渐平衡，位移不再明显增加.

在基质‒层理‒裂缝模型中， 位移起始值与基质‒层理模型相近，但上升速度在初期非常快.这是因为裂缝为CO₂的快速运移提供了通道，CO₂能够迅速通过裂缝在模型中扩散，从而在短时间内对周围岩石产生较大的压力，引起剖面位移的快速上升.随着CO₂的持续注入，位移上升速度开始减慢，但仍然高于前两个模型.这是因为裂缝虽然在初期加速了CO₂的运移，但随着CO₂的扩散，裂缝可能会部分被岩石颗粒等堵塞，或者裂缝周围的压力逐渐平衡，使得二氧化碳的运移速度降低，位移增长变缓.停止注入CO₂后，位移在扩散阶段仍然有一定的增长，但增长幅度逐渐减小并趋于稳定.裂缝在扩散阶段继续影响CO₂的运移，使得二氧化碳能够在更大范围内扩散，但随着压力平衡，位移的变化也逐渐停止.

4.4　CO₂封存稳定性

图16展示了储盖层系统中CO₂体积分数随时间的演化特征.随着CO₂注入过程的持续，系统内CO₂体积分数呈递增趋势；注入至第10天左右，CO₂突破盖层压力阈值进入盖层，盖层中CO₂体积分数开始缓慢上升，但始终维持在较低水平.这一现象表明，仅有少量CO₂能够穿透盖层，大部分CO₂被有效封存在储层内部.注入停止后，受重力分异作用影响，CO₂在盖层下方持续积聚并不断向盖层运移，但与注入阶段相比，盖层中CO₂体积分数的增长速率显著降低.具体到不同地质结构体系：

（1）在仅含基质的体系中，岩石结构致密均一，缺乏优势运移通道.从曲线特征看，其斜率较小且呈平缓上升趋势，表明CO₂向盖层的逃逸速率较低，大部分CO₂被有效限制在储层内部.

（2）当储层存在层理结构时，渗透率呈现显著各向异性——地质体中的层理均平行于层理方向分布，其渗透率相对较高，为CO₂运移提供了顺畅通道，但同时也抑制了CO₂的垂向运移，因此逃逸至盖层的CO₂量较少.

（3）在含裂缝的结构中，裂缝的发育显著改变了储层渗透性能及CO₂运移路径，成为CO₂运移的主导通道，导致CO₂在裂缝处形成高浓度条带.由于模型中裂缝连通性较差且主裂缝倾角较小，CO₂以横向运移为主并在储层中大量聚集，其扩散速率较基质模型及基质‒层理模型更快；但垂直方向上逃逸至盖层的CO₂体积分数未发生明显变化，因此未造成CO₂的大量逃逸.

为定量描述CO₂扩散至盖层的体积分数变化特征，本文引入“逃逸比”指标，定义为盖层中CO₂逃逸体积与总注入体积的比值.

\mathrm{逃逸}\mathrm{比}=\frac{V_{\mathrm{c}\mathrm{a}\mathrm{p}\mathrm{r}\mathrm{o}\mathrm{c}\mathrm{k}}}{V_{\mathrm{i}\mathrm{n}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}}}，(28)

式中：V_caprock是注入25 d后累积迁移并存于盖层中的CO₂总体积，m³；V_injected是注入25 d后累计注入储层的CO₂的总体积，m³.

V_{\mathrm{i}\mathrm{n}\mathrm{j}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}}=q\times t，(29)

式中：q是注入流量，在二维模型中值为注入速率r_in，m/s；t是累计注入时间.

V_{\mathrm{c}\mathrm{a}\mathrm{p}\mathrm{r}\mathrm{o}\mathrm{c}\mathrm{k}}={\iiint }_{{\Omega }_{\mathrm{c}\mathrm{a}\mathrm{p}\mathrm{r}\mathrm{o}\mathrm{c}\mathrm{k}}}^{}{S}_{\mathrm{n}}\left(x,t\right)\cdot \varphi \left(x\right)\mathrm{d}\Omega，(30)

式中：Ω_caprock为盖层的几何区域，S_n（x，t）为时刻t、空间位置x处的CO₂饱和度，无量纲.

图17展示了不同地质条件下CO₂的总体积分数及逃逸至盖层的体积分数分布.计算结果表明：在注入25 d后，基质模型、基质‒层理模型、基质‒层理‒裂缝模型的逃逸比依次为0.279 7、0.266 3、0.265 0，呈现逐步递减规律，这一特征与各模型的渗流特性具有显著关联性.具体而言，层理与裂缝通过改变渗流路径影响CO₂运移——本文所设置的层理为横向展布特征，且裂缝具有连通性差、倾角小的特点，二者共同作用使CO₂在储层内的运移更趋分散，从而降低了向盖层泄漏的风险.

5 结论

本文通过构建天然裂缝性储盖层系统模型，对基质、基质‒层理、基质‒层理‒裂缝三种地质模型中CO₂的运移特征及储层压力的演化过程进行了数值模拟，并引入CO₂向盖层扩散的逃逸比指标，用于定量评估CO₂地质封存的稳定性.主要结论如下：

（1）在单一基质储层模型中，CO₂注入后羽流呈均匀扩散特征；当储层含层理与复杂裂缝时，介质非均质性显著增强，形成高渗透性运移通道网，导致CO₂运移方向与速率发生改变——其优先沿层理与裂缝展布方向运移，并在裂缝端口处发生积聚.停止注入后，CO₂虽沿原路径持续扩散，但受重力分异作用影响，逐渐在盖层下方形成聚集.

（2）注入初期，低粘度CO₂驱替原位咸水，导致储层压力快速抬升；随着注入过程持续进行，CO₂获得更充分的运移空间，储层压力逐渐下降.当储层含层理与裂缝时，CO₂优先沿此类高渗透通道扩散，使注入点附近的CO₂聚集量减少，储层压力相应降低，表明在相同储层条件下层理与裂缝对压力的缓解作用更显著.停止注入后，储层压力骤降，CO₂羽流在盖层下方逐渐积聚，咸水向中下部区域回填，系统表面水压力则缓慢回升.

（3）CO₂注入后，高渗透性的裂缝作为CO₂运移的主要通道，CO₂扩散范围增大，孔隙压力升高，显著降低有效应力，增大位移.而水平层理提供定向变形弱面，促进垂向应变，导致位移增大.

（4）CO₂向盖层的逃逸比计算结果显示：单一基质模型的25 d逃逸比为0.279 7，而含层理与裂缝模型的逃逸比分别降至0.266 3与0.265 0.模型中的裂缝网络连通性较低，这在一定程度上促进了CO₂的横向扩散，同时抑制了其垂向窜流，同时，天然裂缝中较大裂缝以低倾角裂缝为主，更易促进CO₂在储层内部的横向展布，从而降低了逃逸比.因此，在靶区优选和工程设计中，应通过地球物理探测与地质建模精确刻画裂缝倾角的分布规律.

（5）在进行CO₂注入时，应基于地质结构实行差异化注入，对于裂缝发育程度高的区域，应严格控制注入速率与井底压力，利用裂缝的高传导性实现CO₂的大范围扩散，同时防止形成贯穿性流道.对于以基质为主的均质区域，可适当提高注入速率，因其压力积累和位移响应相对缓和.其次，充分利用层理导向作用，将注入井应优先部署在层理发育的层段.利用层理的横向高渗透性，引导CO₂远离井筒横向扩展，增大波及体积，从而降低井周局部压力和垂向逃逸风险.

参考文献

[1] Bai, Y., Bai, B., Xu, W. L., et al., 2022. Pore Characteristics of Shale and Occurrence Mode of Shale Oil in Member 7 of Yanchang Formation in Southern Ordos Basin. Acta Petrolei Sinica, 43(10): 1395-1408 (in Chinese with English abstract).

[2] Bigi, S., Battaglia, M., Alemanni, A., et al., 2013. CO₂ Flow through a Fractured Rock Volume: Insights from Field Data, 3D Fractures Representation and Fluid Flow Modeling. International Journal of Greenhouse Gas Control, 18: 183-199. https://doi.org/10.1016/j.ijggc.2013.07.011

[3] Cai, B. F., Zhang, L., Lei, Y., et al., 2023. A Deeper Understanding of the CO₂ Emission Pathway under China’s Carbon Emission Peak and Carbon Neutrality Goals. Engineering, 9(11): 27-29.

[4] Du, S.H., Shi, Y.M., Guan, P., 2019. Fluid Filling Rule in Intra⁃Granular Pores of Feldspar and Fractal Characteristics: A Case Study on Yanchang Formation Tight Sandstone Reservoir in Ordos Basin. Earth Science, 44(12): 4252-4263 (in Chinese with English abstract).

[5] Fu, J. H., Li, S. X., Niu, X. B., et al., 2020. Geological Characteristics and Exploration of Shale Oil in Chang 7 Member of Triassic Yanchang Formation, Ordos Basin, NW China. Petroleum Exploration and Development, 47(5): 870-883 (in Chinese with English abstract).

[6] Han, L., Shi, X., Ni, H. J., et al., 2025. Fracture Initiation and Propagation Behaviours of Supercritical CO₂ Enhanced Fracturing in Layered Shale of Horizontal Wells. Geoenergy Science and Engineering, 252: 213938. https://doi.org/10.1016/j.geoen.2025.213938

[7] Hyman, J. D., Jiménez⁃Martínez, J., Gable, C. W., et al., 2020. Characterizing the Impact of Fractured Caprock Heterogeneity on Supercritical CO₂ Injection. Transport in Porous Media, 131(3): 935-955. https://doi.org/10.1007/s11242⁃019⁃01372⁃1

[8] Jiang, L. L., Tian, L., Chen, Z. X., et al., 2025. Research on the Micro⁃Pore Structure and Multiscale Fractal Characteristics of Shale under Supercritical CO₂ Action: A Case Study of the Chang 7₃ Submember in the Ordos Basin, China. Journal of Natural Gas Geoscience, 10(3): 159-178. https://doi.org/10.1016/j.jnggs.2025.05.003

[9] Kong, D. H., Wu, F. Q., Saroglou, C., 2020. Automatic Identification and Characterization of Discontinuities in Rock Masses from 3D Point Clouds. Engineering Geology, 265: 105442. https://doi.org/10.1016/j.enggeo.2019.105442

[10] Kong, X. Y., Chen, F. L., Chen, G. Q., 1999. Mathematical Models and Feature Parameters of Non Newtonian Liquid Flows in Porous Media. Journal of University of Science and Technology of China, 29(2): 141-147 (in Chinese with English abstract).

[11] Lee, I. H., Ni, C. F., 2015. Fracture⁃Based Modeling of Complex Flow and CO₂ Migration in Three⁃Dimensional Fractured Rocks. Computers & Geosciences, 81: 64-77. https://doi.org/10.1016/j.cageo.2015.04.012

[12] Li, D. C., Saraji, S., Jiao, Z. S., et al., 2023. An Experimental Study of CO₂ Injection Strategies for Enhanced Oil Recovery and Geological Sequestration in a Fractured Tight Sandstone Reservoir. Geoenergy Science and Engineering, 230: 212166. https://doi.org/10.1016/j.geoen.2023.212166

[13] Li, S., Chen, J. B., Wang, H. Q., et al., 2020. Experimental Study on the Scale Effect of Strength and Deformation of Chang 7 Shale in Ordos Basin. Journal of China Coal Society, 45(12): 4121-4131 (in Chinese with English abstract).

[14] Liao, Z. W., Yang, J. M., Zhong, X. Y., et al., 2024. Review on Research Progress of Carbon Dioxide Geological Sequestration Technology. Chinese Journal of Underground Space and Engineering, 20(S1): 497-507 (in Chinese with English abstract).

[15] Liu, Y. Q., Chen, J. P., Tan, C., et al., 2022. Intelligent Scanning for Optimal Rock Discontinuity Sets Considering Multiple Parameters Based on Manifold Learning Combined with UAV Photogrammetry. Engineering Geology, 309: 106851. https://doi.org/10.1016/j.enggeo.2022.106851

[16] Nordbotten, J. M., Celia, M. A., Bachu, S., 2005. Injection and Storage of CO₂ in Deep Saline Aquifers: Analytical Solution for CO₂ Plume Evolution during Injection. Transport in Porous Media, 58(3): 339-360. https://doi.org/10.1007/s11242⁃004⁃0670⁃9

[17] Oh, J., Kim, K. Y., Han, W. S., et al., 2013. Experimental and Numerical Study on Supercritical CO₂/Brine Transport in a Fractured Rock: Implications of Mass Transfer, Capillary Pressure and Storage Capacity. Advances in Water Resources, 62: 442-453. https://doi.org/10.1016/j.advwatres.2013.03.007

[18] Pavan, T. N. V., Devarapu, S. R., Govindarajan, S. K., 2024. Numerical Investigations on the Performance of Sc⁃CO₂ Sequestration in Heterogeneous Deep Saline Aquifers under Non⁃Isothermal Conditions. Gas Science and Engineering, 130: 205437. https://doi.org/10.1016/j.jgsce.2024.205437

[19] Qin, X. J., Wang, H., Xia, Y. X., et al., 2025. Pore⁃Scale Investigation of Water⁃CO₂⁃Oil Flow in Shale Fractures for Enhanced Displacement Efficiency and CO₂ Sequestration. Engineering Geology, 348: 107969. https://doi.org/10.1016/j.enggeo.2025.107969

[20] Ren, F., Ma, G. W., Wang, Y., et al., 2017. Two⁃Phase Flow Pipe Network Method for Simulation of CO₂ Sequestration in Fractured Saline Aquifers. International Journal of Rock Mechanics and Mining Sciences, 98: 39-53. https://doi.org/10.1016/j.ijrmms.2017.07.010

[21] Rui, Z. H., Zeng, L. B., Dindoruk, B., 2025. Challenges in the Large⁃Scale Deployment of CCUS. Engineering, 44: 17-20. https://doi.org/10.1016/j.eng.2024.11.031

[22] Sheng, D. N., Wang, H. M., Sheng, J. C., et al., 2025. Effect of Random Fracture Network Orientations on Sealing Performance of Caprock in CO₂ Geological Sequestration. Earth Science, 50(1): 349-360 (in Chinese with English abstract).

[23] Su, X. C., Gong, L., Fu, X. F., et al., 2023. Fracture Distribution Characteristics and Effectiveness Evaluation of Tight Sandstone Reservoir of Chang 7 Member in Sanbian Area, Ordos Basin. Earth Science, 48(7): 2601-2613 (in Chinese with English abstract).

[24] Sun, H., Jia, C., Xiong, F., et al., 2024. Numerical Modelling of CO₂ Leakage through Fractured Caprock Using an Extended Numerical Manifold Method. Engineering Analysis with Boundary Elements, 162: 327-336. https://doi.org/10.1016/j.enganabound.2024.02.013

[25] Sun,Z.D., Song, X.M., Feng, G., et al.,2020.Influence of Supercritical, Liquid, and Gaseous CO₂ on Fracture Behavior in Sandstone.Energy Science and Engineering, 8(11): 3788-3804. https://doi.org/10.1002/ese3.736

[26] Wang, H. D., Chen, Y., Ma, G. W., 2020. Effects of Capillary Pressures on Two⁃Phase Flow of Immiscible Carbon Dioxide Enhanced Oil Recovery in Fractured Media. Energy, 190: 116445. https://doi.org/10.1016/j.energy.2019.116445

[27] Wang, W., Liang, Z. Z., Zuo, J., et al., 2025. The Pore Structure Changes and CO₂ Migration Dynamic Characteristics in Tight Sandstone during Supercritical CO₂ Geosequestration: A Case Study in the Chang 7 Layer, Ordos Basin, China. Fuel, 379: 133019. https://doi.org/10.1016/j.fuel.2024.133019

[28] Wang, Y., Li, T., Chen, Y., et al., 2019. Numerical Analysis of Heat Mining and Geological Carbon Sequestration in Supercritical CO₂ Circulating Enhanced Geothermal Systems Inlayed with Complex Discrete Fracture Networks. Energy, 173: 92-108. https://doi.org/10.1016/j.energy.2019.02.055

[29] Yang, S. Q., Hong, W. X., Sun, B. W., et al., 2023. Experimental Study on Triaxial Mechanics and Failure Characteristics of Shale in Different Brine Environments. Chinese Journal of Geotechnical Engineering, 45(11): 2217-2226 (in Chinese with English abstract).

[30] Yang, S. Q., Xu, S. B., Liu, Z., 2022. Experimental Study on Mechanical and Permeability Behaviors of Sandstone under Deep Saline Environments. Chinese Journal of Rock Mechanics and Engineering, 41(2): 292-304 (in Chinese with English abstract).

[31] Yuan, Z., Ren, P. G., Liu, J. H., et al., 2025. Research on Whole Strata Geological Modeling Technology for CO₂ Geological Storage in Salt Water Layer. Earth Science, 50(5): 1987-1998 (in Chinese with English abstract).

[32] Zheng, F. N., Jha, B., Jafarpour, B., 2025. Mitigating Caprock Failure and Leakage Risks through Controlled CO₂ Injection and Coupled Flow⁃Geomechanics⁃Fracturing Simulation. International Journal of Greenhouse Gas Control, 144: 104387. https://doi.org/10.1016/j.ijggc.2025.104387

基金资助

国家自然科学基金委员会面上项目(42077258)