黄金科学技术,2020,28(6):846-858doi:10.11872/j.issn.1005-2518.2020.06.122
矿产勘查与资源评价
2.昆明理工大学城市学院，云南昆明650051
2.CityCollege，KunmingUniversityofScienceandTechnology，Kunming650051，Yunnan，China
收稿日期:2020-07-09修回日期:2020-08-03网络出版日期:2021-01-29
Received:2020-07-09Revised:2020-08-03Online:2021-01-29
作者简介Aboutauthors
由于断层宽度远小于其延伸方向的尺寸，造成数值模型建模困难、计算效率低等问题。把断层概化为无几何厚度的空间曲面，并引用裂隙渗流理论中的立方定律对其进行渗流计算，可以有效降低建模难度。本研究旨在验证断层概化方法的可行性和合理性，并解决断层中的裂隙流与基岩中的达西流之间的流—热耦合问题。采用公式推导得出断层与基岩之间的流—热耦合控制方程，并使用数值模型的计算结果进行了验证分析。结果表明：当不考虑断层内部结构影响时，利用立方定律计算断层渗流的方法是可行的；耦合控制方程对基岩与断层之间的流—热耦合计算合理有效；随温度而变化的流体粘滞性对数值模型计算结果影响显著。
关键词：立方定律;断层;流—热耦合;数值计算;热液型矿床;裂隙渗流
Keywords：cubiclaw;fault;flow-thermalcoupling;numericalcalculation;hydrothermaldeposit;fissureseepage
本文引用格式
陈刚,马玲,龚红胜.基于立方定律的断层流—热耦合数值计算方法[J].黄金科学技术,2020,28(6):846-858doi:10.11872/j.issn.1005-2518.2020.06.122
CHENGang,MALing,GONGHongsheng.NumericalCalculationMethodofFaultFlow-ThermalCouplingBasedonCubicLaw[J].GoldScienceandTechnology,2020,28(6):846-858doi:10.11872/j.issn.1005-2518.2020.06.122
热液型矿床的物理过程非常复杂，构造往往是其主要的热液渗流通道，成矿动力学数值模拟通常会面临以下问题：（1）岩体中裂隙、断层的渗透性远大于岩石基质，但裂隙、断层数量众多，二者的宽度远小于其延伸方向尺寸，造成模型建立困难；（2）从渗流角度而言，岩石基质中的渗流符合达西定律，而裂隙中的渗流多数情况下符合立方定律，二者相互作用、相互影响形成复杂的流—热耦合问题；（3）流体的动力粘滞性具有随温度和压力变化的特点，粘滞性的变化对模拟成矿结果的影响问题。
本研究的目的是验证断层概化为空间曲面方法的可行性和合理性，并解决断层中的裂隙流与基岩中的达西流之间的流—热耦合问题。初步探讨流体粘滞性随温度变化这一性质对数值模型计算结果的影响。解决上述问题，可以降低建模难度，在数值模型中增加大量对热液渗流有影响的断层或裂隙成为可能，使热液渗流数值计算更为精细、准确，并扩展数值方法在成矿动力学数值模拟中的应用范围。
式中：q为流量（m2/s）；g为重力加速度（m/s2）；a为缝隙的隙宽（m）；J为水力梯度；ν为流体的运动粘滞系数（m2/s）。
从式（1）可以看出，通过缝隙的流量大小q与缝隙隙宽a的三次方成正比，即著名的立方定律。为更好地表达裂隙本身的渗透特性，使用渗透率概念κ（m2），其表达式为
式中：f为裂隙表面的粗糙度；κ为裂隙渗透率。
受限于篇幅，本文对于裂隙渗流主要考虑开度和粗糙性2个方面。
三维条件下裂隙内流量方程表示为
裂隙内连续性方程（渗流平衡方程）表示为
图1裂隙流—达西流耦合概念模型
Fig.1Conceptualmodeloffractureflow-Darcyflowcoupling
代入到裂隙的连续性方程式（5）中得出：
由于裂隙隙宽很小，且裂隙内部为连续介质场，所以对一个微小单元格来说，可认为其内部压力一致，表示为
Pup=Pdown=PF（8）
上、下裂隙面压力相等，这一特性造成裂隙面两侧的基岩内的压力在靠近同一裂隙时，有压力逐渐相等的趋势。当裂隙流平衡方程为式（7）时，可以实现裂隙流与达西流的耦合。因为是2种物理场的耦合，数值建模时使用不同维度的网格，裂隙网格相较于基质网格而言要低一维度。当裂隙内部压力和岩石基质与裂隙接触面上的压力一致时，可以完成裂隙与基岩、裂隙与裂隙之间的流量交换计算，计算结果可以反映出裂隙对基岩渗流的影响，以及裂隙交叉后的渗流。
裂隙中的能量守恒方程为
由于裂隙的隙宽很小，所以内部温度是相同的，则有：
Tup=Tdown=TF（12）
式（11）代入式（9）后就是裂隙内能量守恒方程。利用式（8）和式（12）作为基质和裂隙的压力、温度边界条件，可以保证裂隙与岩石基质之间的质量和能量连续性；联立式（6）和式（11）方程作用裂隙流量、能量的输入项，可以保证裂隙内的质量、流量守恒，同时保证了达西流和裂隙流2种物理场的耦合。
流—热耦合是一个复杂的多物理场强耦合过程，本文中的数值模型采用COMSOL软件进行建模和计算，模型中分别使用了达西渗流模块、裂隙流模块和固体传热模块，用于模拟不同的物理过程。为保证模型中裂隙与基质之间的流量、能量的耦合，渗流计算模块中在裂隙上面和下面分别施加式（6）所代表的质量源边界条件，以实现式（7）方程；传热模块中在裂隙上面和下面分别施加式（11）所代表的热能量边界条件。
（2）裂隙宽度为2b，裂隙内水的流速恒定；
（3）裂隙及上下两层岩石初始温度为T0，注入水的温度为Tin；
（4）忽略岩石基质中平行于裂隙方向的热传导；
（5）裂隙周围岩石基质为无限厚。
图3裂隙内部温度变化曲线
Fig.3Temperaturechangecurveinsidethefracture
图4地质模型结构及测点空间位置
（a）模型结构（b）模型中测点分布位置
Fig.4Geologicalmodelstructureandspatiallocationofmeasuringpoints
表1地层参数列表
Table1Listofformationparameters
图5模型施加的温度—压力曲线
Fig.5Temperatureandpressurecurveappliedbythemodel
图6数值模型断层面通量曲线
Fig.6Fluxcurvesoffaultplaneofnumericalmodel
图7第二层基岩中多点温度、压力曲线
Fig.7Multi-pointtemperatureandpressurecurvesinthesecondlayerofbedrock
Fig.8Slicesoftemperaturedistributionindifferenttimemodels（temperatureunit：℃）
将模型中的断层抽象为三维曲面裂隙，其隙宽参数为2mm，考虑裂隙粗糙度因素，其水力隙宽换算后约为1.58mm。这一宽度数值对自然界中的断层而言也不算大，但其对热能、渗流的传导影响却造成较大范围内基岩流—热产生巨大的变化，从侧面反映出断层流—热的影响巨大，尤其对于热液成矿的矿床，断层影响是不能忽略的重要因素。
图9断层通量、热传导通量的比值曲线
Fig.9Ratiocurvesoffluxandheatconductionflux
成矿动力学模拟是了解成矿过程的有效方法之一，成矿过程由于其复杂性涉及多个学科领域。本文以裂隙渗流理论为基础，探讨了断层在成矿动力学模拟过程中的影响，综合分析得出如下结论：
（1）裂隙和岩石基岩中的渗流有着不同的渗流理论，本文得出的裂隙流和达西流耦合控制方程，保证了数值模拟过程中需要保证数值模型计算域内渗流场压力、速度、质量和能量的连续性；通过对数值模型的分析，认为耦合控制方程可以合理、准确地实现裂隙流和达西流的耦合。
（2）由于裂隙和断层的空间形态具有宽度远小于延伸尺寸的特性，本文将断层简化为曲面，基于热—流耦合控制方程，使用数值方法计算了涵盖达西渗流、裂隙流和热传导的成矿动力学过程的数值模型；由于裂隙流内渗流速度快、热传导效率高，对模型最终的温度场和能量传导影响十分显著。
（3）经过数值模型的验证认为，在不考虑断层内部结构时，成矿动力学模拟中将断层概化为无几何厚度的空间曲面，使用裂隙渗流理论进行计算是可行的。
赵义来，刘明亮.复杂形态岩体接触带成矿耦合动力学三维数值模拟：以安庆铜矿为例
ZhaoYilai，LiuMingliang.3D-numericalmodelingofcoupledgeodynamicprocessesandmineralizationatthecontactzonesofcomplexplutons：ExampleformtheAnqingdeposit，AnhuiProvince，China
朱静，陈建平.基于FLAC3D的成矿过程模拟研究现状
ZhuJing，ChenJiangping.ResearchstatusofFLAC3D-basedmineralizationprocesssimulation
刘向冲.构造—流体耦合有限元模拟：以石英脉型钨矿为例
LiuXiangchong.Finite-elementsimulationsofstructure-fluidcoupling：Acasestudyinvein-typetungstendeposits
戴文强，李晓晖，袁峰，等.安庆铜矿床典型矽卡岩矿物形成过程数值模拟
DaiWenqiang，LiXiaohui，YuanFeng，etal.NumericalsimulationofformationprocessoftypicalskarnmineralsinAnqingcopperdeposit
赵崇斌，HobbsbeBE，OrdA.用计算地球科学研究方法探讨地质现象的动力学机制——以断层中等距成矿分布为例
ZhaoChongbin，HobbsBE，OrdA.Investigatingdynamicmechanismsofgeologicalphenomenausingmethodologyofcomputationalgeosciences：Anexampleofequal-distantmineralizationinafault
池国翔，薛春纪.成矿流体动力学的原理、研究方法及应用
ChiGuoxiang，XueChunji.Principles，methodsandapplicationsofhydrodynamicstudiesofmineralization
李瑞红，安平，恽孟河，等.焦家断裂带三维结构模型及其数值模拟：以新城金矿床控矿构造为例
LiRuihong，AnPing，YunMenghe，etal.3Dore-controllingstructuralmodelandnumericalsimulationofJiaojiafaultzone：AcasestudyoftheXinchenggolddeposit
杨永春，刘家军，王学银，等.甘肃滴水山金矿不同构造—岩相带岩石地球化学特征及构造控矿机理探讨
YangYongchun，LiuJiajun，WangXueyin，etal.Geochemicalcharacteristicsandstructuralore-controlmechanismaboutdifferentstructural-lithofacieszonesoftheDishuishangolddepositinGansuProvince
赵少攀，徐书奎，徐宗蛟，等.河南洛宁县龙门店银矿区断裂构造分期及含矿性评价
ZhaoShaopan，XuShukui，XuZongjiao，etal.Stagedivisionandore-bearingevaluationoffaulttectonicintheLongmendiansilverdeposit，Luoning，HenanProvince
JuanesR，SamperJ，MolineroJ.Ageneralandefficientformulationoffracturesandboundaryconditionsinthefiniteelementmethod
陈必光，宋二祥，程晓辉.二维裂隙岩体渗流传热的离散裂隙网络模型数值计算方法
ChenBiguang，SongErxiang，ChengXiaohui.Anumericalmethodfordiscreatefracturenetworkmodelforflowandheattransferintow-dimensinalfracturedrocks
LauwerierHA.Thetransportofheatinanoillayercausedbytheinjectionofhotfluid
PruessK，BodvassonGS.Thermaleffectsofreinjectioningeothermalreservoirswithmajorverticalfractures
ChengAHD，GhassemiA，DetournayE.Integralequationsolutionofheatextractionfromafractureinhotdryrock
赵坚.岩石裂隙中的水流—岩石热传导
ZhaoJian.Experientalstudyofflow-rockheattransferinrockfractures
赵阳升，王瑞凤，胡耀青，等.高温岩体地热开发的块裂介质固流热耦合三维数值模拟
ZhaoYangsheng，WangRuifeng，HuYaoqing.3DnumericalsimulationforcoupledTHMofrockmatrix-fracturedmediainheatextractioninHDR
张树光，李志建，徐义洪，等.裂隙岩体流—热耦合传热的三维数值模拟分析
ZhangShuguang，LiZhijian，XuYihong，etal.Three-dimensionalnumericalsimulationandanalysisoffluid-heatcouplingheat-transferinfracturedrockmass
唐志伟，米倡华，张学峰，等.增强型地热系统固流耦合数值模拟与分析
TangZhiwei，MiChanghua，ZhangXuefeng，etal.Numericalsimulationandanalysisofthecoupledforheat-fluid-solidinenhancedgeothermalsystems
曲占庆，张伟，郭天魁，等.基于COMSOL的储层参数与层理缝对地热产能影响规律研究
QuZhanqing，ZhangWei，GuoTiankui，etal.ResearchontheeffectofgeothermalreservoirandbeddingfracturesongeothermaldeliverabilitybasedonCOMSOL
张伟，孙江，曲占庆，等.高温地热开采热流固耦合模型及综合评价方法
ZhangWei，SunJiang，QuZhanqing，etal.Thermo-hydro-mechanicalcouplingmodelandcomprehensiveevaluationmethodofhightemperaturegeothermalextraction
DavidTS.Anisotropicpermeabilityoffracturedmedia
LouisC.Astudyofgroundwaterflowinjointedrockanditsinfluenceonthestabilityofrockmass
速宝玉，张文捷，盛金昌，等.渗流—化学溶解耦合作用下岩石单裂隙渗透特性研究
SuBaoyu，ZhangWenjie，ShengJinchang，etal.Studyofpermeabilityinsinglefractureundereffectsofcoupledfluidflowandchemicaldissolution
王媛，速宝玉.单裂隙面渗流特性及等效水力隙宽
WangYuan，SuBaoyu.Researchonthebehavioroffluidflowinasinglefractureantitsequivalenthydraulicaperture
刘明亮，周瑞超，赵崇斌.构造应力环境对浅成岩体成矿系统的制约：从安庆月山岩体冷却过程动力学计算模拟结果分析
LiuMingliang，ZhouRuichao，ZhaoChongbin.Constraintsoftectonicstressregimeonmineralizationsystemrelatedtothehypabyssalintrusion：ImplicationfromthecomputationalmodelingexperimentsonthegeodynamicsduringcoolingprocessoftheYuenshanintrusioninAnqingdistrict，China
式中：q为流量（m2/s）；g为重力加速度（m/s2）；a为缝隙的隙宽（m）；J为水力梯度；ν为流体的运动粘滞系数（m2/s）....42004...岩石裂隙中的渗流不同于孔隙介质的渗流，其具有特定的渗流规律，形成一门独立的学科，即岩石水力学［14］.裂隙广泛分布于各类地层中，并构成主要渗流通道；而裂隙周围的岩石基质的渗透性很小，但占有地层绝大部分体积和孔隙［15］.在成矿动力学模拟过程中，对热液的渗流模拟是其中非常重要的一个方面，渗流模拟如果忽略裂隙的渗流影响，其模拟结果可能会存在较大的偏差.裂隙中的流体与基岩之间的流量、热量交换计算十分重要，尤其是在流—热耦合计算过程中.热液渗流对于成矿结果的影响显著，这其中就涉及裂隙与基岩之间的流量、热量交换问题....
流—热耦合是一个复杂的多物理场强耦合过程，本文中的数值模型采用COMSOL软件进行建模和计算，模型中分别使用了达西渗流模块、裂隙流模块和固体传热模块，用于模拟不同的物理过程.为保证模型中裂隙与基质之间的流量、能量的耦合，渗流计算模块中在裂隙上面和下面分别施加式（6）所代表的质量源边界条件，以实现式（7）方程；传热模块中在裂隙上面和下面分别施加式（11）所代表的热能量边界条件....