随着人类对地下围岩工程的开发，工程围岩有时处于高温环境中，如在地下煤炭气化炉[1-2]、隧道灾难性火灾[3-4]、核废料处置库[5-6]、增强型地热开发系统EGS[7-8]、工业废热储存[9]、油页岩储层的原位开采[10]、聚光太阳能发电厂CSP中热能储存系统TES[11-12]，均会涉及经高温作用后围岩热破裂的问题．在地下核废料处置库中，当高温热应力对围岩造成热破裂后，便存在着地下水穿透屏障系统，向人类环境输送放射性核素的潜在风险[5-6]．高温将会影响岩石内部的孔隙结构，较高的温度还会诱发热破裂，显著弱化岩石力学性能，降低围岩稳定性，对地下围岩工程形成安全隐患．因此，开展高温作用后的岩石力学研究，对于正确评价高温对围岩稳定性的影响具有重要意义．目前，在高温作用后岩石的物理与力学性质方面，国内外已取得了大量成果，主要表现在矿物组分变化、微裂纹分布、岩石渗透率、力学强度与弹性模量、泊松比等方面[13-15]．然而，由于岩石内部分布着大量跨尺度孔隙，在外部载荷作用下，岩石的力学行为是细观孔隙的宏观表现形式，因此有必要探究高温作用后每种孔隙与宏观力学参数的内在关联．岩石本构模型能够从数学逻辑角度客观反映出岩石力学温度效应的物理意义．自从1976年Dougill首次将损伤力学引入到岩石类材料以来，关于岩石热损伤的研究已取得了较多成果．刘泉声[16]、邵宏甫等[17]将弹性模量作为评价岩石损伤的基本参数，推导出了岩石损伤本构模型；张连英等[18]利用该模型，展示了大理岩由脆性向延性断裂的转化过程．在此基础上，Xu等[19]则建立了三轴压缩条件下岩石热损伤的统计本构方程．然而，这些研究并未考虑岩石的压密阶段，对于热处理后富含孔隙的多孔岩石，传统本构模型的适用性可能会降低．基于此，本研究基于高温处理后的砂岩力学试验，发现在温度效应影响下，岩石力学性能的变化呈现出双重性；然后，利用分形理论，研究了跨尺度孔隙各自与力学性能的内在关联；最后，建立了考虑峰后应力跌落的本构模型，该模型能够合理反映出高温作用后砂岩应力应变曲线的变化特征．1 试验过程1.1　试件制备与加热试验本研究选择砂岩作为研究材料，为了尽量保证试件的均质性，所有试件均从同一块母岩获取．将岩石预制成直径为50 mm、高为100 mm的圆柱体，岩石上下两个端面的非平整度不超过0.02 mm，以满足国际岩石力学学会(ISRM)对岩石力学试验的要求．然后，测量试件的基本物理参数，例如密度与波速．将部分密度与波速异常的试件舍去，剩余试件进行加热试验．选择可编程高温炉用于加热岩石，温度值分别为150，300，450，600，750，900与1 050 ℃，该温度范围基本涵盖了目前高温岩石力学所涉及的温度[13-15]，每个温度组均选择5个试件．本试验选择足够慢的2 ℃/min加热速率对试件进行加热，随后将砂岩试件置于目标温度下2 h以确保砂岩内部被均匀加热．1.2　岩石核磁共振与力学测试加热实验完成后，首先使用真空饱和装置使砂岩处于饱水状态，再利用核磁共振微观分析系统获取砂岩的孔隙分布数据(图1)．根据Carr-Purcell-Meiboom-Gill(CPMG)序列测试得到信号强度与弛豫时间的衰减曲线，通过数值反演获得各种尺度孔隙的体积占比．须说明的是，在核磁共振测试中，砂岩经历了饱水处理．高温作用后经饱水处理的岩石，与只经历高温处理的岩石相比，两者的力学行为可能存在一定差异．因此，在高温热处理完成后，应将每组中的2个试件进行核磁共振测试，剩余3个试件进行后续的单轴压缩试验．10.13245/j.hust.240216.F001图1试验装置在单轴压缩试验中，在岩石两端部均匀涂抹凡士林以降低端面摩擦，随后利用WHY600型伺服压力试验机进行力学试验．为了尽量模拟出工程围岩在外载荷作用下的变形与破坏，选择加载速率为0.12 mm/min[5]．2 高温对岩石力学性能的双重影响图2为高温作用后砂岩在力学试验中的应力(σ)-应变(ε)曲线．通常来讲，岩石应力应变曲线可分为压密(INS)、线弹性(LES)、塑性(PS)与峰后应力跌落阶段(PPFS)．众所周知，压密阶段与岩石内部的张开型微裂纹的存在有关．砂岩起始的变形主要表现为张开型微裂纹的闭合，较小的外界载荷能够引起较大变形．理想状态下的晶粒膨胀是能够完全恢复的，但岩石属典型非均质材料，晶粒并非呈现简单的规则排布，而是各向异性的复杂排列．高温使晶粒间的接触关系发生差异性变化，砂岩冷却至室温后，通常矿物晶粒间的接触关系只能得到部分恢复[12]．在压密结束后，岩石系统相继进入线弹性与塑性阶段．在超过峰值应力之后，砂岩便进入峰后阶段．从材料力学角度讲，砂岩的脆性与延性互相转变．在150~450 ℃内，力学曲线的起始非线性阶段不明显，砂岩在达到抗压强度后的破坏速率较快，砂岩呈现出较强的脆性破裂特征；反之，在600~900 ℃内，砂岩的起始非线性阶段较明显，破坏速率较慢，具有较强的延性破裂特征．10.13245/j.hust.240216.F002图2高温作用后砂岩的力学曲线根据砂岩峰值应力的变化可知，在25~300 ℃内，岩石强度被强化．从物质构成上讲，岩石不仅包含各种成矿矿物，还有结合水、矿物结构水等．当岩石经150与300 ℃处理后，水的溢出增大了矿物间的摩擦力，从而使岩石力学性能被强化；当温度高于450 ℃时，石英相变与矿物热分解使岩石力学性能被弱化，这便是高温对岩石力学性能的双重影响．3 砂岩跨尺度孔隙的分形特征通过傅里叶变换能够得到砂岩内部各个位置的CPMG T2脉冲序列信号[20]，该信号是不同大小孔隙内1H信号叠加的结果，根据文献[21]可知，R=ρ2FsT2，(1)式中：R为孔隙半径；ρ2为表面驰豫强度；T2为横向弛豫时间；Fs为孔隙几何因子(Fs，球体=3，Fs，圆柱体=2)．由于岩石内部孔隙的空间形状呈各向异性分布，因此选择ρ2Fs=0.01 μm/ms[22]．根据式(1)可知：在多孔材料砂岩中，R与T2成正比，当T2增大时，R变大；反之，R变小．T2谱曲线的纵轴可用来表示某种孔径孔隙含量的高低．根据本研究中砂岩孔隙的分布特点，按照孔径增大的顺序，将砂岩内部不同孔径的孔隙依次称为微孔、小孔、中孔、大孔．由于T2谱曲线是岩石孔隙的基本参数，因此图3仅展示了几种高温作用后砂岩孔隙的分布特征，图中由左向右四种颜色区域分别代表微孔、小孔、中孔与大孔(核磁共振信号强度无量纲)．不难发现，砂岩的T2谱呈现多峰分布，高温将会使砂岩T2谱发生显著变化．详细的孔隙分布讨论见文献[12]．10.13245/j.hust.240216.F003图3砂岩孔隙的T2谱曲线3.1　砂岩孔隙的分形特征根据分形理论可知，物体的数目与测量尺度存在幂律分布，具有相同孔隙半径的孔隙数量(N)与孔隙半径有以下关系，N(R)=aR-D，(2)式中：a为比例常数；D为孔隙结构的分形维数．假设岩石的孔隙空间是由一系列毛管束组成，根据Brooks-Corey模型可知，孔隙半径小于R的累积体积分数Sv为Sv=(R/Rmax)3-D，(3)式中Rmax为最大孔隙半径．进一步可获得表征孔隙结构分形维数的表达式为lgSv=(3-D)lg {T2}s+(D-3)lg {T2,max}s，(4)式中T2，max为横向弛豫时间的最大值．表1为砂岩中四种孔隙对应的分形维数，表中T2k (k=a～h)为分布范围．不难发现，高温能够显著改变砂岩跨尺度孔隙的分布范围．a．对于微孔而言，孔隙的最小孔径整体呈现先增大、后减小的变化．相比于未加热的岩石，经150~450 ℃处理后，微孔的最小孔径变大，表明高温对砂岩微孔隙/微裂纹的影响表现为闭合效应．当温度较高时，严重的热破裂与矿物分解，显著劣化了砂岩微观结构，拓宽了孔隙间距(拓宽效应)．b．对于小孔、中孔而言，150~750 ℃的高温不足以改变这两种孔隙的分布范围；仅当温度为900 ℃时，小孔的最大孔径被压缩，这仍与拓宽效应有关．c．对于大孔而言，随温度的升高，最大孔径先减小(孔径范围收缩)、后增大(孔径范围扩大)．10.13245/j.hust.240216.T001表1砂岩跨尺度孔隙的分形维数温度/℃微孔小孔中孔大孔T2a/msT2b/msDT2c/msT2d/msDT2e/msT2f/msDT2g/msT2h/msD250.320.98-0.291.059.662.6910.3595.482.92102.341 245.892.9941500.420.98-1.211.059.662.6210.3595.482.85102.34943.792.9973000.340.98-0.811.059.662.6410.3595.482.83102.34270.502.9774500.400.98-2.271.059.662.6110.3595.482.86102.34270.502.9766000.130.980.471.059.662.4610.3595.482.75102.34505.262.9847500.070.981.971.059.662.8010.3595.482.58102.34410.272.954900———1.382.25-5.9813.6795.482.10102.341 162.322.980Qin等[23]指出岩石孔隙分形维数的有效范围为2~3，当某种孔隙的分形维数小于2时，表明该种孔隙的分形性质较差或不具有分形特征．由此可知，砂岩内部的微孔不存在分形特征．除去900 ℃组以外，对于小孔而言，其分形特征是显著的．在本研究中，高温作用后砂岩的中孔与大孔，均具有明显的分形特征．不难发现，随着温度的升高，小孔及大孔的分形维数呈现波动变化，中孔的分形维数基本呈现逐渐减小的趋势，即中孔分形维数存在显著的温度效应．3.2　孔隙结构与力学参数的关联性分析由于高温将会改变岩石内部的孔隙结构，因此细观跨尺度孔隙的分布特征应与宏观力学性能存在某种内在关联．图4为小孔、中孔、大孔分形维数与力学参数的关系(因在900 ℃组中，小孔不具有分形特征，故未显示)，图中：U为强度；E为弹性模量．可结合图4与表1获取不同温度处理后每种孔隙的分形维数．根据图中孔隙分维值与力学参数的关系可知：对于孔隙直径为0.01~0.10 μm的小孔，其孔隙结构越均质，力学性能越强，相反，越复杂的孔隙结构，力学性能越差；对于孔隙直径为0.1~1 μm的中孔、与大于1 μm的大孔，其孔隙结构越均质，力学性能越弱，相反，越复杂的孔隙结构，力学性能越强．10.13245/j.hust.240216.F004图4砂岩孔隙分形维数与力学参数的关系4 考虑峰后应力跌落的本构模型上述研究结果表明，高温处理使砂岩的孔隙分布与力学性能发生了显著变化．为了深入了解砂岩力学性质的温度效应，建立温度效应下的砂岩本构模型更能客观反映出力学变化背后的数学逻辑．砂岩损伤演化方程的讨论见文献[24]．引入Drucker-Prager准则作为砂岩单元体的破坏判据，砂岩损伤C的公式如下[25-26]，      C=1-exp-xF0(1-DT)m0(1-DT)=1-exp-α0I1+J2F0(1-DT)m0(1-DT); (5)α0=sinφ9+3sin2φ；(6)I1=σx*+σy*+σz*=σ1*+σ2*+σ3*；(7)J2=16(σ1*-σ2*)2+(σ2*-σ3*)2+(σ3*-σ1*)2, (8)式中：m0与F0为经高温作用后的砂岩韦布尔分布参数；DT为损伤变量；α0为中间参量；φ为内摩擦角；I1为应力第一不变量；J2为应力偏量第二不变量；σi*(i=1，2，3)为有效应力；x为微元体强度．经过推导，可得到在单轴压缩试验中，考虑温度效应的岩石本构方程，σi=Eεi(1-η+ηexp[-(f/F)m])，(9)式中：εi为应变值；η为损伤修正系数；m与F为韦布尔分布参数；f =α0I1+J2．4.1　岩石本构模型的构建岩石经高温处理后，内部将会产生大量跨尺度分布的孔隙．在单轴压缩试验中(无围压条件)，受孔隙闭合过程的影响，应力-应变曲线的压密阶段表现为凸出形态．设dj (j为不同温度值)为模型值与试验值差值的绝对值，图5为传统损伤本构模型与试验值的对比(以温度为600 ℃时为例)，结果显示两者误差较大．这表明上述传统的损伤本构模型无法客观反映多孔砂岩在单轴压缩条件下的力学行为．10.13245/j.hust.240216.F005图5在砂岩单轴压缩试验中传统本构模型与试验值的对比由此可知，研究单轴压缩条件下经高温处理后的砂岩本构模型，必须将岩石的压密阶段与峰后应力衰减阶段作为重点．根据韦布尔分布函数特征可知，压密阶段内应力应变曲线的形状由exp[-(f/F)m]所主导．岩石力学曲线呈现显著的向上凹形特征，这与图5中模型值的形状特点相反．邓华锋等[27]提出了将exp[-(f/F)m]代替为1-exp[-(f/F)m]的方法．本研究基于该方法，对式(9)进行修正．当岩石的压密阶段完成后，岩石将先后进入弹塑性阶段与峰后应力跌落阶段．此时损伤本构形式与原形式不变，只须调整应力与应变的范围即可．同时，对峰后应力跌落阶段增加一个指数衰减函数．至此，高温作用后砂岩应力-应变全曲线的本构模型为σi=     Eεi1-exp-1FEεi3σ1∙σcsinφ3+sin2φ+σcm    (εiεe);     g(n)E(εi-εe)∙1-η+ηexp-1FE(εi-εe)3σ1∙σ1sinφ3+sin2φ+σ1m+σe    (εi≥εe),g(n)=1    (εiεc);(εc/εi)n    (εi≥εc),式中：g(n)为表征峰后应力非线性衰减的函数；n为岩石内部空间中呈各向异性分布的微元体对整体力学响应程度的指标，本研究将n称为峰后应变软化因子；εe为压密结束时的应变；εc为峰值应力所在应变．4.2　分段式统计本构模型的验证目前，在常见的算法效果评估函数中，均方根误差(RMSE，δ)能够用于比较不同数据集间预测误差的精确度．为了确定合理的n，对理论值与试验值做RMSE处理，设峰后应力应变曲线的数据点个数为K，σi为试验值，σi'为模型理论值，RMSE的计算公式为δ=1k∑i=1K(σi-σi')2．表2为砂岩n对RMSE的影响．当RMSE最小时，理论计算值与岩石试验值的误差最小，表明此时的n是合理的．10.13245/j.hust.240216.T002表2不同峰后应变软化因子下的RMSE值温度/℃n02468253.754.174.595.005.381503.572.672.051.811.953008.026.275.084.384.044501.231.852.543.213.836008.268.548.768.959.117500.531.141.882.533.07在确定了合理的n后，便可以构建考虑峰后应力衰减的分段式统计本构模型，图6为在单轴压缩试验中砂岩理论值与模型值的对比结果，限于篇幅，图中仅展示3个温度组．结果表明：在岩石被压缩过程中，分段式统计本构模型呈现出显著的压密阶段、弹性阶段、塑性阶段与峰后跌落阶段．10.13245/j.hust.240216.F006图6分段式本构模型与砂岩应力应变曲线的对比该模型能够客观反映出含大量微观热裂纹与孔隙的岩石在压密阶段内力学曲线的变化特征．同样在弹塑性阶段内，该本构模型也能较精确地表征出岩石的塑性变形能力．由于理论值与试验值的符合性较好(相比于峰值应力σc，峰前di曲线对应的数值较小)，这表明考虑多孔砂岩压密阶段的分段式本构模型是合理的．5 结论a．砂岩力学性能的温度效应，表现为在150~300 ℃内，力学性能被显著强化；温度超过450 ℃后，力学性能被逐渐弱化．随着温度升高，温度对砂岩内部微裂纹的影响先为闭合效应、后为拓宽效应．影响大部分微裂纹呈现闭合或张开的临界温度约在300 ℃附近．b．高温作用后砂岩内部的微孔不具有分形特征，小孔(除了900 ℃组)、中孔与大孔具有显著的分形特征．小孔分形维数与力学参数呈负相关，较大孔径的孔隙(中孔和大孔)与力学参数呈正相关．c．在温度效应因素下，建立了考虑峰后应力跌落与残余强度的砂岩本构模型．鉴于岩石实际峰后应力呈现“摆动式+跳跃式”的衰减过程，引入均方根误差，准确获得了合理的峰后应变软化因子．