康苗業(yè),肖偉華,魯 帆,王義成,侯保燈,張雪蕾
(中國水利水電科學研究院 流域水循環(huán)模擬與調控國家重點實驗室,北京 100038)
河川徑流在氣候變化與人類活動的影響下發(fā)生了不同程度的變化,分析徑流演變規(guī)律對流域水資源的開發(fā)利用具有重要意義[1]。我國北方地區(qū)水資源匱乏,地下水超采嚴重,制約了經(jīng)濟社會的發(fā)展[2]。黃河是我國華北地區(qū)的主要水源,其徑流演變規(guī)律會對流域水資源的利用以及中下游水利工程的建設產(chǎn)生重要影響。目前已有相關學者對黃河流域徑流演變規(guī)律進行了大量研究,何毅[3]對黃河花園口水文站的徑流量進行M-K突變檢驗,發(fā)現(xiàn)徑流量在1986年發(fā)生突變;劉超等[4]對黃河流域上游4個主要干支流進行假設統(tǒng)計檢驗,結果表明這4個干支流的徑流量均有不同程度的減小趨勢;楊志峰等[5]采用M-K突變檢驗法對黃河三花(三門峽—花園口)區(qū)間的徑流量進行檢驗,結果表明該區(qū)間徑流量主要在1979年、1982—1983年發(fā)生突變。
河北省缺水問題突出,尤其是白洋淀曾出現(xiàn)過多次干淀現(xiàn)象,對社會經(jīng)濟發(fā)展與生態(tài)環(huán)境安全產(chǎn)生不利影響,因此對白洋淀進行生態(tài)補水迫在眉睫。引黃入冀補淀工程是大型引調水工程,該工程自河南濮陽渠村引黃口引水,渠線穿越海河向沿線受水區(qū)及白洋淀輸水,這是目前對白洋淀進行長效補水的最為現(xiàn)實措施之一[6]。但在引水初期,受黃河主河槽擺動影響,河勢南移遠離引水口且引水口淤塞,導致渠首閘引水流量僅4 m3/s,無法滿足引水需求,因此在引黃入冀補淀工程運行過程中保障引黃水量的有效供給十分重要?;▓@口水文站位于渠村引黃口上游,是距離該引黃口最近的水文站,同時雄安新區(qū)的建設與發(fā)展對該引黃口的水源供給提出了更高要求,因此本文以1956—2017年花園口水文站年徑流量與引水期徑流量[7]為基礎數(shù)據(jù),檢驗年徑流量、引水期徑流量的變化趨勢與顯著性,分析徑流年內(nèi)分配不均勻性,確定年徑流量、引水期徑流量的突變年份,進而總結黃河花園口水文站的徑流演變規(guī)律,以期為引黃可用徑流量研究提供一定的理論支持。
花園口水文站設立于1938年7月,位于鄭州市北郊,東經(jīng)113°40′15″,北緯34°54′9″,距離河南濮陽渠村引黃口約134 km,是黃河干流的重要控制站以及黃河下游防汛的標準站。該水文站集水面積73萬km2,占黃河流域總面積的97%,是國家級重要水文站[8],在氣候變化和人類活動的影響下,分析長時間序列的實測徑流數(shù)據(jù)有助于精確獲取水文過程的統(tǒng)計規(guī)律[9]。本文以黃河花園口水文站以上流域為研究區(qū)域,選取1956—2017年花園口水文站的實測徑流數(shù)據(jù)對徑流演變規(guī)律進行具體分析。圖1為黃河流域、花園口水文站與引黃入冀補淀工程位置的簡單示意。
圖1 黃河流域、花園口水文站與引黃入冀補淀工程位置示意
本文利用5 a滑動平均法、M-K趨勢檢驗法與Sen’s斜率估計法分析徑流年際變化趨勢與顯著性,采用不均勻系數(shù)Cv分析徑流年內(nèi)分配不均勻性,通過M-K突變檢驗法確定徑流量發(fā)生突變的年份。
(1)5 a滑動平均法。在對長時間序列的徑流進行研究時,很難從序列徑流本身發(fā)現(xiàn)其變化趨勢,采用5 a滑動平均法可以對序列徑流進行光滑處理,從而消除一些偶然變動因素的影響,公式如下[10]:
式中:yj為第j個滑動平均值;c為滑動長度;xi+j-1為時間序列的第i+j-1個數(shù)據(jù)值。
(2)M-K趨勢檢驗法與Sen’s斜率估計法。M-K趨勢檢驗法僅可以對時間序列變化趨勢的顯著性進行檢驗[11],而Sen’s斜率估計法僅可以計算得到時間序列的變化幅度即時間序列的斜率[12],因此將2種方法結合使用,可以更好地分析長時間序列徑流的年際變化。M-K趨勢檢驗法是一種非參數(shù)統(tǒng)計檢驗方法[13],非參數(shù)統(tǒng)計檢驗又稱無分布檢驗,其優(yōu)點是樣本不必遵循某一特定的分布規(guī)律,受異常值的干擾較小,計算簡便[14]。對于數(shù)據(jù)序列{x1,x2,…,xn},n>10,服從標準正態(tài)分布的時間序列統(tǒng)計量Z可由下式計算:
式中:S為趨勢檢驗的統(tǒng)計量;Var(S)為統(tǒng)計量的方差;xi為時間序列的第i個數(shù)據(jù)值;xj為時間序列的第j個數(shù)據(jù)值;t為時間長度。
Z>0表示序列呈上升趨勢,Z<0表示序列呈下降趨勢。對于顯著性水平α,當Z>Z1-α/2或者Z<Zα/2時,表明趨勢通過顯著性檢驗,即變化趨勢顯著。
采用Sen’s斜率估計法計算時間序列的斜率β,可以降低或避免數(shù)據(jù)異?;蛉笔τ嬎憬Y果的影響,公式為
式中:Median()為中值函數(shù)。
當β>0時,序列呈上升趨勢;當β=0時,序列變化趨勢不明顯;當β<0時,序列呈下降趨勢。
不均勻系數(shù)Cv可以定性表示徑流分配的不均勻程度[15-16],公式為
式中:σ為均方差系數(shù);ˉr為月平均徑流量,m3;m為月份;rm為第m個月的徑流量,m3。
M-K突變檢驗法可用于檢驗長時間序列發(fā)生突變的時間、突變次數(shù)及其顯著性[17-21]。對時間序列構成秩序列Sk:
式中:秩序列Sk為第i時間序列值大于j時間序列值個數(shù)的累計數(shù)。
當x1、x2、…、xn相互獨立且連續(xù)分布時,定義統(tǒng)計量[18]:
式中:U Fk為M-K檢驗統(tǒng)計量,服從正態(tài)分布,U F1=0;E(Sk)、Va r(Sk)分別為Sk的均值和方差。
設Uα為顯著性水平為α的統(tǒng)計量,其中U0.05=±1.96、U0.01=±2.57,若|U Fk|>|Uα|,說明時間序列的變化趨勢顯著,并且|U Fk|越大、序列變化趨勢越顯著。當U Fk>0時,表示序列呈上升趨勢;當U Fk<0時,表示序列呈下降趨勢。將時間序列逆序并重復M-K檢驗統(tǒng)計量的計算過程,得到逆序列檢驗統(tǒng)計量U Bk=-U Fk,當統(tǒng)計量序列曲線UFk與統(tǒng)計量逆序列曲線UBk出現(xiàn)交點且交點在顯著性水平臨界線之間時,交點所對應年份即為發(fā)生突變時間。
1956—2017年花園口水文站長時間序列的實測年徑流量的5 a滑動平均曲線、實測年徑流量變化曲線及其趨勢線見圖2。由實測年徑流量變化曲線可以看出,花園口水文站的年徑流量總體呈減小趨勢,減小率為51.16億m3/10 a,變化幅度較大,最大年徑流量出現(xiàn)在1964年,為861.41億m3;最小年徑流量出現(xiàn)在1997年,為142.57億m3;多年平均徑流量為351.43億m3,年徑流量的極值比為6.04,表明年徑流量的豐枯變化比較明顯。
圖2 1956—2017年花園口水文站實測年徑流量變化與5 a滑動平均過程
1956—2017年花園口水文站引水期實測徑流量的5 a滑動平均曲線、引水期實測徑流量變化曲線及其趨勢線見圖3。由引水期實測徑流量變化曲線可以看出,花園口水文站引水期徑流量整體呈減小趨勢,減小率為6.99億m3/10 a,最大徑流量出現(xiàn)在1961年,為149.48億m3;最小徑流量出現(xiàn)在1960年,為24.36億m3;引水期多年平均徑流量為69.04億m3,引水期徑流量的極值比為6.26。與年徑流量的減小率相比,引水期的變化程度較小,表明該時期徑流比較穩(wěn)定。
圖3 1956—2017年花園口水文站引水期實測徑流量變化與5 a滑動平均過程
采用M-K趨勢檢驗法與Sen’s斜率估計法對年徑流量與引水期徑流量進行了趨勢分析和顯著性檢驗。對于年徑流量,計算得出時間序列統(tǒng)計量Z=-10.281 0,|Z|>Z0.01=2.570 0,Sen’s斜率估計值β=-0.172 4<0,說明年徑流量減小趨勢非常顯著;對于引水期徑流量,計算得出時間序列統(tǒng)計量Z=-8.206 8,|Z|>Z0.01=2.570 0,Sen’s斜率估計值β=-0.006 2<0,同樣說明引水期徑流量減小趨勢非常顯著。但與年徑流量相比,引水期徑流量的變化幅度較小,表明引水期徑流更加穩(wěn)定。
1956—2017年花園口水文站月平均徑流量的年內(nèi)分布見圖4,可以看出,7—10月連續(xù)4個月的徑流量較大,范圍為41.72億~51.99億m3,徑流量之和為184.80億m3,占全年徑流量的52.5%,說明年內(nèi)徑流分配比較集中,原因為該時段處于黃河汛期,上游暴雨增多,導致徑流量增大。除此以外,與1、2月份相比,3—6月的徑流量有稍微增大趨勢,這是由于該時段黃河上游冰雪融化,導致徑流量增大。
引黃入冀補淀工程的引水時段在11月至次年2月,由圖4可以看出,這4個月的月平均徑流量在12.67億~25.30億m3之間,相對于其他月份來說徑流量較小,但是非汛期徑流主要由地下水補給、上游水庫調節(jié),含沙量較小,因此該時段內(nèi)的徑流大部分可用于向白洋淀供水。
圖4 1956—2017年花園口水文站月平均徑流量
1956—2017年花園口水文站長時間序列的年徑流量不均勻系數(shù)的計算結果見圖5,可以看出,最大不均勻系數(shù)出現(xiàn)在1992年,為1.190,表明1992年徑流分配最不均勻;最小不均勻系數(shù)出現(xiàn)在1990年,為0.290,表明1990年徑流分配最為均勻。
圖5 1956—2017年花園口水文站年徑流量不均勻系數(shù)
1956—2017年花園口水文站引水期徑流量的不均勻系數(shù)計算結果見圖6,可以看出,最大不均勻系數(shù)出現(xiàn)在1961年,為0.543,表明1961年引水期徑流分配最不均勻;最小不均勻系數(shù)出現(xiàn)在1995年,為0.038,表明1995年引水期徑流分配最為均勻。與年徑流量的不均勻系數(shù)相比,引水期徑流量的不均勻系數(shù)較小,均在0.600以下,說明引水期徑流分配更加均勻。
圖6 1956—2017年花園口水文站引水期徑流量不均勻系數(shù)
利用M-K突變檢驗法分析了顯著性水平α=0.05時1956—2017年花園口水文站長時間序列年徑流量的變化趨勢與突變年份,如圖7所示。1956—1970年,UFk值在0上下波動,說明該時期徑流量變化趨勢不顯著;1970—1991年,UFk<0,且介于0.00與0.05顯著性水平臨界線之間,說明該時期徑流量呈減小趨勢,但是不顯著;1991年以后,UFk值超過臨界線且直線下降,說明1991年以后徑流量呈顯著減小趨勢。UFk、UBk曲線在1986年出現(xiàn)交點,說明年徑流量在1986年發(fā)生突變。
圖7 花園口水文站年徑流量M-K突變檢驗
花園口水文站引水期徑流量的M-K突變檢驗結果見圖8,可以看出,該時段的徑流量變化趨勢與年徑流量的變化趨勢相似,同樣是在1986年發(fā)生突變。
圖8 花園口水文站引水期徑流量M-K突變檢驗
20世紀80年代以前,人類活動對黃河流域的影響相對較小,下墊面變化較緩慢,因此徑流量變化不大;80年代后期開始,人類活動及下墊面的影響逐步體現(xiàn),如退耕還林草、梯田建設等引起綠水增加、藍水減少,同時流域水資源開發(fā)強度加大,消耗量增加,河川徑流量減少,導致黃河徑流量呈現(xiàn)出顯著的衰減趨勢。本文中花園口水文站年徑流量與引水期徑流量均在80年代前期呈現(xiàn)出微弱的減小趨勢,在80年代后期發(fā)生突變,進入90年代后呈現(xiàn)出顯著的減小趨勢,該結論與黃河流域近年來的變化規(guī)律相符合。
通過對1956—2017年花園口水文站長時間序列的年徑流量與引水期徑流量進行分析,可以得出以下結論:
(1)通過62 a的徑流量數(shù)據(jù)來看,年徑流量在142.57億~861.41億m3之間,總體呈減小趨勢,減小率為51.16億m3/10 a;引水期徑流量在24.36億~149.48億m3之間,同樣呈現(xiàn)出減小趨勢,減小率為6.99億m3/10 a。兩者相比可以看出,引水期徑流量的變化率較小,相對穩(wěn)定。
(2)通過對62 a的徑流量進行趨勢分析與顯著性檢驗,可知在此時間序列上年徑流量與引水期徑流量均呈現(xiàn)出顯著減小的趨勢,并且引水期徑流量的減小幅度與顯著性均小于年徑流量。
(3)年內(nèi)徑流分配不均勻,主要集中在7—10月,這4個月的月平均徑流量之和占全年徑流量的52.5%;年徑流量的不均勻系數(shù)在0.290~1.190之間,引水期徑流量的不均勻系數(shù)在0.038~0.543之間,因此引水期的徑流變化比較緩和。
(4)通過M-K突變檢驗法對62 a序列徑流進行檢驗,發(fā)現(xiàn)年徑流量與引水期徑流量均在1991年以后呈現(xiàn)出顯著減小的趨勢,并且均在1986年發(fā)生突變。