郑明国

（1. 水利部黄土高原水土流失过程与控制重点实验室，郑州 450003；2. 广东省生态环境技术研究所广东省农业环境综合治理重点实验室/广东省面源污染防治技术工程中心/华南土壤污染控制与修复国家地方联合工程研究中心，广州 510650；3. 中国科学院地理科学与资源研究所陆地水循环及地表过程实验室，北京 100101）

0 引 言

河流水沙变化会带来诸多的环境问题。从全球角度言，河流水沙及其挟带的养分物质的变化会导致生物地球化学循环过程的相应变化；从区域角度言，河流水沙变化不仅关系到河流本身的健康，而且在水利规划、防洪减灾、水资源利用保护、水土保持及生态环境建设等方面也是须考虑的重要因素。作为全球变化研究的一个部分，气候变化和人类活动影响下的流域水沙过程研究已成为一个全世界广泛关注的重要科学问题[1]。

作为典型的多沙河流，黄河的泥沙问题为中国历代政府所重视。自20 世纪80 年代中期以来，黄河水沙情势发生了巨大变化[2]，引起了社会各界的广泛关注。黄河水沙变化研究已成为中国水科学领域的重大科学问题之一，有众多的科研项目对此问题展开了研究(详见文献[3])。这些项目研究中，基本都设置有水沙变化贡献分割相关研究内容，即定量评估气候因素和各种人类活动对观测的水沙变化的贡献，藉此确定黄河水沙变化的主导因素，从而为治黄决策服务。然而，这些项目的贡献分割结果并不一致，如冉大川等[4]和刘万铨[5]比较了2000 年以前主要的水土保持减水减沙效益计算成果，发现各成果之间差异较大，难以作为宏观决策的依据[5]。

针对黄河水沙变化贡献分割研究，已提出有大量方法，但仍以水文法和水保法为主[3]。水文法首先建立基准期（即非治理期）的降雨和产沙关系，然后将治理期的降雨因子代入，得到仅降雨发生变化时的产沙量，其与基准期产沙的差值即为降雨对产沙变化的贡献，降雨以外的其他贡献全部归为人类活动结果。该方法的计算精度主要取决于所建立的年降雨产沙回归关系的精度。黄河中游年降雨产沙关系通常表现出一定程度的离散，利用该关系式计算年产沙会有一定的误差，但回归方程计算结果为随机变量的条件均值，因此利用该关系式计算多年均值时误差会小很多。这样，如忽略基准期和治理期两个时段的雨强变化，水文法是一种可靠的产沙变化贡献分割方法[6]，刘万铨[5]甚至建议以水文法为基准对水保法结果进行检验和校正。但水文法只能区分人类活动和气候变化的贡献，无法进一步区分不同人类活动，如不同水土保持治理措施的贡献。

与水文法不同，水保法可以区分不同水土保持治理措施对流域产沙变化的贡献。该方法首先根据坡面小区观测结果，确定各水土保持措施的减蚀率，然后统计各措施治理面积，在此基础上得到各项水土保持措施的减沙贡献。水保法存在诸如小区代表性（包括空间、时间以及措施强度代表性等）等问题[6]，但更为严重的是，水保法完全忽略了减沙效益的空间尺度变异，理论基础存有严重缺陷。世界上许多地区的研究表明[7-8]，尽管坡面治理后坡面侵蚀强度变小，但水流进入沟道后会继续冲沙，可能会导致沟道侵蚀增强，使得流域出口产沙量变化小于坡面甚至保持不变。如对巴西一个1.19km2小流域的研究表明，水土保持措施实施后坡面对流域出口的产沙贡献降低，但同时沟道侵蚀的贡献却变大，导致流域出口的洪水平均含沙量治理前后保持不变[8]，同样的现象在黄土高原也有报道[9-10]。以坡面减蚀率直接代替流域出口的减沙率，显然会导致对水土保持措施减沙贡献的高估。

水沙变化贡献分割的精度取决于数据的质量以及方法的合理性两个方面[5,11]。目前，随着遥感等新技术的进步以及投入经费的增大，数据质量问题已有了很大改善，但仍缺乏有效的贡献分割方法[6]，使得难以准确评估黄土高原不同治理措施的减沙效益[6,12-13]。黄土高原水土保持治理已取得了巨大成效，但减沙贡献分割方法的滞后使得治理效果的归因目前存有较大争议，这严重影响了治黄策略的制定。

黄土区的水保治理措施可分为坡面（主要包括植被和梯田措施）和沟道措施（主要指淤地坝）两大类。基于对比试验流域共计15a 的观测资料，郑明国等[10]认为，在流域尺度上，植被、梯田等坡面措施仅通过减少产流来减沙，并不改变流域出口的水沙关系。基于此，本文提出了一种新的黄土区流域减沙贡献分割方法。首先分析了在水沙关系框架下，不同水保措施减沙机制的异同；在此基础上，给出了不同水保措施及降雨的减沙计算公式；最后，通过与水文法计算结果以及实测数据的比较，对该方法进行了验证。

1 研究区概况、数据来源及研究方法

1.1 研究区概况及数据来源

黄河的泥沙主要来自中游的黄土高原。黄土高原为典型温带半干旱气候，降雨量集中在400～600 mm 之间，大部分地区黄土厚度超过100 m，土质疏松，抗蚀性差，在大规模水土保持治理前，许多地区的土壤侵蚀强度超过10 000 t/km2·a。长期的高强度侵蚀使得本地区沟壑密布，地形起伏大，重力侵蚀如崩塌等频繁发生[14]。

陕北是黄土高原淤地坝建设最为集中的区域。本文选择陕北3 条典型河流—延河、清涧河、大理河为研究区（图1），用于验证淤地坝减沙效益计算结果。1970年以前，当地农业活动强烈，水土保持治理措施较少。1970 年以后开展大规模水土保持治理，治理措施主要包括植树、种草、修建梯田和淤地坝等，2000 年以前流域减沙以淤地坝为主，各流域历年治理面积见表1。根据所搜集到的淤地坝拦沙数据，具体验证工作在甘谷驿、延川、绥德、李家河、曹坪5 个水文站点展开（图1 和表2）。

为验证坡面措施减沙效益计算结果，本文也选择山西离石站和甘肃天水站的两对水土保持对比沟为研究区。山西离石站位于黄土丘陵沟壑区第一副区，1956—1970 年期间曾对离石县王家沟流域内两个毗邻小流域—羊道沟和插财主沟进行过对比观测。羊道沟流域面积0.21 km2，完全未经治理，植被稀疏，农地面积占总面积的59.9%。插财主沟流域面积0.19 km2，自1956 年开始治理，治理措施有梯田、地梗、草田轮作、造林以及封禁等，以林草为主，无淤地坝措施，观测后期治理面积为78.3%。由于严禁人为破坏，植被结构已和天然林区近似，林下有枯枝落叶层和地衣苔藓。羊道沟观测期间年均产沙模数为20 811 t/km2，插财主沟为8 503 t/km2，流域治理的减沙率为59.1%。

图1 大理河、清涧河、延河流域监测站点位置 Fig.1 Locations of the monitoring stations in the Dalihe, Qingjianhe, and Yanhe river basins

表1 延河、清涧河及大理河流域不同水土保持措施治理面积 Table 1 Cumulative areas of various soil and water conservation measures in the basins of Yanhe, Qingjianhe and Dalihe River km2

甘肃天水站位于黄土丘陵沟壑区第三副区，2 个对比观测流域—桥子西沟和桥子东沟，空间上相邻，自然特性也相似。桥子西沟流域面积1.09 km2，为非治理流域。桥子东沟为治理流域，流域面积1.36 km2，于1987年开始治理，2006 年以前以造林为主，随后开始修建有淤地坝[16]。在2006 年，桥子西沟坡耕地面积占比36.4%，林地为49.3%，桥子东沟坡耕地仅占4.5%，林地面积比例高达88.8%[17]。1987—2010 年期间，桥子西沟年均产沙模数为5 756 t/km2，桥子东沟为2 005 t/km2[17]，流域治理的减沙率65.2%。

本文研究所使用的降雨、径流和泥沙数据，桥子东沟和桥子西沟为天水站观测资料；羊道沟和插财主沟为山西省水土保持科学研究所观测，摘录自《山西省水土保持科学研究所径流测验资料（1955—1981）》；曹坪站1970 年以前摘录自《黄河流域子洲径流试验站水文实验资料（1959—1969）》，其他来自《中华人民共和国水文年鉴》。桥子东沟和桥子西沟数据年限为1987—2006年，羊道沟和插财主沟为1956—1970 年；延河、清涧河、大理河数据为1950s—1989 年期间观测，各站点观测年限均超过30 年。为方便，以下称1970 年以前时段为P0 时段，1970—1979 年为P1 时段，1980—1989 年为P2 时段，1970—1989 年为P3 时段。本文中“产沙”一词指观测时段内通过河流某一断面的悬移质总量，产沙模数指单位流域面积的产沙量。

为保证数据精度，本文搜集的淤地坝拦沙量数据均为逐坝测量结果，列于表2 中Δ 列。Δ 中也包括水库淤积量，但所占比例一般较小，故下文中称其为淤地坝拦沙量。Δ 计算公式为：Δ=aλV，V 为淤地坝和水库内淤积体体积，λ 为淤泥干容重，取1.35 t/m3。由于Δ 包括了推移质及淤地坝修建前的填土垫方量，因此需按一定系数折减，a 为折减系数。对淤地坝，a 取0.9；对水库，a 取0.95。V 数据的施测时间为表2 中对应时段的终止年份。在1989 年，陕西省水土保持局对陕北淤地坝展开了一次全面普查。延河、清涧河、大理河和李家河的V 数据均为该次普查期间测量结果，数据摘录自惠养瑜等[15]和张胜利等[18]。岔巴沟的V 数据摘录自徐建华等[19]，为黄河水利委员会吴堡水文总站和子洲县水利局测量结果。

表2 中ΔSYCD为根据Δ 计算的对应时段内坝库年均拦沙量。对延河、清涧河、大理河流域，惠养瑜等[15]根据历年淤地坝控制面积和降雨估算了各时段的淤积比例，然后将Δ 分配到P0、P1、P2 各时段，本研究按相同比例进行Δ 在各时段的分配；李家河站的分配比例直接采用绥德站的结果；曹坪站P0 时段的淤积量采用徐建华等[19]的计算结果。下文将以ΔSYCD为标准对淤地坝拦沙量计算结果进行验证，并称之为“实测拦沙量”。

表2 研究区各时段降雨、径流、泥沙数据及坡面措施减水和淤地坝拦沙量 Table 2 Rainfall, runoff , sediment yield, reduced runoff by slope measures, and trapped sediment by check dams for the selected catchments

表2 中坡面措施的年均减水量，即ΔRSM全部为地表径流（surface runoff）组分，无基流组分。惠养瑜等[15]根据具有长期观测记录径流小区资料，确定了不同质量等级的梯田、林地和草地的减水效益，并按丰、平、枯水年进行了校正，最终确定了延河、清涧河、大理河的梯田、林地和草地不同时期的减水量，对该成果数据进行汇总后即得到ΔRSM；李家河流域按大理河的数据估算；岔巴沟数据为徐建华等[19]计算结果。

下文计算过程中假定P0 时段无水保措施，故利用淤地坝拦沙量（即表2 中ΔSYCD）和水保措施减水量（表2中ΔRSM）对P0 时段的SSY、Q 数据进行了还原。不过，由于P0 时段的水保措施很少，还原后SSY 和Q 变化不大，还原与否对最后的计算结果影响很小。

1.2 相关研究方法

1.2.1 基流分割

本文中地表径流数据利用总径流减去基流获取，基流分割采用广泛应用的Lyne and Hollick 滤波方法[20]。文献[21]推荐使用该方法时，滤波参数a 取值为0.98，并连续进行三轮滤波处理（一次前向滤波，一次后向滤波，最后再进行一次前向滤波）。本文使用该方法时，a 取值为0.925，且仅对日径流时间序列数据进行了一次前向滤波处理。黄河中游8 条支流的数据验证结果表明，按此设置的计算结果和Zheng[22]提出的水沙关系法计算结果更为接近。Zheng[22]认为，径流输沙关系中的截距项即为研究时段内的平均基流量。滤波法把所有的低频组分全部归类为基流，因此会系统地高估基流。水沙关系法的计算结果比滤波法更为接近真值，但当径流输沙关系不好时，如枯水年或少沙年份，该方法难以应用。

1.2.2 减水贡献估算

本文计算过程中需要定量评估降雨以及各水保措施对地表径流的影响。地表径流主要产生于汛期，故首先建立各研究区在P0 时段的年汛期降雨（PFa，mm）和年地表径流（Ra，mm）的回归关系，设为Ra=a0PFa+b0。分析表明，与将其他月份（如6-9 月或5-9 月）作为汛期比，采用5-8 月为汛期建立的Ra和PFa回归方程的R2普遍较高，故本文采用5-8 月降雨量为PFa值。设P0 时段的多年平均汛期雨量为PF0（mm/a），流域治理后某一时段Pi 的多年平均汛期雨量为PFi（mm/a），则与P0 时段比，Pi 时段由于降雨变化导致的地表径流变化量ΔRP（mm /a）可按下式计算

各坡面措施的减水量直接参考相关文献结果，详见表 2 中ΔRSM列及1.1 部分说明。在2000 年以前，黄土区的植被改善和梯田修建面积很有限[11]，故也可假定ΔRSM=0，下文也针对此情形进行了相关计算。设R0和Ri分别为P0和Pi 时段的年均地表径流（mm/a），ΔR= Ri-R0。排除ΔRSM和ΔRP后，ΔR 中剩余部分可认为是Pi 时段内为淤地坝的年均减水量ΔRCD（mm/a）：

1.2.3 利用水文法计算人类活动减沙贡献

首先建立P0 时段的年汛期降雨PFa和年产沙模数SSYa（t/km2）的回归关系，设为SSYa=d0PFa+e0。按水文法，降雨变化对产沙变化的贡献为d0（PFi-PF0）。如降雨外的因素全部归于人类活动的结果，则人类活动的减沙贡献ΔSSYH（t/km2·a）可估计为：

式中ΔSSY = SSYi-SSY0，SSY0和SSYi分别为P0 和Pi 时段的年均产沙模数（t/km2·a）。

表3 中给出了研究区各流域P0 时段的PFa-Ra和PFa-SSYa回归关系。其中PFa为面数据，雨量站分布见图 1。在年尺度，各研究流域内不同站点雨量差别不大，因此直接根据算术平均法计算面雨量值。图2 给出了延河各时段的PFa-Ra和PFa-SSYa关系。

表3 研究区各流域P0 时段汛期降雨PFa 与年地表径流Ra、年产沙模数SSYa 的回归关系式 Table 3 Regressions between Ra and PFa, and between annual specific sediment yield (SSYa ) and PFa for the reference period

图2 延河流域甘谷驿站不同时段年汛期降雨量(PFa)与年地表径流量(Ra)、年产沙模数(SSYa)的关系 Fig.2 Relationships between flood-season precipitation (PFa) and annual surface runoff (Ra) and annual specific sediment yield (SSYa) for different periods in the Yanhe basin

2 不同水保措施及降水的减沙机制

流域产沙主要发生在地表径流时段，基流产沙很少。如忽略基流产沙，则任一流域任一时段的水沙关系总可以表示为：SSY=CR，其中C 和R 分别为该时段的地表径流平均含沙量（kg/m3）和地表径流量（mm/a）。在该水沙关系框架下，任一因素对SSY 的影响机制可概括为3 类：1）仅通过影响R 来影响SSY；2）仅通过影响C来影响SSY；3）通过同时影响C 和R 来影响SSY。黄土高原的水保措施主要包括植被、梯田和淤地坝3 类，其中植被和梯田为坡面措施，淤地坝为沟道措施。在以上水沙关系框架下，不同水保措施的减沙机制存有显著差异。

2.1 植被、梯田等坡面措施

对植被等坡面措施减沙机制的多尺度分析表明[9-10]，在坡面尺度，植被等既减少R 值，也降低C 值，导致减沙率显著高于减水率；但在流域尺度，植被等坡面措施减沙机制与坡面显著不同。其原因为是黄土区沟谷重力侵蚀极其发育，而重力侵蚀基本上与植被无关[23]。在黄土区，即使乔木，根系也主要分布在地表以下80 cm 深度范围内[24-27]，而现代侵蚀沟（指切沟和冲沟）的深度一般为几米到十几米甚至几十米，远超植被根系的作用范围。由于植被等不能有效地抑制沟谷重力侵蚀，坡面水下沟后会继续冲刷，这导致流域出口的C 值和水沙关系并没有发生变化，洪水事件（文献[9-10,17]）和年尺度（图3）均如此。因此，在流域尺度上，植被等坡面水保措施仅通过减R 来减沙，使得减沙率和减水率非常接近。目前文献中支持该结论的观测证据主要有：1）插财主沟和羊道沟1956—1970 年期间的对比观测[9,10]；2）桥子东沟和桥子西沟1987—2006 年期间的对比观测[17]；3）甘肃天水站3 个水土保持治理流域（即桥子东沟、吕二沟和罗玉沟）1980s—2010 年期间水沙关系随流域治理的演变[16-17,28-29]。这3 个流域的长时间序列观测数据均表明，如无淤地坝措施，植树、修建梯田等坡面变化均未导致水沙关系的变化。如桥子东沟1987—2006 年期间，流域坡耕地面积占总面积比例从54.2%下降到8.3%，植被和梯田面积比例从40.1%增长到90.5%，下垫面的巨大变化造成了流域降雨-产流和降雨-产沙关系的显著改变，但并没有改变流域的径流-产沙关系，导致单位径流输沙量即C 值变化不大[16]。只有当2006 年后流域内建设了大量淤地坝后，水沙关系才有明显改变，C 值才显著降低[16]。Zhang 等[30]最近的计算也表明，在黄河中游，植被主要通过减少产流来减沙，而径流含沙量的变化主要与淤地坝建设有关。

在黄土区大规模治理前，由于泥沙来源充沛，绝大部分水流都可以达到冲刷限制下的极限含沙量（水流此时不能冲刷，但输沙能力未必达到饱和）[31-33]。其他发生在坡面的人类活动，如修路、耕作等，和植被一样，也不能改变黄土区泥沙来源充沛的特点，因此也应仅仅通过改变地表径流量R 来影响流域产沙。这一特性在大流域应更为突出，因为小流域受坡面过程影响显著，随流域面积增大，坡面过程对水沙输移影响会逐渐减小[9]。在黄土高原，也有许多中大尺度流域，如岔巴沟流域[34]、大理河流域[35,36]、窟野河流域[9]、河龙区间[37]、渭河流域[9]等，报道有流域治理未改变径流/输沙关系或流域治理的减水率和减沙率接近的情况。

美国亚利桑那州Walnut Gulch 试验流域也有和王家沟、天水试验站类似现象[31]。该流域坡面小区的观测表明，灌木覆盖小区的平均含沙量为草地覆盖小区的2 倍以上[38]。但在流域尺度，草地覆盖流域比灌木覆盖流域产沙模数偏小37.5%，产流模数偏小33.2%，两者非常接近，这导致两者的径流平均含沙量也非常接近（1.34% vs. 1.36%）[31]。该流域泥沙来源也非常充沛[38]，因此认为，坡面变化仅通过改变R 来影响流域产沙，该结论可能是泥沙来源充沛地区的一个普遍规律。

图3 羊道沟（非治理流域）和插财主沟（治理流域）年径流深(Ra)与年产沙模数(SSYa)关系的比较 Fig.3 Comparison of the Ra-SSYa relationships between Yangdaogou (non-managed) and the Chacaizhugou (managed) watersheds

2.2 降雨

图4 表明，在羊道沟小流域，降雨对洪水平均含沙量的影响很小。尽管低的含沙量基本发生在小雨量或小雨强情况下，但小雨量或小雨强时也可能出现很高的含沙量。对主要的产沙事件（定义为洪水平均含沙量高于500 kg/m3的事件，其产沙量占总量的86.9%），洪水平均含沙量不仅与雨量和雨强完全不相关（岔巴沟等流域也报导有类似现象[39-40]），而且变异较小，变差系数（CV）为16.9%。其他许多黄土流域也有类似的特点，这导致洪水产沙模数和径流深关系可以用比例函数很好地拟合[41-42]。因此认为，降雨对流域产沙的影响机制与植被相同，也仅通过影响R 来影响产沙。

2.3 淤地坝

淤地坝的减沙机制取决于其是否有良好的排水、排洪设施。无排水设施的淤地坝，如闷葫芦坝，如未发生漫溢和垮坝，会全部拦蓄上游来水来沙。在黄土区大规模治理前，C 值不仅历次洪水之间变化不大，空间上变异也很小[31-32]，如在无定河黄土区，C 值的空间变差系数仅为7%左右[32]。因此可认为，当流域某一范围的来水来沙全部为无排水设施的淤地坝拦截时，流域出口的含沙量并不会改变，无排水设施的淤地坝仅通过减R 来减沙。有良好排水设施（包括溢洪道和放水建筑物）的淤地坝，坝体的壅水和拦截效应会使得大量泥沙落淤，导致排放到下游的洪水含沙量大大变小。此类淤地坝如蓄水，则通过减R 和减C 两种机制减沙，如不蓄水，则仅通过减C来减沙。这几类淤地坝在黄土区流域一般均存在，因此可认为淤地坝既通过减R 来减沙，也通过减C 来减沙。

图4 羊道沟流域洪水事件平均含沙量(Ce) 随降雨量(Pe)和最大30min 雨强(I30)的变化 Fig.4 Mean sediment concentrations of individual runoff events (Ce) varying with rainfall amount (Pe) and maximum 30-min intensity (I30) in the Yangdaogaou watershed

在黄土高原流域治理早期，大部分淤地坝排水设施配置均较差，淤地坝的减沙机制会以减R 为主；后期有良好排水设施的淤地坝数量逐渐变多，因此通过减C 而实现的减沙效应会逐渐增强。图2a 比较了延河流域P0、P1、P2 时段的降雨径流关系，结果表明，相同降雨情况下P1 和P2 时段的R 值显著小于P0 时段，R 值的降低显然与淤地坝的拦蓄作用密切相关；但3 个时段的C 值分别为460、469、341 kg/m3（表2），这表明只有到了P2 时段，淤地坝才通过降低C 值这种方式来减沙，清涧河流域有相同趋势。在大理河流域，尽管C 值在P1时段就已开始降低，但降幅有限，仅在P2 时段才开始显著降低（表2）。1980s 期间新建淤地坝数量较少[43]，出现这种情况的原因很可能是由于黄河中游1977—1978 年大洪水期间，发生了严重的垮坝现象，为防止类似现象再发生，其后对许多淤地坝补建了溢洪道等排水设施[44]。

图5 表明，在P0 时段，和事件尺度（图4）类似，主要输沙年份的地表径流平均含沙量年际之间变异也较小，CV 为15%左右，这些主要输沙年份贡献了约90%的总产沙。假定某一流域在P0 时段的地表径流平均含沙量为C0，根据以上分析可知，除非该流域后期修建有较多的排水良好的淤地坝，则该流域地表径流平均含沙量可始终估计为C0。

图5 P0 时段主要输沙年份的地表径流平均含沙量(Ca) Fig.5 Mean sediment concentrations of surface runoff (Ca) among the major sediment-yield years in P0 period

3 方法的提出

给定某一黄土区流域，假定多年平均降雨、地表径流量、地表径流平均含沙量和产沙模数在无治理期P0 时段的分别为P0、R0、C0、SSY0，流域治理后某一时段Pi分别为Pi、Ri、Ci、SSYi，设ΔP=Pi-P0，ΔR=Ri-R0，ΔC = Ci-C0，ΔSSY=SSYi-SSY0。忽略水土保持治理以外的其他人类活动，并假定流域治理措施仅包括坡面措施和淤地坝，则可认为：ΔR=ΔRP+ΔRSM+ΔRCD，ΔSSY=ΔSSYP+ ΔSSYSM+ΔSSYCD，其中ΔRP、ΔRSM、ΔRCD分别为降雨、坡面、淤地坝对ΔR 的贡献，ΔSSYP、ΔSSYSM和ΔSSYCD分别为降雨、坡面措施和淤地坝对ΔSSY 的贡献。我们的目标为：根据P0 和Pi 时段的降雨、径流和产沙观测资料，确定ΔSSYP、ΔSSYSM和ΔSSYCD。

可用R、C、SSY 作为描述流域水沙状态的3 个变量。在P0 时段，流域状态为（R0，C0，SSY0），由于降雨变化和后期流域治理的实施，在Pi 时段流域状态变为（Ri，Ci，SSYi）。按第三部分分析，降雨、坡面治理和淤地坝均可影响R，但仅淤地坝会影响C 值。流域状态从（R0，C0，SSY0）变为（Ri，Ci，SSYi）可分解为以下4 个步骤：

1）降雨变化ΔP，导致R 变化了ΔRP，由于降雨不改变C 值，故C 值仍等于C0，经过此步骤，SSY 变化了C0ΔRP；

2）坡面措施进一步使R 变化ΔRSM，同上，C 值仍等于C0，经过此步骤，SSY 变化了C0ΔRSM；

3）含沙量为C0的水流到达淤地坝，其中有ΔRCD的径流被淤地坝拦截，经过此步骤，SSY 变化了C0ΔRCD；

4）排出淤地坝到达流域出口的径流量变为Ri，这些径流在汇入淤地坝时含沙量为C0，到达流域出口时为Ci，经过此步骤，SSY 变化值为RiΔC，完成从状态（R0，C0，SSY0）到（Ri、Ci、SSYi）的变化。

以上第一个步骤中SSY 的变化完全由降雨引起，故可认为：

同理，第二个步骤中SSY 的变化可认为是坡面措施的贡献，故有：

第三步骤中SSY 的变化是淤地坝通过减R 导致的贡献，第四个步骤中SSY 的变化是淤地坝通过改变C 值导致的贡献，故：

式（4）—（6）等式右端诸项之和刚好为ΔSSY。ΔSSY的分割公式因此为：

式中ΔSSYH=ΔSSYSM+ΔSSYCD，为人类活动对ΔSSY 的贡献。根据公式（7），ΔSSYH可估算为：

按以上分析，如某一坡面措施m 使得R 变化了ΔRm，则其对ΔSSY 的贡献为C0ΔRm。如欲进一步区分不同坡面措施的贡献，则：

式中ΔSSYm为第m 种坡面措施的贡献，∑ΔSSYm即ΔSSYSM，为各坡面措施贡献之和。

以上公式中，计算结果中的负值表示减沙效应，正值表示增沙效应。式（6）中，如ΔC＞0，这表示淤地坝修建后并没有改变C 值，而非淤地坝具有增沙效应。实际计算时，如出现ΔC＞0，可置ΔC=0，并将SSY 的变化全部归因于R 的变化。

4 方法的验证

本部分将利用以上提出的方法（以下简称水沙关系法），计算各研究流域的减沙效益，并对结果进行验证，表4 给出了延河、清涧河、大理河流域的计算结果。验证方案分两种：1）比较水文法和水沙关系法计算的ΔSSYH是否一致，以此验证公式（8）和（4）；2）通过与观测值比较验证公式（4）—（6）。

分别按水沙关系法和水文法，计算了表2 中各流域各治理时段的ΔSSYH。水文法被认为是一种较可靠的贡献分割方法[5-6]，图6 表明，水文法和水沙关系法计算的结果非常接近。这证明了公式（8）的有效性，也间接证明了公式（4）的有效性。

对某一特定流域，如下垫面无变化，则历年产沙的变化可认为完全是降雨的变化所导致。基于此假设，利用羊道沟的数据继续对公式（4）进行了验证，结果见图 7。ΔSSYP的观测值和计算值非常接近。

图8 利用羊道沟/插财主、桥子东沟/桥子西沟两对对比流域的观测数据对公式（5）进行了验证。研究时段内，治理流域插财主沟和桥子东沟均只有坡面治理措施，而无淤地坝。在此情形下，治理流域和非治理流域产沙模数的差异即为ΔSSYSM的观测值（真实值）。图8 表明，在多年时间尺度，公式（5）计算值对观测值的偏差非常小，两个地点都仅约6%；在年时间尺度，计算精度有所下降，但计算值和观测值的R2仍高于0.8。

表4 水沙关系法计算结果 Table 4 Results of runoff-sediment relationship method

图6 水文法和水沙关系法计算的人类活动减沙量比较 Fig.6 Comparison of the calculated sediment-reduction effects of human activities between the H (hydrology) method and the RSYR (runoff-sediment yield relationship) method

图9 计算了表2 中各流域除P0 外其他时段的淤地坝拦沙量，按忽略（图9a）和考虑（图9b）坡面措施减水量两种方式进行计算，并分别与表2 中给出的实测值进行了比较。当考虑坡面措施减水量时，计算值和实测值非常吻合。然而，即使忽略坡面措施减水量，计算值和实测值也偏差不大，其原因是在研究时段，流域减沙效益以淤地坝为主， 坡面措施贡献有限。根据表4 计算结果，对P1、P2、P3 时段，各流域坡面措施对减沙贡献为8.4%～17.9%，平均为12.4%，而淤地坝为78%～101%，平均为86.9%。图 9a 中两者的偏离有随淤地坝拦沙量变大而变大的趋势，主要原因是淤地坝拦沙量较大的一般对应大的流域，而大流域的坡面措施面积一般也大，其减沙量的绝对值也大，忽略坡面措施的贡献因此在大流域会造成更大的绝对误差。一旦考虑了坡面措施的贡献，这种趋势就完全消失。

图7 利用羊道沟数据验证降雨的减沙效应计算公式（公式（4）） Fig.7 Verifying the calculated effect of rainfall on sediment yield by equation (4) using annual observations of the Yangdaogou (non-managed) watersheds

图8 利用对比观测流域数据验证坡面水保措施的减沙效应计算公式（公式（5）） Fig.8 Verifying the calculated effect of slope measures on sediment yield by equation (5) using the paired-watersheds observations

图9 利用淤地坝拦沙量实测数据验证淤地坝的减沙效应计算公式（公式（6）） Fig.9 Verifying the calculated effects of check dams on sediment yield by equation (6) using the measured values of the trapped sediment by check dams

5 讨 论

本文淤地坝拦沙效益的计算公式，主要在多年时间尺度进行了验证，但也有两次年时间尺度的结果，即曹坪站1977 和1978 年。1977 和1978 年黄河中游发生特大暴雨，估算的曹坪站淤地坝拦沙量与实测值的偏差仅为4.2%和11.1%（表4）。1978 年汛期降雨为487 mm，远超P0 时段多年均值（296 mm），计算结果表明，降雨导致产沙模数增加了23 163 t/km2，淤地坝拦沙导致其减少了22 784 t/km2，两者大致相等，因此曹坪站该年的产沙模数（21 000 t/km2；表2）保持和P0 时段的平均值（21 901 t/km2）接近。但这并不意味着本文方法也能运用于枯水年，因为枯水年含沙量低的原因，基本不是由于淤地坝的拦沙作用，而是因为地表径流流量很小，侵蚀和搬运泥沙的能力较弱的缘故。坡面措施减沙效益计算公式也不适用于枯水年，其原因是枯水年以小洪水或小流量为主，离石和天水站的观测均表明，流域进行坡面治理后，小洪水或小流量情况下含沙量明显变小[10,17]。总之，坡面措施和淤地坝减沙效益计算公式均不适用于枯水年。但由于枯水年对流域产沙贡献很小，故本文方法能适用于多年时间尺度。

如上述，坡面治理可显著降低小洪水或小流量情况下的含沙量，因此植被等坡面措施的减沙机制可更精确地概括为“植被等坡面措施不改变流域主要输沙时段的水沙关系”。当流域治理达到很高程度，以至于地表径流以小洪水或小流量为主时，植被等坡面水保措施的减沙机制也可能发生变化。在黄土高原未治理前，当流量小于某一阈值时，同流量下的含沙量变幅极大，当流量超过该阈值时，含沙量的变幅明显缩小[45]，趋近于冲刷限制下的极限含沙量[31-33]。如坡面治理能将洪峰流量都控制在小于该流量阈值范围内，则洪水的含沙量会有非常显著的降低。但这一阈值非常小（参考文献[41]图2），将主要输沙时段的流量控制在该阈值以下会非常困难。

坡面措施减沙效益的计算公式（公式（5））只在小流域尺度进行了验证。但大的河流，其水沙特性会比小的河流更难以改变，若流域治理未能改变小流域的C 值，则可以推断，流域治理更不可能改变大流域的C 值[9]。因此认为，公式（5）也应该可以很好地运用于大流域的植被减沙效益计算。

2000 年以后，随着退耕还林政策的实施，黄河中游植被覆盖度显著提高，同时，也新修了大量排水良好的淤地坝。本文对淤地坝减沙效益计算公式（公式（6））的验证，未包括2000 年以后的数据，因此需要进一步验证。任何一种治理措施，其减沙效益均通过改变R 或C值实现，因此按公式（6）计算结果应包含淤地坝的全部减沙效益。除通过坝体直接拦沙外，淤地坝尚有其他的减沙途径，如淤积泥沙抬高了侵蚀基准面，从而能稳定坡谷、减少重力侵蚀的发生。退耕政策实施前，沟谷中泥沙淤积厚度可能尚有限，淤地坝的减沙效益主要通过坝体拦沙来实现，故本文方法计算结果和实测拦沙量非常接近。但随着淤积年限的增加，通过稳定坡谷这一途径实现的减沙效益应该会越来越重要，这可能导致公式（6）计算结果显著超过实测拦沙量，两者间的差异可认为是淤地坝通过稳定坡谷而实现的减沙效益。

6 结 论

基于降雨和不同水保措施的减沙机制差异，本文提出了一种新的流域减沙贡献分割方法。验证结果表明，该方法可较好进行降雨、不同坡面措施以及淤地坝等的减沙贡献分割。常规的水保法基于坡面小区观测结果，其计算结果为减蚀率，即水保措施在坡面减少了多少侵蚀。本方法考虑了减沙效益的空间尺度变异，计算结果为减沙率，即水保措施在流域出口减少了多少泥沙。本方法有较好的理论基础，可以有效地估算气候变化和不同流域治理措施或不同人类活动对入黄泥沙的影响，对厘清黄河泥沙锐减的主导驱动因素和未来治黄策略的制定具有重要意义。