南岭山地森林流域退水规律及影响因素

李泽华，周 平，黄远洋，徐 卫，谭兆伟

（1. 广东省科学院广州地理研究所广东南岭森林生态系统国家野外科学观测研究站/广东省地理空间信息技术与应用公共实验室/广东省遥感与地理信息系统应用重点实验室，广州 510070；2. 重庆地质矿产研究院自然资源部重庆典型矿产生态修复野外科学观测研究站/重庆市万盛矿区生态环境保护修复野外科学观测研究站，重庆 401120）

对于包括广东省的中国大部分地区而言，降水的年内分配不均导致年径流中绝大多数在丰水期以洪水形式入海而难以利用（黄楚珩等，2019）。与丰水期径流相比，枯水期径流因为在维持河流生态系统稳定和保障下游供水安全等方面起到更重要的作用，而受到越来越多的关注（Bruijnzeel, 2004；曾松青等，2010）。在枯水期内由于降水持续偏少甚至无降水的时间较长，径流表现为逐渐消退的流域退水过程，该过程中径流成分以地下水出流为主，主要受流域特征和水文地质特性影响（穆文相，2020）。Brutsaert 等（1977）通过对流域概化成具有不透水底板的含水层，根据Dupuit-Boussinesq 方程的解析解和枯水期径流消退特征，将流域水文地质参数与退水过程建立联系，其中提出的利用多事件的蓄泄关系推求消退系数的方法得到广泛应用。随后有学者发现利用多事件的蓄泄关系识别存在较大的不确定性，并转而采用基于独立事件的蓄泄关系分析流域退水过程的变化规律及影响因素（Shaw et al.,2012；高满等，2013）。然而，在以往基于独立事件的蓄泄关系研究中，针对退水过程的影响因素，不同学者得出的结论存在较大差异。有学者将主要影响因素归结为蒸散发（Federer,1973;Cadol et al., 2012;Tashie et al., 2020），也有学者认为土壤水分起着关键作用（Shaw et al., 2012;2013），还有学者提出地形特征和河网结构等主导了退水过程的变化（Biswal et al., 2010, 2014; Patnaik et al.,2015）。可见，有关退水规律和影响因素的结论多有不同，对于地形、植被、土壤和地质结构等下垫面条件更趋复杂的山地森林流域而言，其退水过程在目前更是缺乏完整的机理解释和系统的实证研究，有关森林是否具有“削峰补枯”作用的争论也一直在持续（陈军锋等，2001；Wang et al.,2016; Ellison et al., 2017; Li et al., 2021; Zhao et al.,2022）。

南岭保存有地球同纬度带上最完整的亚热带常绿阔叶林，是中国南方重要的生态屏障带，也是北江和东江等众多河流的源头区，因为被称为“广东水塔”，对粤港澳大湾区的供水安全起到至关重要的作用。南岭山地的地形地貌复杂性和生态系统多样性为代表水文学发展方向的生态水文学提供了重要的研究对象（周国逸等，2018）。根据广东南岭森林生态系统国家野外科学观测研究站降水观测数据显示，广东南岭国家级自然保护区2021年降水量1 231 mm，与2019年降水量2 592 mm 和2020年降水量2 333 mm 相比减少1 000 mm 以上。在此背景下，本研究基于广东南岭森林生态系统国家野外科学观测研究站获取的一手观测数据，采用基于独立事件的蓄泄关系分析方法，揭示南岭典型的山地森林流域退水规律，识别影响退水过程的主要因素，以期为南岭山地流域的枯水期径流模拟及预测、水源涵养功能评估、河流生态系统保护与修复等提供科学依据。

1 材料与方法

1.1 研究区概况

研究区位于广东南岭国家级自然保护区，属亚热带湿润季风气候，集中分布着中国最典型的中亚热带常绿阔叶林，土壤类型以黄壤土为主。依托广东南岭森林生态系统国家野外科学观测研究站，位于保护区乳阳片区内的实验流域（经纬度范围为24°55'12"—24°56'06"N 和112°59'46"—113°-0'32"E，集水面积约0.86 km2，海拔范围为1 115～1 775 m）布设有测流堰（图1），采用6541型浮子水位计监测流域出口断面水位变化，在实验流域内采用基于TDR原理的CS655型土壤水分传感器监测不同深度 （5、10、15、20、40、60、100 cm） 土壤水分体积含量变化，采用CS451型压力式水位计监测地下水位变化；采用IRGASON-NM-BB-IC 型开路式涡度相关系统监测实际蒸散发变化，并同步监测降水、气温等气象要素变化；以上数据监测频率均设定为10～30 min。

图1 研究区示意图Fig.1 Map of the study area

1.2 数据来源

利用在测流堰采集的水位数据，通过堰流水力学公式计算获得实验流域出口河道断面流量数据；对不同深度的土壤水分进行加权平均，获得0～100 cm深度层的土壤水分体积含量数据；根据地下水监测井压力式水位计的布设埋深（14.0 m），将监测的地下水位数据转换为地下水埋深数据；将以上数据及同步监测的气象数据按时间频率24 h 统一整理，最终获得逐日平均流量、土壤水分、地下水埋深、潜在蒸散发、实际蒸散发等的时间序列。南岭实验流域监测数据基本信息见表1所示。

表1 南岭实验流域监测数据基本信息Table 1 Description of the data observed at the Nanling experimental catchment

1.3 研究方法

天然流域河道径流消退过程中，流量Q作为时间t的函数，其导数（即随t的变化率） 与流域的调蓄特性有关：

式中：Q为流量；t为时间；f为反映流域调蓄特性的函数，亦称流域退水斜率曲线。当监测时间频率为Δt，连续观测的流量数据为Qi和Qi+1，式（1）可近似为：

对于具有不透水底板的含水层，在忽略源汇项情况下，潜水水流运动可通过Dupuit-Boussinesq方程描述：

式中：Sy为给水度；K为饱和渗透系数；x为横坐标，表示在含水层剖面的水平位置；h为地下水深度。Brutsaert等 （1977） 基于对 Dupuit-Boussinesq方程的推导发现流域退水斜率曲线呈现幂律函数形式，可表示为：

式中：a和b为消退系数，分别对应于流域退水斜率曲线在双对数图中的截距和斜率。其中，截距a主要反映径流的消退速率，a值越大，消退速率越快；斜率b则是曲线非线性程度的量度，b值增大，曲线凹度越大，意味着高流量时退水加快而低流量时退水趋缓。

以均方根误差RMSE（Root Mean Square Error）作为流域退水斜率曲线的拟合准则：

式中：（-dQ/dt）i为观测的流域退水速率；f（Qi）为基于观测的流量Qi和拟合的流域退水斜率曲线计算的流域退水速率；n为流域退水事件的数据序列长度。基于公式（4）和（5）逐一对退水事件的流域退水斜率曲线进行拟合，分为不对消退系数a和b设置限制（拟合方案1）和对b设置限制（拟合方案2）两种方案。

2 结果与分析

2.1 退水事件提取

采用2019—2021 年逐日流量时间序列筛选退水事件，考虑到本研究的流域集水面积较小（0.86 km2），由降水形成的地表径流可在1～2 d 内快速流出，为排除地表径流对退水过程分析的影响，退水事件以流量拐点（洪峰）后第3 d 作为退水过程的起始点（Brutsaert et al., 1998），同时退水过程中满足流量Q及退水速率-dQ/dt连续递减的要求，即排除当日流量大于前日流量、当日退水速率快于前日退水速率的事件。根据以上原则共筛选出25 个退水事件（图2 和表2），其中2019、2020 和 2021 年分别筛选出 10、6 和 9 个退水事件，这些退水事件中持续日数最短仅4 d，最长达15 d，初始流量Q0大小介于0.002～0.092 m3/s。

图2 2019—2021年退水事件流量过程线Fig 2 Recession events selected during 2019-2021

2.2 退水特征变化

通过在log-log双对数坐标空间中绘制-dQ/dt-Q关系（图3）发现，-dQ/dt与Q基本上呈正比例关系，即较快的退水速率-dQ/dt对应于较大的流量Q，但数据点整体尤其在中至尾部（-dQ/dt较小时）显得较为散乱。通过不同颜色和形状区分不同年份的数据点，可以发现2021年数据点与2019和2020年的数据点相比明显向左偏移，由此分别拟合的流域退水斜率曲线斜率b基本保持不变，但截距a将相应增大，说明相同流量情况下，2021年的退水速率比2019和2020年更快。

图3 退水速率-流量关系（-dQ/dt-Q）Fig.3 Scatter points of the estimated-dQ/dt-Q

拟合方案1不对消退系数a和b设置限制（拟合结果见表2 和图4），可以发现拟合的斜率b值大小介于0.6～3.0，均值为2.0，且大部分集中分布在1.5～2.5区间；在此基础上，为便于比较不同退水事件的流域退水斜率曲线，拟合方案2假定各流域退水斜率曲线具有相同的斜率b（本研究假定b=2.0），可以发现拟合的截距a值大小介于0.6～36.3，均值为14.4（拟合结果见表2）。

图4 斜率b的拟合值及概率密度分布Fig.4 Fitted values and their probability density distribution of slope b

表2 退水事件特征及参数Table 2 Recession events with their characteristics and estimated parameters

通过对比可以发现（图5-a、c），2019年的截距a值在8月前仍然保持在6以下水平，8月后才开始出现增加趋势，最大值为14.2；相比而言，2021年的截距a值增加时间提前至6月，且增幅加大，最大值达到36.3。相应地，2019和2021年的潜在蒸散发PET 分别为2 046 和2 135 mm，差别不大；而2019和2021年的降水量分别为2 592和1 223 mm，2021年降水量不足2019年降水量的一半，可见-dQ/dt与Q的关系即流域退水斜率曲线与来水条件关系密切。

图5 2019—2021年截距a值及逐日累积降水/潜在蒸散发变化Fig.5 Temporal variation of intercept a and daily accumulated precipitation/potential evapotranspiration during 2019-2021

2.3 退水影响因素分析

通过比较发现（图6），土壤水和地下水对降水的响应总体较好，地下水埋深较土壤水分含量变化略有延迟，说明土壤入渗能力较强，降水能迅速补充土壤水和地下水。7 月、8 月中旬至9 月中旬2个时间段连续无雨日数增多，同时潜在蒸散发PET 也相对较大，在此作用下土壤水分显著下降，与之相比由于地下水埋深较深（5.0 m 以上），受蒸散发影响较小，因此地下水下降幅度相对缓和。

图6 2021年逐日降水和潜在蒸散发（a）土壤水分含量和地下水埋深（b）变化Fig.6 Temporal variation of daily precipitation and potential evapotranspiration(a)soil moisture volumetric content and groundwater depth(b)in 2021

为识别退水过程的主要影响因素，将拟合方案2下各退水事件拟合的截距a值分别与退水事件的同期潜在蒸散发PET、实际蒸散发ET、地下水埋深及土壤水分含量作相关分析。从图7可以发现，与截距a值相关程度最高的为地下水埋深（r=0.93*），为正相关；其次分别为土壤水分含量（r=-0.76*）、实际蒸散发ET（r=-0.61）和潜在蒸散发PET（r=-0.32），均呈负相关关系。其中，截距a值与地下水埋深及土壤水分含量的相关系数r通过了99%置信水平下的显著性检验。可见，地下水埋深在流域退水过程中起最主要的控制作用，地下水埋深越深，退水速率越快。

图7 截距a值与退水事件同期潜在蒸散发（a）实际蒸散发（b）地下水埋深（c）土壤水分含量（d）关系散点图（*99%置信水平下显著相关）Fig.7 Scatterplots of intercept a and cocurrent potential evapotranspiration(a)evapotranspiration(b)groundwater depth(c)soil moisture volumetric content(d)(*Significant at p=0.01)

3 讨论

有林地包气带土层比无林地疏松，土壤入渗能力较强，降水和地表径流更容易补充土壤水和地下水，从而增加流域的蓄水量，削减洪峰流量和洪水总量，已有研究对森林“蓄水削峰”作用的认识基本一致（Bruijnzeel, 2004; Ellison et al., 2017），在本研究中也得到了实验流域水文气象长期定位观测数据的进一步证实。然而，枯水期由于降水减少，河川径流需要由流域蓄水（主要为地下水）补给时，对于森林能否保持径流均匀和稳定，即所谓的“放水补枯”作用，学者们的认识存在较大分歧（陈军锋等，2001；Ellison et al., 2017; Zhao et al.,2022）。Li 等（2021）基于气候情景的模拟研究发现，近几十年来中国东江流域上游的大规模植树造林对枯水期径流有较为明显的促进作用。而Zhao等（2022）则发现位于东南亚季风区的泰国湄南河流域在经历森林砍伐后，枯水期径流呈现增加趋势。森林覆盖率的增加或减少对流域枯水期径流乃至总径流造成怎样的影响在目前仍未有定论。对于一个流域而言，其所处的自然气候带、地质地貌条件及受人类活动的干扰等都是影响森林与水的关系的重要因素。不同的自然气候带有着不同的降水、入渗、蒸发和产汇流规律。研究区域和尺度的差异性，以及森林生态系统本身的复杂性都会对森林水文效应的评价造成较大的影响，因此，森林水文效应的区域性对比和跨尺度综合研究是未来亟需加强的方向（白晋华等，2008；Ellison et al.,2017）。

对于本研究所依托的实验流域，流域出口处地下水位长期低于河床位置，可以判断该河段为渗失河段，河道水长期向地下水渗失（地下水补给区）。与以往研究认为影响退水过程的主要因素是蒸散发和土壤水分等结论（Cadol et al., 2012; Shaw et al.,2013;Tashie et al.,2020）不同，本研究表明地下水埋深是影响山地森林流域出口退水过程的最主要因素，由于地下水位低于河床位置，流域出口河段处于地下水补给区，因此地下水埋深越深，河道水向地下水渗失越快，从而导致流域出口退水速率加快。对于山地森林流域而言，河岸带森林植被对水陆生态系统间的物质、能量、信息和生物流动具有明显的边缘效应（白晋华等，2008），由此决定了地表水和地下水之间必然存在紧密而复杂的交换作用，这种作用对山地森林流域的枯水期径流的影响目前仍然知之甚少。在森林的生态水文效应下，山地流域河岸带的水文地质条件、河流-含水层之间的水力梯度、包气带厚度及土壤特性等是如何控制河流与含水层系统的水热交互过程（王平，2018），以及由此对河流水文情势的具体影响机制亟需进一步揭示，从而为山地森林流域的枯水期径流模拟及预测、水源涵养功能评估、河流生态系统保护与修复等提供科学依据。研究区属于山地森林流域，集水面积虽然仅有0.86 km2，但海拔跨度比较大（1 115～1 775 m），剧烈的地形起伏决定了下垫面具有显著的非均质性和各向异性，未来需要根据流域的水文地质条件和土壤特性开展更大范围的长期定位观测，结合温度示踪和水化学方法（马瑞等，2013；朱金峰等，2017），量化分析不同河段地表水与地下水的相互作用，从水陆界面的物质通量角度进一步揭示山地森林流域退水过程的潜在机理。

4 结论

采用广东南岭森林生态系统国家野外科学观测研究站实验流域2019—2021年逐日平均流量、土壤水分、地下水埋深、潜在蒸散发、实际蒸散发等实时监测数据，基于由Dupuit-Boussinesq方程推导的流域退水斜率曲线，对南岭山地森林流域的退水特征进行分析和比较，识别影响流域退水过程的主要因素，得到以下结论：

1）退水速率与流量的关系即流域退水斜率曲线在不同来水条件下变化明显。流域退水斜率曲线在斜率b保持不变的情况下，截距a在枯水年2021年较大，在丰水年2019 年和平水年2020 年较小，表现出流域退水斜率曲线受丰枯变化影响明显，相同流量情况下退水速率在枯水年要比在丰水年更快。

2）流域土壤入渗能力较强，降水能迅速补给土壤水和地下水。0～100 cm土壤水分体积含量和地下水埋深对降水过程有较好响应，相比而言，地下水较土壤水变化略有延迟，说明土壤入渗能力较强，降水能迅速入渗补给土壤水和地下水；相对土壤水而言，地下水由于埋深较深（5.0 m以上），受蒸散发影响较小。

3）地下水埋深是影响流域退水斜率曲线变化的最主要因素。退水事件拟合的截距a值与事件同期潜在蒸散发PET、实际蒸散发ET、土壤水分和地下水埋深的相关分析结果表明，地下水埋深是影响流域退水斜率曲线变化的最主要因素，两者呈显著正相关关系（r=0.95），说明地下水埋深越深，流域退水速率越快。