地质封存二氧化碳沿断层泄漏数值模拟研究

张志雄，谢 健，戚继红，胡立堂，张可霓，吴礼舟

(1.成都理工大学环境与土木工程学院，四川 成都 610059；2.北京师范大学水科学研究院，北京 100875； 3.劳伦斯伯克利国家实验室，美国加州 伯克利 94720)

利用深部咸水层进行CO2地质封存(GCS)被认为是实现大规模减排、缓解温室气体效应最具前景的技术之一[1]。随着近年GCS工程在全球的不断推进和实施，其安全性和泄漏风险愈来愈引起人们重视。由于含水层不似废弃油气田，未有储气“经历”，因此其盖层(隔水层)完整性和地层垂向组合特征是否满足GCS密封性要求未得到充分验证。断层作为影响储层地质条件复杂性的重要不确定性因素，对GCS工程的安全性、封存效果具有至关重要的影响[2]。断层及其破碎带往往是结构复杂的非均质系统[3]，GCS场地如果存在穿透主要储层和盖层的断层，CO2极有可能沿着断层破碎带迁移至浅部地层甚至泄漏至地表。这不仅降低储层的封存能力，而且可能导致上覆含水层污染、激活断层等环境地质效应[4]。

研究地质封存CO2及其储层咸水在注入和封存阶段沿断层迁移泄漏的主要特征和受控因素是极为必要的。然而，目前国内外关于这一课题的研究还十分有限。Pruess[5～7]研究了CO2泄漏过程中的物理化学作用。郑菲等[8]对苏北盆地CO2地质封存泄漏风险的全局敏感性分析，发现与毛细压力有关的断层参数对断层中CO2总量的影响最大。董华松等[9]重点关注了GCS泄漏监测技术的发展。崔振东等[4,10～11]研究了CO2泄漏的风险及后果，指出泄漏可导致浅部淡水资源被污染、生态系统和人类健康受威胁，还有可能激活断层。总体来看，国内学者多数停留在介绍、综述国外研究成果和理论，缺少这一课题的专门定量研究[12]。

本研究运用多相流达西定律，分析储层流体渗流的影响因素。在此基础上利用储层多相流模拟软件TOUGH2/ECO2N建立一系列场地泄漏模型，旨在揭示和阐明GCS流体沿断层泄漏的主要机理或规律，评估断层对GCS泄漏的可能影响，为场地评估和封存效果研究提供理论依据。

1 理论分析

CO2注入到含水层中会引起H2O-NaCl-CO2混合系统多种复杂的物理化学耦合作用，其中最主要的是多组分多相流体的渗流过程，这可用多相流达西定律进行描述[13]。向含水层中注入CO2会引起孔隙压力的抬升，这必然会改变储层流体所受的压力梯度。而系统毛细压力和流体黏度的变化也会导致流体流动性的变化。如果CO2所在含水层被断层切割且断层导水，则储层流体在压力梯度和密度差引起的浮力作用下必然沿断层发生迁移，进而突破盖层进入上覆含水层。CO2与咸水在含水层乃至断层中的迁移扩散的过程[1,14]：

(1)

(2)

(3)

式中：t——注入时间；

Vn——流动系统的体积；

M——单位体积的质量或能量；

K——组分；

F——质量通量；

n——表面单元dΓn的法向矢量；

q——质量或能量源汇项；

φ——孔隙度；

S——饱和度；

β——相态(气相或液相)；

ρ——密度；

X——质量分数；

k——岩体固有渗透率；

kr——相对渗透率；

μ——黏度；

g——重力加速度。

由式(3)可知，在相同的储层环境下，增大CO2注入速率，系统内气相相对渗透率及压力梯度Pβ均增大，这会导致Fβ也增大，即泄漏速率及泄漏量增大；断层发育位置距离注入井越近，相同注入压力条件下断层中的压力梯度Pβ越大，因此泄漏速率及泄漏量也增大；断层倾角越小，泄漏途径增长，因此系统内压力梯度Pβ减小，泄漏速率及泄漏量也因此减小；断层断距的变动可改变断层两盘的连通性，若因断距使两盘的连通性增大，也就相当于k增大，因此泄漏速率及泄漏量也增大；对于断层厚度的改变，断层单位宽度上的Fβ不变，但厚度增大后，泄漏量会增多；增大断层渗透性，显然泄漏速率及泄漏量都增大。

2 模型说明及模拟方案

2.1 地质模型及网格剖分

选取一套砂岩泥岩互层的地层为模型研究场地，包括3层泥岩和3层砂岩(图1)，以砂岩3为目标储层，泥岩2为主盖层。模型水平方向以注入井为中心，径向延伸1 000 m，地层模拟段埋深1 000～1 220 m。基础模型断层发育位置距注入井300 m，断层垂直。另外设置了一组研究倾斜断层的模型，模型水平径向范围为900 m，地层埋深1 000～1 200 m。所有模型均采用矩形网格剖分，网格大小0.5～50 m不等，于注入井和断层附近加密，加密区网格水平间距为0.3～1.0 m。为了更好地刻画断层，将断层带分成渗透性较高的破碎强烈带和渗透性较小的断层破碎影响带，而且在研究断层渗透性的模型中还将断层分成上下1区和2区(图2)。

图1 地质概念模型及网格剖分图Fig.1 Geological model and grid design

图2 断层带分区刻画示意图Fig.2 Schematic diagram of the zonation of faults

2.2 模型方案及参数设置

为揭示可能影响GCS泄漏的断层特征和注入方式，设计了6个方案总共15个模拟进行研究。方案I、II、III、IV、V、VI是分别用于探讨注入速率、断层位置、断层倾角、断层断距、断层带厚度及断层渗透性对泄漏的影响，各模拟方案及参数设置，见表1。为了方便对照，表中基础模型(编号1)有重复。

模型中储层及盖层的参数是参考神华CCS项目的主要储层刘家沟组及上部泥岩层设置的，其它参数则是在储层及盖层的基础上设计的；断层的倾角设置成90°(研究倾角的除外)，这主要是为了方便刻画模型；切穿深部地层的断裂倾角一般较大，因此在方案III中设置了大于45°的三个倾角；对于注入速率3.17 kg/s是参考神华项目的10×104t/a设置。概念模型中各地质单元的水文地质参数具体设置见表2。

2.3 模型边界条件与初始条件

模型右边侧向边界设置为常压边界，以保证含水层水平方向的开放(无限含水层)，设置成常压边界是将边界上的网格体积设置成超大体积(1051m3)实现的。每个模型的左侧边界分成三段，上下两段都设置成无渗流边界，渗透率设成很低(10-20m2)；中段设置成常压边界，保证CO2能正常注进储层，只需把渗透率设置成比储层的渗透率略大。

模型中各地层温度假设为42.5 ℃，模拟过程不考虑非等温效应。模型顶部压强参考神华CCS场地取12.0 MPa，利用重力平衡计算得出模型的初始压力分布。盐度也参考神华的刘家沟组咸水层取0.03%。各地层中的初始CO2质量分数为0。

表1 模型方案及断层参数设置Table 1 Model schemes and associated fault parameters

注：表中V为CO2注入速率；ζ为断层距注入井的距离；α为断层倾角；d为断层断距；λ为断层破碎带厚度；K1x,K1y为断层破碎影响带水平方向渗透率；K1z为断层破碎影响带垂直方向渗透率；K2x,K2y为断层破碎强烈带水平方向渗透率；K2z为断层破碎强烈带垂直方向渗透率。

表2 概念模型的各地质单元物性参数设置一览表Table 2 Physical properties of the geological units in the models

注：表中d1为主模型各地质单元厚度；d2为辅助模型各地质单元厚度；n为孔隙度；ρ为密度；Kx,Ky为水平方向渗透率；Kz为垂直方向渗透率；P0为进气毛细压力；m为孔隙分布指数；Slr为残余液相饱和度；Sgr为残余气相饱和度。

3 模型结果及分析

3.1 泄漏时间

文中以有CO2通过盖层底部断层位置的时刻为开始泄漏时间，是表征CO2在储层中迁移速率的一个变量。文中储层及流体的性质都一致，所以CO2从注入到泄漏的时间长短，主要取决于压力梯度Pβ的大小和断层位置的远近。表3为6组模拟方案的泄漏时间列表，可见注入速率、断层位置及断层倾角对泄漏时间的影响最大，其次是断层带的渗透性，而断层的断距、厚度对其影响不大。注入速率由1.585 kg/s增大到6.34 kg/s，泄漏时间提前了3 853 d；断层位置由500 m减至100 m，泄漏时间由9 821 d提前至277 d；断层倾角由45°增至90°，泄漏时间由132 d提前至91 d，突破模型顶部时间提前了360 d；断层渗透性放大一倍后，虽然开始泄漏时差别不大，但突破模型顶部时间提前了3 782 d；但考虑到断层上下部分渗透性不同时，发现其泄漏时间是由两部分共同决定的，但并不是对等关系。

3.2 泄漏速率

表3 不同模拟方案CO2泄漏时间对比Table 3 Comparison of the time at which CO2 leaks for various modeling schemes

图3 地下水及CO2的泄漏速率与时间的关系Fig.3 Leakage rate of CO2 and groundwater

泄漏速率不同于泄漏时间，它表征的是单位时间泄漏量的多少。图3为各组模型的地下水及CO2泄漏速率曲线。需要注意的是讨论泄漏速率皆为Z向的。注入速率对CO2及地下水的泄漏速率的影响显著(图3a)，注入速率由1.585 kg/s抬高至6.34 kg/s，CO2的泄漏速率加快了0.07～0.075 kg/s，地下水的泄漏速率加快了0.18～0.20 kg/s。断层位置的远近对CO2泄漏速率影响较大，距注入井500 m位置的CO2泄漏速率相对于100 m的加快了0.09～0.16 kg/s；但其对地下水的影响更大，这主要是因为地下水密度更大，而且CO2驱替地下水是一个积累的过程，断层位置越远地下水泄漏量越大。图3c中可见，断层倾角越大泄漏速率偏小，90°倾角相对于45°倾角，CO2泄漏速率小了0.18～0.34 kg/s。图3d中的规律不够显著，应该是断距的大小会影响Z向和R向的泄漏量大小的比例。断层带厚度的影响很显著，厚度越窄速率越快。断层厚度由20 m减为10 m，地下水速率增加了0.065～0.112 kg/s，CO2泄漏速率增加了0.048～0.054 kg/s。图3f四个曲线差别不大，断层渗透性提升一倍，CO2泄漏速率提高了约0.005 kg/s。

3.3 CO2羽体分布及泄漏量

CO2发生泄漏后，泄漏量的大小可以用于评判CO2封存效果及封存安全性。有学者认为，CO2的泄漏量不应该超过注入总量的1‰[9]。图4为六个方案10a或20a后的CO2羽体分布，水平方向为R向，垂直方向为Z向(深度)，注入井位置在图的左侧。泄漏量的大小可大致通过CO2羽分布范围大小来作比较，但为了更加精确的分析，图中右侧有注明储层残余CO2的比例。注入总量是注入速率与注入时间的乘积，基础模型的CO2注入总量为1.99×106t。图中可见CO2羽体在储层中因密度小呈现“漏斗”状，而在上部砂岩层因迁移先后呈现“倒漏斗”状。

综合10a或20a的CO2分布和储层残余CO2发现，注入速率、断层位置对CO2泄漏量影响最大。注入速率(图4a)为1.585kg/s时，储层残余CO2占注入总量的81.97%；速率增大一倍后，残余CO2为60.40%；速率增大3倍后，残余CO2为49.54%，泄漏量增多32.43%。断层位置(图4b)为100 m， CO2泄漏量高达注入量63.39%；而在模型5(断层位置500 m)中还未发生泄漏。模型倾角由45°增大至90°，泄漏量相比增大了1.95%。图4d可见，断层右盘越靠下，泄漏量越大，当断距为40 m时，泄漏量高达注入量的43.46%。断层带厚度越大，泄漏量越大，厚度20 m相对于10 m的情况，泄漏量增大了5.26%(图4e)。断层渗透性越强，泄漏的越多，渗透率加倍后，泄漏量增加了2.11%(图4f)；当考虑到断层上下部分渗透性不同时，发现其泄漏量是由两部分共同决定的，但并不是对等关系。

图4 不同模拟方案的CO2羽分布图Fig.4 Comparison of the CO2 plumes for different model schemes

3.4 结果分析

综合上述6个方案的模拟结果发现，研究注入速率、断层位置的I、II两个方案结果极其显著。在相同的储存环境下，增大CO2注入速率，系统内气相相对渗透率及压力梯度Pβ均增大，这会导致泄漏速率及泄漏量增大；断层发育位置距离注入井越近，相同注入压力条件下断层中的压力梯度Pβ越大，因此泄漏速率及泄漏量也增大。

然而研究断层倾角、断距、厚度及渗透性的III、IV、V、VI 4个方案结果不是很显著。断层倾角越小，泄漏途径增长，因此系统内压力梯度Pβ减小，泄漏速率及泄漏量也因此减小；断层断距的变动可改变断层两盘的接触关系，极可能影响泄漏通道的连通性，若因断距使连通性增大，也就相当于k增大，因此泄漏速率及泄漏量也增大；对于断层厚度的改变，系统内压力梯度Pβ差别不大，导致断层单位宽度上的Fβ变化很小，但厚度增大后，泄漏量会增多；增大断层渗透性，显然泄漏速率及泄漏量都增大。

4 结论

(1)6组模拟方案中，研究注入速率、断层位置的I、II 2个方案结果极其显著。注入速率由1.585 kg/s增大至6.34 kg/s，CO2泄漏将提前3 853 d，泄漏量增多32.43%；断层距注入井的位置由500 m减小至100 m，CO2泄漏将提前9 544 d，泄漏量增多63.39%。

(2)6组模拟方案中，研究断层倾角、断距、厚度及渗透性的III、IV、V、VI 4个方案结果较显著。断层倾角由45°增大至90°，CO2泄漏将提前41 d，泄漏量增多1.95%；断层断距由-20 m增大至40 m，泄漏时间差别不大，但泄漏量增多5.97%；断层厚度由10 m增大到20 m，泄漏时间差别不大，但泄漏量增多5.26%；断层渗透性增大一倍，泄漏时间差别不大，但泄漏量增多2.11%；断层上下部分渗透性不同时， CO2泄漏的影响是由两部分共同决定的，但并不是对等关系。

(3)模型结果分析表明，CO2注入方案、断层性质(发育位置、产状、几何形态、内部结构)、系统内岩层的组合形态对CO2泄漏均有影响。注入速率越大、断层位置越近、断层倾角越大、断层越厚、渗透性越高，CO2泄漏的越早、泄漏量也越大。但断距对CO2的泄漏主要受各岩层厚度及组合形态，与断距大小关系不大。在这些因素中，对CO2泄漏影响最大的是注入速率和断层位置。

[1] 谢健，张可霓，王永胜, 等.鄂尔多斯深部咸水层CO2地质封存效果评价[J]. 岩土力学, 2016, 37(1): 166-174.[XIE J, ZHANG K N, WANG Y S,etal. Performance assessment of CO2geological storage in deep saline aquifers in Ordos Basin, China[J]. Rock and Soil Mechanics, 2016, 37(1): 166-174.(in Chinese)]

[2] Chuanhe Lu, Yunwei Sun, Thomas A,etal. Uncertainty quantification of CO2leakage through a fault with multiphase and nonisothermal effects[J]. Greenhouse Gases, 2012, 2(6): 445-459.

[3] VIALLE S, DRUHAN JL, MAHER K. Multi-phase flow simulation of CO2leakage through a fractured caprock in response to mitigation strategies[J]. International Journal of Greenhouse Gas Control, 2016,44: 11-25.

[4] 崔振东, 刘大安, 曾荣树，等. CO2地质封存工程的潜在地质环境灾害风险及防范措施[J].地质论评, 2011, 57(5): 700-706.[CUI Z D, LIU D A, ZENG R S,etal. Potential geological and environmental risks and its prevention measures for CO2geological storage projects[J]. Geological Review, 2011, 57(5): 700-706.(in Chinese)]

[5] PRUESS K. Integrated modeling of CO2storage and leakage scenarios including transitions between super-and subcritical conditions, and phase change between liquid and gaseous CO2[J]. Greenhouse Gases, 2011, 1(3): 237-247.

[6] PRUESS K. On CO2fluid flow and heat transfer behavior in the subsurface, following leakage from a geologic storage reservoir[J]. International Journal of Environmental Geology, 2008, 54(8): 1677-1686.

[7] PRUESS K. Leakage of CO2from geologic storage: Role of secondary accumulation at shallow depth[J]. Greenhouse Gas Control, 2008, 2(1): 37-46.

[8] 郑菲, 施小清, 吴吉春, 等. 苏北盆地盐城组咸水层CO2地质封存泄漏风险的全局敏感性分析[J]. 高校地质学报, 2012, 18(2): 232-238.[ZHENG F, SHI X Q, WU J C,etal. Global sensitivity analysis of leakage risk for CO2geological sequestration in the saline aquifer of Yancheng formation in Subei Basin[J]. Geological Journal of China Universities, 2012, 18(2): 232-238.(in Chinese)]

[9] 董华松, 黄文辉. CO2捕捉与地质封存及泄漏监测技术现状与进展[J]. 资源与产业, 2010, 12(2): 123-128.[DONG H S, HUANG W H. Researach of CO2capture, geological storage and leakage technologies[J]. Resources & Industries, 2010, 12(2): 123-128.(in Chinese)]

[10] 李毅, 张可霓, 王笑雨. CO2地质封存泄漏对浅层地下水影响的分析评价[J]. 工程勘察, 2014, 42(11): 44-50. [LI Y, ZHANG K N, WANG X Y. Assessment of the impact on shallow groundwater system by leakage of CO2geological storage[J]. Geotechnical Investigation & Surveying, 2014, 42(11): 44-50.(in Chinese)]

[11] 田地, 马欣, 查良松, 等. 地质封存CO2泄漏对近地表陆地生态系统的影响综述[J]. 生态与农村环境学报, 2013, 29(2): 137-145. [TIAN D, MA X, ZHA L S,etal. Review of impact of CO2leakage from geologic storage on near-surface terrestrial ecological system.[J]. Journal of Ecology and Rural Environment, 2013, 29(2): 137-145. (in Chinese)]

[12] 谢健，魏宁，吴礼舟，等. CO2地质封存泄漏研究进展[J].岩土力学, 2017,38(增刊1):181-188. [XIE J, WEI N, WU L Z,etal. Progress in leakage study of geologic CO2 storage[J]. Rock and Soil Mechanics, 2017,38(Sup1):181-188.(in Chinese)]

[13] ZHANG Keni, XU Yaqin, LING Lulu,etal. Numerical investigation for enhancing CO2injectivity in saline aquifers[J]. Energy Procedia, 2013,37: 3347-3354.

[14] Pruess K. ECO2N: A TOUGH2 Fluid Property Module for Mixtures of Water, NaCl, and CO2[R]. Berkeley: Lawrence Berkeley National Laboratory, 2005.