增强型地热储层多场耦合数值模拟研究进展

张炜韬，韩东旭，李敬法，李庭宇，孙东亮，巩 亮，宇 波

(1.中国石油大学(华东)新能源学院，山东 青岛 266580；2.北京石油化工学院机械工程学院，北京 102617；3.中国空气动力研究与发展中心低速空气动力研究所，四川 绵阳 621052)

地热能是产生和储存于地球内部的一种热能，主要来源于地球板块构造运动释放的能量和天然放射性同位素衰变产生热量，可分为浅层地热、水热型地热(地下热水)和干热岩地热资源三种[1].浅层地热能是蕴藏于地表下约几百米处的岩土、地下水、地表水中具有可开发价值的热能，根据纬度不同，温度通常在10 ℃～21 ℃之间，具有分布广、可再生、可就近开发等特点，被广泛用于建筑供暖、制冷、制取生活热水等方面[2].水热型地热能是地下裂隙岩体、孔隙或深部断层中的水被附近高温熔岩或岩浆加热形成热液，渗流出地表或被开发抽取出来的热水，根据储层条件不同，温度通常在环境温度到350 ℃以上，被广泛用于供暖、理疗、种植养殖等方面[1，3]，水热型地热能的地质环境、水文环境、热储介质等都比较复杂[4]，热量提取有限，一般采用注冷水或地埋管方式取热[5，6].干热岩是埋藏于地层3 km～10 km范围内、含天然裂缝的结晶岩热储，不含或仅含少量流体，温度高于180 ℃，是地热能资源的重要赋存形式之一，其热能是在当前技术条件下可以利用的岩体，也是地热资源研究的热点和前沿[7].干热岩的开发利用对全球能源多元化利用、助力我国实现“碳达峰”与“碳中和”有重要价值[8].

干热岩热储层多由脆、硬且岩石力学强度高的基底花岗岩组成，具有低渗透性.在上世纪70年代，美国洛杉矶国家实验室首次提出EGS的概念[9]，目的是通过压裂改变储层渗透率达到提高产热率.该过程主要包括钻井到干热岩储层、压裂产生有效裂缝或使原有裂缝扩展或两者有机结合、注冷水循环取热.在循环取热时，流体在储层中运移、换热不仅改变岩体孔隙压力，也使温度分布发生改变，同时在水力、温度的诱导作用下，岩体固有应力状态也会发生改变，从而发生变形，甚至破裂、崩坍、发生地震等[10]，而在高温、高压环境下，流体与岩体长期作用，还会使岩体壁面的硅质矿物发生溶解/沉淀反应，造成岩体的孔隙结构、渗透率发生改变.可见，EGS取热时，流(Hydraulic：H)、热(Thermal：T)、力(Mechanical：M)、化(Chemical：C)等物理场是相互耦合、随时空演化的过程如图1所示，研究和理解各物理场耦合过程对EGS生产实践非常重要[11].

图1 EGS物理过程示意图：黑色实线为水力裂缝

理论、实验和数值模拟是科学讨论上述问题的三大重要手段.对于实验而言，由于深部干热岩热储层的实际尺寸较大，很难开展有针对性的室内研究，现场许多先进的检测手段对储层中裂缝动态扩展过程、裂缝总面积、裂隙岩体的渗透率等分布难以给出准确合理的估计和解释[12]，理论和数值模拟被认为是研究相关问题的有效途径.EGS储层中存在着复杂的缝网系统，是典型的裂缝型储层.尽管对裂缝型储层的表征在油、气藏工程等领域已有多年的发展，但应用于增强型地热系统中需要在其基础上叠加复杂的THMC多场耦合模型，目前仍待深入发展与完善.本文主要从裂缝型储层表征、多场耦合模型及多场耦合快速求解方法三个方面展开综述，相关内容对未来EGS多场耦合数值模拟研究具有重要的借鉴和参考价值.

1 裂缝型储层表征研究

EGS储层中存在着复杂的缝网系统，是典型的裂缝型储层，天然裂缝或压后形成的水力裂缝是EGS增产关键，准确描述与模拟裂缝系统是精细产量与产能预测的核心.随着致密油气田、碳酸盐岩油气藏、页岩气以及地热储层开发的需要，研究人员开发了多种可以表征裂缝型储层的方法，主要分为：等效介质方法、双重孔隙介质方法、离散裂隙方法以及嵌入式离散模型，各个表征方法特点如表1所示.

表1 裂缝型储层表征方法[13]

等效介质方法不显示描述裂缝，而是将裂缝的孔隙度和渗透率等效为岩体的原生孔隙度和渗透率，采用经典的连续介质理论建模.Long等[14]基于裂隙岩体和均质、各向异性的多孔介质渗流理论，指出裂隙可以被等效成连续介质的两个条件：(1)当试验样本体积变化一个小量时，介质等效渗透率变化不大；(2)当水力梯度方向改变时，等效渗透率张量能够正确地预测通量.Oda[15]发展出了一种裂隙岩体流-固耦合的等效介质模型，并作了两点假设：(1)裂缝占整个弹性连续体很小；(2)每条裂缝可以由两个弹簧连接的平行平板替代.在实际应用中，等效介质模型需考虑两个关键问题：(1)裂缝渗透率的等效方法；(2)保证等效后储层参数的有效性.对于简单问题，该方法能够得到令人满意的结果，但是对于裂缝主导的流动，其误差较大.

双重孔隙介质方法认为任意表征单元内同时存在基岩系统和裂缝系统，岩块中的孔隙是流体的主要储存场所，而裂缝是流体的主要流动通道，两者形成相互独立又相互联系的水动力学系统.表征单元空间任意一点既代表基岩系统，又代表裂缝系统，可同时定义两个不同的孔隙度和渗透率，两个系统间的流体交换与其压力差值密切相关[16].该方法最早由Barenblast等[17]提出，随后Warren等[18]对其进行了简化，简称为BW模型，模型中的岩块可以理想化为层状、圆柱、圆球等形状，不同形状对应的形状因子不同.当储层中存在大量且连通性较好的裂缝时，采用该方法计算量小且效率高.但是，计算过程需要确定基岩系统和裂缝系统间的窜流函数，而窜流函数与形状因子有关，形状因子只在特定的流动阶段保持不变[19]，因此，对于裂缝明显主导的流动问题，该方法往往显得“力不从心”.此外，在化学溶蚀过程中，也不能体现大裂缝对流体流动的主导作用[20].

离散裂隙类方法以显示表征裂缝为核心，对裂缝采用降维处理的策略，该类方法能够避免对裂缝附近进行网格加密，提高计算效率，降低计算内存的消耗，易于表达裂缝与基岩之间的窜流关系，可方便地处理局部强非均质性问题.但是，当裂缝数量较大时，计算量会增大，且随着裂缝分布的复杂性增加，网格划分也将变得困难.根据物理问题不同，离散裂隙类方法发展有离散裂缝模型(Discrete Fracture Model，DFM)、离散裂缝网络模型(Discrete Fracture Network，DFN)以及嵌入式离散裂缝模型(Embedded Discrete Fracture Model，EDFM).

DFM模型由Noorishad等[21]首次提出，随后Baca等[22]对其进一步完善.该方法实施过程中假设裂缝开度较小，忽略流体在开度方向上的梯度变化，采用线单元表征裂缝.网格划分时，保证了裂缝能沿基质网格边界生成，在基岩和裂缝交接处，二者共用同一个节点，基于有限元方法，将基岩和裂缝方程整合到一起进行求解.Aliyu等[23]基于DFM和局部热平衡假设，建立了地热开发过程的热-流耦合模型，并分析了温度对密度、黏度、热导率和比热容的影响，讨论了人为控制因素和天然因素对地热开发过程的影响.Karvounis等[24]在自适应DFM模型的基础上，通过引入核函数来表征不同裂缝之间的流量、能量交换，建立了干热岩开发过程的流动和换热模型.Slough等[25，26]将DFM推广到三维情形，模拟了多相流问题.对于三维问题，裂缝被简化为二维平面，网格划分时，裂缝与基岩网格面边界重合，需要大量的网格来模拟裂缝的影响，当网格较少时，精度很低，甚至得到错误的结果.

DFN模型是目前比较先进的裂缝型储层建模技术，裂缝数据来自于压裂、微地震事件、测井、示踪剂等方面，通过整合各个尺度上的裂缝，利用地质统计学方法，建立综合反映储层实际特性的裂缝模型[27-29].在DFN模型中，将裂缝分为大尺度、中尺度和小尺度三种类型.大尺度裂缝和断层的数据主要来源于微地震检测，因此可确切地知道裂缝和断层的方位、走向和倾角等信息，不需要随机生成.中尺度和小尺度裂缝根据岩心资料、测井数据进行先验分析，粗略统计裂缝的大小、开度、分布、方位等，采用确定和随机方法相结合生成裂缝.然后根据示踪剂及动态模拟曲线，判断裂缝之间的连通情况，通过反复修正模型参数，实现最大限度使裂缝模型逼近真实储层.应用DFN模型可以帮助地质学家、压裂工程师判断裂缝或断层的状态.但是，由于实际裂缝型储层具有很强的非均质性，如何建立各种复杂工况下的DFN模型，并精确地指导生产实践是该模型未来的发展方向.

DFM模型和DFN模型虽然能够较精细地模拟裂缝主导的流动问题，但当裂缝数量较多时，网格生成困难，且当裂缝不沿坐标方向或不垂直于坐标方向时，会引入非结构化网格.为了解决这一问题，Lee等[30-32]和Li等[33]提出了EDFM模型，该模型将裂缝降维嵌入到基岩网格系统中，如对于二维问题，裂缝被看作一维线段，对于三维问题，裂缝被看作二维平面.该方法使基岩网格划分时可以采用结构化网格，而裂缝网格划分时不用考虑基岩网格，如图2所示.基岩系统和裂缝系统、裂缝系统与裂缝系统、及各自系统内的物理量交互，通过交换系数和传递系数进行.为了使EDFM模型处理交叉裂缝通量交换更为方便，Karimi-Fard等[35]借鉴了电路中串并联电阻的等效计算方法，可避免在交叉裂缝处引入控制容积，提高了计算效率.对于单相流动问题，EDFM模型计算的压力场与DFM模型相当.但对于多相流问题，流体的相饱和度在裂缝面两侧具有不连续性，会导致EDFM模型在处理多相流问题时存在困难[36].为此，学者们提出了一些改进方法，如投影基的嵌入式离散裂缝模型(pEDFM)[37，38]、单调的嵌入式离散裂缝模型(mEDFM)[39，40]、三维鲁棒性的嵌入式离散裂缝模型(3D-EDFM)[41]等.

图2 二维嵌入式离散裂缝模型示意图

综上，等效介质方法和双重孔隙介质方法主要应用于裂缝影响不大的流动和换热问题，对于裂缝主导的流动问题，主要应用的是离散裂隙类方法.但由于实际物理问题中裂缝错综复杂，每条裂缝的粗糙度、开度、渗透率及孔隙度等均不相同，因此想要完全准确地应对实际复杂的裂缝型储层物理问题，对于裂缝型储层的精细化表征还有很长的路要走.

2 多场耦合模拟进展

EGS多场耦合模拟是室内试验、现场测试的必要补充，可为实验研究、场地测试、工程应用等提供一些准确、可靠的预测和有价值的参考.因此，其耦合结果应该能够准确地表征EGS生产中的变化，例如，热诱导裂缝宽度改变问题、流体注入导致裂缝剪切滑移问题、矿物溶解/沉淀导致储层的热耗散和渗透率改变问题等[42].本节主要从模拟软件发展、涉及化学场的THMC耦合模型两方面进行叙述.

2.1 模拟软件

在多场耦合模拟中，对岩体在复杂载荷和多场耦合作用下的力学行为，常采用的数值方法有：连续介质力学方法、非连续介质力学方法、连续-非连续介质力学耦合方法.连续介质力学包括有限差分法(FDM)、有限单元法(FEM)和边界元法(BEM)，适用于研究连续体的小变形和损伤断裂；非连续介质力学包括离散元方法(DEM)、DFN，适用于研究非连续体的破坏和运动；连续-非连续介质力学耦合方法包括FEM-BEM耦合、DEM-BEM耦合、FEM-DEM耦合等，适用于研究岩体从连续到非连续的全过程变化[43].对储层中的流动、传热、组分输运、化学反应等物理过程，较多采用有限差分法(FEM)、有限体积法(FVM)、FEM或这些方法的组合和改进.

国内外学者针对裂隙岩体在多场耦合方面倾注了大量心血，一些可用来研究EGS多场耦合作用的商用和开源求解器被陆续研发出来，如表2所示.这些模拟器有的高度专业化，可实现快速、并行的求解，但只能有效地处理特定物理场的耦合问题；有的可实现顺序迭代多场耦合，但是其耦合顺序固定，或者不允许耦合其他物理场；有的可实现部分场的全耦合模拟，但是由于前期耦合开发策略的复杂性，导致其耦合可拓展性受限.此外，多数模拟器仅限于宏观大尺度下的模拟，而微观、介观以及微-宏观之间多尺度的耦合模拟功能比较有限.对于EGS的多场耦合模拟，学者们较多应用的为COMSOL Multiphysics和TOUGH+FLAC3D[45].

表2 地热多场耦合求解器[44]

2.2 耦合模型

由于EGS储层岩体为天然属性材料，其内部结构复杂，具有不连续性、各向异性、不均匀、非弹性等特点.长期受地壳隆起、沉降、滑移等地质构造影响，储层内岩体中的裂隙和孔隙中的流体介质在原地应力场和原地温度场的影响下将发生不同程度变化，如相态转变，这种不确定性是EGS储层的特征之一.如果在模化分析时，考虑所有不确定性因素，将难以用数学语言展开描述，若再同时考虑热、流、力、化各场的变化，将使THMC多场耦合模化分析的难度更大.研究人员对EGS的多场耦合模化分析经历了由简到繁、由易到难，逐步发展的过程.在维度方面，从二维到三维发展；在尺度方面，由宏观、微观分析向宏、微观多尺度发展；在多场耦合方面，经历包括了从热-流(TH)耦合、热-流-化(THC)耦合、热-流-力(THM)耦合到热-流-力-化(THMC)耦合的发展，并已经拓展到注入诱导地震的孔-宏观模型构建[65].建立EGS多场耦合模型时，需要根据具体的物理问题进行合理的假设和简化.例如，美国能源部地热技术办公室通过严格的物性假设、边界条件和初始条件等的设置，提出了七个影响EGS取热过程的基准问题，这些问题被用来评估不同求解器建立多场耦合模化分析的有效性和准确性，各个问题对应的耦合场，如表3所示.

表3 不同耦合过程中研究内容汇总(翻译)[42]

图3 THMC耦合过程[66，69]

由于THMC耦合模型中各偏微分方程具有非线性特点，研究人员也发展出了多种耦合场的求解方法，如全隐式求解、顺序求解、迭代求解等.全隐式求解需要大量计算资源，比较适合研究小而中等类的问题，而顺序求解和迭代求解比较普遍使用.Ogata等[71]将THMC耦合过程分为达西渗流、传热、力学变形、化学反应输运四个求解子过程，采用COMSOL Multiphysics按照顺序求解框架进行，每个时步内完成各场求解后，对裂缝宽度和渗透率进行一次更新.Taron等[72，73]基于显式迭代求解框架，首先在TOUGHREACT中完成THC的耦合计算，输出温度、孔隙压力和孔隙度等关键参数，然后将温度、孔隙压力计算结果通过插值方法，映射到对应的网格单元节点上，再采用FLAC3D完成力平衡迭代计算，每个时步迭代收敛后，计算并更新渗透率变化和化学反应造成的应变项变化.

到目前为止，EGS多场耦合模化的研究工作已历经40多年，取得了一定的进展.除了已被应用和开发的商用和开源求解器外，也建立了不同场景下的THMC耦合模化模型，但是在THMC全耦合模化分析方面的研究工作还比较有限.限于实验和场地测试条件，现有的耦合模型只能依靠有限数据进行修正.Xia等[59]指出，如果缺少建造预期缝网系统的地质-技术支撑，即使是已被数值验证了的EGS取热模型，也很难确定其可投入应用的可行性.综合分析，建立EGS多场耦合模型仍旧是一项长期复杂和挑战性的工作.

3 多场耦合快速求解方法

当计算域裂缝数量较多、网格划分较密时，多场耦合微分方程离散形成的矩阵维数较大，求解不仅消耗大量计算资源，且比较耗时.因此，研究人员发展出多种快速求解方法，比如并行计算、多重网格、MsFV、POD等，本节主要介绍MsFV和POD.

3.1 多尺度有限容积法

除了常规的有限体积法(FVM)、有限单元法(FEM)、边界元法(BEM)等数值解法外，为了提高数值计算效率，Jenny等[74]在FVM基础上发展了一种新的算法，即MsFV.该方法需要一套粗网格和一套细网格，或多层嵌套粗、细网格系统，如图4所示.MsFV通过网格系统，将原来细尺度上的物理问题转换到局部粗尺度上求解，从而降低计算量.在求解中需要构建两组基函数：一组基函数用于构建多尺度传递系数，另一组基函数用于计算细尺度上的速度通量.Jenny通过理论分析和数值算例验证了该方法的计算精度和计算效率，并将其应用到了地下多孔介质单相与多相渗流问题上.虽然MsFV在一定程度上可提高计算速率，但是在高各向异性或大网格纵横比的算例中，MsFV计算精度大幅下降.

图4 多尺度网格[13]

众多学者对MsFV进行了改进，在计算速度方面，Lunati等[75]通过计算局部基函数、计算粗尺度压力解、重构细尺度守恒通量和求解输运方程四个步骤，对MsFV进行了改进，新改进的方法在计算二维各向异性不可压缩流动问题时，可节省1.4倍～2.2倍时间；在计算精度方面，Lunati等[76]指出了精度下降是由于局部问题的边界条件不精确，导致在大网格比的情况下误差被逐步放大，提出了阻尼剪切边界条件用以缓解该情况.Hajibeygi等[77，78]针对多孔介质单相和多相流动问题，提出了迭代多尺度有限体积法(i-MSFV)，并应用线松弛技术(Line-Relaxation，LR)来消除粗网格边界上引入的误差余量，通过两相流算例表明，每10个时间步长更新一次边界函数，计算结果与细网格上的参考结果一致，该方法保证了计算精度，同时提高了计算效率.在稳定性方面，Lunati[79]构造了一个迭代多尺度算法，该算法核心是在多尺度求解中，在方程右端增加一个修正项，以弥补局部边界通量无法传递带来的误差，采用预条件广义最小余量法加快收敛，与LR技术相比，该方法无条件稳定.

目前，多尺度类方法已在流动、换热问题得到应用.Praditia等[80]基于传统MsFV方法，构建了描述裂隙岩体单相流动的热-流耦合模型.宇波等[34]采用两步迭代多尺度求解策略，使得多尺度求解误差余量先后在粗、细尺度上被耗散，通过裂缝型储层开采模拟显示，多尺度解与细尺度解相对偏差为10-3～10-4，同时计算效率提高了4倍～5倍.HosseiniMehr等[81]将MsFV方法应用到三维地热开发过程，并提出了代数动态多级(ADM)求解技术.此外，MφYner等[82]针对多孔介质渗流问题，提出了多尺度约束光滑基(MsRSB)方法，该方法适用于结构网格和非结构网格.同其他多尺度方法，MsRSB依赖底层细网格的粗尺度划分和一组局部延拓算子(多尺度基函数)，其计算采用雅克比迭代.该方法主要有三个优点：(1)粗网格和细网格都具有一般的多面体几何和非结构化拓扑结构，特别是由局部网格重建引入的非结构化单元连接，可以很容易地构造多尺度分区和良好的延拓算子；(2)与现有的多尺度方法相比，该方法具有较高的精度和鲁棒性，并且局部基函数更新时不需要高昂代价；(3)由于MsRSB方法是建立在以单元为中心的FVM之上，因此适用于任何压力分离求解的流动问题.Shah等[83]将MsRSB方法应用到裂缝流动问题，该方法在计算局部基函数时，需要对裂缝进行网格划分.Khayrat等[84]将MsRSB方法应用到了任意非结构孔隙网络流动问题.

虽然传统的MsFV方法在一定程度可加快多孔介质渗流问题计算效率，但是其应对复杂的储层，如高各向异性的储层介质，或粗网格边界存在源项问题时，求解精度会下降.其它诸如i-MSFV、MsRSB等方法，均在尝试拓宽MsFV方法的应用范围，并在使计算保持高精度的基础上，提高算法的效率.目前，多尺度类方法主要在流动、换热领域应用，在力学领域和THMC多场耦合领域仍需要研究.

3.2 POD模型降阶方法

POD是一种数据处理手段，POD与Galerkin投影法相结合可以建立原始高维物理问题的降阶模型，即离散方程的求解自由度大幅降低，从而提高数值计算的速度，由于其取样需要时间，因此适用于大规模多算例的场合，如可靠性评价、方案优化、结构反演等，可通过较少算例直接计算结合POD降阶模型实现其他多数算例的计算.一般问题的POD-Galerkin低价模型的建模过程如下：(1)采用数值模拟技术对原始物理问题进行求解，获取样本数据；(2)通过POD快照法(Snapshot)获取基函数[85]；(3)将原始方程向基函数空间投影，获得原始物理问题低阶模型；(4)低阶模型精度和计算效率验证，如图5所示.

图5 POD低阶模型建模过程[86]

POD与其它求解方法相结合，可以建立计算效率高的低阶模型.Chaturantabut等[87]采用离散经验插值法(DEIM)与POD相结合，对一个二维多孔介质非线性混相不互溶的黏性驱替问题，建立了有限差分离散系统的低阶模型.计算表明，黏性驱替过程的全阶系统求解需要自由度为15000维，而低阶模型只需要40维，计算时间减少了原来的1000多倍.Ghasemi等[88]通过引入辅助变量，将饱和方程转化为二次双线性形式，基于POD进行模型降阶时，无需使用DEIM近似计算非线性项，但却保留了原始模型的特征性，对两相非混相流动计算表明，该方法降低了计算的复杂度.Wang等[89]将POD与Galerkin投影相结合，对不可压缩多孔介质单相流动进行降阶处理，计算结果显示，相对偏差为1.0×10-4%～2.3×10-1%，加速效果为880倍～98 454倍.但是，Wang等[90，91]发现，对于双重孔隙介质内部的气体流动问题，应用POD降阶时，会出现质量不守恒问题.为此分别建立基质和裂缝方程，将两套方程相加后再投影到基质基函数或裂缝基函数.计算结果显示，基质压力计算最大偏差为3.87%，裂缝压力计算最大偏差为2.84%，计算效率可提高720倍.此外，为减少POD构造低阶模型时的计算量，Ushijima等[92]将耗时的系统矩阵投影分解成一组简单的标量矩阵乘法，与标准投影法相比，计算量极大减少.

目前，POD已应用到多数物理问题的快速求解方面.Siade等[93]在研究承压含水层问题时，通过一个简单的指数函数获取样本，计算显示，该方法运行速度约为全阶模型的1000倍.Rizzo等[94]在研究多孔介质被动输运问题时，采用POD减少对流弥散方程的计算量，通过将快照分裂技术与不同时间尺度上的建模策略相耦合，并将全阶模型的解与其降阶模型的解相结合，实现较快地捕获系统的显著特征，该方法不仅能显著降低计算成本，而且可精确地模拟输运过程.Pau等[95]针对热水驱提高原油采收率问题，建立了裂缝型储层多相流高分辨率快速预测降阶模型，该模型可以看作是一种基于回归的降阶技术，其特别之处是通过POD获得奇异向量来解决最小二乘的最小化问题，直接将低分辨率解映射到高分辨率解.计算结果显示，即使循环105次后，POD低阶模型还能准确地再现采油速度和含油饱和度的非均质结构，计算成本至少降低了85%.Li等[86]在研究裂缝型岩体热-流-力三场耦合过程时，采用一种新的POD低阶模型建模技术，将流动、换热和力学方程离散后，再投影到POD基函数空间，计算结果显示，相比于全阶的求解技术，该方法得到的离散矩阵条件数、求解自由度极大降低，计算效率可加速到13.2倍.

综上，POD将系统矩阵从全阶模型空间投影到一个子空间上，子空间的范围跨越全阶模型空间，但维数比全阶模型空间小得多，并且运行速度相当快，其模型误差也可以接受，由于其取样需要时间，因此较适用于大规模多算例的场合，如可靠性评价、方案优化、结构反演等，可通过较少算例直接计算结合POD降阶模型实现其他多数算例的计算.在裂隙型储层多场耦合方面虽然POD降阶模型已经在两场、三场的物理问题上得到证明和应用，但在THMC多物理场耦合领域，该技术还需要进一步研究与发展.

4 结 论

近年来深部地热能的开发利用被广泛关注，国内外学者针对EGS的研究快速发展.查明EGS取热过程中涉及到的流场、温度场、应力场、化学场有助于高效、稳定的热量提取.但是现有对于EGS热流力化各场耦合的研究多集中在两场、三场方面，对于THMC耦合作用下取热过程中的热质传输规律、储层时空演变、取热性能研究并不全面.本文总结了国内外关于裂缝型储层表征、EGS在多场耦合作用下的数值软件和耦合模拟、以及多场耦合快速求解方法的研究现状，分析探讨了EGS在多场耦合作用下研究的重难点，并得出相关结论，以期对我国EGS的研究发展提供参考，具体如下：

(1)EGS储层中存在着复杂的缝网系统，是典型的裂缝型储层，天然裂缝或压后形成的水力裂缝是EGS增产关键，准确刻画和描述裂缝系统是精细产量与产能预测的核心.尽管对裂缝型储层的表征在油、气藏等工程领域已有多年的发展，但应用于增强型地热系统时，需要在其基础上叠加复杂的THMC多场耦合模型，目前仍待深入发展与完善.

(2)截止目前，EGS已历经40多年研发工作，建立了可以描述EGS取热过程不同场景下的TH、THC、THM和THMC耦合模型，一些可用来研究EGS取热过程的商用和开源求解器也被陆续开发和应用.但是，对于THMC全耦合模型的研究工作还比较有限，相关的求解器也具有一定的局限性，缺乏微观、孔观到宏观的多尺度耦合模拟功能.因此要建立一个匹配实际的EGS多场耦合模型和研发一套完善的EGS多场耦合模拟求解器，仍是一项长期复杂和挑战性的工作.

(3)MsFV在一定程度上可加快多孔介质渗流问题的计算效率，但应对复杂的储层，如高各向异性的储层介质，或粗网格边界存在源项问题时，求解精度会下降.一些尝试拓宽MsFV应用范围的方法，如i-MSFV、MsRSB等，不仅保持了计算的高精度，同时还提高了算法的计算效率.目前该方法已应用到流动、换热领域，在力学、多场耦合领域仍在研究.

(4)POD可实现将系统矩阵从全阶模型空间投影到一个子空间上，子空间的范围跨越全阶模型空间.相比于全阶模型，降阶后的离散矩阵条件数、自由度数都极大降低，但计算效率却可以提高数倍，其模型误差也在可接受范围内.POD降阶模型已经在很多物理问题得到应用，但是在THMC多物理场耦合领域的适用性还需要进一步研究与发展.