人工冻结斜井洛河组砂岩解冻后蠕变模型试验研究

任建喜，易 归，陈 旭，曹西太郎

（西安科技大学 建筑与土木工程学院，陕西 西安 710054）

陕西榆神矿区煤矿斜井所穿越的白垩系地层多为富水洛河组砂岩[1]。为了保证斜井井筒安全顺利的穿越富水洛河组砂岩层，需采用人工冻结法施工[2]。人工冻结斜井井筒洛河组砂岩解冻后物理力学性质复杂，内部含有各种裂隙和空洞等微观缺陷[3-4]，在长期应力作用下围岩内部微裂隙发生闭合、扩展直至宏观裂缝产生，导致洛河组砂岩破坏，进而引发井筒井壁破坏，影响煤矿安全生产。目前对于人工冻结斜井解冻后洛河组砂岩的蠕变损伤破坏机理研究还远未成熟，西部浅埋煤层人工冻结斜井解冻后洛河组砂岩蠕变损伤破坏特性研究对人工冻结斜井井筒解冻后灾害治理及长期安全运行具有重要理论意义和潜在的工程应用价值。

软岩蠕变力学特性的研究较早。Li Y S 等[5]、Brijes Mishra 等[6]先后采用砂岩、粉砂岩及页岩为研究对象，进行了一系列的单轴和三轴蠕变试验，得出岩石蠕变开始及岩石蠕变破坏的阈值。国内陈宗基院士等[7]对宜昌砂岩进行了蠕变力学特性的试验研究；梁冰等[8]进行了片麻岩的蠕变力学特性的研究。国内外学者对于冻融循环后及低温饱水冻结状态下的岩石蠕变试验的研究也较多[9-13]。单仁亮等[14]通过对冻结温度为-10 ℃的冻结红砂岩进行了3 种不同应力水平下的单轴压缩蠕变试验，研究表明，冻结红砂岩加速蠕变阶段门槛值为13 h；张强勇等[15]分析得到了温度及应力对于片麻状花岗岩的蠕变力学特性的影响规律；陈国庆等[16]对经历不同冻融循环次数后的石英砂岩进行常规三轴和三轴蠕变试验，分析了冻融作用对岩石在短期受荷和长期受荷下的力学特征。

关于软岩蠕变理论模型方面也取得了丰硕的成果[17-20]。杨秀荣等[21]完成了软岩在不同含水状态下的三轴压缩蠕变试验，提出了1 个对西原模型进行改进的蠕变损伤模型；曹煜等[22]通过把定常黏性元件替换成非定常黏性元件，提出了1 个对Burgers 模型进行改进的蠕变模型；Tang H 等[23]提出了1 种基于变阶分数阶导数和连续损伤力学的四元蠕变模型；周宏伟等[24-25]结合分数阶微积分理论，提出了1 种新的黏性元件（Abel 黏壶），以此构建了1 种全新的非线性分数阶蠕变模型；宋勇军等[26]通过引入可以表现随应力水平及时间变化的损伤变量，结合分数阶微积分理论，提出了1 种四元件的分数阶蠕变损伤模型；苏腾等[27]人通过引入损伤系数来描述岩石加速蠕变阶段的非线性变化，并且提出了1 种变阶分数阶黏弹塑性蠕变模型，通过试验数据拟合验证了模型的适用性。

综上所述，国内外学者对于岩石的蠕变力学特性及蠕变模型有较多的研究成果，但对于解冻后砂岩蠕变力学研究依然较少；另外，部分学者所提蠕变模型能够很好地描述岩石蠕变过程，但由于其过于复杂，参数较多，从而导致适用性降低。为此，依托可可盖煤矿冻结斜井工程，进行解冻后砂岩的单轴压缩蠕变试验，研究解冻后洛河组砂岩的蠕变破坏规律，并以此建立更为简便的分数阶蠕变模型；基于研究成果，为富水洛河组砂岩地层中人工冻结斜井解冻后井筒的长期稳定性评价与煤矿安全生产提供技术保障。

1 试验方案

1）试样的制备。此次试验岩样取自陕西省可可盖煤矿，将砂岩加工成φ50 mm×100 mm 的国际标准圆柱体试样。通过室内的物理力学试验，分析得到洛河组砂岩的物理力学性质指标见表1。

表1 洛河组砂岩的物理力学指标Table 1 Physical and mechanical indexes of sandstone in Luohe Formation

2）试验方案。通过单轴压缩蠕变试验，考察不同冻结温度解冻后洛河组砂岩的蠕变力学特性。具体的试验步骤为：对遴选后的洛河组砂岩试件进行饱水处理，在不同冻结温度下（-5、-10、-15 ℃）饱水冻结24 h，常温饱水解冻24 h。采用加载速率为0.02 MPa/s 的应力控制方式将轴（围）压加至初始静水阶段。待围压恒定后，以0.02 MPa/s 速率应力控制方式加载轴压，蠕变加载应力σ1i′（i=1，2，3，4）分级为单轴压缩试验峰值强度的10%、30%、50%、70%，每隔12 h 后加载1 级，直至试件破坏。四级蠕变荷载分别为0.8、2.4、4.0、5.6 MPa。

3）试验设备。本次蠕变均采用长春朝阳试验仪器有限公司生产的TAW-1000 微机控制高温三轴蠕变试验机。TAW-1000 岩石三轴压缩蠕变试验系统主要参数：轴向最大试验力1 000 kN；最大围压80 MPa。变形测控范围：轴向0～10 mm，径向0～8 mm；位移测控范围：0～100 mm。

2 砂岩单轴压缩蠕变试验结果

2.1 不同冻结温度下解冻后砂岩蠕变量的分析

不同冻结温度下Tf（-5、-10、-15 ℃）解冻后的砂岩岩样进行单轴压缩蠕变试验结果如图1。

图1 不同冻结温度下砂岩分级加载蠕变曲线Fig.1 The creep curves of sandstone under different freezing temperatures

分析图1 可知，不同冻结温度下解冻后砂岩蠕变曲线趋势基本一致，整体上呈现阶梯状上升变化，共进行了4 级加载，前3 级荷载作用下，砂岩经历衰减蠕变阶段与稳定蠕变变形2 个阶段，最后一级荷载作用下，砂岩经历瞬时蠕变变形、稳定蠕变变形及加速蠕变变形3 个阶段。随着冻结温度的不断降低，砂岩每一级荷载作用下蠕变量逐渐增大，且最终的蠕变量也有所增大，主要原因是，冻结温度的降低，洛河组砂岩的内部水分逐渐冻结，导致内部孔隙受“冰劈作用”影响增大，解冻后砂岩内部孔隙裂隙增大，导致低冻结温度下的洛河组砂岩解冻后受同一轴压荷载应力作用下，其应变量有所增大，而由于冻结循环次数仅有1 次，冻结温度对于砂岩内部损伤的影响有限。

为了进一步的了解每一级的蠕变过程，对于每一级荷载作用下的蠕变量进行统计，不同冻结温度时的每级蠕变量见表2。不同冻结温度下的蠕变变形曲线如图2。

图2 不同冻结温度下的蠕变变形曲线Fig.2 Creep deformation curves at different freezing temperatures

表2 不同冻结温度时的每级蠕变量Table 2 Creep per level at different freezing temperatures

分析表2 可知，随着冻结温度的降低，砂岩每一级的蠕变量逐渐增大。砂岩受低水平应力作用下蠕变量随着冻结温度的降低增长较多，由于长期荷载作用下砂岩内部的孔隙逐渐闭合，且冻结造成的孔隙在第1 级荷载作用下闭合明显，导致其蠕变量增长较多。随着荷载级数的增长，蠕变量增长率也逐渐增大，说明随着应力水平的增加，砂岩内部孔隙发育明显，导致其蠕变量增大。

由图2 可知，冻结温度的降低，冻结过程中砂岩内部损伤增大，导致其蠕变变形量的增大，且表现为线性增长的趋势。砂岩最大的蠕变量随着冻结温度的降低逐渐增大，但增大的速率随冻结温度的降低逐渐减小，说明冻结温度为-10 ℃时，每一级荷载对砂岩的蠕变力学性能影响较大，随着冻结温度的进一步降低，其影响逐渐降低。

2.2 不同冻结温度下解冻后砂岩蠕变速率的分析

对砂岩试样的每一级荷载作用下蠕变速率进行分析时，低应力水平作用下，砂岩稳定蠕变阶段其蠕变速率基本为0，因而此处并不作比较分析，主要分析第3 级荷载作用下砂岩的稳定速率k3随冻结温度的变化；第4 级荷载作用下砂岩的稳定速率k4－1及加速速率k4－2随冻结温度的变化，解冻后砂岩蠕变速率随冻结温度的变化曲线如图3。

图3 解冻后砂岩蠕变速率随冻结温度的变化曲线Fig.3 The change curves of sandstone creep rate with freezing temperature after thawing

分析图3 可知，第3 级荷载作用下蠕变速率随冻结温度的降低呈现线性增大的趋势，第4 级荷载作用下蠕变速率随冻结温度的降低呈现二次项线性增大的趋势，冻结温度从-5 ℃降至-10 ℃和-15 ℃时，稳定蠕变速率k4－1分别增大16.86%和33.8%；加速蠕变速率k4－2分别增大4.34%和15.1%，即冻结温度为-15 ℃时，解冻后的砂岩蠕变速率会急速增大，说明冻结温度的降低导致砂岩内部损伤增多，相同荷载作用下更容易破坏。随着应力水平的增加，第4级荷载作用下蠕变速率明显增大。

3 砂岩分数阶蠕变模型

分数阶微积分有多种定义，广泛应用是Riemann-Liouville 提出的分数阶微积分定义。周宏伟等[24-25]和宋勇军等[26]基于分数阶导数理论构建的1 种能够描述物体黏弹性性质的原件Abel 黏壶，可以描述处于流固体之间的材料的非线性应变过程。然而，对于砂岩的流变过程，尤其是加速流变阶段，砂岩黏性系数不可能是恒定不变的，苏腾[27]等将分数阶导数的常阶系数由时间函数所代替成为随时间变化的变阶系数，建立了1 个变阶数的Abel 黏壶元件。通过弹性元件、常阶数的Abel 黏壶元件以及变阶数的黏壶元件相连接便可较好的表达洛河组砂岩蠕变过程。但均过于复杂，对于砂岩材料解冻后蠕变特性，提出1 个参数更少，更为简洁的模型。

对于砂岩蠕变力学特性试验，根据恒定荷载σ0大小的不同，可以将蠕变分为2 种类型，即稳定蠕变和非稳定蠕变。当σ0≤σs（σs为砂岩蠕变应力阀值）时，砂岩的蠕变过程为稳定蠕变，其蠕变速率随时间逐渐减小，最终保持恒定；当σ0>σs时，砂岩的蠕变过程为非稳定蠕变，其蠕变速率随时间的变化为先减小，后逐渐增大，最终导致砂岩的破坏。

3.1 分数阶蠕变模型

3.1.1 稳定蠕变的蠕变模型

根据蠕变试验结果可知，加载的瞬间，砂岩便会产生与时间无关的1 个瞬时的弹性变形，因而可以利用弹性元件来模拟表达，即施加应力σ0与应变εe的关系可以表达为：

式中：ε 为总应变；t 为时间；ηα(t)为稳定蠕变过程中的黏性系数；α（t）为时间函数。

蠕变本构模型示意图如图4。图中，η1、η2为Abel 黏壶元件与变阶数Abel 黏壶元件的黏性系数；α1、α2可为t1～t2时段和t2～t3时段的分数阶阶数；参数ω 是可由试验确定的与材料有关的参数。

图4 蠕变本构模型示意图Fig.4 Schematic diagram of creep constitutive model

根据图4 可知，稳定蠕变的黏弹性应变为t0～t1时间段内，则总应变ε（t）可以表达为：

式中：εt0t1（t）为t0～t1时间段内稳定蠕变的黏弹性应变。

3.1.2 非稳定蠕变模型

对于砂岩的非稳定蠕变，t0时刻与稳定蠕变模型相同，砂岩产生1 个瞬时弹性变形，与时间无关，利用弹性元件模拟表达；t0～t1时段，砂岩呈现黏弹性，应变速率随时间不断减小的阶段被称为初始蠕变阶段，利用Abel 黏壶元件模拟表达；随后的t1～t2时段，砂岩呈现黏塑性，利用变阶Abel 黏壶元件模拟表达。根据叠加原理，砂岩的总应变ε（t）为：

引入Mittag-Leffler 函数：

3.2 分数阶蠕变模型参数确定

蠕变本构模型的参数识别是验证模型的重要一环。本研究中采用Levenberg-Marquardt 最小二乘法，通过1stOpt 数学软件实现，拟合得到的参数包括：η1、α1、η2、α2、w，不同冻结温度时分数阶蠕变模型参数见表3。根据不同冻结温度下的解冻后洛河组砂岩的蠕变试验结果拟合验证砂岩的蠕变模型，蠕变试验数据与理论模型曲线对比如图5，图中ε1t为蠕变应变量。

图5 蠕变试验数据与理论模型曲线对比图Fig.5 Comparison of creep test data and theoretical model curves

表3 不同冻结温度时分数阶蠕变模型参数Table 3 Fractional creep model parameters at different freezing temperatures

分析图5 可知，蠕变试验结果与模型拟合计算得到的结果比较可知，单轴压缩蠕变分数阶模型（式14）能够较好的模拟单轴压缩蠕变试验全过程曲线，尤其是加速蠕变阶段，其拟合效果更好。

分析表3 可知，随着冻结温度的降低，蠕变模型中参数E0、η1、η2、及ω 均有所减小，但其减小的幅度非常小，说明冻结温度的不同，对于解冻后砂岩的蠕变特性有所影响，但其影响程度较弱。

根据蠕变试验结果与模型拟合曲线的对比分析，可以得到非线性分数阶蠕变模型可以较好的模拟解冻后洛河组砂岩的蠕变过程，尤其是加速蠕变阶段R2至均在0.98 以上，说明模拟效果更为明显，且蠕变模型简单，参数较少，使用方便。

4 结 语

冻结温度的降低导致砂岩内部损伤增多，相同荷载作用下更容易破坏。解冻后洛河组砂岩呈现非线性蠕变力学特征，非线性程度与蠕变时间、砂岩所受应力水平及冻结温度均有关系，蠕变加载时间越长，荷载应力水平越大，冻结温度越低，砂岩的非线性程度越高。

基于分数阶微积分理论，构建了更为简洁的分数阶黏壶与含损伤的变阶分数阶黏壶，建立了解冻后砂岩的分数阶一维蠕变模型。通过最小二乘法确定了模型中相应的参数，经过试验数据与拟合曲线的对比分析，非线性分数阶蠕变模型可以较好的模拟解冻后洛河组砂岩的蠕变过程，该蠕变模型具有参数少且使用方便的特点。

基于以上研究成果，为富水洛河组砂岩地层中人工冻结斜井解冻后井筒的长期稳定性评价与煤矿安全生产提供理论支持。避免淹井、井筒漏水、井筒破坏等灾害的发生。