基于D_P准则的三维弹塑性有限元增量计算的有效算法

基于D-P准则的三维弹塑性有限元

增量计算的有效算法

A practical3D ela sto2pla stic incremental method in FEM

ba sed on D-P yield criteria

杨　强,陈　新,周维垣

(清华大学水利系,北京　100084)

摘　要:针对岩土材料常用的D-P准则,提出了一种新的增量分析方法,不用形成弹塑性增量矩阵,直接导出了符合正交流动法则的转移应力的解析解。该方法无论是对小步长还是大步长加载均有良好的收敛性。当采用精细的步长划分时,它就是严格意义上的理想弹塑性增量计算。在大步长情况下,在收敛域内最大载荷低于结构真实的极限承载力;对应的应力场是一个静力容许应力场;同时由于正交流动法则在平均意义下得到满足,收敛域内最大载荷接近结构真实的极限承载力。按此法所得结果接近真解且偏于安全。将整个计算模型装入三维非线性有限元程序TFI NE中,对某拱坝进行了超载分析。

关键词:转移应力;极限载荷;点安全度

中图分类号:T U452 文献标识码:A 文章编号:1000-4548(2002)01-0016-05

作者简介:杨　强(1964-),男,云南人。1988年在清华大学获硕士学位,1996年在奥地利Innsbruck大学获博士学位,现为清华大学水利系教授。主要从事水工结构及岩石力学方面的研究工作。

Y ANG Qiang,CHE N X in,ZH OU Wei2yuan

(Department of Hydraulic Engineering,Tsinghua University,Beijing100084,China)

Abstract:In this paper,focused on popularly used D-P yield criteria in geomaterials,a new incremental method in which the stresses to be trans2 ferred according to normal flow rule are directly derived without forming elasto2plastic increment matrix,was proposed.This method converges for either small load steps or large load steps.When very small load steps are used,the method is equivalent to standard elasto2plastic incremental method.When large load steps are used,the maximum load applied is lower than limit load in structure,the calculated stress field is an static ad2 missible one.As normal flow rule is satisfied in average,the maximum load is close to limit load.The soltion calculated by the method is on the safe side and close to real solution.The method was embedded into a3D nonlinear FEM software named TFINE,and overloading analysis was performed on an arch dam.

K ey words:the stresses to be transferred;normal flow rule;limit load

1　引 言Ξ

岩土材料具有很复杂的本构特性,如各向异性、硬化、软化等,目前描述岩土材料的本构模型非常多。但在实际工程三维有限元计算分析中,尤其是在岩体工程里,大量使用的仍是最简单D-P准则及理想弹塑性分析。其主要原因是参数选取不易。如在二滩高拱坝建设中,做了大量坝肩岩体现场大型抗剪试验,但具体到某一岩级,试验点数仍然很有限,且离散性很大。很难完全依赖试验确定参数,一般都要进行工程类比,对中、小工程工程类比更是参数确定的主要手段。最终一般只能给出岩体的抗剪参数f,c值。在这种情况下,从工程实用角度来说,追求本构关系的精致、完备并无太多实用意义。

相对而言,在岩土工程三维非线性有限元分析里,计算收敛性是一个较大的问题。弹塑性增量计算要采用精细的步长划分,才能确保计算收敛到正确解。在岩土工程,尤其是岩体工程里,荷载量级都很大,如高拱坝对水荷载的极限承载力可达上亿吨,而这对两岸高陡边坡所承受的的自重荷载来说,还只是一个小数,又如高地应力区大型地下洞室、高边坡(如三峡船闸高边坡)开挖过程中的释放荷载量级也十分巨大。若采用精细的步长划分,计算量将很大。岩体地质构造复杂,三维网格划分时经常会有畸形单元。由于地址缺陷或加固措施导致相邻单元材料性质差异过大,再加上高水平的荷载,各种因素交互影响,使得在计算过程中,经常出现局部发散现象,使得增量计算难以进行下去,最终结果可信度低,也难以从计算结果判断何时结构丧失稳定性。而对岩土工程来说,往往更关注结构的稳定性和极限承载力,而非应力和位移分布。

Ξ基金项目:国家自然科学基金资助项目(59879005);清华大学基础研究基金资助项目

收稿日期:2001-04-12

　第24卷　第1期岩　土　工　程　学　报V ol.24　N o.1 2002年 1月Chinese Journal of G eotechnical Engineering Jan.,　2002　

本文针对目前三维非线性有限元分析岩土材料常用的D -P 准则,提出了一种新的增量分析方法,它直接导出了符合正交流动法则的转移应力的解析解,不用形成弹塑性增量矩阵。该方法无论是对小步长还是大步长加载均有良好的收敛性。当采用精细的步长划分时,它就是严格意义上的理想弹塑性增量计算。当采用大步长增量计算,所得到是在平均意义上符合正交流动法则的一个静力容许应力场;在收敛域内,所能施加的最大荷载接近极限载荷的下限。

2　转移应力的解析解

在三维有限元弹塑性增量分析里,D -P 准则由于

简单实用,且和当前地质方面所能提供参数相适应,一直是目前应用最广的岩土材料屈服准则,其形式为

f =αI 1+

J 2-k ≤0

(1)

式中　J 2=

16

[(σ1-σ2)2+(σ2-σ3)2+(σ3-σ1)2],I 1=σ1+σ2+σ3,其中σ1,σ2,σ3为主应力。α和k 可通过拟合莫尔-库仑准则而得。例如,在π平面上,若

D -P 准则为库仑六边形的外接园,则

α=2sin φ3(3-sin φ),k =6c cos φ3(3-sin

φ)(2)式中　

φ和c 为材料的摩擦角和粘聚力。设某一高斯点在某一加载步或迭代步前初始应力为σ0且满足f (σ0)≤0。对某一个加载步或迭代步,由位移法求得该点应变增量为Δ

ε,它对应于弹性试应力σ1=σ1

ij =σ0+D :Δε,这里D 为弹性张量;若f (σ1)>0,则需进行应力调整。若此加载步或迭代步中塑性

应变增量为Δεp

,则调整后的应力为σ=σij =σ0+D :(Δε-Δεp

)=σ1-D :Δεp

。将正交流动法则d εp

=d λ

5f 5σ近似写成增量形式Δεp =Δ

λ5f 5σ,并以σ1确定5f Π5σ的代表值。则由条件

f (σ)=0,σ=σ1-ΔλD :5f 5σσ=σ1

(3)

即可确定调整后应力

σ=σij =(1-n )σ1

ij +p δij (4)

这里

n =

w

μJ 2

,p =-mw +

1

3nI 1,m =α(3λ+2μ

),w =f

3αm +μ

(5)

其中J 2,I 1,f 均由σ1确定。λ,μ为拉梅常数,即λ=E ν(1+ν)(1-2ν),μ=

E

2(1+ν

)(6)式中　E ,ν为杨氏模量和泊松比。容易证明由式(4)

确定的调整后应力σ必在屈服面上,即f (σ)=0。故对每个增量步或迭代步,应力转移值为

Δσ=σ1-σ=n σ1

ij -p δij (7) 求得由式(7)确定的转移应力Δσ后,可采用最普通的常刚度迭代进行应力转移计算。该方法避免了复杂的塑性矩阵的运算,直接求得满足正交流动法则的应力转移值。显然当采用精细的步长划分时,正交流动法则d εp

=d λ(5f Π5σ)得到充分满足,它就是严格意义上的理想弹塑性增量计算。以σ1确定5f Π5σ对计算收敛起到重要作用,此时转移应力的方向指向屈服锥体的中心线,所以调整后应力总能退回到屈服面上来。若以σ0确定5f Π5σ,如果σ1和σ0相差较大,则有可能调整后应力无法退回到屈服面上来。这实际上也是一般弹塑性计算当增量步过大就不收敛的重要原因之一。一般弹塑性计算也不进行校核。本文建议的方法在任何步长下都能确保调整后应力在屈服面上。大量的数值计算表明该方法具有良好的收敛性,所以可以用计算是否收敛来作为整体结构是否失稳的判据。

在大步长情况下求得的应力场与精细的步长划分下的所得是有差异的。我们可以从极限承载力的角度来探讨一下这个问题。在比例加载下,结构所能承受的最大荷载为极限承载力。在精细的步长划分下,理想弹塑性增量计算的最大荷载就是严格意义上的极限承载力。极限分析的下限定理认为:“一个满足平衡条件并且到处都不破坏材料屈服条件的应力场,是一个静力容许应力场;跟静力容许应力场对应的外荷载是极限荷载的下限;最高的下限便是极限荷载”。显然,在大步长情况下,在收敛域内最大载荷低于结构真实的极限承载力;对应的应力场是一个静力容许应力场;同时由于正交流动法则在平均意义下得到满足,收敛域内最大载荷接近结构真实的极限承载力。所以按此法所得的极限承载力接近真解且偏于安全。

3　其他考虑

岩土材料为低抗拉材料,故还应考虑抗拉条件:

σ1≤σt ,σ2≤σt ,σ3≤σt

(8)式中　σt 为材料单轴抗拉强度。在程序流程上,首先判断抗拉条件是否满足,若不满足,则调整应力满足之;其次再按D -P 准则判断调整应力。整个计算模型被装入三维非线性有限元程序TFI NE 里。在TFI NE 里,材料有三种破坏模式Ξ

:

7

1Ξ对理想弹塑性模型而言,屈服即为破坏。

　第1期杨　强,等1基于D -P 准则的三维弹塑性有限元增量计算的有效算法