0 引言

当地下工程的基础向着深部不断推进时，基础所处区域内的深部围岩遭遇的地质条件也更加复杂多变，深部围岩在高围压、高温、高地应力以及强扰动条件下，表现出了复杂的流变特性，而岩石蠕变变形特性以及蠕变变形破坏机理是怎样变化发展的，对于研究材料蠕变破坏机理具有重要的意义。故需要通过对岩石室内进行蠕变试验研究，分析岩石在不同力学状态下蠕变特性，建立一种合适的蠕变本构模型来描述岩石的流变性质，对建立数值力学模型及解决工程实际问题都具有重大的影响^[1-4]。

近些年，学者们对于岩石流变性质、变形破坏机理以及如何构建一个合适的蠕变模型做了大量研究^[5]，主要表现在以下几个方面:张楚旋等^[6]发现围压可以有效地提升能量的耗散，且卸载试验的能量耗散要远远大于加载试验中岩石的能量耗散量。王凯等^[7]通过试验得到在巷道开挖过程中围岩的变形规律，结合数值计算软件对煤岩开采过程中变形特性进行模拟，最终围岩的耗能随卸载速率逐渐增大，且导致围岩进一步变形破坏。雷涛等^[8]通过外部手段对岩石进行不同初始损伤试验，进一步对含有初始损伤的岩石进行卸载试验研究，得出岩石能量耗散与孔隙度均随初始损伤程度的增大而增大。上述试验研究中对于岩石在加载、卸载过程中能量耗散规律进行了深入的研究，验证了通过能量耗散理论可以有效地描述岩石的变形破坏特性。为了进一步研究和预测岩石在不同条件下的变形规律，将结合岩石在蠕变过程中的能量耗散规律，反映岩石在蠕变变形过程中内部结构的劣化情况，从而更好地反映岩石的变形破坏行为。

因此，本研究对取自海棠山隧道围岩进行三轴蠕变试验，研究岩石在不同围压条件下变形劣化特性，然后结合现有最小能量耗散理论推导损伤演化模型，建立一个基于能量耗散理论的岩石非线性损伤蠕变方程，最后通过三轴蠕变试验数据验证模型的正确性与合理性。

1 三轴蠕变试验研究

本研究选取的砂岩取自海棠山隧道围岩，为减小试样的差异性，确保试验的可对比性，挑选的岩块均来自同一断面，海棠山隧道所在地区为低山丘陵地貌，山势陡峭为深埋隧道。本研究设定围压为20 MPa和30 MPa，试验每个轴向应力水平分成4级，采用分级加载的方式以达到逐级加载偏差应力。本次试验采用偏差应力加载方法见表 1。

试样的加工精度按照《岩石试验方法标准》(GB50218—94)^[8]进行，制备成形状为标准圆柱体，直径50 mm，高度100 mm。本次试验采用MTS812.02岩石试验机如图 1所示。

试验方案：首先，在试件两端面均匀涂抹凡士林，最大限度地消除端部效应；在加载时设定位移上限值，保证岩石变形不受外界干扰；加载方式采用荷载控制方式，以0.5 MPa/s的速率均匀施加轴向应力和围压，轴向应力水平每次加载10 MPa直至岩石失稳破坏为止，采样间隔设置为5 s；最终通过输出的试验数据绘制出轴径应变-时间历时曲线。

1.1 蠕变历时曲线分析

根据上述试验方法和试验方案，对海棠山隧道围岩进行不同围压的常规三轴压缩蠕变试验，σ₃=20 MPa和σ₃=30 MPa作用下的分级加载轴向历时蠕变曲线如图 2所示。

由图 2可知，在围压20 MPa作用下，荷载作用4级后发生失稳变形破坏，在围压30 MPa作用下，荷载作用也经历5级后发生失稳变形破坏；瞬时应变和蠕变应变都随着应力水平增大不断增大，且瞬时应变占总应变比例的变化规律呈现先减小后增大趋势。

1.2 轴向蠕变曲线分析

将不同围压条件下所得蠕变试验数据运用Boltzmann叠加原理^[9]进行处理，得到各级轴向蠕变曲线如图 3所示。

由图 3可知，只有当施加在岩石的应力水平较高时，才会出现稳定蠕变和加速蠕变现象，否则蠕变变形只有衰减蠕变变形；以围压30 MPa为例，在σ₁=80 MPa时，蠕变速率从初始值衰减至零，最终的蠕变变形为0.110 2%(全部为衰减蠕变)；在90 MPa和100 MPa作用下，衰减蠕变变形时间大约为8 h左右，岩石的蠕变变形便转变为稳定蠕变变形，蠕变速率从初始值衰减至一个稳定值；在最后一级应力水平110 MPa作用下，蠕变变形出现加速蠕变变形，变形与时间之间的规律呈现出明显的非线性变化趋势。同时，围压的增高不仅增大了试样的破坏应力水平，同时延缓了轴向蠕变变形。

1.3 蠕变试验瞬时变形分析

根据试验数据绘制不同围压下砂岩瞬时变形量与应力水平关系如图 4所示。

岩石在低应力作用下，不考虑岩石初始损伤^[10-12]，瞬时应变阶段的变形可以当做弹性应变，此时岩石内部空隙开始闭合，在施加较大应力水平后，瞬时应变不仅包含了弹性应变，也包含了塑性应变，此时岩石内部空隙在完全闭合后，会产生一定裂隙使得岩石体积由压缩状态转化为膨胀状态，瞬时应变非线性特性更加明显。

2 非线性损伤蠕变方程 2.1 损伤变量

通过引入开关元件来改进西原模型^[13]，则改进后的西原模型，如图 5所示。

改进后西原体一维蠕变方程为^[14-17]：

当 σ<σ_cr时

(1)

当σ_cr≤σ < σ_s时，

(2)

当σ_s≤σ时，

(3)

式中，E₀，E₁为弹性模量；η₁，η₂为黏滞性系数；σ_s为长期强度；σ_cr为开关元件启动应力；ε为应变；σ为应力。

根据应变等效原理^[18-19]将改进后的西原体一维蠕变方程表示如下

(4)

式中，D₁为弹性体劣化的损伤变量；D₂为弹性体和黏弹性体劣化的损伤变量；D₃为改进后整个西原体劣化的损伤变量。

假设不考虑对时间t求导所得应变速率，分别对上式关于时间t求导数得：

(5)

式中，为应变率，为损伤变量关于时间t的一阶导数，i=1，2，3。

根据材料内变量理论和最小耗能原理^[20]，在不考虑材料热耗散能时，其耗能率φ (t)为：

(6)

式中，为由耗散引起的不可逆非线性应变率；σ_i为不同主应力方向的应力。

本研究蠕变试验在常温条件下进行假三轴试验(σ₁>σ₂=σ₃)，屈服函数选取以下形式，即^[21]：

(7)

式中J₂为应力偏量第二不变量。

根据最小耗能原理知材料在损伤破坏过程中，岩石任意时刻在约束条件(屈服函数)下耗能最小，即

(8)

式中λ为拉格朗日比例常数。

联立式(5)~ (8)式可得：

(9)

式中，λ_i为拉格朗日比例常数，i=1，2，3；C₁，C₂和C₃为材料常数，可拟合得到; a=1/E₀+1/E₁× [1-exp(-E₁×t/η₁)]；b=exp (-E₁×t/η₁)/η₁。

根据文献[22]知，在不考虑体积变化的前提下，令S=σ₁-σ₃，将式(1)~ (3)转化为西原体三维状态下的非线性损伤方程，得出岩石的轴向蠕变损伤方程为:

当σ < σ_cr时，

(10)

当σ_cr≤σ < σ_s时，

(11)

当σ_s≤σ时，

(12)

式中 S为应力偏张量。

2.2 损伤蠕变参数确定

根据文献[23-24]可知, 通过等时应力-应变曲线(如图 6所示)可以得出围压30 MPa下砂岩的长期强度为σ_s=90 MPa。

本研究采用Origin软件对试验曲线进行参数辨识(参数见表 2)。由于所施加应力水平较大，远远超过蠕变起始应力门槛值σ_cr，不会存在只有瞬时应变且不产生蠕变变形的情况，故只考虑其他两种情况下的蠕变损伤模型与试验数据的对比。

将表中不同应力水平作用下的蠕变参数值代入到改进后的模型中，得到岩石在不同应力水平作用下模型曲线与试验数据对比如图 7(b)所示。同理，可得围压20 MPa蠕变损伤变量与时间关系的曲线如图 7(a)所示。从试验数据与模型曲线对比可知，结合能量原理和最小能耗原理确定内变量D，可以较好地反映岩石内部变形变化状态，并可以将岩石损伤演化规律以能量变化方式呈现出来，进一步说明岩石损伤贯穿于整个加载试验过程中，也可以通过损伤演化规律来描述岩石能量耗散规律。

将表 2中确定的蠕变参数代入式(9)中，绘制出不同围压条件下，砂岩试样最后一级损伤变量D₃与时间t的变化关系如图 8所示。

由图 8可知，损伤变量随着时间的增加，先呈现出近水平、小幅度增长，在蠕变变形进入加速蠕变后，损伤变量开始急剧增大，数值逐渐趋近于1，这说明了在变形进入加速蠕变变形阶段后，岩石内部空隙首先在软弱区积聚，使得岩石损伤程度进一步加大，内部空隙和裂隙逐渐贯通为宏观裂缝，最终岩石产生失稳破坏。

为了进一步验证本研究建立蠕变损伤模型的适用性，引用文献[25]中花岗岩蠕变试验数据，基于最小二乘法获得蠕变模型参数，然后将计算模型曲线与试验数据对比如图 9所示。

由图 9可知，本研究建立的蠕变损伤模型，对于不同围压作用下花岗岩蠕变特性也有较好地描述，计算曲线和试验数据拟合度较高，说明了该损伤模型不仅适用于描述砂岩的蠕变特性，也适用于描述花岗岩的蠕变特性。这是由于建立的损伤模型是考虑岩石耗散机制的基础上，而能量的耗散机制存在于任何一个自然作用中，故本研究建立的损伤模型可以适用于任何岩石的蠕变变形情况。

3 结论

以阜新海棠山隧道围岩(砂岩)为研究对象，对三轴蠕变作用机制展开研究，构建了一个新型的能量耗散蠕变模型。得出以下主要结论：

(1) 在施加较大应力水平后，瞬时应变不仅包含了弹性应变，也包含了塑性应变，此时岩石内部空隙在完全闭合后，会产生一定裂隙使得岩石体积由压缩状态转化为膨胀状态，瞬时应变非线性特性更加明显，且瞬时变形逐渐由线性转变为非线性，具有较为明显的转折点。

(2) 在相同围压、不同应力下砂岩蠕变试验曲线与模型曲线拟合度较高，说明了基于最小耗能原理来建立非线性蠕变损伤模型，对砂岩蠕变全过程演化规律描述是合适可行的，也充分地说明将岩石作为耗散结构来确定岩石损伤程度以及反映岩石蠕变全过程变形规律是正确的。

(3) 通过能量原理引入内变量D，可以较好地反映岩石内部应力-应变变化状态，也较好地呈现出岩石内部能量耗散具体过程，且该模型对于不同围压作用下花岗岩蠕变特性也有较好地描述，计算曲线和试验数据拟合度较高，说明了该损伤模型的适用性广泛。