田昊瑋， 陳伏龍， 龍愛華， 劉 靜， 海 洋

（1.石河子大學水利建筑工程學院，新疆 石河子 832000；2.中國水利水電科學研究院水資源所，北京 100038；3.天津大學水利工程仿真與安全國家重點實驗室，天津 300072）

艾比湖流域（79°53′～83°53′E，44°02′～45°23′N）位于新疆西北部，是準噶爾盆地西部最低的匯水中心。自20 世紀50 年代以來，由于大規(guī)模的水土開發(fā)，灌區(qū)人口、灌溉面積和引水量大幅度增加，原來有水注入艾比湖的大部分河流先后斷流，導致入湖水量急劇減少，流域徑流不僅要承擔人們生產(chǎn)生活供水任務，還要為維持當?shù)亓己蒙程峁┥鷳B(tài)用水，2019 年生態(tài)環(huán)境用水量僅占總用水量的3.1%，當?shù)厣鐣?jīng)濟-生態(tài)環(huán)境用水矛盾日益增加。博爾塔拉河（簡稱博河）是艾比湖的主要入湖水源，近40 a補給水量達入湖總量的70%，博河水量變化直接關乎艾比湖面積增減及湖周生態(tài)狀況，對艾比湖生態(tài)修復保護發(fā)揮基礎性的支撐作用。博河第一次冰川編目條數(shù)為255 條、面積216.91 km2；第二次冰川編目條數(shù)已經(jīng)退減到243條，面積消融至113.75 km2。兩次冰川編目數(shù)據(jù)相比，博河冰川消失12 條，面積縮減103.15 km2。隨著全球進一步增暖，氣候模式預估全球水循環(huán)進一步增強，水循環(huán)季節(jié)差異增大［1］。博河上游地區(qū)位于高山區(qū)，在氣候變化背景下，山區(qū)冰川消融加劇，導致極端水文事件的發(fā)生概率增多，河流水文過程將會變得更為復雜［2］。因此，開展氣候變化背景下博河徑流變化預測，對預防氣候變化帶來的潛在水資源風險，保障當?shù)厣鐣?jīng)濟發(fā)展與生態(tài)文明建設顯得十分迫切。

目前，針對徑流的研究主要包括尋找數(shù)據(jù)間的最優(yōu)數(shù)學關系為目標的數(shù)據(jù)驅(qū)動模型［3］和以水文學概念為基礎考慮物理過程的水文模型［4］。與數(shù)據(jù)驅(qū)動模型相比，水文模型能真實、有效地模擬流域徑流形成過程，客觀地反映氣候或下墊面變化對流域水循環(huán)的影響程度［5］，是研究水文循環(huán)機理和過程、解決復雜水文問題的有效基礎工具［6］。在眾多水文模型中，SWAT模型具有很高的普適性，在各流域的模擬應用中［7-11］，均取得了較好的效果；在對氣候變化的研究過程中，SWAT 模型成為國內(nèi)外研究人員開展工作最常用的模型之一，如Shukla等［12］、Cao等［13］、Zhao 等［14］分別運用SWAT 模型在Satluj 流域、呼倫湖流域、白山盆地等開展了氣候變化情景對徑流的影響研究。同時，SWAT 模型和其他模塊具有耦合性，楊明智等［15］基于SWAT 模型嵌入水資源配置模塊，實現(xiàn)了“自然-社會”二元水循環(huán)與水資源實時配置一體化模擬功能；Yin等［16］將SWAT模型嵌入冰川模塊，解決了SWAT 模型在高寒山區(qū)模擬精度不高的問題，并定量估算了葉爾羌河流域冰川對徑流的貢獻?？梢钥闯觯琒WAT 模型不僅在廣泛應用于各流域和氣候變化的研究中，還能針對自身薄弱點，因地制宜地耦合相關模塊，提高模型模擬精度。

在博河以往的研究中，多側(cè)重從融雪描述徑流機理［17-18］和水質(zhì)［19］、水化學［20］、降水［21］、地下地表水同位素變化特征［22］等方面的研究，很少關注氣候變化對博河流域水文過程的影響。本研究以博河上游為研究區(qū)，將冰川模塊嵌入SWAT模型中，構(gòu)建能夠刻畫冰川徑流的冰川增強型SWAT 模型，并耦合CMIP5 氣候變化模式中的不同情景，預測未來徑流增減情況，揭示引起徑流變化的主要驅(qū)動因子，以期為變化環(huán)境下的流域水資源管理措施提供輔助支撐。

1 研究區(qū)概況

博河流域（79°53′～83°53′E，44°02′～45°23′N）地處西天山北坡的阿拉套山，全長252 km，近60 a 平均徑流量為5.1×108m3，經(jīng)溫泉縣和博樂市后注入艾比湖，并以溫泉縣和博樂市為界，將博河分為上、中、下游，其中冰川面積為110.3 km2，冰川補給量約為1.05×108m3（圖1）。博河溫泉站以上集水面積為2206 km2，全長97 km。年均氣溫為6.6～7.8 ℃，其中5—9 月的氣溫達到10 ℃以上，最高是7 月；多年平均降水量為288.63 mm，主要集中在4—9 月，6 個月內(nèi)降水量約占全年的80.7%，其中降水量最大的月份是7月?？梢?，夏季是氣溫增高、降水增多的集中時段，最高月份均出現(xiàn)在7月，雨熱同期。全年徑流呈現(xiàn)中間高兩頭低的單峰分布，年內(nèi)徑流主要集中發(fā)生在6—8月，3個月內(nèi)徑流量約占全年的40.62%，其他9 個月份徑流量相對平均，均為0.2×108m3左右。

圖1 艾比湖水系分布Fig.1 Distribution of Ebinur Lake river system

2 數(shù)據(jù)與方法

2.1 研究方法

2.1.1 SWAT模型SWAT 模型是在美國農(nóng)業(yè)部（USDA）的農(nóng)業(yè)研究中心Jeff Arnold 博士的主持下，于1994 年將SWRRB 模型與ROTO 模型整合到一起，開發(fā)出了SWAT 模型。該模型經(jīng)過幾十年的發(fā)展，已經(jīng)被廣泛應用于徑流模擬、泥沙、點源污染等方面［23-25］。在SWAT模擬各計算單元的水循環(huán)過程時，水量平衡法則是根本驅(qū)動力。模型首先根據(jù)數(shù)字高程模型（DEM）數(shù)據(jù)、土壤類型、土地利用類型及坡度，細分流域反映不同作物和土壤類型的蒸散發(fā)差異，然后模擬各水文響應單元的徑流量，最后匯總得到流域的總徑流量，物理意義明確。

式中：SWt為土壤最終含水量（mm）；SW0為第i天的土壤初始含水量（mm）；t為時間（d）；Rday為第i天的降水量（mm）；Qsurf為第i天的地表徑流量（mm）；Ea為第i天的蒸散發(fā)量（mm）；Wseep為第i天從土壤剖面進入包氣帶的水量（mm）；Qgw為第i天回歸流的水量（mm）。

2.1.2 冰川模塊算法SWAT 模型具有較強的徑流演化與溶質(zhì)遷移模擬功能，但對干旱區(qū)內(nèi)陸河流域徑流的重要組分——冰川融水模擬考慮不足［26］。Yin 等［16］根據(jù)冰川融化因子和氣溫的關系，考慮太陽輻射和地形等因素的影響，修正了溫度指數(shù)模型，對冰川水文過程的模擬更具物理意義。

式中：M為冰川日融化量（mm）；FM為冰溫融化因子；Rice為冰輻射融化因子；Ipot為潛在太陽輻射；T為日平均氣溫（℃）；Tmlt,ice為冰川融化溫度閾值（℃）。有關應用表明［27-28］，這一改進對提高干旱內(nèi)陸河徑流模擬具有較好的完善作用。

2.1.3 模型評價本文選定4 個評價指標作為評判模擬結(jié)果好壞的標準，分別是納什系數(shù)（NSE）、均方根誤差與實測值標準差的比值（RSR）、偏差百分比（PBIAS）和決定性系數(shù)（R2），其中R2表示模擬值與實測值的相關程度，取值范圍為0～1，越接近1 表明相關性越好，具體評價標準［29］見表1。各指標計算公式如下：

表1 模型評價標準Tab.1 Model evaluation criteria

2.2 數(shù)據(jù)來源及模型建立

2.2.1 數(shù)據(jù)來源與處理驅(qū)動SWAT 模型的數(shù)據(jù)有DEM 數(shù)據(jù)、土壤數(shù)據(jù)、土地利用數(shù)據(jù)、氣象數(shù)據(jù)等。DEM數(shù)據(jù)來自地理空間數(shù)據(jù)云平臺，空間分辨率為90 m（http://www.gscloud.cn/）。土壤數(shù)據(jù)來自中國科學院南京土壤所發(fā)布的1:1000000 土壤數(shù)據(jù)庫，其中土壤容重（SOL_BD）、飽和導水率（SOL_K）和田間持水量（SOL_AWC）等參數(shù)通過SPAW 模型計算得到，然后將數(shù)據(jù)重分類成12種土壤類型。土地利用數(shù)據(jù)來自地理空間數(shù)據(jù)云2018年30 m 分辨率土地利用數(shù)據(jù)（http://www.gscloud.cn/），經(jīng)過重分類以后共有10種土地利用類型。

水文數(shù)據(jù)來源博爾塔拉水文探勘局，其中包括溫泉站1972—2018 年逐月徑流數(shù)據(jù)。氣象數(shù)據(jù)為流域內(nèi)國家氣象站點1968—2018年的監(jiān)測數(shù)據(jù)，包括降水、平均氣溫、最高氣溫、最低氣溫。用于預測的氣象數(shù)據(jù)來自第五次國際耦合模式比較計劃（CMIP5）的預測數(shù)據(jù)，從中選取8個GCM（Global climate model）數(shù)據(jù)（表2），包括中低濃度排放情景（RCP4.5）和高濃度排放情景（RCP8.5）2種情景模式下2020—2050 年的最高氣溫、最低氣溫、降水數(shù)據(jù)。對這些數(shù)據(jù)首先運用雙線性插值法將各模式的分辨率統(tǒng)一成0.5°×0.5°，然后采用RoMBC（Robust multivariate bias correction）方法對模式輸出數(shù)據(jù)進行偏差校正，該方法能有效降低數(shù)據(jù)的系統(tǒng)誤差［30］，最后根據(jù)多模式等權(quán)重的集合平均方法［31］輸出2種情景模式下降水、氣溫數(shù)據(jù)，分析2種情景模式下對研究區(qū)徑流變化的影響。

表2 CMIP5全球氣候模式基本信息Tab.2 Basic information of CMIP5 global climate models

2.2.2 博河上游SWAT模型構(gòu)建本文以0.8 km2作為集水閾值，將研究區(qū)劃分為19 個子流域，1095 個水文響應單元，在此基礎上輸入歷史氣象數(shù)據(jù)構(gòu)建出SWAT 模型?；谒恼緦崪y月徑流資料，通過調(diào)整參數(shù)獲得歷史最優(yōu)模擬情況，運用模型評價指標評定模擬結(jié)果，最后將CMIP5多模式集合未來情景數(shù)據(jù)帶入率定好的模型中，對未來博河上游徑流進行預測和變化分析。圖2為在構(gòu)建SWAT模型模擬方案時的DEM、子流域劃分、土地利用和土壤類型的空間分布圖。

圖2 數(shù)字高程模型（DEM）、子流域劃分、土地利用類型和土壤類型的空間分布Fig.2 Spatial distributions of DEM,subwatershed division,land use types and soil types

3 結(jié)果與分析

3.1 模型參數(shù)化與效果評估

官方發(fā)布的SWAT模型對冰川消融過程模擬欠缺，因此在模型構(gòu)建時加入冰川模塊［16］。采用swatcup 中的SUFI-2 算法對參數(shù)進行敏感性分析，結(jié)合敏感性分析結(jié)果、各參數(shù)的物理意義以及研究區(qū)特殊的水文過程，采用試錯法逐步調(diào)整關鍵參數(shù)，得到最終率定結(jié)果（表3）。

表3 參數(shù)率定結(jié)果Tab.3 Parameter calibration results

以1968—1971、1972—1996、1997—2018 年分別作為模型的預熱期、率定期、驗證期，對研究區(qū)月徑流進行模擬（圖3、表4）。校準期和驗證期NSE分別為0.81 和0.83，NSE＞0.75，模擬結(jié)果為優(yōu)；校準期和驗證期RSR 分別為0.44 和0.41，RSR＜0.5，模擬結(jié)果為優(yōu)；校準期和驗證期PBIAS 分別為0.26%和-6.84%，|PBIAS|＜10%，模擬結(jié)果為優(yōu)；校準期和驗證期R2分別為0.82和0.86，很接近1，模擬值與實測值相關性很好，模型整體評價為優(yōu)。以上綜合表明，加入冰川模塊的SWAT 模型在博河上游的模擬結(jié)果較為可靠，模型參數(shù)設置合理，可用于該流域的徑流過程模擬和氣候變化的響應研究。

表4 博爾塔拉河上游流域模擬結(jié)果評價Tab.4 Evaluation on simulation results of upper watershed of Bortala River

圖3 1972—2018年溫泉站月平均流量模擬結(jié)果Fig.3 Simulation results of monthly average flow of Wenquan station during 1972—2018

3.2 水文氣象因子變化分析

根據(jù)博河上游溫泉站1968—2018 年實測水文資料及情景模式水文資料，分析其降水、氣溫變化趨勢（圖4）。2 種模式下平均氣溫分別為5.83 ℃和6.93 ℃，較歷史時期4.15 ℃都升溫1.50 ℃以上，升溫趨勢顯著，這也與IPCC第一工作組評估報告［32］預測結(jié)果相符。2 種氣候模式下未來降水量分別為328.26 mm 和302.14 mm，較歷史時期年均降水量288.63 mm 分別增長了13.73%（RCP4.5）和4.70%（RCP8.5），其增速分別為13.21 mm·(10a)-1（RCP4.5）和4.35 mm·(10a)-1（RCP8.5）。RCP4.5降水量的增速明顯大于RCP8.5，其原因可能是增溫導致蒸發(fā)量加劇，進而抵消了水循環(huán)加快使降水增加的積極影響［33］。此外，通過累積距平法發(fā)現(xiàn)，2種氣候模式的氣溫和降水到2050年并未出現(xiàn)明顯拐點，這表明增溫增濕還將繼續(xù)。

圖4 溫泉站水文數(shù)據(jù)及情景模式年際變化趨勢Fig.4 Inter annual variation trends of hydrological data and scenario models of Wenquan station

用2 種氣候模式的降水、氣溫數(shù)據(jù)驅(qū)動率定好的SWAT 模型，分別得到各模式下的模擬徑流。從年際變化趨勢（圖5）來看，在氣溫和降水的綜合影響下，RCP4.5和RCP8.5情景模式中，年均徑流量分別為4.09×108m3和4.36×108m3，較歷史時期年徑流量3.16×108m3分別增加了0.93×108m3和1.2×108m3，增速分別為0.31×108m3·(10a)-1和0.40×108m3·(10a)-1。冰川徑流對總徑流的貢獻率由歷史時期的27.61%提升到32.45%（RCP4.5）和36.99%（RCP8.5），較歷史時期分別增長了4.84%和9.38%，冰川加速消融導致冰川貢獻率占比提升是徑流量增加的主要原因。從年內(nèi)變化趨勢（圖6）來看，RCP4.5 冰川徑流主要是從4 月開始到9 月結(jié)束，而RCP8.5 冰川徑流主要是從3 月開始到9 月結(jié)束，這也說明升溫會導致冰川消融時間提前。2種模式下模擬的月內(nèi)徑流主要集中在7、8 月，7 月達到峰值，峰值流量分別為41.00 m3·s-1（RCP4.5）和50.65 m3·s-1（RCP8.5），這也與冰川徑流峰值重合，7、8 月冰川徑流占當月徑流量65%以上，占全年冰川徑流的60%。值得注意的是，冰川消融峰值與汛期重合，氣候變暖導致冰川加速消融會加劇水系統(tǒng)脆弱性和河川徑流的波動性，使汛期極端水文事件的發(fā)生概率多于歷史時期［2］。

圖5 實測徑流、模擬年徑流、冰川貢獻率年際變化趨勢Fig.5 Inter annual variation trends of measured runoff,simulated annual runoff and glacier contribution rates

圖6 模擬情景徑流、情景冰川徑流、冰川貢獻率年內(nèi)變化趨勢Fig.6 Annual change trend of simulated scenario runoff,scenario glacier runoff and glacier contribution rates

3.3 氣候變化對徑流的影響

將1972—2018 年實測數(shù)據(jù)和2020—2050 年2種氣候模式模擬數(shù)據(jù)的月徑流、月降水、月均溫數(shù)據(jù)按照季節(jié)春（3—5 月）、夏（6—8 月）、秋（9—11月）、冬（12月—翌年2月）劃分，分別采用Pearson相關系數(shù)進行雙側(cè)檢驗（表5、表6、表7）。按照季節(jié)對比分析可以發(fā)現(xiàn)，徑流在春季與氣溫的相關性隨著溫度的上升相關性也隨之上升，這也證實了氣候變暖會導致冰川消融時間提前；徑流在夏季與氣溫均通過0.01 顯著性檢驗，且相關性不斷提高，這也說明夏季增溫會導致冰川消融加速，使博河上游徑流量增加，這也印證了迪麗努爾·阿吉等［34］得到的結(jié)論。秋冬兩季的溫度基本在0 ℃以下，氣溫與徑流的相關性不高，氣溫與降水的相關性較高，秋冬兩季主要通過降水的形式促進冰川的積累。但秋季氣溫與降水的相關性隨著氣溫的升高，相關性也隨之降低，這改變了冰川的積累、消融規(guī)律，導致冰川面積進一步的縮減，冰川調(diào)節(jié)能力下降。

表5 歷史水文數(shù)據(jù)突變對比分析Tab.5 Comparison and analysis of abrupt changes in historical hydrological data

表6 RCP4.5情景模式徑流、氣溫、降水相關系數(shù)Tab.6 Correlation coefficients of runoff,temperature and precipitation in RCP4.5 scenario model

表7 RCP8.5情景模式徑流、氣溫、降水相關系數(shù)Tab.7 Correlation coefficients of runoff,temperature and precipitation in RCP8.5 scenario model

4 討論

博爾塔拉河上游屬于高寒山區(qū)，徑流主要以高山降水和冰川融雪為主。眾多學者在研究博河徑流時，多是考慮降水和融雪對徑流的影響［18,21］，而忽視了高寒山區(qū)冰川消融對徑流的影響，孟現(xiàn)勇等［17］采用SWAT+CMADS 模擬精博河發(fā)現(xiàn)溫泉站的模擬精度低于精河站，猜測是因為溫泉上游高寒山區(qū)冰川對溫泉站模擬結(jié)果影響較大。本文將SWAT模型和冰川模塊耦合，該模塊在眾多的以冰川徑流占比較高的徑流模擬中均有較好的適用性［11,16,27-28］，彌補了SWAT 模型冰川消融過程模擬欠缺的問題，能有效提高模型在高寒山區(qū)的模擬精度。

氣溫變化是導致冰川消融速率改變的主要因素。通過對未來的氣溫分析發(fā)現(xiàn)，未來氣溫平均升高1.5 ℃。有學者預估了在1.5 ℃升溫條件下，亞洲高山區(qū)到21世紀末冰川損失量超過一半以上［35］，這意味著因冰川消融速率加快帶來的徑流量增加是不可持續(xù)的，冰川調(diào)節(jié)能力的下降會增加徑流量變率［2］，給水資源管理者帶來新的挑戰(zhàn)。

本研究仍存在不足之處。研究區(qū)枯水期相較于豐水期擬合較差，導致枯水期模擬效果較差的原因是，博河自西向東流至溫泉時，經(jīng)過斷陷盆地、基地隆起等地質(zhì)結(jié)構(gòu)，導致地表水地下水頻繁轉(zhuǎn)換［36］，頻繁的地表水和地下水轉(zhuǎn)換，使得枯水期的徑流量呈現(xiàn)出無規(guī)則的波動，進而增加模擬枯水期徑流量時的難度。今后研究中考慮將SWAT 模型與MODFLOW 模型等地下水模型進行耦合，進一步提高模擬精度。冰川模塊中的冰川消融因子（冰輻射融化因子）屬于經(jīng)驗參數(shù)，模擬徑流時是以多年徑流過程線為參照物進行調(diào)參的，單純的徑流校正會因降水徑流模擬的不準確導致冰川消融因子的低估和高估，從而影響徑流組分的分割［37］。如何降低因調(diào)參帶來的模型模擬差異，從而做到精準化分析徑流組分是今后值得討論的問題之一。

5 結(jié)論

本文以博爾塔拉河源流區(qū)為靶區(qū)，構(gòu)建嵌入冰川模塊的SWAT模型，在CMIP5不同情景模式下，預測徑流變化趨勢，評估氣候變化對徑流的影響。主要結(jié)論如下：

（1）考慮冰川融化對徑流的影響，構(gòu)建嵌入冰川模塊SWAT 模型。結(jié)果表明，嵌入冰川模塊的SWAT 模 型，全 時 段 的NSE 為0.82，RSR 為0.42，RBIAS 為-3.22%，R2為0.84，模型評價結(jié)果為優(yōu)，證明本次模型的構(gòu)建在該流域有很好的適用性。

（2）在CMIP5的RCP4.5和RCP8.5情景模式中，未來平均氣溫分別為5.83 ℃和6.93 ℃，未來降水量分別為328.26 mm和302.14 mm，兩者較歷史時期均表現(xiàn)出增加的趨勢。用2 種氣象數(shù)據(jù)模擬徑流，得到的年均徑流量歷史時期年均徑流量分別增加了0.93×108m3和1.2×108m3。冰川徑流對總徑流的貢獻率分別增長了4.84%和9.38%，冰川加速消融導致冰川徑流占比增大是徑流增大的主要原因。RCP8.5 比RCP4.5 冰川消融提前了一個月，表明升溫會導致冰川消融時間提前。冰川消融峰值和汛期重合，會加劇了水系統(tǒng)脆弱性和河川徑流的波動性，使汛期極端水文事件的發(fā)生概率多于歷史時期。

（3）按照季節(jié)將歷史和2 種氣候模式的徑流、降水、氣溫進行分類并進行相關性分析發(fā)現(xiàn)，隨著氣溫的升高，冰川消融時間提前，冰川消融加速，冰川的積累時間減少，導致冰川面積進一步的縮減。