一种新的岩石黏弹塑性流变模型

， ， ，

(1.中国地质大学(武汉) 工程学院，武汉 430074; 2.广西地质环境监测总站，广西 桂林 541004)

1 研究背景

岩石是一种具有流变特性的介质材料，主要表现为蠕变、应力松弛和弹性后效。在一定的应力和时间条件下，岩石就会表现出明显的流变特性。岩土工程中，岩石的流变特性与工程的长期稳定性和安全性有着密切的联系。因此，岩石的流变力学特征在岩土工程中得到了广泛的研究。

近几年，随着新理论、新方法的采用，岩石流变模型理论得到了很大的发展并取得了丰硕的成果。黄耀英和郑宏[1]对节理岩体等效损伤流变模型进行了初步研究，假设岩石为各向同性体且各向同性损伤，节理面的法向和切向损伤不同，分别建立了岩石和节理面流变损伤演变函数。徐卫亚和杨圣奇[2]对龙滩水电站无充填节理岩石进行剪切流变试验，提出了非线性黏塑性体，将其与5元件黏弹性剪切流变模型串联，建立了7元件非线性流变模型，并用试验数据验证了所建模型的合理性和正确性。N.D.Cristescu[3-4](1993，1994年)提出了可描述盐岩初始蠕变和稳态蠕变的弹黏塑性模型。宋勇军等[5]将含有分数阶微积分的软体元件与弹簧元件串联，结合1个幂函数黏塑性体，得到了一种新的4元件非线性黏弹塑性流变模型，并给出了该模型的本构方程和蠕变方程，得到了不同应力状态下的蠕变曲线。长期以来，研究者在岩石稳态流变和定常流变方面提出了多种成熟的流变模型。由于受到试验条件的限制，岩石流变模型理论并不是十分完善，特别是在岩石非线性流变方面还有待进一步地研究。因此，本文将在已有流变理论的基础上，同时借鉴徐卫亚、杨圣奇等的研究思路，对岩石非线性流变模型进行研究。

2 岩石黏弹塑性流变模型的研究

2.1 岩石黏塑性流变模型的建立

图1 岩石流变的3种典型曲线

岩石剪切流变会经历减速流变阶段、稳态流变阶段和加速流变阶段，如图1所示。由图可知，当剪应力τ0≤长期抗剪强度τs时，岩石流变处于减速流变和稳态流变阶段，如图中曲线1，2所示；当剪应力τ0>长期抗剪强度τs时，岩石流变将经历减速流变阶段和稳态流变阶段进入到加速流变阶段，如图中曲线3所示。图中t1为岩石从减速流变进入稳态流变的时刻；t2为岩石从稳态流变进入加速流变的时刻，可用回归系数二阶差分法确定；tc为岩石破坏的时刻。

对图1中曲线3的BC段，其变形量的增长不收敛于定值，具有明显的黏塑性特征。由于线性元件组合模型(如Burgers模型、Kelvin模型、Maxwell模型)无法精确地描述岩石在该阶段的流变特征，故需用一个非线性加速元件组合模型来描述其加速流变特征。岩石非线性流变模型有2种建立方法[6-11]:①将新的非线性元件代替原有的线性元件；②采用断裂及损失力学理论、内时理论等新的理论。

图2 岩石黏塑性流变模型

根据第1种方法，本文假定岩石非线性流变变形与时间的幂函数成指数型关系，提出了一种新的非线性黏性元件，并将其与一个塑性元件并联起来，构成了一个新的黏塑性流变模型，如图2所示。

在恒定剪应力τ0作用下，该模型的变形与时间关系曲线如图3所示，相应的流变方程为

(1)

式中：u为变形量(以剪切位移量表示)；η3为岩石流变模型参数；n为流变指数，反映岩石加速流变速率的快慢程度；t为流变时间。

图3 岩石黏塑性流变模型的流变趋势曲线

图3中曲线1，2，3，4分别表示式(1)中流变指数n取值为0.3，0.5，0.7，0.9时岩石黏塑性流变模型的流变趋势曲线。由图可知，曲线1并未出现明显的加速现象，而曲线2，3，4则表现出了明显的加速特征。当n>0.3时，岩石黏塑性流变模型中岩石的变形与速率均随着时间的增加而呈非线性增加，说明该模型可以用来描述岩石的非线性加速流变变形。

进一步分析可知，图3中曲线2以较快的加速度增大，而曲线3、曲线4则先经历缓慢增大的过程而后再迅速增大。这里将曲线2定名为A类流变曲线，曲线3、曲线4定名为B类流变曲线。对比研究徐卫亚和杨圣奇[2]提出的NVPB模型(参照文献[2]中图8)可知，NVPB模型的流变指数n取值大于1时，其变形量以较快的加速度增大，仅可出现A类流变曲线；而本文所提出的岩石黏塑性流变模型则可出现A，B两类流变曲线。通过改变本文所提出的黏塑性流变模型的流变指数，即可得到不同类型的流变曲线，因而该模型具有更广泛的适用性。

2.2 岩石黏弹塑性流变模型的建立

岩石流变是弹性、黏性、塑性、黏弹性和黏塑性等多种变形共存的一个复杂过程，需采用多种元件(线性和非线性元件)的组合来对其进行模拟[6]。对于图1中的曲线1、曲线2，可以采用徐卫亚和杨圣奇[2]提出的5元件黏弹性剪切流变模型进行了较好的拟合；而对于曲线3，则需用能够描述岩石黏弹塑性流变特征的模型来拟合。为了全面描述岩石的流变过程，将本文提出的岩石黏塑性流变模型与文献[2]中的5元件黏弹性剪切流变模型串联起来，建立一种新的岩石黏弹塑性流变模型，其模型如图4所示。

图4 岩石黏弹塑性流变模型

由模型的串联关系可知，建立的岩石流变组合模型应满足以下条件：

τ0=τ1=τ2;

(2)

u=u1+u2。

(3)

式中：τ1，τ2分别为模型的第1部分和第2部分所受的剪应力；u为模型总变形量；u1，u2分别为模型的第1部分和第2部分的变形量(如图4)。

根据τ0与τs之间的关系，所建立的岩石黏弹塑性流变模型的流变方程可分为以下2种情况：

(1) 当τ0≤τs时，岩石黏弹塑性流变模型的第1部分处于禁止状态，因此模型发挥作用的只是第2部分，由文献[2]可知其相应的流变方程为

(4)

式中：u(t)为变形量；G1为瞬时剪切模量；G2和G3均为黏弹性剪切模量；η1和η2均为黏滞系数。

(2) 当τ0>τs时，岩石黏弹塑性流变模型的第1和第2部分均处于启动状态，其相应的流变方程为

(5)

式中η3为黏滞系数。

显然，当τ0≤τs时，可用式(4)来描述图1中的曲线1和曲线2；当τ0>τs时，可用式(5)来描述图1中的曲线3。若将式(4)和式(5)合并成一个式子，则有

(6)

其中H(τ0-τs)为正值函数，即

(7)

2.3 一维通式与三维流变方程

根据2.2节中式(4)所示的黏弹性流变模型，可得到算子形式的通式[12]为

P(D)σ=Q(D)ε。

(8)

对于2.2节中的式(6)，当τ0<τs时，其一维通式与式(8)相同；当τ0≥τs时，若采用式(8)表达一维通式，应将式中的σ用σ-σs来代替，即一维本构方程的通式可表示为

(9)

在三维应力状态下，可将岩石材料内部的应力张量分成球应力张量和偏应力张量。对于球应力张量σm的表达式为

(10)

偏应力张量为

(11)

从而可得

σij=Sij+δijσm。

(12)

同理，也可将应变张量分成球应变张量εm和偏应变张量eij分别为

(13)

(14)

于是可得

εij=eij+δijεm。

(15)

设岩石弹性剪切模量为G、弹性体积模量为K，则二者与弹性模量E及泊松比μ之间的关系为

(16)

基于弹性理论，弹性本构关系的一维形式可表示为

σ=E0ε。

(17)

弹性本构关系的三维张量形式可表示为

Sij=2G0eij, σii=3Kεii。

(18)

式中：Sij为应力偏量；eij为应变偏量；σii为应力;εii为应变。

结合式(18)，对于三维状态下的黏弹性体，其流变方程可参照一维应力状态下的流变方程求得，但方程中的偏应变须用试验时的恒定偏应力Sij来代替，即

(19)

对于黏塑性模型的三维流变方程有

(20)

由式(19)、式(20)可得三维应力状态下的岩石流变方程为：

当Sij≤τs时，与式(19)相同;

当Sij>τs时，

(21)

3 模型参数的确定及应用

笔者曾利用精度为0.001 mm的JQ200型岩石剪切流变仪对西藏邦铺矿区花岗岩进行了岩石剪切试验。试验采用分级连续加载的方式，并保证试验过程中环境温度恒定。现以天然状态的花岗岩流变试验数据为依据，给出本文的岩石黏弹塑性流变模型的参数确定方法。

(1) 当τ0≤τs时，H(τ0-τs)=0，此时岩石黏弹塑性流变模型退化为5元件黏弹性剪切流变模型，即由式(6)转化为式(4)。对于5元件黏弹性剪切流变模型的参数求解方法可参照文献[2]。具体的求解方法是：首先根据试验数据确定瞬时剪切模量G1。由式(4)可知，当t=0时，G1=τ0/u(0)，u(0)为泥岩的瞬时弹性变形量，可从花岗岩流变试验曲线上直接得到，而τ0为已知量，因此可求得G1。然后再确定花岗岩黏弹性剪切模量G2，G3及黏滞系数η1，η2。这里借助Origin软件中自带的功能函数Exp-Assoc对花岗岩的剪切流变试验数据进行拟合，最终可获得花岗岩黏弹性剪切模量G2，G3及黏滞系数η1，η2。

(2) 当τ0>τs时，则有H(τ0-τs)=τ0-τs。此时岩石将出现图1中曲线3的流变状态，求解模型参数的具体步骤如下：

①利用花岗岩在(0，t2)时段内的流变试验数据，并根据求解5元件黏弹性剪切流变模型参数的具体方法确定G1，G2，G3与η1，η25个参数。

②在(t2，tc)时段内，利用得到的5元件黏弹性剪切流变模型求出该段时间内岩石的黏弹性理论解。然后用花岗岩剪切流变试验数据对相应黏弹性理论解作差，将得到n组(Δu，t)的数据。

图5 岩石黏弹塑性流变模型与试验结果的对比

③在Origin中输入得到的n组数据并绘制散点图。由于Origin中没有能够表达公式(1)的内置函数，故需自定义函数进行非线性拟合。在Origin中自定义函数并初始化系数后对n组数据进行非线性回归分析，即可得到岩石黏塑性流变模型的剪切流变参数η3和n。至此，7个剪切流变参数便可全部求得。

图5给出了天然状态下花岗岩出现加速流变的全过程试验结果与岩石黏弹塑性流变模型的对比曲线，求得的模型参数如表1所示。从图中可以看出，所建立的岩石黏弹塑性流变模型与花岗岩剪切流变试验结果吻合得比较理想，这表明了该模型的具有一定的正确性和合理性，因此可以用它来描述岩石的加速流变阶段。

表1 岩石黏弹塑性流变模型参数值

4 结 语

(1) 假定岩石非线性流变变形与时间的幂函数成指数型关系，提出了一种新的非线性黏性元件，将其与一个塑性元件并联起来，构成了一个新的黏塑性流变模型，并给出了模型的流变方程。

(2) 将新的岩石黏塑性流变模型与5元件黏弹性剪切流变模型串联起来，建立了一种新的岩石黏弹塑性流变模型。同时，基于花岗岩剪切流变试验，给出了模型参数的求解方法。

(3) 利用建立的岩石黏弹塑性流变模型对天然状态下的花岗岩流变全程曲线进行拟合，得到了模型的流变参数。流变模型曲线与试验结果的比较显示：本文所建立的岩石黏弹塑性流变模型具有一定的正确性和合理性。

参考文献：

[1] 黄耀英，郑 宏.节理岩体等效损伤流变模型初步研究[J].岩土力学，2011，32(12)：3566-3570.(HUANG Yao-ying, ZHENG Hong. Preliminary Study of Equivalent Damage Rheological Model for Jointed Rock[J]. Rock and Soil Mechanics, 2011，32(12)：3566-3570. (in Chinese))

[2] 徐卫亚，杨圣奇.节理岩石剪切流变特性试验与模型研究[J].岩石力学与工程学报，2005，24(增2)：5536-5542.(XU Wei-ya, YANG Sheng-qi. Experiment and Modeling Investigation on Shear Rheological Property of Joint Rock[J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(Sup.2): 5536-5542. (in Chinese))

[3] CRISTESCU N D. A General Constitutive Equation for Transient and Stationary Creep of Rock Salt[J]. International Journal of Rock Mechanics and Mining Sciences & Geo-mechanics Abstracts, 1993, 30(2):125-140.

[4] CRISTESCU N D. A Procedure to Determine Non-associated Constitutive Equations for Geo-materials[J]. International Journal of Plasticity，1994，10(2):103-131.

[5] 宋勇军，雷胜友，韩铁林.一种新的岩石非线性黏弹塑性流变模型[J].岩土力学，2012，33(7)： 2076-2080.(SONG Yong-jun, LEI Sheng-you, HAN Tie-lin. A New Nonlinear Viscoelasto-plastic Rheological Model for Rocks[J]. Rock and Soil Mechanics, 2012,33(7): 2076- 2080. (in Chinese))

[6] 杨圣奇，徐卫亚，杨松林.龙滩水电站泥板岩剪切流变力学特性研究[J].岩土力学，2007，28(5)：895-902.(YANG Sheng-qi, XU Wei-ya, YANG Song-lin. Investigation on Shear Rheological Mechanical Properties of Shale in Longtan Hydropower Project[J].Rock and Soil Mechanics, 2007, 28(5):895-902. (in Chinese))

[7] 徐卫亚，杨圣奇，谢守益，等.绿片岩三轴流变力学特性的研究(II)：模型分析[J].岩土力学，2005，26(5)：693-698.(XU Wei-ya, YANG Sheng-qi, XIE Shou-yi,etal. Investigation on Triaxial Rheological Mechanical Properties of Greenschist Specimen (II): Model Analysis[J]. Rock and Soil Mechanics, 2005, 26(5): 693-698. (in Chinese))

[8] 徐 辉，胡 斌，唐辉明，等.饱水砂岩的剪切流变特性试验及模型研究[J].岩石力学与工程学报，2010，29(增1):2776-2781.(XU Hui, HU Bin, TANG Hui-ming,etal. Experiment and Model Research on Shear Rheological properties of Saturated Sandstone[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(Sup.1): 2776-2781. (in Chinese))

[9] 胡 斌，蒋海飞，胡新丽，等.紫红色泥岩剪切流变力学特性分析[J].岩石力学与工程学报，2012，31(增1)：2796- 2802.(HU Bin, JIANG Hai-fei, HU Xin-li,etal. Analysis of Shear Rheological Mechanical Properties of Fuchsia Mudstone[J]. Chinese Journal of Rock Mechanics and Engineering, 2012, 31(Sup.1): 2796- 2802. (in Chinese))

[10] 邓荣贵，周德培，张倬元，等.一种新的岩石流变模型[J].岩石力学与工程学报，2001，20(6)：780-784.(DENG Rong-gui, ZHOU De-pei, ZHANG Zhuo-yuan,etal. A New Rheological Model for Rocks[J]. Chinese Journal of Rock Mechanics and Engineering, 2001，20(6)：780-784. (in Chinese))

[11] 杨圣奇，倪红梅，于世海.一种岩石非线性流变模型[J].河海大学学报(自然科学版)，2007，35(4)：388-392.(YANG Sheng-qi, NI Hong-mei, YU Shi-hai. A Kind of Nonlinear Rheological Model for Rock[J]. Journal of Hohai University(Natural Sciences), 2007,35(4):388-392. (in Chinese))

[12] 王芝银，李云鹏.岩体流变理论及其数值模拟[M].北京：科学出版社，2008.(WANG Zhi-yin, LI Yun-peng. Rock Rheological Theory and Numerical Simulation[M]. Beijing: Science Press, 2008. (in Chinese))