岩体冻融损伤的温度场-渗流场-应力场耦合控制方程

岩体冻融损伤的温度场-渗流场-应力场耦合控制方程

1冻结岩体内部温度和渗流方程岩体的冷冻和解冻破坏包括复杂的温度场、渗流场和应力场的耦合。造成岩体工程冻害的主要原因是岩体裂(孔)隙中的水随温度变化发生相变,产生冻胀融缩效应。裂隙水结冰时,会产生约9%的体积膨胀,从而产生冻胀荷载;冰融化时,水会渗入新扩展的裂隙。反复冻融循环导致岩体损伤。冻结缘是岩体中相变最为活跃区域,也是温度场–渗流场–应力场(THM)耦合作用强烈的区域。数十年来,不少科研工作者对冻结岩体多场耦合问题进行了研究。R.D.Miller提出,冻结缘带位于冻结锋面和最暖冰透镜之间,并认为冻结缘处水分不断凝结成冰是造成冻胀的主要原因。J.P.G.Loch和B.D.Kay认为,冻结冰(水)压与渗透性及温度相关,可假定符合热动力学平衡条件,忽略岩石基质的变形,得出水–冰热动力学平衡公式:其中,式中:pi为冰压力,pw为水压力,ρi为冰的密度,ρw为水密度,π为渗透压力,L为相变潜热,T为温度,T0为冰点,Cw为水的热容,R为常数。在国内,杨更社和张全胜通过试验证明温度梯度是水分迁移的主要驱动力。温度梯度越大,则水分场越快达到重新分布状态;且孔隙率越高,冻结过程越长。张玉军考虑水冰潜热的影响,并对THM耦合模型及二维有限元程序做了改进,对岩土体冻融问题进行了分析。徐光苗等推导了冻结岩体的质量守恒方程、平衡方程及能量守恒方程,指出岩体应变会因冰体的冻胀产生膨胀效应,并研究了冻结温度下冰与岩石的膨胀耦合关系。国内外学者对岩体低温多场耦合问题的研究取得了一定成果,但以往的研究多将冻结岩体视为等效连续多孔介质,很少考虑裂隙的影响,对冻结率的研究也不多见。冻结过程会引起岩体渗透性的改变,是冻岩多场耦合的重要特征之一,但目前关于此问题研究也十分鲜见。本文在前人研究的基础上,立足冻结岩体的特殊性,首先,分析了冻岩THM耦合特征和相变物态方程,并基于能量守恒原理得出冻结率的表达式;然后,基于双重孔隙介质理论,推导了裂隙岩体低温THM耦合条件下的应力平衡方程、水–冰系统连续性方程及能量守恒方程;最后,通过一个含裂隙隧道低温多场耦合算例,研究了其温度场、应力场和孔隙压力等的分布规律。运用等效热膨胀系数法模拟夹冰(含水)裂隙的冻胀效应,并考虑了冻结过程对渗透系数的影响。以期为裂隙岩体低温多场耦合的研究提供一定的参考。2冻岩低温thm耦合-冰质体系相互作用机制冻结缘的概念是由R.D.Miller提出的,最初应用于冻土的研究。J.S.Walder和B.Hallet提出冻岩内冻结缘的水分迁移及透镜体增长模式。冻结缘的THM耦合是个十分复杂的过程,归纳起来,其相互作用机制可用图1所示关系来表示。如图1所示,冻岩低温THM耦合与常温耦合的重要区别在于水–冰相变的参与。当温度高于冰点时,低温THM耦合机制与常温相同。但当温度低于冰点并产生冰晶时,THM耦合的机制会发生变化。水结成冰引起相变潜热释放,会影响温度场,同时,体积膨胀会产生附加应力场。冰晶还会导致裂隙岩体渗透系数降低。常温下,地下水的存在形式为完全流体状态,而在寒区工程岩体中,当温度降至冰点以下时,自由水结冰导致多孔介质的渗流性降低。此外,相变还会引起岩体其他物理力学性质的变化。3可萃取过程和冻结率3.1冰乳化体积增量fsv岩体裂隙水在一定的温度和外界压力条件下,液、固相比例不再变化时即达平衡状态。此时冻结率到达稳定值,称此状态为相平衡状态,假定满足Clapeyron方程:式中:ΔH为焓变增量(J),ΔfusH为单位质量的冰融化焓增量(J·kg-1),ΔfusV为单位质量的冰融化体积增量(m3·kg-1),P为外界压力(Pa)。由式(2)积分可得式中:k为常数。式(3)即为外界压力与冰点的函数关系。设0℃(273.2K)和标准压力(1.01325×105Pa)下水和冰的密度分别为999.9kg/m3和916.8kg/m3,冰的融化焓ΔfusH为333.5kJ/kg,则单位质量冰的融化体积增量为将式(4)代入式(3),得标准压力为1.01325×105Pa时,水的冰点为273.15K,代入式(5)得k=273.16。运用Matlab绘制0.1~20.0MPa压力区间内冰点与压力的关系曲线,如图2所示,由图可知,冰点–压力关系曲线十分接近直线。3.2融冰系统及潜热系统概述设在封闭裂隙内水的冰点为T0,未冻结时水的初始体积为V0,初始水压为P0,降温时冰晶出现,当温度降至T时,冻结率为uTi,冻结压力为PT,降温过程中水–冰系统的温差为ΔT=T-T0,冻结过程如图3所示。假定:(1)冻胀力的增长是一个准静态过程;(2)相变释放的潜热全部被水–冰系统吸收,忽略传递至外界的潜热。由初始状态S0(P0,T0)降温至稳态S2(P2,T2),可等效为以下2个路径叠加(见图4):(1)S0→S1(路径1):由冰点T0等压降温至T0+dt,充分实现热扩散,最终完全冻结,至稳态S1;(2)S1→S2(路径2):由稳态S1绝热等温压缩至状态S2。根据热力学第一定律,能量增加值ΔU等于吸收的热量Q与外界做功W之和,即其中,式中:ΔV为体积增量;p(V)为以体积V为参变量的压力函数,其值随体积V的变化而变化。依据假定(2),可将路径2视为绝热等温压缩,此过程热力学能不变,所有的焓变等于外界做功,内能增量为w,全部用于增加相变潜热,导致部分冰融化,从而储存能量。设融化的冰体积为Vtr,称为融冰,储存潜热为LVtr,其中L为潜热系数。由内能的计算公式:可得从而可得状态S2的冻结率ui为式中:ΔP为压力增量。4充分反渗裂隙的影响以往关于冻岩THM耦合的研究多采用等效多孔介质的方法,不能充分反映裂隙的影响。本文尝试采用双重孔隙介质理论对冻结岩体低温THM耦合问题进行研究,根据质量守恒定律、能量守恒定律和静力平衡原理及相应的物性方程,推导裂隙岩体低温THM耦合方程。4.1裂隙介质的应变–应力增量关系式及柔度矩阵将裂隙岩体视为裂隙介质和岩块组成的双重孔隙介质系统。水(冰)主要赋存在裂隙介质中。局部坐标系(x′oy′)及整体坐标系(xoy)下的夹冰(含水)裂隙如图5所示。采用双重孔隙介质模型对冻结裂隙岩体进行研究。参照张玉军对裂隙介质的描述,考虑平面应力问题,设岩体内分布有n组裂隙,应变–应力增量关系式可表示为其中,式中:C为总的柔度矩阵;C0为孔隙介质(岩块)的柔度矩阵;Cfi为整体坐标系下第i组裂隙的柔度矩阵;αs为岩石骨架的线性热膨胀系数;dT为温度变化量;δij为Kronecker符号;dεij为应变增量张量;σij为应力张量;μ为泊松比;C′fi为局部坐标系下第i组裂隙的柔度矩阵;Li为坐标转换矩阵;Kni,Ksi,li,Ai分别为裂隙的法向刚度、切向刚度、连通率和间距。冰与水共同充填孔隙,保持良好的联通性,可近似认为pw=pi=p。增量形式应力的计算公式为式中:σijs为孔隙介质应力张量,κ为冻胀传压系数。4.2裂隙介质网络模型水分迁移的驱动力主要包括:冻胀压力引起的水势梯度、初始应力,分凝势及温度梯度。如图5所示,冻结压力势导致水分背离冻结缘迁移,而分凝势和温度梯度导致未冻水向冻结缘迁移。笔者认为,冻结条件下水的连续性方程可表示为对于一般的孔隙介质而言,水流矢量可表示如下:式中:为水流矢量;k为岩石介质的固有渗透系数张量;kw为水的相对渗透系数;μw为水的黏滞系数;∇Pw为水势梯度;klt为热流耦合系数张量,表征热梯度对水的流速影响;∇T为温度梯度。对于双重孔隙介质模型,裂隙的渗透系数远大于基质孔隙介质。假定单裂隙的渗流遵从立方定律:其中,式中:kf为裂隙的渗透系数;q为单宽流量,其方向由裂隙分布决定;b为张开度;Jf为沿裂隙方向的水力梯度;ν为水的动黏滞系数;γw为水的容重;ξ为与裂隙面粗糙度相关的系数,光滑平直时取ξ=1/12;ψ为冻结率–渗流耦合系数,因相变冻结而导致的裂隙渗透性降低。多组裂隙时,假定不同方向裂隙组在裂隙网络中互相连通且各方向渗流互不干扰,则裂隙介质岩体的渗流方程可表示为式中:Kf为渗透系数张量,m为岩体内裂隙组数,Si为第i组裂隙的平均密度,bi为第i组裂隙的平均隙宽,I为单位矢量,αi为第i组裂隙的平均方向余弦矢量。考虑外界应力及冻结过程的影响,采用渗透系数张量修正式为式中:Kf0为初始渗透系数张量;Df为裂隙分布的分维数,一般有0<Df≤2,Df=0时表示岩体内无裂隙;ψ(u)为冻结率–渗透率耦合系数,0<ψ(u)≤1时考虑其为冻结率的函数,当u=0时ψ(u)=1,当u→1时ψ(u)→0;σ为应力。笔者认为,关于引起裂隙水迁移的水力梯度,应主要考虑冻胀压力、分凝势及温度梯度的影响,此时的计算公式为式中:φT为分凝势,与冻结速度、介质特征和含水率等因素相关。4.3岩体热传导模型裂隙岩体的热传递方式主要包括热传导、相变潜热和裂隙水渗流迁移热。定义2种系统传热方式:静态传热和动态传热。前者是不引起宏观物质运动而产生的热迁移方式,后者是由宏观物质交换引起的热转移。冻结岩体热平衡方程可表示为式中:Cv为系统等效容积热容;Cw为水的容积比热;λ为岩体等效导热系数,冻岩和未冻岩取值不同。岩体的体积热容表达式为岩体的密度表达式为式中:θs,θw,θi分别为单位体积内岩石基质、水和冰的体积比;Cs,Cw,Ci分别为岩石基质、水和冰的体积热容;ρs为岩石基质。温度梯度引起的热传递符合傅里叶定律,即式中:λ为热传导系数(W·m-1·℃-1),qT为温度梯度引起的热流。由水渗流作用产生的热量变化为在冻融岩体内,水冻结放热,融化吸热,水–冰相变产生的潜热为5温度耦合单元法冻害问题严重影响寒区隧道的正常运行。冻害常造成衬砌风化、拱顶开裂、侧壁挤出等破坏,对围岩的长期稳定构成潜在威胁。一般的隧道围岩经长期的地质构造作用而富含节理裂隙等弱面结构,尤其在断层带内,裂隙更为发育。本文尝试建立含裂隙的隧道模型,并模拟低温环境下的THM耦合过程。本文采用ANSYS建立含裂隙的隧道模型,并划分网格,之后导入FLAC3D中进行计算分析,计算流程如图6所示。建立的隧道计算模型如图7所示,计算模型的尺寸为20m×60m×40m(长×宽×高),隧洞毛断面拱高8m,边墙宽12m,衬砌厚度0.4m。裂隙内含水初始水压0.5MPa。裂隙开度为2mm,表面设置接触面Interface单元(见图8),为透水界面。按平面应变问题考虑,隧洞走向前、后界面固定位移。计算中模型顶部施加10MPa初始地应力,水平地应力为5MPa+γH(渐变面力荷载),其中,γ为岩石容重,H为单元体距顶面距离。初始岩体温度为0℃,裂隙水温度为1℃,隧洞内气温为-20℃。裂隙介质和岩石基质分别定义不同的渗透系数、热传导系数、比热容等物理力学参数值,裂隙的渗透系数远大于岩石。低温多场耦合的重要特征之一是相变,耦合过程中,相变潜热的主要作用体现在降温热传递过程中水–冰系统的温度较未发生相变时偏高,可理解为潜热导致冰–水系统温度上升,其增加量计算公式如下:其中,式中:Cwi为冰–水系统的体积热容。计算中,对于冻结率u,FLAC3D没有自带的内置变量,需通过Fish函数定义。为简化计算,定义其在冰点(本文取-2℃)~-20℃区间内近似呈线性分布,T=-2℃时取u=0,T=-20℃时取u=1。同时,相变导致围岩的渗透系数变化,设围岩单元的渗透系数与冻结率线性相关。定义潜热系数变量L,根据式(25),(26),计算中调用该函数对温度进行调整,同时采用Fish遍历单元法对不同Group单元渗透系数进行调整。参考刘泉声等的研究结果,采用等效热膨胀系数方法模拟冻胀过程,设水的冰点为-2℃,设置-2℃~-20℃范围内夹冰(含水)单元的热膨胀系数为负值,温度降低时体积增加,与常温的热胀冷缩效应相反。冻结约10h后,隧道围岩的温度场、孔隙压力及最大主应力分布分别如图9~11所示。受裂隙的影响,温度场、应力场及孔隙压力都呈现明显的不连续性,裂隙构造的影响得以显现。由图9可知,假定冰点为-2℃时,10h后冻结深度约0.5m。最大主应力主要以压应力为主,但部分裂隙尖端附近存在明显的拉应力集中区,且裂隙面两侧的压力分布等值线断续密集。计算结束时隧道围岩渗透系数如图12所示,可见,渗透系数呈现不均匀性。因冻结作用影响,隧洞内壁近场围岩渗透性明显降低。6基于anasas的双裂隙介质模型冻岩THM耦合是个十分复杂的过程。冻结率问题很早就引起大家的关注,但对冻结率的理论研究很少。本文基于能量守恒原理推导出冻结率表达式,可为深入研究冻胀效应提供一定的参考。以往关于冻岩问题的研究多将岩体视为等效连续多孔介质,本文结合前人关于双重孔隙介质理论的研究成果,将岩体视为岩块和裂隙组成的双重孔隙介质,推导出裂隙岩体低温THM耦合控制方程。算例分析中运用ANSYS建立含裂隙隧道模型并划分网格,之后导入FLAC3D进行运算,模拟低温环境下的含裂隙隧洞的THM耦合过程,计算结果中裂隙的影响得以显现,这是一次有意义的尝试。水–冰相变对渗流场的一个重要影响是渗透系数的变化。水结冰后,冰的填充密实作