软岩非定常分数阶蠕变模型研究

王红娟，邓辉, 唐梓豪，李万才

(成都理工大学 地质灾害防治与地质环境保护国家重点实验室，四川成都610059)

流变学的研究较为复杂，其中岩土体蠕变是流变的重要研究内容之一。蠕变被定义为在恒应力作用下应变与时间的关系。软弱工程岩体的蠕变特性对工程的影响有着很重要的影响[1]，加速蠕变阶段的描述既是重点，又是难点。

现如今，有大量的蠕变模型被提出用以描述岩石蠕变的变形特征。按照构建模型的方法，可分为元件和经验模型。赵宝云[2]等对砂岩短时三轴蠕变模型的黏性元件进行了非线性改进，与Kelvin模型串联，得到了一种新的黏弹塑性蠕变模型，该模型可正确表征蠕变4个阶段的特征。最后，基于BFUS非线性优化算法对该模型参数进行了识别，同时将模型拟合曲线与改进西原模型曲线进行了比较分析，验证了非线性黏弹塑性蠕变模型的正确性；郭佳奇[3]等用FC元件取代整数阶微积分Kelvin-Voigt流变模型中的牛顿体元件，形成了基于分数阶微积分的蠕变模型，该模型在保证了拟合精度的情况下减少了模型参数；陈家瑞[4]等在分数阶微积分理论的基础上，基于三参量流变模型，提出了一个新的流变模型，结合数学软件1stOpt拟合参数，最终得到的模型能较好地描述破碎泥岩的流变规律；杨文东[5]等提出了由瞬弹性Hooke体，黏弹塑性村山体和非线性黏弹塑性体串联而成的岩石非线性黏弹塑性流变模型，该流变模型可以描述岩石流变的减速、等速特别是加速蠕变阶段，且通过反演，该模型的可靠性也得到证实；张明[6]等就锦屏一级水电站的大理岩进行流变模型的选择，提出了Hooke-Kelvin模型，该模型较Burgers模型更适合表达大理岩的蠕变特征；徐卫亚[7]等提出了由黏弹塑性体和五元件线性黏弹性模型串联的河海模型，该模型能较好地反应岩石加速流变的性质。以上大多只通过与某一种岩石的试验数据进行拟合来验证模型正确性，不具有广泛代表性。

综合以上考虑，本文基于分数阶微积分，提出了一种能正确表征加速蠕变阶段的蠕变方程，以西原模型为基础，引入了加速流变起始时间点[8]，采用软体元件表示非线性阶段，提出了一种参数可变的非线性黏弹性蠕变模型，给出了该模型的本构方程。并使用不同类型岩石在不同围压下的蠕变试验数据对该模型进行了验证，结果表明，试验数据与预测曲线都能较好地吻合，验证了模型的正确性。

1 分数阶微积分定义简介

本文提出的蠕变模型需要借助于分数阶微积分理论，故对分数阶微积分做如下简单介绍。

分数阶微积分有多种不同定义。在岩土工程学科中，Riemann-Liouville定义应用最广。以下对R-L定义进行简单介绍[9]。对于函数f(t)的β阶积分定义为

(1)

式中，d为微分算子。分数阶微分则定义为:

(2)

式中，β>0，且n-1<β≤n(n为正整数)。

Γ(β)为Gamma函数，其定义为：

(3)

2 基于分数阶微积分建立的软体元件

理想弹性体的应力应变关系满足Hook定律，即σ(t)～ε(t)，其蠕变形式表达为：

(4)

式中σ——常应力；E——弹性模量[10]。

(5)

式中η——黏性系数。

(6)

式中η——黏性系数；β——求导阶数。

当应力为常数时，即σ保持不变。元件将描述流变行为的蠕变，对式(6)两边进行分数阶积分，根据分数阶R-L积分定义，可得：

(7)

当β=0或β=1时，即Γ(β+1)=Γ(1)或Γ(β+1)=Γ(2) 时，式(7)退化为Hook体和Newton 体。

3 分数阶流变模型的建立

3.1 线性黏弹性流变模型

为得到一个应用范围广、适用性强且能更大程度模拟岩石蠕变过程的模型，本文建立了一个由弹簧体、开尔文体和一个可反应加速流变起始时间的黏塑性元件串联组成的蠕变模型，见图1。

在线性流变模型中，当σ<σs时，即i和ii共同作用时[11]，应变有如下关系：

εi,ii=εe+εve

(8)

式中εe——胡克体中弹簧的弹性应变；εve——黏弹性应变；εi,ii——瞬时蠕变和稳态蠕变的总和。Hook体中，应力应变有如下关系：

σ=EH·εe

(9)

式中，EH为胡克体中弹簧的弹性模量。Newton/Kelvin体中，应力应变关系如下：

(10)

εe1为开尔文体中弹簧弹性应变，γ意义如上，由并联模型性质可知：

εe1=εve=εii

(11)

所以得到：

(12)

即

(13)

对式(13)进行Laplace变换得：

(14)

式中，E(s)为复变量s的函数，解得：

(15)

对式(15)进行拉普拉斯逆变换得：

εve(t)=

(16)

3.2 非线性黏塑性元件

李良泉[8]假定破坏荷载和破坏时间、破坏荷载和加速流变历时存在一定的关系，可用指数形式来表示为：

(17)

式中t——破坏时间；σ0——岩石瞬时强度；σ∞——长期强度；Δt——加速流变历时；α和β均为常系数，则可得岩石加速流变的起始时间为：

(18)

徐卫亚[7]就绿片岩三轴流变试验曲线，提出了一个可反映加速流变特性的非线性黏塑性体。本文引用其形式，同时引入一个新的参数Δt，可得非线性阶段蠕变方程为

(19)

即：

(20)

3.3 考虑参数可变的流变模型

现有的蠕变模型大多是由元件模型经过串并联组合而成，且力学参数很多为常参数，但在实际情况下，如黏滞系数和弹性模量等的取值都是时间的函数。丁志坤[12]根据试验，提出了EH是定常而EK是时间t的函数，并给出了EK的表达式为：

EK(t)=p1+p2ep3t

(21)

式中，p1、p2、p3为常参数，本文考虑了弹性模量为非定常的情况，综合三部分应变，可得到岩石非线性黏弹塑性本构模型，该模型的本构方程为：

ε(t)=

(22)

4 模型验证

本文所提出的五元件蠕变模型，结合黄海峰等[13-15]对不同种类岩石的蠕变试验数据对模型进行拟合。红砂的应力大小和长期强度分别为55、23.56 MPa，岩盐的应力大小和长期强度分别为45、37.8 MPa，泥岩的应力大小和长期强度分别为40、36.21 MPa，由于所施加应力水平均大于长期强度，故需考虑加速阶段。此时待拟合参数为11个，这是由于本文考虑了变参数[15]的影响，虽然拟合参数较多，但也可通过拟合工具快速获取。本文利用数学软件1stOpt对式(22)中的参数进行拟合，所得结果见表1。经计算得知：红砂岩、盐岩和泥岩的相关系数R2的值分别为0.97、0.96和0.99。试验值与拟合曲线对比见图2，可见试验值和拟合曲线符合很好，这是由于采用了分数阶进行建模[17-19]。表1和图2仅仅展示了围压100 MPa下的拟合效果，为了进一步检验该模型，采用200 MPa下的红砂岩蠕变数据模拟，结果见图3。可见，新模型能够较好地反映软岩的蠕变特性。

表1 非线性蠕变模型参数

通过对比曲线可以发现，对于3种不同的岩石来说，瞬时蠕变阶段和稳态阶段都能较好地吻合，非线性阶段的试验曲线也能和预测曲线大致重合，证明了该模型的适用性和正确性。

5 结论

本文依据元件模型，在元件本构方程基础上引入非定常参数，使模型能更准确描述岩石蠕变全过程，通过红砂岩、岩盐以及泥岩的蠕变试验数据对该模型的验证，得到以下结论。

a) 通过对红砂岩，岩盐以及泥岩试验数据的拟合，发现该模型能较好地拟合岩石流变过程，加速阶段同样吻合，具有合理性，同时也证明了模型的适用性。

b) 所建模型能较好地模拟软岩的蠕变全过程，便于使用。在今后的研究中，将考虑把模型嵌入数值模拟软件之中，以期将其很好地应用于工程实际之中。