湖泊科学   2019, Vol. 31 Issue (4): 1120-1131.  DOI: 10.18307/2019.0407.
0

研究论文

引用本文 [复制中英文]

包为民, 侯露, 沈丹丹, 倪用鑫, 黄土高原大理河流域水沙耦合模型应用研究. 湖泊科学, 2019, 31(4): 1120-1131. DOI: 10.18307/2019.0407.
[复制中文]
BAO Weimin, HOU Lu, SHEN Dandan, NI Yongxin. Application of flow-sedimentation coupled model in Dali River Basin of Loess Plateau. Journal of Lake Sciences, 2019, 31(4): 1120-1131. DOI: 10.18307/2019.0407.
[复制英文]

基金项目

国家重点基础研究发展计划“973”项目(2016YFC0402703)和中央级公益性科研院所基本科研业务费专项(HKY-JBYW-2017-12)联合资助

作者简介

包为民(1956~), 男, 博士, 教授; E-mail:wmbao@163.com

通信作者

侯露, E-mail:583974810@qq.com

文章历史

2018-09-10 收稿
2018-12-07 收修改稿

码上扫一扫

黄土高原大理河流域水沙耦合模型应用研究
包为民1 , 侯露1 , 沈丹丹1 , 倪用鑫2     
(1: 河海大学水文水资源学院, 南京 210098)
(2: 黄河水利委员会黄河水利科学研究院, 郑州 450003)
摘要:水沙模型是定量描述水沙关系及水沙规律的重要工具,现阶段国内外对于水沙模型的研究大都为基于某个典型流域的经验统计模型或基于流域大量基础资料的物理模型,极大限制了其使用范围及模拟精度.本文建立了结构与参数均具有物理意义的流域水沙耦合物理概念模型,其优点是物理概念清楚,模拟精度高,实用性强,易于深入研究泥沙基本规律.该模型将概念性水文模型和泥沙模型耦合,提出水流挟沙能力和土壤抗侵蚀能力概念,用对数曲线近似描述流域土壤抗侵蚀能力的空间变异性,在拜格诺河道水流悬移质泥沙公式基础上建立概念性沟蚀产沙公式,参照水流汇集相似性建立泥沙汇集演算公式.选取黄河中游大理河流域4个不同流域尺度的实际流域对模型进行应用检验,模拟结果表明,该模型的水流泥沙两部分均有很高的模拟精度,可以很好地模拟黄河中游地区不同流域尺度和年际尺度上的洪水过程和泥沙产生及输移过程,充分证明了该模型结构、参数和计算方法的合理性,可用于定量分析流域内各项水保措施的减水减沙效益及水沙关系变化趋势,对该模型的推广应用可做进一步分析研究.
关键词黄土高原    大理河流域    概念模型    耦合模拟    水流    泥沙    
Application of flow-sedimentation coupled model in Dali River Basin of Loess Plateau
BAO Weimin1 , HOU Lu1 , SHEN Dandan1 , NI Yongxin2     
(1: College of Hydrology and Water Resources, Hohai University, Nanjing 210098, P. R. China)
(2: Yellow River Conservancy Commission, Yellow River Institute of Hydraulic Research, Zhengzhou 450003, P. R. China)
Abstract: The flow-sedimentation coupled model is a crucial tool for quantitatively describing the relationship and the law of flow and sedimentation. Studies on flow-sedimentation model at home and abroad are almost empirical statistical models of a typical basin or physics-based model based on considerable description of the watershed at this stage, which greatly limits its usage and simulation accuracy. A conceptual flow-sedimentation coupled simulation model for basins, which combines the conceptual hydrological model with the sediment model, is presented with physical structures and parameters in this paper. Its clear physical concept, high simulation precision and strong practicability make it easy to study the basic law of sediment. We put forward the concept of hydraulic erosion capacity and soil erosion resistance capacity to estimate the sediment yields on the slopes and use logarithmic curve to approximately describe the spatial variability of watershed soil erosion resistance capacity. The calculus formula for sediment yield in gully and sediment concentration is respectively established based on the suspended sediment formula of the Bagnold's River channel and by referring to the similarity of flow. The proposed model is tested in four watersheds of different scales in the Dali River Basin and the results indicate that the model can provide a good simulation of runoff and sediment generation as well as transportation processes at both watershed scale and inter-annual time scale in the middle reaches of the Yellow River, and thus further confirm that the model is reasonable in structure, parameters and calculation methods. The model can be used to quantitatively analyze the benefits of flow and sediment reduction and how their relationship changed due to various water conservation measures implemented in the basin, and can be further analyzed and studied for its popularization and application.
Keywords: Loess Plateau    Dali River Basin    conceptual model    coupled simulation    runoff    sediment    

黄河泥沙问题是治黄的根本问题.近年来,黄土高原地区土地资源利用的变化、水资源供需矛盾的突出以及各种水利工程的建设,加上由降水和温度因素引起的自然气候变化,使得流域内来水来沙量锐减,产汇流、产汇沙过程和水沙关系均产生了巨大变化[1-4].如何定量分析人为与自然因素对流域水沙关系的影响是流域综合管理的前提,长期以来备受国内外专家的高度重视[5-8].但是,由于水沙过程所涉及的物理机理和泥沙在流域上变化的复杂性,外加流域内各项水保措施对泥沙的作用,水沙条件及水沙系列资料受全球气候条件变化、人类活动、流域尺度等影响造成非周期性变化,国内外的水沙模型常为针对于某个典型流域或是资料充足的较小试验流域建立的经验统计模型.通用土壤流失方程(USLE)[9]在国外统计模型中最具影响力,该模型考虑了土壤可蚀性因子、水保措施因子、作物覆盖因子等,但仍为经验性土壤侵蚀预报模型. LISEM、WEPP模型[10-11]等在描述径流、土壤侵蚀产沙量方面基于一定的物理公式如质量守恒和动量守恒方程,由于模型输入需要大量的试验流域基础资料,应用于不同流域尺度或无资料地区的水沙过程模拟精度低,通用性差.综观目前的经验统计模型和基于一定物理过程的土壤侵蚀模型,各模型考虑的因素不尽相同,但对不同流域尺度的泥沙规律没有得到很好的体现,模型模拟精度及适用范围受到极大限制.

水沙模型的发展趋势是能同时进行洪水预报及泥沙冲淤计算,且具有可移植性[12].因此,本文基于流域自然地理特性及流域产汇流、产汇沙过程的主要物理机制建立结构和参数均具有物理意义的流域概念性水沙耦合模型,将概念性水文模型与泥沙模型耦合,模拟水流和泥沙分量.不断优化坡面、沟道产沙计算中的参数结构,在汇沙计算中引入河道动力非线性冲刷系数,突出水流对泥沙的作用,以物理成因机制为基础,研究泥沙的基本规律.经黄土地区不同流域尺度的实际流域和前期大量的试验流域验证[13-14],该模型可应用于不同流域尺度的水沙模拟,模拟精度高,结构简单,各部分均具有一定的物理意义,可以很好地模拟黄土高原地区流域的水沙规律,易于推广应用.

本文采用黄河中游大理河流域的王茂沟、韭园沟两个小流域和岔巴沟、大理河两个大流域的多年实测水文泥沙资料对建立的模型进行应用检验,利用函数曲面参数率定方法[15]率定水流、泥沙两部分参数,分析水流泥沙两部分模拟精度结果的主要影响因素.黄河流域特别是中游多沙粗沙地区水土流失十分严重,在该地区进行概念性水沙耦合模型的研究与应用不仅有利于研究干旱半干旱地区的水沙过程,而且有利于定量分析水土流失情况,定量评价气候、下垫面变化及人类活动对水沙关系的影响,准确掌握流域水沙情况及变化趋势,为流域综合治理规划、水土保持、合理利用水沙资源等提供必要的基础数据[16].水沙模型作为一项重要的非工程措施,可成为研究流域产汇流、产汇沙过程的重要研究手段及合理利用水沙资源、定量评价水土保持效益的有效工具.

1 流域水沙耦合模型 1.1 模型模拟内容与流程

概念性水沙耦合模型,主要分为产流计算、汇流计算、产沙计算和汇沙计算,据水沙在流域上产生、运动各个环节机制的差异及理论、应用方便和模拟效果综合考虑,流域水流概念模型中产流计算部分采用垂向混合产流计算,汇流计算中山坡汇流采用线性水库、沟道汇流采用马斯京根法进行演算;概念性泥沙演算模型[17-18]中流域产沙包括面蚀和沟蚀产沙,坡面侵蚀产沙量利用水流挟沙能力和土壤抗侵蚀能力平衡方程推求,沟蚀产沙量利用基于拜格诺河道水流悬移质泥沙公式概念化的概念性沟蚀产沙公式进行计算,并依据水流汇集规律相似性,将泥沙汇集分为坡面汇沙和沟道汇沙.将这些机制的成因关系和概念模拟方法有机结合起来,组成一个较为完整的流域概念性水沙耦合模型[13],模型结构如图 1.

图 1 水沙物理概念模型结构图 Fig.1 Water and sediment conceptual model diagram
1.2 模型计算原理与方法 1.2.1 垂向混合产流计算

包为民[19]从物理成因下渗方程出发,根据土壤含水量对毛管水压力的影响关系,把忽略地面滞水深条件下的格林-安普特下渗公式改进为:

$ F M=F C(1+K F(W M-W) / W M) $ (1)

实际下渗FA计算公式为:

$ F A=\left\{\begin{array}{l}{F M-F M(1-P / F M M)^{1+B F}, P <F M M} \\ {F M, P \geqslant F M M}\end{array}\right. $ (2)

式中, FMM为流域最大下渗能力(mm/min);BF为流域下渗能力分布曲线指数;KF为渗透系数,反映土壤含水量对下渗的影响;FM为流域平均下渗能力(mm/min);P为降雨量(mm);WM为土壤田间持水量(mm);W为土壤含水量(mm);FC为饱和条件下的下渗率,即稳定下渗率(mm/min).

根据流域特性垂向混合产流将流域下渗能力分布曲线和流域蓄水容量分布曲线垂向组合,结构图如图 2.垂向混合产流计算中地面径流RS取决于降雨强度和初始土壤含水量,地下径流RR取决于土壤缺水量和下渗水量FA,分别为超渗产流和蓄满产流计算模式.计算方程为:

$ R S=\left\{\begin{array}{l}{0, P E \leqslant F A} \\ {P E-F A, P>F A}\end{array}\right. $ (3)
$ R R=\left\{\begin{array}{l}{F A-W M+W+W M\left(1-\frac{F A+W^{\prime}}{W_{\operatorname{mm}}}\right)^{1+B}, F A+W^{\prime} <W_{\operatorname{mm}}} \\ {F A-W M+W, F A+W^{\prime} \geqslant W_{\operatorname{mm}}}\end{array}\right. $ (4)
$ R=R S+R R $ (5)
图 2 垂向混合法结构图 (F为流域上某点的下渗能力(mm/min);αF为下渗小于某个定值F的面积比例;ΔW为土壤缺水量(mm)) Fig.2 The structure of vertical mixing method

式中,PE为净雨量(mm),B为土壤蓄水容量分布曲线指数,W′为流域上某一点的土壤含水量(mm),Wmm为流域最大蓄水容量(mm).

1.2.2 汇流计算

(1) 水源划分.水源划分采用自由水箱结构划分法,该方法考虑了包气带的垂向调蓄作用.按混合产流模型计算出的总径流量R,先进入自由水箱调蓄,再划分水源.自由水箱结构如图 3所示.

图 3 自由水箱结构图 Fig.3 The structure of the free water reservoir

自由水箱划分水源的具体计算式为:

$ F R=R / P E $ (6)
$ R S=\left\{\begin{array}{l}{R+(S+S M) \cdot F R, P E+A u \geqslant S M M} \\ {R+\left\{S-S M+S M\left[1-\frac{(P E+A u)}{S M M}\right]^{E X+1}\right\} \cdot F R, P E+A u <S M M}\end{array}\right. $ (7)
$ R I=K I \cdot S \cdot F R $ (8)
$ R G=K G \cdot S \cdot F R $ (9)

式中,RSRIRG分别为地面径流、壤中流和地下径流(mm),FR为产流面积比或径流系数,S为初始自由水蓄量(mm),KIKG分别为壤中流和地下径流出流系数,SM为流域平均蓄水深(mm),Au为相应平均蓄水深的最大蓄水深(mm),SMM为流域最大蓄水径流深(mm),EX为反映蓄水深流域分布特征的参数.

(2) 山坡汇流.根据水量平衡和蓄泄关系[20],可得山坡汇流演算式为:

$ Q S_{2}=C S \cdot Q S_{1}+(1-C S)\left(R_{1}+R_{2}\right) / 2 \cdot U $ (10)
$ Q I_{2}=C I \cdot Q I_{1}+(1-C I)\left(R_{1}+R_{2}\right) / 2 \cdot U $ (11)
$ Q G_{2}=C G \cdot Q G_{1}+(1-C G)\left(R_{1}+R_{2}\right) / 2 \cdot U $ (12)

式中, QSQIQG分别为地面径流、壤中流、地下径流出流量(m3/s);CSCICG分别为地面径流、壤中流、地下径流线性水库消退系数;R1R2分别为时段初、末的径流量(m);U为单位换算系数,$U=\frac{F S\left(\mathrm{m}^{2}\right)}{3.6 \Delta t(\mathrm{s})}$.

(3) 河道汇流.运用马斯京根河道演算法[21],可得河道汇流演算式为:

$ Q_{2}=C_{0} \cdot I_{2}+C_{1} \cdot I_{1}+C_{2} \cdot Q_{1} $ (13)
$ C_{0}=\frac{0.5 \Delta t-K E \cdot X E}{0.5 \Delta t+K E-K E \cdot X E}, \quad C_{1}=\frac{0.5 \Delta t+K E \cdot X E}{0.5 \Delta t+K E-K E \cdot X E}, \quad C_{2}=\frac{-0.5 \Delta t+K E-K E \cdot X E}{0.5 \Delta t+K E-K E \cdot X E} $ (14)

式中, I1I2分别为时段始、末的河段入流量(m3/s), Q1Q2分别为时段始、末的河段出流量(m3/s), C0C1C2为马斯京根洪水演算方法的演算系数;KE为马斯京根法河段传播时间(min);XE为马斯京根法流量比重系数.

1.2.3 产沙计算

(1) 坡面侵蚀产沙.定义流域坡面水流挟沙能力为当流域坡面上各处均有足够供水流挟带的疏松土壤时,水流能带走的该流域泥沙量.全流域平均的坡面水流挟沙能力SC为:

$ S C=A C M \cdot R \cdot A $ (15)
$ A C M=C M \cdot\left(P / P_{\max }\right) $ (16)

式中, R为全流域坡面平均的径流深(m);ACM为坡面水力侵蚀能力系数(kg/m3);A为全流域坡面面积(m2);P为流域时段平均面降雨量(m);Pmax为流域最大降雨量(m);CM为全流域坡面水流最大可能含沙量(kg/m3),主要取决于坡面水流因素的变化幅度.

然而,对于实际流域,水流除挟带疏松土壤外,还要剥离一些非疏松的土壤颗粒.在剥离过程中,由于流域上各点产生的抗动力不同,所以流域实际的侵蚀量与水流挟沙能力也不同.流域土壤抗侵蚀能力的空间变化分布曲线如图 4所示. 图 4中,REC为土壤抗侵蚀能力(kg);REM为流域平均抗侵蚀能力(kg);αR为抗侵蚀能力小于REC的面积比;REMM为流域最大抗侵蚀能力(kg);α0为抗侵蚀能力为零的面积比.

图 4 土壤抗侵蚀能力分布曲线图 Fig.4 Distribution curve of soil erosion resistance

由此可知,流域坡面实际产沙应是水流挟沙能力与土壤抗侵蚀能力之差的积分[13].则流域坡面实际产沙量SS可表示为:

$ S S=\left\{\begin{array}{l}{S C-R E M\left[1-\left(1-\frac{S C}{R E M M}\right)^{B S+1}\right], S C <R E M M} \\ {S C-R E M, S C \geqslant R E M M}\end{array}\right. $ (17)
$ R E M=\int_{\alpha_{0}}^{1} R E C \mathrm{d} \alpha=\frac{1-\alpha_{0}}{B S+1} $ (18)

式中, BS为抗侵蚀能力分布曲线指数.

(2) 沟蚀产沙.沟道产沙含沙量随沟道水力条件的变化而变化,据拜格诺的河道水流悬移质泥沙公式概念化为:

$ S G=C G M \cdot \frac{u}{u M} \cdot Q $ (19)

式中, SG为沟道侵蚀产沙速率(kg/s), CGM为平均沟道产沙浓度(kg/m3), u为断面平均流速(m/s), uM为断面时空平均流速(m/s), Q为沟道水流流量(m3/s).

(3) 流域产沙.根据上述分析可知,流域次洪产沙量(S)为次洪坡面侵蚀产沙量和沟蚀产沙量之和,即:

$ S=\sum\limits_{i=1}^{N} S S_{i}+\sum\limits_{i=1}^{N} S G_{i} \Delta T $ (20)

式中, N为次洪水时段数;ΔT为时段长(min).

1.2.4 汇沙计算

泥沙汇集概念模拟,忽略水沙间的相互作用,可把水沙分离模拟.流域上任一控制元的泥沙质量平衡可表示为:

$ \frac{\mathrm{d} W S}{\mathrm{d} t}=I-O+q $ (21)

式中, IO分别为输入和输出控制元的泥沙速率(kg/s), WS为控制元内随水流运动的泥沙量(kg), q为控制元内的冲淤速率(kg/s).式(21)在泥沙汇集模拟应用中,Iq为已知量,WSO为未知量,但OWS密切相关.考虑坡面控制元时有关系[18]:

$ W S=K S[X S \cdot I+(1-X S) \cdot O] $ (22)

式中, XS为泥沙比重系数,或为沙峰随洪水波运动过程中的坦化系数;KS为泥沙颗粒平均运动时间(min).

研究河段中由于水流的非稳定性,水流和泥沙条件沿程也会发生变化,河段的平均冲淤速率表达式为[16]

$ q=\alpha \cdot \zeta\left(C I_{1}-C_{1}\right) \cdot I_{q}+(1-\alpha) \zeta \cdot\left(C I_{2}-C_{2}\right) \cdot Q $ (23)

式中, α为上下游断面冲淤速率比重系数;CI1CI2分别为上、下断面的水流挟沙能力;C1C2分别为上、下断面的实际含沙量;Iq为上断面流量;Q为下断面流量;ζ为冲淤系数.

将式(22)、(23)代入(21)再求差分得坡面泥沙演算式为:

$ O_{t}=C G \cdot Q_{1}+(1-C G)\left(I_{t}+I_{t-1}\right) / 2 $ (24)

式中, CG为消退系数, Ot为断面出口坡面泥沙速率(kg/s).

泥沙随水流运动而运动,泥沙汇集与水流汇集有许多相似之处,类似地可推得沟道泥沙演算式为:

$ S_{t}=C S_{0} \cdot O_{t}+C S_{1} \cdot Q_{t-1}+C S_{2} \cdot S_{t-1}+e+f $ (25)
$ \left\{\begin{array}{l}{C S_{0}=(0.5 \Delta t-K S G \cdot X S-0.5 \alpha \cdot \zeta \cdot \Delta t) / R R} \\ {C S_{1}=(0.5 \Delta t+K S G \cdot X S-0.5 \alpha \cdot \zeta \cdot \Delta t) / R R} \\ {C S_{2}=1-C S_{0}-C S_{1}} \\ {e=\zeta\left[0.5 \alpha\left(I_{1} \cdot C I_{1}+I_{2} \cdot C I_{2}\right)+0.5(1-\alpha) \cdot\left(Q_{1}+Q_{2}\right)\right] / R R} \\ {f=\zeta\left[K X D \cdot 0.5(I+Q) \cdot 0.5\left(C_{1}+C_{2}\right) \cdot\left(\frac{0.5(I+Q)}{B Q M}\right)^{a .5}\right] / R R} \\ {R R=0.5 \Delta t-K S G \cdot X S+K S G+0.5 \zeta \cdot \Delta t(1-\alpha)}\end{array}\right. $ (26)

式中, S为河道出口输沙速率(kg/s), KXD为河道动力非线性冲刷系数, BQM为断面平均流量(m3/s), KSG为泥沙颗粒在河段中的平均传播时间(min), XS为比重系数.

2 模型检验

本文选用黄土高原大理河流域中王茂沟、韭园沟两个小流域和岔巴沟、绥德两个大流域的实测水文泥沙资料,对以上建立的流域水沙耦合物理概念模型进行模拟分析检验.

2.1 流域概况和资料

王茂沟流域(37°34′~37°36′N, 110°20′~110°22′E)是无定河的一条二级支沟,流域面积5.97 km2.多年平均降水量513 mm,降水主要集中在7-9月,多暴雨,降水量年内分配不均,年际变化率大,汛期降水量占年降水量的70 %以上.流域水土流失严重,侵蚀方式主要为水力侵蚀和重力侵蚀,95 %的泥沙流失量集中在汛期,未实施各项水保措施前流域年平均侵蚀模数为1.8万t/(km2·a)[22].

韭园沟流域(37°33′~37°38′N, 110°16′~110°26′E)是无定河中游一个分支流域,流域面积70.7 km2,多暴雨,多年平均降水量约为500 mm,降水多集中在汛期的6-9月,年内分配极不均匀,年际变化大,汛期降水量占年总降水量的70 %以上.流域内梁峁起伏,沟壑纵横,地形破碎,治理前流域多年平均侵蚀模数为1.8万~1.99万t/(km2 ·a)[23].

岔巴沟流域位于黄河支流无定河水系上,属于黄土丘陵沟壑第一副区,流域面积187 km2,曹坪是岔巴沟流域的控制出口水文站.流域内干旱少雨,多年平均降雨量为447 mm,降水年内分配极不均匀,6-9月降水占全年降水的80 %以上,且降雨多以短历时强降雨出现,暴雨强度大,洪水含沙量高,流域平均侵蚀模数为2.22万t/(km2·a), 年内径流和输沙量变化十分明显[24].

大理河是无定河最大的一级支流,流域面积3906 km2,绥德站是大理河流域的出口控制站.流域植被稀疏,地形破碎,气候干燥,水土流失严重,治理前流域多年平均输沙模数为1.6万t/(km2·a)[24].多年平均降水量为439.5 mm,年内分配不均,降雨大都为暴雨形式,强度大,历时短,5-9月降水量占全年总降水的80 %以上.大理河流域内的雨量站和流量站站点分布如图 5所示.

图 5 大理河流域站点分布 Fig.5 Distribution of the sites in Dali River Basin

模型采用的各流域实测资料统计见表 1.

表 1 各流域资料情况统计 Tab. 1 Statistical table of data in various basins
2.2 参数率定原理及结果

本文建立的水沙模型采用函数曲面率定方法[15]对产流、汇流、产沙、汇沙四部分参数进行率定,各部分参数相互独立,按照先水流后泥沙的顺序进行参数率定.

流域产沙模型的应用需预先确定沟道平均流速,坡面水流和CMBSREMMCGMACM 5个参数.其中产沙层次中的参数CMCGM,用直接法根据其物理概念信息通过观测资料直接分析确定[25],其余参数用水沙模型进行率定.

α0为土壤抗侵蚀能力为零的面积比,与流域前期气候条件有着密切关系.当前期气候干燥时,α0值随坡面疏松土壤的增多而增大.因此,在应用中把α0作为时变参数,其余作为时不变参数,经研究[17],该值可以用洪前土壤含水量来反映,可表示为:

$ \alpha_{0}=A W[1-(1-W / W M) / W B]^{B W} $ (27)

式中, AWBW均为常系数, WWM分别为洪前实际土壤含水量和土壤含水量上限(mm), WB= $\frac{1}{n} \sum\limits_{i=1}^{n}\left(1-\frac{W_{i}}{W M}\right)$, n为洪水次数.

沟蚀产沙中流速中断面平均流速(u)用以下经验关系确定:

$ u=A V\{[\ln (Q+1)] / L Q B\}^{B V} $ (28)

式中, LQB为ln(Q+1)的时间平均值(m2/s);AV为ln(Q+1)=LQB时的沟道平均流速(m/s);BV为常参数.

式(28)亦可以写为:

$ \frac{u}{u M}=[\ln (Q+1) / L Q B]^{B V} $ (29)

各流域参数率定结果见表 2表 3.

表 2 各流域水流参数率定结果 Tab. 2 Calibration results of flow parameters in various basins
表 3 各流域泥沙参数率定结果 Tab. 3 Calibration results of sediment parameters in various basins
2.3 模拟结果分析

产流和产沙分别以次洪径流深和产沙量为标准,比较模型计算值与实测值;汇流和汇沙分别以流量、输沙率过程为标准,比较模型模拟过程与实际过程的吻合程度.模型的有效性指标根据模型计算结果与实测之差来估计, 使用的定量化指标为:

$ E=\sum\limits_{i}\left|X_{i}-X C_{i}\right| / \sum\limits_{i} X_{i} $ (30)
$ D C Q=1-\sum\left(Q_{i}-Q C_{i}\right)^{2} / \sum\left(Q_{i}-\overline{Q}\right)^{2} $ (31)

式中, XXC分别为评估量实测值和模型计算值, E为相对误差, DCQ为流量过程确定性系数, $\overline{Q}$为实测流量的平均值(m3/s),QiQCi分别为第i时刻的实测流量和计算流量.

根据模拟精度评定规则,次洪径流深的合格标准以实测值的20 %作为许可误差,次洪产沙量的合格标准以相对误差小于30 %为准.以岔巴沟流域为例,表 4表 5分别列出各场次洪的水流泥沙模拟结果及各项精度统计,表内合格度为“√”表示合格,“×”表示不合格.

表 4 岔巴沟流域水流模拟 Tab. 4 Runoff simulation in Chabagou drainage area
表 5 岔巴沟流域泥沙模拟 Tab. 5 Sediment simulation in Chabagou drainage area

表 4中可以看出,各场次洪的模型计算径流深与实测值的相对误差均在20 %内,合格率为100 %,实测与计算洪峰值相差不大,但计算洪峰值偏小,主要是降雨均化所致.除3场洪水外,其余次洪的峰现时差均在一个时段内,30场洪水的确定性系数整体较高,平均确定性系数达0.713,说明将该模型应用于岔巴沟流域的水流模拟有较高精度.

进一步分析个别场次洪确定性系数不高的原因,主要是受雨量资料时段均化及实测资料观测误差的影响,由于四个流域均处于干旱半干旱地区,地面径流主要为超渗产流,计算时段受不同流域尺度的流域汇流时间影响需要6 min或更短,而雨量观测资料时段最短1 h,多为3 h、8 h或更长,经过线性内插产生的时段均化给产流时间分布和计算流量过程带来较大误差.且部分次洪存在实测流量资料缺测、漏测或误测等现象,应用模型模拟时对实测资料进行插补展延,导致实测流量过程存在水平段、斜线段或者存在复峰形式,汇流过程模拟不佳,进而影响确定性系数.

泥沙模拟中,大部分次洪的计算产沙量与实测值相对误差均小于30 %,合格率较高(表 5).产沙量模拟主要受坡面水力侵蚀参数ACM和冲淤系数ζ的影响,从式(16)、(23)可知,这两个参数均受水力要素影响,表明产沙量与水流在运动中的变化幅度有关.因此,为了减少水流模拟存在的误差带入到泥沙模拟中,泥沙模拟采用的时段产流量已经过实测径流深校正,且沟蚀产沙计算中的沟道流量采用实测值,进而提高泥沙模拟精度.

表 5可知,1973年后的实测产沙量明显较1973年前少,主要是由于岔巴沟流域在1975年后建成大量的坝库工程并发挥显著蓄水拦沙效益,流域下垫面变化剧烈,不同时期的水沙关系变化趋势如图 6所示,从图中可以看出天然期1960-1973年与人类干预期1977-2009年的水沙关系线明显不同,流量、输沙率较天然期显著减小,而概念性模型受一定的概化影响,在泥沙研究中只考虑主要影响因素,故模型中参数多为时不变参数,因此利用多年实测资料率定出的同一参数集会使计算产沙量在人类干预期系统偏大.

图 6 不同时期岔巴沟流域水沙关系 Fig.6 Relationship of water and sand of the Chabagou River Basin during different periods

为了进一步说明模拟结果精度,选取岔巴沟流域2场典型洪水(以19710723、19770811号次洪为例)展示水流和泥沙模拟过程,其流量和输沙率过程见图 7图 8,从图上可以看出水流泥沙过程模拟与实测过程一致,主要体现在洪量、沙量、峰形和峰值的拟合上.总体来说,模型模拟的水流和泥沙两部分结果均具有很高的精度,说明该模型是合理可行的.

图 7 19710723号次洪流量过程和输沙率过程 Fig.7 Flow and sediment transport rate process of No. 19710723 flood
图 8 19770811号次洪流量过程和输沙率过程 Fig.8 Flow and sediment transport rate process of No. 19770811 flood

概念性水沙耦合模型应用于4个流域的水流和泥沙模拟结果及各项精度统计可以看出,模型应用于不同流域尺度的水流、泥沙模拟均具有一定的模拟精度,可同时进行洪水预报及泥沙冲淤计算(表 6).水流模拟合格率均大于95 %,计算径流深与实测值相对误差小,且洪峰拟合较好,但确定性系数偏低,主要是由于降雨资料时段均化和实测资料的观测误差影响.其中,王茂沟、韭园沟流域的模拟结果相对于岔巴沟、绥德流域模拟结果较差的原因是前2个流域面积较小,且次洪大部分为小洪水,试验场的实测资料观测误差较大.泥沙模拟的合格率整体较水流偏低,产沙量计算值与实测值相差小但平均相对误差多为负值,表明计算值较实测值偏大,主要是由于流域内各项水保措施对泥沙影响较水流影响大,减水减沙效益显著,在人类干预期的产沙量明显较天然期少,导致模型采用长期实测资料率定出的参数会使模拟结果分阶段地产生系统偏差,进而影响模拟精度.

表 6 各流域水流泥沙模拟精度统计 Tab. 6 Statistical accuracy statistics of water and sediment in various basins
3 结论

流域水流基本规律是研究泥沙规律的基础,经验统计模型物理意义不明确,难以揭示泥沙的自然机制,且精度低,适用范围受限.本文建立的概念性水沙耦合模型将水流模块与泥沙模块耦合,各部分结构、公式及参数均有清楚的物理意义,可以模拟水流和泥沙分量,模拟流域出口断面的流量过程和输沙率过程.泥沙规律研究参照水流模型产、汇流研究方法,提出泥沙概念性模拟模型,根据水流对颗粒拖曳力和颗粒对水流抗动力的物理概念,提出水流挟沙能力和土壤抗侵蚀能力概念,并建立流域内土壤抗侵蚀能力分布曲线描述空间变异性,进而推求流域的坡面侵蚀产沙量,沟蚀产沙量采用基于拜格诺河道水流悬移质泥沙公式建立的概念性沟蚀产沙公式推求.

将模型应用于黄河中游大理河流域不同流域尺度的4个实际流域,利用多年实测水文泥沙资料并应用函数曲面参数率定方法对模型各部分参数进行率定,对模型进行应用模拟检验,分析结果表明:该模型具有很高的模拟精度,水流泥沙模拟合格率高,流量、输沙率过程与实测过程一致,符合水流与泥沙的运动规律,可以很好地模拟不同流域尺度的水沙过程.其中,大流域的模拟精度较小流域模拟精度高,主要是由于较小试验流域的实测资料观测误差较大,且雨量资料时段均化给产流时间分布及流量过程模拟带来较大影响.该模型可用于整个黄土地区大中小流域长系列的水沙模拟问题,利于深入研究黄河流域的泥沙规律,以满足水土保持措施效益分析和黄河中游入黄水沙变化原因及趋势分析等问题研究的大尺度、长系列模拟的需要,为人类活动和气候变化对流域水沙变化的影响提供基础数据.

通过概念性模拟研究建立概念性水沙耦合模型,其研究方法的主要特点是假设和概化,突出自然规律的主要影响因素,忽略次要因素,概化为与实际规律较为接近的关系,从而便于研究.如流域产沙计算公式均建立在流域下垫面条件相对均匀的前提下,即忽略地形、土质、人类活动等因素对产沙的影响,着重考虑气候因素.黄河流域内大量水利工程的修建及植树造林等水保措施的开展导致流域下垫面条件及水沙关系变化剧烈,由于模型中水流和泥沙参数大都为时不变参数,使得人类活动成为制约水沙模拟的关键因素,但模型输入因素受小范围水土保持措施的影响不大,模拟较短年际尺度上的水沙过程时精度比模拟较长年际尺度时更高.可进一步改进模型参数结构,提高模型应用于人类活动影响和气候变化年份下的长系列水沙模拟精度,提高模型通用性,同时分离评估减水型措施引起的减水、减沙和植被型措施引起的减水减沙效果,对该模型的推广应用可做进一步分析研究.

4 参考文献

[1]
Yao WY, Gao YJ, An CH et al. Analysis of trend of runoff and sediment load in upper and middle reaches of Yellow River at century scale. Advances in Science and Technology of Water Resources, 2015, 35(5): 112-120. [姚文艺, 高亚军, 安催花等. 百年尺度黄河上中游水沙变化趋势分析. 水利水电科技进展, 2015, 35(5): 112-120.]
[2]
Liu C, He Y, Liu A. Key drivers of changes in sediment loads of rivers. Advances in Science and Technology of Water Resources, 2017, 37(1): 1-7. [刘成, 何耘, 刘桉. 河流输沙量变化的主要驱动因素. 水利水电科技进展, 2017, 37(1): 1-7.]
[3]
Yao W, Xiao P, Shen Z et al. Analysis of the contribution of multiple factors to the recent decrease in discharge and sediment yield in the Yellow River Basin, China. Journal of Geographical Sciences, 2016, 26(9): 1289-1304. DOI:10.1007/s11442-016-1227-7
[4]
Shi FC, Zhang R. Cause analysis and recognitions on the recent sharp decreasing of the Yellow River water and sediment amount. Yellow River, 2013, 35(7): 1-3. [史辅成, 张冉. 近期黄河水沙量锐减的原因分析及认识. 人民黄河, 2013, 35(7): 1-3. DOI:10.3969/j.issn.1000-1379.2013.07.001]
[5]
Beven KJ, Kirkby MJ. A physically based variable contributing area model of basin hydrology. Hydrological Science Bulletin, 1979, 24(1): 43-69. DOI:10.1080/02626667909491834
[6]
Wigmosta MS, Vail LW, Lettenmaier DP. A distributed hydrology vegetation model for complex terrain. Water Resources Research, 1994, 30(6): 1665-1679. DOI:10.1029/94WR00436
[7]
Si W, Bao WM, Jiang P et al. A semi-physical sediment yield model for estimation of suspended sediment in Loess Region. International Journal of Sediment Research, 2017, 32(1): 12-19. DOI:10.1016/j.ijsrc.2015.10.002
[8]
Li WJ, Wang XK, Li DX et al. A physically-based distributed watershed water erosion prediction model. Journal of Hydraulic Engineering, 2012, 43(3): 264-274. [李文杰, 王兴奎, 李丹勋等. 基于物理过程的分布式流域水沙预报模型. 水利学报, 2012, 43(3): 264-274.]
[9]
Wischmeier WH, Smith DD. Predicting rainfall erosion losses-a guide to conservation planning. Agic Handbook, 1978, 537.
[10]
De Roo APJ, Jetten VG. Calibrating and validating the LISEM model for two data sets from the Netherlands and South Africa. Catena, 1999, 37(3/4): 477-493.
[11]
Laflen JM, Lane LJ, Foster GR. WEPP: A new generation of erosion prediction technology. Journal of Soil & Water Conservation, 1991, 46(1): 34-38.
[12]
Wan YY, Jin ZW, Huang RY. Outlook and prospect on sedimentation model study. South-to-North Water Transfers and Water Science & Technology, 2006, 4(1): 48-51. [万远扬, 金中武, 黄仁勇. 泥沙模型研究述评与前景展望. 南水北调与水利科技, 2006, 4(1): 48-51. DOI:10.3969/j.issn.1672-1683.2006.01.014]
[13]
Bao WM. A conceptual flow-sedimentation coupled simulation model for small basins. Geographical Research, 1995, 14(2): 27-34. [包为民. 小流域水沙耦合模拟概念模型. 地理研究, 1995, 14(2): 27-34. DOI:10.3321/j.issn:1000-0585.1995.02.004]
[14]
Bao WM, Chen YT. A conceptual flow-sedimentation coupled simulation model for large basins. Advances in Water Science, 1994, 5(4): 287-292. [包为民, 陈耀庭. 中大流域水沙耦合模拟物理概念模型. 水科学进展, 1994, 5(4): 287-292. DOI:10.3321/j.issn:1001-6791.1994.04.005]
[15]
Bao WM, Zhang XQ, Zhao LP. Parameter estimation method based on parameter function surface. Science China Technological Sciences, 2013, 56(6): 1485-1498. DOI:10.1007/s11431-013-5224-3
[16]
Tang LQ, Chen GX. Application of physical conceptual model in evaluating the benefits of water and soil conservation. Journal of Hydraulic Engineering, 1998, 29(9): 62-65. [汤立群, 陈国祥. 物理概念模型在水保效益评价中的应用. 水利学报, 1998, 29(9): 62-65. DOI:10.3321/j.issn:0559-9350.1998.09.012]
[17]
Bao WM. A conceptual modelling study for small basin sediment yields in Loess regions. Advances in Water Science, 1993, 4(1): 44-50. [包为民. 黄土地区小流域产沙概念性模拟研究. 水科学进展, 1993, 4(1): 44-50. DOI:10.3321/j.issn:1001-6791.1993.01.007]
[18]
Bao WM. A tentative study of conceptual sedimentation routing model. Journal of Hohai University, 1990, 18(6): 24-29. [包为民. 概念性汇沙模型初探. 河海大学学报, 1990, 18(6): 24-29. DOI:10.3321/j.issn:1000-1980.1990.06.004]
[19]
Bao WM. Improvement and application of the Green-Ampt infiltration curve. Yellow River, 1993(9): 1-3. [包为民. 格林——安普特下渗曲线的改进应用. 人民黄河, 1993(9): 1-3.]
[20]
Zhao RJ. Hydrological simulation of watershed. Beijing: Water Resources and Electric Power Press, 1984: 55-105. [赵人俊. 流域水文模拟. 北京: 水利电力出版社, 1984: 55-105.]
[21]
Chen JQ. Application of Muskingum methods in flood forecasting. China Rural Water and Hydropower, 2005(5): 33-35. [陈金庆. 马斯京根流量演算法在洪水预报中的应用. 中国农村水利水电, 2005(5): 33-35. DOI:10.3969/j.issn.1007-2284.2005.05.013]
[22]
Gao HD, Li ZB, Li P et al. Concept and calculation methods of erosion control degree: A case study of the Wangmaogou Watershed. Science of Soil and Water Conservation, 2013, 11(1): 17-24. [高海东, 李占斌, 李鹏等. 流域侵蚀控制度的概念与计算方法——以王茂沟流域为例. 中国水土保持科学, 2013, 11(1): 17-24. DOI:10.3969/j.issn.1672-3007.2013.01.003]
[23]
Wang D, Li ZB, Li P et al. Evaluation of overall distribution of check dam system in the Jiuyuangou Watershed. Research of Soil and Water Conservation, 2016, 23(5): 49-55. [王丹, 李占斌, 李鹏等. 韭园沟流域淤地坝坝系布局评价. 水土保持研究, 2016, 23(5): 49-55.]
[24]
Qi JY, Cai QG, Cai L et al. Scale effect of runoff and sediment reduction effects of soil and water conservation measures in Chabagou, Dalihe and Wudinghe Basins. Progress in Geography, 2011, 30(1): 95-102. [綦俊谕, 蔡强国, 蔡乐等. 岔巴沟、大理河与无定河水土保持减水减沙作用的尺度效应. 地理科学进展, 2011, 30(1): 95-102. DOI:10.11820/dlkxjz.2011.01.012]
[25]
Bao WM. Automatic calibration of Xinan Jiang model parameters. Journal of Hohai University, 1986, 14(4): 22-30. [包为民. 新安江模型的参数自动率定. 河海大学学报, 1986, 14(4): 22-30. DOI:10.3321/j.issn:1000-1980.1986.04.003]