基于温度示踪的河流潜流带地下水流速研究

潘维艳，杨姗姗，刘俊锋，钱秀红，徐征和

（1.济南大学水利与环境学院，山东 济南 250022；2.山东省海河淮河小清河流域水利管理服务中心，山东 济南 250014）

0 引 言

潜流带是流域系统动态水文过程与生物地球化学反应的关键带，潜流带内地表水与地下水的相互作用直接影响水在河床沉积带内的滞留时间、沉积带营养物质的通量、河床孔隙水的酸碱性及氧化还原条件［1，2］。在河岸潜流带内，地表水通过河岸沉积层与地下水发生物质和能量交换，生物活动强烈。地表水入渗过程中，水分运移过程能够携带污染物进入地下水中，而潜流带对地表水中的污染物具有重要的截污作用，因此，潜流带对地下水水量和水质变化过程发挥着重要作用［1，3］。潜流带也已成为近年来国内外学者研究的热点［4-6］。因此，定量分析潜流带水量与流速的时空变化特征是充分认识地表水-地下水交互作用规律的关键问题。

潜流带中时刻伴随着水热交换，温度作为载体，可以直观反映地表水和地下水流动过程及其交换关系，地表水与地下水交互作用的研究已经从水文学和水动力学法发展到温度示踪法［7，8］，温度示踪法检测成本低、灵敏度高、数据获得直接、稳定且无污染，在示踪地下水流速方面具有明显优势［9-11］。近年来，温度监测手段的进步和温度场研究的深入，利用温度示踪法进行流速定量计算的研究越来越多，如Caissie 等［12］运用温度示踪法评估并分析了季节性冰盖河流的地下水的流速和方向，张文兵等［13］通过温度示踪的潜流交换通量解析模型对比研究发现，浅层潜流带中潜流交换频繁，且深度越大，交换作用减弱。乔晓英等［14］将温度示踪应用于湖水与地下水交互作用，发现基于温度示踪的数值模拟法结果更理想。董林垚等［15］利用温度示踪法测算地下水流速发现地质体扩散率和地表温度是影响结果准确性的主要因子。温度是研究河流等地表水与地下水交换的较为有效的方法之一，但不同研究带的水文地质特征和水动力条件的不同对地下水流速有不同影响，目前，基于温度示踪法探讨河岸潜流带的研究成果还相对较少，仍需加强。

因此，本文选取再生水补给河湖的典型河道为研究区，定期监测地表水、地下水的水位和水温，利用水温和水位数据资料分析潜流带温度场的变化规律及其影响因素，并利用水动力法和温度示踪法计算流速，进行对比分析，探讨采用温度示踪法用于潜流带的适用性和合理性。再生水回灌区潜流带渗漏研究是再生水安全、高效回灌的关键，可以为再生水资源配置及地下水环境保护提供理论依据和决策支持。

1 研究区概况

本研究的实验场地选择在北京市朝阳区奥林匹克森林公园内的清洋河河岸（40°00′～40°02′N，116°22～116°24′E）（图1）。研究区多年平均气温11～12 °C，平均水面蒸发量在1 200 mm 左右，平均降雨量约600 mm，降雨年内分配不均，85%以上主要集中在汛期。清洋河河道补水全部来自再生水、降雨和径流，再生水补水期为每年的3-11月，年供给水量约1.5×106m3。清洋河是人工河道，采样分析发现研究区河岸带的土壤属于壤土，河床底部有防渗处理。

图1 试验区位置图Fig.1 Location of study area

2 材料与方法

2.1 试验布置与数据采集

在河岸带选取一处典型断面，沿垂直河水方向，在断面上布置1 处河水温度观测点，在距离河岸位置2.5 和5 m 处分别布置 2 组观测井（W1 和W2），每眼井分别在距离地面0.5、1、1.5、2、2.5、3、4、5、6、7、8、9 m 处安装热电偶热线，以观测包气带∕饱和带各深度处的温度，并在2 眼井附近布设水位探头（HOBO U20L-01，美国）。利用数据采集器（CR3000，美国）采集并记录各观测点温度，采集频率设置为4 次∕h，自2013年9月25日至2014年1月20日进行了连续观测。试验布置和现场试验图如图2所示。采用土钻法每月采集一次包气带土壤含水量，采样深度为0.5、1、1.5和2 m，采集时间为2013年9月-2014年1月。

图2 试验监测井布置和试验现场图Fig.2 Test monitoring well layout and test site diagram

2.2 入渗速率计算方法

采用的热运移扩散方程为［16-18］：

式中：T为温度，℃；t为时间，s；q为流体的流速，m∕s；Cw为水的体积热熔，J∕（m3·℃）；C是饱和沉积物体积热容，J（m3·℃）；z是深度，m；ke为饱和介质有效热扩散系数，m2∕s。

Hatch 等［16］提出了基于上述方程的流速解。Hatch 通过振幅或相位差计算的流速解析解分别为：

式中：q为垂向向下流速，m∕s；Ar为下侧传感器数据振幅与上侧之比；Δz为传感器间距离，m；v为温度锋速度，m∕s；Δt为信号传播的时间差，s；P为温度信号波动周期，s；α定义为：α=ke为等效热扩散系数；ke定义为：β|Vf|；λ0为水流的热传导系数；β为热弥散度；Vf为温度前端运移速度，m∕s；H定义为H=C∕λ0，在上述解析解中通过振幅法计算的解析解包含流速的大小和方向，通过相位计算得到的解析解只有流速大小。

本文利用Hatch 相位和振幅两种方法对温度时间序列进行分析并计算不同深度处的水流流速。计算过程中所需的参数主要通过试验和相关文献料获得，参数见表1。采用VFLUX2.0程序针对上述方法进行流速计算。VFLUX2.0 是一种基于matlab平台，利用热传输方程计算饱和多孔介质中一维垂直流体流动（渗流通量）的程序［19］。文中选择温度波动相对强烈的时段为例，进行流速计算，计算时段为2013年11月7日至11月17日，共计240 h，Hatch位相法计算所用的温度采用频率为0.5 次∕h。水动力学法基于达西定律进行计算［20］，根据多点原位取样发现研究区河岸潜流带质地比较均一，通过变水头入渗试验测得平均渗透系数，然后按照达西定律计算监测井附近的流速［21］。

表1 相关计算参数Tab.1 Related calculation parameters

3 结果与讨论

3.1 河水深度与地下水埋深变化特征

试验期间河水深及地下水井W1 和W2 的地下水位埋深随时间的变化如图3所示。河水水深波动范围为1.0～1.8 m，W1的地下水埋深波动范围为在1.5～3.0m，W2 处的地下水埋深波动范围为2.5～5.5 m。在观测期内，W1和W2处地下水埋深整体上呈增大趋势。其中，W1和W2处地下水埋深在2016年10月5日至15日呈先减小后增大的趋势，在2017年1月3日至1月23日地下水埋深呈现逐渐增大趋势，W1 处地下水埋深从2.2 m 增大到2.8 m，W2 处地下水埋深从3.5 m 增大到5.1 m。从图3中可看出，W1和W2的地下水埋深随时间变化呈现相似的变化规律，地下水埋深的变化主要取决于河水位的波动，即河水水位减小导致地下水水位抬升，并呈现一定的滞后性。结果表明，在观测期间，河水位持续高于地下水位，河道再生水持续补给地下水。

图3 地下水埋深与河水水深随时间变化曲线Fig.3 Temporal variation of groundwater table and river table

3.2 河岸带温度场时空分布规律

本次试验对W1 和W2 两处监测井进行了持续近4 个月的温度观测，由于设备故障等原因造成部分时段数据出现缺失（如图4、5 中两处缺口部分）。根据采集的温度数据显示，2 眼井中的温度数据变化规律和趋势一致，仅在数据上有差异，选择数据更加完善的W2 数据进行分析。图4是W2 中不同各深度处温度随时间变化的曲线，其中，0.5、1、1.5和2 m观测点长期处于地下水位以上，属于包气带层的地温，6、7、8 和9 m 四个观测点长期处于地下水位以下，为地下水温度。

由图4可知，河岸包气带地温在时间上，2016年9月-2017年1月，不同深度的包气带地温随时间逐渐减小，温度波动范围分别为23.4～3.4 ℃（0.5 m），24.5～6.0 ℃（1 m），24.8～7.8 ℃（1.5 m）和25.0～10.2 ℃（2 m），且整体上呈现线性下降趋势（各R2均大于0.96）。

河岸包气带温度在空间上，包气带平均地温随着深度的增加呈现逐渐增加的趋势，平均温度从0.5 m 处的10.8 ℃增加到2 m 处的16.9 ℃。观测期内，气温波动剧烈（图4），包气带地温（2 m 以内）随气温波动变化明显，两者呈显著正相关关系（R2=0.87～0.89，图5）。结果表明，气温是影响包气带土壤温度的主要因素。在4 个月的观测期内（9月-次年1月），除个别时间点外，地温均高于气温，其中，表层地温最接近气温，随着深度增大，地温与气温的差异越大，地温和气温的平均差异在5～11 ℃范围内。

图4 W2中不同深度处温度随时间变化规律Fig.4 Variation of temperature at different depths with time in W2

图5 W2处包气带地温与气温之间的相关性Fig.5 Correlation between ground and air temperature in the vadose zone in W2

包气带不同深度处的土壤含水量与包气带地温之间的相关性分析如图6所示，结果发现，在0.5～2.0 m 范围包气带土壤质量含水率与包气带地温之间呈现显著的正相关关系（R2=0.680 4），包气带地温随着含水量的增加而增大，在观测期内，当含水率在17%～26%内变化时，相应深度处的包气带温度在1.9～25.0 ℃范围内变化，可见，土壤含水量是影响包气带温度的因素之一。有研究指出降水也是影响包气带温度的重要因素之一，由于试验期间研究区降水量较少，故未考虑其影响。

图6 土壤质量含水率与包气带地温之间的关系Fig.6 Relationship between soil mass moisture content and ground temperature in vadose zone

由图7可知，河岸饱和带地下水温度呈现出与包气带地温不同的时空变化特征。从时间上看，与包气带地温变化规律不同，饱和带地下水温度整体上呈先现逐渐增加后逐渐下降的变化趋势，但不同深度呈温度变化的拐点时间略有差异，其中，6、7 和8 m 三处的地下水温度自2016年9月25 至10月23 呈缓慢增加趋势，随后至次年1月中旬，温度缓慢下降；9 m 处，地下水温度自2016年9月25日至12月14日，呈缓慢增加趋势，随后至1月中旬又缓慢下降，不同时间节点温度波动范围分别为24.0～24.2～16.9 ℃（6 m），21.8～22.3～17.6 ℃（7 m），19.0～20.1～18.0 ℃（8 m），16.2～18.8～18.1 ℃（9 m）。4 处观测点处的地下水平均温度波动范围为17.4～21.4 ℃，非饱和带平均地温波动范围为6.8～24.4 ℃，可见，地下水温度随时间的波动性比非饱和有所减弱。

图7 W2处地下水温度随时间的变化过程Fig.7 Variation of groundwater temperature with time in W2

空间上，从包气带到过渡带到饱和带，不同深度处平均温度分别为10.8 ℃（0.5 m）、13.7 ℃（1 m）、15.4 ℃（1.5 m）、16.9 ℃（2 m）、18.2 ℃（3 m）、20.0 ℃（4 m）、21.1 ℃（5 m）、20.9 ℃（6 m），20.3 ℃（7 m），19.2 ℃（8 m），18.1 ℃（9 m）。在不同季节呈现出不同的变化趋势，在秋季，温度随着深度增加呈现先增后减的趋势，表现为上冷中暖下冷的分层现象；在冬季，温度随着深度增加呈现逐渐增加的趋势，表现为上冷下暖的分层现象。在观测期内，地下水平均温度随深度变化规律与包气带平均地温相反，即地下水平均温度呈随深度的增加而下降的趋势。

3.3 水动力学法与温度示踪法计算结果对比

水动力学法和Hatch 相位法、振幅法计算得到的流速结果如图8所示。水动力学法计算得到的流速值变化范围为6.74～7.74×10-5m∕s，最大值出现在第102 h，流速整体变化幅度不大，这主要是因为该河段在试验期间河水水位波动不大。Hatch 相位法计算的流速值变化范围为6.34～8.36×10-5m∕s，最大值出现在114 h。

图8 水动力学法和温度示踪法结果Fig.8 Results of hydrodynamic method and temperature tracer method

从数值上来看，水动力学法的计算结果与Hatch 相位法的计算结果相似。但是，与水动力学法相比，Hatch 相位法对计算的流速值变化波动范围和变化幅度更大，且波动趋势呈现一定的滞后性，这是因为温度在多孔介质中的传到和扩散速度要滞后于地下流速［21］。Hatch 振幅法计算的流速值变化范围为0.61～1.01×10-5m∕s，相位法计算的流速值要明显大于振幅法，相差接近一个数量级，两者的变化趋势相似。振幅法计算结果为正值，表明河岸带的渗透流速均为垂向向下，即河水补给地下水［22］，这与实际情况相符。对比河水位、地下水位与流速的关系还发现，当河水与地下水水位差较大时，对应的流速值相对较大，结果表明，当河水补地下水时，河水与地下水的水位差对研究河段河水入渗率的变化有一定影响。根据水动力学法的计算结果可以发现，Hatch 相位法计算流速的准确性要高于Hatch 振幅法，这与李英玉等［9］的研究结果相似，这主要取决于两种计算方法解析解表达式的差异。因此，后文的比较分析均采用Hatch相位法计算的流速值结果。

3.4 河岸带不同深度处流速对比

用Hatch 相位法计算得到的河岸带不同深度处的q1流速（距地面2.5 m）、q2流速（距地面3 m）和q3流速（距地面4 m）如图9所示。由图9可知，3 个深度处的流速变化范围为5.5～8.3×10-5m∕s，最大值出现在浅层q1流速处，最小值出现在深层q3流速处。不同深度流速的时程变化趋势在大部分时段相似，但不同深度处流速值存在一定差异，总体来看，流速值随着深度的增加呈减小趋势，这是因为随着深度增大，水压消散和水头损失增大，从而导致流速减小。

图9 不同深度处流速随时间的变化曲线Fig.9 Variation curves of flow velocity with time at different depths

对比q1流速处和q3流速处的变化趋势发现，虽然不同深度处温度的监测频率一样，但温度示踪法对浅层流速变化的刻画描述要比深层好，这是由热传导和热扩散过程中温度波动信号的衰减导致的［21］。因此，温度示踪法用于潜流层地下水流速的计算能够较好的反映有水位波动引起的浅层地下水流速时程变化趋势，并均有较好的准确性。

4 结 论

（1）从包气带到过渡带到饱和带，温度在不同秋季和冬季呈现不同的变化特征，在秋季，温度随着深度增加呈现先增后减的趋势，表现为上冷中暖下冷的分层现象；在冬季，温度随着深度增加呈现逐渐增加的趋势，表现为上冷下暖的分层现象。在观测期内，地下水温度随深度增加呈下降趋势，这与河岸包气带地温变化规律相反，且地下水温度随时间的波动性比包气带有所减弱。气温和土壤含水量与包气带地温成正相关关系。

（2）分别采用水动力学法和Hatch 相位法、振幅法进行地下水流速计算，结果发现，水动力学法计算得到的流速值变化范围为6.74～7.74×10-5m∕s，Hatch 相位法计算得到的流速值变化范围为6.34～8.36×10-5m∕s，Hatch 振幅法计算得到的流速值变化范围为0.61～1.01×10-5m∕s，在计算时段内，3 种方法计算得到的流速的变化趋势在大致相似，温度示踪法计算结果呈现一定的滞后性。对比发现，对于研究区来说，Hatch 相位法用于河岸带潜流层地下水流速的计算更合理更准确。

（3）比较不同深度处的地下水流速发现，深度越大，地下水流速越小，且温度示踪法对流速的刻画描述效果越差。