刘丽君，范雪枫，吴震宇

(1.四川大学 水力学与山区河流开发保护国家重点实验室 水利水电学院，四川 成都 610065；2.昆明勘测设计研究院，云南 昆明 650216)

我国西南山区大多数水闸建在覆盖层地基上。这些水闸的地基在材料组成、渗流特性、防渗方式和渗流监测方面具有一定的相似性。山区河流的河床覆盖层普遍较为深厚，结构和级配变化大，局部存在架空现象，透水性强，渗透稳定性差。渗压监测数据可以真实反映闸基渗流性态，但缺点是测点数量有限，只能观测局部渗流情况。采用有限元法等数值模拟技术进行渗流分析，可以比较真实地模拟地基材料组成、各种材料的渗流特性、防渗设施和渗流边界条件等，得到整个计算模型区域内的渗流场，因此，与基于监测数据方法相比，数值模拟更加全面地分析地基渗流状态[1- 10]。渗透系数等材料渗流参数的取值对地基渗流数值模拟结果影响较大，如果参数选取不合理，则渗流计算结果与实际渗流情况会存在较大差异。根据渗流监测数据率定材料的渗流参数(即渗流参数反演分析)，然后采用率定后的参数进行渗流数值模拟，可使数值模拟结果更接近实际情况。

某水闸位于四川岷江干流，主要功能为发电，地基主要由漂卵石层和砂壤土层相间组成，闸基渗流安全是该工程重要问题之一。在闸基防渗墙后钻孔埋设了8个渗压计，用于监测闸基渗压变化情况。本文首先根据该水闸2010—2016年期间的渗压监测数据，分析闸基渗压变化规律及其原因。然后根据渗压监测数据，采用基于响应面法的参数反演分析方法，对覆盖层材料的渗透系数进行率定。最后采用率定后的渗透系数，对枯期和汛期的闸基渗流场进行有限元数值模拟，验证了渗压变化的原因，并对闸基整体渗流情况进行了分析，获得了监测数据未能反映的闸基渗透坡降情况，为闸基渗流安全评估提供了参考。

1 闸基渗压监测数据分析

1.1 工程概况

水闸于20世纪70年代建成。该水闸从左岸至右岸由取水口拦污栅闸、拦河闸、挡水坝等建筑物组成。拦河闸横河向长度为156m，顺河向长度为25m，最大闸高21.4m，设有6孔泄洪冲沙闸。水库总库容93万m3，正常蓄水位945.00m，设计洪水位945.60m，校核洪水位948.87m。闸基覆盖层最大深度45m，采用上游混凝土铺盖和悬挂式防渗墙防渗，防渗墙位于闸室上游20m，厚度0.7m，最大深度27m。闸室上游铺盖长75～95m，闸室下游护坦长56m。闸址处两岸岸坡较缓，临河坡高300～500m，河谷较开阔，河床宽80～100m。闸址处河床覆盖层主要由漂卵石层和砂壤土层相间组成，厚度约40～62m。1#闸孔顺河向剖面如图1所示。

图1 1#闸孔顺河向剖面图

为监测闸基渗压，以及上游铺盖和悬挂式防渗墙的防渗效果，在防渗墙后沿坝轴线方向钻孔埋设了两排渗压计，其中上游侧布置了5支渗压计(编号为UP1～UP5)，下游侧布置了3支渗压计(编号为UP6～UP8)，总共8支渗压计。闸基渗流监测布置如图2所示。

图2 闸基渗压监测平面布置图

1.2 闸基渗压监测数据分析

闸基渗压分析采用2010年1月—2016年12月期间的监测数据。为了便于和库水位进行比较，渗压测值已转化为渗压水头高程。各测点渗压水头高程随时间变化规律一致，受篇幅限制，只给出了UP2测点的渗压测值过程线，如图3所示。

图3 闸基渗压实测过程线

从图3可以看出，渗压测值呈明显周期性变化规律，即汛期(5—10月)水位高，枯期水位低。各渗压计的最大测值基本出现在7月，在938.34～945.15m之间(受2013年7月9号泥石流影响，下游水位壅高超945m高程，各渗压计测值也相应偏高)；最小测值基本出现在12月～次年4月，在935.07～939.34m之间。渗压测值与库水位的相关性弱，如2015—2016年渗压测值的年变幅为2.74～4.25m，而相应库水位的年变幅仅为1.09～1.64m，因此，库水位变化不是引起渗压测值周期性变化的主要原因。

在枯期，拦河闸关闭，通过左岸取水口引水发电，下游河道基本无水。在汛期，开闸泄洪，下游河道水位较高。地基采用上游混凝土铺盖和悬挂式防渗墙防渗，防渗墙位于闸室上游20m。从地基防渗布置来看，其主要作用是控制上游库水向下游渗透，而不能有效控制下游河水向上游渗透。所以，渗压测值的周期性变化估计是汛、枯期下游河道水位大幅变化而行成较大水位差变化所引起的。

采用靠近防渗墙的一排渗压计(UP1～UP5)的枯期测值(上、下游水位差较大)计算闸基渗压折减系数

α=(HSP-HRB)/(HR-HRB)

(1)

式中，HSP—渗压水头高程；HR—库水位；HRB—下游河床高程。

计算渗压折减系数基本在0.1～0.3之间，表明闸基防渗系统的防渗效果较好。

由于没有下游河道水位测值，故渗压测值周期性变化的原因还需进一步验证。此外，渗压计的数量有限，只能观测局部渗流情况，且只能获得渗压测值，无法确定闸基渗透坡降。因此，采用有限元法进行闸基渗流分析。为了使有限元计算结果更接近实际情况，先根据渗压监测数据率定闸基覆盖层材料的渗透系数。

2 闸基覆盖层渗透系数反演分析

2.1 闸基渗流有限元模型

采用ANSYS软件[11]建立1#闸孔和2#闸孔的闸基平面渗流有限元模型。模型包含闸底板、上游铺盖、下游护坦(含排水孔)、防渗墙、闸基覆盖层及下部花岗岩，如图4—5所示。护坦中的排水孔通过增大相应位置处的单元的渗透系数的方式模拟。有限元模型的模拟范围为：自上游铺盖末端向上游延伸75m，自下游护坦末端向下游延伸80m，自覆盖层底部向下延伸50.4m。为了对比计算渗压和实测渗压，根据渗压计的埋设位置，在有限元模型中相应位置布置结点(1#闸孔闸基中的渗压计为UP2和UP6，2#闸孔闸基中的渗压计为UP3和UP7)。

由于渗流分析和热分析的基本方程相似，因此采用ANSYS软件的热分析模块计算闸基渗流。采用PLANE55热分析单元离散有限元模型，1#闸孔闸基的有限元模型划分为6922个单元和7095个结点，2#闸孔闸基的有限元模型划分为6354个单元和6527个节点。模型上下游边界及底部边界为不透水边界。模型顶部边界水头取枯期和汛期的上、下游水位，上游水位采用实测值，枯期下游按无水考虑，汛期的下游水位根据水力学计算结果确定。

2.2 渗透系数反演分析方法

混凝土闸底板、混凝土铺盖、混凝土护坦、混凝土防渗墙和覆盖层下部花岗岩的渗透系数均有试验值。因为难以取到闸基覆盖层原状材料，以及材料颗粒尺寸大难以进行室内试验，所以，闸基覆盖层材料渗透系数的设计值是通过现场勘察并结合工程经验确定的，具有不确定性。而材料的渗透系数对闸基渗流计算结果影响较大，因此，根据渗压监测数据对覆盖层材料的渗透系数进行反演分析。

渗透系数的反演分析即寻找一组渗透系数使有限元计算的闸基渗压与渗压计测到的渗压最吻合。需要采用渗压的计算值和实测值建立目标函数，求解一个以渗透系数为未知变量的优化问题。采用有限元法计算渗压，既耗时又无法建立目标函数，给优化问题的求解造成困难。本文采用渗压计算值的响应面方程(即以渗压计算值为因变量和需要率定的渗透系数为自变量的函数)代替有限元计算，可以建立目标函数，提高参数反演分析的效率。

2.2.1渗压计算值的响应面方程

本文通过构造如下响应面方程[12- 13]，代替耗时的有限元计算，快速确定不同渗透系数对应的渗压计算值，极大提高参数反演分析的效率。

(2)

对需要率定的渗透系数进行组合。然后对每一组渗透系数，采用有限元法进行闸基渗流计算，提取渗压计埋设处的渗压计算值。最后通过回归分析确定响应面方程(2)的系数。

本文采用正交试验设计[14- 15]进行渗透系数的组合。以三种材料的渗透系数组合为例进行说明。每种材料的渗透系数取3个值，即μ、μ(1-kδ)、μ(1+kδ)，其中μ为渗透系数均值，δ为渗透系数变异系数，k为系数。然后根据正交试验表，对3种材料的9个渗透系数进行组合，得到9个试验点，如图6中立方体网格中的黑色圆点所示。由图6可见，立方体网格共有6个平面，每个平面上都有三个试验点，每条直线上都有1个试验点，并且这些试验点是均匀分散的，能很好地代表大量试验点的情况。

图6 正交设计示意图

2.2.2渗透系数反演分析的目标函数与求解

利用渗压计算值的响应面方程(1)和渗压的实测值建立如下渗透系数反演分析的目标函数：

(3)

(4)

渗透系数的反演分析就是在各种材料的渗透系数的上、下限范围内搜寻一组渗透系数使式(3)表达的目标函数值最小。因此，渗透系数的反演分析可以转化为带约束条件的非线性规划问题。本文利用LINGO软件[16]求解式(3)—(4)表达的非线性规划问题。

2.3 覆盖层渗透系数反演分析结果

根据UP2、UP3、UP6和UP7四个渗压计对2015年汛期和枯期的渗压监测数据，采用2.2节所述方法对漂卵石、卵砾石和砂壤土三种材料的渗透系数进行反演分析。需要率定渗透系数的材料有三种，故渗压计算值的响应面式(2)共有10个系数。采用回归分析确定这些系数需要10个以上的试验点。因此，进行两次正交试验设计，共得到18个试验点。第一次正交试验设计每种材料的渗透系数取μ、μ(1-δ)、μ(1+δ)3个值，第二次取μ、μ(1-3δ)、μ(1+3δ)3个值，其中μ采用渗透系数的设计值，见表1，δ取0.3[17]。采用SPSS软件[18]的正交试验设计模块对漂卵石、卵砾石和砂壤土三种材料的渗透系数进行组合。然后对每一组渗透系数，采用图4—5所示的有限元模型进行汛期和枯期两种工况下的闸基渗流计算，提取UP2、UP3、UP6和UP7四个渗压计埋设处的渗压计算值。最后利用SPSS软件的回归分析模块确定汛期和枯期工况下UP2、UP3、UP6和UP7四个测点处渗压计算值的响应面方程(总共8个)系数。8个响应面方程的回归分析的复相关系数均大于0.9，同时采用响应面方程计算的渗压和有限元计算结果的差别很小，最大相对误差小于0.34%，大部分相对误差小于0.1%。

表1 材料渗透系数的设计值和率定值

利用渗压计算值响应面方程和渗压实测值，按照式(2)—(3)建立渗透系数反演分析目标函数，渗透系数的上限设为μ(1+3δ)，下限设为μ(1-3δ)。采用LINGO软件搜索漂卵石、卵砾石和砂壤土渗透系数的最优组合，使目标函数值最小，结果见表1。由表1可见，漂卵石和卵砾石渗透系数的率定值与设计值有较大差别。采用覆盖层材料渗透系数的率定值进行闸基渗流计算。图7显示了UP2、UP3、UP6及UP7四个测点处的渗压计算值和实测值的对比。可以看出，渗压的计算值和实测值总体比较接近。计算值与实测值的差异可能是渗压计的测量误差和有限元模型误差(没有采用三维模型，此外，用于建模的闸基地质剖面是根据有限钻孔数据绘制的，与实际地层情况可能存在差别)造成的。

图7 汛期和枯期渗压计算值与实测值的对比

3 闸基渗流有限元分析

3.1 上下游水位变化对闸基渗压的影响分析

根据1.2节的分析，汛、枯期下游河道水位大幅变化是引起渗压计测值周期性变化的主要原因。本节采用有限元数值模拟，研究上下游水位变化对闸基渗压的影响，验证渗压计测值周期性变化的原因。首先将上游水位固定在944.5m，将下游水位从936m以0.5m的间隔升高至938.5m，计算UP2、UP3、UP6及UP7四个测点处的渗压，分析下游水位变化对渗压的影响。然后将下游水位固定在936m，将上游水位从942.5m以0.5m的间隔升高至945m，计算各测点处的渗压，分析上游水位变化对渗压的影响。测点处渗压随上下游水位的变化情况如图8所示，可以看出渗压随上游或下游水位的升高而增大，并呈线性关系。渗压随下游水位变化的斜率明显大于随上游水位变化的斜率，上游水位升高2.5m引起的渗压增大值不超过1m，而下游水位升高2.5m引起的渗压增大值为2～2.5m。因此，说明渗压计测值的变化主要受下游水位变化影响。

图8 上下游水位对渗压的影响

3.2 闸基渗流特性分析

采用图4—5所示的有限元模型(其中覆盖层材料渗透系数采用率定值)进行汛期(上游水位944.60m，下游水位938.50m)和枯期(上游水位944.59m，下游水位936.00m)两种工况下的闸基渗流计算。由于枯期和汛期两种工况下闸基渗压分布规律相似，只是量值不同。故仅给出枯期时1#闸孔的闸基渗流等势线和渗透坡降，分别如图9—10所示。

图9 枯期1#闸孔闸基渗流等势线

图10 枯期1#闸孔闸基渗透坡降

由图9—10可以看出，闸基渗压从上游至下游逐渐降低，且在防渗墙两侧渗压下降较快，说明防渗墙发挥了降低闸基渗压的作用。闸基三个部位的渗透坡降相对较大，即防渗墙底端附近、混凝土铺盖上游端附近和护坦上游排水孔附近，其中以防渗墙底端附近的渗透坡降最大。由于防渗墙靠上游布置，不能控制下游河水向上游渗透，因此，汛期下游水位大幅升高后，防渗墙后的渗压比枯期时明显增大。由于枯期上下游水位差较大，因此，相应的闸基渗透坡降也较大。上述三个部位在汛期的最大渗透坡降分别为0.302、0.129和0.148，枯期分别为0.425、0.181和0.209。除上述几个部位外，闸基其余部分的渗透坡降都很小，闸基渗流总体比较稳定。由有限元计算结果可见，防渗墙底端、混凝土铺盖上游端和护坦上游排水孔附近的闸基土体发生渗透破坏的风险较大。特别是在枯期运行时，需要加强巡视检查，观察是否存在闸前水库旋涡、护坦排水孔中析出细颗粒等异常现象。到目前为止，该水闸闸基无异常渗流现象。

4 结论

与基于监测数据的方法相比，数值模拟可以分析地基整体渗流状态。根据渗压监测数据率定材料渗透系数，再采用率定后的渗透系数进行渗流计算，可使计算结果更接近实际情况。

本文提出基于响应面法的渗透系数反演分析方法。通过构造渗压计算值的响应面方程，代替耗时的有限元计算，而且可以建立显式表达的目标函数，因此，能够提高渗透系数反演分析的效率。

结合监测数据和数值模拟对水闸的地基渗流进行了分析。监测数据表明该水闸地基渗压呈汛期高、枯期低的周期性变化规律，原因是汛、枯期下游河道水位大幅变化。根据渗压监测数据，对闸基覆盖层材料渗透系数进行了反演分析。采用率定后的渗透系数，对枯期和汛期的闸基渗流场进行了有限元数值模拟，验证了渗压周期性变化的原因，并揭示了闸基中渗透坡降较大的部位。本文提出的方法可用于类似水闸地基的渗流分析与安全评估。