曾亚武，黎 玲，熊 俊，曹 源

(武汉大学土木建筑工程学院，武汉 430072)

1 研究背景

岩石是一种特殊的地质材料，其内部包含有大量的微裂纹。诸多实验研究表明，岩石内部微缺陷的萌生、发育、扩展及至连接成核的一系列过程形成了岩石在外荷载作用下的宏观应力－应变过程。当岩石材料变形到接近或超过峰值强度时，材料会出现局部的应变软化带，原来的均匀变形模式会被一种局限在某个很狭窄的带状区域内的严重不均匀的变形模式所代替，即产生很大的应变梯度，这种现象即所谓的应变局部化现象［1］。这种现象通常最先发生在材料内某些有限区域内，伴随着局部化带的发展而进一步扩展，它将直接导致材料的承载能力下降，也就必然导致岩土工程材料的承载力下降，故常常可以把它看作是材料局部破坏和结构失稳的一种先兆。

近年来，岩石应变局部化现象的研究已是国内外岩土界的研究热点之一。至今，国内外相关领域的研究者围绕应变局部化带问题展开了一系列工作，很多重要理论如局部化分叉理论、广义空隙压力理论、复合体理论、Cosserat理论、非局部理论、梯度塑性理论等在应变局部化研究中得到应用和发展。在这些理论方法中，梯度塑性理论引入非局部化的思想，将软化参数的梯度项引入到材料的屈服模式中，考虑了一点周围介质点的应力状态对该点的影响，解决了应变局部化分析中的网格依赖性问题。从梯度塑性理论的角度出发，对岩石应变局部化现象进行相关研究将有着重要的应用价值和工程背景;而且，基于梯度塑性理论的岩石应变局部化研究问题在学术上具有很大的挑战性，长远看来，它也必将推动弹塑性力学理论、有限元计算理论等的发展。

2 梯度塑性理论的分类

目前，对梯度塑性理论的研究可分为2类:应变梯度塑性模型与内变量梯度塑性模型。

2.1 应变梯度模型

在材料每一个质点上的变形不仅取决于通常意义的应变(与位移一阶梯度相关)，还取决于应变梯度(与位移二阶梯度相关)。一般情况下，这类模型也可以考虑二阶或高阶的应变梯度项。

2.2 内变量梯度模型

引进与某个热力学力相共轭的内变量梯度，热力学力只出现在内变量的演化方程中，在平衡方程中不出现。这类模型不改变平衡方程，只须修改本构方程，给数值计算带来了很大的便利。在弹性范围内，内变量不会变化而始终保持为初始值(通常为零)。因此，内变量梯度模型的初始响应由标准的弹性力学理论决定，而只有在非弹性阶段才会有梯度效应。

2.3 两者之间的差别

2类梯度塑性模型的根本区别在于，应变梯度模型中的应变梯度被看作是附加的可观察状态变量，这类变量与平衡方程中的高阶应力共轭;而内变量梯度模型中的内变量梯度与某些耗散热力学力共轭，它们作为内变量出现在演化方程中，但是不会出现在平衡方程中。因此，内变量梯度塑性理论仅仅修正塑性模型本构方程的基本形式，而运动学方程和平衡方程仍然保持不变。从热力学角度来讲，可以说内变量梯度塑性理论仅仅提高了潜在的热力势和耗散能，而应变梯度理论既需要对内力功也需要对外力功进行概括。

当然，某些模型也可能会带有混合的性质。例如，Zervos等［2］提出的梯度塑性模型，既可以被理解为应变梯度塑性模型，也可以被认为是包含内变量的二阶梯度的软化流动塑性模型。

3 应变局部化的梯度塑性理论研究

3.1 梯度塑性理论一般本构关系

经典的塑性理论认为，对于各向同性材料，材料的屈服函数可写为

式中κ为硬化参数(这里的硬化是包含软化在内的一个广义概念)，它可采用塑性功或等效塑性应变等内变量表示。梯度塑性理论在此基础上引入了非局部化思想，认为一点的屈服不仅与自身的应力状态和硬化参数相关，还与相邻点的硬化参数有关。它将硬化参数的梯度项引入到材料的屈服函数中，使一点的屈服极限还受相邻区域硬化参数的影响，这样其影响区域的大小将由所给定材料的内部特征长度决定。

唐天国等［3］基于经典塑性理论，结合梯度塑性理论的特点，以非局部硬化参数代替传统屈服函数中的局部硬化参数，推导出了梯度塑性理论的一般本构关系

3.2 内变量梯度塑性模型

由(2)式可知，在引入硬化参数的梯度项之后，一致性条件变为一个关于塑性乘子的偏微分方程，无法写出塑性乘子在局部意义上的显式表达式，进而也无法写出单元切线刚度矩阵的显式表达式。为了使计算结果接近真实的数值解，必须使用适当的正则化方法以防止应变变成一系列的零值。Lasry和Belytschko［4］在应变－位移关系中引入了位移的高阶梯度项;Fleck和Hutchinson［5］将高阶应变梯度作为附加的、可观测的状态变量;Aifantis［6］，Muhlhaus等［7］在非弹性本构方程中引入了描述非弹性变形的内变量的某种梯度，在材料模型中引入了作为正则化机制的内部长度参数。

考虑到岩石材料的特殊性，以塑性体积应变为内变量，在Drucker-Prager屈服模式下引入其梯度项，建立内变量梯度塑性模型。

3.2.1 塑性体积应变

按照经典弹塑性理论，在应力空间中岩石的屈服函数可以表示为

式中:I1是主应力第一不变量;J2是偏应力第二不变量;θσ是罗德角。按照复合求导法则可得

式中:δij是 Kronecker符号;sij为偏应力;tij是偏张量，且

其中:sim，smj也均为偏应力。

使用关联流动法则，可知塑性应变增量为

式中:dλ是1个非负的比例因子(表示塑性应变增量的大小)，上式括号中第1项反映塑性体积应变增量，而第2，3项反映塑性偏应变增量，即

对于岩石材料而言，屈服后的加载过程中，由于塑性变形不断增加，引起加载面不断改变，硬化程度取决于塑性加载历史，使用内变量来描述，于是，加载面的改变取决于应力状态和内变量。而在其变形过程中塑性体积变形的影响非常显著，应用塑性体积应变来定义内变量更利于反映体积应力与塑性体积应变的关系。

3.2.2 岩石应变局部化增量本构关系

材料内总应变与弹性应变和塑性应变之间存在以下关系:

在本模型中，选用考虑体积力对塑性变形产生影响的Drucker-Prager屈服函数:

由式(8)得到此时塑性偏应变增量为

为了简便，假定岩石在单轴压缩状态下，其本构关系为双线性，即峰值强度之前为线弹性，超过峰值强度后为线性应变软化，降模量为一个常值λ'，如图 1 所示。用，σ)分别代替单轴压缩下应力-塑性应变关系式中的塑性应变和应力，得到 σ()的表达式为

图1 单轴压缩线性软化模型Fig.1 Linear softening model in uniaxial compression

根据梯度塑性模型，在岩石各向同性假设下，可得如下关系:

由于l为材料的内部特征长度，于是式(10)变为

采用相关流动法则，有

而由式(11)有

结合式(15)、式(16)与式(17)，可得岩石应变局部化的增量本构关系为

式中:D为弹性矩阵。

3.3 平面应变下单向受力岩石应变局部化分析

3.3.1 应变局部化带宽度公式推导

假定岩石在单轴压缩状态下，压应力大小为σ0。由Drucker-Prager屈服函数得

设由上式得到的各应力分量流动方向分别为

考虑平面应变情况，将平面应变情况弹性矩阵D代入式(18)，得应力率为:

式中:G为剪切模量;μ为泊松比。

由一致性条件及屈服条件(15)并结合式(21)化简得

式中:

考虑到应变局部化发生过程中产生的显著体积变化，在单轴压缩平面应变情况下，考虑应变局部化区域内产生的体积扩容，设扩容系数为ξ，并有如下关系存在

将式(24)代入式(20)第二项，得到

由于平面应变条件下σ22的变化很小，可以认为，于是由式(25)可得

假设岩石试样在外荷载作用下形成的应变局部化带宽度为w，倾角为θ，非应变局部化带区域仍为弹性变形。引入一个垂直于应变局部化带的坐标系XoY，如图2所示。则该坐标系与现有xoy坐标系之间的转换关系为:

图2 坐标转换示意图Fig.2 Coordinate conversion

假设应变局部化带内参数仅沿带的Y方向发生变化，而在X方向是均匀不变化的，则:

将式(29)代入式(26)，得

将式(30)代入式(29)，得

将式(31)代入式(17)，并简化得

下面将讨论方程(32)的解，分以下2种情况讨论。

求解式(34)，得

由于b≠0，要使式(35)有解，则w=0，这与实际明显矛盾，故U＜0的情形下无解。

(2)当U＞0时，令U=u2，同上得方程的解为

于是有wu/2=π，由此可以得到应变局部化带的宽度表达式为

3.3.2 局部化带参数对局部化带宽度的影响

由式(38)可看出，岩石试件应变局部化带的宽度是随着应变局部化过程的发展不断变化的。又易知

于是由式(19)可得:

将式(40)分别代入式(23)与式(26)中，并化简可得:

于是式(32)中的U为

由式(41)及式(42)可看出:在U的表达式中，E，λ'，l，G，μ，α，ξ均为常数，θ为应变局部化带倾角。当U取最大值时，应变局部化带宽度w最小，应变局部化开始发生，此时的θ表示初始局部化带倾角。当U取最小值时，w最大，岩石已经发生破坏，此时的θ表示岩石破坏时对应的倾角。

4 结论

本文引入了梯度塑性理论的一般本构关系，将塑性变形梯度引入材料的本构方程中，说明梯度塑性理论中具有能够表征材料内部特征长度尺寸的长度量纲。考虑到岩石材料在受到外荷载作用时内部体积变形的特殊性，尝试采用以塑性体积应变作为内变量，建立内变量梯度塑性模型进行岩石应变局部化的研究，并从经典塑性理论下的屈服函数一般形式出发推导了塑性体积应变的表达式。在此基础上，还推导了岩石应变局部化增量本构关系，并在Drucker-Prager屈服函数下写出了塑性体积应变的具体形式，将其梯度项引入岩石材料的屈服关系中对其应变局部化问题进行研究。在平面应变条件下，分析了岩石试样受单轴压缩作用应变局部化的发生发展情况，得到其应变局部化带宽度的一般表达式，并对其进行分析讨论。

［1］王学滨，潘一山.地质灾害中的应变局部化现象［J］.地质灾害与环境保护，2001，12(4):1－5.(WANG Xue-bin，PAN Yi-shan.Strain Localization Phenomenon in Geological Disaster［J］.Journal of Geological Hazards and Environment Preservation，2001，12(4):1 － 5.(in Chinese))

［2］ZERVOS A，PAPANASTASIOU P，VARDOULAKIS I.A Finite Element Displacement Formulation for Gradient Elastoplasticity［J］.International Journal for Numerical Methods in Engineering，2001，50(6):1369 －1388.

［3］唐天国，谢新生，段云岭.高地应力场岩体深部卸荷断裂带的梯度项塑性算法分析［J］.水力学报，2011，12(42):1432 － 1437.(TANG Tian-guo，XIE Xin-sheng，DUAN Yun-ling.Analysis of Gradient-Plasticity Algorithm for Deep Unload Rock in High Ground Stress Field［J］.Journal of Hydraulic Engineering，2011，12(42):1432 －1437.(in Chinese))

［4］LASRY D，BELYTSCHKO T.Localization Limiters in Transient Problems［J］.International Journal of Solids and Structures，1988，24(6):581 －597.

［5］FLECK N A，HUTCHINSON J W.Strain Gradient Plasticity［J］.Advances in Applied Mechanics，1997，33:295－361.

［6］AIFANTIS E C.On the Microstructural Origin of Certain Inelastic Models［J］.Journal of Engineering Materials and Technology，1984，106(4):326 －330.

［7］MUHLHAUS H B，ALFANTIS E C.A Variational Principle for Gradient Plasticity［ J］.International Journal of Solids and Structures，1991，28(7):845－857.