非饱和黄土边坡降雨数值模拟研究

张绍献

（广东工业大学 土木与交通工程学院，广东 广州510006）

1 概述

目前国内外均有学者在不断进行降雨入渗方面的研究，而研究的重点主要放在对不同深度非饱和黄土的体积含水率的变化规律的研究上。Liang-tong ZHAN 等 对西安地区黄土的储水能力进行了现场试验研究，结果表明掺碎石黄土地表能提高土体41%的储水能力[1]。L. Z. Wu & A. P. S. Selvadurai 通过有限元软件COMSOL Multiphysics 对降雨入渗下二维非饱和多孔介质进行了渗流- 变形耦合分析，模拟出了水分入渗与土体变形随着时间的变化过程[2]。Chang-Liang Zhang 等 在我国甘肃省正宁县对该地区黄土进行降雨入渗研究，现场监测结果表明，地表2m 范围内的含水率容易被降雨影响[3]。许方领等针对降雨入渗作用下非饱和土边坡的变形和稳定问题, 运用ABAQUS 有限元软件，对边坡渗流场, 应力场和位移场进行分析；结果表明：降雨入渗引起边坡饱和度,等效应力和位移增大，边坡表层响应早于坡体内部[4]。

2 黄土非饱和入渗过程分析

单元土体在水的入渗过程中会经历3 中状态：（1）天然状态；（2）非饱和状态；（3）饱和入渗状态。单元土体中水分的状态可以用饱和度S 来描述：

在一维降雨入渗深度公式中，通过Boltzmann 变换，可得到对应入渗深度的变换公式为：

取前2 项作为近似，则有

对z 求导，可得

由达西定律，有

比较式（4）和（5）可得

式中：λ—与土体特性和降雨入渗边界条件相关的参数；a2—与渗透系数k 相关的特性参数，a1—与土体的扩散性有关的参数。

3 边坡入渗数值模拟

本文数值模拟采用有限元软件abaqus 对土体内部水分的移动的进行模拟。

模型的尺寸与网格划分如图1 所示，降雨荷载施加在整个坡面上，边界孔压条件设置为随着高度线性增加的静水孔压边界，地下水位设置在边坡的底面即BC 处；边坡土体设置为均质同一黄土，该黄土的基本力学参数见表1。

图1 黄土边坡尺寸

表1 土体参数

4 模拟结果

本文首先分析了降雨后边坡的孔压场变化，然后对降雨前后的边坡进行强度折减，并比较安全系数的变化。模拟结果如下：

图2 降雨24h 后边坡孔压场分布

图3 降雨后折减边坡的等效塑性区分布

图4 降雨前折减边坡的等效塑性区分布

图5 降雨后边坡安全系数随水平位移U 的关系曲线

图6 降雨前边坡安全系数随水平位移U 的关系曲线

5 结论

本文利用abaqus 有限元软件先计算天然黄土边坡的渗流场，然后基于摩尔库伦准则，利用强度折减法对比降雨前后的边坡稳定性，通过对比分析得到以下结论：

5.1 降雨前初始孔压呈线性分布，并且随着深度的增加而不断增大，孔压的最小值在边坡的下边界，为0kPa，最大值位于坡顶，数值为-250kPa。在降雨15.5h 后，整个坡面的孔压场迅速减少，平均孔压降低为61.25kPa，降雨24h 后，坡面的平均孔压降低为39.56kPa，降雨影响深度随着降雨时长的增加而增加。

5.2 在降雨前边坡的最大塑性应变为11.51，安全系数取降雨前安全系数随位移变化曲线的拐点，即降雨前边坡的安全系数为3.28。在降雨后边坡的最大塑形应变为11.74，安全系数取降雨后安全系数随位移变化曲线的拐点，即降雨24h 后边坡安全系数为2.46。对比可知降雨后边坡的塑性应变增大，安全系数明显降低，这是因为降雨后边坡的基质吸力降低，导致边坡的稳定性明显降低，更容易发生滑坡事故。