（一）研究现状
地下多相、多组分流体运移数值模拟是在质量和能量守恒的基础上，建立的多相流体运动以及反映地球化学运移扩散的数学模型，通过离散建立大量的线形或非线形方程组，然后利用计算机计算求解，再通过图像显示模拟结果，达到对工程问题和物理问题乃至相关其他问题研究的目的。CO2地质封存数值模拟就是利用计算机模拟的方法，来解决CO2进入地质封存系统后运移、转化、水-岩-气之间的相互反应、CO2泄漏对浅部含水层影响及诱发的储盖层物性变化等一系列问题，从而指导CO2地质封存工程的实施。
目前，国内外已开展的关于CO2地质封存数值模拟的研究工作包括以下几个方面：
1.超临界CO2-水多相流体运动模拟
Pruess等（2003）模拟了均质各向同性咸水含水层中以恒定流量灌注CO2条件下，灌注井井周非等温径向流情况。当忽略重力和惯性力效应时，模拟结果中存在相似变量ζ＝R2/t（其中，R为径向流动距离，t为时间），CO2饱和度、溶解CO2质量分数、沉淀盐的体积分数及流体压力都是相似变量的函数。这与O' Sullivan（1981）及Doughty和Pruess（1992）的结果一致。两相流的模拟考虑了CO2和水的相对渗透率及毛细管力作用问题（Van Genuchten，1980），考虑了流体密度、黏度和CO2溶解性随压力、温度和含盐度的变化，以及盐的沉淀导致含水层渗透率的减小等因素。
Doughty和Pruess（2004）利用Fro咸水含水层封存CO2监测资料，反推了CO2灌注后发生的物理和化学过程。他们采用TOUGH2数值模拟软件对两相（液、气）三组分（CO2、水和溶解NaCl）系统进行模拟。考虑15MPa和65℃条件下，超临界CO2在咸水中为非混溶流体，并能部分溶解于咸水的情况，分析了多相流系统边界设定的影响及相对渗透率取值问题，即模拟中对侧向边界的设置为均开（或均闭），结果导致压力的模拟结果与实际相比过低（或过高）。研究表明，由于上覆断层对CO2的封堵作用，侧向边界对CO2弥散羽的影响不大。模拟结果还显示，相对渗透率函数对CO2弥散羽的演化有很强的影响。如何确定一个合适的相对渗透率以表征CO2灌注咸水含水层的变化，仍是亟待解决的问题。Doughty和Pruess模拟了两种CO2残余饱和度条件下的CO2羽扩展随时间的变化，发现存在较大差异。残余饱和度较大的情况下，CO2羽表现出紧缩状，在浮力作用下运移较慢；相反，在残余饱和度较小的情况下，CO2羽流弥散很快，溶解性显著提高。
2.多组分反应地球化学运移模拟
水-砂岩-CO2相互作用往往形成一系列次生矿物，或次生矿物组合。Worden et al.（2006）通过岩石学和CO2灌注长石砂岩的地球化学模拟表明，北海Magnus油田上侏罗统浊积亚长石砂岩中的铁白云石、高岭石和石英可能具有成因联系。其中，铁白云石中的碳来自有机成因的CO2。Watson et al.（2004）通过CO2气与CH4气储集砂岩的比较岩石学研究，证实澳大利亚Otway盆地Ladbroke Grove CO2气储集砂岩中与CO2气灌注有关的次生矿物组合为铁白云石-高岭石-次生石英。
Xu et al.（2005）利用一维砂岩-页岩系统模型模拟了储层中灌注的CO2与矿物发生的化学反应过程，以及对储层环境的影响。模拟显示，在砂岩环境下，CO2主要被方解石所固定，而方解石的沉淀导致孔隙度减小，进而导致渗透率相应减小。10万年间，砂岩封存能力达到90kg/m3的封存能力，这些被矿物固定下来的CO2可以永久封存。Zwingmann等运用地球化学模拟软件EQ3/6进行的水-矿物-CO2相互作用模拟也表明，若将CO2灌注到日本本州岛中北部新潟盆地更新世灰爪组砂岩，CO2以溶于水和形成碳酸盐矿物两种形式封存，其中后者封存量最大为21.3mol/kgH2O，可达总封存量的90％，形成的碳酸盐矿物中也出现了片钠铝石。
3.耦合岩石力学模拟
从目前发表的论文及各国研究计划的综合报告上看，在CO2咸水含水层封存研究方面，对于CO2运移机制的分析和模拟很少考虑应力场的耦合作用。事实上，CO2灌注压力和超临界CO2的浮力作用将改变地层应力状态，即CO2在上浮运移和侧向扩散过程中，孔隙压力可能会对原始裂隙和断裂产生影响；CO2在咸水含水层中的长时期（千年级尺度以上）的封存，将改变含水层的地球化学状态，CO2-咸水-含水层矿物的化学作用将可能导致岩体力学和水力学性质发生变化。
日本因位于4大板块交界处与环太平洋构造带中，活断层密集发育，地震频繁，地应力分布复杂，在CO2地质封存评价方面，非常重视CO2地质封存的力学稳定性研究（李琦等，2002；李小春等，2003）。李琦等（2002；2004；2006）提出了一个考虑初始地应力场、灌注压力、CO2浮力及含热传导作用的热-水-力（THM）耦合模拟框架，考虑盖层底部附近存在不同倾角断层的二维平面应变地质封存问题。采用有限元算法，对灌注CO2流体对断层稳定性的影响进行模拟分析。计算结果表明，为了避免断层位移需要特别注意对灌注压力的控制，因为CO2灌注压力对断层滑动的影响远大于CO2羽流浮力带来的影响。停止灌注CO2后，CO2羽流的上升则成为应力场扰动的主要因素。
（二）主要软件介绍
近年来，计算机模拟技术在许多研究领域得到了广泛的应用，开发出了许多优秀的模拟软件和程序。同样，可用于研究CO2地质封存的数值模拟软件也很多，主要有PHREEQC、GEM、ECLIPSE、TOUGHREACT、PetroMod、MUFTE-UG和NUFT等，它们都有各自的特点和适用性。在进行数值模拟之前，需对这些数值模拟软件进行评价分析，选择适用于所研究问题的模拟软件。现对国际上常用的几款软件简介如下。
1.PHREEQC
PHREEQC是一款用于计算多种低温水文地球化学反应的计算机软件。以离子缔合水模型为基础，PHREEQC可完成以下任务：(1)计算物质形成种类与矿物的溶解饱和指数；(2)模拟地球化学反演过程；(3)计算批反应与一维运移反应。另外，与多组分溶质-运移模型耦合的PHREEQC可生成PHAST，一个用于模拟地下水流系统的三维反应-运移模拟器。但由于PHREEQC是在单相水流的基础上建立的模型，因而不能模拟超临界CO2-水的两相流运动。
PHREEQC最简单的应用就是计算溶液中各种化学物质的分布，以及溶液中矿物与气体的饱和状态。反演模拟功能可推导和量化在流动过程中，能够反应化学物质变化的化学反应方程。PHREEQC可处理的反应方程包括建立矿物、表面配合物、阳离子交换剂、土壤溶液、气体组分单位分压、给定压力或给定体积气相间平衡的物质运移反应。在模拟这些均衡反应的同时，PHREEQC还可以模拟动力化学与生物反应，以及模拟从简单的线性衰变（代谢物降解或放射性衰变）到复杂的依赖于溶液化学组成和微生物数量确定的反应速度。这些反应处理功能可在批反应模拟或一维对流、弥散、反应型运移模拟中使用。
2.GEM
GEM v.2009.13（Nghiem et al.，2004）是一款用来模拟利用CO2和酸性气体提高石油采收率的模拟器，该模拟器完全耦合了地球化学组成状态方程。GEM采用一步求解法进行状态方程的求解。GEM可以用来模拟：对流和弥散流体、油（或超临界CO2）、气和咸水间的平衡、水相物种间的化学平衡，以及矿物的动态溶解和沉淀。该模拟器采用自适应的隐式离散技术利用一维、二维或者三维模型来模拟孔隙介质中溶质的运移。油相和气相用一个状态方程来模拟，气体在水相的溶解度采用亨利定律模型来计算。水通过蒸发进入到气相、盖层的穿透、热效应和裂隙的封闭作用也可以利用GEM来模拟。
3.ECLIPSE
ECLIPSE是一个并行化的可以模拟黑油、组分、热采等问题的成熟软件。1994年，胜利石油管理局引进了ECLIPSE油藏数值模拟串行软件，广泛开展了从油藏到气藏，从常用油田到特殊油气田、从常规模拟研究到特殊模拟研究等多方面的应用。主要模块有主模型、黑油、组分、热采、流线法、运行平台和ECLIPSE Office等。
ECLIPSE是一个商业软件，在使用中其内核部分是封闭的，使用者只能将其作为一个“黑箱”来操作。其不足之处有：不可能免费的获得和随意地使用和修改；无法耦合最前沿的地质流体热力学模型；无法加入更多影响因素来研究具体问题。因此，ECLIPSE不适宜用于前沿科学研究。
4.TOUGH2/TOUGHREACT
TOUGH2是Transport of Unsaturated Groundwater and Heat（非饱和地下水流及热流传输）的英文缩写，是一个模拟一维、二维和三维孔隙或裂隙介质中，多相流、多组分及非等温的水流及热量运移的数值模拟程序。TOUGH2使用积分有限差分（Integral Finite Differences，IFD）（图3-8）的方法来解决多相流动和多组分化学运移模拟中的空间离散化问题（Pruess et al.，1999s；Xu et al.，2004）。为了满足大规模计算需要，Zhang et al.（2008）开发了TOUGH2的平行计算版本，即TOUGH2-MP。
该方法在对地质介质的离散化上是比较灵活的，允许使用不规则的网格，十分适合对多区域非均质系统和裂隙岩石系统中流体流动、运移和水岩相互作用的模拟。而对于规则的网格剖分，积分有限差分方法相当于传统的有限差分法。其中，对于任意区域Vn，它的质量（对于水、气体和其他化学组分）和能量（对于热）守恒方程可以用积分的方式（式3-5）表达：
图3-8 积分有限差分法中的空间离散化和几何参数数据构成图
中国二氧化碳地质封存选址指南研究
式中：下角标n为表示一个单元格；下角标m为表示和单元格n互相连接的网格m；Δt为时间步长；Mn为单元格n的平均质量或能量密度；Anm为单元网格n和m交界的面段；Fnm为通过面积为Anm的质量或能量通量；qn为单元格n内单位体积的平均源汇率。
许天福等（1998）在TOUGH2的框架基础之上，增加了多组分溶质运移和地球化学反应的模拟功能，形成了一套较为完善的可变饱和地质介质中非等温多相流体反应地球化学运移模拟软件——TOUGHREACT。该软件不仅包括了TOUGH2的全部功能，而且适用于不同温度、压力、水饱和度、离子强度、pH值和氧化还原电位（Eh）等水文地质和地球化学条件下的热-物理-化学过程。还可以应用于一维、二维或三维非均质（物理和化学的）孔隙或裂隙介质中的相关数值模拟研究。在理论上可以容纳任意数量的以固相、液相或气相存在的化学组分（但是在实际模拟中会受到计算能力和计算时间等硬件条件的限制），并且考虑了一系列化学平衡反应，如溶液中的配合反应、气体的溶解或脱溶、离子吸附作用、阳离子交换及受平衡控制或反应动力学控制的矿物溶解或沉淀反应等。可以说TOUGHREACT、是TOUGH2的升级版，近年来在世界范围内CO2地质封存研究和工程实践中得到了广泛的应用。
除包含TOUGH2所有的功能外，TOUGHREACT还可以应用于一系列的反应性流体和地球化学迁移问题。比如：(1)伴随Kd线性吸附和放射性衰变的污染物迁移问题；(2)在周围环境条件下，自然界中地下水的化学演变；(3)核废料处置地点评估；(4)深部岩层的沉积成岩作用；(5)CO2地质处置。多相流体运动，多组分反应地球化学，各种封存形式封存量以及随时间空间变化；(6)矿物沉积（如表生铜矿富集）；(7)自然和补给环境下热水系统中的矿物变化。
通过最近几年相关研究者的不懈努力，TOUGHREACT在实际应用中得到了进一步完善和提高，增加了部分新功能，如水相内部反应动力学和生物降解作用，改进了矿物-水反应表面积计算方法，以及气-水反应中气的活度系数的修正等。
5.PetroMod
由德国IES（Integrated Exploration System）公司研究开发的PetroMod多组分、多相态的多维含油气系统模拟软件综合平台已被世界石油业所公认。该软件融入了断层活动性、盐丘上涌和刺穿、火山岩的侵入、气体扩散效应、油气水三相运移和油气吸附模型等相关技术。
该模拟软件平台推出和采用的油气运移组合模拟算法（Hybird）是当今最先进的油气运移模拟算法，既可以保证模拟的精度，又可以极大地提高模拟的运算速度。其中的PetroFlow3D用于油气运移、聚集、圈闭和散失等情况的模拟，同时PetroCharge Express为我们提供了基于图件的油气运移和圈闭模拟的快速分析工具。
6.MUFTE-UG
MUFTE-UG是MUFTE和UG.MUFTE的结合。MUFTE即多相流（Muliphase Flow）、运移（Transport）和能量（Energy）模型。该软件包主要包括物理模型概念和孔隙裂隙介质中等温和非等温多相多组分流动和运移过程的离散方法（Helmig，1997；Helmig et al.，1998）。它能对裂隙孔隙介质进行离散性描述（Dietrich et al.，2005）。UG是非结构性网格（Unstructured Grid）的缩写，它提供的数据结构能快速解算以平行、自适应多网格法为基础的离散型偏微分方程。具有模块化结构的MUFTE-UG很容易解决各种有特殊要求的问题。
模块化结构的MUFTE-UG具有许多不同的环境与技术应用。例如，在环境应用领域，MUFTE-UG能够模拟如下两个问题。
（1）NAPL（非液相流体）向饱和与非饱和土壤的渗流。优化改进的修复技术在MUFTE中具有广泛的研究和发展空间。
（2）地下CO2的消散。CO2以高温高压灌注地表以下几百米的地层中，MUFTE-UG可用于非均质含水层（对流和弥散运移）中羽状体演化评价，伴随温度效应（由于膨胀和压缩）和组分间相互溶解（卤水和CO2）。
7.NUFT
NUFT（Nonisothermal Unsaturated-Saturated Flowand Transport model）是一套用来解决在多孔介质中多相、多组分非等温流动和溶质运移过程中地下污染物运移的数值解法器。此软件利用简单的代码来利用通用的实用程序和输入文件的格式。最近，此代码在Unix和DOS系统下运行成功。
该程序利用一套完整的有限差分空间离散法求解平衡方程组。每一个时间步长内利用Newton-Raphson方法求解非线性方程组，而在每一步迭代过程中利用直接解法和预共轭梯度法求解线性方程组。该模型可以解决一、二和三维水流及溶质运移问题。将来该模型会耦合进毛细滞后、非正交网格离散、有限单元剖分和固体非线性等温吸附等功能。
（三）研究方法
通常情况下，CO2地质封存数值模拟包括以下主要过程。
（1）建立概念模型：根据各种方法获取的实际资料来概化和建立CO2地质封存概念模型，包括边界范围、地层或储盖层高程、储盖层确定、参数及分区、源汇项、主要物理化学过程以及模型维度（一维、二维和三维）。
（2）建立数学模型：建立一套描述深部咸水层中多相流动和多组分反应性溶质运移的偏微分方程组，包括初始条件和边界条件问题。
（3）模型离散化：把概念模型中的各种信息通过网格剖分进行离散，形成大量的网格单元，然后通过有限差分、有限单元和积分有限差分等方法转化成单元的质量和能量守恒方程组，再用多种方法将非线性方程组线性化，形成线性代数方程组，然后求解方程组。
（4）模型识别和校正：根据模型计算结果和实际监测数据进行对比拟合，适度合理调整参数，使模型能够综合反映实际情况。在历史拟合过程中出现较大误差，应重新检查概念模型，修正概念模型。对所建模型进行参数敏感性分析，对于较敏感的参数应该慎重选取，甚至需要做大量的试验来确定。
（5）模型预测：建立了可靠的模型后，便可以进行模拟预测。
数值模拟的关键是地质模型概化、计算精度和计算速度。由于计算的精度取决于离散的程度，而离散的程度又决定了计算的速度，这是一对矛盾，要根据解决问题的需要来选择离散化的程度和计算速度。
CO2在储层中的运移、溶解以及与围岩的化学反应形成了一个多相、多组分的反应体系，涉及的主要数学方程有超临界CO2-水的两相流体运动控制方程、溶质运移控制方程和化学反应方程等。建立数值模型时，通常采用有限差分法、有限元法和积分有限差分法等。
由于实际应用时多采用已有的数值模拟软件对CO2地质封存的某一过程进行模拟，不涉及软件的开发及程序代码的编写，只需根据研究的需要选择合适的软件进行模拟预测，而软件一旦选定，数学模型和数值模型基本上已经确定。