第25卷第3期岩石力学与工程学报V ol.25 No.3 2006年3月Chinese Journal of Rock Mechanics and Engineering March，2006
岩石非线性黏弹塑性流变模型(河海模型)
及其应用
徐卫亚，杨圣奇，褚卫江
(河海大学 岩土工程研究所，江苏 南京 210098)
摘要：通过将作者提出的非线性黏塑性体与五元件线性黏弹性模型串联，建立一个新的岩石非线性黏弹塑性流
变模型(河海模型)，该流变模型可以充分反映岩石的加速流变特性。推导了岩石在恒应力与恒应变情况下的流
变方程，从理论上对岩石非线性蠕变和松弛特性进行分析，同时利用在岩石流变伺服仪上获得的蠕变全程试验
曲线，对提出的非线性黏弹塑性流变模型进行充分验证，试验曲线与非线性流变模型的比较，显示了所建流变
模型的正确性与合理性。基于提出的岩石非线性黏弹塑性流变模型，推导该非线性流变模型在三维条件下的中
心差分格式，采用FLAC3D所提供的二次开发程序接口，研制岩石非线性流变数值程序，同时通过一个单轴压
缩算例，进一步验证提出的非线性流变模型的正确性与合理性。采用研制的非线性流变数值程序，对锦屏一级
水电站坝基岩石工程进行三维流变数值模拟，分析结果为锦屏一级水电站坝基岩石工程的长期稳定与安全性提
供了合理建议与评价。
关键词：岩石力学；非线性黏弹塑性流变模型；蠕变特性；松弛特性；程序研制；工程应用
中图分类号：TU 485+.3 文献标识码：A 文章编号：1000–6915(2006)03–0433–15 NONLINEAR VISCOELASTO-PLASTIC RHEOLOGICAL MODEL (HOHAI MODEL) OF ROCK AND ITS ENGINEERING APPLICATION
XU Wei-ya，YANG Sheng-qi，CHU Wei-jiang
(Institute of Geotechnical Engineering，Hohai University，Nanjing，Jiangsu210098，China)
Abstract：By connecting the nonlinear viscoplastic body(NVPB) model put forward by the authors and five-component viscoplastic model in series，a new seven-component nonlinear viscoelasto-plastic rheological model of rock(Hohai model) is proposed，which can reflect fully the accelerative rheological property of rock material. The rheological equations of rock at the constant stress and strain are deduced. The nonlinear creep property and relaxation property of rock are carried out based on the theoretical analysis and investigation. At the same time，on the basis of complete creep curves gained by using rock servo-controlling rheology equipment，the identification of the proposed seven-component nonlinear viscoelasto-plastic rheological model of rock is validated. The comparison between nonlinear rheological model and experimental curve shows that the proposed nonlinear viscoelasto-plastic rheological model is available and reasonable. Based on the proposed seven-component nonlinear viscoelasto-plastic rheological model of rock，the central difference equation in three-dimensional condition of this nonlinear rheological model is deduced. By adopting the secondary development routine interface of software—fast Lagrangian analysis of continua (FLAC3D)，the numerical analysis code of nonlinear rheological model is investigated and developed. Through a uniaxial compression numerical
收稿日期：2005–01–15；修回日期：2005 - 04–05
基金项目：国家重点基础研究发展规划(973)项目(2002CB412707)；国家自然科学基金重点项目(50539110)
作者简介：徐卫亚(1962–)，男，博士，1982年毕业于华东水利学院工程地质专业，现任教授、博士生导师，主要从事岩石力学与工程方面的教学与
研究工作。E-mail：wyxu@https://www.doczj.com/doc/ee18262795.html,。
? 434 ? 岩石力学与工程学报 2006年
example，the correctness and rationality of the proposed nonlinear rheological model is validated further. In the end，by adopting the developed nonlinear rheological numerical code，the three-dimensional simulation is conducted for the dam foundation of Jinping First Stage Hydropower Station，and the analysis results bring forward reasonable suggestion and evaluation for the long-term stability and safety of dam foundation of Jinping First Stage Hydropower Station.
Key words：rock mechanics；nonlinear viscoelasto-plastic rheological model；creep property；relaxation property；code development；engineering application
1 引言
岩石流变模型是岩石流变力学理论研究中重要的组成部分[1～6]。岩石流变模型理论至今还不很成熟，许多重大岩石工程均为岩石流变模型理论的研究带来了严峻的挑战。当前岩石流变模型理论尤其是能反映岩石加速流变特性的模型理论仍是岩石流变力学研究中的热点和难点问题之一。然而由于岩石流变试验设备的限制，对岩石流变模型理论仍有许多值得探讨的地方。重大岩体工程建设均迫切需要了解岩石流变特性，以促进工程建设的顺利进行，并确保岩体工程在长期运营过程中的安全与稳定性，这就需要对岩石流变模型理论进行进一步的深入研究。
目前建立岩石非线性流变模型的方法主要有如下两种[7～12]：一是采用非线性流变元件代替常规的线性流变元件，建立能够描述岩石加速流变阶段的流变本构模型；二是采用内时理论、损伤断裂力学等新的理论，建立岩石流变本构模型。这两种方法建立的流变本构模型均能较好地描述岩石的加速流变阶段。本文采用第一种方法建立新的岩石非线性黏弹塑性流变模型。近年来，岩石力学中发展了一些非线性流变模型理论，较有代表性的有：韦立德等[12]根据岩石黏聚力在流变中的作用提出了一个新的SO非线性元件模型，建立了新的一维黏弹塑性本构模型；金丰年和浦奎英[13]基于试验结果，结合传统线性黏弹性模型的分析，提出了非线性黏弹性模型；邓荣贵等[14]根据岩石加速蠕变阶段的力学特性，提出了一种非牛顿流体黏滞阻尼元件，将该阻尼元件与描述岩石减速蠕变和等速蠕变特性的传统模型结合，构成了新的综合流变力学模型；曹树刚等[15]采用非牛顿体黏性元件构成五元件的改进西原正夫模型，探讨了与时间有关的软岩一维和三维本构方程和蠕变方程；陈沅江等[16]提出了蠕变体和裂隙塑性体两种非线性元件，并将它们和描述衰减蠕变特性的开尔文体及描述瞬时弹性的虎克体相结合，建立了一种可描述软岩的新的复合流变力学模型；张向东等[17]基于泥岩的三轴蠕变试验结果，建立了泥岩的非线性蠕变方程，并以此分析了围岩的应力场和位移场；王来贵等[18]以曹树刚等[15]改进的西原正夫模型为基础，利用岩石全程应力–应变曲线与蠕变方程中参数的对应关系，建立了参数非线性蠕变模型。由于岩石非线性流变元件模型有助于从概念上认识变形的弹性分量和塑性分量，且表达式通常能直接描述蠕变与松弛，所以许多岩石力学研究工作者用非线性流变元件模型来解释岩石的各种力学特性。当前岩石非线性流变元件模型仍是目前岩石流变力学理论研究中的一个重要课题。
鉴于此，本文提出一个新的岩石非线性黏弹塑性流变模型，并从理论上对其非线性蠕变与松弛特性进行分析，同时研制岩石非线性流变数值分析程序，以促进岩石非线性流变模型理论在工程实践中的应用。
2 岩石非线性黏弹塑性流变模型
笔者曾对绿片岩三轴流变试验结果进行分析研究[19，20]，提出了一个新的非线性黏塑性体(nonlinear viscoplastic body，NVPB)，该模型不仅能将西原模型中的塑性体和黏性体并联组合包括进去，而且能充分反映岩石加速流变特性，所以可以用于描述岩石流变的全过程曲线。岩石流变过程往往是弹性、黏性、塑性、黏弹性和黏塑性等多种变形共存的一个复杂过程，因而需要采用多种元件(线性和非线性元件)的复合来对其进行模拟。通过将提出的非线性黏塑性体(NVPB模型)与五元件线性黏弹性流变模型串联起来，可以建立一个新的岩石非线性黏弹塑性流变模型(如图1所示)，该流变模型可以充分反映岩石的加速流变特性。其流变模型满足如下条件：
当仅有1，2和3部分参与流变时，此时流变模型为五元件线性黏弹性模型，其相应的状态方程为
第25卷 第3期 徐卫亚等. 岩石非线性黏弹塑性流变模型(河海模型)及其应用 ? 435·
图1 岩石七元件非线性黏弹塑性流变特性分析示意图 Fig.1 Analysis sketch figure of seven-component nonlinear
viscoplasto-plastic rheological property of rock
?
?
?
?
?
????
++====+=+==32132132333212221
11εεεεσσσσεηεσεηεσεσ&&E E E (1) 但当1，2，3和4部分均参与流变时，此时流变模型为七元件非线性黏弹塑性模型，其相应的状态方程为
?
?
?
?
?
?
??
???+++=====+=+=+==?432143211
43S 432333212221
11)/(εεεεεσσσσσεησσεηεσεηεσεσn nt E E E &&& (2) 式(1)，(2)中：σ和ε分别为模型总的应力和应变；
1σ，2σ，3σ和4σ分别为1，2，3和4部分的应力；
1ε，2ε，3ε和4ε分别为1，2，3和4部分的应变；E 1，E 2，E 3，1η，2η和3η分别为材料的弹性、黏性
和塑性参数；n 为流变指数。
根据式(1)可得到五元件黏弹性流变模型的本构方程为 ????++++=+????????+
+21211
1212132231
21
ηηηηησεηηεηηεE E E E E E E &&&&& σηηση????????+++211313221213
)E E E E E E E E E & (3) 根据式(2)可得到七元件黏弹塑性流变模型的本构方程为 ????++++=+????????+
+21211
121213
2231
21ηηηηησεηηεηηεE E E E E E E &&&&& ++???
?????+++?????σησηηση&&31
211313221213n nt E E E E E E E E E +?????????++????3S 1231
23
S 2)
()
()1(ησσηηησσn n nt E E t n n ∫
??t nt E E n d )
(3
S 12
132ησσηη (4)
下面依据式(1)～(4)来推导岩石非线性黏弹塑性蠕变和松弛方程，并对岩石非线性黏弹塑性蠕变和松弛特性进行分析。
3 岩石非线性黏弹塑性流变特性
3.1 非线性黏弹塑性蠕变特性
在t =
0时加0σ≤S σ，并保持不变，将0σ代
入式(1)，对其进行Laplace 变换，并注意如下初始 条件：
??
??
???=====)
3 2 1(0)0()0(0)0()0(，，i i i i i σσ
εε
&&&&&& (5) 则
)
()
(~320
210
10
E S S E S S S
E ++
++
=ησησσε
(6)
式中：ε
~为ε的拉氏变换，S 为拉氏变换中的复变 量。对式(6)进行拉氏逆变换，可得
??
??
???+?????
?
?+
=
??E
E E E E 231
2e 1e
130201
ηησσσε (7) 式(7)即为五元件黏弹性流变模型的蠕变方程。两边对时间t 分别进行一次和二次求导，可得
t
E E 23
12e e 2
010ησησε
??+=& (8) E t E E E 23
12e e 2
2
032102ηηησησε????=&& (9) 由式(8)，(9)可知，ε
&＞0而ε&&＜0。此时模型有瞬时弹性变形及蠕变，且随着时间的不断增加，蠕变速率逐渐减小，故蠕变模型反映的是稳态蠕变，即应力0σ施加后，
模型产生瞬时弹性变形和衰减蠕变变形，最终趋于稳定，为稳态蠕变。
在=t 0时加0σ＞S σ，并保持不变，将0σ代入 式(2)，对其进行拉氏变换，并注意到式(5)的初始条件，可得
? 436 ? 岩石力学与工程学报 2006年
1
30320
210
10
!
)()
()
(~+?+
++
++
=n S S n E S S E S S S
E ησσησησσε
(10) 对式(10)进行拉氏逆变换，可得
n t E
t
E
t E E E 3S 030201
231
2e 1e
1ησσσσσε?+???
??????+??????
?
??+
=
?? (11)
式(11)即为七元件非线性黏弹塑性流变模型的蠕变方程。两边对时间t 分别进行一次和二次求导，可得
1
3
S 02010)(e e 2312????++=n t E
t E t n ησσησησε& (12) 23
S 02
2
32
1
2)
1()(e
e
2
31
2?????+
?
?=n t
E t
E
t n n E E ησσησησε
ηη&& (13)
可见当1=n 时，有ε
&＞0而ε&&＜0，此时模型只能反映岩石流变全过程曲线的初期流变和稳态流变 阶段，而不能反映岩石的加速流变阶段。而当n ＞1时，此时应变的一阶导数恒大于0，但应变的二 阶导数可以小于0、等于0或大于0，恰好对应着岩石流变全过程曲线的三阶段，如图2(a)所示。流变模型此时存在初期蠕变、稳态蠕变和加速蠕变，在 加速蠕变阶段，模型的蠕变随时间而迅速增长，蠕变增长的程度随着流变指数n 的变化而变化。模型 的蠕变特征曲线如图2(b )所示。 3.2 非线性黏弹塑性松弛特性
对于五元件黏弹性流变模型，在=t 0时加0ε，并保持不变，将0εε=代入式(1)，对其进行拉氏变换，并注意初始条件式(5)，可得
3
22110~~~E S E S E S ++++=ησησσε (14) 移项得
=σ
~ =
++++++++0211321322132211)]
()())([()
)((εηηηηηηE S E E S E E S E S S E S E S E =??
????+++????????++++???
???+????????++0212
1323123112122132122321εηηηηηηηηE E E E E E S E E E E S S E E S E E S E
))(())(())((210213
21210122312
101P S P S S E E E P S P S E E E P S P S S E +++++????????++
++εηηεηηε (15)
(a)
(b)
图2 岩石非线性黏弹塑性流变模型的蠕变特征曲线 Fig.2 Creep characteristic curves of seven-component
nonlinear viscoelasto-plastic rheological model of rock
其中，
??
???
??????
????
++=+++=??=?+=
0 0 2
42
4212
132312
311212221ηηηηE E E E E E B E E E E A B A A P B
A A P (16)
由式(16)中的后两式可知：1P ＞0且2P ＞0。 对式(15)中等号右边三项分别进行拉氏逆变换，有
01121
2210
11)e e (1))((L 12εεE P P P P P S P S S E t P t P ?????=?
?????++ (17)
＞ ＞
ε
，ε
.
ε
第25卷 第3期 徐卫亚等. 岩石非线性黏弹塑性流变模型(河海模型)及其应用 ? 437·
=????
?
???????++????????+?))((L 210122
311P S P S E E E εηη 01122
312)e e (121εηηE E E P P t P t
P ???????????+? (18) =????
????
???++?))((L 210213211
P S P S S E E E εηη 0213
2112122
112e 1e 111εηηE E E P P P P P P t P t P ????????????????+?? (19) 五元件黏弹性流变模型的松弛方程为 ??????
??+????????+?=
?t P P P P E E E P E E E P E 2e )(122213212112232
21ηηηησ 02121321121213
2111122321
11e )(εηηηηηη??
??
?+
?+????????+??P P E E E P P P E E E P E E E P E t P (20) 很显然，松弛曲线为下降的指数曲线。 当0=t 时，应力为
??????
??+????????+?=
)(122213212112232
210P P P E E E P E E E P E ηηηησ 02121321121213211
112232
11)(εηηηηηη?????
?+?+????????+?P P E E E P P P E E E P E E E P E (21)
当t →∞时，应力为
02
132313
2102121321lim εεηησσE E E E E E E E E P P E E E t ++===∞→∞(22)
对于七元件非线性黏弹塑性流变模型，在=t 0时加0ε，并保持不变，由于NVPB 模型无松弛特性，任何时刻应力都等于长期强度，因此岩石七元件非线性黏弹塑性流变模型的松弛特性求解与五元件线性黏弹性流变模型基本等同。经过推导发现，当
t →∞趋于无穷时，岩石七元件非线性黏弹塑性流变模型的应力松弛到岩石长期强度，即
S lim σσσ==∞
→∞t (23)
岩石七元件非线性黏弹塑性流变模型的松弛特征如图3所示，由图可见，松弛为不完全松弛，这符合岩石材料本身的特性。
图3 岩石非线性黏弹塑性流变模型的松弛特征曲线 Fig.3 Relaxation characteristic curve of nonlinear viscoelasto-plastic rheological model of rock
4 岩石非线性流变模型的试验验证
为了验证提出的岩石非线性黏弹塑性流变模型的正确性与合理性，利用岩石全自动流变伺服仪对绿片岩进行了三轴流变试验研究[19]，图4，5分别是绿片岩在围压15 MPa ，
应力水平为100 MPa 条件下的轴向与侧向流变全程试验曲线。由图可见，流变全程曲线均具有明显的加速流变阶段BC 。考 虑到该流变试验曲线是在等围压三轴压缩状态下得到的结果，下面首先根据前述推导得到的岩石一维
图4 岩石典型的轴向流变全程曲线 Fig.4 Typical axial rheology complete curve of rock
σσσ
轴向应变/ 10
－2
0.0
0.30.6
0.9 1.2 1.5
1.8
时间/ h
? 438 ? 岩石力学与工程学报 2006年
图5 岩石典型的侧向流变全程曲线
Fig.5 Typical lateral rheology complete curve of rock
非线性黏弹塑性蠕变方程，进一步推导岩石三维应力状态下的蠕变方程，然后提出岩石非线性流变模型参数的确定方法。
4.1 非线性黏弹塑性蠕变三维形式
通常岩石在三维应力状态下，其内部的应力张量可分解为球应力张量m σ和偏应力张量ij S ，其表达式分别如下：
??
?
?
?
??
?=?==++=kk ij ij ij ij ij kk S σδσσδσσσσσσ31
31
)(31
m 321m (24)
则有
m σδσij ij ij S += (25)
一般而言，球应力张量m σ只能改变物体体积，而不能改变其形状；而偏应力张量ij S 只能引起形状变化而不引起体积的变化。因此，同样可以将应变张量分解成球应变张量m ε和偏应变张量ij e ，其表达式分别如下：
??
?
?
?
??
?=?==++=kk ij ij ij ij ij kk e εδεεδεεεεεε31
31)(31m 321m (26)
则有
m εδεij ij ij e += (27)
对于五元件黏弹性流变模型而言，总应变是由弹性应变和黏弹性应变构成，即
321εεεε++= (28)
设岩石剪切模量为G ，体积模量为K ，则
)
1(2ν+=E
G (29)
)
21(3ν?=
E
K (30)
式中：E 和ν分别为岩石弹性模量和泊松比。
对于三维应力状态下的虎克体，则
m m 3εσK = (31)
ij ij Ge S 2= (32)
这样对于三维应力状态下的黏弹性体，其蠕变方程可比照单轴应力状态下的蠕变方程，只要将0σ换成0)(ij S (0)(ij S 为岩石试验时恒定偏应力)，即
???
?????
?+???
?????
?=
??t E
ij t E
ij ij E S E S e 2312e 12)(e 12)(3020ηη (33) 而非线性黏塑性体(NVPB)的蠕变方程为
?
?
?
???==) )(( )() )(( 0S 03S 0S 0σησσij n
ij ij ij ij S t S e S e (34) 根据叠加原理，结合前面各式，便可以得到三维应力状态下岩石的蠕变方程：
????????????????+???
??????+?????
????+
=?????
????+?????????+
=????)
)(( )( e 12)(e 12)(2) )((e 12)( e 12)(2S 03S 030201m S 030201m
2312
23
12σησσσσηηη
ηij n
ij t E ij t E ij ij ij t E ij t E
ij ij S t S E S E S E e S E S E S E e (35)
据此可以得到等围压三轴压缩时岩石在恒定应
力条件下的蠕变方程，即 ?
?
?????
????????
???+??????????+???
??
?????+?+?=????
???????+???
???????+
?+?=????) (3e 13e 1339)(2) (e 13e 1339)(2S
313S 3133123113131S 313312
3
11313123
122312σσσησσσσσσσσσσσεσσσσσσσσσσσεηη
η
ηn t G G
t G
G
t G G G K G G G K (36)
4.2 非线性黏弹塑性流变模型参数确定及验证
确定岩石流变模型参数主要方法有：回归反演法、最小二乘法以及流变曲线分解法[21]等，鉴于回
侧向应变/ 10－
2
0.0 0.3
0.6
0.9 1.2 1.5
1.8
时间/ h
≤
≤ ＞ ＞ ＞≤
第25卷第3期徐卫亚等. 岩石非线性黏弹塑性流变模型(河海模型)及其应用? 439·
归反演法与最小二乘法对于非线性流变模型效果不
是很理想，本文采用流变曲线分解法，提出岩石非
线性流变模型参数确定方法。设a1=9G1K/(3K+2G1)，
下面讨论轴向偏应力差大于岩石长期强度的岩石非
线性黏弹塑性流变模型参数的计算方法。对图4和
5所示的黏弹塑性流变全程曲线而言，首先利用岩
石在(0，t2)(t2指从稳态流变向加速流变过渡的起始
时刻)时段内的流变试验数据，采用直接迭代法，对
式(36)中的五元件黏弹性模型进行求算，这样可求
得5个流变参数，即a1，G2，G3，1η和2η，具体求
算方法如下：根据试验得到的n组数据(ε，t)，首
先给定一组初始近似值(01a，02G，03G，01η，02η)，
以此为基础进行迭代，这样反复计算，直到满足所
需要的精度，则此时对应的流变参数值即为所求。
然后利用得到的五元件黏弹性流变模型，求出在(t2，
t c)时段内岩石应变的黏弹性理论解，对图4所示的
岩石轴向流变全程曲线而言，t2 =1.2 h，而对图5所
示的岩石侧向流变全程曲线而言，t2 =1.31 h；t c指
流变破坏时间，为1.66 h，通过找出在(t2，t c)时段
内的试验值与黏弹性理论解差值的m对应变数据，
再通过最小二乘法对m对应变数据与相对应的时间
t进行幂律回归分析，即可求得NVPB模型的流变
参数
3
η和n。至此，7个岩石流变参数已全部得到，
代入到式(36)中的第二项，即可得到常规三轴压缩
状态下的岩石七元件非线性黏弹塑性流变模型。采
用上述方法求解得到的以轴向应变和侧向应变辨识
的岩石七元件非线性黏弹塑性模型的流变参数，分
别如表1和2所示。根据辨识得到的a1值以及
表1 以岩石轴向应变辨识的非线性黏弹塑性流变参数
Table 1 Nonlinear viscoelasto-plastic rheological parameters
of rock in accordance with axial strain curve
a1 /GPa
G2
/GPa
G3
/GPa
η1
/(GPa·d)
η2
/(GPa·d)
η3
/(GPa·d)
n
25.8 238.4 32.3 0.87 11.0 694.4 12.673
表2 以岩石侧向应变辨识的非线性黏弹塑性流变参数Table 2 Nonlinear viscoelasto-plastic rheological
parameters of rock in accordance with lateral
strain curve
a1 /GPa
G2
/GPa
G3
/GPa
η1
/(GPa·d)
η2
/(GPa·d)
η3
/(GPa·d)
n
64.4 96.8 96.8 7.43 7.2 578.7 15.990式(29)与(30)，通过试验首先确定出岩石的泊松比ν值，然后就可以确定出K和G1的值。若=
ν0.4，则以轴向应变辨识的K和G1的值分别为45.87和9.83 GPa，而以侧向应变辨识的K和G1的值分别为114.48和24.53 GPa。
图6给出了得到的出现加速流变的全过程试验结果与提出的岩石非线性黏弹塑性流变模型的对比。由图可见，建立的岩石非线性黏弹塑性流变模型与轴向以及侧向流变试验结果吻合得相当理想，这表明本文所建岩石七元件非线性黏弹塑性流变模型的正确性与合理性。
图6 岩石非线性黏弹塑性流变模型计算结果与试验结果的比较
Fig.6 Comparisons between results of the proposed rheology model of rock specimen and experimental result
5 岩石非线性流变数值程序研制
近年来随着重大岩石工程的大量兴建，迫切需轴
向
应
变/
1
－
2
0.38
0.39
0.40
0.41
0.42
0.43
0.44
时间/ h
0.00.30.60.9 1.2 1.5 1.8
非线性黏弹塑性流变模型
轴向应变试验结果
侧
向
应
变/
1
－
2
0.15
0.19
0.23
0.27
0.31
0.35
时间/ h
0.00.30.60.9 1.2 1.5 1.8
侧向应变试验结果
非线性黏弹塑性流变模型
? 440 ? 岩石力学与工程学报 2006年
要研制能够真实反映岩石流变特性的非线性流变数值程序。本文采用FLAC 3D 所提供的二次开发程序接口，研制岩石非线性流变数值程序，这样不仅能保证计算效率，而且程序编写代码相对较小。 5.1 非线性流变数值程序研制的基本原理
前文已经给出了提出的岩石非线性黏弹塑性流变模型的一维形式与非线性蠕变的三维形式，在岩石非线性流变程序研制过程中，需给出该非线性流变模型在三维条件下的中心差分形式，这里采用应变偏量与球应变张量形式来处理。
根据图1可知，岩石总应变有如下形式：
4
321ij ij ij ij ij εεεεε+++= (37)
将式(37)转成偏量速率可得
4
321ij ij ij ij ij e e e e e &&&&&+++= (38) 其中，
ij x kk
x ij
x
ij
e δεε3
?
= (=x 1，2，3，4) (39) 对于1部分来说，有
1
12ij ij e G S = (40)
对于2部分来说，有
222
122ij
ij
ij e G e
S +=&η (41) 对于3部分来说，有
333222ij ij ij e G e S +=&η (42)
式(37)～(42)中：G i (=i 1，2，3)为剪切模量，可由 式(29)确定。
由作者以前的研究[20]可知，一维条件下NVPB 模型的应变速率形式为
η
σσε
)(S 01?=?H nt n & (43) 将式(43)转成三维条件下，则应变速率为
ij
n ij
g
F H nt σηε
??=?314
2)(& (44) 式中：F 为屈服函数，g 为塑性势函数。
这样NVPB 模型，即图1中的4部分，应变偏量速率可以写为
ij kk ij n ij
e
g F H nt e δση32)(4314
&&?
??=? (45) 其中，
???
???????+??+??=?3322113142)(σσσηg g g
F H nt e n kk
& (46) 通常在塑性力学中一般假定塑性变形不引起体积改变，因而体积应力只由弹性体积应变产生，则
3
32211kk kk kk m e K e
K e K &&&&++=σ (47) 式中：K i (=i 1，2，3)为岩石体积模量，可由式(30)确定。
将式(38)写成增量形式，有
4
321ΔΔΔΔΔij ij ij ij ij e e e e e +++= (48)
式(40)的差分形式为
1
1
2ΔΔG S e ij
ij =
(49) 将式(41)写成中心差分形式，有
t e G e t S ij ij ij Δ22Δ222
1+=η (50)
将式(42)写成中心差分形式，有
t e G e t S ij ij ij Δ22Δ3332+=η (51)
其中，
2O
N ij
ij ij S S S +=
(52) 2
O N ij
ij ij e e e +=
(53) 式(52)，(53)中：上标“N ”表示一个流变计算时间
步新的量值，而“O ”表示计算时间步中旧的量值。
将式(52)，(53)代入式(49)可得
??
????++=)(4Δ1O N 1O 2N
2ij
ij ij ij
S S t De C e η，， (54) 其中，
t G C Δ211
2
η+= (55) t G D Δ211
2
η?
= (56) 同理，将式(52)，(53)代入式(50)可得
??????+Δ+=)(41O N 2O 3N
3ij ij ij
ij
S S t Ne M
e η，， (57) 其中，
t G M Δ2123
η+
= (58) t G
N Δ212
3η?= (59)
类似可以得到2部分Kelvin 体新的球应变为
??
????++=)(6Δ1O m N m 2O 2m N
2m
σσK t Je I e ，， (60) 其中，
t K I Δ211
2
η+= (61) t K J Δ?
=1
2
21η (62)
第25卷 第3期 徐卫亚等. 岩石非线性黏弹塑性流变模型(河海模型)及其应用 ? 441·
3部分Kelvin 体新的球应变为 ??
????+Δ+=)(61O m N m 3O 3m N
3m
σσK t Te R e ，， (63)
其中，
t K R Δ+=2
3
21η (64) t K T Δ?
=2
3
21η (65) 将式(51)，(54)与式(57)代入式(47)中，可得
=N ij S
??????????????????????+Δ?ΔO 2O 2O
4111，，ij ij ij ij ij e M N e D C bS e e a (66) 其中，
???
?????Δ+Δ+=
2114421
ηηC t C t G a (67) ???????
?Δ+Δ?=2114421
ηηC t C t G b (68) 对于整个非线性黏弹塑性流变模型，球应力张量为
3
32211O m N m kk kk kk e K e K e K Δ+Δ+Δ+=σσ (69)
式(66)，(69)是编写代码的关键所在，采用这两个式子，可以在一个计算时间步中根据前一步的一些量值得到新的应力值。
5.2 非线性流变数值程序研制的关键技术
FLAC 3D (Version 2.1)是由面向对象的语言C++编写而成的，因而程序接口表现出高效性、稳定性和简洁性。软件中所有的本构模型(包括流变模型)都是以动态链接库文件(后缀名为.dll)的形式提供给用户，同时软件支持用户编写自定义的本构模型，即通过在VC6.0的开发环境下将自定义的本构模型代码编译成动态链接库文件，来实现将自定义的本构模型嵌入到FLAC
3D
软件中。在FLAC
3D
求解过
程中，主程序自动识别并加载用户指定的本构模型的动态链接库文件。FLAC
3D
采用面向对象技术将
自定义本构模型的全部信息封装在抽象类class
constitutive model 中，二次开发过程实质上就是由 抽象类派生出用户所需具体本构类的过程。
基于节4.1推导得到的非线性流变模型在三维条件下的中心差分形式，建立的岩石非线性流变数值程序研制流程图如图7所示。
图7 岩石非线性流变数值程序研制流程图 Fig.7 Flow chart of rock nonlinear rheology numerical code
在编译动态链接库文件时需要用到stensor.h ，axes.h 和conmodel.h 共3个头文件，其中stensor.h 文件中有如下两函数在进入屈服阶段后会经常用到，需引起重视：函数EXPORT bool Resoltopris(double *
pdMin ，double *pdMid ，double *pdMax ，Axes* aDir ，unsigned dim ，int* icase ，double *sdif ，double *psdif)的功能是求出张量的主应力和主轴；函数EXPORT
void Resoltoglob(const double & pdMin ，const double & pdMid ，const double & pdMax ，Axes& aDir ，const double & dMin ，const double & dMid ，const double & dMax ，unsigned dim ，int* icase ，double *sdif ，double * psdif ，bool fast)的功能是由张量的主值和主轴求出全局坐标(在实际问题中是计算网格模型的整体坐标系)下的张量。屈服函数往往写成主应力或是主 应力不变量的形式，因而经常需要用第一个函数求出主应力，然后由流动法则和相关方程得到新的主应力后，再用第二个函数将新的主应力转换成全局
? 442 ? 岩石力学与工程学报 2006年
坐标下的应力张量。
Conmodel.h文件中的抽象类ConstitutiveModel 是程序研制的关键所在，研制本构模型的过程也就是继承ConstitutiveModel类，并重载其中几个关键函数的过程，其中最关键的两个函数是Virtual const char *Initialize(unsigned uDim，State *pst)和Virtual const char *Run(unsigned uDim，Staste *pst)。Initialize函数通常用于检查FISH语言输入的材料参数是否正确，为了提高程序执行效率，可以对此函数计算过程中用到的一些中间参数进行初始化。对于提出的岩石非线性流变模型，9个材料参数G1，G2，G3，K1，K2，K3，1η，2η与3η均不能小于0，此外由于采用的是Mohr-Couloumb屈服准则，所以必须在Initialize函数中检验摩擦角、黏聚力、膨胀角和抗拉强度是否大于0。FLAC3D在求解时会在每一个计算时间步内对每一个单元的子单元调用Run 函数。本构方程是通过重载ConstitutiveModel函数来实现的，具体地说就是根据子单元State数据类型所提供的量值来得到新的应力值。
在派生ConstitutiveModel类时，本文岩石非线性流变数值程序研制有如下几点需引起注意：
(1) 在派生ConstitutiveModel类时，除了申明几个必要的计算参数变量外，还需要额外地申明14个双精度浮点数类型的变量，这14个变量分别用来存放一个计算时间步内2和3部分Kelvin体新的偏应变和球应变。
(2) 在FLAC3D计算过程中常常需要暂停求解，为了使中断的求解继续进行，需要重载EXPORT const char *SaveRestore(ModelSaveObject *mso)这个虚函数。当FISH命令SA VE或是RESTORE时，此函数将被调用，函数主要功能是保存或取出计算过程中所需要的参数值。
(3) 进入塑性变形后，应力、应变呈非线性关系，由应变增量值来求应力值通常需要迭代求解，此时提高程序执行效率显得尤为重要。为了简化一个时间步内的计算量，在重载Virtual const char * Run (unsigned uDim，State *pst)函数时，应尽量避免大规模的循环语句。
5.3算例考证
通过笔者研制的非线性流变数值程序，基于一个单轴压缩算例，对提出的岩石非线性流变模型进行了考证。试件尺寸为10 cm (X向)×20 cm(Y向)×1 cm (Z向)，共划分200个单元，462个节点。模型底部在Y方向约束，顶部施加一个100 MPa的分布压力，如图8所示。基于研制的非线性流变数值程序，分别采用广义Kelvin模型、西原模型以及七元件非线性黏弹塑性流变模型进行非线性流变数值分析。同时与FLAC3D中自带的广义Kelvin模型(三元件模型)(即图1中的1与2部分串联组合)进行比较，以验证非线性流变数值程序研制的正确性与合理性。
图8 单轴压缩算例示意图
Fig.8 Sketch figure of uniaxial compression example
表3分别给出了由研制的非线性流变数值程序退化的广义Kelvin模型与FLAC3D中自带的广义Kelvin模型计算得到的节点A(如图8所示)竖直向(Y 向)位移的比较，数值算例采用的计算参数如下：K1=30 GPa，G1=40 GPa，G2=50 GPa，=
1
η100
表3 研制程序与FLAC3D计算得到的竖直向位移比较Table 3 Displacement comparison gained by developed code FLAC3D
节点位移/(10－4m) 节点位移/(10－4m) 时间
/d 研制程序FLAC3D
时间
/d 研制程序FLAC3D
0 2.406 6 2.406 660 3.730 1 3.730 1
1 2.531 1 2.531 180 3.730 5 3.730 5
2 2.82
3 1 2.823 1100 3.730 9 3.730 9
3 3.210 8 3.210 8150 3.731 8 3.731 8
4 3.576 4 3.576 4200 3.732 6 3.732 6
5 3.623 5 3.623 5250 3.733 4 3.733 4
10 3.720 2 3.720 2300 3.734 2 3.734 2
20 3.728 4 3.728 4350 3.735 0 3.735 0
40 3.729 6 3.729 6400 3.735 9 3.735 9
X
Z
1)
P
X
第25卷 第3期 徐卫亚等. 岩石非线性黏弹塑性流变模型(河海模型)及其应用 ? 443·
GPa ·d 。由表可见，采用广义Kelvin 模型，分别基于研制程序与FLAC
3D
计算得到的节点流变位移完
全吻合。图9给出了由研制程序计算得到的节点A 竖直向位移(以向下的位移为正)流变曲线。由图可见，随着时间的增加，节点位移逐渐增加，但当时间增加至10 d 以后时，便基本趋于稳定，不再明显增加，这符合广义Kelvin 流变模型自身的特性，也充分表明了研制的非线性流变数值程序黏弹部分的正确性与合理性。
图9 由研制程序计算得到的节点A 竖直向位移曲线 Fig.9 Displacement curve of node A in vertical direction
gained by generalized Kelvin model of code
为了验证研制的非线性流变数值程序黏塑部分的正确性，以岩石七元件非线性黏弹塑性流变模型退化为西原模型的数值程序，分别依次对模型施加
30，40，50，60，70与80 MPa 六级荷载，采用的数值计算参数如下：K 1=30 GPa ，G 1=40 GPa ，
G 2=50 GPa ，=1η100 GPa ·d ，η2=1 GPa ·d ，=c 15 MPa 和=?45°。图10给出了计算得到的节点A 竖直向位移(以向下的位移为正)流变曲线，曲线附近的数字为施加的应力水平值。由图10可知，低应力水平时，节点位移表现为典型的黏弹性质，但当应力水平高于屈服应力后，节点位移出现加速流变，表现为塑性特征，数值结果与笔者所得流变试验结果[19]在定性上相一致。
基于七元件非线性黏弹塑性流变模型，采用研制的非线性流变数值程序，分别依次对模型施加
30，40，50，60，70与80 MPa 六级荷载，通过对模型进行流变数值分析，由此得到的节点A 竖直向位移流变曲线如图11所示。
采用的计算参数如下所示：
K 1=30 GPa ，G 1=40 GPa ，G 2=50 GPa ，G 3=60 GPa ，η1=100 GPa ·d ，η2=150 GPa ·d ，η3=200
图10 基于西原模型得到的节点A 竖直向位移曲线 Fig.10 Displacement curve of node A in vertical direction
gained by Kormanura model
图11 基于研制的非线性流变数值程序计算得到的节点A
竖直向位移曲线
Fig.11 Vertical displacement of node A gained by developed
nonlinear rheological numerical code
GPa ·d ，=c 15 MPa ，=?45°以及不同的n 值。由图可见，低应力水平时，节点位移亦表现为典型的黏弹性质，但当应力水平达到80 MPa 时，节点位移表现为塑性特征，而且流变指数对岩石的黏弹特性几乎没有影响。
为了清楚地看出流变指数对岩石塑性特征的影响规律，对80 MPa 应力水平下的流变曲线进行了放大处理，如图12所示。由图12可见，随着流变指数的逐渐增加，岩石流变全程曲线加速段特征越来越明显，这充分反映了岩石的非线性流变特性，数值分析结果与笔者所得到的流变理论研究结果相吻合[20]。通过上述分析可以看出，基于本文提出的七元件非线性黏弹塑性流变模型，研制的非线性流
0.6
1.11.6
2.12.6
3.13.6竖向位移/(10－
4 m )
050
100 150 200
250
时间/d
50MPa
60 MPa
70 MPa 80MPa
40MPa
30MPa
2.4
2.8
3.23.6
4.0－50
50150250350450
时间/d
竖向位移/(10－
4 m )
0.0
1.0
2.0
3.0
4.05.0
竖向位移/(10－
4 m )
100
200300400500
时间/d
? 444 ? 岩石力学与工程学报 2006年
图12 NVPB 模型流变指数对流变曲线的影响 Fig.12 Effect of rheological index of NVPB model on
rheological curve of specimen
变数值程序能够合理反映岩石材料的非线性黏弹塑性流变特性。
6 工程应用
基于研制的岩石非线性流变数值程序，对锦屏一级水电站坝基岩石工程进行了三维流变数值模拟分析，力图为锦屏一级水电站坝基岩石工程的长期稳定与安全性提供合理建议与评价。 6.1 工程概况
锦屏一级水电站位于雅砻江中上游普斯罗沟与手爬沟之间，水能经济指标十分优越。大坝为混凝土双曲拱坝，坝高305 m ，它是目前拱坝坝高居世界第一的特大型工程。该水电站以发电为主，兼有蓄能、蓄洪和拦沙等作用。水库总库容77.6亿m 3
，调节库容为49.1亿m 3
，装机容量3 600 MW ，年发电量166.2亿kW ·h 。坝基地处青藏高原向四川盆 地过渡的斜坡地带，地貌区为强烈切割的中～高山区。普斯罗沟坝段出露的地层主要为三叠系中上统杂谷脑组及第四系冲洪积、崩坡积和残积物。坝基岩体出露的地层主要是大理岩和绿片岩，由于经受了多期多方位的强烈构造作用和区域变质作用，岩体成分及结构十分复杂，既存在原生的层面和节理面，又有经后期改造形成的构造裂隙、断层错动 带、挤压带及深拉裂缝，同时受风化、卸荷等表生改造作用影响，从而使得局部岩体变得松散，岩体质量下降，力学强度降低。 6.2 数值模型
三维数值模型的研究范围确定为：顺河向取
1 000 m ，约3.3倍坝高，其中向上游延伸约1.5倍坝高，向下游延伸约1.8倍坝高；横河向取1 300 m ，约为4.1倍坝高；模型最大高程为
2 400 m(高程从
1 300 m 起)，约为3.6倍坝高。
根据地质剖面图，自下游向上游依次绘制了8个剖面，按坝基范围较密、上下游延伸范围内较稀的方式布置。根据拱坝各控制高程的水平剖面建立坝体，然后利用基岩和坝体的各层剖面，建立起拱坝、拱座和高边坡的三维联合数值模型，基岩内的
断层采用实体单元模拟，共生成29 972个单元，
分为变质砂板岩、粗晶大理岩、细晶大理岩、绿片岩、断层以及坝体6种材料。数值模型上、下游边界顺河向约束，左、右岸边界横河向约束，基础底部沿铅垂向约束。考虑自重以及坝体所受到的水荷载。
6.3 弹塑性计算参数
为了与流变计算结果进行比较，利用建立的三维拱坝和高边坡联合数值模型，进行了弹塑性数值计算。Drucker-Prager 强度准则简单实用，一直是目前应用最广的岩土屈服准则，其形式为
k I J F ?+=12 α (70)
式中：I 1和J 2分别为第一应力不变量和应力偏量第二不变量；α和k 为与材料的黏聚力c 和内摩擦角
?有关的参数，其表达式分别为
)
sin 3(3sin 6??α?=
(71)
3.23.43.63.8
4.04.24.44.6400.0 400.4 400.8 401.2
401.6
时间/d
节点位移/(10－
4 m )
第25卷 第3期 徐卫亚等. 岩石非线性黏弹塑性流变模型(河海模型)及其应用 ? 445·
)
sin 3(3cos 6???=
c k (72)
计算时坝基岩体采用D-P 强度准则，而坝体采用修正的Mohr-Coulomb 强度准则，其形式如下所示：
??σσN C N F 231+?= (73)
式中：1σ和3σ分别为第一与第三主应力；?N 为参数，与内摩擦角?有关，其形式为
?
?
?sin 1sin 1?+=
N (74)
弹塑性计算参数如表4
[22]
所示。
表4 弹塑性计算参数[22]
Table 4 Elastoplastic parameters of numerical simulation [22]
岩石 E /GPa
泊松比ν ρ/(kg ·m －
3)
c /MPa ? /(°)
变质砂板岩 28 0.200 2 750 1.00 46.1 粗晶大理岩 25 0.250 2 750 1.05 51.3 细晶大理岩 28 0.250 2 800 1.34 63.4 绿片岩 12 0.250 2 850 0.70 37.2 断层 2 0.275 2 600 0.02 19.0 拱坝
24 0.167 2 400 5.00 59.5
6.4 流变计算参数
流变计算时采用的是提出的七元件非线性黏弹塑性流变模型，根据实际工程地质条件、坝基岩体工程地质模型以及实际工程背景，在试验、监测反馈与综合类比的基础上，采用的流变计算参数如 表5[22]所示。
表5 流变计算参数
[22]
Table 5 Rheological parameters of numerical simulation
[22]
岩石 E 1 /GPa
E 2
/GPa η1 /(GPa ·
d)E 3 /GPa η2 /(GPa ·d)η3 /(GPa ·d)n
变质砂板岩 28 600 600 600 600 100
0.8粗晶大理岩 25 900 300 900 300 100 0.8细晶大理岩 28 800 600 800 600 100 0.8绿片岩 12 400 500 400 500 150 0.6断层 2 20 40 20 40 200 0.4拱坝
24 500 300 500 300 200
0.8
6.5 计算结果分析
以高程1 710 m 处坝肩左岸岩体关键点2454为例，对比坝肩岩体弹塑性与流变计算结果。计算表明，坝肩岩体该关键点处弹塑性顺河向位移为4.92
mm ，而考虑流变的黏弹塑性顺河向位移为 5.27 mm ，增加了7.11%。弹塑性横河向位移为1.01 mm ，而考虑流变的黏弹塑性横河向位移为1.65 mm ，增加了63.4%。弹塑性竖直向位移为7.01 mm ，而考虑流变的黏弹塑性竖直向位移(0=t )为7.68 mm ，增加了9.56%。这表明不考虑流变特性影响计算得到的弹塑性数值结果偏于保守，会造成工程某种程度上的安全隐患，对此应该有足够清楚的认识。
为了清楚地看出正常蓄水位工况下坝肩岩体的流变规律，进一步分析高程1 710 m 处坝肩左岸岩体关键点2454的流变变形。正常蓄水位工况下坝肩左岸岩体关键点的竖直向位移(以向下方向为正)随时间的变化曲线分别如表6与图14所示。由图可 见，在开始发生流变的0～120 d 内，坝肩岩体关键点位移随着时间逐渐增加，但是在达到120 d 以后，随着时间的继续增长，关键点的位移便趋于稳定，不再有明显变化。因此120 d 可以看作是锦屏一级水电站坝肩岩体流变趋向于稳定的临界值。
表6 关键节点岩体竖直向位移
Table 6 Vertical displacements of key points in rockmass
时间/d
位移/mm
时间/d
位移/mm
0 7.679 8 180 14.279 5 7.620 1 240 14.279 10 8.053 5 300 14.279 20 11.691 0 360 14.279 40 13.252 0 480 14.279 80 14.146 0 600 14.280
120 14.279 0 720 14.280
图14 关键节点岩体流变曲线
Fig.14 Rheological curves of key points in rock mass
时间/d
节点位移/m m
6810
12
14
160
200
400600800
0.00.10.2
0.3
0.4
节点位移速率/(m m ·d －
1)
? 446 ? 岩石力学与工程学报 2006年
7 结论
(1) 基于提出的岩石非线性黏弹塑性流变模型(河海模型)，推导了岩石在恒应力与恒应变情况下的流变方程，从理论上对岩石的非线性蠕变特性和松弛特性进行了分析。
(2) 西原模型是河海模型在流变指数为1时的一个特例，提出的河海模型可以很好地描述岩石的非线性加速蠕变变形。
(3) 利用流变全程试验曲线，对提出的岩石非线性黏弹塑性流变模型(河海模型)进行了成功辨识，试验结果与非线性流变模型的比较，显示了所建的河海模型的正确性与合理性。
(4) 推导了岩石非线性黏弹塑性流变模型(河海模型)在三维条件下的中心差分格式，采用FLAC3D 所提供的二次开发程序接口，研制了岩石非线性流变数值程序。
(5) 采用研制的非线性流变数值程序，通过一个单轴压缩算例，并结合锦屏一级水电站坝基岩石工程三维流变数值分析，验证了提出的河海模型与研制的非线性流变数值程序的正确性与合理性。
参考文献(References)：
[1] Li Y S，Xia C C. Time-dependent tests on intact rocks in uniaxial
compression[J]. Int. J. Rock Mech. Min. Sci.，2000，37(3)：467–475.
[2] Yang C H，Daemen J J K，Yin J H. Experimental investigation of creep
behavior of salt rock[J]. Int. J. Rock Mech. Min. Sci.，1999，36(2)：
233–242.
[3] Maranini E，Yamaguchi T. A non-associated viscoplastic model for the
behaviour of granite in triaxial compression[J]. Mechanics of Materials，2001，33(5)：283–293.
[4] Maranini E，Brignoli M. Creep behaviour of a weak rock：
experimental characterization[J]. Int. J. Rock Mech. Min. Sci.，1999，36(1)：127–138.
[5] Shao J F，Zhu Q Z，Su K. Modeling of creep in rock materials in terms
of material degradation[J]. Computers and Geotechnics，2003，30(7)：549–555.
[6] 朱珍德，徐卫亚. 岩体黏弹性本构模型辨识及其工程应用[J]. 岩石
力学与工程学报，2002，21(11)：1 605–1 609.(Zhu Zhende，Xu Weiya. Visco-elastic constitutive model recognition of rock mass and its engineering application[J]. Chinese Journal of Rock Mechanics and Engineering，2002，21(11)：1 605–1 609.(in Chinese))
[7] 陈沅江，潘长良，曹平，等. 基于内时理论的软岩流变本构模
型[J]. 中国有色金属学报，2003，13(3)：735–742.(Chen Yuanjiang，Pan Changliang，Cao Ping，et al. Endochronic rheological constitutive
model of soft rock[J]. Chinese Journal of Nonferrous Metals，2003，13(3)：735–742.(in Chinese))
[8] 金丰年，范华林. 岩石的非线性流变损伤模型及其应用研究[J]. 解
放军理工大学学报，2000，1(3)：1–5.(Jin Fengnian，Fan Hualin.
Study on nonlinear rheology damage property of rock[J]. Journal of PLA University of Science and Technology，2000，1(3)：1–5.(in Chinese))
[9] 曹树刚，鲜学福. 煤岩蠕变损伤特性的试验研究[J]. 岩石力学与工
程学报，2001，20(6)：817–821.(Cao Shugang，Xian Xuefu. Testing study on the characteristics of creep and damage of coal and other rocks[J]. Chinese Journal of Rock Mechanics and Engineering，2001，20(6)：817–821.(in Chinese))
[10] 孙钧. 岩土材料流变及其工程应用[M]. 北京：中国建筑工业出
版社，1999.(Sun Jun. Rheology of Geomaterial and Its Engineering Application[M]. Beijing：China Architecture and Building Press，1999.(in Chinese))
[11] 巫德斌，徐卫亚，朱珍德，等. 泥板岩流变试验与黏弹性本构模型
研究[J]. 岩石力学与工程学报，2004，23(8)：1 242–1 246.(Wu Debin，Xu Weiya，Zhu Zhende，et al. Study on rheological test and viscoelastic constitutive model for shale[J]. Chinese Journal of Rock Mechanics and Engineering，2004，23(8)：1 242–1 246.(in Chinese)) [12] 韦立德，徐卫亚，朱珍德，等. 岩石黏弹塑性模型的研究[J]. 岩土
力学，2002，23(5)：583–586.(Wei Lide，Xu Weiya，Zhu Zhende，et al. Study on constitutive model of viscoelastoplastic rock[J]. Rock and Soil Mechanics，2002，23(5)：583–586.(in Chinese))
[13] 金丰年，浦奎英. 关于黏弹性模型的讨论[J]. 岩石力学与工程学
报，1995，14(4)：335–361.(Jin Fengnian，Pu Kuiying. Discussion on visco-elastic models[J]. Chinese Journal of Rock Mechanics and Engineering，1995，14(4)：335–361.(in Chinese))
[14] 邓荣贵，周德培，张倬元，等. 一种新的岩石流变模型[J]. 岩石力
学与工程学报，2001，20(6)：780–784.(Deng Ronggui，Zhou Depei，Zhang Zhuoyuan，et al. A new rheological model for rocks[J]. Chinese Journal of Rock Mechanics and Engineering，2001，20(6)：780–784.(in Chinese))
[15] 曹树刚，边金，李鹏. 岩石蠕变本构关系及改进的西原正夫
模型[J]. 岩石力学与工程学报，2002，21(5)：632–634.(Cao Shugang，Bian Jin，Li Peng. Rheologic constitutive relationship of rocks and a modified model[J]. Chinese Journal of Rock Mechanics and Engineering，2002，21(5)：632–634.(in Chinese))
[16] 陈沅江，潘长良，曹平，等. 软岩流变的一种新力学模型[J]. 岩
土力学，2003，24(2)：209–214.(Chen Yuanjiang，Pan Changliang，Cao Ping，et al. A new mechanical model for soft rock rheology[J].
Rock and Soil Mechanics，2003，24(2)：209–214.(in Chinese)) [17] 张向东，李永靖，张树光，等. 软岩蠕变理论及其工程应用[J]. 岩
石力学与工程学报，2004，23(10)：1 635–1 639.(Zhang Xiangdong，Li Yongjing，Zhang Shuguang，et al. Creep theory of soft rock and its engineering application[J]. Chinese Journal of Rock Mechanics and Engineering，2004，23(10)：1 635–1 639.(in Chinese))
第25卷 第3期 徐卫亚等. 岩石非线性黏弹塑性流变模型(河海模型)及其应用 ? 447·
[18] 王来贵，何 峰，刘向峰，等. 岩石试件非线性蠕变模型及其稳定
性分析[J]. 岩石力学与工程学报，2004，23(10)：1 640–1 642.(Wang Laigui ，He Feng ，Liu Xiangfeng ，et al. Nonlinear creep model and stability analysis of rock[J]. Chinese Journal of Rock Mechanics and Engineering ，2004，23(10)：1 640–1 642.(in Chinese))
[19] 徐卫亚，杨圣奇，杨松林，等. 绿片岩三轴流变力学特性的研究(I)：
试验结果[J]. 岩土力学，2005，26(4)：531–537.(Xu Weiya ，Yang Shengqi ，Yang Songlin ，et al. Investigation on triaxial rheological mechanical properties of greenschist specimen (I)：experimental results[J]. Rock and Soil Mechanics ，2005，26(4)：531–537.(in Chinese))
[20] 徐卫亚，杨圣奇，谢守益，等. 绿片岩三轴流变力学特性的研究(II)：
模型分析[J]. 岩土力学，2005，26(5)：693–698. (Xu Weiya ，Yang
Shengqi ，Xie Shouyi ，et al. Investigation on triaxial rheological mechanical properties of greenschist specimen (II)：model analysis[J]. Rock and Soil Mechanics ，2005，26(5)：693–698.(in Chinese)) [21] Boukharov G N ，Chanda M W ，Boukharov N G . The three processes of
brittle crystalline rock creep[J]. Int. J. Rock Mech. Min. Sci. and Geomech. Abstr.，1995，32(4)：325–335.
[22] 徐卫亚，杨圣奇，谢守益. 锦屏一级水电站坝基岩石流变特性及参
数研究报告[R]. 南京：河海大学岩土工程研究所，2004.(Xu Weiya ，Yang Shengqi ，Xie Shouyi. Study report on the rheological properties and parameters of rock in dam base of Jinping First Stage Hydropower Station[R]. Nanjing ：Institute of Geotechnical Engineering ，Hohai University ，2004.(in Chinese))
2003～2004年中国相关学科11种期刊来源指标和被引指标
来源指标
被引指标 来源
平均
平均
基金
总被影响 即年
他引 引用 被引
引用
期刊分类 期刊 名称
年份 文献量引文数 作者数论文比引频次
因子
指标 总引比 刊数 半衰期半衰期
2003 389(1)
11.18(2) 3.08(1)0.51(4)1888(1) 1.407(1)0.203(1)0.47(8) 246(3) 3.16(11)6.33(11)岩石力学
与工程学报 2004 801(1)13.09(2) 3.10(4)0.60(4)3057(1) 1.650(1)0.164(1)0.48(8) 304(3) 3.28(10)6.73(9)2003 234(4)7.39(7) 2.88(5)0.48(5)478(8)0.421(9)0.072(8)0.86(4) 156(8) 3.93(8)7.67(5)岩土力学
2004 438(2)11.82(4) 3.00(6)0.58(5)1004(4)0.848(3)0.116(3)0.69(7) 202(8) 3.73(9)7.88(2)2003 120(9)
12.06(1) 2.83(6)
0.57(2)
568(7)
0.48(8)
0.066(9)0.93(1) 228(4) 6.73(1)9.57(2)
力学 力学学报
2004 112(10) 11.91(3) 2.92(10)0.67(2)
648(8)0.674(6)0.062(5)0.94(1) 244(5) 6.06(2)9.25(1)2003 185(5)7.39(7) 2.80(7)0.42(7)1705(2)0.766(2)0.129(3)0.92(2) 343(1) 5.59(3)7.80(4)岩土工程学报
2004 177(5)
8.97(6) 2.98(7)0.47(8)2499(2) 1.058(2)0.129(2)0.94(1) 345(2) 5.85(3)7.39(7)2003 237(3)
8.58(3) 2.68(9)
0.45(6)
729(4)0.728(3)0.100(5)0.93(1) 213(5) 4.61(6)7.29(7)建筑科学 土木工程学报
2004 231(4)
8.58(7) 2.76(11)0.53(7)992(5)0.621(8)0.060(6)0.94(1) 268(4) 4.66(6)7.49(4)2003 168(6)7.83(5) 3.02(3)0.38(8)595(6)0.512(5)0.130(2)0.72(6) 171(7) 5.23(4)7.07(8)水文地质工程地质 2004 167(7) 6.59(10) 3.04(5)
0.40(9)
654(7)0.643(7)0.059(7)0.88(5) 221(6) 5.63(4)7.83(3)2003 79(11)
7.15(9) 2.83(5)0.26(10)306(9)0.487(7)0.075(7)0.69(7) 112(9) 3.85(9)6.43(10)地 质学 工程地质学报 2004 79(11)16.69(1) 2.97(8)0.82(1)499(9)0.778(4)0.037(10)0.74(6) 153(9)
2003 140(8)7.22(8) 3.05(2)0.63(1)670(5)0.505(6)0.085(6)0.92(2) 199(6) 4.69(5) 6.72(9)煤炭学报
2004 172(6)7.37(9) 3.28(2)
0.55(6)
866(6)0.581(9)0.023(11)0.90(4) 217(7) 5.23(5)6.42(10)
2003 167(7)7.58(6) 3.05(2)0.24(11)232(10)0.303(10)0.035(11)0.85(5) 110(10) 3.67(10)7.63(6)矿业工程
矿冶工程
2004 161(8)8.22(8) 3.36(1)0.22(11)264(10)0.297(10)0.055(8)0.90(4) 123(10) 4.10(8)7.40(6)2003 267(2)8.29(4) 2.99(4)0.54(3)1291(3)0.541(4)0.123(4)0.88(3) 295(2) 6.01(2)9.11(3)水利学报
2004 256(3)
9.69(5) 3.12(3)
0.63(3)1841(3)0.776(5)0.097(4)0.92(3) 376(1) 6.08(1)7.44(5)2003 107(10) 5.46(10) 2.71(8)
0.30(9)
160(11)0.248(11)0.046(10)0.85(5) 74(11) 4.12(7)9.63(1)
水利
工程 长江科学院院报 2004 115(9) 5.37(11) 2.96(9)0.32(10)227(11)0.285(11)0.043(9)0.93(2) 107(11) 4.50(7)5.88(11)
注：( )中数字为排名顺序。
(摘编自《中国期刊引证研究报告》(2005年9月出版))
岩土类材料弹塑性力学模型及本构方程 TYYGROUP system office room 【TYYUA16H-TYY-TYYYUA8Q8-
岩土类材料的弹塑性力学模型及本构方程 摘要：本文主要结合岩土类材料的特性，开展研究其在受力变形过程中的弹性及塑性变形的特点，描述简化的力学模型特征及对应的适用条件，同时在分析研究其弹塑性力学模型的基础上，探究了关于岩土类介质材料的各种本构模型，如M-C、D-P、Cam、D-C、L-D及节理材料模型等，分析对应使用条件，特点及公式，从而推广到不同的材料本构模型的研究，为弹塑性理论更好的延伸发展做一定的参考性。 关键词：岩土类材料，弹塑性力学模型，本构方程 不同的固体材料，力学性质各不相同。即便是同一种固体材料，在不同的物理环境和受力状态中，所测得的反映其力学性质的应力应变曲线也各不相同。尽管材料力学性质复杂多变，但仍是有规律可循的，也就是说可将各种反映材料力学性质的应力应变曲线，进行分析归类并加以总结，从而提出相应的变形体力学模型。 第一章岩土类材料 地质工程或采掘工程中的岩土、煤炭、土壤，结构工程中的混凝土、石料，以及工业陶瓷等，将这些材料统称为岩土材料。 岩土塑性力学与传统塑性力学的区别在于岩土类材料和金属材料具有不同的力学特性。岩土类材料是颗粒组成的多相体，而金属材料是人工形成的晶体材料。正是由于不同的材料特性决定了岩土类材料和金属材料的不同性质。归纳起来，岩土材料有3点基本特性：1.摩擦特性。2.多相特性。3.双强度特性。另外岩土还有其特殊的力学性质：1.岩土的压硬性，2.岩土材料的等压屈服特性与剪胀性，3.岩土材料的硬化与软化特性。4.土体的塑性变形依赖于应力路径。 对于岩土类等固体材料往往在受力变形的过程中，产生的弹性及塑性变形具备相应的特点，物体本身的结构以及所加外力的荷载、环境和温度等因素作用，常使得固体物体在变形过程中具备如下的特点。 固体材料弹性变形具有以下特点：(1)弹性变形是可逆的。物体在变形过程中，外力所做的功以能量（应变能)的形式贮存在物体内，当卸载时，弹性应变能将全部释放出来，物体的变形得以完全恢复；(2)无论材料是处于单向应力状态，还是复杂应力状态，在线弹性变形阶段，应力和应变成线性比例关系；(3)对材料加载或卸载，其应力应变曲线路径相同。因此，应力与应变是一一对应的关系。 固体材料的塑性变形具有以下特点：(l)塑性变形不可恢复，所以外力功不可逆。塑性变形的产生过程，必定要消耗能量(称耗散能或形变功)；(2)在塑性变形阶段，应力和应变关系是非线性的。因此，不能应用叠加原理。又因为加载与卸载的规律不同，应力与应变也不再存在一一对应的关系，也即应力与相应的应变不能唯一地确定，而应当考虑到加载的路径(即加载历史)；(3)当受力固体产生塑性变形时，将同时存在有产生弹性变形的弹性区域和产生塑性变形的塑性区域。并且随着载荷的变化，两区域的分界面也会产生变化。 第二章弹塑性力学中常用的简化力学模型 对于不同的材料，不同的应用领域，可以采用不同的变形体模型。在确定力学模型时，要特别注意使所选取的力学模型必须符合材料的实际情况，这是非常重要的，因为只有这样才能使计算结果反映结构或构件中的真实应力及应
第七章 粘弹塑性模型的基本概念 7 . 1 引言 为了描述土体应力一应变关系受时间的影响，需要采用与时间有关的类模型（如粘弹胜模酬、粘塑性模型，粘弹塑隆模型）来描述土的性状。 弹性、塑性和粘性是连续介质的三种基本性质，各在定条件F 独自反映材料本构关系的一个方面的特性。理想弹性模型、理想塑胜模型（或称刚塑性模型）和理想粘性模型是反映这三种性质的理想模型，通常称为简单模型。实际工程材料的本构关系可以用这些简单模型的各种组合来构成。 理想弹性模型又称虎克弹性模型，通常用理想弹簧表示（图7-1( a )）。其本构方程为虎克定律。一维条件下，如单轴压缩和纯剪清况下，表达式分别为： （7.1.1） （7.1.2） 式中E —— 弹性模量、 G——剪切模量。 剪切模量与弹性模量和泊松比的关系如下式所示： （7.1.3） 式中 ——泊松比。 三维条件下本构方程可表示为下述形式： （7.1.4） 式中 K——体积弹性模量。
（a） （b） 图7-1 理想弹性模型 体积弹性模量与弹性模量和泊松比的关系如下式所示： （7.1.6） 理想粘性模型又称牛顿粘滞体模型。通常用一粘壶（或称阻尼器）表示（图7-2 ( a ) ）。粘壶内充满粘滞液体和一个可移动的活
塞。活塞在粘滞液体中的移动速度与所受阻力成正比关系，反映了粘性介质内一点的应力与该点处应变速率成正比例关系的性质。一维条件如单轴压缩或纯剪情况下，表达式分别为： （7.1.7） （7.1.8） 式中 、 ——粘滞系数。 由上两式可以看出，从数学表达的形式上与理想弹性体单轴压缩和纯剪时的本构方程相类似。 与理想弹性体的方程相对应，类似式7.1.3，存在下述关系： （7.1.9） 式中 ——粘性应变速率的横向比值。
常用弹塑性材料模型下表列出了ANSYS/LS-DYNA材料模型以及相应的LS-DYNA命令 ANSYS Material Model LS-DYNA Command LS-DYNA MAT # Example Isotropic Elastic*MAT_ELASTIC1Yes Bilinear Isotropic Plasticity *MAT_PLASTIC_KINEMATIC 3 Yes Bilinear Kinematic *MAT_PLASTIC_KINEMATIC 3 Yes Plastic Kinematic *MAT_PLASTIC_KINEMATIC 3 Yes Piecewise Linear Plasticity *MAT_PIECEWISE_LINEAR_PLASTICITY24 Yes Rigid *MAT_RIGID 20 Yes 7.2.1.1各向同性弹性模型 各向同性弹性模型。使用MP命令输入所需参数： MP，DENS—密度 MP，EX—弹性模量 MP，NUXY—泊松比 此部分例题参看B.2.1，Isotropic Elastic Example:High Carbon Steel。 B.2.1. Isotropic Elastic Example: High Carbon Steel MP,ex,1,210e9 ! Pa MP,nuxy,1,.29 ! No units MP,dens,1,7850 ! kg/m3 7.2.3.1 双线性各向同性模型 使用两种斜率（弹性和塑性）来表示材料应力应变行为的经典双线性各向同性硬化模型 （与应变率无关）。仅可在一个温度条件下定义应力应变特性。（也有温度相关的本构模型； 参看Temperature Dependent Bilinear Isotropic Model）。用MP命令输入弹性模量（Exx）， 泊松比（NUXY）和密度（DENS），程序用EX和NUXY值计算体积模量（K）。用TB和TBDATA 命令的1和2项输入屈服强度和切线模量： TB，BISO TBDATA，1， Y σ（屈服应力） TBDATA，2， tan E（切线模量） 例题参看B.2.7,Bilinear Isotropic Plasticity Example：Nickel Alloy。
第七章粘弹塑性模型的基本概念 7 . 1 引言 为了描述土体应力一应变关系受时间的影响，需要采用与时间有关的类模 型（如粘弹胜模酬、粘塑性模型，粘弹塑隆模型）来描述土的性状。 弹性、塑性和粘性是连续介质的三种基本性质，各在定条件F独自反映材料本构关系的一个方面的特性。理想弹性模型、理想塑胜模型（或称刚塑性模型）和理想粘性模型是反映这三种性质的理想模型，通常称为简单模型。实际工程材料的本构关系可以用这些简单模型的各种组合来构成。 理想弹性模型又称虎克弹性模型，通常用理想弹簧表示（图 7-1（ a ））。其本构方程为虎克定律。一维条件下，如单轴压缩和纯剪清况下，表达式分别为： E （7.1.1） G （7.1.2）式中E——弹性模量、 G剪切模量。 剪切模量与弹性模量和泊松比的关系如下式所示： G E—（7.1.3） 2 1 式中——泊松比。 三维条件下本构方程可表示为下述形式： (7.1.4) m K 式中K ——体积弹性模量。
(a) (b) 图7-1理想弹性模型 体积弹性模量与弹性模量和泊松比的关系如下式所示: (7.1.6) 理想粘性模型又称牛顿粘滞体模型。通常用一粘壶(或称阻尼器)表示(图7-2 ( a ))。粘壶内充满粘滞液体和一个可移动的活塞。活塞在粘滞液体中的移动速度与所受阻力成正比关系，反映了粘性介质内一点的应力与该点处应变速率成正比例关系的性质。一维条件如单轴压缩或纯剪情况下，表达式分别为： & (7.1.7) & (7.1.8) 式中、——粘滞系数由上两式可以看出，从数学表达的形式上与理想弹性体单轴压缩和纯剪时的本构方程相类似。 与理想弹性体的方程相对应，类似式7.1.3，存在下述关系： (7.1.9)
第七章-粘弹塑性模型的基本概念
第七章 粘弹塑性模型的基本概念 7 . 1 引言 为了描述土体应力一应变关系受时间的影响，需要采用与时间有关的类模型（如粘弹胜模酬、粘塑性模型，粘弹塑隆模型）来描述土的性状。 弹性、塑性和粘性是连续介质的三种基本性质，各在定条件F 独自反映材料本构关系的一个方面的特性。理想弹性模型、理想塑胜模型（或称刚塑性模型）和理想粘性模型是反映这三种性质的理想模型，通常称为简单模型。实际工程材料的本构关系可以用这些简单模型的各种组合来构成。 理想弹性模型又称虎克弹性模型，通常用理想弹簧表示（图7-1( a )）。其本构方程为虎克定律。一维条件下，如单轴压缩和纯剪清况下，表达式分别为： E σε= （7.1.1） G τγ= （7.1.2） 式中E —— 弹性模量、 G ——剪切模量。 剪切模量与弹性模量和泊松比的关系如下式所示： () 21E G ν=+ （7.1.3） 式中 ν ——泊松比。 三维条件下本构方程可表示为下述形式： m K νσε= （7.1.4） 式中 K ——体积弹性模量。
（a ） （b ） 图7-1 理想弹性模型 体积弹性模量与弹性模量和泊松比的关系如下式所示： () 312E K ν=- （7.1.6） 理想粘性模型又称牛顿粘滞体模型。通常用一粘壶（或称阻尼器）表示（图7-2 ( a ) ）。粘壶内充满粘滞液体和一个可移动的活塞。活塞在粘滞液体中的移动速度与所受阻力成正比关系，反映了粘性介质内一点的应力与该点处应变速率成正比例关系的性质。一维条件如单轴压缩或纯剪情况下，表达式分别为： σ?ε=& （7.1.7） τηγ=& （7.1.8） 式中 ?、η ——粘滞系数。 由上两式可以看出，从数学表达的形式上与理想弹性体单轴压缩和纯剪时的本构方程相类似。 与理想弹性体的方程相对应，类似式7.1.3，存在下述关系： ()*21? ην=+ （7.1.9）
《软土地基》课程论文 学院建工学院 姓名王洋 学号
软土本构模型综述 1 引言 土体具有复杂的变形特征，如剪胀性、各向异性、受应力路径影响等。土体变形的这种复杂性是在复杂受力状态下表现出来的。复杂应力状态存在 6 个应力分量，也有 6 个应变分量。其间的关系是一种多因素物理量与多因素物理量之间的关系，不能由试验直接建立。须在简化条件的试验基础上，做某些假定及合乎规律的推理，从而提出某种计算方法，把应力应变关系推广到复杂应力状态。这种计算方法叫本构模型。 1.1 土的本构模型 发展到现在，土的本构模型数目众多，大致可以分为以下几大类: ( 1) 非线性模型; ( 2) 弹塑性模型; ( 3) 粘弹塑性模型; ( 4) 结构性模型。 对于软土而言，比较适用的一般为弹塑性模型。弹塑性模型是把总的变形分成弹性变形和塑性变形两部分，用虎克定律计算弹性变形部分，用塑性理论来解塑性变形部分。 1.2 变形假定 对于塑性变形，要作三方面的假定: ( 1) 破坏准则和屈服准则; ( 2) 硬化准则; ( 3) 流动法则。 不同的弹塑性模型，这三个假定的具体形式也不同。最常用的弹塑性模型为剑桥模型及其扩展模型。 2 剑桥模型与修正剑桥模型 1958 年，Roscoe 等发现了散粒体材料在孔隙比-平均有效应力-剪应力的三维空间里存在状态面的事实，1963 年，提出了著名的剑桥模型，1968 年，
形成了以状态面理论为基础的剑桥模型的完整理论体系。 Roscoe 等人将“帽子”屈服准则、正交流动准则和加工硬化规律系统地应用于Cam 模型之中，并提出了临界状态线、状态边界面、弹性墙等一系列物理概念，构成了第一个比较完整的土塑性模型。剑桥模型又被称为临界状态模型，是一个非常经典的弹塑性模型，它是第一个全面考虑重塑正常固结或弱超固结粘土的压硬性和剪胀性的模型，标志着土的本构理论发展新阶段的开始。 1968 年，Roscoe 等人在剑桥模型的基础上提出了修正剑桥模型，将原来的屈服面在p'，q 平面上修正为椭圆，并认为在状态边界面内土体变形是完全弹性的。在状态边界面内，增加的剪应力虽不产生塑性体积变形，但可产生塑性剪切变形。修正剑桥模型是一种“帽子”型模型，在许多情况下能更好地反映土的变形特性。修正剑桥模型至今仍在工程中广泛应用，是因为它具有很多优点: 形式简单，模型参数少，参数确定方法简单( 只需常规三轴试验即可) ，参数有明确的物理意义，能够很好的反映重塑正常固结或弱超固结粘土的压硬性和剪缩性，因此修正剑桥模型是土力学中比较成熟而且应用广泛的弹塑性本构模型。同时，修正剑桥模型也有一定的局限性: 屈服面只是塑性体积应变的等值面，只采用塑性体积应变作硬化参量，因而没有充分考虑剪切变形; 只能反映土体剪缩，不能反映土体剪胀; 没有考虑土的结构性这一根本内在因素的影响; 假定的弹性墙内加载仍会产生塑性变形等。修正剑桥模型对实际情况进行了一系列假定: ①屈服只与应力球量p 和应力偏量q 两个应力分量有关，与第三应力不变量无关; ②采用塑性体应变硬化规律，以为硬化参数; ③假定塑性变形符合相关联的流动法则，即g( σ) = f( σ) ; ④假定变形消耗的功，即塑性功为: 剑桥模型是当前在土力学领域内应用最广的模型之一，其主要特点有: 基本概念明确; 较好地适宜于正常固结粘土和弱超固结粘土; 仅有3个参数，都可以通过常规三轴试验求出，在岩土工程实际工作中便于推广; 考虑了岩土材料静水压力屈服特性、剪缩性和压硬性。王清等分析了修正剑桥模型的应力应变关系，以其为基础引进了接触单元和杆单元，运用修正合格模型，用有限元程序模拟了
第!"卷第#期应用力学学报$%&’!"(%’# +,-’)**! )**!年!)月!"#$%&%’()*$+,(-+..,#%/0%!"+$#!& !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !文章编号：!***.#/0/（)**!）*#.**01.*2 一种适合橡胶类材料的非线性粘弹性本构模型" 安群力危银涛杨挺青 （华中科技大学武汉#0**1#） 摘要 借助非线性流变模型建立大变形情况非线性粘弹材料的本构关系，考虑到大多数橡胶类材料具有的几乎不可压缩性，以及体积响应和剪切响应的流变性能不同，将 变形梯度乘法分解为等容部分和体积变形部分，给出了一种适合橡胶类材料的非线 性粘弹性本构模型，并模拟了粘滞效应。对于极快或极慢的过程，该模型退化为橡胶 弹性理论；在小变形情况下退化为经典的广义3456,&&粘弹性材料。模型与热力学 第二定律相容，适合于大规模数值分析。 关键词：橡胶；粘弹性；有限变形；本构关系 中图分类号：70#2；8900文献标识码：: !引言 在橡胶结构的设计与分析中因橡胶类材料力学性能的复杂性使得数值方法起着越来越重要的作用［!］。目前，应用数值方法时缺乏适于大规模计算用的本构关系，本构模型成为解决问题的关键［);!*］。构造粘弹材料的本构模型，一种方法是从连续介质力学本构理论的基本原理出发，经过简化而得到［!*;!)］。另外一种常用的方法是基于内变量理论，借助于连续介质热力学和流变模型来确定材料的本构模型［#;/，!0;!<］。在通常的内变量理论中，自由能的构造、内变量的选取及演化方程的确定有一定的困难。 本文利用非线性流变模型，认为总应力等于弹性应力与非弹性应力的和，通过平衡应变能函数表述其演化方程，绕过了通常内变量理论的困难，在参考位形内建立了以=>%&4.?>@-AA%BB 应力和C@,,D应变表示的大变形非线性粘弹性本构关系，给出了一种适合橡胶类材料的非线性粘弹性本构模型，物理意义简明。在一定条件下模型可以退化为相应的弹性或线粘弹性模型，讨论了材料的粘滞现象。 "基金项目：国家自然科学基金资助项目（!/<0)*0*）来稿日期：)***.*#.*0修回日期：)***.!!.!1 万方数据 第一作者简介：安群力，男，!/<"年生，博士，华中科技大学力学系；研究方向：粘弹塑性理论及其应用E
沥青混合料粘弹塑性本构模型的实验研究沥青混凝土路面是近年来高速公路广泛采用的一种结构形式,随着公路运输量日益增长和运输向重型方向发展,路面破坏日趋严重。进行沥青混合料本构模型的研究,对掌握路面变形规律,预测路面结构永久变形大小,预防和抑制路面损害具有十分重要的意义。 文章针对沥青混合料单轴压缩、蠕变和恢复等力学特性,在实验基础上,结合理论和数值拟合分析,建立了沥青混合料不同形式的粘弹塑性本构模型,提出了模型参数确定方法,讨论了加载应力和环境温度对混合料力学行为的影响,并将模型预测结果与实验结果进行了比较,最后还初步分析了集料级配对沥青混合料力学行为的影响。主要内容包括:(1)提出并建立了沥青砂微分型粘弹塑性本构模型。 依据沥青砂蠕变特性,将总变形分解为粘弹性、粘塑性二种分量,采用Burgers模型描述粘弹性变形,采用滑块与粘壶并联模型描述粘塑性变形,然后加以组合,提出了基于二变形分量的粘弹塑性本构模型;进一步细分,将总变形分解为粘弹性、粘塑性和弹塑性三种分量,分别采用不同子模型描述上述分量,然后组合这些子模型,提出了基于三变形分量的粘弹塑性本构模型。基于较优模型,利用实验数据建立了参数与环境温度和加载应力的函数表达式,通过模型预测与实验结果的比较,证实模型可以较好地描述沥青砂三个蠕变阶段的变形特点。 (2)提出并建立了沥青砂、沥青混合料积分型粘弹塑性本构模型。将总变形分解为粘弹性和粘塑性变形,分别采用Schapery非线性模型描述粘弹性变形,采用Uzan模型描述粘塑形变形,提出了改进的Schapery积分模型,建立了积分型的非线性粘弹塑性本构关系,提出了非线性参数的实验确定方法,分别采用蠕变回
弹塑性分析 在这一册中,我们将详细地介绍由于塑性变性引起的非线性问题--弹塑性分析,我们的介绍人为以下几个方面: ?什么是塑性 ?塑性理论简介 ?ANSYS程序中所用的性选项 ?怎样使用塑性 ?塑性分析练习题 什么是塑性 塑性是一种在某种给定载荷下,材料产生永久变形的材料特性,对大多的工程材料来说,当其应力低于比例极限时,应力一应变关系是线性的。另外,大多数材料在其应力低于屈服点时，表现为弹性行为，也就是说，当移走载荷时，其应变也完全消失。 由于屈服点和比例极限相差很小，因此在ANSYS程序中，假定它们相同。在应力一应变的曲线中，低于屈服点的叫作弹性部分，超过屈服点的叫作塑性部分，也叫作应变强化部分。塑性分析中考虑了塑性区域的材料特性。 路径相关性： 即然塑性是不可恢复的，那么这种问题的就与加载历史有关，这类非线性问题叫作与路径相关的或非保守的非线性。 路径相关性是指对一种给定的边界条件，可能有多个正确的解—内部的应力，应变分布—存在，为了得到真正正确的结果，我们必须按照系统真正经历的加载过程加载。 率相关性： 塑性应变的大小可能是加载速度快慢的函数，如果塑性应变的大小与时间有关，这种塑性叫作率无关性塑性，相反，与应变率有关的性叫作率相关的塑性。 大多的材料都有某种程度上的率相关性，但在大多数静力分析所经历的应变率范围，两者的应力－应变曲线差别不大，所以在一般的分析中，我们变为是与率无关的。 工程应力，应变与真实的应力、应变： 塑性材料的数据一般以拉伸的应力—应变曲线形式给出。材料数据可能是工程应力 ）。（P A0）与工程应变（?l l0），也可能是真实应力（P/A）与真实应变（n L l l() 0大应变的塑性分析一般采用真实的应力，应变数据而小应变分析一般采用工程的应力、应变数据。 什么时候激活塑性： 当材料中的应力超过屈服点时，塑性被激活（也就是说，有塑性应变发生）。而屈服应力本身可能是下列某个参数的函数。 ?温度 ?应变率 ?以前的应变历史 ?侧限压力 ?其它参数 塑性理论介绍 在这一章中，我们将依次介绍塑性的三个主要方面： ?屈服准则 ?流动准则 ?强化准则 屈服准则： 对单向受拉试件，我们可以通过简单的比较轴向应力与材料的屈服应力来决定是否有塑性变形发生，然而，对于一般的应力状态，是否到达屈服点并不是明显的。 屈服准则是一个可以用来与单轴测试的屈服应力相比较的应力状态的标量表示。因此，知道了应力状态和屈服准则，程序就能确定是否有塑性应变产生。
土体动本构模型的研究现状 土体实际动本构关系是极其复杂的,它在不同的荷载条件、土性条件及排水条件下表现出极不相同的动本构特性. 要建立一个能适用于各种不同条件的动本构模型的普遍形式是不切实际的,其切实的方法是对于不同的工程问题,应该根据土体的不同要求和具体条件,有选择地舍弃部分次要因素,保留所有主要因素,建立一个能反映实际情况的动本构模型. 目前,具体建立的动本构模型已达数十个,大致可分为两大类,即粘弹性模型和弹塑性模型.曲线模型,均属于等效线性模型[2 ] 。Masing 类模型以曲线Hardin Drnevich 或Ram2berg Osgood 曲线等为骨干,改用瞬时剪切模量代替前面的平均剪切模量。为使这类动本构模型更接近实测的动应力应变曲线,很多学者做了大量的工作,以使其能够描述不规则循环荷载作用下土的动本构关系[3 ] 。Iwan 用一系列具有不同屈服水平的理想弹塑性元件来描述土的动本构关系,它分串联型和并联型2 种构成方式。串联型和并联型的伊万模型所描述的动应力应变特性基本上一致,只是前者以应变为自变量,后者以应力为自变量[4 ] 。郑大同在伊万模型的基础上,提出了一个新物理模型,该模型的骨架曲线可为加工硬化状,也可为加工软化状,骨架曲线与滞回曲线的2 个分支既可相同,也可不同[5 ] 。一般的粘弹性模型不能计算永久变形(残余变 形) ,在主要为弹性变形的情况下比较合适。但实际上,土在往复荷载作用下还会因土粒相互滑移,形成新的排列而产生不可恢复的永久变形。为此,Mar2tin 等人根据等应变反复单剪试验结果,提出了循环荷载作用下永久体积应变的增量公式[6 ] 。后来,日本学者八木、大冈和石桥等分别由等应力动单剪试验及扭剪试验各自提出了计算永久体积应变增量的经验公式。国内的姜朴、徐亦敏、娄炎根据动三轴试验应变与破坏振次的关系式。沈珠江[7 ] 对等价粘 弹性模型进行了较全面的研究,认为一个完整的粘弹性模型应该包含4 个经验公式: (1) 平均剪切模量; (2) 阻尼比; (3) 永久体积应变增量和永久剪切应变增量; (4) 当饱和土体处于完全不排水或部分排水条件下,还需给出孔隙水压力增长和消散模型。粘弹性理论是目前应用中的主流,但存在多方面的不足,如不能考虑应变软化,不能考虑应力路径的影响,不能考虑土的各向异性以及大应变时误差大等,但它是试验结果的归纳,形式上直观简单,经过处理改进后,结合有限元程序,就可以计算出循环荷载作用下土工构造物的孔隙水压力和永久变形的 平均发展过程。 211 粘弹性理论 人们早在生产实践中认识到土体的应力—应变关系是非线性的,但实际工程中常用线性理论对这种非线性关系进行简化。自Seed 提出用等价线性方法近似考虑土的非线性以来,粘弹性理论已有了较大的发展。在土体的动力反应分析中,常用的粘弹性理论有等效线性模型和曼辛型非线性模型2 大类。前者把土体视为粘弹性材料,不寻求滞回曲线(即描述卸载与再加载时应力应变规律的曲线) 的具体数学表达式,而是给出等效弹性模量和等效阻尼比随剪应变幅值和有效应力状态变化的表达式,即以G 和λ作为它的动力特性指标引入实际计算;后者则根据不同的加载条件、卸载和再加载条件直接给出动应力应变的表达式。在给出初始加载条件下的动应力应变关系式(骨干曲线方程) 后,再利用曼辛二倍法得出卸荷和再加荷条件下的动应力应变关系,以构成滞回曲线方程[1 ] 。Hardin Drnevich 模型、Ramberg Osgood 模型、双线性模型及一些组合 基于阻尼的地震循环荷载作用下黏土非线性模型 尚守平刘方成王海东 ( 湖南大学, 湖南长沙410082) 摘要: 提出一种基于阻尼比的黏土动应力应变模型, 通过在滞回曲线中显示地引入代表阻尼比大小的形状系数,使得理论滞回曲线真实地反应土体的滞回阻尼性能。首先推导在等幅对称
第27卷 第3期 岩 土 工 程 学 报 Vol.27 No.3 2005年 3月 Chinese Journal of Geotechnical Engineering Mar., 2005 考虑强度非线性的堆石料弹塑性本构模型研究Study on elasto-plastic constitutive model of rockfills with nonlinear strength characteristics 刘萌成1，高玉峰2，黄晓明1 （1.东南大学 交通学院, 江苏 南京 210096；2.河海大学 岩土工程研究所,江苏 南京 210098） 摘 要：根据堆石料大型三轴试验成果，提出了非线性强度包线在q—p坐标面上的幂函数表达式。将此非线性强度准则模式应用于剪切屈服面与体积屈服面表达方程的推导，从而建立了一个可反映堆石料实际工程特性的新型弹塑性本构模型。应用自编计算程序将该模型与殷氏模型用于不同应力路径条件下堆石料试验的数值模拟研究，对比分析表明：本文模型计算结果与试验实测结果更为一致，能较好地描述各种加载路径下堆石料的强度与变形性状。 关键词：堆石料；大型三轴试验；非线性强度模式；双屈服面本构模型；数值模拟 中图分类号：TU 521文献标识码：A文章编号：1000–4548(2005)03–0294–05 作者简介：刘萌成(1972–)，男，江西樟树人，博士，主要从事路基与路面材料工程特性方面的研究工作。 LIU Meng-cheng1, GAO Yu-feng2, HUANG Xiao-ming1 (1.Transportation College in Southeast University, Nanjing 210096，China 2. Research Institute of Geotechnical Engineering in Hohai University，Nanjing 210098，China) Abstract: Based on the results of large-scale triaxial test of rockfills, an exponential formulation of the nonlinear strength envelope is presented. When the nonlinear strength criterion is introduced to determine the shear and volume yield surfaces, a new elasto-plastic constitutive model is herein developed for reproducing practical engineering behaviors of the rockfills. The numerical simulations for the various stress-path test are performed in a developed program with the proposed model and Yin’s model. It is shown that the proposed model can realistically describe the strength and deformation characteristics of rockfills subjected to some complex loading conditions. Key words: rockfills; large scale triaxial testing; nonlinear strength model; two-surface constitutive model; numerical simulation 0 引 言 材料的本构关系是反映材料力学性状的数学表达式。自20世纪50年代以来，岩土介质本构关系研究随工程实践的进行而不断深入发展。目前，Duncan-Chang模型与MCC（修正剑桥）模型已为人们所普遍接受和认可。为了更为合理地反映岩土材料的复杂工程特性，国内外许多学者提出了剪切与体积屈服面相结合的双屈服面模型[1，2]。其中，沈珠江与殷宗泽提出的双屈服面弹塑性本构模型兼具MCC模型和Duncan-Chang模型的优点，能较好地反映土体非线性、剪胀性、变形的应力路径相关性等复杂特性[3，4]。 堆石料是指山体开采所得岩块经一定程度粉碎而得到的岩石碎块类集合体，属无粘性混合料。通常认为，堆石料特性研究是在砂的特性研究基础上展开的[5]，堆石料变形特性的数学表述大多是在砂土本构模型基础上修改而得到的[6]，因此细粒土本构模型亦被推广用于堆石料变形特性研究。然而，大量研究表明，堆石料具有粒径大、沉陷变形小、低压剪胀及颗粒破碎明显的特性，使其强度变形特性更为复杂而与细粒土性质存在较大差别[5,7,8]。国内外众多学者认为：低 压下的剪胀与高压下的破碎使得抗剪强度包线呈现为非线性。为此，Vesic&Clough、De Mello、Indraratna 与张启岳等人提出了反映堆石料强度包线非线性的相应函数表达式[7~11]。鉴于此，为反映堆石料强度特性对变形的显著影响，有必要引入非线性强度准则，重新推导得到新的双屈服面本构模型表达式，以便较为真实地描述堆石料在各种加载路径下的力学性状。 1 强度非线性研究 采用河海大学岩土工程研究所和四川大学岩土仪器研究所联合研制而成的YS—30型微机控制应力路径大型三轴试验仪进行堆石料大型三轴试验研究[12]。在宜兴抽水蓄能电站筑坝材料中，3-10-15号试样（以下简称试样），岩性为西梅园料场石英砂岩，试样基本性质指标如表1所示。图1给出了该试样在不同围压下各向等压固结排水（CID）剪切试验下的摩尔圆与抗剪强度包线拟合曲线，由图可见：在围压较低(kPa 300 3 < σ)时，内摩擦角增大，强度包线明显弯曲且通过原点,显示出堆石料在低围压下就具有较明─────── 收稿日期:2004–03–11
使用纤维模型做桥梁的动力弹塑性分析 北京迈达斯技术有限公司 2004.12
目 录 1.概要 2.纤维单元的特性 3.桥梁资料 4.建立结构模型 5.定义纤维单元 6.结构的非线性特性 7.定义时程分析数据 8.运行结构分析 9.定义分析结果函数 10.查看分析结果
使用纤维单元做预应力桥梁的动力弹塑性分析 1. 概要 纤维单元是将梁单元截面分割为许多只有轴向变形的纤维的模型，使用纤维模型时 可利用纤维材料的应力-应变关系和截面应变的分布形状假定较为准确地截面的弯 矩-曲率关系，特别是可以考虑轴力引起的中和轴的变化。但是因为使用了几种理 想化的骨架曲线(skeleton curve)计算反复荷载作用下梁的响应，所以与实际构件 的真实响应还是有些误差。 MIDAS/Civil中的纤维模型使用了下面的几个假定。 ?截面的变形维持平截面并与构件轴线垂直。 ?不考虑钢筋与混凝土之间的滑移(bond-slip)。 ?梁单元截面形心的连线为直线 通过下面例题，介绍使用纤维单元做动力弹塑性分析的步骤。 因为本例题说明侧重于利用纤维单元做动力弹塑性分析的介绍，所以省略了前期建 模的过程，并认为用户已经熟练掌握了MIDAS/Civil的建模方法。 本例题模型为三维预应力梁桥的实际桥梁模型(2002年11月建，韩国)，但为了说明 上的便利，进行了一些简化处理，最后的结果有可能与实际设计稍有差异。 MIDAS/Civil中使用纤维单元做动力弹塑性分析的步骤如下: 1.定义纤维模型的材料特性 2.定义纤维模型的截面特性 3.定义并分配构件的非弹性铰特性 4.输入动力弹塑性时程分析数据 5.运行分析 6.定义分析结果函数 7.检查并验算分析结果 1