第４５卷第８期 ２０１２年８月 天津大学学报 Ｊｏｕｒｎａｌ ｏｆ Ｔｉａｎｊｉｎ Ｕｎｉｖｅｒｓｉｔｙ Ｖｂ１．４５ ＮＯ．８ Ａｕｇ．２０１２ 降雨入渗条件下边坡湿陷变形的数值模拟 王柳江，刘斯宏，李 卓，白福清 （河海大学水利水电学院，南京２１００９８） 摘要：降雨入渗是边坡发生失稳滑动破坏的主要因素之一 为研究降雨入渗条件下边坡湿陷变形的基本特性，采 用修正的巴萨罗那模型（ＢＢＭ）和非饱和多孔介质力学理论，建立了基于弹垫｝生分析的非饱和流固耦合计算模型，开 发相应的有限元程序，进行降雨入渗下非饱和土边坡渗流场和应力场耦合的数值分析，研究了降雨持时、降雨强度 以及饱和渗透系数对边坡渗流和变形的影响．结果表明：该计算模型能够很好地反映边坡的湿陷变形特性，最大湿 陷变形发生在坡面拐角及坡顶位置．降雨强度和土体饱和渗透系数的大小对边坡内孔压和位移影响显著：在降雨强 度小于饱和渗透系数时，坡内孔压及位移随降雨强度的增大而增大；当饱和渗透系数小于且接近降雨强度时，坡内 孔压和位移增量明显增大；而当饱和渗透系数远大于或小于降雨强度时，坡内孔压及位移变化量则较小． 关键词：降雨入渗；非饱和土边坡；湿陷变形；流固耦合；修正ＢＢＭ模型 中图分类号：ＴＵ４４３ 文献标志码：Ａ 文章编号：０４９３ ２１３７（２０１２）０８—０７１４—０９ Ｎｕｍｅｒｉｃａｌ Ｓｉｍｕｌａｔｉｏｎ ｏｎ Ｗｅｔｔｉｎｇ Ｄｅｆｏｒｍａｔｉｏｎ ｏｆ Ｕｎｓａｔｕｒａｔｅｄ Ｓｌｏｐｅ Ｕｎｄｅｒ Ｒａｉｎｆａｌｌ Ｉｎｆｉｌｔｒａｔｉｏｎ ＷＡＮＧ Ｌｉｕ－ｊｉａｎｇ，ＬＩＵ Ｓｉ－ｈｏｎｇ，ＬＩ Ｚｈｕｏ，ＢＡＩ Ｆｕ－ｑｉｎｇ （Ｃｏｌｌｅｇｅ ｏｆ Ｗａｔｅｒ Ｃｏｎｓｅｒｖａｎｃｙ ａｎｄ Ｈｙｄｒｏｐｏｗｅｒ Ｅｎｇｉｎｅｅｒｉｎｇ，Ｈｏｈａｉ Ｕｎｉｖｅｒｓｉｔｙ，Ｎａｎｊ ｉｎｇ ２ １ ００９８，Ｃｈｉｎａ） Ａｂｓｔｒａｃｔ：Ｒａｉｎｆａｌｌ ｉｎｆｉｌｔｒａｔｉｏｎ ｉｓ ｏｎｅ ｏｆ ｔｈｅ ｒｅａｓｏｎｓ ｆｏｒ ｓｌｏｐｅ ｆａｉｌｕｒｅｓ．Ｉｎ ｏｒｄｅｒ ｔｏ ｉｎｖｅｓｔｉｇａｔｅ ｔｈｅ ｗｅｔｔｉｎｇ ｄｅｆｏｒｍａｔｉｏｎ ｏｆ ｓｌｏｐｅ ｕｎｄｅｒ ｒａｉｎｆａｌｌ ｉｎｆｉｌｔｒａｔｉｏｎ，ａｎ ｕｎｓａｔｕｒａｔｅｄ ｆｌｏｗ—ｄｅｆｏｒｍａｔｉｏｎ ｃｏｕｐｌｉｎｇ ｍｏｄｅｌ ｗａｓ ｐｒｅｓｅｎｔｅｄ ｂａｓｅｄ ｏｎ ｔｈｅ ｉｍ— ｐｒｏｖｅｄ ｂａｓｉｃ Ｂａｒｃｅｌｏｎａ ｍｏｄｅｌ ａｎｄ ｔｈｅ ｍｅｃｈａｎｉｃａｌ ｔｈｅｏｒｙ ｏｆ ｕｎｓａｔｕｒａｔｅｄ ｐｏｒｏｕｓ ｍｅｄｉａ．Ｎｕｍｅｒｉｃａｌ ｓｉｍｕｌａｔｉｏｎｓ ｗｅｒｅ ｃｏｎｄｕｃｔｅｄ ｔｏ ｅｘａｍｉｎｅ ｔｈｅ ｓｅｅｐａｇｅ ｆｉｅｌｄ ａｎｄ ｗｅｔｔｉｎｇ ｆｉｅｌｄ ｏｆ ｕｎｓａｔｕｒａｔｅｄ ｓｌｏｐｅ ｄｕｒｉｎｇ ｒａｉｎｆａｌｌ ｉｎｆｉｌｔｒａｔｉｏｎ．Ｔｈｅ ｉｎｆｌｕ－ ｅｎｃｅｓ ｏｆ ｒａｉｎｆａｌｌ ｄｕｒａｔｉｏｎ，ｒａｉｎｆａｌｌ ｉｎｔｅｎｓｉｔｙ ａｎｄ ｓｏｉｌ ｓａｔｕｒａｔｅｄ ｐｅｒｍｅａｂｉｌｉｔｙ ｏｎ ｓｅｅｐａｇｅ ｆｉｅｌｄ ａｎｄ ｄｅｆｏｒｍａｔｉｏｎ ｗｅｒｅ ｓｔｕｄｉｅｄ．Ｔｈｅ ｓｉｍｕｌａｔｅｄ ｒｅｓｕｌｔｓ ｓｈｏｗ ｔｈａｔ ｔｈｅ ｐｒｅｓｅｎｔｅｄ ｍｏｄｅｌ ｉｓ ｒｅａｓｏｎａｂｌｅ ｉｎ ｒｅｆｌｅｃｔｉｎｇ ｔｈｅ ａｃｔｕａｌ ｂｅｈａｖｉｏｒ ｏｆ ｓｌｏｐｅ ｕｎｄｅｒ ｒａｉｎｆａｌｌ ｉｎｆｉｌｔｒａｔｉｏｎ，ａｎｄ ｔｈａｔ ｔｈｅ ｍａｘｉｍｕｍ ｗｅｔｔｉｎｇ ｄｅｆｏｒｍａｔｉｏｎ ＯＣＣＵＬＴＳ ａｔ ｂｏｔｈ ｔｈｅ ｃｏｍｅｒ ｏｆ ｓｌｏｐｅ ｓｕｒｆａｃｅ ａｎｄ ｔｈｅ ｔｏｐ ｏｆ ｓｌｏｐｅ．Ｔｈｅ ｒｅｌａｔｉｏｎｓｈｉｐ ｂｅｔｗｅｅｎ ｒａｉｎｆａｌｌ ｉｎｔｅｎｓｉｔｙ ａｎｄ ｓａｔｕｒａｔｅｄ ｐｅｒｍｅａｂｉｌｉｔｙ ａｌｓｏ ｈａｓ ａ ｇｒｅａｔ ｉｎｆｌｕｅｎｃｅ ｏｎ ｔｈｅ ｖａｒｉａｔｉｏｎ ｏｆ ｐｏｒｅ－ｗａｔｅｒ ｐｒｅｓｓｕｒｅ ａｎｄ ｄｉｓｐｌａｃｅｍｅｎｔｓ．Ｔｈｅ ｐｏｒｅ—ｗａｔｅｒ ｐｒｅｓｓｕｒｅ ａｎｄ ｄｉｓｐｌａｃｅｍｅｎｔｓ ｉｎｃｒｅａｓｅ ｗｉｔｈ ｒａｉｎｆａｌｌ ｉｎｔｅｎｓｉｔｙ ｗｈｅｎ ｔｈｅ ｒａｉｎｆａｌｌ ｉｎｔｅｎｓｉｔｙ ｉｓ ｇｒｅａｔｅｒ ｔｈａｎ ｔｈｅ ｓａｔｕｒａｔｅｄ ｐｅｒｍｅａｂｉｌｉｔｙ．Ｔｈｅ ｐｏｒｅ—ｗａｔｅｒ ｐｒｅｓｓｕｒｅ ａｎｄ ｄｉｓｐｌａｃｅ— ｍｅｎｔｓ ｉｎｃｒｅａｓｅ ｍｏｒｅ ｏｂｖｉｏｕｓｌｙ ｗｈｅｎ ｓａｔｕｒａｔｅｄ ｐｅｒｍｅａｂｉｌｉｙ ａｐｐｒｏａｃｈｅｓ ｒａｉｎｆａｌｌｔ ｉｎｔｅｎｓｉｔｙ ｔｈａｎ ｔｈｏｓｅ ｗｈｅｎ ｔｈｅ ｆｏｒｍｅｒ ｉｓ ｆａｒ ｇｒｅａｔｅｒ ｏｒ ｌｅｓｓ ｔｈａｎ ｔｈｅ ｌａｔｔｅｒ． Ｋｅｙｗｏｒｄｓ：ｒａｉｎｆａｌｌ ｉｎｆｉｌｔｒａｔｉｏｎ；ｕｎｓａｔｕｒａｔｅｄ ｓｌｏｐｅ；ｗｅｔｔｉｎｇ ｄｅｆｏｒｍａｔｉｏｎ；ｆｌｏｗ—ｄｅｆｏｒｍａｔｉｏｎ ｃｏｕｐｌｉｎｇ；ｉｍｐｒｏｖｅｄ ＢＢＭ 降雨入渗是引起边坡产生滑动破坏的主要因素 之一，且非饱和土边坡的变形和稳定也一直是岩土工 大的经济损失和人员伤亡．因此，深入研究降雨引起 斜坡失稳的规律并建立定量的分析模型对于滑坡的 程界研究的热点和难点之一．由于我国分布有广泛 的土质边坡，每年降雨引起的滑坡都对国家造成了重 收稿日期： ２０１１．０３—２２；修回日期：２０１１．０５．１２． 预防有重要的指导意义． 目前，研究降雨入渗对边坡稳定性影响的手段主 基金项目： 汀苏省基础研究计划自然科学基金资助项目（ＢＫ２００９３４４）． 作者简介： 王柳江（１９８５一 ），男，博士研究生，ｌｉｕｊ ｗ＠ｙａｈｏｏ．ｃａ 通讯作者： 刘斯宏，ｓｉｈｏｎｇｌｉｕ＠ｈｈｕ．ｅｄｕ．ｃｎ． ２０１２年８月 王柳江等：降雨入渗条件下边坡湿陷变形的数值模拟 ・７１５・ 要有２类：①采用原位监测，通过采集大量的监测数 据建立相应的统计模型；②建立系统的数学物理计算 模型对降雨入渗下的边坡进行定量分析．迄今为止， 人们主要采用第２种方法对边坡稳定性进行分 析．其中，通过渗流计算和极限强度平衡法相结合分 析降雨人渗下边坡稳定性的方法被广泛使用＿１４Ｊ，然 而该法通常根据非饱和土的强度特性来反映边坡的 失稳机制而更能直观反映边坡失稳的变形规律却无 法得到．由于降雨入渗下的边坡变形实质上是典型 的非饱和土固结问题，因此，它的解决必须基于非饱 和土的渗流和应力耦合理论 曲Ｊ．Ｃｈｏ等 开发了水、 气、固三相耦合程序，且对降雨人渗下的均质土坡进 行了分析，但没有对边坡变形进行过深入研究．徐晗 等 Ｊ建立了一个考虑水力渗透系数特征曲线、土水特 征曲线和修正Ｍｏｈｒ－Ｃｏｕｌｏｍｂ破坏准则的非饱和土 流固耦合计算模型，然而所采用的非饱和土本构模型 并不是真正意义上的弹塑性模型．因此，在非饱和土 多孔介质力学理论的基础上采用非饱和土弹塑性模 型对降雨人渗下边坡的湿陷变形进行分析具有重要 意义．为了更准确地反映非饱和土的变形特性，很多 学者在弹塑性理论框架内建立了非饱和土本构模型， Ａｌｏｎｓｏ等 。Ｉ＿提出的ＢＢＭ模型目前已被广泛应用． 但该模型也存在一些缺陷：如通过ＢＢＭ模型中的 Ｌｏａｄ．Ｃｏｌｌａｐｓｅ（ＬＣ）屈服函数计算得到的湿化变形随 着围压的增大而增大，这与大量的试验结果不符；模 型采用净应力作为应力状态变量，无法很好地反映非 饱和土的水力．力学耦合特性，且在数值计算中运用 困难．为此，在ＢＢＭ模型的基础上对其修正进行建 模的思路已经得到了很多学者的关注ｌＪ Ｊ． 笔者在非饱和土多孔介质力学理论的基础上，采 用能够考虑吸力影响以及非饱和土水力．力学耦合特 性的修正ＢＢＭ模型，开发相应的有限元程序，模拟 了降雨入渗下边坡的湿陷变形，研究了降雨持时、降 雨强度以及饱和渗透系数对边坡渗流和变形的影响． １ 非饱和土弹塑性本构模型 ＢＢＭ模型最早由西班牙学者Ａｌｏｎｓｏ提出，该模 型能通过吸力变化来反映非饱和土干湿循环下的胀 缩特性．但其仅是基于特定试验结果的本构模型，尚 有一些地方不够完善，例如：ＬＣ曲线在描述高围压 下土体的湿化变形时与试验结果不符；模型中参考应 力Ｐ。的物理概念不明确，且当Ｐ 等于饱和状态下土 体的前期固结应力Ｐ：时，ＬＣ曲线在Ｐ Ｓ平面上为一 条竖直线，此时无法反映吸力变化引起的土体塑性体 应变；模型采用净应力作为应力状态变量，在描述非 饱和土的饱和一非饱和状态过渡时存在困难．为了使 其更好地适用于非饱和土的流固耦合分析，这里从应 力状态变量、压缩线斜率和ＬＣ屈服面３方面对其进 行了修正． １．１修正ＢＢＭ模型简介 修正后的ＢＢＭ模型描述如下． 屈服函数采用了修正剑桥模型的形式，即 ｆ＝ｑ 一Ｍ （Ｐ ＋ ）（ｐ：一Ｐ ）＝０ （１） 式中： 为平均主应力；ｑ为广义剪应力；Ｐ：为随吸力 增大的附加黏聚力；Ｐ：为非饱和土的前期最大固结 应力，其值随吸力的变化而变化，在Ｐ—Ｓ平面上绘制 出来即为ＬＣ屈服线． 与初始ＢＢＭ模型不同，本文中的ＬＣ曲线方程 是在高应力状态下推导得到的，当围压增大到一定值 时，所有吸力下的压缩曲线交于一点，此时土体吸力 降低不引起变形，则ＬＣ曲线函数为 ，、 ， 、—善一 ｐ，－ ｆ堕 ｆＰ （２） ｎ ｐａＩ／ 式中：Ｐ 为不同吸力下土体压缩曲线交点所对应的应 力；Ｐ：为饱和状态下土体的前期最大固结应力； （Ｏ）为饱和状态下土体的压缩系数；】；（为土体回弹系 数；Ｓ为吸力；Ｐ ｔ为大气压； 为围压等于１ ＭＰａ时 土体干化或湿化时由吸力变化引起的体积压缩系数； ２（ｓ）为非饱和状态下吸力等于Ｓ时对应的土体压缩 曲线斜率．根据Ｓｕｎ等【１ 】试验结果可知在高围压下， 非饱和土体的压缩曲线斜率随吸力的增大而增大，其 关系式为 （ ）： （０）＋ （３） Ｐａ１＋Ｓ 式中 为试验参数，由吸力趋向于无穷大时的土体压 缩曲线斜率与吸力等于０时的土体压缩曲线斜率相 减即可得到．图１为根据式（２）和式（３）在Ｐ—Ｓ平面上 绘制的ＬＣ曲线，其中参数 （０）＝０．２， ＝０．０３， Ｐ ＝１ ＭＰａ， ＝０．００２保持不变，变化参数ｐ 和 得到． 为更好地描述非饱和土模型中水力一力学的相互 作用关系，使模型更适用于有限元流固耦合分析，该 模型采用了平均骨架应力 作为应力状态变量，即 ＝ 一 （１一Ｓ） 一 “ （４） 式中： 为Ｋｒｏｎｅｃｋｅｒ符号；Ｓ为饱和度；“。和”ｒ 分别为孔隙气压力和孑Ｌ隙水压力，通常岩土工程中认为 天 津 大 土体与大气相通，则孔隙气压为大气压，“ ＝０． （ｂ）Ｐ：＝Ｏ．１ＭＰａ 图１不同参数下的ＬＣ曲线 Ｆｉｇ．１ ＬＣ ｃｕｒｖｅｓ ａｔ ｄｉｆｅｒｅｎｔ ｖａｌｕｅｓ ｏｆ ｐａｒａｍｅｔｅｒｓ Ｐｏ ａｎｄ丑 １－２应力．应变关系式推导 有限元计算过程中需采用应力一应变的增量关 系，因此进行有限元计算前需将屈服函数先进行变 换．首先，根据相关联流动法则得到塑性应变增量为 ｄｅ￣Ｐ＝ＡＯｆ （５） 。ｕ式中 为比例系数，可根据屈服函数的连续性条件获 得．将式（１）写成 ＝ ｄ ＋ ＋善（６） ｓ 其次，根据式（２）可得 誓 （７） 由于塑性体应变作为硬化参数，则屈服面上塑性 体应变增量相等，可采用饱和状态下土体的应力状态 及其应力增量计算得到，即 ｄ ： （１＋Ｐ） （８） 然后，将式（５）代入式（８），得 蒡 学 学 报 第４５卷第８期 非饱和土的应力一应变增量关系可表示为 ｄｏ－ ：Ｄ。ｄｅ ＋Ｗ。ｄｓ （１０） 式中：Ｄ。为弹性刚度矩阵；Ｗ。为与吸力对应的弹性 刚度向量；ｄｅ。为弹性应变增量，可表示为 一 ￣ｆ）Ｔ ＤＯｄｅ＋１ ３ｆｂ Ｐ＇ｏ，＋ 。警＋（ ］Ｔ ｗｅ］ （ ］ 。。 一 ａｆ ａｐ； －。 堕ａａ＇ｆ。 一 Ｃ ａａ＇） － （ ］ 。。 ａｆ 一瑟瑟 誓 丽 Ｄ ｆ 。：一 塑 ＋ 塑 ＋ｆ 、１Ｔ ｅ １ ｆ ］ Ｄ。 一 亟 望 Ａ Ｐ＋Ｓｒ ：０ （１６） ２０１２年８月 王柳江等：降雨入渗条件下边坡湿陷变形的数值模拟 ｒ／Ｋｗ－Ｃ（ｓ）］警 鲁一 ，，— 、 ３边坡降雨入渗变形的有限元模拟 ｓｋ （ ） ，Ｏｈ ｆ：０ （１７） Ｕ ｉ＼ Ｕ ｉ Ｊ 式中：刀为孔隙率；Ｋ 为水的体积压缩模量；ｈ为水 头；Ｃ（ｓ）为容水度，Ｃ＝ａ ／ ，０ 为体积含水量； 为土骨架单元体应变；露 ｓ为饱和渗透系数张量； ｋｒ（ｓ）为相对渗透系数；７ｗ为水的容重． ２．２土一水特征曲线模型 土体土一水特征曲线指基质吸力和饱和度之间的 关系，代表土体孔隙系统的持水特性．下文对文献 ［１ ５］中的边坡进行数值模拟，根据土一水特征实测值曲 线拟合情况，其结果与文献［７］介绍的土．水特征曲线 模型一致，则本文采用的饱和度与吸力之间的关系为 Ｓｒ＝Ｓｒ。＋（Ｓｍ—Ｓｒ。）口 ／（ａｓ＋６ｓ ） （１８） 式中： 为残余水饱和度； 为最大饱和度；ａ 、ｂｓ 和ｃｓ为参数．下文计算中土体水力渗透系数与吸力 之间也满足文献［７］中的水力渗透特性曲线函数，即 ｋｒ＝口 Ｉ（ａ ＋ｂｗＳ ） （１９） 式中：ｋｒ为相对渗透系数；ａ 、ｂｗ和ｃ 为拟合参 数．图２（ａ）和（ｂ）为文献［１５］中土．水特征及水力渗透 系数的实测值采用式（１ ８）和式（１９）的拟合情况． （ａ）土 水特征曲线 （ｂ）水力渗透特性曲线 图２土．水特征曲线和水力渗透特性曲线 Ｆｉｇ．２ Ｓｏｉｌ－ｗａｔｅｒ ｃｈａｒａｃｔｅｒｉｓｔｉｃ ｃｕｒｖｅ ａｎｄ ｈｙｄｒａｕｌｉｃ ｃｏｎｄｕｃ ｔｉｖｉｔｙ ｃｈａｒａｃｔｅｒｉｓｔｉｃ ｃｕｒｖｅ 为研究降雨入渗下非饱和土边坡的变形特性，根 据上文介绍的计算原理和应力．应变本构关系，开发 了有限元计算程序，对一降雨边坡进行了数值模 拟 ． ３．１ 有限元计算模型及条件参数 图３为原型边坡断面的有限元网格，图中黑点代 表监测点，在计算结果分析中应用．其中ａ～ｅ在坡 面以下５—１０ ｍ，而厂和ｇ分布较深，在坡面以下 ５０～１００ ｍ．其边界条件如下：底边界为位移固定和 不透水边界；左右两侧为水平向位移约束及零流量边 界；边坡表面为降雨人渗边界，人渗量的确定将在下 节介绍．初始条件包括初始应力和初始孔压，其中初 始应力根据边坡自重确定，而初始孔压的确定采用如 下方法：设定右侧水位为３５０ １ＴＩ作为常水头边界，模 拟未降雨条件下边坡渗流至稳定，将计算得到的孔压 和饱和度作为初始值，图４为初始孔压分布情况．根 据文献［１５］论述，边坡土体为砾类土，容重 ＝１６．７ ｋＮ／ｍ ，孔隙比Ｐ＝０．６２，比重 ＝２．５１，饱和 渗透系数ｋｓ＝４￣１０～ｍ／ｓ，弹性模量Ｅ＝１０ＭＰａ，泊松 比 ：０．３，黏聚力Ｃ ４．０ｋＰａ，内摩擦角 －－３６．９。， 其中土体的土一水特征曲线及水力渗透特性曲线如图 ２（ａ）和（ｂ）所示．由于文献［１５］并未针对修正ＢＢＭ模 型专门进行参数测定，因此本文通过已知参数变换以 及工程类比的方法近似获取修正ＢＢＭ模型参数．首 先，根据弹性模量 与回弹系数 间的关系式为 Ｅ： ！ 二 ！ （２０） 式中Ｐ为围压，通常在１００ １５０ｋＰａ之间．将已知参 数代入式（２０）即可得到 ，且 比 小，一般 ＝ｆ０．１～０．２１ ，则 可近似得到． 在Ｐ． 平面上临界状态线的斜率 可采用内摩 擦角表示，即 Ｍ：鱼 ３一ｓｉｎ （２１） 根据工程类比的方法，其余参数参考文献［１６１， 则得到的修正ＢＢＭ模型计算参数见表１．本文共模 拟了３０ ｈ的降雨过程，平均降雨量为１．２ ｍ／ｄ，在本 文中将该情况设为基本工况．为研究降雨强度和土 体饱和渗透系数对边坡渗流和变形的影响，本文根据 基本工况，通过改变降雨强度和饱和渗透系数，分别 ・７１８・ 天 津 大 学 学 报 第４５卷第８期 计算了饱和渗透系数保持４．０×１０ ｍ／ｓ不变，降雨 强度为０．６ ｒｎ／ｄ、１．８ ｍ／ｄ、２．４ ｍ／ｄ的情况；以及降雨强 度保持１．２ ｒｎ／ｄ不变，饱和渗透系数为１．０× １０ ｍ／ｓ、５．０×ｌ０～ｍ／ｓ、１．０×１０一 ｍ／ｓ的情况． ０ ２００ ４００ ６００ ０ ∞ ０ 加 ｘ／ｒｎ Ｏ ∞ Ｏ 加 Ｏ ｍ ０Ｉ， 图３有限元计算网格及其监测点 Ｆｉｇ．３ ＦＥＭ ｍｅｓｈ ａｎｄ ｓｅｌｅｃｔｅｄ ｎｏｄｅｓ ｆｏｒ ｒｅｓｕｌｔ ａｎａｌｙｓｉｓ 图４初始孑Ｌ压分布等值线（单位：ｋＰａ） Ｆｉｇ．４ Ｃｏｎｔｏｕｒ ｏｆ ｉｎｉｉｔａｌ ｐｏｒｅ－ｗａｔｅｒ ｐｒｅｓｓｕｒｅ ｉｎ ｓｌｏｐｅ （ｕｎｉｔ：ｋＰａ） 表１修正ＢＢＭ模型参数 Ｔａｂ．１ Ｐａｒａｍｅｔｅｒｓ ｏｆ ｉｍｐｒｏｖｅｄ ＢＢＭ ｃｏｎｓｔｉｔｕｔｉｖｅ ｍｏｄｅｌ ｌＭ ２（０） ｐ ／ｋＰａ Ｐ ／ＭＰａ Ｉ１．５ ０．１ Ｏ．Ｏ１７ Ｏ．１５ ０．００１ ５ １ＯＯ １．６５ ３．２降雨入渗边界的处理 降雨入渗水流在土体非饱和区是变化的，属于水 分通过岩土体包气带运移的两相流问题，且雨水入渗 量通常受降雨强度、降雨持时、岩土初始含水量以及 人渗面几何特征等因素的影响，因此降雨人渗是一个 复杂的边界非线性渗流问题．在渗流计算中，人渗面 通常简化为边界条件，当降雨强度大于人渗能力产生 地表径流时，为水头边界；当降雨强度小于饱和渗透 系数时，为流量边界．本文采用了以积水点为判断标 准的降雨入渗计算方法【１ ，其地表入渗能力通过地 表水头和入渗率进行判断．当地表节点非饱和时，以 地表节点为临界水头（ｈ＝０）时的人渗量为人渗能力 来判断是否由流量边界转化为水头边界；饱和时，以 地表节点的人渗率为人渗能力来判断是否由水头边 界转化为流量边界． ４计算结果及其分析 图５和图６为不同时刻边坡水平位移和沉降等 值线．由图可见，水平位移主要发生在地下水位线附 近及其以上的非饱和区，其中坡脚位置土体的水平位 移向右，而地下水位线以上的向左，最大值出现在边 ０ ∞ Ｏ 如 Ｏ 加 ０ ｍ 《≥ ｘ／ｍ （ｂ）ｔ：３０ ｈ 图５不同时刻水平位移分布等值线（单位：ｃｍ） Ｆｉｇ．５ Ｃｏｎｔｏｕｒｓ ｏｆ ｈｏｒｉｚｏｎｔａｌ ｄｉｓｐｌａｃｅｍｅｎｔ ｉｎ ｓｌｏｐｅ ａｔ ｄｉｆｅｒｅｎｔ ｉｔｍｅ（ｕｎｉｔ：ｃｍ） ０ ２００ ４００ ６００ ｘ／ｍ （ｂ）ｔ＝３０ ｈ 图６不同时刻竖向位移分布等值线（单位：ｃｍ） Ｆｉｇ．６ Ｃｏｎｔｏｕｒｓ ｏｆ ｖｅｒｔｉｃａｌ ｄｉｓｐｌａｃｅｍｅｎｔ ｉｎ ｓｌｏｐｅ ａｔ ｄｉｆｅｒｅｎｔ ｔｉｍｅ（ｕｎｉｔ：ｃｍ） Ｏ 如 ２０１２年８月 王柳江等：降雨入渗条件下边坡湿陷变形的数值模拟 坡坡面拐角处，且边坡顺坡向的水平位移随着降雨持 续而增大；同时，降雨导致了边坡内出现了沉降和隆 起，沉降主要发生在地下水位以上非饱和区，而地下 水位以下饱和区的竖向位移朝上．其中，沉降主要是 由非饱和区吸力减小导致土体强度降低软化引起，而 隆起则是由饱和区土体孔压增大导致有效应力减小 而引起的土体回弹．由图可知，最大湿陷变形发生在 坡顶及边坡坡面拐角处，且沉降值及隆起量均随降雨 持续而增大．因此，根据边坡的变形分布可知坡面拐 角和坡顶是决定该边坡稳定性的主要部位，由降雨引 起的边坡滑动面极有可能贯穿这２个部位． 图７为ａ～ｇ点的孑Ｌ压增量随时问变化情况．由 图７可知：降雨初期，ａ ｅ的孔压增量基本一致，随着 降雨持续，孔压增速逐渐减小，且与位置有关，其中 ａ、ｂ和ｃ的减小量依次增大，而ｄ和ｅ基本相同；边 坡浅层ｂ和ｄ的孔压增速明显大于深层土体厂和ｇ． 由此可知，上述点的孔压变化规律取决于土体所处的 位置．对于边坡表层土体，由于初始吸力随高度增大 而增大，在降雨条件下，吸力越低的表层土体能越快 达到饱和状态，且当表层土体趋向于饱和时，由于雨 水人渗率逐渐减小，则表面土体的孔压增大速度也相 应减小；对于深层土体单元，由于降雨入渗是雨水由 表及里流动的过程，因此，深层土体的孔压增大与消 散过程较表层土体有明显的滞后性．图８为监测点塑 性体应变的变化规律，而降雨３０ ｈ后的塑性剪应变等 值线分布如图９所示．由图可见，表层土体单元的塑 性体应变随高度的增大而减小，随着降雨的持续而增 大，塑性剪应变最大值位于坡脚地下水溢出部位，说 明塑性变形不仅与位置有关，还与降雨持续时间有 关，且由降雨引起的边坡塑性破坏区从坡脚往坡顶方 向发展．图１０为ｂ—ｅ在降雨过程中的应力路径方 向．由图１０可见，其应力路径方向朝左上方，这主要 与降雨的人渗有关．可以解释如下：降雨人渗导致表 层土体的孔压增大而使平均应力Ｐ 减小，同时，土体 图７监测点孔压增量与时间关系 Ｆｉｇ．７ Ｐｏｒｅ—ｗａｔｅｒ ｐｒｅｓｓｕｒｅ ｉｎｃｒｅｍｅｎ￣ｏｆ ｓｅｌｅｃｔｅｄ ｎｏｄｅｓ ｖｅｒｓｕｓ ｔｉｍｅ 模量降低以及土体自重增大又使广义剪应力ｇ 增大． 图８监测点塑性体应变与时间关系 Ｆｉｇ．８ Ｐｌａｓｔｉｃ ｖｏｌｕｍｅｔ￣ｃ ｓｔｒａｉｎｓ ｏｆ ｓｅｌｅｃｔｅｄ ｎｏｄｅｓ ｖｅｒｓｕｓ ｔｉｍｅ 图９ 降雨结束后塑性剪应变分布等值线（单位：％） Ｆｉｇ．９ Ｃｏｎｔｏｕｒ ｏｆ ｐｌａｓｔｉｃ ｓｈｅａｒ ｓｔｒａｉｎ ａｔ ｔｈｅ ｅｎｄ ｏｆ ｒａｉｎ— ｆａ１１（ｕｎｉｔ：％） 图１０监测点的应力路径 Ｆｉｇ．１０ Ｓｔｒｅｓｓ ｐａｔｈ ｏｆ ｓｅｌｅｃｔｅｄ ｎｏｄｅｓ 图１１～图１３为不同降雨强度下监测点孔压、沉 降和水平位移随降雨持时的变化．由于本文计算采 用的降雨强度均小于土体的饱和渗透系数（４．０× １０ ｒｎ／ｓ），因此雨水基本能够完全入渗．由图可知， 该边坡内孔压、沉降和水平位移随降雨强度的增大而 增大，且当降雨强度为２．４ｍ／ｄ时，降雨持续３０ｈ 后，深层土体ｇ点出现了正孔压，且该位置土体沉降 出现了回弹，这与降雨引起的地下水位上升且有效应 力减小有关． 图１４～图１６为不同饱和渗透系数时监测点孔 压、沉降和水平位移随降雨持时的变化．设降雨强度 保持１．２ ｍ／ｄ，改变土体饱和渗透系数．结果表明饱 ・７２０・ 天 津 大 学 学 报 第４５卷第８期 和渗透系数对边坡内孔压和变形的影响较大．由图 蛊 幽１± ｌ４可见，对于坡面浅层土体，饱和渗透系数越小，吸 力降低速度越快；而深层土体则是在饱和渗透系数接 近降雨强度时吸力降低速度较快，而当饱和渗透系数 时间／ｈ （ａ）ｂ点 时间／ｈ （ｂ）ｇ点 图１１不同降雨强度下监测点子Ｌ压与时间关系 Ｆｉｇ．１１ Ｐｏｒｅ－ｗａｔｅｒ ｐｒｅｓｓｕｒｅ ｖｅｒｓｕｓ ｔｉｍｅ ｗｉｔｈ ｄｉｆｆｅｒｅｎｔ ｒａｉｎ— ｆａｌｌ ｉｎｔｅｎｓｉｔｙ 昌 世 蜉 时间，ｈ （ａ）ｂ点 量 世 蜉 时间／ｈ （ｂ）ｇ点 图１２不同降雨强度下监测点沉降与时间关系 Ｆｉｇ．１２ Ｓｅｔｔｌｅｍｅｎｔ ｖｅｒｓｕｓ ｔｉｍｅ ｗｉｔｌ１ ｄｉｆｆｅｒｅｎｔ ｒａｉｎｆａｌｌ ｉｎｔｅｎｓｉｔｙ 日 、出１± 量 １时间，ｈ （ａ）６点 鲁 饕 （ｂ）ｇ点 图１３不同降雨强度下监测点水平位移与时间关系 Ｆｉｇ．１３ Ｈｏｒｉｚｏｎｔａｌ ｄｉｓｐｌａｃｅｍｅｎｔ ｖｅｒｓｕｓ ｔｉｍｅ ｗｉｔｈ ｄｉｆｆｅｒｅｎｔ ｒａｉｎｆａｌｌ ｉｎｔｅｎｓｉｔｙ 时间，ｈ （ａ）ｂ点 时间，ｈ （ｂ）ｇ点 图１４不同饱和渗透系数时监测点子Ｌ压与时间关系 Ｆｉｇ．１４ Ｐｏｒｅ—ｗａｔｅｒ ｐｒｅｓｓｕｒｅ ｖｅｒｓｕｓ ｔｉｍｅ ｗｉｔｈ ｄｉｆｅｒｅｎｔ ｓａｔｕ— ｒａｔｅｄ ｐｅｒｍｅａｂｉｌｉｔｙ 远大于或远小于降雨强度时，吸力减小速度较慢，该 现象与降雨强度和土体渗透系数之间的关系有关．通 常降雨初期边坡表层土体的入渗量等于降雨量，所 ２０１２年８月 王柳江等：降雨入渗条件下边坡湿陷变形的数值模拟 ・７２１・ （ａ）ｂ点 （ｂ）ｇ点 图１５不同饱和渗透系数时监测点沉降与时间关系 Ｆｉｇ．１５ Ｓｅｔｔｌｅｍｅｎｔ ｖｅｒｓｕｓ ｔｉｍｅ ｗｉｔｈ ｄｉｆｅｒｅｎｔ ｓａｔｕｒａｔｅｄ ｐｅｒ— ｍｅａｂｉｌｉｔｙ （ａ）ｂ点 （ｂ）ｇ点 图１６不同饱和渗透系数时监测点水平位移与时间关系 Ｆｉｇ．１６ Ｈｏｒｉｚｏｎｔａｌ ｄｉｓｐｌａｃｅｍｅｎｔ ｖｅｒｓｕｓ ｔｉｍｅ ｗｉＵＩ ｄｉｆｅｒｅｎｔ ｓａｔｕｒａｔｅｄ ｐｅｒｍｅａｂｉｌｉｔｙ 以根据渗流连续性条件以及达西定律可知：当表层土 体具有相同入渗量时，渗透系数越小，孔压增量越 大；对于深层土体，当渗透系数远小于降雨强度时， 由于人渗速度缓慢，导致孔压增速明显减慢，而当饱 和渗透系数远大于降雨强度时，由于土体透水性极 强，孔压消散速度同样较快，则降雨引起的孔压增量 较小．由图１５和图１６可见，当土体的饱和渗透系数 为１．０×１０～ ５．０×１０ ｍ／ｓ时，边坡的变形远大于饱 和渗透系数为４．０×１０ 和１．０×１０ ｍ／ｓ的变形，说 明当降雨强度大于且接近饱和渗透系数时，边坡发生 滑动破坏的可能性较大．通过不同饱和渗透系数下 的耦合计算表明：当边坡的渗透系数远大于或远小于 降雨强度时，边坡失稳的可能性较小，这与边坡加固 工程中设置坡内排水系统或表面隔水措施的原理 相似． ５结论 （１）降雨人渗时，最大水平位移位于坡面拐角 处，而最大湿陷变形位于坡面拐角和坡顶位置，说明 坡面拐角和坡顶是决定边坡稳定性的主要部位．则 在降雨过程中，这２个部位极有可能率先发生塑性破 坏而诱发局部滑动，因此，对该部位进行防／排水或加 固处理对边坡的稳定性具有重要意义． （２）降雨期间边坡内孔压的变化与土体位置有 关，边坡底部孔压增量较顶部小，深层较表层小．同 时，在初始应力状态接近的情况下，土体的塑性变形 与初始吸力有关，初始吸力越大，塑性变形越小，且 降雨引起的塑性破坏区从坡脚往坡顶发展． （３）在饱和渗透系数大于降雨强度的条件下，边 坡内孔压和变形的发展与降雨持时和降雨强度成 正比． （４）降雨条件下，边坡内渗流和变形的变化受土 体饱和渗透系数的影响显著，当土体的饱和渗透系数 小于且接近降雨强度时，边坡表层土体位移较大；而 当饱和渗透系数远大于或小于降雨强度时，边坡表层 土体的位移较小． （５）通过对该边坡实例的分析，说明了本文建立 的计算模型以及开发的程序能够客观地反映非饱 和土边坡在降雨入渗下的湿陷变形特性，为降雨人渗 下边坡稳定性的评价及其滑坡预测提供了新的 手段． 参考文献： ［１］吴宏伟，陈守义，庞宇威．雨水入渗对非饱和土坡稳 定性影响的参数研究［Ｊ］．岩土力学，１９９９，２０（１）： １．１４． ・７２２・ 天 津 大 学 学 报 第４５卷第８期 Ｎｇ Ｗａｎｇｗａｉ Ｃｈａｒｌｅｓ，Ｃｈｅｎ Ｓｈｏｕｙｉ，Ｐａｎｇ Ｙｕｅｗａｉ．Ｐａ— ｒａｍｅｔｒｉｃ ｓｔｕｄｙ ｏｆ ｅｆｆｅｃｔｓ ｏｆ ｒａｉｎ ｉｎｆｉｌｔｒａｔｉｏｎ ｏｎ ｕｎｓａｔｕ－ ｒａｔｅｄ ｓｌｏｐｅｓ［Ｊ１．Ｒｏｃｋ ａｎｄ Ｓｏｉｌ Ｍｅｃｈａｎｉｃｓ，１９９９， ２０（１）：１－１４（ｉｎ Ｃｈｉｎｅｓｅ）． ［２］Ｃｈｅｎ Ｈ，Ｃｈｅｎ Ｒ Ｈ，Ｙｕ Ｆ Ｃ，ｅｔ ａ１．Ｔｈｅ ｉｎｓｐｅｃｔｉｏｎ ｏｆ ｔｒｉｇｇｅｒｉｎｇ ｍｅｃｈａｎｉｓｍ ｆｏｒ ａ ｈａｚａｒｄｏｕｓ ｍｕｄｆｌｏｗ ｉｎ ａｎ ｕｒ－ ｂａｎｉｚｅｄ ｔｅｒｒｉｔｏｒｙ［Ｊ］．Ｅｎｖｉｒｏｎｍｅｎｔａｌ Ｇｅｏｌｏｇｙ，２００４， ４５（７）：８９９．９０６． １ ３ Ｊ Ｎｇ Ｃ Ｗ Ｗ，Ｐａｎｇ Ｙ ｗ．Ｉｎｆｌｕｅｎｃｅ ｏｆ ｓｔｒｅｓｓ ｓｔａｔｅ ｏｎ ｓｏｉｌ— ｗａｔｅｒ ｃｈａｒａｃｔｅｒｉｓｔｉｃｓ ａｎｄ ｓｌｏｐｅ ｓｔａｂｉｌｉｔｙ［Ｊ］．Ｊｏｕｒｎａｌ ｏｆ Ｇｅｏｔｅｃｈｎｉｃａｌ ａｎｄ Ｇｅｏｅｎｖｉｒｏｎｍｅｎｔａｌ Ｅｎｇｉｎｅｅｒｉｎｇ， ２０００．１２６（２）：１５７－１６６． ［４］吴长富，朱向荣，尹小涛，等．强降雨条件下土质边 坡瞬态稳定性分析［Ｊ］．岩土力学，２００８，２９（２）： ３８６．３９１． Ｗｕ Ｃｈａｎｇｆｕ，Ｚｈｕ Ｘｉａｎｇｒｏｎｇ，Ｙｉｎ Ｘｉａｏｔａｏ，ｅｔ ａ１． Ａｎａｌｙｓｉｓ ｏｆ ｓｏｉｌ ｓｌｏｐｅ’ｓ ｔｒａｎｓｉｅｎｔ ｓｔａｂｉｌｉｙｔ ｕｎｄｅｒ ｉｎｔｅｎｓｉｖｅ ｒａｉｎｆａｌｌ［Ｊ］．Ｒｏｃｋ ａｎｄ Ｓｏｉｌ Ｍｅｃｈａｎｉｃｓ，２００８，２９（２）： ３８６．３９１（ｉｎ Ｃｈｉｎｅｓｅ）． ［５］ 沈珠江．非饱和土简化固结理论及其应用［Ｊ］．水利水 运工程学报，２００３（４）：１－６． Ｓｈｅｎ Ｚｈｕｊｉａｎｇ．Ｓｉｍｐｌｉｉｆｅｄ ｃｏｎｓｏｌｉｄａｔｉｏｎ ｔｈｅｏｒｙ ｆｏｒ ｕｎ— ｓａｔｕｒａｔｅｄ ｓｏｉｌｓ ａｎｄ ｉｔｓ ａｐｐｌｉｃａｔｉｏｎ［Ｊ］．Ｈｙｄｒｏ Ｓｃｉｅｎｃｅ ａｎｄ Ｅｎｇｉｎｅｅｒｉｎｇ，２００３（４）：１－６（ｉｎ Ｃｈｉｎｅｓｅ）． ［６］ 陈铁林，沈珠江，杨代泉．湿陷性黄土渠道浸水变形 的数值模拟［Ｊ］．水利学报，２００５，３６（３）：３０９—３１３． Ｃｈｅｎ Ｔｉｅｌｉｎ，Ｓｈｅｎ Ｚｈｕｊｉａｎｇ，Ｙａｎｇ Ｄａｉｑｕａｎ．Ｎｕｍｅｒｉｃａｌ ｓｉｍｕｌａｔｉｏｎ ｏｎ ｗｅｔｔｉｎｇ ｄｅｆｏｒｍａｔｉｏｎ ｏｆ ｃｏｌｌａｐｓｉｂｌｅ ｌｏｅｓｓ ｄｉｔｃｈ［Ｊ］．Ｊｏｕｒｎａｌ ｏｆ Ｈｙｄｒａｕｌｉｃ Ｅｎｇｉｎｅｅｒｉｎｇ，２００５， ３６（３）：３０９—３ １３（ｉｎ Ｃｈｉｎｅｓｅ）． ［７］Ｃｈｏ Ｓ Ｅ，Ｌｅｅ Ｓ Ｒ．Ｉｎｓｔａｂｉｌｉｙｔ ｏｆｕｎｓａｔｕｒａｔｅｄ ｓｏｉｌ ｓｌｏｐｅｓ ｄｕｅ ｔｏ ｉｎｆｉｌｔｒａｔｉｏｎ［Ｊ］．Ｃｏｍｐｕｔｅｒｓ ａｎｄ Ｇｅｏｔｅｃｈｎｉｃｓ， ２００１，２８（３）：１８５—２０８． ［８］徐晗，朱以文，蔡元奇，等．降雨入渗条件下非饱 和土边坡稳定分析［Ｊ］．岩土力学，２００５，２６（１２）： １９５７．１９６２． Ｘｕ Ｈａｎ，Ｚｈｕ Ｙｉｗｅｎ，Ｃａｉ Ｙｕａｎｑｉ，ｅｔ ａ１．Ｓｔａｂｉｌｉｔｙ ｎａａｌｙｓｉｓ ｏｆ ｕｎｓａｔｕｒａｔｅｄ ｓｏｉｌ ｓｌｏｐｅｓ ｕｎｄｅｒ ｒａｉｎｆａｌｌ ｉｎｆｉｌｒｔａ ｔｉｏｎ［Ｊ］．Ｒｏｃｋ ａｎｄ Ｓｏｉｌ Ｍｅｃｈａｎｉｃｓ，２００５，２６（１２）： １９５７—１９６２（ｉｎ Ｃｈｉｎｅｓｅ）． 【９ Ｊ Ａｌｏｎｓｏ Ｅ Ｅ，Ｇｅｎｓ Ａ，Ｊｏｓａ Ａ．Ａ ｃｏｎｓｔｉｔｕｔｉｖｅ ｍｏｄｅｌ ｏｆｒ ｐａｒｔｉａｌｌｙ ｓａｔｕｒａｔｅｄ ｓｏｉｌｓ［Ｊ］．Ｇｅｏｔｅｃｈｎｉｑｕｅ，１ ９９０．４０ （３）：４０５—４３０． １ １０ ｊ Ｂｏｌｚｏｎ Ｇ，Ｓｃｈｒｅｌｆｅｒ Ｂ Ａ，Ｚｉｅｎｋｉｅｗｉｃｚ Ｏ Ｃ．Ｅｌａｓｔｏｐｌａｓ． ｔｉｃ ｓｏｉｌ ｃｏｎｓｔｉｔｕｔｉｖｅ ｌａｗｓ ｇｅｎｅｒａｌｉｚｅｄ ｔｏ ｐａｒｔｉａｌｌｙ ｓａｔｕｒａｔｅｄ ｓｔａｔｅｓ［Ｊ］．Ｇｅｏｔｈｅｃｎｉｑｕｅ，１ ９９６，４６（２）：２７９．２８９． １ １１ ｊ Ｓｈｅｎｇ Ｄａｉｃｈａｏ，Ｆｒｅｄｌｕｎｄ Ｄ Ｇ，Ｇｅｎｓ Ａ．Ａ ｎｅｗ ｍｏｄｅｌ— ｉｎｇ ａｐｐｒｏａｃｈ ｆｏｒ ｕｎｓａｔｕｒａｔｅｄ ｓｏｉｌｓ ｕｓｉｎｇ ｉｎｄｅｐｅｎｄｅｎｔ ｓｔｒｅｓｓ ｖａｒｉａｂｌｅｓ ］．Ｃａｎ Ｇｅｏｔｅｃｈ Ｊ，２００８，４５：５１１． ５３４． １ １ ２ Ｊ Ｗｈｅｅｌｅｒ Ｓ Ｊ，Ｓｉｖａｋｕｍａｒ Ｖ Ａｎ ｅｌａｓｔｏ—ｐｌａｓｔｉｃ ｃｒｉｔｉｃａｌ ｓｔａｔｅ ｆｒａｍｅｗｏｒｋ ｏｆｒ ｕｎｓａｔｕｒａｔｅｄ ｓｏｉｌｓ［Ｊ３． ＧｅｏｔｅｃｈｎｉｑＨｅ，１９９５．４５（１）：３５－５３． ［１３］缪林昌．非饱和土的本构模型研究［Ｊｊ．岩土力学， ２００７，２８（５）：８５５—８６０． Ｍｉａｏ Ｌｉｎｃｈａｎｇ．Ｒｅｓｅａｒｃｈ ｏｆ ｃｏｎｓｔｉｔｕｔｉｖｅ ｍｏｄｅｌ ｏｆ ｕｎ— ｓａｔｕｒａｔｅｄ ｓｏｉｌ［Ｊ］．Ｒｏｃｋ ａｎｄ Ｓｏｉｌ Ｍｅｃｈａｎｉｃｓ，２００７， ２８（５）：８５５．８６０（ｉｎ Ｃｈｉｎｅｓｅ）． １ １４ ｊ Ｓｕｎ Ｄｅ’ａｎ，Ｓｈｅｎｇ Ｄａｉｃｈａｏ，Ｘｕ Ｙｏｎｇｆｕ．Ｃｏｌｌａｐｓｅ ｂｅ— ｈａｖｉｏｕｒ ｏｆ ｎｕｓａｔｕｒａｔｅｄ ｃｏｍｐａｃｔｅｄ ｓｏｉｌ ｗｉｔｈ ｄｉｆｅｒｅｎｔ ｉｎｉ— ｔｉａｌ ｄｅｎｓｉｔｉｅｓ［Ｊ］．Ｃａｎ Ｇｅｏｔｅｃｈ Ｊ，２００７，４４：６７３—６８６． ［１５ ｊ ＣｈｅｎＲＨ，ＣｈｅｎＨＰ，ＣｈｅｎＫ Ｓ，ｅｔ ａ１．Ｓｉｍｕｌａｔｉｏｎ ｏｆａ ｓｌｏｐｅ ｆａｉｌｕｒｅ ｉｎｄｕｃｅｄ ｂｙ ｒａｉｎｆａｌｌ ｉｎｆｉｌｒｔａｔｉｏｎ［Ｊ３．Ｅ，２ｖｉｒｏｎ— ｍｅｎｔａｌ Ｇｅｏｌｏｇｙ，２００９，５８（５）：９４３—９５２． ［１６］Ｓｕｎ Ｄ Ａ，Ｓｈｅｎｇ Ｄ，Ｌｉ Ｘ，ｅｔ ａ１．Ｅｌａｓｔｏｐｌａｓｔｉｃ ｍｏｄｅｌ— ｌｉｎｇ ｏｆ ｈｙｄｒａｕｌｉｃ ａｎｄ ｓｔｒｅｓｓ・・ｓｔｒａｉｎ ｂｅｈａｖｉｏｕｒ ｏｆ ｕｎｓａｔｕ・－ ｒａｔｅｄ ｓｏｉｌｓ ｕｎｄｅｒ ｕｎｄｒａｉｎｅｄ ｃｏｎｄｉｔｉｏｎｓ［Ｊ］．Ｃｏｍｐｕｔｅｒｓ ａｎｄ Ｇｅｏｔｅｃｈｎｉｃｓ．２００８，３５（６）：８４５．８５２． ［１７］朱伟，程南军，陈学东，等．浅谈非饱和渗流的几 个基本问题［Ｊ］．岩土工程学报，２００６，２８（２）：２３５． ２４０． Ｚｈｕ Ｗｅｉ，Ｃｈｅｎｇ Ｎａｎｊｕｎ，Ｃｈｅｎ Ｘｕｅｄｏｎｇ，ｅｔ ａ１．Ｓｏｍｅ ｆｕｎｄａｍｅｎｔａｌ ｐｒｏｂｌｅｍｓ ｏｆ ｕｎｓａｔｕｒａｔｅｄ ｓｅｅｐａｇｅ［Ｊ］．Ｃｈｉ— ｎｅｓｅ Ｊｏｕｒｎａｌ ｏｆ Ｇｅｏｔｅｃｈｎｉｃａｌ Ｅｎｇｉｎｅｅｒｉｎｇ，２００６， ２８（２）：２３５—２４０（ｉｎ Ｃｈｉｎｅｓｅ）．
