基于温度场-渗流场-应力场耦合的双重介质流固热耦合分析

基于温度场-渗流场-应力场耦合的双重介质流固热耦合分析

1裂隙岩体耦合对岩体（尤其是断裂岩体）的温度场耦合研究具有重要的理论意义，并具有广泛的实践背景。高温干热岩体的开发、核废料贮存库的安全评估、环境工程中污染源的迁移及深部开采等都迫切需要研究岩体的温度场–渗流场–应力场三场耦合作用。从20世纪80年代以来,对该问题的研究取得了重要进展。J.Noorisha等提出了饱和裂隙岩体的THM耦合基本方程式,N.Barton和S.C.Bandis针对工程岩体对三场耦合进行初步探讨;L.Jing等给出了相对系统的三场耦合模型,但对模型简化和实用方面的研究还不够深入。近十几年来,提出了许多耦合模型,如单纯的热传导模型,水流动模型、水热耦合模型、固流热耦合模型;考虑岩体介质性态,有均质模型、裂隙介质模型等。相比之下较全面的模型是T.Kohl等提出的固流热耦合数学模型,T.A.Sprecker等提出并完善了裂隙岩体的温度场、渗流场与应力场耦合的有限元解法。当前,随着能源储库、核废料贮存库及干热岩体地热开发等大型岩土工程的出现,三场耦合研究成果已逐步服务于工程,用于评价工程的长期稳定性、安全性[9~12]。国内学者李宁等建立了一些耦合模型,但大多是基于多孔连续介质或离散介质假设。本文采用双重介质(拟连续介质+裂隙介质)提出研究裂隙岩体温度场–渗流场–应力场耦合数学模型,开发了三场耦合数值模拟软件,通过与渗流场–应力场耦合,应力场–温度场耦合进行数值比较可知,进行双重介质温度场–渗流场–应力场耦合分析是十分必要的。2裂隙岩体膜渗流规律模型裂隙岩体是岩土工程中广泛遇到的一类复杂岩体,是由随机分布的裂隙和被裂隙切割的岩块组成的不连续介质。裂隙岩体温度场–渗流场–应力场耦合作用是一个相对复杂的研究问题,主要表现在裂隙岩体赋存地质环境的各个组成部分,即渗流场、应力场与温度场都随时间、空间而发生变化;与此同时,以上各个组成部分之间的耦合作用处于一种复杂动态变化过程中。传统的离散介质假设(忽略岩块孔隙系统的透水性)和连续介质假设(将裂隙中的渗流平均到岩体的渗流中去)往往很难真实地反映裂隙岩体的温度场—渗流场–应力场耦合作用,而双重介质(拟连续介质+裂隙介质)三场耦合模型在一定程度上能真实地模拟裂隙岩体渗流、传热过程,对大型裂隙(断层、节理、人工裂缝等)可以查明,按裂隙介质对待;对被这些大中型裂隙切割而成的含众多低序次的小裂隙介质则不可能查明,认为是不确定的,采用拟连续介质模型模拟。然后根据两类介质接触处水头、温度、位移相等来建立裂隙岩体的三场耦合模型。本文建立的三场耦合模型基于如下假设:(1)岩体是由拟连续介质、裂隙介质组成的结构体,拟连续介质为均质同性的弹塑性体,裂缝介质服从节理单元模型。(2)渗流规律服从非线性达西定律。(3)不考虑水的汽化,岩体被单相水饱和。(4)热质在固相、液相介质中传递方式以传导和强迫对流为主,各介质的比热容及热传导系数不随温度和压力而变化。(5)介质的变形以孔隙变形为主。2.1裂隙介质拟连续介质及kfi式中:e为体积变形;βT,βP分别为流体热膨胀系数和压缩系数。对于拟连续介质:对于裂隙介质:式(3)~(7)中,kij,kfi分别拟连续介质和裂隙介质渗透系数;Θe为有效体积应力;b为裂隙宽度;v(T)为运动黏滞系数;φ为孔隙率;qi为比流量;p为水压;si为切向坐标。2.2耦合的能量平衡方程根据热力学第一定律:外界传入系统的能量与内部热源产生的能量之和等于物质内能的增量与对外做功之和,建立分别考虑岩体(固相)与温度场、液相与温度场耦合的能量平衡方程,该方程考虑了液体热对流(强迫对流为主)和固相热传导而导致的能量迁移。对于拟连续介质,有裂缝中主要介质为水,能量的迁移是通过裂隙水流的热传导、热对流(强迫对流为主)及岩体裂隙表面和裂隙水流之间的热传导的双过程实现的,即式(8),(9)中:λr,λw分别为岩体和水的热传导系数;ρr,ρw分别为岩体和水的密度;Cpr,Cpw分别为岩体和水的热容系数;Tr,Tw分别为岩体和水的温度;Trb为岩体裂缝边缘的温度。2.3温度场–渗流场—静力平衡控制方程对于拟连续介质,静力平衡控制方程为对于裂隙介质,静力平衡控制方程为式(10)~(12)中:λ,µ均为拉梅常数;uj,ji,ui,jj均为岩体位移;F为体力;β为岩体热膨胀系数;σn′,σs′分别为裂缝法向和切向有效应力;kn,ks分别为裂缝法向与切向刚度;εn,εs分别为裂缝法向与切向变形;α为等效孔隙压系数。式(10)~(12)这3组控制方程,再辅以初始和边界条件,就构成了基于双重介质假设的裂隙岩体温度场–渗流场–应力场耦合数学模型。这一模型具有如下特点:流体渗流方程中考虑了岩体的变形和温度对渗透系数与水份迁移的耦合影响。能量守恒控制方程(式(8))中左边第1项为温度变化引起的能量变化;第2项为热传导引起的能量变化,第3项为热对流引起的能量变化,第4,5项分别表示热–液耦合和热–固耦合项。固体变形方程式(10)中(αp,)i,βT,i描述渗流和温度对变形的耦合作用,如果是等温过程,能量方程不必要,在渗流方程和应力本构方程中略去含温度T的各项,则方程组退化为流–固耦合渗流的控制方程。如果只考虑热–固耦合,则在静力平衡方程中没有(αp,)i项。在能量守恒控制方程式(8)中去掉左边第3项热对流引起的能量变化,第4项热–液耦合项即可。图1给出了温度场–渗流场–应力场耦合模型(以拟连续介质为例忽略变形生热)。3双介质三场耦合模型的数值分析3.1高温岩体计算模型双重介质渗流场–应力场–温度场耦合非常复杂,基本求解策略是将不同介质的固体变形、能量方程和渗流看作独立的子系统,根据两类介质接触处水头、温度、位移相等,温度、应力、水头在同一介质内及不同介质间随时间传递,而将各子系统耦合求解。对时间采用变时步差分方法,对空间采用有限元方法离散求解,对拟连续介质采用空间等参单元,裂缝采用八节点等参节理单元,开发双重介质三场耦合模型的三维有限元程序。下面用本文开发的三维有限元计算程序对高温干热岩体地热开发系统进行多场耦合数值模拟分析,探讨多物理场之间的耦合效应。图2给出了地下2200m的高温岩体的计算模型(人工控制致裂形成众多近似平行的人工裂缝构成了水渗流的裂隙网络通道,模型底部为固定约束,四周取不透热边界)及有限元网格图,(x,y)坐标系下注入井坐标为(37.5m,37.5m);生产井坐标为(162.5m,162.5m),计算参数见表1。本文具体模拟了3种工况:工况1:温度场–渗流场–应力场耦合;工况2:渗流场–应力场耦合(不考虑热应力、对流项对温度的影响和温度变化引起水压力的变化);工况3:温度场–应力场耦合(不考虑流体渗流)。选取裂隙面、过注入井和生产井的纵剖对角面及裂隙面上靠近注入井的5428#节点、靠近生产井的5598#节点,研究3a内岩体应力、温度、水压分布和裂隙单元的张开度变化及介质渗透系数的变化。3.2围岩应力分析三场耦合作用下,水压和温度变化梯度产生附加节点载荷,注入的低温水使得注入井附近围岩温度下降,围岩发生收缩变形,引起岩体应力上升,为协调围岩收缩变形,裂隙单元应力降低;水压变化梯度产生附加节点载荷,导致研究区域的有效应力低于原岩应力。三场耦合作用下,岩体应力总体上要小于原岩应力,并随渗流方向应力逐步上升。图3,4分别给出了三场耦合下,纵剖对角面和裂隙面上σz分布图,拟连续介质内应力为40.5~52.5MPa,而注入井附近裂隙单元应力只有23.0MPa左右。(2)裂隙水土流失区的黄土壤压射在三场耦合作用下,水在基岩(拟连续介质)中渗流速度远不如裂隙介质渗流快,主干裂隙系统的强导水性使水压沿此通道迅速上升,在非稳定渗流阶段,裂隙水压要高于拟连续介质内的水压(见图5),主干裂隙显著地控制着研究区域的水压分布。由于固定井口压力的影响,使得井口附近区域裂隙水压梯度大,中间区域水压梯度少,表现为图6中井口附近区域水压等值线相对密集,中间区域水压等值线相对稀疏。(3)储层流变质岩渗流场的耦合图7,8分别给出了三场耦合下纵剖对角面和裂隙面上温度分布图。流体渗流过程也是能量的传递过程,能量在渗流过程中通过热对流和热传导,从岩体中传到了流体中,能量被注入水提取,这体现了能量方程中的φCpwρwkiip,iiTw,i,λrTr,ii(拟连续岩石介质)和(ρwCpwTwkfip,i)i,2λr(Trb-Tw)/b(裂隙介质)项,在拟连续岩石介质中由于渗透系数远低于裂隙介质的渗透系数,孔隙水流的热对流速度远不如裂隙水流的热对流速度,使裂隙介质的温度比它毗邻的拟连续岩石介质的温度要低(见图7)。这体现了渗流场对温度场的耦合。由于主干裂隙的存在使水头与温度分布呈现出各向异性,因而这是单纯连续介质耦合模型无法描述的。4双介质三场耦合效应分析(1)渗流场及其应力场耦合分析图9~12为不同耦合工况下围岩应力分布情况。在工况3中,温度场只对注入井附近围岩应力产生显著影响,该区域温度降低梯度大,产生较大的热应力,注入井附近80m区域内围岩应力明显上升,高于原岩应力0%~12%,裂隙介质为协调岩体的收缩变形,出现了应力降低区,在注入井附近裂隙单元应力仅为45MPa(见图9,10);考虑工况2时,渗流场对围岩应力的作用是全局性的,渗流场对应力场的耦合效应通过修正太沙基有效应力来实现,岩体有效应力随渗流方向逐步上升(见图11,12)。在三场耦合分析下,存在复杂的耦合路径如渗流场和温度场之间存在热对流耦合作用(见图1)。因此,三场耦合分析下围岩应力不是图9,11的简单叠加。在注入井附近区域内,温度场、渗流场及应力场的三场耦合作用是显著的,而靠近生产井区域,渗流场–应力场耦合作用控制着此区域内应力分布。由此可见,一定空间、时间域内,往往存在起主导控制作用的双场耦合系统。对于裂隙介质来说渗流场、温度场对应力场的耦合作用可归结为裂隙水压、温度梯度产生附加节点载荷导致裂隙张开度的变化。从图13中可以看出:裂缝张开度随时间的延长总体呈增加趋势,在t=1a内,裂缝张开度急骤增加,注入井附近5428#裂隙单元三场耦合分析,裂隙张开度变化最大∆b=0.0356cm,其次是渗流场–应力场耦合分析,∆b=0.0288cm,t=1a后裂隙张开度趋于稳定;而应力场–温度场耦合分析时,t=5100h该裂隙单元张开度变化仅为∆b=0.0085cm,在5100h后该裂隙单元张开度趋于稳定。这显示了渗流场对裂隙应力场的强耦合效应。三场耦合分析下,生产井附近5598#裂隙单元水压、温度梯度产生的附加载荷相对较小,在t=1a内,该裂隙单元张开度变化量∆b=0.0192cm,1a后,随着围岩应力的调整,该区域内裂隙单元的应力略有增加致使单元张开度稍有下降,t=3a时,∆b=0.0185cm。(2)裂隙单元分布应力场对渗流场的耦合作用主要体现在拟连续岩石介质的渗透系数是有效体积应力的负指数函数;裂隙介质的渗透系数是裂隙变形的函数。取z=60m的水平岩体面,在耦合工况1下,该面上的渗透系数分布如图14所示,渗透系数从低应力区的3.4×10-4m/s到高应力区的0.7×10-4m/s变化;应力场、温度场变化导致裂隙张开度和流体黏性的改变,在耦合工况1,2下裂隙渗透系数总体上呈时间相依的增加趋势。在t=10079h时对5428#裂隙单元进行三场耦合分析,裂隙渗透系数为0.87m/s,渗流场–应力场耦合分析裂隙渗透系数仅为0.65m/s,为三场耦合分析所得裂隙渗透系数的74.7%(见图15)。这是由于未考虑温度场耦合作用所致,可见裂隙渗流场对温度场的变化较为敏感;这显示出了温度场对裂隙渗流场的较强耦合作用。5流固热三种耦合作用机理(1)双重介质温度场–渗流场–应力场全耦合数学模型包含了众多耦合项的作用,并定义了不同介质的本构关系及耦合变量;该模型考虑了不同介质间的水量和热量交换,体现了主干裂隙的优先渗流和强热对流,具有较好的拟真性和操作性。(2)双重介质三场耦合路径非常复杂,本文将固体变形、能量方程和渗流看作独立的子系统,根据两类介质接触处水头、温度、节点位移相等的原则,采用有限元和有限差分法相结合,将具有耦合效应的参数值在各子系统相互传递耦合求解。(3)建立了高温干热岩体地热