丁倩倩, 刘友存, 焦克勤, 边晓辉, 刘 燕
(1.江西理工大学 建筑与测绘工程学院, 江西 赣州 341000; 2.江西理工大学 资源与环境工程学院,江西 赣州 341000; 3.中国科学院 西北生态环境资源研究院, 甘肃 兰州 730000)
河流的水沙变化不仅受流域内降水、下垫面和河床演变等自然因素的影响,亦是水利工程和河道改造等人为因素的反映[1]。如何定量区分自然因素和人为因素对河流水沙变化的影响,已成为国内外河流水沙研究的热点。同时,水沙作为天然水体中所携带的固体颗粒,不仅对河道形态、洪水特征、河道生态和河口沉积地貌等有着很大的影响[2-3],而且是水利和防洪工程、水资源合理利用和保护以及水土保持等方面必须考虑的重要因素[4-5]。自1753年布朗姆斯(Brahms)提出水沙起动流速与重量的1/6次方成正比关系以来[6],国际上很多学者[7-8]对河流水沙变化及其影响因素以及河流输沙变化所引起的河口地区的冲淤演变和环境变化等进行了观测研究。1950s以来,随着我国对河流水沙治理的日益重视,使得水沙学科体系逐步得到完善,并在非均匀悬移质输沙[9]、泥沙运动统计理论[10]、水流挟沙能力[11-12]、水库淤积和河流模拟[13-14]等方面取得了一批重要成果。近年来,围绕黄河小浪底[15]和长江三峡[16]工程的水沙问题以及黄河水沙调控[17]等重大工程问题的研究也取得了引人注目的进展。然而,针对江南典型流域输沙过程的时间序列特征及成因分析,以及进一步探讨河流输沙与气候变化相互关系的研究较少,尤其是定量地分析流域的径流和输沙时空变化特征,并探讨水沙过程对自然因素与人类活动的响应研究更少。为此,本研究以赣江重要的支流桃江流域为研究对象,运用数理统计分析(双累积曲线和Mann-Kendall法)和3S技术等方法,依据流域内过去36 a逐日实测降水、径流和输沙资料,分析了桃江流域输沙过程的时间序列特征及其与气候变化的相互关系,进一步探讨桃江河流输沙与气候和人类活动的关系,以期为桃江乃至赣江流域的社会、经济与生态环境的可持续发展,特别是流域内土壤侵蚀治理、土地利用和功能调整等提供科学的依据。
桃江流域位于24°28′—25°53′N,114°11′—115°19′E,流域面积约8 440 km2。桃江流域主要由中低山和丘陵组成,海拔从南至北逐渐降低[17]。流域上游和东侧边缘山地海拔较高,呈弧形,最高海拔为1 383 m;中下游地势虽有起伏,但较为和缓,平均海拔约300 m,而与贡江交汇处海拔低于200 m。桃江流域属于亚热带湿润季风气候区,气候温润,光照充足,雨量丰沛。年均降水量为1 580.6 mm,降水年内分配不均匀,3—8月(汛期)约占全年降水量的74%左右。降水量年际差别较大,丰水年约为枯水年的2倍。桃江的丰水期、平水期和枯水期特征很明显,而径流的年内和年际变化受降水控制。
选取了桃江流域控制站—居龙滩水文站36 a的逐日径流量和输沙量实测资料;以及流域内27个站点的逐日降水量实测资料,运用泰森多边形的插值方法由各站点降水量推算出整个流域的降水量,从而得到1980—2015年桃江流域的逐日降水量。
2.2.1 水文法 水文法是用流域水沙的观测资料分析水土保持措施减水减沙作用的一种方法[19]。它以治理前(基准期)的实测水沙资料,建立降水产水产沙量数学模型,然后将治理后的降水因子代入所建模型,计算出相当于治理前的产水产沙量,再与治理后的实测水沙量进行对比,其差值即为经过治理后减少的水沙量。如果将治理前的实测水沙量视为天然产水产沙量,那么,根据治理后的降水因子用产流产沙模型计算出的产流产沙量,就相当于治理后降水条件下所产生的天然产流产沙量,两个时段天然产流产沙量的差值即为降水变化对产流产沙的影响量。相应地,如果将模型计算的天然产流产沙量与同一时段实测的水沙量相减,即可视为人类活动对产流产沙的影响量。
2.2.2 双累积曲线法 双累积曲线法是检验两个参数间关系一致性及其变化的常用方法。Searcy等[20]认为,如果双积累曲线的斜率发生突变(Break)则意味着2个变量之间的比例常数发生了改变或者其对应的累积值的比可能根本就不是常数。若接受2个变量累积值之间直线斜率已发生改变,那么斜率发生突变点所对应的年份就是两个变量累积关系出现突变的时间。
设有2个变量X(参考变量或基准变量)和Y(被检验变量),在N年的观测期内,有观测值Xi和Yi,其中i=1,2,3,…,N。首先对变量X和变量Y按年序计算各自的累积值,得到新的逐年累积序列Xi′和Yi′,其中i=1,2,3,…,N,即:
(1)
(2)
然后,在直角坐标系中绘制两个变量所对应点累积值的关系曲线。绘制的曲线图一般以被检验的变量为纵坐标(即Y轴)、参考变量或基准变量为横坐标(即X轴)。
2.2.3 Mann-Kendall法(M-K法) M-K法[33]是用来评估水文气候要素的时间序列趋势的检验方法。
设一时间序列为:X1,X2,X3, …,Xn,构造一序列mi,mi表示Xi>Xj(1≤j≤i)的样本累积数。定义dk:
(3)
其中:dk的均值、方差定义为:
(4)
(5)
在时间序列随机独立的假定下,定义统计量:
(6)
式中:UFk——标准正态分布,给定显著性水平a,查正态分布表得到临界值t,当UFk>t,表明序列存在一个显著的增长或减少趋势,绘制时间序列UFk的曲线L1,通过显著性检验可得知其是否具有显著性趋势。将时间序列X按逆序排列,再重复上述的计算过程,并将计算值乘以-1,得到UBk,UBk的时间序列曲线用L2表示,当曲线L2超出置信区间,即表示存在明显的变化趋势。若L1和L2在置信区间内存在交点,则此点为突变点的开始。
2.2.4 连续小波变换(CWT)法 如果函数Ψ(t)满足以下容许条件:
(7)
则称Ψ(t)为一个基本小波[21],并定义如下的积分变换:
(8)
式中:f(t)——Ψ(t)为基本的积分连续小波变换;a——尺度因子,表示与频率相关的伸缩;b——时间平移因子。
3.1.1 水沙年际变化 依据桃江流域1980—2015年的年均径流量和输沙率的实测数据,在1980—2015年的时间序列上,径流量和输沙率总体上呈降低趋势,但径流量的降低趋势不是很明显,只是略有下降。输沙率的下降趋势则较为明显,降低幅度比较大。在输沙变化过程中,1991和2004年达到极小值, 同期的径流量也为极小值。从2007—2015年输沙量呈显著下降趋势,而同期的径流量并没有明显的减小。另外,年平均输沙率通过了0.01显著性水平的相关性检验,而年平均径流量则未通过显著性相关检验。
3.1.2 水沙年内变化 依据桃江流域1980—2015年的年内月平均径流量和输沙率的实测数据,桃江流域的径流量在上半年不断增加,6月达到最大值,下半年则不断减小。总体上,径流量在3—6月较大,为丰水期。输沙率在3—6月远远大于其他月份,其中4月的输沙率达到最大值为94.50 kg/s,5—6月稍有递减,分别为92.83和87.58 kg/s。
3.1.3 水沙变化的突变分析 采用M-K法对桃江流域逐年径流量和输沙率做M-K突变检验,结果如图1所示。由图1可以看出,径流量的UFk统计量在1997—1998年,2001—2003年为正值,其余年份均为负值,表明径流量呈下降趋势。但是UFk统计量值除了1989和1991年超过了显著性水平为 的临界值,其余年份的统计量值均未超过该显著性水平的临界值,故径流量的下降趋势并不显著。从图2总体来看,径流量的UFk统计量值波动变化较多。在1991年之前,虽然其值上下波动过3次,但总体为下降趋势,其变化值域没有超过0。在1991—1998年UFk呈上升趋势,而在1997—1998年超过了0值,表明径流量为增加,紧接着在1998—2002年UFk呈先下降后上升。从2002年开始UFk开始呈明显的下降,但其值并未超过显著水平 的临界值,说明径流量虽有下降,但并不显著。正逆序列的M-K统计量值UFk和UBk在2002年后出现了7个交叉点(分别为2003—2004年,2005年,2007年,2009—2010年,2010—2011年,2011—2012年和2012—2013年)。结合桃江流域水利工程情况(表2)可知,五渡港水库于2004—2005年加固维修,其库容量为3.33×107m3,桃江电站于2004年建成投入使用,其库容量为3.71×107m3,而居龙滩电站于2007年建成投入运行,其库容量为7.76×107m3,3个水电站都属于大型水库,其拦水拦沙作用比较明显,故径流量在上述几个年份出现突变是合理的。采用M-K法对桃江流域逐年输沙率进行的突变检验,由图1可以看出,输沙率的UFk统计量值除了1996年为正值,其余年份均为负值,表明输沙率在整个时间序列中以下降为主。在2007年输沙率的UFk统计量值超过了显著性水平 的临界值,其下降趋势更为显著。正逆序列的M-K统计量值UFk和UBk在2005—2006年间出现一个交叉点,且交叉点位于显著性水平临界值的范围之内,因此,可以认为2005年即为输沙率的突变年份。
图1 桃江流域径流量和输沙率的M-K突变检验
总体看来,径流量和输沙率的M-K统计量值绝大部分为负值,故径流量和输沙率均呈现减小趋势。在1980—1982年径流量和输沙率的M-K统计量值UFk不仅为负值,而且呈现递减趋势,到了1983年径流量和输沙率的M-K统计量值UFk又上升到1980年的位置。1983—1991年径流量的M-K统计量值UFk呈现波动减小,同期的输沙率M-K统计量值UFk同样呈现波动减小,但在其变化过程中波动起伏很小,波动起伏的时间节点基本与径流量M-K统计量值UFk的时间节点相吻合。1991—1998年径流量的M-K统计量值UFk呈递增趋势,同期的输沙率M-K统计量值UFk同样呈现上升趋势,但其上升的速度较缓慢,且上升幅度没有前者大。从1999年开始,径流量和输沙率的M-K统计量值UFk均呈现下降趋势,前者的下降幅度相对较小,且起伏波动较多。后者则呈现递减趋势,1999—2006年虽有较小幅度的起伏,但总体下降幅度较大。从整个时间序列来看,径流量和输沙率的M-K统计量值UFk的变化趋势基本上是同步的,只是在1999年前,前者的变化幅度大,起伏波动较明显;而在1999年后,后者的减小幅度要比前者大很多。
3.1.4 水沙变化周期分析 图2中细实线表示影响锥(the cone of influence, COI),在该曲线以外的能量谱由于受到边界效应的影响而不予考虑,粗实线表示通过95%置信水平显著性检验的临界值。桃江流域1980—2015年的逐月降水量、径流量和输沙量在时频域中都存在不同尺度且通过0.05置信水平检验的显著周期(图2)。降水量在1980—1990年和1991—2015年均在时频域低能区存在9—15个月的周期。径流量在1980—1990年间、1991—1997年和2004—2013年在时频域低能区分别存在6—16个月,10—16个月和10—15个月的周期。输沙量在1981—1990年和1991—1997年在时频域低能区存在10—16个月的周期。
3.2.1 降水 桃江流域在1980—2015年时段内逐年降水量如图3所示。从图3可以看出,在1983年达到最大值为2 007.9 mm。从1983—1989年为递减波动趋势,在1989年达到最小值为1 277.2 mm,1990年又上升到1 704.4 mm,而到1991年再次降到整个时间序列的最小值为1 106.9 mm。1992年上升到1 912.4 mm,在1992—2002年时段内除1999年为1 404.6 mm外,其余年份均达到1 500 mm以上。到2003和2004年降水量再次降低,且均低于1 200 mm。2005—2015年时段内呈现波动变化,但在1 500 mm上下波动,所有值均大于1 200 mm。从1980—2015年的降水量变化过程来看,呈现微弱的递减趋势,其变化趋势可用线性拟合(图3),但年平均降水量没有通过显著性相关检验。
图2 桃江流域连续小波变换能量谱
桃江流域的年平均降水量为1 580.63 mm,由流域降水量的年内分布情况可知(图3),3—8月为丰水期,占全年降水量的73.7%;10—1月为枯水期,占全年降水量的12.8%;2和9月为平水期,占全年降水量的13.5%。
由1980—2015年桃江流域月降水量的时间序列单累积曲线可以看出(图4a),桃江流域的降水量随着时间其累积值与时间呈线性相关关系,其线性拟合方程的拟合优度R2达到0.999 7,说明桃江流域的降水量随着时间的变化并没有显著性的增加或减少。由1980—2015年桃江流域月径流量的时间序列单累积曲线可知(图4b),桃江流域径流量总体减小趋势不明显,但存在小幅度的下降,结合图1,可将2003年确定为径流量的一个分界点,并以此求出径流量的减少量。第一阶段(1980—2002年)的线性拟合的拟合优度R2为0.998 3,第二阶段(2003—2015年)的R2为0.996 9,第二阶段相对第一阶段的拟合曲线其斜率有所减小,说明径流量在第二阶段相对第一阶段有所减少。从2004年开始,年均径流量减少7.26×1010m3。由流域1980—2015年的逐月输沙率的时间序列单累积曲线可知(图4c),输沙率在2005年发生了突变,其输沙量明显存在两个阶段的变化。第一阶段(1980—2004年)的线性拟合的拟合优度R2为0.995 7,第二阶段(2005—2015年)的R2为0.896 2。第二阶段的线性拟合曲线斜率明显小于第一阶段,也就是说,从2005年开始每年输沙率的累积值减小,平均每年的输沙量减少量为8.82×107t。由桃江流域1980—2015年降水量和径流量的时间序列双累积曲线可知(图4d),双累积曲线在2003年发生突变,第一阶段(1980—2002年)拟合曲线的拟合优度R2为0.999 6,第二阶段(2003—2015年)的R2为0.999,第二阶段拟合曲线的斜率要小于第一阶段的。根据“水文法”计算得到(表1),2003—2015年桃江流域人类活动造成的年均径流量较1980—2002年的年均径流量减少了5.77×1010m3,占年均总径流量减少量的76.77%;其次因降水量减少引起的径流量减少了1.49×1010m3,占年均径流量减少量的23.23%。人类活动和降水量变化对径流量的影响之比约为0.8∶0.2。
表1 桃江流域减水减沙量
图3 桃江流域降水量年际和年内变化
图4 降水量、径流量和输沙率单累积和双累积曲线
由桃江流域1980—2015年降水量和输沙率的时间序列双累积曲线可知(图4e),双累积曲线在2005年发生突变,可将曲线分为两个阶段。第一阶段(1980—2004年)的线性拟合曲线的拟合优度R2为0.998 3,第二阶段(2005—2015年)的R2为0.993 1。由表1可以看出,人类活动造成桃江流域第二阶段的输沙量比第一阶段年均减少量为8.44×107t,占总体输沙量减少量的95.7%;其次为降水量减少造成的输沙量减少量为3.82×106t,占总体年均输沙量减少量的4.3%。人类活动和降水量变化对输沙量的影响之比约为0.95∶0.05。由流域1980—2015年径流量和输沙量的时间序列双累积曲线曲线可知(图4f),整个趋势与降水量和输沙率的时间序列双累积曲线图相一致,在2005年发生突变。1980—2004年的线性拟合曲线的拟合优度R2为0.997 3,2005—2015年的R2为0.9174,且第二阶段的拟合曲线斜率小于第一阶段的斜率。由径流量和输沙量间的线性关系可知,2005—2015年年均径流产生的输沙量减少量为6.82×107t。
由图4可知,输沙量的时间序列单累积值从2005年开始呈现阶段性的跃升变化。而降水量和输沙率、径流量和输沙率的双累积图也呈现同样的变化,故这些线性曲线的拟合优度均不是高。结合桃江流域输沙率的实测资料可知,输沙率在2005年6月,2006年7月,2007年6月,2008年6月,2010年6月,2011年5月,2012年3月,2013年5月,2014年5月和2015年5月均为极大值,分别为167,140,163,66.6,79.9,57.2,48.6,50.8,52.5和52.6 kg/s, 最小输沙率为0.37 kg/s,除去这些极大值之后的月平均输沙率为4.66 kg/s,输沙率的极大值要比平常值高出几十倍甚至数百倍。故而,输沙率的累积值会出现阶段性的跃升变化,进而导致降水量和输沙率、径流量和输沙率的双累积曲线呈现同样的变化。
对输沙率做归一化处理,剔除极值的影响,重新做时间序列单累积图以及与降水量和径流量的双累积曲线,得到这些线性曲线的拟合优度均在0.99以上(图5)。
图5 桃江流域输沙率归一化处理的相关性分析
3.2.2 水利工程 在人类活动对河道水沙变化的影响因素中,水土保持措施和水利工程建设是影响河道径流和输沙变化的主要因素。本文主要讨论水利工程措施对河道来水来沙量的影响。
桃江流域的水利工程至2015年共有9项,其中水电站和水库分别为3座(居龙滩、龙头滩和桃江电站)和6座(龙井、上迳、白兰、五渡港、走马垅和中村水库),其主要特征参数表2。而对1980—2015年间径流量和输沙量有影响的水利工程,主要是1980年后建成运行的龙头滩、居龙滩和桃江电站。由前面分析得到的径流量和输沙量分别在2003和2005年开始发生突变。而3座电站(龙头滩、桃江和居龙滩电站)分别于2000,2004和2007年建成。由于龙头滩电站的库容较小(1.38×107m3),故对桃江流域的径流量和输沙量影响也较小;居龙滩和桃江电站的库容较大,尤其是居龙滩电站的库容达到7.76×107m3,控制面积约为7 739 km2,占桃江流域总面积的91.7%,故对桃江流域的径流量和输沙量的影响也较大。综合来看,电站建成时间与径流量和输沙量的突变时间基本一致。
从表1—2可知,水库单位库容的减水量和减沙量分别为4.49×107m3/a和6.57 t/a。运用趋势前推法可知,1980年之前建的走马垅、白兰、五渡港、上迳、龙井和中村水库对桃江流域的减水减沙量分别为4.73×1010m3/a和6.91×107t/a。由此可见,水库库容的大小与其拦水拦沙的能力成正相关关系。
表2 桃江流域水利工程情况
(1) 桃江流域1980—2015年的年均输沙率呈现明显的下降趋势,而年均径流量下降则不明显趋势;年内月均径流量在3—6月较大为丰水期,年内月均输沙率在4月达到全年最大值为94.50 kg/s,5—6月稍有递减。
(2) 降水量在1980—1990和1991—2015年的时频域低能区存在9—15个月的周期;径流量在1980—1990,1991—1997和2004—2013年的时频域低能区分别存在6—16个月、10—16个月和10—15个月的周期;输沙量在1981—1990和1991—1997年的时频域低能区均存在10—16个月的周期。
(3) 桃江流域1980—2015年的径流量和输沙率的M-K突变检验中得知,输沙率在2005年发生突变,并呈现明显的下降趋势;径流量则呈现不明显的下降趋势,2003年为突变开始点。桃江流域的径流量从2004年开始,年均径流量减少量为7.26×1010m3。输沙量从2005年开始,年均减少量为8.82×107t。
(4) 2003—2015年桃江流域人类活动造成的年均径流量较1980—2002年减少了5.77×1010m3,占年均径流量的76.8%;而因降水量减少引起的径流量减少量为1.49×1010m3,占年均径流量的23.2%。 水利工程造成2005—2015年输沙量比1980—2004年年均减少量为8.44×107t,占总输沙量减少量的95.7%;而降水量减少造成输沙量减少量为3.82×106t,占年均输沙量减少量的4.3%。
(5) 在桃江流域人类活动对水沙变化的影响因素中,只考虑水利工程建设的水土保持治理措施。桃江流域的水库单位库容的减水量和减沙量分别为4.49×106m3/a和6 569 t/a。1980年前所建的走马垅、白兰、五渡港、上迳、龙井和中村水库总体对桃江流域的减水量和减沙量分别为4.73×1010m3/a和6.91×107t/a。