大兴庄水库降水作用下大坝稳定性分析

王起芝

(山东省临沂市费县水利局，山东 临沂 273400)

我国水库数量位居世界榜首，多数修建于20世纪50—70年代，但多数水库由于施工处理不当、填筑标准低等问题，长期处于“带病”运行的状态，时常出现坝体裂缝和坝坡滑移、溃坝等事件[1- 2]。库水位变化是当前引发大坝破坏的重要因素之一，频繁的库水位变化会降低坝体内部的抗剪强度，及时开展水库水位变化后大坝渗流状态的研究十分重要[3]。岑威钧等[4]利用刚体极限平衡理论，通过考虑土工膜后渗流对坝坡稳定的作用，对其进行了数值模拟和分析；唐书童等[5]根据水文监测数据，对城西水库堤防渗漏的安全性和今后高水位作业情况下的防渗措施进行了研究；贺芳丁等[6]对水库的蓄水安全性评价与水库大坝的安全性评价进行了比较；尹芳芳等[7]从多方面对磨刀坑水库大坝安全性进行了全面复核和评价，针对大坝的结构渗流稳定等进行了分析和论证；罗江锋[8]参考一个除险加固的工程实例，考虑渗流场和应力场间的耦合，对土石坝加固前、后的坝体受力作用进行计算分析；孙玮玮等[9]针对黄河源电站的实际情况，采用数值模拟方法，对黄河源发电站2种不同的渗流场进行了数值模拟和计算，对其渗透安全性进行了综合评估；简鸿福等[10]根据试验资料，对堤体各个岩石和土体的力学性能进行了反演，以确保堤防渗漏的安全性；岳庆河等[11]建立老岚水库土石坝渗流典型断面有限元模型，对不同上下游水位工况及假设可能出现的不利工况进行数值模拟，。

文章以大兴庄水库为例，采用GEO-STUDIO软件并结合非饱和土性质，主要探讨了大坝的坝基土的导水率、水库水位下降率、水库水位下降等因素对大坝的作用，研究成果可为相关工程提供参考。

1 工程概况

大兴庄水库坝址在薛庄镇大兴庄村北200m，位于沂河水系的薛庄河，该水库的控制流域为2.8km2，库容22.6万m3，兴利库容7.3万m3，死库容0.08万m3。大兴庄水库坝体是一种均质坝，坝体宽4.5m，高10m。上下游坡比均为1∶2.0。大兴庄水库地处季风型、半湿润地区，是一种具有显著季节性的温带大陆性气候。春天多风，降水稀少，夏天多雨，秋天有丰富的阳光，冬天又冷又干燥。本区年平均降水量为852.5mm，降水年际差异大，造成了干旱、洪涝年的不均衡；年平均日照时间为2529.4h，平均日晒百分比为59%；最大冻结深度为28cm；年平均蒸发面积达1605mm。该区年均最高风速为3.4m/s，最大风速为17.0m/s。风向最多的是东风。

2 数值模型及计算参数

本文建立的大坝数值模型结构如图1所示。

图1 大坝数值模型(单位：m)

大坝横断面的尺寸、材料区域和土壤性质根据现场试验取得，在土体非饱和状态下，土体含水量特性曲线是土体结构研究中的一个重要环节。因为实测的SWCC不能获得大坝的具体用材，因此本文的心墙计算参数根据勘测数据进行估算取得，大坝和沙滤层是由其颗粒大小来估计的。本次采用GEO-STUDIO软件进行分析，采用瞬态渗透法测定了水库在下游的孔隙水压，并对上游坝体进行了稳定计算。2种方法均在自动耦合模式下进行。计算参数见表1。

表1 数值计算参数

2.1 渗流方程

通过SEEP/W软件进行了渗透数值计算，分析了大坝在水力作用下的二维稳定渗流场和水压作用下的瞬时渗透特性。应用有限单元法，对渗透问题进行了数值模拟，得到了由公式(1)表示的广义的控制微分方程：

(1)

式中，H—总水头，m；kx—X方向的导水率，cm/s；ky—Y方向的导水率，cm/s；Q—边界通量，m3/s；θ—体积含水量，无量纲；t—时间，s。

坝体有限单元网格的剖分情况如图2所示。

图2 坝体的有限元网格剖分状况

网格剖分分别由4和3个结点组成的四边形与三角形单元。假定坝顶的最高水位(19.2m)是ABC与DE的界限的一个关键条件，即一个常量的总水头，它的数值为19.2和9.2m。在此基础上，通过稳定的渗透计算，得到了初步的孔隙水压力分布。上游的边界线是以总水位的形式来确定的，这样的界限使得使用者可以根据时间来确定总水位的改变。在瞬时渗透计算中，随着水库水位的降低，上游边坡面的水压由19.2m直线降低到9.2m。在水头降低过程中，坝体由饱和向非饱和转变，因而，各物料的饱和导水系数被称为吸力的函数。通过对坝体在降水过程中的瞬态渗透进行数值模拟，得到了全网格内孔隙水压随着水位的降低而发生的动态变化。

2.2 坝体稳定分析方程

为了准确分析出水位下降时，坝体上游边坡的稳定性。本文利用SLOPE/W中已有的GLE法对边坡失稳进行了数值仿真。通过计算各个斜坡的潜在失效面，得到了在不同的水位变化条件下，最大安全因子的最优解。在这种情况下，由于大坝的水位降低，出现了黏土心墙坝、沙滤层等坝体不饱和的情况，故本文将材料的不饱和问题和剪切强度也考虑在内，并给出了相应的计算公式(2)。利用抗剪强度参数和SWCC来估算材料的非饱和抗剪强度：

(2)

式中，τ—饱和或非饱和土的抗剪强度，kPa；φ′—有效摩擦角，(°)；(ua-uw)—土壤破坏面内的基质吸力，kPa；Ua—非饱和土壤中的孔隙空气压力，kPa；Uw—孔隙水压力，kPa。

公式(2)可以用于瞬态渗透计算得到孔隙水压力变化，可应用于上游坝坡在各种水位下的稳定性分析。

3 数值结果分析

3.1 导水率对边坡稳定性的影响

为分析土体饱和导水率(ksat)对上游边坡稳定性造成的影响，本文对不同水位下降情况下出现的4种不同ksat边坡安全系数变化情况进行分析，如图3所示。

图3 不同ksat下土坝上游边坡的安全系数随水位变化

由图3可知，对于给定ksat下，随着大坝的水位降低，大坝的上游边坡面稳定性逐渐降低，在库容降低到1/2～2/3时为最低。这主要是因为大坝的蓄水量减少和大坝暂时不排水造成的，而且随着孔隙水压力的增大，土体剪切强度也会随之减弱，从而降低了安全系数。在蓄水过程中，在蓄水量由原来的1/2降低到2/3的情况下，稳定性由最低值逐渐增大，直到水位不再下降。坝体材料的抗剪强度因其孔隙水的消失而增大。这对斜坡体系的稳定起到了很好的作用。当水库的水位降低时，从图3中可以看出，土体的饱和导水率(ksat)对上游边坡稳定性有显著影响。在水位下降的早期阶段，堤坝的上游边坡采用低饱和导水率建造边坡比高导水率更稳定。然而，随着水位的进一步降低，边坡稳定性会随着ksat值升高而改变。值得注意的是，饱和带(潜水表面以下的区域)在稳定状态下，高ksat坝体的渗流条件大于低ksat坝体的渗流条件。由于高ksat土壤允许水更深地渗透到大坝中，饱和区会逐渐扩大，这是坝体在水位下降早期阶段，安全系数降低的主要原因，低ksat土壤中较小饱和区域的孔隙水压力可能比高ksat土壤中较大饱和区域的孔隙水压力消散得更快，因此，超孔隙水压力较小的边坡比超孔隙水压力较大或饱和区域较大的边坡更稳定。

上游坝坡稳定性在饱和与非饱和状态下其抗剪强度的变化情况如图4所示。

图4 饱和与非饱和抗剪强度下边坡安全系数变化趋势

由图4可知，在水位下降期间，通过计算抗剪强度，得出的结果比非饱和抗剪强度的结果稍好。因为非饱和物料可以使水库的水流更深入地渗透到堤坝中，从而使堤坝更加湿润，这是导致非饱和物料堤坝安全性下降的重要因素。相对于不饱和物料坝的高湿度区域，在饱和物料坝体的小湿度区域，其孔隙水压力的消散速度要快，所以在低孔隙水压或在低湿区的水坝要比在高水压条件下的水坝更稳定。但是，应该指出，大坝的饱和状态和现实中的不一样，比如黏土心墙、路堤和沙滤等，在大坝的水位高度到达大坝顶部时，也不会完全饱和。

3.2 水位下降率对边坡稳定性的影响

本节对水位下降率对上游边坡稳定性造成的影响进行分析，分别为8、6、2和0.5m/d，如图5所示。

图5 水位下降率对上游边坡稳定性的影响

由图5可知，当水库的水位下降率减小时，水库的上坡面稳定性明显提高。由公式(2)可知，由于坝顶的高度下降速率减小，导致大坝中的孔隙水压力通过坝内饱和区渗透的时间延长，因而提高了大坝的抗剪能力。当孔隙水压力消散时，可以使物料内部的孔隙水压力下降到负数(吸力)，因此边坡稳定性增加。相比之下，更高的水位下降率或快速下降(例如，8m/d)，坝体没有足够的时间进行孔隙水压力渗透，因此在较高的收缩率下，材料的抗剪强度增加小于在较低的收缩率下材料的抗剪强度增加。

4 结论

本文以大兴庄水库为实例，采用GEO-STUDIO程序分析了大坝的饱和导水率、水库水位下降率、水库水位下降等因素对大坝上游边坡稳定性作用。研究表明当库容降至总水位2/3左右时，大坝的上游坡体稳定性减小，而当水库蓄水量继续下降时，大坝上游坡体的安全性逐步增大。此外，当水库的蓄水量降低时，高渗透率土壤的土坝稳定性要大于低渗透性土坝同时当水库的水位下降率比较低时，坝体的上游边坡稳定性相对较高。因此，建议在进行水库大坝稳定性设计时，应当考虑库水位下降高度以及水位下降速率对大坝稳定性的影响，从而提高坝体稳定性。