基于SWAT模型的赤水河流域徑流年內分配特征及其對降水的響應研究

劉家威, 蔡 宏, 鄭婷婷, 唐 敏

(貴州大學 礦業(yè)學院, 貴州 貴陽 550025)

赤水河是長江上游的重要支流，也是經過貴州省的主要河流之一。赤水河流域降水和徑流的年內分配和年際變化情況，深刻地影響著流域內的水文水資源特征、生態(tài)環(huán)境以及當地居民的生活質量與生命財產安全。河川徑流是由流域內氣象要素和下墊面共同作用的產物[1]，對徑流影響最大的因素包括降水、氣溫和人類活動等[2-4]。其中，氣溫主要是從兩個方面對徑流產生影響，一是影響蒸散發(fā)，二是影響融雪。赤水河發(fā)源于云南省鎮(zhèn)雄縣板橋鎮(zhèn)長槽村滮水巖，發(fā)源地以林地為主，沒有顯著的雪山[5]。因此相較黑河[4]等流域內存在雪山的河流，赤水河徑流量受溫度的影響僅需考慮蒸散發(fā)。另外，赤水河是少數沒有修建大型水壩的河流，且該地區(qū)長期注重生態(tài)保護，人類活動對赤水河流域徑流量的影響較小。綜上所述，相比其他河流，赤水河徑流受降水影響更加顯著，因此研究該流域徑流對降水的響應特征顯得很有必要。

國內學者已經做了很多關于區(qū)域降水/徑流變化特征的研究。Zhang等[6]把候降水量看作矢量，提出一種度量降水年內非均勻分配的方法，很好地反映年總降水量年內非均勻分配特性。劉賢趙等[7]將這種方法推廣到月尺度并用到徑流年內分配的研究上，研究發(fā)現，月徑流集中度比徑流年內不均勻系數具有更高的分辨能力和敏感性，用集中度和集中期能夠充分反映徑流年內分配的不均勻性，具有客觀準確、概念明確、適用性強的優(yōu)點。但國內學者在利用集中度和集中期進行相關問題的研究時，仍存在兩個問題： ①相關研究僅單獨對降水或徑流的集中度與集中期進行討論，少有人同時對兩者進行分析，并以此探討徑流對降水的響應情況，這使得降水與徑流間年內分配的特征關系難以被很好地解釋； ②在對徑流的相關研究中，主要基于水文站點多年的實測數據進行。由于數據收集難度較大，往往只能以一兩個站點的實測數據來反映流域整體狀況，無法對流域內各支流的情況進行詳細探究，所得結論也是各區(qū)域綜合影響下的結果。這使得對各支流流域徑流的年內分配特點及其對降水年內分配的響應特征描述不夠細致。

SWAT模型在大尺度流域水文過程的模擬預測方面有著廣泛的應用，故本文以赤水河流域中上游為研究區(qū)，利用多源數據，建立SWAT模型，從而完成對流域內各子流域降水及徑流的模擬。以SWAT模型輸出結果，結合降水和徑流的集中度與集中期，研究流域降水與徑流的分配和變化情況以及徑流對降水變化的響應關系，以期為赤水河流域水文預報、河流水質與流域生態(tài)保護及防洪治澇等提供有效的參考。

1 數據來源

1.1 研究區(qū)概況

赤水河是長江上游的主要支流，發(fā)源自云南省鎮(zhèn)雄縣，自東向西流經云南、貴州和四川3省，至合江城東匯入長江。地理位置：27°15′2″—28°25′59″N，104°44′19″—106°58′36″E。本文以赤水河流域中上游為研究區(qū)，以丙安古鎮(zhèn)作為赤水河中上游和下游的分界點。赤水河中上游河長約258.23 km，流域集水面積約為11 610 km2，巖溶十分發(fā)育，兩岸陡峭，多險灘急流，地表坡降大，用地類型以林地、草地和耕地為主。研究區(qū)內以典型的中亞熱帶濕潤氣候為主，氣溫較高，降水主要表現形式為降雨，且集中在夏季，冬季較少，年內分配極不均勻。研究區(qū)地勢特征見圖1。

圖1 赤水河流域地勢

1.2 數據來源

本次試驗中所使用的數據包括赤水河流域中上游的數字高程模型(DEM)，2015年土地利用數據、土壤數據、研究區(qū)內部及周圍5個氣象站點的日值氣象數據以及二郎壩站和茅臺站兩個水文站點的實測徑流數據。其中數字高程模型(DEM)、土壤類型分布圖以及土地利用類型分布圖定義的投影坐標系均為WGS_1984_UTM_Zone_48 N，各項數據詳細情況及來源見表1。

表1 數據詳情及來源

2 研究方法

2.1 SWAT模型構建

SWAT(soil and water assessment tool)模型近年來已得到了廣泛應用。國內部分學者在不同流域成功建立了SWAT模型并完成了對流域徑流的模擬[8-9]。水文響應單元(HRU)是SWAT模型的基本計算單元，由土地利用、土壤類型和坡度共同定義，通過閾值劃分[10]。SWAT模型可結合給定的氣象數據、水文數據以及管理措施等信息，單獨計算每個HRU的內部循環(huán)，并在子流域進行累計匯總，再通過河網對子流域進行有機連接以模擬地表匯流。

2.1.1 土地利用數據庫 本研究所使用的土地利用數據來自中科院資源環(huán)境科學數據中心，研究區(qū)土地利用類型共被劃分為了18個二級類，與SWAT模型內采用的由美國地質調查局指定的分類系統(tǒng)不一致[11]。故根據研究區(qū)實際情況，綜合考慮模型自帶的土地利用屬性數據庫，將研究區(qū)土地利用類型重分類為6類，分別為：林地、耕地、草地、水域、建設用地及未利用地。

2.1.2 土壤數據庫 土壤數據庫中所需土壤相關的物理屬性參數較多，主要存在3種收集途徑。其中大多數參數可以直接在HWSD的屬性數據庫中查詢獲?。煌寥缹咏Y構(TEXTURE)、土壤濕密度(SOL_BD1)、土壤可利用有效水(SOL_AWC1)以及飽和水力傳導系數(SOL_K1)可利用由美國農業(yè)部與華盛頓州立大學開發(fā)的SPAW軟件計算得到；土壤侵蝕力因子(USLE_K1)則可利用現今被普遍認可的由美國學者Williams提出的計算方法進行求算。

2.1.3 氣象數據庫 氣象發(fā)生器用于當站點氣象數據缺測時對缺失值進行填補[12]。氣象數據庫應采用盡量長時間序列的氣象數據進行構建，以使得其中的各項參數盡可能地趨近于真實值。本次試驗共收集到了研究區(qū)內及周圍共5個氣象站2000—2020年逐日的降水、最高氣溫、最低氣溫、平均氣溫、平均氣壓、平均風速、相對濕度以及日照時數的實測數據用以構建氣象發(fā)生器。SWAT模型氣象發(fā)生器參數計算工作量大，本研究運用北京師范大學數字流域實驗室楊霞等[13]研究集成的各個參數公式建立的計算程序SWATWeather進行相關參數的計算。

2.1.4 子流域與水文響應單元的劃分 已有研究發(fā)現子流域劃分存在一個合理的水平[14]，子流域劃分個數不宜過多或過少[15]。本文參考史曉亮等[16]的研究成果，設置最小集水面積為研究區(qū)總面積的2%(250 km2)，將研究區(qū)劃分為38個子流域，劃分結果見圖2。考慮到下墊面因素對流域水文過程的影響顯著，故本文結合研究區(qū)實際情況和研究目的，分別設定土地利用、土壤與坡度的最小面積百分比閾值為5%，10%與5%，將研究區(qū)劃分為776個HRU。

注：阿拉伯數字1—38為子流域編號。下同。

2.2 模型評價方法

本研究以二郎壩站和茅臺站的實測逐月徑流數據為參照，利用SWAT-CUP工具中的SUFI-2算法對SWAT模型輸出結果進行率定及驗證。由于研究區(qū)涉及區(qū)域較大，赤水河中上游除干流外還存在有較多支流，其中尤以桐梓河(赤水河一級支流)涉及的范圍最廣。在劃分的38個子流域中，桐梓河共涉及13個子流域。考慮到赤水河中上游干流與桐梓河所涉及區(qū)域的各項參數可能存在差異，故本文采用先支流后干流的率定驗證方法，先以二郎壩站實測徑流數據對桐梓河徑流模擬結果進行率定驗證，再以茅臺站實測徑流數據對赤水河流域中上游干流徑流模擬結果進行率定驗證。選用決定系數(R2)及納什效率系數(Ens)作為模型模擬精度的評定指標，其計算公式如下：

(1)

(2)

式中：R2表示確定性系數； Ens表示納什效率系數；Qs,i為第i個模擬值；Qsa表示模擬值的平均值；Qo,i為第i個實際觀測值；Qoa為實際觀測值的平均值；n為值的個數。R2值越接近于1，表明模擬值越接近實測值，模擬效果越好； Ens用于評價模型模擬值與實測值的擬合程度，其值越接近1，說明模擬值與實測值擬合程度越高[17]。一般來說R2≥0.6，Ens≥0.5同時達到時，就說明結果達到模擬精度要求[18-20]。

2.3 流域降水總量與徑流總量的計算

值得注意的是，SWAT模型輸出降水量模擬結果的方式是：分別輸出各子流域降水量(mm)。因此在計算流域降水總量(m3)時，不能簡單地將模型輸出的各子流域降水量相加，而是需將模型輸出的各子流域降水量進行單位換算后，乘以對應子流域的面積，得到各子流域的降水總量，最后再相加得到流域的降水總量。類似地，SWAT模型輸出的徑流量模擬結果單位是(m3/s)，需乘以對應的時長，才能與流域降水總量單位統(tǒng)一。

2.4 集中度和集中期的計算

為方便研究流域降水和徑流的時間變化規(guī)律，本文引入降水和徑流的集中度與集中期兩個量化指標。其中，降水集中度(PCD)是一個能反映降水量在各個月份集中程度的無量綱數，取值范圍為0～1，PCD值越接近于1，說明其年降水量越集中于某一時間段；反之越接近于0，說明其各時間段的降水量分配較為均勻。降水集中期(PCP)則用來描述降水的集中時段[21]。徑流集中度(RCD)與徑流集中期(RCP)同理。

集中度、集中期是以矢量圓均等分圓周角為辨別區(qū)間，角度指示降水集中在某一時間段，集中程度的量化可以從集中度的矢量模數量值體現[22]。當采用月為計算時段時，每個月的天數是不同的，因此，必須做一定程度的概化處理，即不考慮月大、月小，均視為同一個時段長[2]。降水集中度與集中期計算公式如下：

(3)

(4)

(5)

(6)

式中：PCDi，PCPi分別為研究區(qū)第i年的降水集中度和集中期；Ri為研究區(qū)在第i年內的降水總量；rij為研究區(qū)第i年第j月的降水量；θj為第j月對應的方位角(一年的方位角為一個圓周，即360°)。

3 結果與分析

3.1 模型參數率定驗證

本次研究以2007—2009年作為模型的預熱期，將收集到的兩個水文站10 a實測逐月徑流數據分為兩個部分，其中2010—2017年部分用于模型率定，2018—2019年部分用于模型驗證。模型率定驗證的結果見表2。水文站逐月徑流量實測值與模擬值對比見圖3。從表2和圖3中可以看出模型模擬結果較好，基本達到精度要求。

表2 亦水河流域二水文站徑流模擬結果評價

圖3 二郎壩站和茅臺站的徑流實測值與模擬值

3.2 徑流年內年際分配情況及對降水的響應

基于模型模擬結果計算得到研究區(qū)2010—2019年月降水和月徑流總量，結果見圖4。從圖4可知，總體上赤水河流域中上游月徑流總量受月降水總量影響顯著，兩者變化趨勢基本一致，且年內降水與徑流分配極不均勻。研究時段內月降水量與月徑流量最大值均出現在2014年7月，該月降水量達3.09×109m3，約占全年降水總量的21.1%；月徑流量達1.74×109m3，約占全年徑流總量的31.5%。而月降水量與月徑流量最小值均出現在2010年2月，該月平均降水量為9.95×107m3，平均徑流量為1.52×107m3，分別占全年平均降水量與徑流量的0.93%與0.55%。各年降水與徑流主要集中在夏季的6—8月。就降水而言，6—8月降水總量占多年年均降水總量的45.4%，最大降水月的降水量是最小降水月降水量的9.8倍。而6—8月徑流總量占多年年均徑流總量的52.0%，最大徑流月的徑流量是最小徑流月徑流量的14倍。由此可見，相對降水而言，徑流的年內變化更為劇烈，不穩(wěn)定性更強，其年內分配的不均勻性更為突出。

圖4 赤水河流域中上游多年月降水/徑流總量變化

3.3 研究區(qū)降水/徑流集中度

2010—2019年赤水河流域中上游降水集中度與徑流集中度計算結果見圖5。由圖5可知，研究區(qū)徑流集中度與降水集中度多年變化趨勢基本一致，降水集中度對徑流集中度影響較為顯著；同時，徑流集中度顯著高于降水集中度。

為探究以上規(guī)律的普適性，本文計算了研究區(qū)13個支流流域2010—2019年年均徑流/降水集中度，其空間分布情況見圖6(通過自然斷點法劃分等級)。其中，西南部各支流降水/徑流集中度較大，西北部各支流降水/徑流集中度較小，徑流與降水集中度空間分布高度一致，進一步證明研究區(qū)徑流集中度受降水影響顯著。13個支流流域2010—2019年共130組降水/徑流集中度計算結果見圖7。其中徑流集中度大于降水集中度的有119組，徑流集中度小于降水集中度的僅有11組。

圖6 赤水河流域各支流多年平均降水/徑流集中度空間分布特征

圖7 赤水河流域各支流流域多年降水/徑流集中度

經計算，赤水河流域中上游年平均徑流系數約為0.31，這意味著流域降水中大部分通過蒸散發(fā)或被下墊面截留等方式損失，少部分轉化為徑流。圖8為赤水河流域中上游2010—2019年各年旱季(1—4月及11，12月)與雨季(5—10月)徑流系數。其中，各年雨季徑流系數均大于同年旱季徑流系數，這使得徑流集中度往往大于降水集中度。

圖8 赤水河流域多年旱季與雨季徑流系數

已有研究發(fā)現，土壤、土地利用及流域地質條件的變化對徑流系數影響較小，而降水量變化對徑流系數的影響較大[23]。根據水量平衡原理，在任意時段內，流域降水量等于產生的徑流量、蒸散發(fā)量及蓄水變化量之和。一方面，當流域發(fā)生降雨事件時，降雨在地表主要有兩個分流過程(由于降雨時間短暫，降雨蒸發(fā)量忽略不計)，即入滲和徑流[24]。當降雨強度小于或等于土壤入滲率時不產流，并以降雨強度向下入滲；當降雨強度大于入滲率時，則形成超滲產流；當土壤最大蓄水量小于累積入滲量時，則形成蓄滿產流[25]。因此，當降雨強度大于入滲強度時，降雨強度越大，則入滲造成的損失占比越低，尤其是當累積入滲量大于土壤最大蓄水量后，入滲作用將不再顯著，此時降水損失比例最小，大部分降水將直接轉化為徑流。另一方面，由于潛在蒸散發(fā)為在充分供水條件下的區(qū)域蒸散發(fā)能力，當未滿足充分供水條件時，區(qū)域的實際蒸散發(fā)小于潛在蒸散發(fā)，此時實際蒸散發(fā)隨供水量增加而增加；當滿足充分供水條件時，則實際蒸散發(fā)等于潛在蒸散發(fā)，并不再隨降水量增加而繼續(xù)增加。

由此可知，當降水量較小時，由于入滲和蒸散作用，其變化對徑流量的影響并不顯著，此時徑流系數較小。當降水量足夠大時，入滲量和蒸散發(fā)量不再隨降水量增加而增加，此時降水量越大，則徑流系數越大。因此研究區(qū)雨季徑流系數通常大于同年旱季徑流系數。Norbiato D等[26]在對意大利北部14個集水區(qū)進行研究時也發(fā)現，徑流系數隨降水總量的增加而增加，且徑流系數的空間變異性可以被降水量較好地解釋。

進一步觀察11組徑流集中度小于降水集中度(RCD≤PCD)的數據發(fā)現，當降水集中度較大時(PCD>0.3)，PCD與RCD的差值均在0.03以內；僅當降水集中度較小時(PCD<0.3)，才會出現RCD顯著小于PCD的情況。這是由于當降水集中度較小時，該地區(qū)降水的年內分配比較均勻且整體處于較低水平。此時，年內降水的小幅變化對徑流集中度的影響極小，而氣溫、人類活動等因素對徑流集中度的影響也變得更為顯著。因此當降水集中度較小時，徑流集中度不再完全以降水集中度為主導，可能出現徑流集中度顯著大于降水集中度的情況。

3.4 研究區(qū)降水/徑流集中期

2010—2019年赤水河流域中上游降水/徑流集中期計算結果見圖9。研究區(qū)降水與徑流集中期多年變化趨勢基本一致，前者對后者影響顯著。研究區(qū)多年平均降水集中期與徑流集中期分別為179和188，其對應日期分別約為6月29日和7月8日，相差約9 d。由于絕大部分降水需流經下墊面才能匯入河道，使得流域徑流對降水變化的響應存在滯后效應[27]，故徑流集中期往往大于降水集中期。

圖9 赤水河流域多年降水/徑流集中期

針對圖9中研究區(qū)2013年PCP

圖10 赤水河流域各支流2013年降水/徑流集中期

對比以上支流流域2013年逐月降水和徑流，結果如圖11所示。各支流流域年內降水量均存在兩個峰值，分別為5月和8月；而徑流量僅在6月有一個峰值，部分支流徑流量在8月有小幅上升，但并不明顯，這使得以上支流流域徑流集中期小于同年降水集中期。

圖11 赤水河流域各支流2013年逐月降水和徑流

對比兩個降水峰值，各支流5月降水量均大于8月，6月降水量也維持在較高水平，而9月降水量大幅下降。經計算，5，6月各支流平均降水總量是8，9月的1.73倍，而5，6月各支流平均徑流總量是8，9月徑流總量的4.97倍。結合前述分析，說明5，6月降水總量足夠大，使得降水損失遠低于8，9月，而降水的損失比例在不同降水量條件下的變化，是各支流徑流量在兩個降水峰值期間有著不同響應的原因。

4 結 論

(1) 本研究實現了SWAT模型對赤水河流域中上游月徑流量的有效模擬。模型率定驗證時采用了先支流后干流的方式，茅臺站與二郎壩站在率定期的R2與Ens均大于0.83，驗證期的R2與Ens均大于0.69，說明SWAT模型能很好地模擬赤水河流域中上游的徑流量。

(2) 赤水河流域中上游徑流量變化趨勢與降水量變化趨勢基本一致，徑流量受降水量影響顯著。降水和徑流年內分配均主要集中于6—8月，但相對于降水，徑流年內分配的不均勻程度更高。

(3) 研究區(qū)徑流集中度的時空分布受降水集中度影響顯著。同時，由于入滲和蒸散作用的影響，徑流系數在旱季和雨季有著不同的表現。旱季徑流系數低，而雨季徑流系數高，因此流域徑流集中度往往大于降水集中度。但當降水集中度較低時(PCD<0.3)，徑流集中度不再完全以降水集中度為主導。

(4) 研究區(qū)徑流集中期受降水集中期影響顯著。同時，流域徑流對降水變化的響應存在滯后性，徑流集中期往往大于降水集中期。但由于入滲和蒸散作用的影響，徑流對短期較小幅度的降水量增加并不敏感。故短期較小幅度的降水量增加會對降水集中期產生較為顯著的影響，但對徑流集中期影響有限。