工程科学学报，第37卷，第4期：407-413,2015年4月 Chinese Journal of Engineering,Vol.37,No.4:407-413,April 2015 DOI:10.13374/j.issn2095-9389.2015.04.002:http://journals.ustb.edu.cn 考虑二维降雨入渗的非饱和土边坡瞬态体积含水率 分析 付建新2》，宋卫东12)四，杜建华12) 1)北京科技大学金属矿山高效开采与安全教有部重点实验室，北京1000832)北京科技大学土木与环境工程学院，北京100083 ☒通信作者：E-mail:songwd@usth.cdu.cm 摘要降雨是非饱和土边坡失稳的关键因素之一，非饱和土体的强度与体积含水率有密切的关系.考虑二维降雨入渗条 件，基于非饱和土达西渗流定理及质量守恒定律建立了非饱和土的二维降雨入渗模型方程.基于交替隐式有限差分法，采用 MATLAB编制了计算程序，对方程进行了数值计算，研究了不同降雨条件、边坡内部不同位置的瞬态体积含水率分布.研究 表明，由于考虑了水平方向的入渗，边坡内随着深度和水平距离的增加，体积含水率都发生了变化，在距离边坡表面1m左右 时，体积含水率变化率变小，边坡表层的体积含水率变化远远大于边坡深部.初始体积含水率越小，边坡体积含水率变化越 大，这对边坡的稳定性不利. 关键词非饱和土边坡；边坡稳定性分析：水含量：数学模型：交替方向隐式法：有限差分法 分类号TU42 Transient volume water content analysis of unsaturated soil slopes considering two-dimensional rainfall infiltration FU Jian-xin,SONG Wei-dong,DU Jian-hua 1)Key Laboratory of High-Efficient Mining and Safety of Metal Mines (Ministry of Education),University of Science and Technology Beijing,Beijing 100083.China 2)School of Civil and Environmental Engineering.University of Science and Technology Beijing,Beijing 100083,China Corresponding author,E-mail:songwd@ustb.edu.cn ABSTRACT Rainfall is one of the key factors of unsaturated soil slope instability and the strength of unsaturated soils is closely related to volume water content.Considering two-dimensional rainfall infiltration,on the basis of the Darcy seepage theorem and mass conservation law of unsaturated soils,this paper constructed a two-dimensional infiltration model of unsaturated soils and coded a MATLAB program based on an alternating-direction-implicit (ADI)finite difference method which is used to numerical calculations to study the transient volume water content distribution at different slope locations under different rainfall conditions.The results show that the volume water content changes with the horizontal distance and depth.Moreover,the volume water content at the slope surface changes far greater than that at the deep part.The smaller the initial volume water content is,the greater the change in volume water content of the slope is,which is disadvantageous to slope stability. KEY WORDS unsaturated soil slopes;slope stability analysis;water content:mathematical models:alternating direction implicit method;finite difference method 非饱和土是一种典型的多孔介质四，由于多孔介质中空裂隙的大小、形状复杂，因此地下水在土体中的 收稿日期：2013-10-29 基金项目：国家自然科学基金资助项目(5137403)：教有部科技发展中心资助项目(20120006110022)
工程科学学报，第 37 卷，第 4 期: 407--413，2015 年 4 月 Chinese Journal of Engineering，Vol． 37，No． 4: 407--413，April 2015 DOI: 10． 13374 /j． issn2095--9389． 2015． 04． 002; http: / /journals． ustb． edu． cn 考虑二维降雨入渗的非饱和土边坡瞬态体积含水率 分析 付建新1，2) ，宋卫东1，2) ，杜建华1，2) 1) 北京科技大学金属矿山高效开采与安全教育部重点实验室，北京 100083 2) 北京科技大学土木与环境工程学院，北京 100083 通信作者: E-mail: songwd@ ustb． edu． cn 摘 要 降雨是非饱和土边坡失稳的关键因素之一，非饱和土体的强度与体积含水率有密切的关系． 考虑二维降雨入渗条 件，基于非饱和土达西渗流定理及质量守恒定律建立了非饱和土的二维降雨入渗模型方程． 基于交替隐式有限差分法，采用 MATLAB 编制了计算程序，对方程进行了数值计算，研究了不同降雨条件、边坡内部不同位置的瞬态体积含水率分布． 研究 表明，由于考虑了水平方向的入渗，边坡内随着深度和水平距离的增加，体积含水率都发生了变化，在距离边坡表面 1 m 左右 时，体积含水率变化率变小，边坡表层的体积含水率变化远远大于边坡深部． 初始体积含水率越小，边坡体积含水率变化越 大，这对边坡的稳定性不利． 关键词 非饱和土边坡; 边坡稳定性分析; 水含量; 数学模型; 交替方向隐式法; 有限差分法 分类号 TU42 Transient volume water content analysis of unsaturated soil slopes considering two-dimensional rainfall infiltration FU Jian-xin1，2) ，SONG Wei-dong1，2) ，DU Jian-hua1，2) 1) Key Laboratory of High-Efficient Mining and Safety of Metal Mines ( Ministry of Education) ，University of Science and Technology Beijing，Beijing 100083，China 2) School of Civil and Environmental Engineering，University of Science and Technology Beijing，Beijing 100083，China Corresponding author，E-mail: songwd@ ustb． edu． cn ABSTＲACT Ｒainfall is one of the key factors of unsaturated soil slope instability and the strength of unsaturated soils is closely related to volume water content． Considering two-dimensional rainfall infiltration，on the basis of the Darcy seepage theorem and mass conservation law of unsaturated soils，this paper constructed a two-dimensional infiltration model of unsaturated soils and coded a MATLAB program based on an alternating-direction-implicit ( ADI) finite difference method which is used to numerical calculations to study the transient volume water content distribution at different slope locations under different rainfall conditions． The results show that the volume water content changes with the horizontal distance and depth． Moreover，the volume water content at the slope surface changes far greater than that at the deep part． The smaller the initial volume water content is，the greater the change in volume water content of the slope is，which is disadvantageous to slope stability． KEY WOＲDS unsaturated soil slopes; slope stability analysis; water content; mathematical models; alternating direction implicit method; finite difference method 收稿日期: 2013--10--29 基金项目: 国家自然科学基金资助项目( 5137403) ; 教育部科技发展中心资助项目( 20120006110022) 非饱和土是一种典型的多孔介质［1］，由于多孔介 质中空裂隙的大小、形状复杂，因此地下水在土体中的
·408· 工程科学学报，第37卷，第4期 运动规律复杂，有些地方甚至不连续.可见，研究地下 式中，0为体积含水率，K()和K(山)分别是以土体 水渗流问题不能直接考察水质点的运动轨迹 积含水率和基质势为自变量的非饱和土导水率，业为 非饱和土是一种复杂的三相土，降雨入渗实质上 土水势. 是水分在土体饱气带中的运动6，由于非饱和土水 对于饱和土来说，导水率K为常数，而非饱和土 分渗流问题的复杂性，现有的研究往往只考虑了一维 的导水率K是土的基质势或者体积含水率的函数，并 入渗情况下土体体积含水率的变化规律0，即只考 随着体积含水率和基质势的减小而减小.将式(2)代 虑了水分的垂直入渗，而忽略了水平入渗的影响. 入式(4)中，假设重力势只与深度有关，可得 李宁和许健聪四、詹良通等四采用Fourier积分 q=-K()V(wm±z). (5) 变换分别对不同降雨强度情况下，斜长边坡的降雨入 式中“±”由z来确定，z向上为正时取“+”，z向下为 渗解析解进行了推导，并建立降雨入渗模型，这种方法 正时取“-” 不能直观反映土体含水量随雨水入渗的变化情况，也 1.3非饱和土边坡二维入渗运动方程 没有考虑到土体的渗透系数和扩散率等水分运动参数 非饱和土边坡在二维降雨入渗时包括两个过程， 随含水量的变化而对入渗过程的影响，只是以饱和的 一个是降水从地表垂直向下进入边坡的垂直入渗问 渗透系数进行研究，这与雨水非饱和土中入渗的实际 题，一个则是侧向入渗问题，如图1所示.斜坡在降雨 过程不符. 过程中，边坡表层单元A会同时有垂直通量V和水平 本文综合考虑降雨过程中，水分的垂直和水平入 通量V 渗的相互影响，研究二维入渗情况下土体瞬态体积含 水率的分布规律.不考虑空气压力、土骨架变形、温度 降雨入渗 及水中的溶质等因素，基于水分基本运动方程建立降 雨入渗模型，在此基础上进行网格化，基于交替隐式差 分法，采用MATBLAB编制计算程序，最后以工程实例 垂直通量V 进行验算，研究在不同条件下的土体瞬态体积含水率 在二维降雨入渗条件下的分布情况. 水平通量V 1非饱和土边坡二维入渗模型 1.1非饱和土的土水势 图1边坡二维降雨入渗示意图 非饱和土是由固体颗粒、水及空气构成的三相体 Fig.I Schematie diagram of two-dimensional infiltration into a slope 土体中水的流动是由土水势所决定的.为了建立非饱 取单元A为研究对象，单元体的入渗情况如图2 和土体的水分运动方程，必须要考虑土水势的影响. 所示，单元体水平和竖直分别为dx和dz,另一方向为 一般而言，土水势是指土中某点的水分状态与标 单位长度1，图中所示方向为正. 准状态的势能差团，即物理点体积元内单位数量的土 中水分具有的势能，由下式表示： 中=中。+中。+中m+.+ (1) 式中，山为总土水势，山，为重力势，中。为压力势，中m为 基质势，山，为溶质势，山，为温度势 对于非饱和土，一般不考虑压力势、溶质势及温度 势，因此非饱和土的土水势如下式所示： 少=中m土中 (2) 式中，基质势与土体体积含水率密切相关，随着土体的 饱和程度的增加而降低. 图2单元体 1.2非饱和土达西定律 Fig.2 Element Richards首先将达西定律引入了非饱和土水分运 假设土体单元A的大小不随时间改变，则时间 动规律的研究中，非饱和土中水分运动的方程可表 内，通过单元体的水量为 示为 (3) (ay+ay drd-dr.. (6) q=-K(山)7山 dx dz 或 另一方面，假设单元体内体积含水率为，则单元 q=-K(o)7山. (4) 体内水量为xdz,在单位时间d:内单元体水量变化
工程科学学报，第 37 卷，第 4 期 运动规律复杂，有些地方甚至不连续． 可见，研究地下 水渗流问题不能直接考察水质点的运动轨迹［2--5］． 非饱和土是一种复杂的三相土，降雨入渗实质上 是水分在土体饱气带中的运动［6--8］，由于非饱和土水 分渗流问题的复杂性，现有的研究往往只考虑了一维 入渗情况下土体体积含水率的变化规律［9--10］，即只考 虑了水分的垂直入渗，而忽略了水平入渗的影响． 李宁和许健聪［11］、詹良通等［12］采用 Fourier 积分 变换分别对不同降雨强度情况下，斜长边坡的降雨入 渗解析解进行了推导，并建立降雨入渗模型，这种方法 不能直观反映土体含水量随雨水入渗的变化情况，也 没有考虑到土体的渗透系数和扩散率等水分运动参数 随含水量的变化而对入渗过程的影响，只是以饱和的 渗透系数进行研究，这与雨水非饱和土中入渗的实际 过程不符． 本文综合考虑降雨过程中，水分的垂直和水平入 渗的相互影响，研究二维入渗情况下土体瞬态体积含 水率的分布规律． 不考虑空气压力、土骨架变形、温度 及水中的溶质等因素，基于水分基本运动方程建立降 雨入渗模型，在此基础上进行网格化，基于交替隐式差 分法，采用 MATBLAB 编制计算程序，最后以工程实例 进行验算，研究在不同条件下的土体瞬态体积含水率 在二维降雨入渗条件下的分布情况． 1 非饱和土边坡二维入渗模型 1. 1 非饱和土的土水势 非饱和土是由固体颗粒、水及空气构成的三相体． 土体中水的流动是由土水势所决定的． 为了建立非饱 和土体的水分运动方程，必须要考虑土水势的影响． 一般而言，土水势是指土中某点的水分状态与标 准状态的势能差［13］，即物理点体积元内单位数量的土 中水分具有的势能，由下式表示: ψ = ψg + ψp + ψm + ψs + ψt ． ( 1) 式中，ψ 为总土水势，ψg 为重力势，ψp 为压力势，ψm 为 基质势，ψs 为溶质势，ψt 为温度势． 对于非饱和土，一般不考虑压力势、溶质势及温度 势，因此非饱和土的土水势如下式所示: ψ = ψm ± ψg ． ( 2) 式中，基质势与土体体积含水率密切相关，随着土体的 饱和程度的增加而降低． 1. 2 非饱和土达西定律 Ｒichards 首先将达西定律引入了非饱和土水分运 动规律的研究中［14］，非饱和土中水分运动的方程可表 示为 q = － K( ψm ) Δ ψ ( 3) 或 q = － K( θ) Δ ψ． ( 4) 式中，θ 为体积含水率，K( θ) 和 K( ψm ) 分别是以土体 积含水率和基质势为自变量的非饱和土导水率，ψ 为 土水势． 对于饱和土来说，导水率 K 为常数，而非饱和土 的导水率 K 是土的基质势或者体积含水率的函数，并 随着体积含水率和基质势的减小而减小． 将式( 2) 代 入式( 4) 中，假设重力势只与深度有关，可得 q = － K( θ) Δ ( ψm ± z) ． ( 5) 式中，“± ”由 z 来确定，z 向上为正时取“+ ”，z 向下为 正时取“－ ”． 1. 3 非饱和土边坡二维入渗运动方程 非饱和土边坡在二维降雨入渗时包括两个过程， 一个是降水从地表垂直向下进入边坡的垂直入渗问 题，一个则是侧向入渗问题，如图 1 所示． 斜坡在降雨 过程中，边坡表层单元 A 会同时有垂直通量 Vwy和水平 通量 Vwx ． 图 1 边坡二维降雨入渗示意图 Fig． 1 Schematic diagram of two-dimensional infiltration into a slope 取单元 A 为研究对象，单元体的入渗情况如图 2 所示，单元体水平和竖直分别为 dx 和 dz，另一方向为 单位长度 1，图中所示方向为正． 图 2 单元体 Fig． 2 Element 假设土体单元 A 的大小不随时间改变，则 dt 时间 内， ( 通过单元体的水量为 Vwx x + Vwz ) z dxdzdt． ( 6) 另一方面，假设单元体内体积含水率为 θ，则单元 体内水量为 θdxdz，在单位时间 dt 内单元体水量变化 · 804 ·
付建新等：考虑二维降雨入渗的非饱和土边坡瞬态体积含水率分析 ·409 率为 形成了地表积水或者径流，地面体积含水率可接近饱 a0 dxdzdt. (7) 和体积含水率日，的值并保持不变，此时土体表面处边 at 界条件为： 由质量守恒可得 )()coa (8) R(t）cosa,z=0,0t 架k@普]+[ko警】冰0 式中，t.为降雨强度超过土的入渗能力的时刻，R()为 降雨强度，0，为饱和体积含水率. (9) (2)下边界条件.假设地下水埋藏较深，土壤深 D()和K()之间的关系为阿 处含水量恒定且分布均匀，即 D(o)=k(0) (10) 0。=0(x,z）,z=L,t>0. C(6) 式中，。为初始体积含水率，L为计算区域的下边界. 式中C(0)为比水容重，定义为 c-总 (3)右边界条件.为了简便的说明问题，假定远 (11 离边坡边界处，没有渗流，即： 则式(9)变为 D(0)0=0,x=L4:D(6)0=0,z=L 翌.(]+o]- ak.(0) 式中，L,为计算区域的右边界坐标，L,为计算区域的 下边界坐标. (12) (4)初始条件. 式(12)即为考虑二维降雨入渗的非饱和土边坡 入渗方程.式中：D(o)为水力扩散系数，cm2·min:0 0l4.=8:(x,z）,z>0,l=0 2.3二维降雨入渗的数值解法 为体积含水率；K(0)为土的导水率，cm·minl;t为时 由上述分析可知，二维降雨入渗过程是一个二阶 间，min. 偏微分方程，采用有限差分法进行求解.为保证精度 2二维降雨入渗模型数值求解 和计算效率，采用交替隐式差分法(ADI). 交替隐式差分法的具体做法是：先在时段，k+ 2.1基本假设 1]内，对x方向取中心隐式差分格式，而z方向则取显 对于二维渗流来说，水的运动方程是一个二阶偏 式差分格式，可列出差分方程并对其求解，得到时段末 微分方程，为求解方便，作如下假设： (t=k+1)的体积含水率分布：然后在时段k+1,k+ (1)假设土体为匀质土体，仅考虑扩散系数D和 2]内对：方向取中心差分分格式，而x方向取显示差 导水系数K随含水量的变化： 分，可列出另一个差分方程，得到时段末(t=k+2)的 (2)为研究方便，在初始状态时，假设土体中各处 体积含水率分布.如此交替进行，便可求出体积含水 的体积含水率处处相等： 率随降雨时间的变化 (3)假设地下水位不变，且埋藏很深，在有限的降 (1)x方向隐式及z方向显式.式(12)中各项导 雨时间内，降雨对地下水位影响不大； 数的差分近似式为： (4)降雨影响深度最大为10m,降雨历时最大为 0_- 24h; dt (13) (5)导水率表达式取K(0)=a ,扩散率取 是{po]- D(8)=b（a. ,系数通过室内实验或原位实验 D)-D( ,(14) △r2 求得： (6)不考虑温度、蒸发等影响因素 是po2]= 2.2边界条件 Dn（t1-)-D-n(,-成-，(15) (1)上边界条件.降雨强度为R()且已知，当降 43 雨开始时，土的入渗能力较大，降雨强度小于土的入渗 aK(0)i-K- (16) 能力，土的实际入渗率就是降雨强度，经历了一定的时 2△z 间后，土的入渗能力减小，降雨强度大于土的入渗率， 式中，Di2K等是D(in)、K(d)等的简写
付建新等: 考虑二维降雨入渗的非饱和土边坡瞬态体积含水率分析 率为 θ t dxdzdt． ( 7) 由质量守恒可得 θ t = Vwx x + Vwz z ． ( 8) 将式( 4) 代入式( 8) ，可得单元体二维渗流方程为 θ t = [ x Kx ( θ) ψm ] x + [ z Kz ( θ) ψm ] z － Kz ( θ) z ． ( 9) D( θ) 和 K( θ) 之间的关系为［15］ D( θ) = K( θ) C( θ) ． ( 10) 式中 C( θ) 为比水容重，定义为 C( θ) = dθ dψm ． ( 11) 则式( 9) 变为 θ t = [ x Dx ( θ) θ ] x + [ z Dz ( θ) θ ] z － Kz ( θ) z ． ( 12) 式( 12) 即为考虑二维降雨入渗的非饱和土边坡 入渗方程． 式中: D( θ) 为水力扩散系数，cm2 ·min － 1 ; θ 为体积含水率; K( θ) 为土的导水率，cm·min － 1 ; t 为时 间，min． 2 二维降雨入渗模型数值求解 2. 1 基本假设 对于二维渗流来说，水的运动方程是一个二阶偏 微分方程，为求解方便，作如下假设: ( 1) 假设土体为匀质土体，仅考虑扩散系数 D 和 导水系数 K 随含水量的变化; ( 2) 为研究方便，在初始状态时，假设土体中各处 的体积含水率处处相等; ( 3) 假设地下水位不变，且埋藏很深，在有限的降 雨时间内，降雨对地下水位影响不大; ( 4) 降雨影响深度最大为 10 m，降雨历时最大为 24 h; ( 5) 导水率表达式取 K( θ) = a1 ( θ θ ) s b1 ，扩散率取 D( θ) = b1 ( θ θ ) s b2 ，系数通过室内实验或原位实验 求得; ( 6) 不考虑温度、蒸发等影响因素． 2. 2 边界条件 ( 1) 上边界条件． 降雨强度为 Ｒ( t) 且已知，当降 雨开始时，土的入渗能力较大，降雨强度小于土的入渗 能力，土的实际入渗率就是降雨强度，经历了一定的时 间后，土的入渗能力减小，降雨强度大于土的入渗率， 形成了地表积水或者径流，地面体积含水率可接近饱 和体积含水率 θs 的值并保持不变，此时土体表面处边 界条件为: － D( θ ( ) θ x sinα + θ z cosα ) + K( θ) cosα = Ｒ( t) cosα，z = 0，0 ＜ t ＜ ta ; θ = θs，z = 0，t ＞ ta { ． 式中，ta 为降雨强度超过土的入渗能力的时刻，Ｒ( t) 为 降雨强度，θs 为饱和体积含水率． ( 2) 下边界条件． 假设地下水埋藏较深，土壤深 处含水量恒定且分布均匀，即 θ0 = θi ( x，z) ，z = L，t ＞ 0． 式中，θ0 为初始体积含水率，L 为计算区域的下边界． ( 3) 右边界条件． 为了简便的说明问题，假定远 离边坡边界处，没有渗流，即: D( θ) θ x = 0，x = Lx ; D( θ) θ z = 0，z = Lz ． 式中，Lx 为计算区域的右边界坐标，Lz 为计算区域的 下边界坐标． ( 4) 初始条件． θ | t = t0 = θi ( x，z) ，z ＞ 0，t = 0． 2. 3 二维降雨入渗的数值解法 由上述分析可知，二维降雨入渗过程是一个二阶 偏微分方程，采用有限差分法进行求解． 为保证精度 和计算效率，采用交替隐式差分法( ADI) ． 交替隐式差分法的具体做法是: 先在时段［k，k + 1］内，对 x 方向取中心隐式差分格式，而 z 方向则取显 式差分格式，可列出差分方程并对其求解，得到时段末 ( t = k + 1) 的体积含水率分布; 然后在时段［k + 1，k + 2］内对 z 方向取中心差分分格式，而 x 方向取显示差 分，可列出另一个差分方程，得到时段末( t = k + 2) 的 体积含水率分布． 如此交替进行，便可求出体积含水 率随降雨时间的变化． ( 1) x 方向隐式及 z 方向显式． 式( 12) 中各项导 数的差分近似式为: θ t = θ k + 1 ij － θ k ij Δt ， ( 13) [ x D( θ) θ ] x = Dk + 1 i + 1 /2，j ( θ k + 1 i + 1，j － θ k + 1 i，j ) － Dk + 1 i － 1 /2，j ( θ k + 1 i，j － θ k + 1 i － 1，j ) Δx 2 ，( 14) [ z D( θ) θ ] z = Dk i，j + 1 /2 ( θ k i，j + 1 － θ k i，j ) － Dk i，j － 1 /2 ( θ k i，j － θ k i，j － 1 ) Δz 2 ， ( 15) K( θ) z = Kk i，j + 1 － Kk i，j － 1 2Δz ． ( 16) 式中，Dk + 1 i + 1 /2，j 、Kk i，j + 1等是 D( θ k + 1 i + 1 /2，j ) 、K( θ k i，j + 1 ) 等的简写． · 904 ·
410 工程科学学报，第37卷，第4期 rDn+n(K-K）. (30) 454 7s △t 25=2△ 等式左边涉及的D均为z方向扩散系数的值，等 式右边涉及的D均为x方向扩散系数的值，K为z方 式(12)可变为： 向的值.同理可以得到如下的三对角矩阵 agc=y (17) 「buC 0 61 ag=rDiing hu (18) a12b12 C12 始 =-+n(Dn+D] he (19) =nDn (20) 01,-2 b1,-2 C1,u-2 h与=-rD-ln,-l+2(D-n+Din）-1]- h1,a-2 0 01,m-l b-LJ nDn+r,（K1-K-i. Lhi.-1 (21) 利用“追赶法”解出每个时段末的体积含水率。对 等式左边涉及的D均为x方向扩散系数的值，等 于i=2这种情况，原理是相同的，只是用的体积含水 式右边涉及的D均为z方向扩散系数的值，K为z方 率是上一个节点求出的体积含水率.以此类推，就可 向的值. 以求出时段末z方向各个节点的体积含水率.这样， 假设上边界条件是第一类边界条件，即体积含水 所有节点的时段末体积含水率都已经求出. 率已知，当=1时对x方向的各个节点写出差分方 当降雨强度大于入渗率时，根据原边界方程补写 程.很明显，我们可以将式(17)写成一个三对角矩阵. 差分方程，如式中的第一个公式，在边界结点j=0处 利用“追赶法”解出每个时段末的体积含水率.对 于=2这种情况原理是相同的，只是用的体积含水率 列差分方程，驰取向前差分，可得差分方程如下： z 是上一个节点求出的体积含水率.这样以此类推，就 可以求出时段末x方向各个节点的体积含水率. :"-:-Ra (31) △z Tou cu 0「1 式中，D。和K分别为z方向的参数D()和 a21 K(:),R1P表示本时段的平均供水强度.式(31) 可写成 0e-2,1 ba-2.1 Cm-2,1 ha-2.1 bo0o+cadi=ho (32) Gn-l,l B- Lhe-1.1 式中，bn=-DAz,cn=-bo,hn=-R1n+K. (2):方向隐式及x方向显式.原方程式(12)中 将上式作为方程组的第一个方程即可组成三对角 各项导数的差分近似式为： 矩阵，用同样的办法可以求得各节点的体积含水率. 0- 基于上述算法，用MATLAB编制了计算程序，可 (22) 迅速的计算出土体内瞬态体积含水率分布. [oo]= 3工程案例分析 De(du-）-2（i-2,(23) 3.1基本参数 △x2 利用本文的程序，对某路基边坡在降雨过程中土 [po0]- 体内部瞬态体积含水率的分布情况（路基边坡坡 脚55)进行如下设定： D ()-D ( 42 ,(24) (1)边坡土体中的水分运动参数为 D,=D,=10.4(0/0.)6,K,=K=0.02(0/8.)28 ak(0) -出 (25) (33) dz 2△z (2)初始体积含水率为0.28，饱和体积含水率 式中，D广nK等是D(n）、K()等的简写. 为0.45： 原方程可变为： (3)计算下边界深度取3m,右边界水平距离取距 ab+c=hg (26) 边坡边缘1.0m,计算步长为5cm,时间步长为1min. ay=rDi (27) 3.2不同工况边坡瞬态体积含水率分析 6,=-0+(Dn+Dn)], (28) (1)设初始体积含水率为0.28，保持降雨强度不 c=aD (29) 变，计算不同时间土体的瞬时体积含水率与深度的关 h=-Di+(Di-i+Di)-1]0j- 系，降雨强度分别选取200、300、400和500mmd,每
工程科学学报，第 37 卷，第 4 期 令 r1 = Δt Δx 2，r2 = Δt Δz 2，r3 = Δt 2Δz ． 式( 12) 可变为: aijθ k + 1 i － 1，j + bijθ k + 1 i，j + cijθ k + 1 i + 1，j = hij， ( 17) aij = r1Dk + 1 i － 1 /2，j ， ( 18) bij = －［1 + r1 ( Dk + 1 i + 1 /2，j + Dk + 1 i － 1 /2，j ) ］， ( 19) cij = r1Dk + 1 i + 1 /2，j ， ( 20) hij = － r2Dk i，j － 1 /2 θ k i，j － 1 +［r2 ( Dk i，j － 1 /2 + Dk i，j + 1 /2 ) － 1］θ k i，j － r2Dk i，j + 1 /2 θ k i，j + 1 + r3 ( Kk i，j + 1 － Kk i，j － 1 ) ． ( 21) 等式左边涉及的 D 均为 x 方向扩散系数的值，等 式右边涉及的 D 均为 z 方向扩散系数的值，K 为 z 方 向的值． 假设上边界条件是第一类边界条件，即体积含水 率已知，当 j = 1 时对 x 方向的各个节点写出差分方 程． 很明显，我们可以将式( 17) 写成一个三对角矩阵． 利用“追赶法”解出每个时段末的体积含水率． 对 于 j = 2 这种情况原理是相同的，只是用的体积含水率 是上一个节点求出的体积含水率． 这样以此类推，就 可以求出时段末 x 方向各个节点的体积含水率． b11 c11 0 a21 b21 c21 … … … an － 2，1 bn － 2，1 cn － 2，1 0 an － 1，1 bn － 1， 1 θ k + 1 11 θ k + 1 21 θ k + 1 n － 2，1 θ k + 1 n － 1， 1 = h11 h21 hn － 2，1 hn － 1， 1 ． ( 2) z 方向隐式及 x 方向显式． 原方程式( 12) 中 各项导数的差分近似式为: θ t = θ k + 1 ij － θ k ij Δt ， ( 22) [ x D( θ) θ ] x = Dk i + 1 /2，j ( θ k i + 1，j － θ k i，j ) － Dk i － 1 /2，j ( θ k i，j － θ k i － 1，j ) Δx 2 ， ( 23) [ z D( θ) θ ] z = Dk + 1 i，j + 1 /2 ( θ k + 1 i，j + 1 － θ k + 1 i，j ) － Dk + 1 i，j － 1 /2 ( θ k + 1 i，j － θ k + 1 i，j － 1 ) Δz 2 ，( 24) K( θ) z = Kk + 1 i，j + 1 － Kk + 1 i，j － 1 2Δz ． ( 25) 式中，Dk i + 1 /2，j 、Kk + 1 i，j + 1等是 D( θ k i + 1 /2，j ) 、K( θ k + 1 i，j + 1 ) 等的简写． 原方程可变为: aijθ k + 1 i，j － 1 + bijθ k + 1 i，j + cijθ k + 1 i，j + 1 = hij， ( 26) aij = r2Dk + 1 i，j － 1， ( 27) bij = －［1 + r2 ( Dk + 1 i，j + 1 /2 + Dk + 1 i，j － 1 /2 ) ］， ( 28) cij = r2Dk + 1 i，j + 1 /2， ( 29) hij = － r1Dk i － 1 /2，jθ k i － 1，j +［r1 ( Dk i － 1 /2，j + Dk i + 1 /2，j ) － 1］θ k i，j － r1Dk i + 1 /2，jθ k i + 1，j + r3 ( Kk + 1 i，j + 1 － Kk － 1 i，j － 1 ) ． ( 30) 等式左边涉及的 D 均为 z 方向扩散系数的值，等 式右边涉及的 D 均为 x 方向扩散系数的值，K 为 z 方 向的值． 同理可以得到如下的三对角矩阵． b11 c11 0 a12 b12 c12 … … … a1，n － 2 b1，n － 2 c1，n － 2 0 a1，n － 1 b1，n － 1 θ k + 1 11 θ k + 1 12 θ k + 1 1，n － 2 θ k + 1 1，n － 1 = h11 h12 h1，n － 2 h1，n － 1 ． 利用“追赶法”解出每个时段末的体积含水率． 对 于 i = 2 这种情况，原理是相同的，只是用的体积含水 率是上一个节点求出的体积含水率． 以此类推，就可 以求出时段末 z 方向各个节点的体积含水率． 这样， 所有节点的时段末体积含水率都已经求出． 当降雨强度大于入渗率时，根据原边界方程补写 差分方程，如式中的第一个公式，在边界结点 j = 0 处 列差分方程，θ z 取向前差分，可得差分方程如下: Dk + 1 i，0 θ k + 1 i，1 － θ k + 1 i，0 Δz － Kk + 1 i，0 = － Ｒk + 1 /2 ． ( 31) 式中，Dk + 1 i，0 和 Kk + 1 i，0 分别 为 z 方 向 的 参 数 D( θ k + 1 i，0 ) 和 K( θ k + 1 i，0 ) ，Ｒk + 1 /2 表示本时段的平均供水强度． 式( 31) 可写成 bi0 θ k + 1 i，0 + ci0 θ k + 1 i，1 = hi0 ． ( 32) 式中，bi0 = － Dk + 1 i，0 /Δz，ci0 = － bi0，hi0 = － Ｒk + 1 /2 + Kk + 1 i，0 ． 将上式作为方程组的第一个方程即可组成三对角 矩阵，用同样的办法可以求得各节点的体积含水率． 基于上述算法，用 MATLAB 编制了计算程序，可 迅速的计算出土体内瞬态体积含水率分布． 3 工程案例分析 3. 1 基本参数 利用本文的程序，对某路基边坡在降雨过程中土 体内部 瞬 态 体 积 含 水 率 的 分 布 情 况 ( 路 基 边 坡 坡 脚 55°) 进行如下设定: ( 1) 边坡土体中的水分运动参数为 Dx = Dz = 10. 4( θ /θs) 3. 6，Kx = Kz = 0. 02( θ /θs) 2. 8 ; ( 33) ( 2) 初始体积含水率为 0. 28，饱和体积含水率 为 0. 45; ( 3) 计算下边界深度取 3 m，右边界水平距离取距 边坡边缘 1. 0 m，计算步长为 5 cm，时间步长为 1 min． 3. 2 不同工况边坡瞬态体积含水率分析 ( 1) 设初始体积含水率为 0. 28，保持降雨强度不 变，计算不同时间土体的瞬时体积含水率与深度的关 系，降雨强度分别选取 200、300、400 和 500 mm·d － 1，每 · 014 ·
付建新等：考虑二维降雨入渗的非饱和土边坡瞬态体积含水率分析 411 个降雨强度下分别取120、240、480、720、960和 1440min进行计算 如图3~图6所示，当降雨强度为300mmd时， 0.46 0.44 降雨时间分别为240、480、960及1440mim时，边坡土 0.42 体内不同位置体积含水率分布图 0.36 是0.34 0.32 0.40 0.30 0.38 0.28 0.36 0.34 050100150200250300350050100150200250300350 深度cm 水平距离em 0.32 未影响区域 0.30 图6降雨时间为1440min时土体体积含水率分布 0.2 00 Fig.6 Distribution of volume water content after 1440 min rainfall 050n00 响区域越来越小，影响深度从1.0m增大到2.5m左 深度em 50300 500 50100150200250300350 水平距离/cm 右，水平影响距离从0.5m增大到1.5m.另外，随着深 度和水平距离增加，曲面上网格逐渐变小，说明降雨对 图3降雨时间为240mim时土体体积含水率分布 Fig.3 Distribution of volume water content after 240 min rainfall 土体体积含水率的影响也越小，由于在边坡表面考虑 了水平方向的入渗，因此表面的体积含水率变化比较 大.对于其他降雨强度下的体积含水率分布规律是相 0.46 同的，这里不再赘述 0.44 图7和图8为降雨时间分别为960min和1440min 0.42 30.40 时边坡土体的体积含水率的等值线分布图.最底部的 0.38 如0.36 一条是体积含水率为0.30的等值线，每条等值线之间 0.34 相隔0.02. 0.32 0.30 0.38 0 0286 0.36 -0.34 0.32 50 1001504h0050300 深度cm 50050100150200250300350 -0.30 水平距离cm 100 图4降雨时间为480min时土体体积含水率分布 150 Fig.4 Distribution of volume water content after 480 min rainfall 200 250 0.46 0 50100150 200 250 300cm 044 0.42 040 图7降雨时间为960min时体积含水率等值线分布图 038 0.36 Fig.7 Contour maps of volume water content after960 min rainfall 0.34 0.32 020 0.28 从各个降雨时间土坡内部的体积含水率等值线可以 50 00 看出：随着降雨时间变长，等值线在不断的下移，并且变得 15 深度em 200 350050100150200230300350 分散，说明降雨影响的深度在不断加深：在斜坡段部分，等 25 3 水平距离/em 值线是向下倾斜，说明了表层的体积含水率大于土体内 部的体积含水率：随着水平距离的变大，等值线变为趋于 图5降雨时间为960min时土体体积含水率分布 平行的直线，说明水平入渗的影响在逐渐变小 Fig.5 Distribution of volume water content after 960 min rainfall (2)边坡初始体积含水率对土体内部瞬态体积含 由图可知，体积含水率是随着深度的增加而逐渐 水率也有很大的影响.选取降雨强度为450mm·d', 减小，逐渐接近土体的初始体积含水率，随着降雨时间 初始体积含水率分别选0.15、0.25、0.30及0.35，降雨 的延长，影响深度及水平距离不断增加，即图3中未影 时间选取480min,距边坡边缘1.0m处作垂直剖面，研
付建新等: 考虑二维降雨入渗的非饱和土边坡瞬态体积含水率分析 个降 雨 强 度 下 分 别 取 120、240、480、720、960 和 1440 min 进行计算． 如图 3 ～ 图 6 所示，当降雨强度为 300 mm·d － 1时， 降雨时间分别为 240、480、960 及 1440 min 时，边坡土 体内不同位置体积含水率分布图． 图 3 降雨时间为 240 min 时土体体积含水率分布 Fig． 3 Distribution of volume water content after 240 min rainfall 图 4 降雨时间为 480 min 时土体体积含水率分布 Fig． 4 Distribution of volume water content after 480 min rainfall 图 5 降雨时间为 960 min 时土体体积含水率分布 Fig． 5 Distribution of volume water content after 960 min rainfall 由图可知，体积含水率是随着深度的增加而逐渐 减小，逐渐接近土体的初始体积含水率，随着降雨时间 的延长，影响深度及水平距离不断增加，即图 3 中未影 图 6 降雨时间为 1440 min 时土体体积含水率分布 Fig． 6 Distribution of volume water content after 1440 min rainfall 响区域越来越小，影响深度从 1. 0 m 增大到 2. 5 m 左 右，水平影响距离从 0. 5 m 增大到 1. 5 m． 另外，随着深 度和水平距离增加，曲面上网格逐渐变小，说明降雨对 土体体积含水率的影响也越小，由于在边坡表面考虑 了水平方向的入渗，因此表面的体积含水率变化比较 大． 对于其他降雨强度下的体积含水率分布规律是相 同的，这里不再赘述． 图7 和图8 为降雨时间分别为960 min 和1440 min 时边坡土体的体积含水率的等值线分布图． 最底部的 一条是体积含水率为 0. 30 的等值线，每条等值线之间 相隔 0. 02． 图 7 降雨时间为 960 min 时体积含水率等值线分布图 Fig． 7 Contour maps of volume water content after 960 min rainfall 从各个降雨时间土坡内部的体积含水率等值线可以 看出: 随着降雨时间变长，等值线在不断的下移，并且变得 分散，说明降雨影响的深度在不断加深; 在斜坡段部分，等 值线是向下倾斜，说明了表层的体积含水率大于土体内 部的体积含水率; 随着水平距离的变大，等值线变为趋于 平行的直线，说明水平入渗的影响在逐渐变小． ( 2) 边坡初始体积含水率对土体内部瞬态体积含 水率也有很大的影响． 选取降雨强度为 450 mm·d － 1， 初始体积含水率分别选 0. 15、0. 25、0. 30 及 0. 35，降雨 时间选取 480 min，距边坡边缘 1. 0 m 处作垂直剖面，研 · 114 ·
·412… 工程科学学报，第37卷，第4期 70 0.46 0.40 0.380.36 0.44 -40em ◆-80cm -0.34 50 0.42 120em 0.32 --160 cm 100 0.40P 0.30 0.38 150 0.36 0.34 200 0.32 250 50100150 200 0.30 250 300cm 0280 50 100 150 200 250 300 图8降雨时间为1440min时体积含水率等值线分布图 深度lcm Fig.8 Contour maps of volume water content after 960 min rainfall 图10距边坡边缘不同水平距离时体积含水率随深度的变化 究剖面上各点体积含水率.图9显示了不同初始体积 Fig.10 Change of volumetric water content at each point on the 含水率时边坡内部体积含水率随深度的变化规律。 vertical section of horizontal distance with depth 0.50 由图可知，随着距离边坡边缘水平距离越小，不同 0.15 0.45 0-0.25 垂直剖面上各点体积含水率随深度的分布差别较大， -0.10 当水平距离大于120cm后，不同位置处垂直剖面上各 0.40 -0.35 点土体体积含水率随深度的变化趋于一致，说明随着 0.35 水平距离的变大，水平入渗对土体体积含水率的分布 0.30 影响在减小，当大于一定值后水平入渗影响可以忽略。 0.25 因此，在降雨过程中，要特别注意边坡表面体积含水率 0.20 的变化情况. 随着距边坡表面水平距离的不同，土体内的体积 0.15 含水率也是不同的.降雨强度为300mm·d,降雨时 0.106 50 100 150200 250 300 间为960min时，不同深度的土体体积含水率随着水平 深度/cm 距离的变化规律如图11所示 图9不同初始体积含水率时体积含水率随深度的变化 Fig.9 Change of volume water content with depth 0.42 0.40 由图可知，土体的初始体积含水率对降雨过程中土 体体积含水率的分布有着很大的影响.初始体积含水率 0.38 -50cm 越小，表层的体积含水率的增长幅度就越大，但影响深度 装036 4100cm 6-150cm 越小，随着深度增加，体积含水率变化速率也越大 .-200cm 非饱和土的抗剪强度与土体的体积含水率息息 相关，较小的初始体积含水率虽然开始抗剪强度较 0.32 大，但受到降雨影响时，下降幅度较大，这对于工程 0.30 实际是非常不利的：较高的初始体积含水率虽然降 0.28 雨影响的深度较大，但是体积含水率的变化幅度较 50 100 150 200250 300 深度/em 小，土体的抗剪强度不会发生突变，对工程比较有 图11距边坡表面不同水平距离时体积含水率随水平距离的变化 利，便于采取措施. Fig.11 Change of volume water content with horizontal distance (3)由于考虑了降雨的水平入渗，距离边坡面不 同水平距离的土体体积含水率是不同的.分别选取距 由图可知：随着水平距离的增加，不同深度的水平 离边坡边缘水平距离为40、80、120及160cm的四个位 剖面上各点的体积含水率变化具有相似的规律，即逐 置，在四个位置处作垂直剖面，研究四个剖面上各点的 渐减小，而当水平距离大于1m后各点体积含水率趋 体积含水率随深度的变化规律.降雨强度取300mm· 于恒定.不同深度的水平剖面上各点的体积含水率， d,降雨时间为960min,得到四个位置体积含水率随 随着深度的增加逐渐减小，且各个剖面之间的体积含 深度的变化规律如图10所示. 水率差值越来越小
工程科学学报，第 37 卷，第 4 期 图 8 降雨时间为 1440 min 时体积含水率等值线分布图 Fig． 8 Contour maps of volume water content after 960 min rainfall 究剖面上各点体积含水率． 图 9 显示了不同初始体积 含水率时边坡内部体积含水率随深度的变化规律． 图 9 不同初始体积含水率时体积含水率随深度的变化 Fig． 9 Change of volume water content with depth 由图可知，土体的初始体积含水率对降雨过程中土 体体积含水率的分布有着很大的影响． 初始体积含水率 越小，表层的体积含水率的增长幅度就越大，但影响深度 越小，随着深度增加，体积含水率变化速率也越大． 非饱和土的抗剪强度与土体的体积含水率息息 相关，较小的初始体积含水率虽然开始抗剪强度较 大，但受到降雨影响时，下降幅度较大，这对于工程 实际是非常不利的; 较高的初始体积含水率虽然降 雨影响的深度较大，但是体积含水率的变化幅度较 小，土体的 抗 剪 强 度 不 会 发 生 突 变，对 工 程 比 较 有 利，便于采取措施． ( 3) 由于考虑了降雨的水平入渗，距离边坡面不 同水平距离的土体体积含水率是不同的． 分别选取距 离边坡边缘水平距离为 40、80、120 及 160 cm 的四个位 置，在四个位置处作垂直剖面，研究四个剖面上各点的 体积含水率随深度的变化规律． 降雨强度取 300 mm· d － 1，降雨时间为 960 min，得到四个位置体积含水率随 深度的变化规律如图 10 所示． 图 10 距边坡边缘不同水平距离时体积含水率随深度的变化 Fig． 10 Change of volumetric water content at each point on the vertical section of horizontal distance with depth 由图可知，随着距离边坡边缘水平距离越小，不同 垂直剖面上各点体积含水率随深度的分布差别较大， 当水平距离大于 120 cm 后，不同位置处垂直剖面上各 点土体体积含水率随深度的变化趋于一致，说明随着 水平距离的变大，水平入渗对土体体积含水率的分布 影响在减小，当大于一定值后水平入渗影响可以忽略． 因此，在降雨过程中，要特别注意边坡表面体积含水率 的变化情况． 随着距边坡表面水平距离的不同，土体内的体积 含水率也是不同的． 降雨强度为 300 mm·d － 1，降雨时 间为 960 min 时，不同深度的土体体积含水率随着水平 距离的变化规律如图 11 所示． 图11 距边坡表面不同水平距离时体积含水率随水平距离的变化 Fig． 11 Change of volume water content with horizontal distance 由图可知: 随着水平距离的增加，不同深度的水平 剖面上各点的体积含水率变化具有相似的规律，即逐 渐减小，而当水平距离大于 1 m 后各点体积含水率趋 于恒定． 不同深度的水平剖面上各点的体积含水率， 随着深度的增加逐渐减小，且各个剖面之间的体积含 水率差值越来越小． · 214 ·
付建新等：考虑二维降雨入渗的非饱和土边坡瞬态体积含水率分析 ·413 (林鸿州，于玉贞，李广信，等.降雨特性对土质边坡失稳的 4结论 影响.岩石力学与工程学报，2009,28(1)：198) (1)一维降雨入渗分析与实际过程是有差距的. Tong L Y,Jian W X.Stability analysis of unsaturated slope con- sidering rainfall infiltration process.J Wuhan Unir Technol, 本文建立了二维降雨入渗模型，并采用交替隐式算法 2012,34(11):78 进行了计算，结果与实际情况相符 (童龙云，简文星.考虑降雨入渗过程的非饱和边坡稳定性分 (2)在相同的降雨强度下，降雨时间越长，降雨入 析.武汉理工大学学报，2012,34(11)：78) 渗影响的范围越大.由于考虑了水平入渗的影响，因 [8]Xu H,Zhu Y W,Cai Y Q,et al.Stability analysis of unsaturated 此距离边坡表面不同水平距离的土体体积含水率是不 soil slopes under rainfall infiltration.Rock Soil Mech,2005,26 同的，这也是传统一维降雨入渗分析无法得到的结果. (12):1957 (徐晗，朱以文，蔡元奇，等。降雨入渗条件下非饱和土边坡 (3)土体的初始体积含水率对降雨过程中土体体 稳定分析.岩土力学，2005,26(12)：1957) 积含水率分布有着很大的影响.土体的初始体积含水 9]Lei Z D,Yang S X.Numerical method of one-dimensional flow 率越小，降雨初期表层土体的体积含水率增长越快，但 through unsaturated soils.Acta Pedol Sin,1982,19(2):141 影响深度越小：随着深度的增加，土体的体积含水率降 (雷志栋，杨诗秀.非饱和土壤一维流动的数值计算.土壤学 低的速率也越大.因此，在降雨过程中，较干燥的边坡 报，1982,19(2)：141) 稳定性变化更加刷烈，需要特别注意 [10]Tan X,Chen S X,Yang M.Saturated-unsaturated seepage anal- (4)由二维降雨入渗的计算结果可知，距离边坡 ysis of slope under rainfall.Rock Soil Mech,2003,24(3):381 (谭新，陈善雄，杨明.降雨条件下土坡饱和一非饱和渗流分 表面越近，土体各处的体积含水率越大，距边坡表面一 析.岩土力学，2003,24(3)：381) 定水平距离范围内，同一垂直剖面上体积含水率随深 [11]Li N.Xu JC.Analytical solutions for rainfall infiltration into ho- 度的分布差别较大，但大于一定距离后随深度的变化 mogenous infinite slopes.Chin J Geotech Eng,2012,34(12): 趋于相同. 2325 (李宁，许健聪.无限长均质斜坡降雨入渗解析解.岩土工程 参考文献 学报，2012,34(12)：2325) [1]Chen Z Y.Soil Slope Stability Analysis.Beijing:People's Press, [12]Zhan L T,Jia G W,Chen Y M,et al.Analytical solution for 2003 rainfall infiltration into infinite long slopes considering properties (陈祖煜.土质边坡稳定性分析.北京：人民出版社，2003) of unsaturated soil.Chin J Geotech Eng,2010.32(8):1214 Dewiest R JM.Flow Porous Media.New York:Academic Press, (詹良通，贾官伟，陈云敏，等。考虑土体非饱和特性的无限 1969 长斜坡降雨入渗解析解.岩土工程学报，2010,32(8)： B] Neuman S P.Galerkin approach to saturated-unsaturated flow in 1214) porous media.Finite Elem Fluids,1974,1:201 [13]Lei Z D,Yang S X,Xie S C.Soil Water Dynamics.Beijing: 4]Wu B.Sun DA.Study of liquefaction characteristics of unsaturat- Tsinghua University Press,1988 ed silt.Rock Soil Mech,2013,32(2)411 (雷志栋，杨诗秀，谢森传.土壤水动力学.北京：清华大学 (吴波，孙德安.非饱和粉土的液化特性研究.岩土力学， 出版社，1988) 2013,32(2):411) [14]Fredlund D G,Rahardjo H.Soil Mechanics for Unsaturated Soil. 5]Sheng DC,Yang C.Discussion of fundamental principles in un- New York:John Wiley Sons Inc,1993 saturated soil mechanics.Chin J Geotech Eng,2012,34(3):438 [15]Hu K L,Li B G,Chen Y.Effect of spatial variability of surface (盛岱超，杨超.关于非饱和土本构研究的几个基本规律的探 soil saturated hydraulic conductivity on water drainage at field 讨.岩土工程学报，2012,34(3)：438) scale.J Hydraul Eng,2006,37(10):1217 6]Lin HZ,Yu YZ,Li C X,et al.Influence of rainfall characteris- (胡克林，李保国，陈妍.表层土壤饱和导水率的空间变异 tics on soil slope failure.J Rock Mech Eng,2009,28(1):198 对农田水分渗漏的影响.水利学报，2006,37(10)：1217)
付建新等: 考虑二维降雨入渗的非饱和土边坡瞬态体积含水率分析 4 结论 ( 1) 一维降雨入渗分析与实际过程是有差距的． 本文建立了二维降雨入渗模型，并采用交替隐式算法 进行了计算，结果与实际情况相符． ( 2) 在相同的降雨强度下，降雨时间越长，降雨入 渗影响的范围越大． 由于考虑了水平入渗的影响，因 此距离边坡表面不同水平距离的土体体积含水率是不 同的，这也是传统一维降雨入渗分析无法得到的结果． ( 3) 土体的初始体积含水率对降雨过程中土体体 积含水率分布有着很大的影响． 土体的初始体积含水 率越小，降雨初期表层土体的体积含水率增长越快，但 影响深度越小; 随着深度的增加，土体的体积含水率降 低的速率也越大． 因此，在降雨过程中，较干燥的边坡 稳定性变化更加剧烈，需要特别注意． ( 4) 由二维降雨入渗的计算结果可知，距离边坡 表面越近，土体各处的体积含水率越大，距边坡表面一 定水平距离范围内，同一垂直剖面上体积含水率随深 度的分布差别较大，但大于一定距离后随深度的变化 趋于相同． 参 考 文 献 ［1］ Chen Z Y． Soil Slope Stability Analysis． Beijing: People's Press， 2003 ( 陈祖煜． 土质边坡稳定性分析． 北京: 人民出版社，2003) ［2］ Dewiest Ｒ J M． Flow Porous Media． New York: Academic Press， 1969 ［3］ Neuman S P． Galerkin approach to saturated-unsaturated flow in porous media． Finite Elem Fluids，1974，1: 201 ［4］ Wu B，Sun D A． Study of liquefaction characteristics of unsaturated silt． Ｒock Soil Mech，2013，32( 2) : 411 ( 吴波，孙德安． 非饱和粉土的液化特性研究． 岩 土 力 学， 2013，32( 2) : 411) ［5］ Sheng D C，Yang C． Discussion of fundamental principles in unsaturated soil mechanics． Chin J Geotech Eng，2012，34( 3) : 438 ( 盛岱超，杨超． 关于非饱和土本构研究的几个基本规律的探 讨． 岩土工程学报，2012，34( 3) : 438) ［6］ Lin H Z，Yu Y Z，Li G X，et al． Influence of rainfall characteristics on soil slope failure． J Ｒock Mech Eng，2009，28( 1) : 198 ( 林鸿州，于玉贞，李广信，等． 降雨特性对土质边坡失稳的 影响． 岩石力学与工程学报，2009，28( 1) : 198) ［7］ Tong L Y，Jian W X． Stability analysis of unsaturated slope considering rainfall infiltration process． J Wuhan Univ Technol， 2012，34( 11) : 78 ( 童龙云，简文星． 考虑降雨入渗过程的非饱和边坡稳定性分 析． 武汉理工大学学报，2012，34( 11) : 78) ［8］ Xu H，Zhu Y W，Cai Y Q，et al． Stability analysis of unsaturated soil slopes under rainfall infiltration． Ｒock Soil Mech，2005，26 ( 12) : 1957 ( 徐晗，朱以文，蔡元奇，等． 降雨入渗条件下非饱和土边坡 稳定分析． 岩土力学，2005，26( 12) : 1957) ［9］ Lei Z D，Yang S X． Numerical method of one-dimensional flow through unsaturated soils． Acta Pedol Sin，1982，19( 2) : 141 ( 雷志栋，杨诗秀． 非饱和土壤一维流动的数值计算． 土壤学 报，1982，19( 2) : 141) ［10］ Tan X，Chen S X，Yang M． Saturated-unsaturated seepage analysis of slope under rainfall． Ｒock Soil Mech，2003，24( 3) : 381 ( 谭新，陈善雄，杨明． 降雨条件下土坡饱和--非饱和渗流分 析． 岩土力学，2003，24( 3) : 381) ［11］ Li N，Xu J C． Analytical solutions for rainfall infiltration into homogenous infinite slopes． Chin J Geotech Eng，2012，34( 12) : 2325 ( 李宁，许健聪． 无限长均质斜坡降雨入渗解析解． 岩土工程 学报，2012，34( 12) : 2325) ［12］ Zhan L T，Jia G W，Chen Y M，et al． Analytical solution for rainfall infiltration into infinite long slopes considering properties of unsaturated soil． Chin J Geotech Eng，2010，32( 8) : 1214 ( 詹良通，贾官伟，陈云敏，等． 考虑土体非饱和特性的无限 长斜坡降 雨 入 渗 解 析 解． 岩 土 工 程 学 报，2010，32 ( 8 ) : 1214) ［13］ Lei Z D，Yang S X，Xie S C． Soil Water Dynamics． Beijing: Tsinghua University Press，1988 ( 雷志栋，杨诗秀，谢森传． 土壤水动力学． 北京: 清华大学 出版社，1988) ［14］ Fredlund D G，Ｒahardjo H． Soil Mechanics for Unsaturated Soil． New York: John Wiley ＆ Sons Inc，1993 ［15］ Hu K L，Li B G，Chen Y． Effect of spatial variability of surface soil saturated hydraulic conductivity on water drainage at field scale． J Hydraul Eng，2006，37( 10) : 1217 ( 胡克林，李保国，陈妍． 表层土壤饱和导水率的空间变异 对农田水分渗漏的影响． 水利学报，2006，37( 10) : 1217) · 314 ·