（岩土工程专业论文）黄土边坡稳定性分析的有限元方法.pdf

西安建筑科技大学硕士论文 黄土边坡稳定性分析的有限元方法 专业：岩土工程 硕士生：刘钟森 指导老师：姚仰平教授 摘要 黄土地区是滑坡泥石流多发地区之一，黄土滑坡( 崩塌) 是黄土高原特有的 常见的地质灾害。由于黄土的水敏性，当降雨时，由于雨水入渗导致土体强度指 标大幅降低，从而造成了大量的滑坡和泥石流。据了解，黄土地区发生的滑坡有 一半以上是由于降雨形成的。这给黄土地区的人民生命财产安全造成了极大的威 胁和损失，随着我国经济建设的飞速发展和西部大开发的实施，黄土滑坡受到了 更多专家和学者的关注。 由于黄土强度指标随着含水量的增加而降低，那么在降雨入渗的过程中，土 坡的稳定性与土坡中土体的含水量有直接的关系，因此研究黄土边坡稳定性与土 体含水量之间的关系显得十分必要。 目前主要的土坡稳定性分析方法还主要是传统的极限平衡分析方法，如圆弧 滑动法、条分法等方法，而这些方法往往需要大量的假定，如j a n b u 法、b i s h o p 法、f e l l e n i u s 法等、而且极限平衡分析法首先要假定土体的滑动面。近年来随着 有限单元法作为一种应用广泛的数值分析方法运用于土坡稳定性分析中，其在计 算上精确、灵活和假定少的特点使其相对于其它方法表现出了其独特的优越性。 有限单元法不但可以计算出土坡中的应力、应变分布，找出土坡中的薄弱区域， 而且与强度折减法相结合还可以得到土坡潜在破坏面的位置和形状，因此有限单 元法是一种比较理想的土坡稳定性分析方法。 本文的主要目的是自己编制出基于修正剑桥模型的有限单元程序，并利用该 程序模拟土坡的降雨入渗过程，并得出土坡中的应力、应变分布情况。利用有限 单元强度折减法得出不同坡角不同入渗深度的土坡安全系数及土坡破坏面形状， 并进一步分析降雨条件下土坡的破坏机理，为实际工程提供一定的科学依据。 关键词：滑坡；边坡稳定性；黄土：有限单元法；渗流 西安建筑科技大学硕士论文 t h ef i n i t ee l e m e n tm e t h o do nl o e s s s l o p es t a b i l i t ya n a l y s e s s p e c i a l t y ：g e o t e c h n i c a le n g i n e e r i n g g r a d u a t e ：l i uz h o n gs e n i n s t r u c t o r ：p r o f e s s o ry a oy a n g p i n g a b s t r a c t l o e s sa r e ai so n eo ft h ep l a c e sw h i c hh a v em o r el a n d s l i d e sa n dm u d r o c k f l o wp r o b l e m s l o e s s l a n d s l i d e si st h es p e c i f i cg e o l o g yd i s a s t e rc o m m o ni nt h el o e s s p l a t e a u b e c a u s eo ft h es e n s i t i v i t yt o t h ew a t e r , t h es h e a rs t r e n g t ho fl o e s sw i l lb e r e d u c e ds h a r p l yw h e nr a i n sf o rr a i n - w a t e rh a v ei n f i l t r a t e di nl o e s s t h u sc a u s eal o to f l a n d s l i d e sa n dm u d - r o c kf l o w a sk n o w n ，o n em o r eh a l fl a n d s l i d e sh a p p e n e do n l y b e c a u s eo fr a i n f a l l i tt a k e sg r e a tt h r e a ta n dd a m a g et ot h ep e o p l esl i f ea n dp r o p e r t yi n t h el o e s sr e g i o n a so u re c o n o m i cc o n s t r u c t i o nd e v e l o pr a p i d l ya n dt h ee n f o r c e m e n to f t h ep o l i c yo fw e s t e r nd e v e l o p m e n t ，t h ep r o b l e mo fl a n d l a n d s l i d eh a sa t t r a c t e dm o r ea n d m o r ee x p e r t sa n ds c h o l a r s b e c a u s eo ft h es h e a rs t r e n g t ho fl o e s sw i l lb er e d u c e da st h em o i s t u r ec o n t e n t i n c r e a s e ，t h em o i s t u r ec o n t e n tb e c o m e sad i r e c t f a c t o rt ot h el o e s s - s l o p es t a b i l i t ya s i n f i l t r a t i n g s ot h er e s e a r c ho nt h er e l a t i o no fm o i s t u r ec o n t e n ta n dt h el o e s s s l o p e s t a b i l i t yb e c o m e se x t r e m e l yn e c e s s a r y n o w a d a y s ，t h em a j o r i t yo fs l o p es t a b i l i t ya n a l y s e s t a k ei n p r a c t i c es t i l l u s e t r a d i t i o n a ll i m i te q u i l i b r i u ma p p r o a c h e ss u c ha sc i r c u l a ra r ca n a l y s i s ，s l i c em e t h o d sa n d s oo n t h e s em e t h o d sn e e dal o to fa s s u m p t i o n ss u c ha sj a n b u s ，b i s h o pa n df e l l e n l u s ， b e s i d e st h ef a i l u r es u r f a c em u s tb ek n o w nf i r s t l y r e c e n t l y , t h ef i n i t ee l e m e n tm e t h o di s u s e di ns l o p es t a b i l i t ya n a l y s e sa saw i d e l yu s e dn u m e r i c a la n a l y s e sm e t h o d i th a sm o r e a d v a n t a g et h a no t h e rs l o p es t a b i l i t ym e t h o d sf o ri t sa c c u r a c y , v e r s a t i l i t ya n df e w e r a s s u m p t i o n s w i t ht h ef m i t ee l e m e n tm e t h o d s ，w ec a l lc a l c u l a t et h es t r e s sa n ds t r a i na n d f i n dt h ew e a kz o n ei n s l o p e ，b e s i d e sw ec a l lg e tt h ep o s i t i o na n ds h a p eo fp o t e n t i a l f a i l u r es u r f a c eo fs l o p ec o m b i n ew i t ht h es t r e n g t hr e d u c e dm e t h o d s ot h ef m i t ee l e m e n t m e t h o di sam o r ei d e a lm e t h o di ns l o p es t a b i l i t ya n a l y s e s t h em a i no b j e c t i v eo ft h i sp a p e ri st ow o r ko u tt h ef i n i t ee l e m e n tp r o g r a mw h i c h m o d e lt h ep r o c e s so fi n f i l t r a t i n go fs l o p eb a s e do nt h em o d i f i e dc a m c l a ym o d e l ，t h u s i i 西安建筑科技大学硕士论文 w ec a ng e tt h ed i s t r i b u t i o no fs t r e s sa n ds t r a i nw i t h i ns l o p e ；g e tt h es h a p eo ff a i l u r e s u r f a c eo fs l o p ea n dt h es a f e t yf a c t o ra s s o c i a t e dw i t hd i f f e r e n ts l o p ea n g l ea n dd i f f e r e n t d e p t ho fi n f i l t r a t i o nw i t hf i n i t ee l e m e n ts t r e n g t hr e d u c e dm e t h o d ，t h u sf u r t h e ra n a l y z e t h ef a i l u r em e c h a n i s mo fs l o p ew h e nr a i n sa n dp r o v i d es o m ec e r t a i ns c i e n t i f i cr e f e r e n c e f o rp r a c t i c a le l l e , i n e e r i n g k e y w o r d ：l a n d s l i d e ；s l o p es t a b i l i t y ；l o e s s ；f i n i t ee l e m e n tm e t h o d ；r a i n f a l l ； i n f i l t r a t i o n 声明 本人郑重声明我羼呈交的论文是我个人在导师指导下进行的研究工作及取得 ，的研究成果。尽我所知，除了文中特别加以标注和致谢的地方外，论文中不包含 其他人已经发表或撰写过的研究成果，也不包含本人或其他人在其它单位已申请 学位或为其他用途使用过的成果。与我一同工作的同志对本研究所做的所有贡献 均已在论文中做了明确的说明并表示了致谢。 申请学位论文与资料若有不实之处，本人承担一切相关责任。 论文作者签名： 刘钟客 日期：矽口多- 歹9 。 关于论文使用授权的说明 本人完全了解西安建筑科技大学有关保留、使用学位论文的规定，即：学校 j 有权保留送交论文的复印件，允许论文被查阅和借阅；学校可以公布论文的全部 内容和部分内容，可以采用影印、缩印或者其它复制手段保存论文。 ( 保密的论文在论文解密后应遵守此规定) 论文作者签名：刘李中狲导师签名：蔓牡亿p 半 日期：印一6 z 。 注：请将此页附在论文首页。 西安建筑科技大学硕士论文 1绪论 1 1黄土高边坡稳定性分析的意义 1 1 1 黄土的分布状况 世界上黄土分布面积达1 3 0 0 万平方公里，约占陆地总面积的9 _ 3 。主要分 布于中纬度干旱、半干旱地区。广泛分布于大陆内部、温带荒漠和半荒漠地区的 边缘。我国黄土分布面积，根据中国科学院地质研究所和北京大学地理系所编制 的中国黄土分布图，为6 5 3 2 8 0 平方公里，占世界黄土分布总面积的4 9 ，主要分 布在北纬3 3 4 7 度，在这区域内，一般气候干燥，降雨量小，蒸发量大，属于干旱、 半干旱气候类型。黄土高原是黄土分布的主要地区。 1 1 2 黄土特性 1 ) 黄土的化学性质 黄土的化学成分主要为二氧化硅，倍半氧化物和碳酸盐类，二氧化硅存在于 由粗颗粒到胶粒的各级粒组中，为主要成分。钙镁呈现固态或液态存在于黄土中， 为主要胶结物。能被盐酸溶解的氧化铝和氧化铁在土中以镁铝胶体的形式存在， 一部分能被盐酸溶解的二氧化硅在土中以无定型硅酸盐的形式存在，都属于胶结 物质。 其中黄土中的水溶盐与土的湿化，收缩，膨胀和透水性关系密切，并影响土 的粘性和强度。由于水溶盐在黄土含水量小的时候是胶结物质，而当含水量增大 时被水溶解，从而土体强度降低，力学性能下降，土体有可能由稳定状态进入不 稳定状态。黄土中有机质含量为0 0 2 一2 0 0 ，平均为0 。6 4 。在各级粒组中的 含量随粒径减小而增加。有机物持水性强，表面能大，常能与二价钙离子相结合 而产生凝聚现象。其存在形态或聚集于大孔壁上，或分散于粘粒中。当其呈现分 散分布时，则成为土的胶结成分。而一旦受水浸湿时会吸收大量水分而使土体崩 解。 2 ) 黄土的水敏性 黄土的水敏性是指其对水的特殊敏感性，即天然低湿度下具有明显高强度和 低压缩性的黄土，在一旦浸水甚至增湿时会发生强度大幅度骤降和变形大幅度突 增的特性。它们在定量上的不可忽视性和在定性上的急速发展性，是黄土变形影 响其上建筑物稳定性的两大突出问题。 西安建筑科技大学硕士论文 ( 1 ) 黄土的水敏性来源于黄土特殊的地域与生成条件造成的粉粒性、富盐性、 大孔性、欠压密性、非饱和性和各向异性，即低的湿度，低的密度及脆弱的架空结 构。这些方面在实际上的差异使黄土反映出颇不相同的水敏性，即颇不相同的湿 陷性，湿陷敏感性和湿剪强度的弱变性。 ( 2 ) 当黄土的含水量增大时，强度的降低与变形的增大是同一个本质反映出 的不同侧面。对承受一定剪应力达到变形稳定的黄土试样逐级增湿( 即增大水荷载) 时，剪切变形将连续逐级增大。如果将此所加的剪应力视为对应于破坏增湿级的 湿剪强度，则它与对此含水量级下的试样逐级增加剪应力所求得的相应抗剪强度 一样，也符合库伦公式，而数值上可能的不同反映出应力路径的影响 6 1 。因此，黄 土的抗剪强度也主要由黄土的含水量来决定，随着含水量的增大，内摩擦角，尤其 是粘聚力均逐渐降低i l l ，( 当黄土的天然含水量低于塑限时，水分变化对强度的影 响最大，试验表明当含水量由7 8 增加到1 8 2 时内摩擦角和内聚力都降低k ) 表明了在黄土土体稳定分析中按土的实际含水量状态分别采用不同强度的十分必 要性。当前的困难在于正确预计土体内原始的或水份入渗时的含水量分布状态。 ( 3 ) 黄土的变形和强度与含水量的变化有非常密切关系，因此，可在不同的 初始含水量下研究黄土的变形和强度，其变形的差值即为增湿的变形增量，强度的 差值即为增湿的强度增量。它们不同于一般土类之处，仅在于这种差值对含水量的 变化具有更加明显的可变性。利用不同含水量对应的模量值应该可以综合反映此 时黄土的变形与强度特性，求得变形的差值来计算湿陷变形以及由双线法求湿陷 系数的方法早己得到了较普遍的应用。 3 ) 黄土的结构性 ( 1 ) 黄土的结构可视为一个由单粒、集粒或凝块等骨架单元共同形成的空间 结构体系。它的单元形态( 单粒的矿物碎屑与集粒或凝块) 确定了力的传递性能和土 的变形性质，它的连接方式( 点接触、面接触) 确定了土的结构强度，它的排列方式 佚孔隙、架空孔隙，粒间孔隙) 确定了土的稳定性。单粒点接触，架空孔隙占优势 的结构，湿陷性大；集粒或凝块面接触，小粒间孔隙占优势的结构，湿陷性小1 3 8 1 。 ( 2 ) 黄土的结构性问题在其结构联结没有遭到破坏以前表现为它维持结构可 稳性的能力，它和颗粒联结的特性与稳定性有关，在结构联结遭到破坏以后表现 为结构可变性的能力，它和颗粒的排列特性与均匀性有关。对它们所反映的结构 特性做出直接的量测是一件困难很大的事情，但使土的结构性遭到破坏对其结构 可稳性和结构可变性作出定量的量测，并按它们在结构性中的作用找出合理的参 数，以反映土的综合结构势，是比较容易的。对黄土的原状土，它的重塑土和饱 和土所作的压缩曲线为基础求取综合结构势这一个土结构性参数的方法，具有自 己的合理性、灵敏性、稳定性和简易性1 3 , 4 1 。 2 西安建筑科技大学硕士论文 1 1 3 研究黄土高边坡的必要性 滑坡是我国最为严重的地质灾害之一，这类灾害在我国分布广泛，突发性强， 造成的危害十分严重。据不完全统计，今百年来死于滑坡和泥石流的人数至少在 万人以上，其中1 9 4 9 年以后的死亡人数在5 5 0 0 人以上。据初步估计，我国每年 因滑坡泥石流造成的经济损失至少达4 0 5 0 亿元，而在黄土地区的滑坡非常严重。 例如1 9 8 3 年3 月7 日的甘肃省东乡族自治县洒勒黄土滑坡摧毁了3 个村庄，造成 了2 2 0 人死亡，2 2 人重伤，埋没牲口4 0 0 余头，粮食6 5 万蚝，在国内外产生了 很大影响。根据陕西省六区四市的统计，共有滑坡4 6 6 4 处，威胁着1 6 万多人的 生命财产安全。由此可见黄土滑坡的危害是不能忽视的1 2 】。 我国现在正在实施西部大开发战略，而我国的黄土大部分都分布在中西部地 区，那么关于黄土滑坡问题是在发展国家西部的过程中所不能回避的问题。如在 修高速公路，铁路时遇到的路基护坡问题，和道路两旁的天然土坡安全性评价问 题，及降雨导致黄土滑坡对人民生命财产安全构成威胁问题。想要解决这些问题 都必须对黄土边坡稳定及滑坡机理做相应的研究。 1 2 边坡稳定性分析方法及研究的现状 1 2 1 边坡稳定分析方法综述 边坡稳定分析主要包括两个步骤：一是计算已知滑动面上的安全系数，常用 的方法有：极限平衡法、极限分析法、有限元法、可靠度法和边坡专家系统方法。 二是对许多可能的滑动面，搜索对英最小安全系数的临界滑动面。 1 ) 极限平衡法 极限平衡法是土坡稳定分析中最早采用且最为广泛的方法，主要方法有： f e l l e n i u s 法、b i s h o p 法、j a n b u 法、m o r g e n s t e m & p r i c e 法、s p e n c e r 法、s a r m a 法 等，并以条分法最为常用。 条分法是在1 9 1 6 年瑞典人p e t e t s e n 提出的，以后经过f e l l e n i u s 、t a y l o r 等人的 不断改进。他们假定土坡稳定问题是个平面应变问题，滑裂面是个圆柱面，计算 中不考虑土条之间的作用力，土坡稳定的安全系数是用滑裂面上全部抗滑力矩与 滑动力矩之比来定义的。四十年代以后，随着土力学学科的不断发展，也有不少 学者致力于条分法的改进工作。他们努力大致有两个方面：其一是着重探索最危 险滑弧位置的规律，制作数表，曲线，以减少计算工作量；其二是对基本假定作 些修改和补充，提出新的计算方法，使之更加符合实际情况。其中b i s h o p 等提出 的关于安全系数定义的改变，对条分法的发展起了非常重要的作用。和一般建筑 材料的强度安全系数相似，b i s h o p 等将土坡安全系数f s 定义为沿整个滑裂面的抗 西安建筑科技大学硕士论文 剪强度f ，与实际产生的剪应力f 之比，即 见。里 f 这不仅使安全系数的物理意义更加明确，而且使用范围更广泛，为以后非圆 弧滑动分析及土条分界面上条间力的各种考虑方式提供了有利条件。 条分法中根据已知条件列出的方程个数比未知数个数少，因此必须建立新的 条件方程。也就是引入适当的假定以减少未知数的个数，或增加方程的个数。 极限平衡法的不足之处是( 1 ) 各种极限平衡法有一个共同假定，那就是假定 土体是理想塑性材料。这一假定显然不符合土是弹塑性材料的事实。( 2 ) 把土体 作为一个剐体，按极限平衡的原则进行力的分析，假设土体破坏是个瞬间过程， 而不是连续过程。没有考虑土体本身的应力应变关系，假设土体破坏时各个部 分同时破坏。研究表明二体破坏是先局部破坏，而后才达到整体破坏，所以这一 假定也不符合土体的破坏模式。( 3 ) 计算土坡稳定性的安全系数时，假定一特定 的滑动面。由于这一假定导致安全系数是否可靠取决于假定的滑动面是否与真实 破坏滑动面相符合。因此搜索最优滑动面成了非常关键的因素。 2 ) 塑性极限分析法 1 9 5 2 年d r u c k e r 和p r a g e r 提出塑性极限分析法 2 7 1 ，其最大优点是考虑了材 料应力一应变关系，并利用极限状态时自重和外荷载所做的功等于滑裂面上阻力 所消耗的功为条件，结合塑性极限分析的上下限定理求得边坡极限荷载与安全系 数。由于塑性极限分析所得解答为波动于某一范围的模糊极值解“”，孙君实提出 滑动机构的概念，并证明了给定滑动机构的耗散功能定理，即滑坡极限分析的极 大值定理，将该定理“模糊化”得到给定滑裂面安全系数的模糊极大值定理：再 把极大的概念模糊化为滑体内力状态的模糊状态条件并构造模糊函数和模糊约束 条件，提出安全系数的模糊解集和最小模糊解集概念，从而建立土坡稳定分析“极 大中极小”问题的模糊极值理论，使长期以来条分法研究在假定多余未知函数方 面存在的随意性问题得到较好的解决1 3 9 1 。但极限分析法，对于复杂荷载、环境条 件等变化很难考虑：对于渗流等因素也较难反映。 3 ) 有限元法 有限元方法在土坡分析中的应用有，1 9 7 4 年s m i t h & h o b b s 用弹塑性模型得出 了机：0 边坡的分析结果，而且其分析所得的安全系数与t a y l o r s ( 1 9 3 7 ) 图表方 法确定的安全系数相当符合。z i e n k i e w i e z ，e t ，a 1 ( 1 9 7 5 ) 在考虑了c 、的情况下 对土坡进行分析，分析结果与条分法所得的结果对比相当一致 3 3 1 。 d v g r i f f i t h s ( 1 9 8 0 ) 把有限元法扩展到了土体几何大变形领域，且扩展到了各种 类型不同的土质中。国内学者中龚晓南，谢康和，雷晓燕等也都编制了弹塑性有 限元程序，对有限元在岩土工程中的应用做出了一定贡献。 4 西安建筑科技大学硕士论文 有限元方法有以下优点： ( 1 ) 有限元法是把土坡当成变形体，按照土体的弹塑性变形特性，计算土坡 内的应力分布，能算出土坡内的应力场和位移场分布，计算土体的整体抗滑稳定性。 如果进行逐步非线性分析，还可了解土坡的逐步破坏机理，跟踪土坡内塑性区的 开展情况。且有限元法可以确定破坏区的位置和破坏区的扩展情况。由此可以得 到临界滑裂面的位置。可见有限元不需要极限平衡分析方法中的各项力学假定和 滑裂面假定。 ( 2 ) 由于土体是不均匀材料，有限元法可以通过划分网格把土体离散化，然 后根据具体情况确定每个单元体的物理力学指标。从而更加真实的反映土体的真 实情况。 但有限元分析不能直接与稳定性评价指标建立关系，需要定义合适的安全系 数，使之在计算时能方便的利用有限元分析的结果。d u n c a n ( 1 9 9 6 ) 【1 7 】指出边坡安 全系数可以定义为使边坡刚好达到临界破坏状态时，对土的剪切强度进行折减的 程度，即定义安全系数是土的实际剪切强度与临界破坏时折减后的剪切强度的比 值。 国内外许多学者，如u g a i ( 1 9 8 9 ) 【1 8 1 ， u g a i & l e s h c h i n s k y ( 1 9 9 5 ) 1 9 】， ( i r j f f i t h s & l a n e ( 1 9 9 9 ) 3 2 1 ，d a w s o n & r o t h ( 1 9 9 9 ) i t l l ，m a n z a r i & n o u r ( 2 0 0 0 ) ， 宋二祥都对此做了进一步的研究。但强度折减有限元法只考虑了强度指标，没有 考虑土体其它参数，；参数折减并不等于强度折减，折减后的应力场并非实际的应 力场。 有限元分析方法已得到广泛应用，但有限元计算结果不能直接与稳定性评价 指标建立联系，需要根据应力场、位移场和材料的破坏准则，制定合理的安全系 数 4 ) 其它方法 随着计算机科学的飞速发展，计算技术不断进步，新兴学科的不断涌现，近 几年在边坡稳定性分析领域又出现了不少新的分析方法，如可靠度法，人工智能 方法，临界滑动面的搜索法等。但是这些方法目前还都不是很成熟，没有统一的 评判标准，还没有到实质应用阶段，还只是处在研究阶段 1 3 本文的研究方法及所做的工作 本文从黄土的水敏性特点着手，运用弹塑性有限单元法对考虑降雨入渗时的 二 二坡进行分析，得出土坡中的应力分布情况及土坡稳定性的安全系数，并得出土 坡的破坏模式。具体思路如下几点： ( 1 ) 建立黄土的强度参数与含水量的关系函数。 ( 2 ) 对不同的坡角土坡进行降雨入渗模拟，得出不同入渗深度时的土体应力， 5 西安建筑科技大学硕士论文 应变分布图及其土体单元网格变形图，找出土体的潜在破坏形式。并用强度折减 法得出土坡的安全系数。 ( 3 ) 总结出不同坡角在不同入渗深度时的土坡破坏规律，为边坡工程提供一 定的科学依据。 本文所做的主要工作如下： ( 1 ) 通过收集有关资料，得出黄土的强度参数与含水量的对应规律，并进一 步建立内摩擦角与含水量之间的函数关系和粘聚力与含水量的函数关系。 ( 2 ) 编制出有限单元网格自动划分程序。该程序简单实用，只需输入坡底长、 坡高、坡角及单元网格的行数和列数即可实现单元网格自动划分，并自动实现单 元序号编码和结点序号编码。 ( 3 ) 用f o r t r a n 9 0 语言编制出基于修正剑桥模型的弹塑性有限元程序。该程 序所采用的单元形式为等参数单元，并引入邓肯对剑桥模型的修正，使其可以更 好的应用于黄土。 ( 4 ) 利用编制的有限单元程序对坡角为3 0 、4 5 、6 0 度的土坡进行有限元分 析，并用m a t l a b 对计算结果进行可视化处理，得到土体网格变形图、土体中平均 主应力分布云图、剪应力分布云图及剪应力与平均主应力比分布云图j 并通过对 图形的分析，找出土坡中的潜在破坏区。 ( 5 ) 结合降雨时土坡中的入渗规律，用有限元强度折减法分别对坡角为3 0 、 4 5 、6 0 度的土坡进行分析，逐次对入渗深度分别为3 、6 、9 、1 2 、1 5 m 情况下的 土坡进行模拟。得到土体的网格变形图，及土坡的安全系数，并通过网格变形图 可以清楚的得出潜在破坏面的形状和位置。 ( 6 ) 结合有限元分析，探讨黄土土坡的破坏模式和破坏机理。 6 西安建筑科技大学硕士论文 2 降雨入渗及其对土体强度的影响 要建立一个能考虑降雨因素的土坡稳定性预测预报模型，首先必须考虑水分 的入渗，并考虑因水分入渗而引起的土的物理力学性状的变化，土体强度参数随 之变化，从而引起土坡稳定的改变。 对于这种情况，考虑土坡稳定性的基本思路是：运用非饱和土壤水分运动的 研究方法，求解给定入渗和蒸发边界条件下的斜坡土体的瞬态含水率分布：并假 定土体抗剪强度与饱和度之间存在一定的函数关系，据此，将瞬态含水率分布换 算为斜坡土体瞬态抗剪强度参数分布，在此基础上使用坡稳定性分析方法计算土 坡瞬态安全因素【7 1 。 2 1 二维非稳态渗流 岩土工程实践中的一个发展方向是能够进行非稳态类型问题的分析，在边坡 渗流中由于降雨入渗和蒸发会造成几乎是不断变化的水流条件。渗流的瞬态分析 在很大程度上受非饱和区条件的影响。非耦合二维非稳态渗流公式分别满足液相 和气相的连续方程。不考虑不同流体流动与土结构平衡条件之间的相互作用，考 虑平面应变问题，假定流体的流动仅在两个方向上发生，在第三个方向上的流动 忽略不计。 各向同性是指土的渗透系数不随方向而变化。这就是说，在土体中的某一点 上x 和y 方向的渗透系数相等。然而，非饱和土的渗透系数可能随空间变化( 即非 均匀性) ，这取决于土体中所考虑点的基质吸力的大小。 1 ) 二维非稳态渗流的液相偏微分方程如下： 式中：u 。孔隙水压力； u 。孔隙气压力； c 。( 1 一历f 册，) ( 新； c ) ，与液相偏微分方程有关的相互作用常数； m ，与法向正应力p - - d , 。) 变化有关的水的体积变化系数； 盯总应力： 坍；与基质吸力 。一) 变化有关的水的体积变化系数； 7 鳖妙川 吆 。 鳖甜 监妙， 监砂 一k监矿 w 竖弘坠打 生k 监舻 w 监孔 t 监以 西安建筑科技大学硕士论文 c 了七。( 风g m ；) ，关于液相的在x 和y 方向的固结系数； 足。液相在工和) ，方向的渗透系数( 对各向同性土) 。在土体内的某一点， 渗透系数是基底吸力的函数； p 。水的质量密度； g 重力加速度； c 。l 佃；，与重力项有关的系数。 2 ) 二维非稳态渗流的气相偏微分方程如下： 鲁= 一e 鲁軎唔等等+ 争- i 每等等泣2 ， 式中：e 与气相偏微分方程有关的相互作用常数： 坍：与法向应力( 盯一) 变化有关的空气体积变化系数： 垅；与基质吸力( 一u w ) 变化有关的空气体积变化系数； 5 饱和度： n 孑l 隙率： u 绝对空隙气压力： c ：在x 和y 方向，关于气相的固结系数： d 二在x 和y 方向( 对各向同性土) 关于气相的传导系数。 3 ) 降雨入渗时土坡水分分布规律 根据文献 1 0 考虑受降雨影响，初始体积含水率0 0 8 ，降雨强度 j = 2 5 x 1 0 c m i m i n ，持续降雨1 6 0 0 m i n 后土坡水分分布的土坡的体积含水率分布 图( 图2 1 ) 。 可以看出土坡中的等体积含水率线基本平行于土坡的自由表面，且当达到一 定的深度时，含水量不随深度的改变而变化。根据总结统计出在i i 断面处的体积 含水率与断面深度的函数关系为 拈一0 4 。坩哆，( 等) 3 + 0 4 9 9 0 ( h 0 ，( 聋) 2 - 0 3 8 8 0 ( h 0 ，( 鸶) 3 瑚眠，虹h 日一日( o ) h o l ( 2 3 ) 式中：h o 体积含水率变化临界深度 0 体积含水率 向断面深度 8 西安建筑科技大学硕士论文 矿 c l1 l ! 竺竺 i 图2 1 【1 0 】总雨量为4 0 c m 时的土坡含水率分布图 o 1 目瑚 世 聪3 0 0 4 0 0 体职台水事 oo 10 ，30 30 4 图2 2 【1 0 l 总雨量为4 0 c m 时的i i 剖面含水率沿深度分布 2 2 土体抗剪强度与含水量的关系 黄土的抗剪强度除与土的颗粒组成、矿物成分、粘粒和可溶盐含量等有关外， 主要取决于含水量和密实度，含水量越低，密实度越高则抗剪强度越大。天然含 水量低的黄土，由于存在架空结构，密度低，因而原始内聚力较小，而加固内聚 力较大。受水浸湿后产生胶溶作用以至于加固内聚力减弱甚至消失，强度降低， 那么在天然黄土边坡中会造成由于降雨入渗而产生滑坡的现象。 试验表明，当天然笛水量低于塑限时，水分变化对强度影响最大。当含水量 超过塑限时，抗剪强度降低幅度小，而超过饱和含水量时，抗剪强度变化不大。 根据文献 1 2 试验资料，黄土强度指标与含水量的关系如下。 9 西安建筑科技大学硕士论文 表2 1黄土强度与含水量的关系 含水量 w 6 78 21 0 21 3 11 5 o1 7 11 9 o2 1 12 3 o2 5 1 摩擦角 2 2 9 42 2 7 92 2 6 12 2 4 72 2 3 02 2 0 82 1 8 92 1 5 42 1 1 82 0 8 0 凝聚，j ck p a 1 l o9 07 05 2 4 2 3 42 82 31 9 1 6 ( 1 ) 根据统计结果发现，粘聚力与含水量的关系呈指数函数关系，通过相关 分析，得到粘聚力与含水量间的相关关系式为： c = 2 0 8 6 2 e - o 1 ”( 2 4 ) 式中：c 粘聚力，k p a w 含水量， 图2 3 粘聚力与含水量关系图图2 4 内摩擦角与含水量关系图 ( 2 ) 根据统计结果发现，摩擦角与含水量的关系呈一次线性函数关系，通过相 关分析，得到摩擦角与含水量间的相关关系式为： = - 0 1 l o w + 2 3 7 9 3( 2 5 ) 式中：舻摩擦角 w 含水量， 1 0 西安建筑科技大学硕士论文 3 有限单元强度折减方法 3 1有限单元方法介绍 有限元法的基本思路是将连续的求解区域离散为一组有限个、且按一定方式 相互连接在一起的单元的组合体。由于单元能按不同的联结方式进行组合，且单 元本身又可以有不同的形状，因此可以模型化几何形状复杂的求解域。有限单元 作为数值分析方法的另一个重要特点是利用在每一个单元内假设的近似函数来分 片地表示全求解区域上待求的未知场函数。单元内的近似函数通常由未知场函数 或及其导数在单元的各个结点的数值和其插值函数来表达。这样一来，一个问题 的有限元分析中未知场函数或及其导数在各个结点上的数值就成为新的未知量， 从而使一个无限自由度问题变成离散的有限自由度的问题。一经求解出这些未知 量，就可以通过插值函数计算出各个单元内场函数的近似值，从而得到整个求解 区域上的近似解。显然随着单元数目的增加，也即单元尺寸的缩小，或者随着单 元自由度的增加及插值函数精度的提高，解的近似程度将不断改进。如果单元是 满足收敛要求的，近似解最后将收敛于精确解。 从应用数学角度来看，有限元基本思想的提出，可以追溯到c o u r a n t 在1 9 4 3 年的工作，他第一次尝试应用定义在三角形区域上的分片连续函数和最小位能原 理相结合，来求解s t v c n a n t 扭转问题。一些应用数学家、物理学家和工程师由于 各种原因都涉足过有限单元的概念。但是直到1 9 6 0 年以后，随着电子数值计算机 的广泛应用和发展，有限单元的发展速度才显著加快。 现代有限单元方法的第一个成功的尝试，是将刚架位移法推广应用于弹性力 学平面问题，这是t u r n e r ，c l o u g h 等人在分析飞机结构时于1 9 5 6 年得到的结果。 他们第一次给出了用三角形单元求得平面应力问题的正确解答。三角形单元的单 元特性是由弹性理论方程来确定的，采用的是直接刚度法。他们的研究工作打开 了利用计算机求解复杂平面弹性问题的新局面。1 9 6 0 年c l o u g h 进一步处理了平面 弹性问题，并第一次提出“有限单元法”的名称，使人们开始认识了有限单元法 的功效。 三十多年来，有限单元法的理论和应用都得到迅速的、持续不断的发展。 从确定单元特性和建立求解方程的理论基础和途径来说，正如上面所提到的， t u r n e r ，c l o u g h 等人开始提出有限单元法时是利用冈4 度法，它源于结构分析的网h 度 法，所以只能处理一些比较简单的实际问题。1 9 6 3 1 9 6 4 年，b e s s e l i n g ，m e l o s h 和 j o n e s 等人证明了有限单元法是基于变分原理的里兹( r i t z ) 法的另一种形式，从 而使里兹法分析的所有理论基础都适用于有限单元法，确认了有限单元法是处理 连续介质问题的一种普遍方法。利用变分原理建立有限元方程和经典的里兹法的 西安建筑科技大学硕士论文 主要区别是有限单元法假设的近似函数不是在全求解区域而是在单元上规定的， 而且实现不要求满足任何边界条件，因此它可以用来处理和复杂的连续介质问题。 从6 0 年代后期开始，进一步利用加权余量法来确定单元特性和建立有限元求解方 程。有限单元中所利用的主要是伽辽金( g a l e r k i n ) 法。它可以用于已经知道问题 的微分方程和边界条件、但是变分的泛函尚未找到或者根本不存在的情况，因而 进一步扩大了有限单元方法的应用领域。 有限单元法的应用已经由弹性力学平面问题扩展到空间问题、板壳问题，由 静力平衡问题扩展到稳定问题、动力问题和波动问题。分析对象从弹性材料扩展 到塑性、粘弹性、粘塑性、和复合材料等。随着现代力学计算数学和计算机技术 等学科的发展，有限单元法作为一个具有巩固理论基础和广泛应用效力的数值分 析工具，必将有广阔的应用前景。 在力学分析中，从选择基本未知量的角度来看，有限单元法可以分为位移法、 力法和混合法三类，它们分别是以结点位移、结点力或部分结点力做为未知量。 就目前而言，位移法在有限单元法中应用最广。现以位移法为例介绍有限单元方 法。 3 1 1 结构离散化 把结构离散化为有限元的第一步，这相当于用一个具有有限自由度数目的系 统代替具有无限自由度系统。离散后的单元体对于二维有三角形、四边形等。对 于三维的有四面体、六面体等。 在把离散化的过程中，要把结点进行编号，以便在合成总体刚度矩阵。由于 刚度矩阵【k 】往往有数百阶甚至数千阶，而且噼】中的非零元素集中在主对角线两 侧，考虑到矩阵是对称矩阵，所以用半宽带存储形式。这就要求相邻结点的序号 相错尽量的小。 3 1 2 位移函数 有限元方法的基本思路是分段逼近，即把感兴趣的区域分为许多小区域，然 后对每一个子区域用简单函数近似求解，最后得到复杂问题的解。因此，最必要 的步骤是为每个单元的解选择一个简单函数用以表示单元内位移场，这一函数 就就是位移函数。位移函数一般用多项式表示，用符号【l 表示。 3 、t 3 单元剐度矩阵及单元结点荷载 有限单元法的核心部分就是建立单元刚度矩阵，有了单元刚度矩阵就可以得 到平衡方程组，就可以进一步求解结点位移，各个单元的应力应变等。在单元刚 度矩阵的建立过程中要用到能量原理、单元位移、单元应力一应变关系。 西安建筑科技大学硕士论文 1 ) 能量原理 由于在各种能量原理中虚位移原理应用最方便，本文就以虚位移原理为例介 绍能量原理。 所谓的虚位移原理可以是任何无限小的位移，它在结构内部必须是连续的，在 结构的边界上必须满足运动学边界条件。如果物体受到f t ，f ：，f ，等的作用。 记忙 = ，t ，：，3 ，。r 在这些外力作用下，物体的应力为 盯) = q ，q 吒，勺k ，吃】7 ( 3 1 ) 假设物体发生了虚位移，在外力作用处与各个外力相应方向上的虚位移为 d i ，6 ；，6 知 6 + ) 。 茸，西，一】7 ( 3 2 ) 由虚位移所产生的虚应变为 ，s ；，砖t y 二 ( 3 3 ) 在产生虚位移时，外力已经作用于物体，而且在虚位移过程中，外力保持不变， 因此，外力在虚位移t - 所做的虚功是 6 y = f + f z + f ，6 ；+ 6 ) 7 f ， q + q s ；+ 呸s ：+ b r 二+ k y 二+ 吃r 二一 s ) ro ) 6 u = 删s 7 料妣纰 7 f 一删s 丫 a d x d y d z 单元结点位移用 唾) 表示结点f 的位移，对于空间问题有 叫甜 ( 3 4 ) ( 3 5 ) ( 3 6 ) 那么在虚位 ( 3 7 ) ( 3 8 ) 曲安建筑科技大学硕士论文 式中心，v i ，w f ，分别是结点f 沿x ，y ，z 方向的位移分量。用 6 。表示单元e 全部结点位 移所构成的向量： 纠。一h 4 ，屯，】7 ( 3 9 ) 单元内任点( x ，y ，z ) 的位移可以表示如下： h = 坩一m m ，以，】 ( 3 1 0 ) 式中m 等为形函数。 3 ) 单元应变与应力 当单元内任一点的位移已知时，通过适当的运算可求出单元内任一点的应变， 一般可表示如下： 斜= 曰 。 式中 8 为形变转换矩阵，对于四边形型平面应变单元其形式为 b ； 根据虎克定理，单元应力可表示如下： 仃) - z ) l , 式中 d 为弹性矩阵，将应变代入上式得 ( 3 1 1 ) ( 3 1 2 ) d r 一【d j 曰】p ra 吲p ( 3 1 3 ) 式中 s = d 【b 1 是应力矩阵。 4 ) 结点力与单元刚度矩阵 为了能用结构力学方法求解连续介质的应力，以作用于单元结点上的等效集 中力代替分布于单元边5 上的应力，并称为结点力，用 e 1 表示f 点的结点力。结 点力的个数和方向必须与结点位移 4 ) 保持一致，对于空间应力问题 1 4 竖秒塑扛 盟l 馨 盟砂 盟秒盟毽 丛缸 盟妙 坐妙堕越 盟缺 墼吵 盟秒盟毽 坐缸 盟砂 。鼢 式中配， k ，嘶分别为沿着工，y ，2 方向作用于结点f 的结点力。用仁