劉 寧,張 霞,祝雪萍,趙志懷,武鵬林
(太原理工大學(xué)水利科學(xué)與工程學(xué)院,山西太原030024)
碧流河流域水資源是供應(yīng)城市用水的主要水源,徑流則是更直接的水源。研究評價SWAT模型對該流域的徑流模擬情況,結(jié)合流域水文要素變化以及利用自動率定法率定參數(shù)并較準確的預(yù)測徑流、模擬徑流意義重大,將為進一步模擬未來氣候情景下的徑流奠定基礎(chǔ)。
基于物理基礎(chǔ)的SWAT模型在研究流域水循環(huán)的影響方面作出了積極貢獻,可以模擬流域內(nèi)發(fā)生的各種物理過程,被應(yīng)運于很多方面[1-2]。在非點源污染方面,Li等研究SWAT模型在江西省漣水流域進行多尺度的泥沙模擬研究以及氮和磷的輸移過程[3]。在徑流模擬方面,姚蘇紅等[4]以內(nèi)蒙古閃電河流域為例,應(yīng)用SWAT模型進行徑流模擬,結(jié)果顯示模擬結(jié)果受參數(shù)影響很大。張余慶[5]等學(xué)者以及宋增芳[6]等學(xué)者,均利用SWAT結(jié)合SUFI-2分別對修水流域、石羊河流域進行徑流模擬,有效控制了參數(shù)對模擬結(jié)果的影響。曹明亮[7]以碧流河流域為例,利用SWAT模型率定法對多站點進行研究,驗證了模型的適用性;但在參數(shù)敏感性分析、模擬結(jié)果評價以及不確定分析方面使用SWAT模型自帶的調(diào)參功能進行參數(shù)校準,該方法不僅耗時長、分析率定的參數(shù)較少、對經(jīng)驗要求高且易受人為因素影響。
故本文將SWAT模型應(yīng)用于碧流河流域徑流模擬。在分析流域水文要素變化的基礎(chǔ)上,結(jié)合自動率定工具SWAT-CUP,利用其中的SUFI-2算法對參數(shù)自動率定、敏感性排序,并對模擬結(jié)果分析評價,彌補了SWAT模型自帶校準功能的不足,最后多角度分析模擬結(jié)果。
碧流河屬于北溫帶濕潤氣候區(qū),6月~9月為汛期,尤其集中在7月~8月多年平均降雨量為742.8 mm,多年平均徑流深為294.5 mm,多年多年平均蒸發(fā)量為1 404.4 mm,流域面積2 814 km2,壩址以上控制流域面積2 085 km2。1975年興建的碧流河水庫是流域內(nèi)的主要水庫(總庫容9.34億m3),主要功能為向大連市供水。流域內(nèi)以林地和耕地為主,以草甸土、棕壤性土為主要土壤類型。2001年建成了控制流域面積為313 km2的玉石水庫。
流域內(nèi)有小石硼、桂云花、天益等9個雨量站,有繭場、碧流河、玉石水庫3個水文站。其中,碧流河水文站作為徑流模擬站點(見圖1)。
圖1 碧流河流域示意
SWAT模型是流域尺度的分布式動態(tài)水文模型[2],以水量平衡方程為問題研究的驅(qū)動力,綜合考慮氣候、下墊面等因素,模擬水文過程。該模型的分布式運算特點如下:先根據(jù)DEM及實際水網(wǎng)把流域劃分為多個子流域;而后對子流域內(nèi)的水文響應(yīng)單元進行聚類分析,在此基礎(chǔ)上獲得不同時間尺度的模擬結(jié)果。
SWAT-CUP中內(nèi)嵌了5種算法。李倩楠[8]等對5種不確定性評價方法進行了比較,得出SUFI-2方法運算效率高,運行次數(shù)少,卻能得到較好的結(jié)果,故本研究選擇此法。
SUFI-2是利用拉丁超立方體抽樣法獲得參數(shù)值帶入模型進行模擬[9]。其敏感性分析通過Global以及One-at-a-time敏感性分析實現(xiàn),其原理可參見SWAT-CUP說明手冊。此外,該算法通過P-factor(95%不確定性區(qū)間內(nèi)的觀測數(shù)據(jù)的百分數(shù))和R-factor(95PPU上下限的平均距離與標準偏差的比值)評價參數(shù)不確定性程度,P-factor的理論值區(qū)間為0~100%、R-factor的理論值區(qū)間為0~+∞。當(dāng)P-factor=1,且R-factor=0時,為完全理想的模擬結(jié)果;故通常以P-factor→1,且R-factor→0的程度來判斷模擬效果。但P-factor的增大會導(dǎo)致R-factor也增大,因此要找到二者的相對最佳值。
SUFI-2算法:首先,需定義較大的參數(shù)范圍(結(jié)合文獻分析以及SUFI-2推薦的參數(shù)范圍確定);然后,按照SUFI-2算法運行要求,設(shè)定各參數(shù)范圍、模擬次數(shù)及實測值等文件,并選擇目標函數(shù);最后,由程序進行多次迭代,依據(jù)迭代結(jié)果納西效率系數(shù)(Ens)、相關(guān)系數(shù)(R2)、相對偏差(PBIAS)、P-factor和R-factor的值)調(diào)整參數(shù)范圍,直至獲得理想的迭代結(jié)果。
前述迭代結(jié)果顯示的5個指標即為模型的評價指標。其中模型的不確定性評價用:P-factor、R-factor;PBIAS、Ens和R23個指標用于評價模擬精度,計算公式如下:
(1)
(2)
(3)
當(dāng)PBIAS<25%,Ens>0.5,R2>0.6時,模型的模擬結(jié)果可信;Ens≥0.65時,認為模擬結(jié)果很好[10-11]。P-factor>0.5且R-factor<1.5時,模擬不確定性可接受[12];當(dāng)P-factor>0.5且P-factor>0.7時,模擬的不確定性較小[13]。
水文要素主要包括降雨、徑流、溫度、等水文循環(huán)相關(guān)的時間序列。其中,降雨是徑流的直接來源,而徑流則是降雨最直觀的輸出形式。為了通過歷史趨勢了解水文要素的總體趨勢,參照文獻[14]給出的1958年~2011年的分析結(jié)果,對研究區(qū)間1978年~2011年的降水及徑流變化趨勢進行了分析(見圖2)。其分析可得出相似的結(jié)論:降雨和徑流都為下降趨勢,徑流的下降趨勢更明顯。此外,與多年平均徑流量對比可知,豐水年在1979年、1985年及1995年出現(xiàn)、枯水年在2002年、2011年出現(xiàn);同時,1995年之前,偏豐水年多于偏枯水年,而1995年之后偏枯水年多于偏豐水年,同樣驗證了徑流的下降趨勢,其下降的原因可能與當(dāng)下氣候全球變暖導(dǎo)致的降水、氣溫等氣象要素有關(guān),有待下一步研究。
圖2 碧流河水庫流域1958年~2011年降雨徑流系列變化
SWAT模型需要輸入的數(shù)據(jù)和來源如表1所示。氣象數(shù)據(jù)包括相對濕度、溫度、風(fēng)速、太陽輻射等。
表1 模型數(shù)據(jù)
本文使用的DEM數(shù)據(jù)為SRTM 90 m精度柵格型數(shù)據(jù),2 000 s 1∶10萬的土地利用分布圖,土壤數(shù)據(jù)選擇1∶100萬數(shù)字化土壤圖件。將各空間數(shù)據(jù)投影類型統(tǒng)一轉(zhuǎn)換為WGS_1984_UTM_ZONE_51N,并對應(yīng)于SWAT模型的代碼,重新分類(見表2、3)。
表2 碧流河水庫上游土地利用類型重分類結(jié)果
表3 碧流河水庫上游土壤類型
屬性數(shù)據(jù)的處理有:土壤屬性數(shù)據(jù)庫的建立和天氣發(fā)生器相關(guān)參數(shù)的計算。
由于土壤分類系統(tǒng)記錄的土壤屬性的參數(shù)與SWAT模型輸入要求的土壤屬性標準不同,在建立土壤屬性數(shù)據(jù)庫前要將國際制標準轉(zhuǎn)換為美制標準。對土壤水特性參數(shù)的獲得方式為:借助軟件SPAW以及直接應(yīng)用。
SWAT模型天氣發(fā)生器的參數(shù)計算由于工作量較大,運用計算程序SwatWeather.exe(北京師范大學(xué)建)進行,將氣象數(shù)據(jù)整理后導(dǎo)入程序,得到模型需要的數(shù)據(jù)[15]。
對碧流河水庫流域構(gòu)建SWAT模型,流域共劃分為29個子流域,碧流河水文站位于第21號子流域,限于當(dāng)前碧流河水文站缺乏日徑流數(shù)據(jù),對該水文站進行月模擬后輸出結(jié)果,并進一步在SWAT-CUP中按算法要求建立SUFI-2工程、輸入相應(yīng)的要求文件得到模擬評價結(jié)果,進而評價分析。
SWAT模型能夠相對準確地模擬徑流的關(guān)鍵在于將眾多的SWAT模型中與徑流有關(guān)的參數(shù)選出,去除對模擬結(jié)果影響較小的參數(shù),進行敏感性分析、率定。本文通過前述SWAT-CUP的敏感性分析方法最終選出14個最為敏感的參數(shù)校準,部分參數(shù)含義具體如下:CN2為平均濕度時的SCS曲線數(shù),值越大,徑流總量越大;SOL_AWC為田間持水量越大,可下滲水量越少,滲透時間越短;GW_DELAY用于計算每天流入地下水層的回灌水量,與水平面的深度和滲漏區(qū)及地下水區(qū)地質(zhì)組成的水力特性有關(guān);GW_REVAP值越大,再蒸發(fā)量越多,當(dāng)?shù)扔?時,蒸發(fā)量等于潛蒸發(fā)量;CH_K2:默認為0,值越大,河道輸移損失越大;REVAPMN:只有當(dāng)淺層地下水層的含水量超過閾值aqshthr,rvp時,再蒸發(fā)過程才發(fā)生。CH_N2:值越大,匯流流速越小;RCHRG_DP值越大,回灌水流入深層地下水層的比例越大;SOL_K為值越大,由土壤層流入下一土壤層的水量就越多,同時產(chǎn)生的壤中流也越多。
表4 參數(shù)敏感性及參數(shù)率定分析結(jié)果
注:①R、V為兩種修改參數(shù)的方法,分別是乘以某值(百分百浮動)、賦值;②t-stat越大表示參數(shù)越敏感,p-value值表示t值的顯著程度,越接近0值越顯著。
按照相應(yīng)的修改方案設(shè)定參數(shù)初始范圍,經(jīng)過多次迭代(每次迭代次數(shù)設(shè)置為200次)逐步滿足評價指標,得到參數(shù)敏感性分析以及率定結(jié)果如圖3及表4所示。
圖3 參數(shù)敏感性分析(P-Value及t-stat)結(jié)果
圖3及表4呈現(xiàn)率定后參數(shù)的最終范圍、最佳值及敏感性判定原則,排在第一位敏感的參數(shù)為:CN2,其P值最大,t檢驗值最小,同樣的判斷原理,接下來依次為:SOL_AWS、GW_REVAP、GW_REVAP、GWQMN等,其對模擬結(jié)果影響逐漸減小。
按照碧流河降雨和徑流變化趨勢分析結(jié)果,將1978年~2011年模擬期的前兩年作為預(yù)熱期,1981年~1995年作為率定期,1996年~2011年作為模型的驗證期,建立模型并獲得模擬結(jié)果。利用前述3個評價指標對模擬結(jié)果進行評價。模型在不同時期的月值模擬結(jié)果如圖4、5所示,評價指標見表5。
模擬結(jié)果顯示:①碧流河水文站的模擬值與觀測值變化趨勢一致性良好、流量過程擬合良好,率定期、驗證期R2分別為0.92、0.91;Ens分別為0.91、0.90;PBIAS分別均控制在5%以內(nèi),且全部落在較小的不確定區(qū)間內(nèi)(95PPU)。②但對峰值的模擬效果較差,模擬值普遍小于實測值,可能由于未考慮融雪影響而致。③驗證期(1995年~2011年)的模擬效果略低于率定期(1980年~1995年)的模擬效果,可能與SWAT模型對豐、枯水年的模擬精度不同有關(guān),偏枯水年模擬精度低于偏豐水年模擬精度[16]。
再從碧流河月徑流分布來分析,如圖6給出了兩個時段(1981年~1995年、1996年~2011年)多年平均月徑流變化,6月~9月是碧流河豐水期,8月達到徑流峰值,之后開始下降,10月至次年5月為枯水期,豐枯時段區(qū)分明顯,這與碧流河年內(nèi)降水時間分布有關(guān)。一年降水峰值主要集中在主汛期(7月~8月附近),圖中同樣可見模型對峰值模擬效果較差。
圖4 月徑流率定期模擬結(jié)果(1981年~1995年)
圖5 月徑流驗證期模擬結(jié)果(1996年~2011年)
時期P-factorR-factorR2EnsPBIAS/%率定期0.880.610.920.91-3.80驗證期0.910.610.910.901.30
圖6 觀測和模擬的多年平均月徑流變化
水文過程包含著各種不確定性因素。對于水文模型,通常是運用大量的概化或經(jīng)驗公式對實際流域的水文循環(huán)過程的一個簡化,各方面的不確定性往往相互作用而影響模型輸出結(jié)果的不確定性。
本文不確定性結(jié)果顯示,率定期P-factor為0.88、R-factor為0.61,依據(jù)判別條件,滿足P-factor>0.7、R-factor<1不確定性較小的條件,模型模擬結(jié)果較適用。驗證期與率定期基本一致,而略高的P-factor表明:90%的觀測值落在95%置信水平的不確定性區(qū)間內(nèi),模擬的不確定性相對更小。但從圖4、5可看出,豐水期模擬的不確定性區(qū)間依然較大,而峰值模擬的不確定性更大。
SWAT模型在碧流河流域的徑流模擬中在前人手動率定較少參數(shù)的基礎(chǔ)上,延長模擬時期,利用自動率定SUFI-2算法提升了模擬效果,并驗證了相對較小的不確定性,為下一步的氣候模式與SWAT模型耦合筑牢了基礎(chǔ);同時也發(fā)現(xiàn),該模型對碧流河的峰值模擬效果不佳、不確定性只是綜合考慮,應(yīng)在后續(xù)的研究中深入綜合考慮更多氣候變化與人類活動等因素。另外限于水文資料的限制,僅利用一個碧流河水文站序列進行分析,代表性方面有一定欠缺,需要進一步獲取更多水文站數(shù)據(jù)深入研究;此外,敏感性分析得到的14個參數(shù)其各參數(shù)對結(jié)果的影響貢獻率多大,有待繼續(xù)探討。