深基坑内渗流破坏的数值模拟研究

李伟沛

(中山市小榄镇水利所，广东 中山 528400)

近年来随着城市的发展，越来越多的深地下室、地下停车场或隧道被修建，这些工程的基础就是深基坑[1]。深基坑的设计通常受到墙周围水流的影响[2]。地下水位下降引起的渗流现象影响了墙的整体稳定性和底部开挖的稳定性，在开挖底部时可能发生体积隆起或管涌[3]。当管涌发生在挖掘层时，破坏尤为严重[4]。

前人研究提出了几种计算土体抗渗透破坏稳定性的方法[5]，但是应用这些方法标记的深基坑中也会发生破坏，因此需要更精确的分析来阐明破坏机理[6]。本文为研究渗流对基坑底部稳定性的影响，更好地阐述渗流破坏现象，采用FLAC2D(快速拉格朗日连续分析)软件进行了数值分析。通过一系列的数值试验，发现土体和界面的特性对破坏机制、破坏形态和破坏水头损失有显著影响。

1 数值模拟程序

本文对围堰坝底部抗砂土渗透破坏稳定性进行了数值模拟研究。在均质各向同性半无限土体中，贯入深度等于板桩的深度D，低于水头H，如图1中所示。

图1 围堰剖面图

在研究破坏时需要考虑的参数：开挖宽度、围堰宽度、挡土墙厚度、渗透性土壤的厚度、墙的柔性、支柱刚度、墙的平移和旋转等[7]。主要的破坏类型包括两种：(1)流土，在向上的渗透水流作用下，表面局部范围内的土体或颗粒同时发生悬浮、移动的现象[8]；(2)管涌，在渗透水流作用下，水中的细颗粒在粗颗粒形成的孔隙中移动，以致流失，最终导致土体内形成贯通的渗流通道[9]。下文通过表的形式，用临界水头损失H/D的取值表示破坏的判断方式。本文研究的主要目的是检验使用有限差分或有限元方法数值分析是否能够正确描述由于向上渗流进入围堰从而出现的各种破坏机制。将由支柱固定的墙看作无限刚性墙，这表示在破坏机制上假定固定壁效应的存在。在经典方法中，仅从液压观点考虑墙壁的破坏，本文的研究可以与其形成对比。

通过计算机程序FLAC2D进行分析，土壤中的变化由FLAC2D程序中完全弹性—塑性非关联莫尔—库仑模型进行模拟。对于γsat/γw=2(γsat为土壤的容重，γw为水的容重)；弹性体积模量K=30 MPa；剪切模量G=11.25 MPa。

2 建模过程

在墙壁粗糙的情况下，需要对土壤和墙壁之间的界面进行模拟分析。材料特性，特别是指定给接口的刚度，取决于界面的使用方式。在土—结构相互作用的情况下，认为界面与周围土体相比是刚性的，但它可能会滑动，并且会因荷载而打开。厚接头更适合模拟墙体与土体界面的摩擦行为。

图2中所示的界面模型被用于模拟由库仑定律定义的土壤/墙壁接触。接口任一侧的触点在性质上类似于离散元方法中使用的接口。

图2 界面模型

图2中S代表滑块；T代表抗拉强度；Kn代表法向刚度；Ks代表剪切刚度；LN代表网格点N作用的长度；LM代表网格点M作用的长度。切线方向上的滑块代表库仑剪切强度标准。法向刚度代表正常接触的极限强度。界面的摩擦角为δ，内聚力c=0 kPa，法向刚度Kn=1.0 GPa/m，剪切刚度Ks=1.0 GPa/m。选择这些Kn和Ks值来近似表达墙与格栅刚性连接时的结果。

在稳定塑性流动条件下的坍塌条件预测是数值模型难以准确模拟的一个问题。国内外学者使用粘塑性或弹塑性算法，已经清楚地解释了平面应变承载力系数N与网格几何的关系，并表明随着土壤膨胀角的减小，该系数随之减小。因此，为了便于后续的计算，对分析方案进行初步模拟，用以测试网格尺寸、单元尺寸、边界条件和静止土压力系数K0的影响。

结果表明，粗网格的破坏机制难以描述，细网格的破坏机制能清晰展现。图3展示了由位移场表示的破坏机制，以及粗网格(40×20，水平×垂直)和细网格(80×40)单元在稳态塑性流动下最大剪切应变率的相应分布。弹性地基参数和静止土压力系数Ko在实际范围内的变化对临界水压损失没有任何显著影响。

图3 由位移场表示的破坏机制

图4为细网格(80×40)的边界条件。在变形和流动梯度集中的壁面附近，网格尺寸很小。为了使边界影响最小化，网格距墙的长度和网格的深度分别位于壁厚的6倍和5倍。从边界条件的一般性原则出发，假定底边界是固定的，右边和左边横向边界固定在水平方向。

图4 网格和边界条件

板桩墙由结构梁单元模拟，梁单元通过连接在其两侧的界面与土网格相连。因此，这堵墙起着不透水构件的作用。为了识别被动土压力为零以及由管涌引起的渗透破坏的极限情况，采用了以下三个模拟步骤：

(1)假设材料是弹性的，计算其地应力。假定有效水平应力与静止有效垂直应力之比为0.5，采用函数建立初始孔隙压力和有效应力。

(2)在图5所示的水力边界条件下，假设地下水位位于地表，即上边界条件，该假定是在地下水配置中将其作为不耦合地下水流动和力学计算原理上设定的，由于图中左右边界和下边界的材料不透水，故左右边界和下边界条件水头为0。并根据上述条件计算孔隙水压力的分布范围。

(3)对前一步建立的孔隙压力分布进行力学响应研究。

在计算被动土压力时，需要注意模拟程序没有施加墙的水平位移。当被动压力消失时，不容易得到H/D的精确值。在不同水头损失H/D、不同控制参数-土壤摩擦、土壤膨胀和界面摩擦的情况下，本文只研究了底部抗土壤渗透破坏的稳定性。

图5 水力边界条件

3 结果与分析

对不同的土壤摩擦角进行了数值研究，φ=20°、25°、30°、35°和40°的结果如表1所示。作为分类，流土破坏的结果用一个或两个星号表示，管涌破坏的结果用三个星号表示。ψ/φ表示膨胀材料比，反映了土颗粒之间竖向应变大小与摩擦力的关系；δ/φ是土壁壁厚与内摩擦角的比值，代表了壁厚与土的摩阻性质的关系。

表1 不同调节参数δ/φ、ψ/φ、φ的临界水头损失H/D

表1中：φ代表土的内摩擦角；δ代表土壁的壁厚；ψ代表土的膨胀角。表1中的结果表明，除了土壤摩擦角φ≥40°、膨胀材料比ψ/φ>1/2和δ/φ≥2/3的情况下，抗渗透破坏的底部稳定性始终与流土现象有关。可以看出，管涌开始于H/D=2.94时。该数值非常接近理论临界水头损失值H/D=2.95，当H/D小于2.94时，发生的都是流土现象。在这种情况下，出口水力梯度达到临界水力梯度值。图6显示了φ=40°、ψ/φ=1、δ/φ=2/3和H/D=3情况下位移场的破坏机理和相应的最大值分布情况。

图6 φ=40°、ψ/φ=1、δ/φ=2/3、H/D=3.16的破坏机理和相应的最大值分布情况

图7为φ=35°、ψ/φ=0、δ/φ=2/3和H/D=2.92情况下的破坏机理。

图7 φ=35°、ψ/φ=0、δ/φ=2/3和H/D=3的破坏机理和相应的最大值分布情况

然而，对于膨胀材料ψ/φ≥1/2，可获得一个三角形土柱，图8为φ=35°、ψ/φ=1/2、δ/φ=2/3、H/D=2.97情况下的破坏机理。模拟结果表明，土的摩擦角对破坏机理形状有重要影响

图8 φ=35°、ψ/φ=1/2、δ/φ=2/3和H/D=3的破坏机理和相应的最大值分布情况

从表1中可以看出，对于δ/φ=0，ψ/φ=1且φ在20°～40°之间变化时，H/D的临界值在2.64～2.93范围内。根据表1给出的结论，在该情况下发生的破坏主要为流土破坏。当φ=40°，ψ/φ=1，δ/φ=0～1时，H/D的临界值在2.93～3.16范围内，参照表1在该临界范围内破坏为管涌破坏。

这种现象可以解释为：当邻近墙的土壤因破坏而发生管涌时，在棱柱体垂直面上产生的剪切力阻止了矩形棱柱体的上升。因此，出现了一个三角形破坏棱柱。这符合极限分析理论框架中的运动容许机制。

当φ和ψ较大时，埋置处墙的剪力和水平土膨胀延缓了三棱柱的破坏。因此，在三棱柱上升之前，体积上升成为临界条件。部分由土壤支撑的柔性墙也可能出现这种情况。

4 结 论

防止渗流破坏的临界水力损失取决于土壤摩擦角和界面土/壁摩擦。管涌只出现在φ≥35°、ψ/φ≥1/2、外壁δ/φ≥2/3的密实膨胀砂土中。在这种情况下，出口坡度等于临界水力坡度，而在其他情况下，只会发生流土现象。模拟结果表明，土的膨胀角对破坏机理形态有显著影响。对于膨胀性材料，会发生三棱柱体的隆起破坏，而对于其他情况则是矩形棱柱体。本文通过FLAC2D程序估算了渗流破坏的情况，结果表明程序表现良好，模拟结果精确。此类方法可推广应用到类似工程的破坏风险评估中。