南沟门水库渗流计算分析

庞耀立，董政良，张晓飞

（1.陕西省水务清洁能源集团有限公司，陕西 西安 710082；2.西安理工大学水利水电学院，陕西 西安 710048）

1 工程概况

南沟门水库枢纽工程位于陕西省延安市黄陵县境内，水库坝址位于葫芦河河口上游约3 km 处的寨头河村南沟门附近。南沟门水利枢纽工程建设任务为工业和城乡生活供水，兼顾灌溉和发电等综合利用。南沟门水库坝址断面多年天然径流量为1.31 亿m3，100 年一遇洪峰流量为1130 m3／s，三日洪量为3410 万m3；5000 年一遇洪峰流量为2910 m3／s，三日洪量为6430 万m3。南沟门水库属Ⅱ等大（二）型工程，坝设计防洪标准为100 年一遇，校核洪水标准为5000 年一遇。工程区地震基本烈度为Ⅵ度。大坝是均质土坝，坝顶高程为852.0 m，大坝坝顶宽度为10 m，在其坝顶上游侧设置了1.2 m 高的防浪墙，墙顶高程为853.2 m，大坝最底点开挖高程789.0 m，最大坝高为63.0 m，坝顶总长为504.4 m；坝体上游边坡1∶2.75 和1∶3，在高程830 m 处设宽3.0 m 马道；下游坝坡比均为1∶2.5，在高程835 m 和818 m 处设宽2.0 m 的马道，在800.5 m 处设棱体排水平台，最大坝底宽为351 m。坝体上游坝坡在高程812.5 m 以上采用厚0.4 m 的干砌石护坡，下游坝坡采用浆砌石网格内植草皮护坡。库区地处陕北黄土高原强烈侵蚀区，河流切割基岩深达50 余米，两岸沟谷纵横，塬、梁、峁黄土地貌景观明显，相对高差200 m～300 m。河流阶地发育不对称，凹岸岩石裸露，地形陡竣；凸岸阶地发育，第四系松散堆积覆盖广泛。库盆周边岩体完整，无大的构造现象，库岸基岩顶板及地下水位高于水库正常水位，无渗漏条件，不存在永久性渗漏问题[1］。

2 计算模型及参数

2.1 模型基本假定

土石坝渗流场与温度场之间是相互影响、相互作用的，一方面，温度影响水的密度和动力粘度，进而影响土的渗透系数，热量的运移又会导致土中水分的迁移；另一方面，渗流作用引起土体中的热量发生转移，进而影响温度场的空间分布。两场耦合作用非常复杂，为了简化求解同时又能体现关键要素，且能保证计算结果可靠的情况下，做以下基本假设：

1）土石坝属于多孔连续介质，土石坝内的渗流是连续介质渗流；

2）水流运动过程中只受重力作用；

3）温度不会导致水体发生相变；

4）不考虑因介质摩擦产生的热量；

5）水流符合达西定律。

2.2 渗流场与温度场耦合数学模型

不考虑流体与多孔介质之间热平衡所需时间，即认为两者之间的热传递是瞬时完成的，忽略热辐射、化学反应等引起的热消耗，且流体密度仅和温度及压强有关，得到三维多孔介质地下水的水-热运移耦合方程[2］：

温度场：

式中：QG为温度场的源汇项；Γ1和Γ2分别表示第一类边界和第二类边界；Φ0为第一类边界上定温度值；n 为第二类边界上的外法线方向；Q0为第二类边界上的单宽热流量；T0为初始温度分布；t 为时间；Ω 表示研究区域。

渗流场：

式中：Qs为渗流场的源汇项；P1为第一类边界上水头高度；P0为初始水头分布；q1为第二类边界上的通量；其他变量如前所述。

2.3 边界条件及参数

2.3.1 边界条件

渗流场方面，坝体上游定水头为 852 m，边界为上游坡面，坝体下游定水头789 m，边界为下游坡面，其余边界均为零通量边界。

温度场方面，当地温度随时间周期性变化的傅里叶级数函数如下：

式中：T0为平均温度；k 代表傅里叶级数的阶数；Ak和Bk代表傅里叶级数的系数；L 为周期[3］。

这里采用陕西省延安市的气温数据（数据来源于中国气象数据网站）并根据傅里叶函数来取5 阶傅里叶级数就可以满足精度要求。拟合得气温回归模型系数为：T0=16.37、A1=-10.28、A2=-0.5417、A3=0.2667、A4=-0.2917、A5=-0.08822、B1=-7.133、B2=-0.2165、B3=0.3、B4=-0.0433、B5=0.08344。将系数带入上式即可求得温度随时间周期性规律，温度为16.37 ℃。

2.3.2 模型参数

由于本文采用COMSOL 软件进行数值模拟，而COMSOL中自带计算传热模块是用热量运移方程进行计算的，默认为导热系数模型为体积平均模型，因此用PDE 模块来代替COMSOL 自带的传热模型进行计算，可以灵活修改导热系数经验模型，并通过修改相关的传热模块来实现不同导热系数条件下的土石坝流热耦合模型求解。模型渗流场及温度场的计算参数[4］，见表1，模型构图见图1[5］。

表1 土石坝流-热耦合模型计算参数

续表2

图1 模型结构示意图

3 计算成果分析

生成的结果见图2，可以发现上游坝体基本处于饱和状态，帷幕灌浆区域轴线附近地带的坝体中间部位处于饱和与非饱和的交界状态，下游坝体则是大部分处于非饱和的状态。

由图2 可知由于帷幕灌浆区域渗透性比较低，等压线相对密集，因此地下水渗流的流速极其缓慢，坝体浸润线（图中零压力线）穿过帷幕灌浆区时急剧降低。该区域前坝体在库水的作用下，其压力水头分布较为平缓。

图2 2019 年6 月典型时刻孔隙水压力水头分布图

在下游区域，坝体浸润线逐渐降低变缓，与排水褥垫层的上表面高程基本一致，说明穿过帷幕灌浆区域的库水主要通过排水褥垫层排出坝体，结合下游坡脚处的排水棱体起到了很好地防渗作用。因此，帷幕灌浆能有效降低坝体浸润线，穿越帷幕灌浆这部分渗漏库水能通过排水褥垫层和排水棱体排出，帷幕灌浆和排水棱体能有效保证土石坝的运行安全，计算结果合理。

4 结论

堤坝流热耦合理论仍然是一个年轻的研究领域，模型和方法的发展至关重要。本文通过理论分析堤坝温度场与渗流场相互影响，建立了土石坝饱和-非饱和流热耦合的数学模型，并选取南沟门水库拦河坝进行土石堤坝的渗漏研究，采用COMSOL Multiphysics 建立有限元模型来分析堤坝温度场和渗流场的变化规律，并通过修改相关的传热模块来实现不同导热系数条件下的土石坝流热耦合模型求解。发现渗流未达到稳定时，坝体浸润线在经过帷幕灌浆区后依然急剧降低，下游浸润线的高程与排水褥垫层基本相同，非饱和带饱和度由坝体浸润线至坝体外表面逐渐减小。结果可为今后陕西省土石坝工程的防渗提供参考依据。