張程鵬,張鳳娥,2,耿新新,2,王 偉,冀俊杰
(1.中國地質科學院水文地質環境地質研究所,石家莊050061;2.自然資源部地下水科學與工程重點實驗室,石家莊050061;3.貴州省畢節市水文水資源局,貴州畢節551700)
中國西南巖溶地區干旱、洪澇災害交替頻發,生態環境脆弱[1],巖溶旱澇災害也是西南地區石漠化治理的障礙之一[2]。除異常氣象條件之外,巖溶地區人類活動也是造成旱澇災害的原因之一。隨著人口經濟的發展,人類開墾荒山、植樹造林、填洼造田,改變了原始的土地利用方式,對流域水資源的形成、分布、轉化產生了極為深遠的影響[3]。因此,研究流域土地利用變化的水文響應,進而分析旱澇災害的產生原因,具有一定的科學意義和應用前景。
大部分學者對土地利用變化引起的流域徑流變化的研究是基于模型展開的。ZHANG[4]等利用MODIS-LAI 校準的SWAT模型研究了澳大利亞北約翰斯頓河流域土地利用變化對水平衡的影響,表明在降雨量、坡度和土壤質地相同的條件下,常綠林的產流量一般小于草地和城市用地,城市地表徑流大而側向徑流和地下水較少;LUO[5]等利用FLUS-SWAT模型研究了中國柘溪河流域城市擴張土地利用變化的水文響應,表明森林和農業用地的減少將導致地表徑流增加64.86%,地下水流量減少9.05%;ZHOU 等[6]基于了長江三角洲城市化水文響應,表明城市化對年產水量的影響較小,但對地表徑流、洪峰流量和洪量的影響顯著;王博威等[7]基于SWAT 模型分析了潘家口水庫土地徑流量減少的原因,耕地、草地向林地的轉化導致產流量減少,蒸發量增加;黎云云等[8]研究表明渭河流域土地利用變化導致徑流減少,支流減少程度大于干流;郝振純等[9]研究表明海河流域林地的增加和草地的減少會降低汛期徑流量以及最大月徑流量。
目前,大多數研究都聚焦在徑流量變化上,研究土地利用變化引起流域調蓄能力變化還不夠深入,且對典型巖溶流域土地利用變化的水文響應研究較少。因此,本文基于SWAT 模型從流域整體和子流域角度進行分析,評估典型巖溶流域土地利用變化引起流域流量、產水量、調蓄能力的變化,以期為決策者提供減輕巖溶旱澇、防治石漠化的技術支持。
研究區(105°5′E~105°16′E,27°15′N~27°22′N)位于貴州省畢節市七星關區(圖1),屬倒天河流域上游部分,區內主要分布村鎮,倒天河下游位于畢節城區內。倒天河發源于畢節市黃泥大婁山脈南麓,屬烏江水系支流白甫河(落腳河)上游。

圖1 研究區位置及區內水文要素Fig.1 Location of study area and hydrological elements
該區屬亞熱帶季風濕潤氣候,冬無嚴寒,夏無酷暑,季風氣候比較明顯,降雨量較為充沛。1951-2017年多年平均降水量為1 134 mm,多年平均氣溫為12.6 ℃,多年平均相對濕度為81.5%。研究區地勢西南高,東北低,海拔1 552~2 199 m,流域面積約105 km2,年均徑流量0.43 億m3,多年平均徑流深409.5 mm,汛期為每年6-9月。
收集數據主要包括數字高程模型(DEM 數據),土地利用/覆蓋數據(Land use 數據),土壤數據(Soil 數據),詳細信息見表1。

表1 地理空間數據介紹Tab.1 Data of geospatial
氣象數據主要包括研究區多年日值降水量、日值高低溫、日值相對濕度、日值太陽輻射量、日值風速,水文數據主要包括多年逐日平均徑流量。氣象站、水文站分布見圖1,詳細信息介紹見表2。其中,畢節站不在研究區內,但數據缺失時亦可使用,但效果不佳。SWAT 模型需要保證輸入氣象數據的起始日期和結束日期相同,其中缺測值可以用-99 代替,因此,綜合氣象和水文數據的時間序列,本文采用1990—2006年的氣象數據和流量數據來運行SWAT模型。

表2 氣象水文數據介紹Tab.2 Data of meteorological and hydrological
數據屬性庫主要包括Land Use屬性庫、Soil屬性庫、坡度。
(1)Land Use 屬性庫。原始Land Use 數據自帶二級分類數據,按模型要求重分類后添加對應屬性,見表3。

表3 土地利用屬性數據Tab.3 Attribute data of Land Use
(2)Soil 屬性庫。模型需要手動添加Soil 屬性數據,根據原始數據自帶的屬性表,利用SPAW 軟件的SWC 模塊計算新的Soil屬性數據,搭建了適用的土壤庫。
(3)坡度。研究區坡度按照自然間斷點分級法(Jenks)分為3級,即0~8.6,8.6~18.2,18.2~72。
SWAT 模型由美國農業部(USDA)農業研究機構的工程師J G Arnold開發[11]。其原理為水量平衡方程[12],具體如下:

式中:?SW、P、Q、ET、DP、QR分別表示土壤含水量、降水、地表徑流、實際蒸散發、深層下滲和淺層回歸流。
基于ArcGIS平臺的ArcSWAT按照輸入的DEM 數據自動提取河網,生成子流域,再按照不同的LandUse 數據、Soil 數據、坡地將子流域劃分為HRU(水文響應單元),本文劃分閾值為5%、10%、10%,通過輸入降水量、溫度、太陽輻射量、相對濕度、風速等信息,計算每個HRU的產匯流情況。
通過SWAT-CUP 軟件[13,14]不斷改變相關參數,使得模型輸出的河道流量水文信息與實測的流量數據不斷擬合,以期得到最佳擬合效果時的最佳參數和最佳參數范圍。本文選取R2(相關系數)和NSE(Nash-Sutcliffe efficiency coefficient,納什效率系數)兩個目標函數來評價擬合效果,其計算方式為:


式中:Q表示流量,m3/s;m和s分別表示觀測值和模擬值;i表示觀測或模擬的次數。
目標函數的值越接近1,表示模擬值與實際值差別越小,模擬效果越好。目標函數置信區間見表4。

表4 目標函數置信區間Tab.4 Confidence interval of objective function
為模擬真實土地利用情況,1990-1997年使用1995年土地利用,1998-2002年使用2000年土地利用,2003-2006年使用2005年土地利用。以月步長運轉模型,將1991-1992年作為預熱期,1993-1997年作為校準期,1998-2002年作為驗證期,2003-2006年作為驗證期。
(1)流域劃分結果。模型劃定的研究區面積為104.970 km2,研究區被分為55 個子流域,其中出水口在12 號子流域,1995年土地利用下HRU 個數為454,2000年土地利用下HRU個數為456,2005年土地利用下HRU個數為431。
(2)參數率定結果。根據大量文獻及學者建議,初步選取12 個相關參數進行率定,在1 000 次模擬后,選擇較敏感的8 個參數作為選定參數,利用SWAT-CUP 軟件進行調參,經每輪1 000 次模擬,5 輪迭代后,得到p-factor和r-factor合理的值(p>0.75,r<1),得到最佳參數范圍和最佳參數。
(3)模型評價。流量模擬結果見圖2。校準期(1993-1997年)R2=0.95,NSE=0.95,p-facto=0.78,r-factor=0.42,實測平均流量1.568,模擬平均流量1.564,誤差0.3%;驗證期(1998-2002年)R2=0.95,NSE=0.95,實測平均流量1.373,模擬平均流量1.447,誤差5.3%;驗證期(2003-2006年)R2=0.90,NSE=0.77,實測平均流量1.067,模擬平均流量1.226,誤差14.9%。對比公式(2),(3)和表4,表明模擬結果良好,參數率定后的SWAT 可以用來評價該研究區水文情況。
研究區主要分布零散的村鎮,土地利用類型主要有3種:耕地、林地、草地。為研究土地利用變化對該流域的水文響應,保持1990-2006年氣象輸入條件不變,設置5種情形。

表5 敏感參數率定結果Tab.5 Calibration results of sensitive parameters

圖2 流量模擬結果Fig.2 Flow simulation results

表6 不同土地利用情景Tab.6 Land use scenarios
由土地利用變化分布圖(圖3)、土地利用類型占比圖(圖4)以及土地利用轉移矩陣(表7、表8)分析得出,在20世紀末到21世紀初,研究區的土地利用變化趨勢為:
(1)耕地面積逐漸增加,林地面積逐漸增加,草地面積逐漸減小;
(2)1995年到2000年,土地利用類型變化不大;2000年到2005年,土地利用類型變化劇烈,主要表現為草地轉化為耕地和林地。
將5 種情景下的土地利用數據代入模型,得到研究區流量數據,因不同情景下流量差異不大,做出了不同情景流量差值圖(圖5),可以更加直觀的比較不同情景的流量差異。不同情景流量差值圖表明:
(1)不同情景下的流量差異主要體現在雨季(汛期),旱季(非汛期)差異不明顯;
(2)雨季差值為正、旱季差值為負,表明研究區調蓄能力減弱;雨季差值為負、旱季差值為正,表明研究區調蓄能力增強;
(3)情景2-情景1 反映了1995-2000年土地利用變化的水文響應,情景3-情景2 反映了2000-2005年土地利用變化的水文響應,兩種土地利用變化均導致流域流量增加,調蓄能力減弱,并且后者比前者的影響更加明顯;表明研究區土地利用類型由草地轉為耕地會導致流域流量增加,調蓄能力減弱;
(4)情景4-情景3、情景5-情景3 反映了在2005年土地利用的基礎上按照國家政策退耕還草、還林后流域流量變化,情景5-情景4 反映了耕地還草、還林兩種方式之間的差異;表明退耕還草、還林均會導致流域流量減小,調蓄能力增強,且還草與還林兩種方式導致流量變化幾乎一致。


圖3 1995-2000土地利用分布圖及土地利用變化圖Fig.3 Land use distribution map and land use change map from 1995 to 2000

圖4 不同時期土地利用類型占比Fig.4 Proportion of Land use types in different periods

表7 1995-2000年土地利用轉移矩陣 km2Tab.7 Land use transfer matrix from 1995 to 2000

表8 2000-2005年土地利用轉移矩陣 km2Tab.8 Land use transfer matrix from 2000 to 2005
流域產水量的計算方法為:

式中:WYLD表示子流域產水量;SURQ表示地表徑流向河道輸水量;LATQ表示側向補給量;GWQ表示地下水向河道輸水量;TLOSS表示徑流向河道傳輸中的傳輸損耗;PA表示抽水損耗。上述變量單位均為mm。
產水量表征一個流域匯水、輸水、產水能力的重要指標,通過分析不同土地利用下子流域產水量變化特征,分析土地利用類型變化對流域產水能力的影響。
在1995-2000年土地利用類型變化中,重點研究27、40、46號子流域產水量變化;在2000-2005年土地利用類型變化中,重點研究19、22、50、51號子流域產水量變化。子流域產水量差值圖表明。
(1)1995-2000年土地利用變化下(圖6)40、46 號(林地→耕地)子流域產水量差異較大,27 號(草地→林地)子流域產水量差異較小;
(2)2000-2005年土地利用變化下(圖7)19 號(草地→耕地)子流域產水量差異最大、22 號(草地→耕地)子流域產水量差異較大,50、51號(草地→林地)子流域產水量差異較小;
(3)林地向耕地、草地向耕地的土地利用類型轉化引起子流域產水量的變化明顯,草地向林地的轉化引起子流域產水量的變化不明顯;

圖5 不同土地利用情景下流量差值圖Fig.5 Flow difference under different Land use scenarios

圖6 1995-2000年土地利用變化下子流域產水量差值圖Fig.6 Water yield difference of subbasin under Land use changes from 1995 to 2000

圖7 2000-2005年土地利用變化下子流域產水量差值圖Fig.7 Water yield difference of subbasin under Land use changes from 2000 to 2005
(4)上述土地利用變化方式下均表現為雨季差值為正、旱季差值為負,表明子流域相較于土地利用方式變化前調蓄能力減弱;林地向耕地、草地向耕地、草地向林地的土地利用轉化均會導致研究區調蓄能力變弱,草地向耕地的土地利用轉化對調蓄能力的改變尤為明顯。
為了研究土地利用方式變化在典型巖溶流域產生的水文響應,構建了倒天河流域的SWAT模型,并利用該模型探討了不同土地利用下流域流量的變化規律和子流域產水量的變化規律及其調蓄能力分析,其具體結論如下。
(1)參數率定后的SWAT 模型在校準期(1993-1997年)R2=0.95、NSE=0.95;驗證期(1998-2002年)R2=0.95、NSE=0.95;驗證期(2003-2006年)R2=0.90、NSE=0.77。構建的SWAT 模型在模擬倒天河流域月流量過程中模擬效果良好,模型適用于該研究區;
(2)3期土地利用數據表明:研究區耕地和林地面積逐漸增加,草地面積逐漸減小,土地利用類型變化主要為由草地轉化為耕地和林地;
(3)研究區5 種情景下土地利用變化引起流域流量變化結果表明:草地轉化為耕地會導致流域流量增加,流域調蓄能力減弱;
(4)不同土地利用方式下子流域產水量結果表明:林地向耕地、草地向耕地的土地利用類型轉化引起子流域產水量的變化明顯,草地向林地的轉化引起子流域產水量的變化不明顯;上述三種土地利用方式的轉化均會減弱子流域的調蓄能力,其中,草地向耕地的轉化對調蓄能力的減弱尤為明顯。
綜上所述,本項工作研究了巖溶流域土地利用變化引起流域產流量、調蓄能力變化的規律,明確了巖溶區土地利用變化產生的水文響應,推動了巖溶區土地利用變化研究進展。據此可以說明“退耕還林還草”措施可以提高巖溶流域調蓄能力,從而降低下游洪澇災害發生的風險;還有利于提高巖溶流域林、草覆蓋率,進而減少石漠化的發生。
盡管如此,本項工作還有不足,未對巖溶流域調蓄能力進行定量評價,沒考慮氣候變化對水文響應的影響。因此,將氣候變化和土地利用變化對流域調蓄能力的定量研究是今后應突破的重點。 □