弹粘塑性3D_FEM在地基蠕变沉降预测中的应用.pdf

收稿日期1999 - 04 - 07 作者简介安关峰1970 - ,男,河北秦皇岛人,博士生. 弹粘塑性3DΟFEM在地基蠕变 沉降预测中的应用 安关峰1,高大钊2 1. 中国地质大学 工程学院,湖北武汉 430074 ; 2.同济大学 地下建筑与结构工程系,上海 200092 摘要采用适用于上海饱和软土的粘弹塑性模型,研制了三维有限元软件.应用该软件对上海华盛路住宅大楼 地基沉降进行了预测,并与实测值进行了对比.结果表明,采用该本构模型和分步施加荷载方法进行地基蠕变沉 降预测是可靠的. 关键词粘弹塑性模型;流变;线性;三维有限元 中图分类号 TU 470.3 文献标识码 A 文章编号 0253 - 374X200102 - 0195 - 05 3DΟFEM Application to the Prediction of Creep Settlement of Soft Clay Considering Elastic2visco2plastic Constitution AN Guan2feng1, GAO Da2zhao2 1. Engineering School ,China University of Geoscience ,Wuhan 430074 ,China; 2.Department of Geotechnical Engineering , Tongji University ,Shanghai 200092 ,China Abstract The elastic2visco2plastic constitution has been employed in programming 3DΟFEM for the soft satu2 rated clay in Shanghai. The program also has been applied to predicting the foundation settlement of residence building at Huasheng road. The comparison between the predicted value and measured value proves that it is reliable to adopt the constitution and step loading in the creep settlement. Key words elastic2visco2plastic constitution model ; rheologic; linear ; three dimensional finite element 土不仅具有弹性和塑性性能,而且还具有粘滞性,即流变性,这已为人们所认识.我国广大沿海地区饱 和软弱土层触变性高、 压缩性大而又强度低,具有明显流变特性.加深对软粘土流变特性的认识,应用有限 元计算方法对实际岩土工程的长期安全进行分析具有重要的现实意义. 本文应用研制的三维粘弹塑性有限元计算软件对上海华盛路高层住宅大楼的地基沉降进行了预测, 得到了比较理想的效果. 1 弹粘塑性模型 据文献[1]对上海地区饱和软粘土的研究可知上海饱和软粘土流变以粘塑性为主,粘弹性是次要的, 可简化为弹粘塑性体.故本文采用线性弹粘塑性理论模型对上海软土地基在荷载作用下的蠕变沉降进行 研究见图1 . 在处理非线性连续体问题时,通常假定,可把总应变ε分为弹性分量εe和粘塑性分量εvp两部分,总 应变率可表示为 εεeεvp1 式中 “ ” 代表对于时间的微分. 第29卷第2期 2001年2月 同 济 大 学 学 报 JOURNAL OF TONG J I UNIVERSITY Vol. 29 No. 2 Feb. 2001 1. 1 弹性应变率 依赖于弹性应变率的总应力率按下式计算 σ Dεe2 式中D为弹性矩阵. 1. 2 粘塑性应变 根据标量形式的屈服条件 F σ,εvp - F003 图1 弹粘塑性模型 Fig. 1 Elastic2visco2plastic constitution 确定是否出现粘塑性状态.式中F0是单向 屈服应力,它本身可以是强化参数κ的函 数.这里假定,只有当FF0时才出现粘塑 性流动.本文的有限元程序含有四个常用屈 服准则,即Tresca准则、Von Mises准则、 Mohr2Coulomb准则、Drucker2Prager准则. 采用粘塑性应变速率只由当前应力确 定这一方案,作为确定粘塑性应变的法则, 因此有 εvp f σ4 由以下粘塑性流动法则给出了式 4的一个显式,即 εvpγ 9Q 9σ 5 式中Q Q σ,εvp,κ是 “塑性” 势;γ是控制塑性流动率的流动参数;Φ x对于x 0是正的单调递增函 数;符号 〈 〉 表示 当 x 0〈Φ x〉Φ x x≤0 〈Φ x〉0 6 采用相关联的塑性流动法则,这里有FQ,因而式1～5可变为 εvpγ 〈Φ F〉 9F 9σ γ 〈Φ〉 α7 式中的α为流动矢量.Φ有以下两种形式 Φ F e M F- F0 F0 -1 Φ F F -F0 F0 N 8 式中M和N是确定的常数,对岩土材料可取为1. 0. 1. 3 粘塑性应变增量 由式7表示的应变率,采用隐式的时间步进方案,则可由式9确定在时间间隔内所产生的应变增量为 Δ ε n vpΔtn [ 1-Θε n vpΘ ε n1 vp ]9 当Θ 0时,为 “全显式法”或前向差分法 ; 当Θ 1时,为 “全隐式法” 或后向差分法 ; Θ 1/2时, 为 “隐式梯形法”. 为求出式9中的 ε n 1 vp ,可以用有限Taylor级数展开式,写成 ε n1 vp ε n vp [ H n ]Δ σ n 10 式中 [ Hn] 9εvp 9σ n [ Hnσ n ] 11 而Δ σn是在时间间隔Δtn tn 1-tn内产生的应力改变量.于是式9可重新写为 691 同 济 大 学 学 报第29卷 Δ ε n vp ε n vpΔtn [ C n ]Δ σ n 12 式中[ Cn] Θ Δtn[ Hn] 13 式11定义矩阵H,它的特征值确定了可用于显式积分法的极限步长Δtn.矩阵H与应力值有关,其 具体计算可参考有关文献. 1. 4 应力增量 用式2的增量形式得到 Δ σ n [ D ]Δ ε n e [ D ]Δ ε n -Δ ε n vp 14 或者用位移增量来表示总的应变增量 Δ ε n [ B ]Δdn15 并由式12将Δ ε n vp代入式14变为 Δ σ n [ D ∧ n ][ B ]Δdn-ε n vpΔtn 16 式中 [ D ∧ n ] I [ D ][ Cn ] -1[ D ] [ D ]-1 [ Cn ] -1 17 在式15和式16中,[B]为应变矩阵.当用相关联的粘塑性法则时,[ Cn]为对称矩阵;而对于不相关 联的情况,[ Cn]是非对称的,进行分析时需要求解非对称方程. 1. 5 平衡方程求解 在任一瞬时tn都要满足平衡方程 ∫Ω[ B ] TσndΩ f n 018 式中f n 是由作用的面力、 体力、 热荷载等所产生的等效结点荷载矢量.在时间增量过程中,必须满足由式 18的增量形式所给出的平衡方程,即 ∫Ω[ B ] T Δ σ ndΩ Δf n ≠019 式中Δf n 表示载荷在时间间隔Δtn内的变化.工程中大量问题的载荷增量,是按不连续的时间间隔施加 的.除了第一次时间间隔有增量变化外,对于其他时间间隔Δf n 0. 用式12可以算出在时间步长Δtn的位移增量为 Δd [ Kn t] -1ΔVn 20 ΔVn∫ Ω[ B ] T[ D ∧ n ]ε n vpΔtndΩΔf n 21 式中[ Knt]是具有下列形式的切线刚度矩阵 [ Knt] ∫ Ω[ B ] T[ D ∧ n ][ B ]dΩ22 而Δ Vn称为伪载荷增量.当把位移矢量回代到式11时,给出应力增量,因此有 σn1σ n Δ σ n dn1 dnΔdn 23 Δ ε n vp [ B ] nΔdn - [ D ]- 1Δ σn 24 于是有 ε n1 vp ε n vpΔ ε n vp 25 用应变率的检查方法,可以检查是否已达到了静态或稳态条件.即在每一时间间隔要计算 εvp,直到这 个量小到容许的数值时就可中止该时间步长内的迭代. 1. 6 时间步长的选择 为了保证时间步进法在数值计算上是稳定的,以及保证在任一步骤中解的精度,有必要对时间步长加 以限制.本文采用下述经验法确定. 每个时间间隔的时间步长可以是常数也可以是变数.若用变步长时,时间步长的大小是用一个因子来 791 第2期安关峰,等弹粘塑性3DΟFEM在地基蠕变沉降预测中的应用 控制的,它限定最大的等效粘塑性应变增量在总有效应变中所占的部分,因此 Δ ε n vp 2 3 ε n ijvp ε n ijvpΔtn≤τ ε n 26 上述方法中的另一中方案则是按照 ε n ii vpΔt n≤τ ε n ii 27 来限定时间步长的.其中是ε n ii总应变的第一不变量,而 ε n ijvp是粘塑性应变率的第一不变量.对于等参元, 所有应变值只在高斯积分点上进行计算,因此,必须对每一高斯积分点算出满足式26或27 的Δ tn,并 把最小值用于分析.对于显式时间步进法时间增量参数τ的值,在0.01 τ 0.15的范围内限能得到精 确的结果.对于隐式法τ在10以内时结果是稳定的. 当用变时间步进法时,还可采用另一种有用的方法,即按照Δtn 1≤kΔtn方法来限制任意两区间之 间时间步长的变化.式中k为给定的常数,经验表明取k 1. 5值是合适的. 2 软件的应用 为了应用所研制的三维有限元计算软件,本文采用文献[2]提供的资料对地基沉降进行预测来检验该 软件的可靠性. 2. 1 工况 上海华盛路高层住宅大楼高12层,有地下室两层,采用天然地基与补偿式箱形基础.基础平面尺寸为 57. 6 m14. 3 m ,基底标高- 6. 45 m ,基底附加压力5. 0 t/ m2.基底下面压缩层平均模量为3 600 kPa. 2. 2 计算几何模型 考虑基础附加应力的水平影响范围为2B B为基础宽度的区域,故水平方向影响范围为288 m 71. 5 m.但考虑到基础水平向的对称形式,故计算时依据对称性取上述区域的四分之一,即144 m35. 75 m.在竖直方向上,参照文献[3]中对地基压缩层厚度的推测,考虑地表25 m内以淤泥质土为主的区域,所 以计算范围取144 m35. 75 m25. 0 m. 2. 3 基本物理参数 土体的计算参数如下[2 ,4]E 3 600 kPa ;μ 0. 4 ;ρ 1. 8 kg/ m3;C 15. 2 kPa ;Φ 17. 7;γ 0. 000 05/ kPad .E ,ρ, C ,Ф值由常规实验实测,μ取经验值,γ由单轴蠕变实验测得. 基础计算参数 E 30 000 MPa ;μ 0.2;ρ 2.4 kg/ m3. 2. 4 计算结果与分析 本次剖分采用八节点六面体等参单元,得1 300个节点,972个单元.其中基础单元24个,土体单元 948个. 图2 基底沉降曲面图 Fig. 2 Settlement surface of foundation bottom 2. 4. 1 应力分布 考虑到施工时建筑荷载逐渐施加, 故计算时采用分级加载方式施加荷载 P 10、20、30、40、50 kPa .在荷载作 用下基底沉降曲面由图2所示,基底 - 6. 45 m水平面上的垂直方向应力 分布由图3所示.由图2可看到在基础 范围内存在一个沉降盆,沉降盆的大小 计算所得范围约为40 m30 m. 计算得最大沉降点在基础中心,在 此范围外该水平面略有隆起.计算得最 大隆起值为1. 355 171 mm.由图3可清 楚得看到基础下垂向应力分布特点.在 891 同 济 大 学 学 报第29卷 基础的角点处垂向应力最大.在X方向做基底水平面垂向应力曲线见图 4 和基底中心点垂向应力随时 间变化曲线见图5 .由于垂向应力量值较大,故图4不能明显反映出应力曲线随时间变化特点,但是它 能反映出应力随时间发生调整.图5可表明基础底面中心点应力随时间而增加的特点. 图3 基底 - 6. 45 m水平面上的垂直方向应力分布 Fig. 3 Vertical stress distribution on the horizontal level of bottom 图4 基底垂向应力随时间变化情况 Fig. 4 Vertical bottom stress2time 2. 4. 2 理论预测值和实测值比较 将基底中心点沉降随时间变化的理预测值和实测值绘制成曲线见图6.由图6可见由于在计算时考虑了施工 荷载的分布施加,计算所得曲线实测曲线比较一致.说明采用本文的本构模型和分步施加荷载方法是可靠的. 图5 基底中心点应力随时间变化 Fig. 5 Foundation bottem center stress2time 图6 基底中心沉降实测曲线与预测曲线 Fig. 6 The comparison between measured cureve and calculated curve on the bottom centor 3 结论 通过以上计算分析和实测对比得出如下结论 1在进行地基蠕变沉降计算时,为反映荷载的施加对沉降的影响,需考虑分步施加荷载. 2合理选用反映软土流变特性的本构关系是进行地基蠕变沉降预测的基础. 参考文献 [1] 谢 宁,孙 钧.上海地区饱和软粘土流变特性[J ].同济大学学报,1996 ,243 233 - 237. [2] 史玉成.上海地区软土流变特性理论分析与工程应用研究[D].上海同济大学,1990. [3] 国家建委建筑科学研究院地基所.上海华盛路高层住宅箱形基础测试报告[R].上海同济大学,1977. [4] 李希元.土体三维非线性流变属性及其在深大基坑开挖工程中的应用研究[D].上海同济大学,1996. 991 第2期安关峰,等弹粘塑性3DΟFEM在地基蠕变沉降预测中的应用