朱熠明, 藍(lán)云龍, 周祖昊, 陳賽男, 蔡靜雅, 劉佳嘉
(1.中交(天津)生態(tài)環(huán)保設(shè)計研究院有限公司, 天津 300461;2.中國水利水電科學(xué)研究院流域水循環(huán)模擬與調(diào)控國家重點實驗室, 北京100038;3.華北水利水電大學(xué), 鄭州 450046; 4.黃河水利委員會西寧水文水資源勘測局, 西寧810008)
黃河上游來水來沙分布不均、水沙異源的問題對黃河流域開發(fā)治理影響較大[1-2]。針對黃河上游水沙變化規(guī)律,學(xué)者們從多角度進行了研究,侯素珍[3]、張世軍[4]、羅春紅[5]、田小靖[6]等采用趨勢性和突變性方法對黃河上游石嘴山、頭道拐等水文站的實測輸沙量進行規(guī)律和發(fā)展趨勢分析,站點斷面年輸沙量持續(xù)性銳減,于1968年和1986年輸沙量發(fā)生突變,呈顯著性減少趨勢;郭彥[7]、許文龍[8]、蘇曉慧[9]、王秀杰[10]等采用小波等周期性分析方法,發(fā)現(xiàn)黃河上游唐乃亥、石嘴山、頭道拐等站點月輸沙量的周期性在70—90年代間開始減弱并于2000年后消失,水沙序列的豐枯變化趨勢不一致,水沙變化具有多時間尺度特性,大中尺度振蕩嵌套較小尺度的周期振蕩;姚文藝[11]、申冠卿[12]等探究黃河上游龍羊峽、劉家峽水庫運行后引起的流量削峰作用與調(diào)節(jié)過程對河道輸沙過程、沖淤演變產(chǎn)生的影響;冉大川等[13]對黃河上游內(nèi)蒙古河段構(gòu)建降雨產(chǎn)流產(chǎn)沙力經(jīng)驗?zāi)P陀嬎憬邓退4胧υ摵佣萎a(chǎn)沙量的影響,降水是影響該河段產(chǎn)沙量減少的重要因素。
目前的研究多是針對黃河上游干流水文斷面或水庫開展,且多是利用實測輸沙資料分析序列的趨勢性、突變性、周期性等規(guī)律,而對于產(chǎn)沙量較大的支流水沙演變規(guī)律研究很少,且很少使用分布式水沙模型模擬探究水沙過程變化的機理。支流的氣候條件、水文條件、地形地貌條件、植被條件、土壤條件與黃河干流不同,其變化特征和時空分布呈現(xiàn)明顯的差異性。本文以黃河上游大夏河流域為研究區(qū),分析該流域的水沙演變規(guī)律、變異特征及驅(qū)動機制。
大夏河發(fā)源于青海省黃南藏族自治州同仁縣南部和甘肅省夏河縣西南部,河源海拔4 221 m,干流河道全長203 km,流域面積7 154 km2,于甘肅省臨夏縣塔張村注入劉家峽水庫,大夏河流域地理位置102°02′—103°23′E,34°51′—35°48′N。折橋水文站為大夏河流域出口控制站,控制流域面積6 967 km2。
研究區(qū)氣象數(shù)據(jù)來源于中國氣象數(shù)據(jù)網(wǎng);地形數(shù)據(jù)采用的高程數(shù)據(jù)為SRTM90;土壤類型數(shù)據(jù)來源于全國第二次土壤普查和《中國土種志》;土地利用信息包括經(jīng)國家相關(guān)部門批準(zhǔn)的1986年、1996年、2000年、2005年、2010年、2015年6個時段的1∶100 000土地利用圖;水土保持建設(shè)信息根據(jù)全流域各縣《水利統(tǒng)計年鑒》公布數(shù)據(jù),并與土地利用圖的疊加提取獲得;植被信息由遙感數(shù)據(jù)反演,主要包含葉面積指數(shù)、植被覆蓋度指數(shù),其中植被覆蓋度根據(jù)歸一化植被指數(shù)(NDVI)計算得出,NDVI來源于兩種數(shù)據(jù)源分別為:8 km精度的GIMMS AVHRR數(shù)據(jù)和1 km精度的MOD13A2數(shù)據(jù);葉面積指數(shù)來源于兩種數(shù)據(jù)源分別為:8 km精度的GlobMap LAI數(shù)據(jù)和1 km精度的MOD15A2數(shù)據(jù);輸沙量數(shù)據(jù)來源于黃河水利委員會提供的折橋水文站觀測數(shù)值。
水循環(huán)表征要素的時空演變規(guī)律主要采用統(tǒng)計學(xué)方法。使用Mann-Kendall非參數(shù)檢驗法[14],對大夏河流域多年實測輸沙序列要素的時間序列進行顯著性檢驗,定量反映大夏河流域變化趨勢和突變情況。
水循環(huán)影響要素的成因機理解析方法主要采用水文模型法及歸因分析法?;赪EP-L模型[15-17]構(gòu)建的WEP-SED模型[18-19]模擬流域的產(chǎn)輸沙過程。圖1為WEP-SED模型坡面及溝道的產(chǎn)輸沙過程,按照土壤侵蝕和泥沙輸移的一般步驟,逐流域進行“坡道—溝道—河道”三級流域輸沙計算。為驗證模型在研究區(qū)的模擬效率,選取Nash-Sutcliffe效率系數(shù)(NASH)、相對誤差(RE)及判定系數(shù)(R2)評價模型模擬效果[20]。一般來說,Re<20%,R2>0.5,Nash>0.5,認(rèn)為模型模擬是可靠的。
圖1 WEP-SED模型坡面及溝道產(chǎn)沙過程
對折橋站模擬水沙過程進行率定和驗證,表1為該站月尺度水沙過程模擬結(jié)果評價,圖2為該站水沙過程模擬與實測結(jié)果對比。經(jīng)過率定和驗證后的水沙過程模擬均滿足模擬評價的精度要求,為大夏河流域水沙過程變化分析和歸因分析提供條件。
通過多因素歸因分析方法進行歸因計算[21],根據(jù)各因素計算率公式,得到各因素分解量如下:
(1)
(2)
(3)
式中:δxi表示因素xi的貢獻量;ei,j是對應(yīng)第i個因素在第j個情景下的權(quán)重系數(shù)(變化期ei,j=1,基準(zhǔn)期ei,j=-1);Sj是第j個情景的模型模擬結(jié)果值;n為影響因素個數(shù);A為所有因素影響貢獻量之和;βi表示在考慮n個影響因素前提下,第i個因素占總變化的貢獻率。分析不同時期各因素對水沙過程影響,采用多因素歸因分析方法對水沙過程進行歸因分析,包括氣溫、降水和土地利用3個影響因素,情景設(shè)置見表2。
表1 1956-2016年折橋水文站月均水沙過程模擬結(jié)果評價
圖2 1956-2016年折橋水文站水沙過程模擬與實測結(jié)果對比
表2 歸因分析因素情景
圖3為大夏河流域?qū)崪y年徑流量和年輸沙量的變化情況,最大年徑流量為24.35億m3(1967年),最小年徑流量為3.85億m3(1991年),最大年徑流量是最小年徑流量的6.3倍;最大年輸沙量為1 391.48萬t(1967年),最小年輸沙量為44.71萬t(2009年),最大年輸沙量是最小年輸沙量的31.12倍;徑流和輸沙變化差異極大。根據(jù)Mann-Kendall趨勢分析法分析,年徑流量整體上呈減少趨勢,變化率為-0.85億m3/10 a,1994年后年徑流量下降趨勢顯著;年輸沙量整體上呈減少趨勢,變化率為-66.25萬t/10 a,1990年后年輸沙量下降趨勢顯著。
圖3 1956-2016年大夏河流域?qū)崪y徑流和輸沙量年際變化情況
圖4為大夏河流域?qū)崪y年徑流量和年輸沙量Mann-Kendall突變性檢驗結(jié)果,對檢驗結(jié)果進行分析,在給定顯著性水平為0.05(U0.05=1.96)的條件下,對于年徑流量,標(biāo)準(zhǔn)正態(tài)分布序列UF和反序列UB在1984—1985年出現(xiàn)交點,年徑流量發(fā)生突變,下降趨勢顯著;對于年輸沙量,UF和UB曲線在1987—1988年出現(xiàn)交點,該段時間年輸沙量發(fā)生突變,且下降趨勢顯著。綜上所述,年徑流量和年輸沙量均在1980s的中期發(fā)生突變。后續(xù)分析中,按照年輸沙量出現(xiàn)的拐點分段,以1956—1987年為基準(zhǔn)期,1988—2018年為變化期,對水沙過程變化規(guī)律進行分析。
圖4 1956-2016年大夏河流域?qū)崪y徑流和輸沙量Mann-Kendall檢驗
圖5為基準(zhǔn)期和變化期的實測水沙過程年內(nèi)變化情況。對于輸沙年內(nèi)分布過程,折橋站的年內(nèi)水沙過程曲線呈明顯“單峰型”,但河道徑流過程與輸沙過程不完全同步。
最大月徑流量出現(xiàn)在9月,而最大月輸沙率出現(xiàn)在8月;徑流量減幅最大的月份也出現(xiàn)在9月,而輸沙率減幅最大的月份也出現(xiàn)在8月。由此推斷,大夏河產(chǎn)流產(chǎn)沙機制存在較大的不同。
圖5 1956-2016年大夏河流域基準(zhǔn)期和變化期實測水沙過程年內(nèi)變化情況
通過計算,折橋站的多年平均輸沙量282萬t,主要集中在汛期的7—9月,占年輸沙量的79%。變化期(1988—2018年)相較于基準(zhǔn)期(1956—1987年)的徑流量和輸沙量均在減小,變化期的年徑流量較基準(zhǔn)期減少31%,輸沙量減小57%,輸沙量的減小幅度是徑流量的1.8倍;在輸沙量最大的7—9月,輸沙量月均減小幅度為56%,徑流量月均減小幅度為32%,與全年基本一致;在輸沙量最大的8月,輸沙量月均減小幅度為58%,徑流量月均減小幅度為33%,與全年基本一致;分析在1—6月和10—12月的輸沙量變化,相較于基準(zhǔn)期,變化期的徑流量和輸沙量同汛期一樣均在減少,該9個月的輸沙量減少52%,徑流量減少29%,與全年基本一致。
表3和圖6為大夏河流域在基準(zhǔn)期與變化期各影響因素變化情況。變化期內(nèi)的年均氣溫和月均氣溫較基準(zhǔn)期均呈增加趨勢,流域多年平均氣溫1.91℃,變化期較基準(zhǔn)期增加0.80℃;其中,7月、8月的月均氣溫為年內(nèi)最大溫度,多年平均月均氣溫分別為11.88℃,11.38℃,變化期內(nèi)月均氣溫較基準(zhǔn)期分別增加0.90℃,0.76℃;1月、12月為年內(nèi)最低月均氣溫,多年平均月均氣溫分別為-9.63℃,-8.04℃,變化期內(nèi)平均月均氣溫較基準(zhǔn)期分別增加0.96℃,0.96℃。
表3 1956-2018年大夏河流域基準(zhǔn)期與變化期各影響因素對比分析
圖6 1956-2018年大夏河流域氣溫、降水年內(nèi)分布情況
大夏河流域多年平均年降水量525.33 mm。分析降水、徑流、輸沙年內(nèi)過程,發(fā)現(xiàn)最大月降水量出現(xiàn)在7—8月,略早于最大月輸沙量(8月),早于最大月徑流量(9月)??梢酝茢?,輸沙量不僅與河道徑流量有關(guān)系,還跟坡面的產(chǎn)流、產(chǎn)沙量有關(guān)系,8月降水量最大,推斷是8月輸沙量最大的原因,在后文模型模擬中將對此進行分析。
變化期的多年平均降水量較基準(zhǔn)期變化幅度很小,兩個時期的年均降水變化幅度幾乎為0,7月和8月為年內(nèi)最大降水月,占年降水量的40%,分析這兩個月發(fā)現(xiàn),變化期8月的降水較基準(zhǔn)期減少11.32 mm,下降幅度10%;變化期7月降水較基準(zhǔn)期增加0.65 mm,增加幅度只有1%。
對于土地利用,大夏河流域林地、草地、裸地、農(nóng)田面積占比分別為30.6%,45.2%,17.8%,5.3%,余下的1.1%為水域和城鎮(zhèn)地表面積,所占比例較小。變化期相較基準(zhǔn)期,林地減少面積與草地增加面積的變化幅度較大且變化量相當(dāng),裸地、農(nóng)田面積變化總量較小,見表3。
3.3.1 水沙年內(nèi)分布規(guī)律成因分析 流域內(nèi)蒸散發(fā)量的變化受氣溫、降水等因素影響,在前文3.2節(jié)分析得到,7—8月的多年平均月均氣溫為年內(nèi)最大月均氣溫,多年平均月降水量在7月出現(xiàn)峰值,通過對大夏河流域的蒸散發(fā)量進行模擬,分析大夏河流域在年內(nèi)蒸散發(fā)量變化過程,大夏河流域年均氣溫蒸散發(fā)量474.06 mm,月蒸散發(fā)量在年內(nèi)隨著溫度的變化呈明顯的先增后減的趨勢,在7月達到年內(nèi)峰值,其中,5—9月的蒸散發(fā)量分別為67.7 mm,80.5 mm,92.7 mm,87.5 mm,60.9 mm,共計約占全年總蒸散發(fā)量的82%。
通過分析,7月之前的降水量較小,7月降水量較6月增加31.0%,達到年內(nèi)峰值。由于7月前的土壤較為干燥,含水率較低,土壤吸水率較大,又因為7月的蒸散發(fā)量較6月增加13.2%達到年內(nèi)峰值,在土壤吸水和蒸散發(fā)大的雙重影響下,徑流量在7月并未達到年內(nèi)峰值,但較6月增幅35.5%;由于7月較大的降水強度,坡面上的泥沙顆粒受降水形成的坡面水流作用,較容易發(fā)生泥沙起動,表層土壤中起動功率較低的泥沙顆粒較容易發(fā)生遷移,可較快完成產(chǎn)匯沙過程,導(dǎo)致7月的月輸沙量增加較多,輸沙量較6月的增加幅度為68.4%。
進入8月后,8月蒸散發(fā)量較7月減少5.9%,降水減少2.5%,蒸散發(fā)的減少幅度高于降水減少幅度,導(dǎo)致8月的徑流量較7月增加23.0%;由于降水過程的持續(xù),8月的土壤含水率較7月高,泥沙顆粒在水流作用下更容易發(fā)生遷移,故年內(nèi)月輸沙量的峰值出現(xiàn)在8月,輸沙量較7月增加42.0%。9月的蒸散發(fā)量較8月減少43.7%,降水的幅度減少為32.2%,蒸發(fā)的減少幅度大于降水,因為9月的實際降水量依然較大,且土壤含水率較7月高,徑流在9月出現(xiàn)年內(nèi)峰值;隨著9月降水量的減少,少部分泥沙起動功率較高的顆粒較難完成泥沙顆粒的遷移,故9月的月輸沙量較8月下降明顯,下降幅度為124.3%。10月后,隨著降水、蒸發(fā)的快速下降,徑流量和輸沙量下降明顯,減少幅度分別為24.0%和124.3%。
3.3.2 各項因子貢獻分析 圖7為大夏河流域徑流過程在氣溫、降水和土地利用等因素影響下的水沙過程年內(nèi)影響量變化情況,表4為該流域在各因素影響下的水沙過程歸因分析。通過分析,降水因素對大夏河流域減水、減沙的影響最大。前文中分析了降水量下降對徑流和輸沙過程的影響方式,尤其在降水量較大的7—9月,徑流量和輸沙量下降明顯。綜合來看,降水對減水、減沙的貢獻率分別為-77.0%和-90.0%。
圖7 各因素對1956-2018年大夏河流域水沙過程的影響
表4 1956-2018年大夏河流域產(chǎn)輸沙減少歸因分析
分析氣溫對徑流和輸沙過程的影響,變化期內(nèi)5—10月的月均氣溫較基準(zhǔn)期增加,該時段的氣溫增加導(dǎo)致流域內(nèi)蒸散發(fā)加大,因此變化期的河道水量較基準(zhǔn)期減少,河道水流挾沙力下降,從而輸沙量下降;其他月的氣溫對徑流、輸沙的影響較小,對年內(nèi)變化總體影響不大。綜合來看,氣溫對水沙過程的影響以減水、減沙為主,貢獻率分別為-29.7%和-16.5%為主。
分析土地利用因素對水沙過程影響,大夏河流域內(nèi)土地利用主要以林地、草地和裸地為主,農(nóng)田面積較小,由于受人類活動的影響,變化期較基準(zhǔn)期的林地面積減少量與草地面積增加量變化相當(dāng),由于草地的保水固沙能力較林地弱,尤其是植被生長旺盛的5—10月,植被的面積變化導(dǎo)致該區(qū)域的綜合保水固沙能力減弱,從而土地利用變化對徑流量和輸沙量的影響呈微弱的增水、增沙作用,貢獻率分別為6.7%,6.6%。
本文采用統(tǒng)計學(xué)方法分析了大夏河流域?qū)崪y輸沙量的變化特征,采用WEP-SED模型模擬該流域水沙過程并進行歸因分析,探究氣溫、降水和土地利用變化對水沙過程的影響,得到如下結(jié)論:
(1) 通過Mann-Kendall趨勢分析,大夏河流域1957—2011年的年徑流量呈減少趨勢,年平均下降速率為0.85億m3/10 a;多年輸沙量呈減少趨勢,平均下降速率約為66.3萬t/10 a。采用Mann-Kendall突變檢驗,大夏河流域的年徑流量在1984—1985年發(fā)生突變后呈顯著下降趨勢,年輸沙量在1987—1988年發(fā)生突變后顯著下降趨勢,輸沙量的減小幅度是徑流量的1.8倍。
(2) 分析大夏河流域氣溫、降水和土地利用的在基準(zhǔn)期和變化期的變化情況。變化期的年內(nèi)氣溫較基準(zhǔn)期增加幅度為0.35~1.28℃,增加幅度較大;兩時期的降水年均變化為幾乎為0,其中,在降水量最大的7月,變化期較基準(zhǔn)期的變化量基本不變,但降水量次大的8月降水減少較多,減少了10%,其余月的降水量均有所增加;對于土地利用,林地面積減少量與草地面積增加量相當(dāng),裸地、農(nóng)田面積變化總量較小。
(3) 采用多因素分析法對大夏河流域水沙過程歸因分析,發(fā)現(xiàn)降水對折橋站徑流量和輸沙量減少的貢獻最大,貢獻率分別為-77.0%和-90.0%;氣溫對徑流和輸沙量變化的貢獻率分別為-29.7%和-16.5%;土地利用變化對徑流和輸沙量的貢獻率分別為6.7%和6.6%,對水沙產(chǎn)沙起到微弱的增加作用。
雖然本文使用的WEP-SED模型在大夏河流域具有較好的適用性,但由于對水沙過程認(rèn)識的局限性,模型模擬結(jié)果存在一定的誤差。模型中所使用的氣象數(shù)據(jù)、地形數(shù)據(jù)、土地利用類型、植被數(shù)據(jù)等由于測量和統(tǒng)計誤差,對模型模擬結(jié)果的準(zhǔn)確性和可靠性也有一定影響。下一步研究中,將側(cè)重于上述不確定性問題的研究以提高量化結(jié)果。
致謝:感謝黃河水文水資源科學(xué)研究院張學(xué)成教高和李東教高、西安理工大學(xué)的李鵬教授、西北農(nóng)林科技大學(xué)的趙廣舉教授、黃河水利科學(xué)研究院的夏潤亮教高、中國水利水電科學(xué)研究院的張曉明教高對本研究的大力支持。