库区洲滩潜流带温度示踪流速计算方法

姬雨雨， 施文卿， 陈求稳， 安 磊

(1.河海大学 水利水电学院， 江苏 南京 210098； 2.南京水利科学研究院， 江苏 南京 210024)

1 研究背景

潜流带是河流地表水和地下水动态交换的饱和过渡区域，是河流中溶质交换、生物新陈代谢与繁殖的重要场所，在维持河流生态系统健康和稳定过程中发挥着重要作用[1-3]。为应对全球能源需求和气候变暖，世界各国推进水电开发[4]，大坝运行导致水位异常波动和泥沙冲淤平衡的改变，在库区和坝下形成众多河湾洲滩潜流带。水库昼夜发电引起水位波动，促使洲滩内部潜流交换和物质交换频繁发生，已逐渐成为国内外生态水文研究的热点区域[5-8]。潜流交换对潜流带内物化性质(溶解性气体、温度、pH、盐度、营养盐等)的空间分布和河流水质具有重要影响[6,9-11]。因此，量化潜流带的潜流交换流速和交换量成为了解其生物地球化学过程的重要基础和关键。

在潜流交换和太阳辐射的双重作用下，潜流带内的温度分布具有时空差异性，该温度差异可作为天然的示踪剂，用于表征地下水的流动。温度示踪法具有无污染、易于观测和准确度较高的优点，在计算潜流交换流速和交换量方面的应用广泛。Alexander等[12]利用连续监测的动态温度数据，估算水文地质参数并量化潜流交换量，Constant等[13]利用温度示踪确定了潜流交换的深度范围，Halloran等[14]利用温度示踪法计算了土体含水饱和度。然而，国内在潜流带研究中对于温度示踪法的应用仍处于起步阶段。

因此，本文以澜沧江漫湾库区洲滩潜流带为研究对象，水库运行期间实时监测洲滩内部水位和水温的变化过程，利用温度示踪法计算潜流交换流速和交换量，丰富温度示踪法在库区河流潜流带中的应用。

2 研究区域概况

选取云南省澜沧江漫湾库区的一处河湾洲滩为研究对象(图1)，漫湾水库位于云南省临沧市境内，建成于1993年，为澜沧江上第一座以发电为单一目标的百万千瓦级电站。河湾目标洲滩(24°43′44″N, 100°23′5″E)位于坝前30 km处，为多年泥沙淤积形成(图2)，洲滩长、宽和淤积深度分别为182、90和60 m。试验区具有两个特点：(1)因水库昼夜发电，库区水位波动频繁，洲滩内部潜流交换活跃；(2)试验区地处云贵高原，太阳辐射强烈且四季气温差异较小。

图1洲滩位置图2试验区沉积物层

3 研究方法

3.1 水位和水温监测

沿库区至洲滩方向布置W1 ～ W4四口水位和水温监测井，距离洲滩边缘分别为0.5、10.5、20.5、35.5 m(图3)。监测井中安放长度均为4 m的盲管，侧面每隔15 cm钻有直径为2 cm的小孔，并外覆以尼龙滤网，以实现井管内腔与外界水流通畅交换。监测井底部安装有水位自动记录仪(U2000101, OneSetHoBo, USA)，用于监测库水位昼夜波动，其精度和适用范围分别为0.5 cm和0 ～ 9 m，沿监测井底部向上每隔60 cm布置一个水温自动记录仪(U22001, OneSetHoBo, USA)，用于监测库水位波动期间地下水温随时间发生的波动，利用水准仪确定监测井的安装高程以保证具有共同的基准面。监测自2016年12月11日开始至2016年12月16日终止，时间间隔为10 min。

3.2 温度示踪流速计算方法

地下浅层温度场受气象和水力条件影响表现出周期性波动，热量作为一种天然示踪剂，可用于计算地下水流速的大小和方向，热量在潜流带中的运移数学模型：

(1)

式中:T为温度，℃；t为时间，s；z为沉积物深度，m；q为水流沿深度方向在孔隙介质中的流速(向下为正)，m/s；ρw为水的密度，kg/m3；ρ为含水介质的密度，kg/m3；cw为水的比热容，J/(kg·℃)；c为含水介质的比热容，J/(kg·℃)；λe为含水介质的热传导系数，J/(s·cm· ℃)。

Hatch和Keery分别得到了热传导方程的解析解，但Keery法只考虑了相位滞后，其解析解仅是流速大小，无法判定方向。Hatch法则综合考虑了波动振幅的衰减、相位滞后和热弥散度，其解析解可判定流速的大小和方向，且具有更高的精确度[15-16]，因此本文采用Hatch法计算流速。Hatch法一维热传导方程的求解公式为：

(2)

图3 研究区域监测井布置图

本渗流模型采用VFLUX计算包计算，该解析解适用于均质多孔介质。通过参考相关资料和室内实验，确定计算过程中所需的有效热传导系数、热弥散度、水比热容、土比热容及洲滩沉积物饱和层孔隙率参数(表1)。

4 结果与讨论

4.1 水位和水温变化

洲滩监测井的水温波动数据规律表现一致，仅在大小上存在差异。W1监测井靠近库区，且W1靠近岸边没有草甸覆盖，受水库水位波动和太阳辐射的双重作用更加明显，因此选择W1监测井水位和水温数据进行分析。监测周期内研究区域的气温和河流水温分别在12～30℃ 和15～16℃之间，洲滩内部水温在垂向上表现出“上凉下暖”。在涨水的0～33 h内，受水库持续蓄水和昼夜发电用水影响，监测井水位呈现持续波动升高趋势(图4(a))，因低温地表水渗入洲滩内部，表层、中层、底层水温均随水位上升而下降，分别从开始时刻的17.9、18.2、18.3℃下降至33 h的17.2、17.8、18.1℃(图4(b))，该时段W1水位高于W2水位，表现为库区持续补给洲滩，在33 h后，由于水库持续泄水，W1和W2监测井水位随之下降(图4(a))，低温地表水渗入量逐渐减小，底层水温在第40 h后开始逐渐上升(图4(b))。受太阳辐射和底层水温传导扩散的双重影响，表层和中层水温分别在第40 h和50 h开始上升(图4(b))。在63 h后，W2水位高于W1水位(图4(a))，表现为洲滩补给库区，低温地表水逐渐潜出，使得表层和中层水温在63 h后均呈较快上升趋势，最终分别稳定在18.4℃和18.5℃(图4(b))。

4.2 水动力学方法

多点原位深层采样发现洲滩为河流细颗粒泥沙淤积形成，质地比较均一，因此通过变水头入渗实验测得平均渗透系数为8.64×102m/d，并按照达西定律计算W1监测井附近的潜流交换流速和交换量(图5)。在0～63 h，库区向洲滩补给，流速呈现先增加后减小趋势，从开始时刻的7.21×10-5m/s 上升为8 h的13.49×10-5m/s，之后呈现递减趋势，在63 h处降至0 m/s，监测井间单位宽度上的补给总量为1.85 m3。在63～120 h，洲滩向库区排泄，流速呈先增加后稳定趋势，从63 h的0增加到115 h的5.15 m/s，监测井间单位宽度上的排泄总量为0.76 m3。

4.3 温度示踪法流速计算

采用VFLUX计算包中的Hatch法计算潜流交换流速和交换量。振幅法计算结果表明，0～66 h为库区补给洲滩阶段，流速呈先增大后减小趋势，从开始时刻的3.32×10-6m/s逐渐增大至30 h的 12.81×10-6m/s，之后减小至81 h的0，单位宽度上的补给总量为1.92×10-1m3。81～115 h为洲滩补给库区阶段，流速呈逐渐增大的趋势，从81 h的0逐渐增加至115 h的2.89 m/s，单位宽度上的排泄总量为0.21×10-1m3(图6(a))。相位法计算结果表明，0 ～ 81 h为库区补给洲滩阶段，流速呈先增大后减小趋势，从开始时刻的3.32×10-5m/s逐渐增大至30 h的 12.92×10-5m/s，之后减小至66 h的0 ，单位宽度上的补给总量为1.87 m3。66～115 h为洲滩补给库区阶段，流速呈先增大后缓慢减小至稳定的趋势，从66 h的0逐渐增加至80 h的5.87×10-5m/s，之后缓慢减小并稳定在115 h的3.92×10-5m/s，单位宽度上的排泄总量为0.79 m3(图6(b))。在库区水位涨落周期内，振幅法和相位法均呈现出库区与洲滩间先补给后排泄的过程，但二者潜流交换流速和交换量的数量级均相差10倍，对比水动力学法计算结果，Hatch相位法的结果准确度更高，所以之后的对比均采用Hatch相位法。

4.4 水动力学方法和温度示踪法计算结果对比

将水动力学法和Hatch相位法计算结果进行对比(图7)。库区补给洲滩阶段，在0～8 h，水动力学法计算的流速逐渐增加，在8 h时达到峰值流速13.49×10-5m/s，Hatch相位法计算的流速在0～30 h呈现逐渐增加的趋势，在30 h时达到峰值流速12.92×10-5m/s；洲滩补给库区阶段，水动力学法计算的流速在80 h时达到峰值流速5.52×10-5m/s，Hatch相位法计算的流速在83 h达到峰值流速5.87×10-5m/s。

可见，利用温度示踪法计算因库区水位波动引起的潜流交换时，流速的峰值存在一定的延迟，这是因为温度在多孔介质中的扩散和传导速度滞后于地下水流速。

对比流速曲线，水动力学计算的流速曲线波动明显，这主要是因为水库发电引起频繁性的水位波动，且具有较为密集的数据采集频率0.1 次/min，而温度示踪Hatch相位法计算的潜流交换流速波动性较小，主要是由于温度示踪法有数据重取样的环节，其数据采集频率为0.5 次/h，数据量相对较少，因此温度示踪计算的流速曲线波动性较小。总体来看，Hatch相位法得到的流速范围为-5.87×10-5～12.92×10-5m/s，单位宽度上潜流交换总量为2.66 m3，水动力学法计算的流速范围为-5.52×10-5～13.49×10-5m/s，单位宽度上潜流交换总量为2.61 m3，二者较为接近。因此，温度示踪法作为一种新兴的计算潜流交换流速和交换量的方法，能够较好地反映水位波动引起潜流交换的流速趋势，并具有较高的准确性。

表1 相关计算参数

图4库区洲滩水位和水温变化特征图5水动力学法流速计算结果

图6 Hatch振幅法和相位法计算得到的流速曲线图7水动力学方法和温度示踪法计算结果对比

4.5 不同深度流速对比

为比较洲滩不同深度处的流速大小，利用Hatch相位法计算了地面下1.1、1.7 m深度处的流速，1.7 m深度处的流速值q2比1.1 m处的流速值q1小(图8)，主要是因为水位波动引起的潜流交换在潜入和潜出洲滩的过程中存在水压消散和水头损失，即深度越深，地下水流速越小。另外，利用温度示踪法计算洲滩不同深度处的流速时，对浅层流速变化过程的刻画要比深层好，这是因为温度波动信号在沿地层深度方向的传导和扩散过程中存在信号衰减。

图8 不同深度地下水流速对比

5 结 论

(1)温度示踪法计算结果中，Hatch相位法和振幅法得到的流速曲线趋势基本一致，仅在水库落水的第80 h相差较大。对比计算结果，水动力学法计算的潜流交换流速和交换总量分别为-5.52×10-5～13.49×10-5m/s和2.61 m3，Hatch相位法计算的潜流交换流速和交换总量分别为-5.87×10-5～12.92×10-5m/s和2.66 m3，Hatch相位法具有更高的准确性。

(2)在库区补给洲滩和洲滩补给库区两个阶段，利用水动力学法计算得到的流速分别在8和79 h达到正负峰值，温度示踪法计算的流速分别在30和83 h达到正负峰值。温度示踪法计算流速达到峰值的时间存在一定滞后。

(3)温度示踪法计算结果表明洲滩上层的流速比下层大，这说明离地越深，传导的温度波动信号越衰减，流速越小。