中高渗倾斜地层与水平地层中CO2地质封存的差异性对比

王福刚，郭 兵，杨永智，汪 芳，田海龙

(1. 吉林大学 地下水资源与环境教育部重点实验室，吉林 长春 130021；2. 中国石油勘探开发研究院，北京 100083)

0 引 言

工业革命以来，大量化石燃料燃烧产生的CO2排放到大气后加剧了温室效应，由此引发了一系列的环境问题[1-2]。因此，第六届国际气候变化委员会(IPCC)的气候报告强调，务必在21世纪中叶将全球温度增幅控制在1.5 ℃以内[3]。控制CO2排放被认为是限制气温增长的关键[4]，而CO2地质封存(CGS)被认为是一种直接有效减少CO2排放的手段，且具有较强的可操作性[5-7]。

CO2地质封存技术较早在欧洲一些国家得到应用，例如北海的斯莱普钠项目、英国石油公司在阿尔及利亚开展的萨拉赫项目等[8]。2010年，中国神华集团在鄂尔多斯盆地开展了CO2捕获和封存示范项目[9-11]，实际注入规模达到了10×104t·年-1。此外，在开展CO2地质封存项目的同时，众多国内外学者进行了基于实际场地的数值模拟，旨在为CO2大规模注入提供技术上的支撑[12-15]。

目前，国内外关于CO2地质封存的研究多基于水平地层展开。然而，由于构造活动和成岩作用，地层往往不是水平的，所以开展倾斜地层对CO2地质封存效率及安全性影响的研究具有实际意义。已开展的倾斜地层相关研究中，Goater等以北海项目为例，评价了有效渗透率、渗透率非均质性及倾斜地层(地层倾角为0.27°～3.00°)对CO2地质封存效率的影响[16]；Hesse等基于数值模拟构建了一个简单的模型来评估二维承压倾斜含水层中的CO2运移[17]；Pruess等利用数值模拟研究了倾斜地层中CO2羽状分布特征[18]；靖晶等以鄂尔多斯盆地石千峰组作为研究对象，研究了低孔低渗倾斜地层(孔隙度为12%～15%，渗透率为1.0～10.0 mD)中注入压力、地层倾角对CO2注入量的影响[19]。目前，对于中高渗倾斜地层的相关研究还不多见，系统性研究亟待加强。由于中高渗地层流动条件较好，是进行CO2地质封存的适宜地层，所以针对该类型倾斜地层进行CO2地质封存效率及安全性的研究是必要的。

新疆准噶尔盆地阜康地区能源丰富，周边分布有华能、鲁能、鸿基焦化、松迪焦化、中泰化学等大型发电厂和化工厂。这些企业在能源开采与利用过程中将会产生大量CO2，在该地区进行CO2地质封存工程能够有效缓解区域内CO2过度排放的状况。本文以新疆准噶尔盆地阜康凹陷CO2地质封存示范场地东沟组地层为目标地层，基于场地地层特征及其环境条件，采用数值模拟方法，研究倾斜地层对CO2地质封存效率及安全性的影响，为该地区进行CO2地质封存提供科学依据。

1 区域地质概况

新疆阿勒泰地区位于中国第二大沙漠——古尔班通古特沙漠的东部，西部为沙漠腹地，东北部为阿尔泰山脉，南部为博格达山，东部为北塔山[20-22]。由于该地区东沟组地层的孔隙度和渗透率较高，并且地层埋深为2 000～2 286 m，温度和压力等条件都比较适合CO2地质封存，所以选取东沟组地层作为目标地层。研究区位于准噶尔盆地阜康凹陷，东沟组地层在整个研究区内均有分布，上覆古近系地层，与下伏白垩系连木沁组与胜金口组地层呈角度不整合接触，地层起伏平缓，倾角最大为6°，倾向SW[23]。东沟组地层主要为河流三角洲相沉积，岩性主要为黄灰色、紫红色、灰黄色泥质粉砂岩，紫红色粉砂岩，灰色细砂岩，黄灰色、褐红色、紫红色泥岩，砂质泥岩，棕褐色泥岩。

此次研究的目标注入井为阜康凹陷中部的一口预探井[24](图1)，该井地面海拔和井深分别为463.83 m和5 140.00 m，位于阜康市东北方向30 km处。根据搜集资料，该井3个含水层射孔层位均为砂岩地层，自下往上分别为:第一层白垩系连木沁组和胜金口组，井段为2 392.00～2 407.00 m；第二层白垩系东沟组，井段为2 246.50～2 265.00 m；第三层白垩系东沟组，井段为2 036.50～2 066.20 m。其中，第三层地层厚度较大，埋深位置适中，孔渗条件较好，适宜CO2注入，选取其作为本次研究CO2注入的目标地层。

图件引自文献[24]图1 新疆准噶尔盆地阜康凹陷目标注入井地理位置Fig.1 Location of Target Injection Well in Fukang Sag of Junggar Basin, Xinjiang

2 数值模型构建

2.1 概念模型

三维地层概念模型能够更加真实地反映实际地层流体流动情况[12]，因此，本次研究将实际地层进行三维概化[图2(a)]。根据实际的勘察资料，目标地层主要为粉砂岩、细砂岩，中间夹杂着泥质粉砂岩和砂质泥岩，顶部盖层和底部基岩以泥岩为主。实际目标地层垂向深度为2 036.80～2 066.20 m，厚度共计29.40 m。为了减小侧向边界对模拟结果的影响，模型中x方向的长度和y方向的宽度分别设置为29.0 km和14.5 km。为了便于与水平地层的注入效果进行对比，在倾斜地层中将目标注入井垂直于地层布设。注入井的位置位于三维地层概念模型的中间部位[图2(b)]，规定注入井左侧和右侧区域分别为上半部分区域和下半部分区域。鉴于模型的计算精度和计算机的计算能力，模型xy平面上采取不等距网格剖分，x、y方向上剖分网格数量分别为51和26，并对注入井附近的网格进行加密剖分，距离注入井越远，网格尺寸越大，网格大小为1～1 000 m。注入井附近的网格设置为1 m，井筒为0.3 m；垂向上，根据地层孔渗变化共剖分了48层，模型共计63 648个网格。为了减小压力扩散到边界对模拟结果造成的影响，侧向边界设置为恒温恒压边界。鉴于顶部盖层和底部基岩岩层厚度比较大，孔隙度和渗透率都很小，因此，将顶板和底板设置为零流量边界条件。CO2初始饱和度视为0，初始CO2含量(质量分数，下同)为0。初始温度和压力条件严格根据目标注入井的实测资料设定，不同深度其值不同。盐度来自于实测地层水样的测定结果，其值为0.043。地层流体流动符合达西定律，其温度和压力分别为59.40 ℃和20.51 MPa。

图2 三维地层概念模型和xz平面Fig.2 3D Stratigraphic Conceptual Model and xz Plane

2.2 数学模型与参数确定

描述CO2注入过程的气液相质量守恒方程为

(1)

式中:F为质量或者能量通量；M为质量或能量；q为单位时间源汇；K为组分，其中K为w、g、s(w、g、s分别表示水相、气相和固相)；t为时间。

模型的初始条件包括温度、压力、盐度和气相饱和度。具体数学表达式为

(2)

式中:Ss为盐度，Sg为气相饱和度；ρw为水的密度；P为压力；x为长度；y为宽度；z为埋深；g0为重力加速度；T为温度;Ti为常数。

在CO2注入到深部地层以后，CO2在目标注入井附近聚集导致地层局部高压，驱动着CO2侧向运移；与此同时，压力也会迅速在地层中向外扩散。模型边界条件为

(3)

式中：Pi为常数。

本次研究将地层概化为48个子地层，同一地层内为均质各向同性，各子地层孔隙度、渗透率条件根据实测值进行设置(图3)。

图3 孔隙度和渗透率随埋深的变化Fig.3 Changes of Pore and Permeability with Buried Depth

根据测井资料中模型的温度和压力，可进一步得出毛细压力和相对渗透率。毛细压力计算采用Van Genuchten模型[25]。其表达式为

Pcap=-P0((S*)-1/λ-1)1-λ

(4)

S*=(Sl-Slr)/(Sls-Slr)

(5)

式中:Pcap为毛细压力；P0为毛细排驱压力(突破压力)；Sl为液相饱和度；Slr为残余液相饱和度；Sls为忽略依附在岩石固体颗粒表面(不可流动的)水分子后的液相饱和度，一般取值为1.0；λ为通过实验获得的参数，一般使用经验值0.46；S*为中间变量。

毛细压力随着岩体孔隙度和渗透率的改变而改变。这种因孔隙度、渗透率变化导致的毛细压力变化可以利用Leverett J-function模型[26](J(Sw))来定量描述。其表达式为

(6)

式中:Sw为水饱和度；k为渗透率；φ为孔隙度；φ为表面张力；θ为气相和液相之间的接触角。

液相相对渗透率(krl)模型采用Van Genuchten模型。其表达式为

(7)

气相相对渗透率(krg)模型选用Corey模型[27]。其表达式为

(8)

(9)

目标地层岩石热传导系数、比热容、密度等参数均通过实验测试获取，具体数值见表1。

2.3 模拟方案和模拟软件

根据示范场地的实际资料，白垩系东沟组地层倾角的最大值为6°。本次模拟设计了4种方案，设置的地层倾角分别为0°、2°、4°和6°。目标注入井与地层采用4 MPa的注入压差进行目标地层的全层位定压注入，连续注入20年，停注后模拟80年，总模拟时间为100年。

表1 模型中参数取值Tab.1 Parameter Values in the Model

图4 注入20年后与初始时刻不同地层倾角下压力差分布Fig.4 Pressure Difference Distributions in Different Stratigraphic Dip Angles After Injection for 20 Years and at Initial Moment

图5 100年后不同地层倾角下压力分布Fig.5 Pressure Distributions in Different Stratigraphic Dip Angles After 100 Years

图6 注入20年后不同地层倾角下CO2饱和度分布Fig.6 Distributions of CO2 Saturation in Different Stratigraphic Dip Angles After Injection for 20 Years

本次研究采用的数值模拟软件为TOUGH2中的ECO2N模块[28]。TOUGH2软件是一个可以用于研究三维孔隙和裂隙介质中多相、多组分非等温流体流动过程的数值模拟程序，主要应用于地热封存工程、饱和或非饱和介质中的流体流动和组分迁移、核废料的隔离研究、环境评估与修复以及CO2地质封存等[3]。

3 结果分析

3.1 压力空间分布与演化

CO2的注入势必会引起地层压力的升高。如图4所示，注入20年后，近井处地层压力出现了明显的上升，随着与井筒间距离的增大，压力变化逐渐变小。水平地层中，注入20年后，压力侧向影响范围扩展到了8 km之外，但并未影响到模型边界。倾斜地层中，注入井两侧压力不再呈对称分布，并且压力分布的差异性随着地层倾角的增大逐渐增大。这是因为CO2的密度小于水，倾斜地层中CO2会在浮力的作用下逐渐向浅部运移(模型的上半部分区域)。同样受地层倾角的影响，在浮力作用下，CO2在倾斜地层中的扩散速度将会更快，从而使得停注后压力更快地消散。如图5所示，停注80年后，水平地层中注入井附近依然存在明显的压力聚集区。相比之下，倾斜地层中压力分布在停注80年后基本恢复正常，并未出现明显的压力聚集区。

图7 100年后不同地层倾角下CO2饱和度分布Fig.7 Distributions of CO2 Saturation in Different Stratigraphic Dip Angles After 100 Years

3.2 CO2空间分布与运移距离

图6为定压注入20年后不同方案中超临界态CO2在地层中的分布情况。由模拟结果可以看出，CO2饱和度的分布模式与压力分布模式相似，但是影响范围要明显小于压力影响范围。注入20年后，水平地层中CO2最大侧向运移距离约为2 860 m；而地层倾角分别为2°、4°和6°时的倾斜地层中，CO2最大侧向运移距离分别为3 060、3 100和3 200 m。图7为停注80年后不同方案中超临界态CO2在地层中的分布情况。由图7可见：随着时间延长，不同倾角地层中CO2侧向运移距离差异显著；地层倾角越大，CO2在地层中的最大侧向运移距离差异越显著。例如，停注80年后，当地层倾角由0°增加到6°时，CO2最大侧向运移距离由3 050 m升高到了7 200 m。倾斜地层中，CO2在浮力作用下快速向浅部运移，虽然在一定程度上有助于压力的消散，促进了CO2的持续注入，但同时这也将显著增加CO2向浅部含水层甚至是大气泄露的风险。

图8展示了整个运行期间(注入20年+停注80年)倾斜地层中CO2向浅部运移过程(由于模拟期间CO2并未穿透上覆盖层，所以认为水平地层中CO2向浅部的运移距离为0 m)。如图8所示，停注80年后，地层倾角分别为2°、4°和6°的倾斜地层中CO2向浅部最大的垂向运移距离(即侧向运移距离在垂向上的投影距离)分别为470.38、595.81和752.60 m(注入层位初始埋深为2 040 m)。由此可见，随着地层倾角的增大，CO2向浅部泄露的风险显著升高。当倾斜地层出露地表时，随着时间延长，CO2可能沿着倾斜地层泄露到大气中。

3.3 CO2封存量随时间的变化

CO2在深部地层的存在形式主要包括超临界态、溶解态和固态(生成矿物沉淀相)，其中超临界态在短时间尺度内占主导地位[29]。停注之后，CO2逐渐溶解，超临界态CO2的比例逐渐降低，溶解态CO2缓慢增加。赵宁宁等研究结果表明，固态CO2形成并占据主导地位所需时间较长，通常需要数千年以上[30]。由于本文研究时间为100年，故不考虑矿物固碳作用的影响。

由表2可知，注入20年后，地层倾角分别为2°、4°和6°的倾斜地层相比于水平地层，CO2的总封存量分别提高了13 294、29 016和65 512 t，总封存量增幅分别为0.274%、0.598%和1.351%。由于在不对称地层压力和浮力的共同作用下，中高渗倾斜地层更有利于注入压力的消散，从而地层倾角增大提升了CO2的总注入量，这一结果与前人基于鄂尔多斯CO2地质封存示范工程中的低孔低渗地层倾角变化对封存量的影响规律[19，31]正好相反。由此可见，中高渗倾斜地层和低渗倾斜地层的流体运移规律存在很大的差异性。不同地层倾角下CO2注入速率随时间的演化关系造成这一现象。从图9(a)可以看出，在注入初期(前8年)，注入井附近空隙空间未被CO2充分填注，在压差一定的情况下，水平地层单位时间储存CO2能力更大，表现为水平地层CO2的注入速率大于倾斜地层，因此，水平地层与倾斜地层CO2注入速率差值为正值；之后，随着压力积累，注入井附近空隙空间被充分填注，由于倾斜地层中压力消散更快，更有利于CO2的持续注入，所以倾斜地层中CO2的注入速率逐渐超过水平地层。在注入速率的差异性变化前提下，CO2的总注入量表现出非单调性变化[图9(b)]。在注入17.5年之前，水平地层CO2的总注入量大于倾斜地层；在第17.5年时，各方案CO2的总注入量基本相等；之后随着注入时间增加，地层倾角越大的倾斜地层CO2的总注入量越大，表现出地层倾角对CO2的总注入量具有促进作用。

表2 注入20年后水平地层和倾斜地层中CO2总封存量Tab.2 CO2 Total Storage Amount of Horizontal and Sloping Strata After Injection for 20 Years

图8 不同地层倾角下CO2最大垂向运移距离变化Fig.8 Variations of CO2 Maximum Vertical Migration Distances in Different Stratigraphic Dip Angles

图9 不同地层倾角下CO2注入速率差值及总封存量差值Fig.9 Differences of CO2 Injection Rate and Total Storage Amount in Different Stratigraphic Dip Angles

4 结 语

(1)在地层压力分布与传递方面，CO2注入引起地层压力的升高，近井处地层压力最先上升，然后压力逐渐向远离井的方向传递；倾斜地层由于受地层倾角的影响，地层压力分布呈不对称形态，压力场的传播范围远大于水平地层。

(2)在运移距离方面，地层倾角的存在会促进CO2向浅部运移；倾角越大，CO2向浅部含水层及大气中泄露的风险也越大。

(3)在封存量方面，中高渗地层在长时间持续注入条件下，地层倾角越大，地层压力消散越快，CO2总封存量越大。这一研究结果与前人基于低渗倾斜地层倾角变化对封存量的影响取得的认识相反。

(4)中高渗倾斜地层虽然有利于地层压力的消散和注入量提高，但同时也显著增大了CO2向浅部泄漏的风险。因此，在CO2地质封存选址过程中，应综合考虑倾斜地层对封存效率及安全性的影响。

(5)本次研究模拟时间较短，未考虑水岩气反应导致的岩石溶蚀与沉淀对地层物性特征的改变，相关研究有待进一步展开。