张文兵，沈振中，陈官运，张琬琳，魏舒萌，汪千敢

(1.河海大学 水文水资源与水利工程科学国家重点实验室，江苏 南京 210098； 2.河海大学 水利水电学院，江苏 南京 210098； 3.中建中东有限责任公司，迪拜 63932； 4.巨野县水务局，山东 菏泽 274900)

1 研究背景

根据我国“十三五”期间能源结构调整的总体趋势和要求，水资源作为一种高效且清洁的能源，将继续在我国的能源发展战略中发挥不可或缺的作用。随着对水资源的不断开发与利用，在复杂地质条件下修库筑坝成为高效开发利用水电能源不可避免的挑战。常见的库址区复杂地质条件主要包括深厚覆盖层[1-2]、地质断层[3]和岩溶[4]等类型。岩溶地貌在我国分布较为广泛，特别是在西南地区，一些水库大坝工程在建设过程中需要穿越岩溶发育区[5]。溶洞是岩溶地层中最为常见的产物，其存在将影响坝基及库岸边坡的渗透稳定性。因此，对相关工程进行渗流分析和渗控效果评价具有重要的工程参考价值。

目前，对于复杂地质条件下的库区渗流研究多集中在地质断层和深厚覆盖层上的坝体、坝基及两岸坝肩的渗流特性分析[6-10]。Jiang等[11]根据四川广安龙滩水库所处的地质环境，构建了岩溶发育地质条件下的水库地下水渗流模型，并对库底岩溶渗漏通道的涌水量进行了模拟；Mohammadi等[12]以伊朗Khersan大坝为案例进行研究，提出了一种研究岩溶区坝址渗漏的方法；柳昭星等[13]针对某矿区岩溶发育地层突水灾害帷幕灌浆防渗工程，开展了岩溶发育地层帷幕灌浆材料的性能及适用性研究；刘胜[14]基于同位素-水文地球化学方法对岩溶区水电工程灌浆帷幕的可靠性进行了综合评价研究。尽管有关岩溶发育地质条件下的工程渗流研究成果颇为丰硕，但针对岩溶发育地质条件下的库区渗流特性研究鲜见报道。

本文结合位于华阳河干流上某岩溶地质发育的水库大坝工程，基于等效连续介质模型，采用饱和-非饱和渗流计算理论，依托此枢纽工程的水文地质条件和工程设计资料，对库区主要建筑物、复杂地质以及防渗结构进行精细化模拟，建立水库运行期的三维渗流分析有限元整体模型，开展三维渗流计算。并在此基础上，结合防渗设计方案，对坝体、坝基和库岸的渗透稳定性等关键性问题进行分析，综合评价所采取的防渗型式的渗控效果。研究成果可为类似地质条件下的水库大坝工程防渗系统设计及优化布置提供参考。

2 饱和-非饱和渗流计算理论

2.1 基本微分方程

土坝及其库岸边坡中的水体运动属于饱和-非饱和渗流问题，其研究难点不仅在于非饱和渗流是一个多因素(固-液-气三相的体积比、空气压力、土体骨架变形和温度等)相互耦合作用的过程，而且还在于寻求非饱和土体的渗透系数[15]。与利用达西定律求解饱和土体渗流所不同的是，非饱和土体的渗透系数不再是一个常量，而是随着含水量变化的函数[16]。

目前，对于饱和-非饱和土体中的渗流问题主要采用多孔介质饱和-非饱和渗流理论来进行计算，其基本方程如下[17-18]：

(1)

土壤水力函数选用van Genuchten模型进行描述[19]：

(2)

θ(hc)=θr+(θs-θr)[1+(αhc)n]-m

(3)

=α(θr-θs)(n-1)(αhc)n-1[1+

(αhc)n]1/n-2

(4)

式中：θr为残余含水率，m3/m3；θs为饱和含水率，m3/m3；α为水分特征曲线进气值的倒数，m-1；n为水分特征曲线坡度的指示参数，通过拟合土壤水分特征曲线得到；m=1-1/n。

2.2 定解条件

饱和-非饱和渗流微分方程的定解条件包括初始条件和边界条件：

(1)初始条件

hc(xi,0)=hc(xi,t0) (i=1,2,3)

(5)

(2)边界条件

hc(xi,t)|Γ1=hc1(xi,t)

(6)

(7)

(8)

(9)

式中：ni为边界面外法线方向余弦；t0为初始时刻，s；hc1为已知水头，m；qn为已知流量，m3/s；qr(t)为降雨入渗流量，m3/s；hc(t0)为t0时刻渗流场水头，m；Γ1为水头边界；Γ2为流量边界；Γ3为饱和出逸边界；Γ4为降雨入渗边界。

2.3 模型离散化

在空间域内采用Galerkin有限元法，在时间域内采用隐式有限差分格式求解公式(1)。利用分部积分并考虑边界条件，则公式(1)可以表示为：

(10)

(11)

(12)

(13)

式中：NE为单元总数；Nn和Nm分别为节点n和m的形函数；t为时间变量，s；Δt为时间增量，s； (hc)t和(hc)t+Δt分别为在t和t+Δt时刻的节点压力水头，m；Ω为计算域；s为计算域边界。

3 计算模型

3.1 工程概况

某水库位于安徽省长江南岸一级支流水阳江的最大一级支流华阳河中上游，是华阳河干流上的骨干控制性工程。该工程是具有防洪、供水、灌溉和发电等综合效益的中型水库工程。水库规划为黏土心墙坝，大坝坝顶高程150.0 m、最大坝高29.0 m，正常蓄水位145.6 m、相应的库容约5 600×104m3，设计洪水位146.4 m、校核洪水位148.07 m。主要永久性水工建筑物大坝、溢洪道、输水隧洞等为3级建筑物，次要建筑物电站厂房等为4级，临时建筑物为5级。水库枢纽工程平面布置见图1。

图1 水库枢纽工程平面布置及断面规划图

库区的复杂地质条件是筑坝时可能面临的重要工程问题。对该水库工程库区地质的勘察结果表明，该库址区存在如下几点工程地质问题：(1)左坝肩岩溶及其渗漏问题；(2)河床砂卵砾石层的渗漏问题；(3)右岸带状山梁的渗漏与库岸再造问题。为此，在库首段左岸山体采用防渗帷幕对岩溶区进行防渗处理，防渗帷幕沿主坝轴线向左侧山体延伸。由于主坝坝基为厚5～13 m的砂卵石层，为解决其防渗问题，沿坝轴线设80 cm厚混凝土防渗墙，防渗墙顶部插入黏土心墙，墙顶高程127.0 m，墙底嵌入风化基岩0.5～1 m，最大墙深15 m。左侧坝下设防渗帷幕，深入相对隔水层以下5 m，帷幕线路沿坝轴线向山体方向延伸。副坝右侧山梁山体较单薄，其砾质土覆盖层下的薄层C2h组灰岩具中等透水性，采用帷幕灌浆防渗处理，灌浆轴线长度约为1 800 m，帷幕防渗标准按q≤5 Lu控制。当帷幕遇岩溶空洞区时，先用混凝土将岩溶空洞填塞，再进行帷幕灌浆，因而在后续的模型计算中，岩溶空洞区将按照混凝土材料进行设定。表1给出了库址区岩溶空洞位置信息。

表1 库址区岩溶空洞位置信息统计表

3.2 有限元模型

根据水库大坝、副坝及左右岸采取的防渗型式，结合工程的实际情况，构建水库运行期的三维渗流分析有限元模型，以期分析正常蓄水位工况下的坝体、坝基及左右岸坝肩山体的渗流分布规律，判断防渗设计方案的有效性以及坝体、坝基和防渗帷幕的渗透稳定性，综合评价所采取的防渗型式的防渗效果。有限元模型计算坐标系定义如下：原点o取自大地坐标(x,y)=(3 398 121.687, 379 443.924)；x方向为垂直河流方向，平行于坝轴线，以指向左岸为正；y方向为垂直坝轴线方向，以指向上游库区为正；z轴垂直向上，以高程为坐标。另外，为了简化模型计算，所构建的有限元模型上游边界为防渗线上游200 m，左岸边界为左岸主坝灌浆帷幕终点以左1.5倍坝高，右岸边界为右岸主坝灌浆帷幕端头以右200 m，下游边界为下游坝脚以外150 m，顶高程取实际地形，底高程截至0。有限元模型截取范围如图1所示。

综合分析库址区计算区域内地形地貌、地层岩性及岩溶等特征，建立库址区三维有限元渗流计算模型超单元网格。根据建筑物的布置、地质分层及计算要求等信息，选取13个控制断面，其中1、4、6、10～12断面存在空洞，2、3、5、7～9、13断面无空洞(见图1)，加密细分后生成的三维有限元网格结点总数为21 437个，单元总数为22 136个。图2为库址区三维有限元模型及网格剖分示意图。

图2 三维有限元模型及网格剖分

3.3 计算参数及边界条件

根据坝体各区材料特性、主要结构面的概化情况和各地层的钻孔注水试验成果，将库区材料按渗透性的不同划分为17个参数区，见表2。

本文中的模型边界主要可分为水头边界、出渗边界以及不透水边界3种。水头边界为库址区水位线以下的河床及地下水位截取边界，其中左、右岸截取边界天然地下水位最大高程为128.3和140.0 m，并由山体高处向河道逐渐降低，最终与河道水位衔接，根据已知资料，取下游河道水位高程为124.1 m；出渗边界为水位线以上的边界，即与大气相接触的所有边界；不透水边界包括模型上、下游截取边界及模型底边界。

表2 库区各材料渗透系数

4 结果分析与评价

4.1 库址区平面渗流特性分析

因库址区左、右岸及副坝处均存在岩溶空洞，其存在将会改变库址区地下水等值线分布，即在岩溶空洞位置处的地下水等值线会出现突降现象。因此，岩溶空洞可能会成为地下水富集区或库水渗漏通道，进而对整个水库工程产生不利的影响。为此，在水库建设过程中采用混凝土砂浆对岩溶空洞先进行填塞，再进行灌浆处理。图3给出了经防渗处理后的库址区地下水位等值线图。

图3 库址区地下水位平面等值线图(单位：m)

从图3中可以看出，对库址区的岩溶发育地质进行防渗处理后，库址区渗流场的分布规律明确，库水由水库通过坝体、副坝和两岸山体渗向下游。因左右岸地形、地质条件不对称，因此，左右岸岩体位势分布并不完全一致。从防渗线路处地下水等值线的稀疏程度来看，相较于右岸山体采用的混凝土防渗墙和防渗帷幕组合防渗，左岸的单独灌浆帷幕防渗效果稍差，但在左岸灌浆帷幕端出现较为明显的绕渗现象，表明左岸帷幕起到了有效阻渗作用。此外，沿着防渗线路区域，地下水等值线分布较为密集，而在防渗线路的上游和下游地下水等值线分布较为稀疏；并且，各岩溶空洞断面处的地下水等值线分布均较为密集，未出现稀疏分布的情况，表明该工程针对岩溶发育地质条件所采取的防渗体系能够起到有效的阻渗作用。

4.2 库址区剖面渗流特性分析

为研究库址区局部区域防渗效果，图4分别给出右岸、副坝、主坝和左岸段部分断面地下水等势线图，并在表3中给出了各断面处浸润线的最高位置及削减水头百分率。由图4(a)～4(d)可以看出，右岸山体中的浸润线在由混凝土防渗墙和灌浆帷幕形成的组合防渗处上下游形成了突降，削减水头约为17.20 m，占总水头的80.00%。断面4为右岸山体存在局部溶洞的断面，图4(c)中显示的地下水等势线分布规律明确，地下水等势线除在组合防渗处发生突降外，在断面的其他部位未出现明显的突降变化，表明对右岸山体岩溶空洞的防渗处理措施得当，且具有较好的防渗效果。图4(e)、4(f)为副坝处横剖面地下水等势线图，等势线在黏土心墙、混凝土墙和灌浆帷幕组合防渗处上下游形成了突降，削减水头12.90 m，占总水头的60.00%，阻渗作用明显。与断面4类似，断面6处的地下水等势线除在组合防渗处存在突降外，其余部位均未发生明显异常，因此对于副坝部位的岩溶空洞防渗处理效果显著。图4(g)、4(h)显示了主坝位置处的横剖面地下水等势线图，从图中可以看出，浸润线在由黏土心墙、混凝土防渗墙和灌浆帷幕形成的组合防渗处的上、下游形成了突降，并且在组合防渗体系下游坝体内的浸润线较为平缓，防渗体的联合阻渗作用共削减水头16.12 m，占总水头的74.98%，防渗效果明显。图4(i)、4(j)为截取的左岸山体处地下水等势线图，相较于混凝土防渗墙、黏土心墙和灌浆帷幕组合防渗体系，由灌浆帷幕单独形成的防渗处浸润线突降变化不明显，但削减水头13.15 m，占总水头的61.16%，并且在溶洞断面11处的地下水等值线分布规律正常，表明对左岸溶洞群的防渗处理方式合理有效，能够起到有效控渗作用。

图4 库址区各断面地下水等势线图(单位：m)

表3 库址区各断面浸润面的最高位置及削减水头百分率

4.3 渗透坡降及渗漏量分析

表4给出了库址区各部位的最大平均渗透坡降和渗透流量。由表4中可以看出，在正常蓄水位工况下，主坝防渗部位的渗透坡降较大，其最大平均渗透坡降出现在主坝黏土心墙起始部位地下水表面附近，达到10.31，表明在组合防渗体上下游水头差较大时，防渗体系能够起到有效阻渗作用；相比较而言，左岸的防渗体系的渗透坡降较小，其最大平均渗透坡降为0.13；副坝防渗部位的最大平均渗透坡降为9.68，出现在副坝黏土心墙起始部位地下水表面附近；右岸山体防渗部位的最大平均渗透坡降为10.75，出现在右岸山体混凝土防渗墙起始部位地下水表面附近。从表4中计算所得到的渗透流量来看，库水主要通过坝基、副坝以及左岸的绕坝渗流至下游，库址区整体渗透流量为8.636 L/s，渗漏量较小。由此可见，该工程所采取的防渗设计效果明显，库区渗漏量能够得到有效控制。

表4 库址区各部位的最大平均渗透坡降和渗透流量

5 结 论

本文采用饱和-非饱和渗流计算理论，结合等效连续介质三维有限元模型对某岩溶发育地质条件的库址区三维渗流场进行分析，得出如下结论：

(1)在当前的防渗体系下，库址区各位置地下水位能够得到控制，对于溶洞的防渗处理措施可以起到有效阻渗作用。

(2)库址区各部位的渗透坡降总体较小，能够满足渗透稳定性要求，防渗设计方案可以满足库址区控渗要求，其设计在技术上是合理的。

(3)该工程总渗漏量为8.636 L/s，库水主要通过坝基、副坝及左岸绕渗至下游河道，在现有的防渗体系下，库址区渗漏量能够得到有效控制。