武海霞，李清雪，孫玉壯，冀雅珍

(1.河北省水資源高效利用工程技術研究中心，河北工程大學，河北 邯鄲 056038；2.河北省煤炭資源綜合開發(fā)與利用協(xié)同創(chuàng)新中心, 河北省資源勘探研究重點實驗室, 河北工程大學，河北 邯鄲 056038；3.山西省水利水電科學研究院，太原 030002)

0 引 言

目前，SWAT分布式水文模型在流域水文、生態(tài)和環(huán)境等領域得到廣泛應用。然而，模型變量和參數(shù)具有很大的隨機性和時空變異性，且這些變量和參數(shù)的變化會造成模擬結(jié)果較大的差異[1]。流域離散是分布式水文模型的基礎，子流域的數(shù)量和面積對模擬結(jié)果影響差異極為顯著，子流域離散方案產(chǎn)生的模擬偏差，可能會超過模型參數(shù)率定所產(chǎn)生的偏差，而且在地形復雜、坡度較大的流域，不同子流域數(shù)量使模擬精度存在不確定性[2]。關于子流域劃分對流域模擬結(jié)果的影響，國內(nèi)外學者進行了相關研究。一些學者認為子流域劃分對流域產(chǎn)流模擬影響較小,而對產(chǎn)沙模擬影響較大[3,4]。郝芳華等和張雪松等的研究表明，當子流域劃分數(shù)量達到一定水平時，增加子流域劃分對模型輸出結(jié)果的影響較小[5,6]。部分學者通過進一步分析，認為在不同集水面積閾值條件下其模擬結(jié)果存在閾值適用范圍的上下限。Tripathi M和Jha M認為對流量、泥沙、NO3-N 和無機磷模擬的最適子流域面積占流域總面積的比例分別為4%、3%、2%和5%[7,8]。同時，基于SWAT模型的子流域劃分方案對模擬結(jié)果的影響和研究尺度有一定的關系。Soni等的研究表明SWAT模型的產(chǎn)流預測對子流域和水文響應單元的數(shù)量具有一定響應，當研究的尺度太大時，流域參數(shù)的集總會影響徑流曲線CN(Curve Number)值的分布，從而導致預測的不準確[9]。雖然國內(nèi)外學者針對分布式水文模型子流域劃分方案進行了相關研究，但由于子流域的數(shù)量和面積對模擬結(jié)果影響差異極為顯著，模擬結(jié)果的不確定因素比較復雜，對于如何準確而有效地確定子流域劃分數(shù)量并沒有統(tǒng)一的結(jié)論和標準，而且在不同流域，結(jié)論并不完全一致。本次以大尺度復雜山區(qū)漳河上游流域為研究對象，分析了在不同水平年流域劃分方案對SWAT模型模擬結(jié)果的影響，對準確模擬該地區(qū)的水沙情況奠定一定基礎，可為該地區(qū)水資源水環(huán)境的進一步研究提供技術支撐。

1 流域概況及模型模擬

1.1 流域概況

本次研究選擇漳河上游流域作為研究對象,漳河上游位于海河流域西南部,具體位置見圖1。流域面積約17 730 km2，占海河流域總面積的6%。流域?qū)贉貛Т箨懶约撅L型氣候，四季分明。全年夏短冬長，多年平均氣溫7.4～10.3 ℃，無霜期170 d左右。流域多年平均降水量為576.9 mm，降水量總的分布趨勢是南部大于北部，迎風坡多于背風坡和盆地中心。流域內(nèi)支流眾多，水系呈扇形分布，分東西兩大支流。東支為清漳，該地區(qū)為石質(zhì)山區(qū)，山高谷深，巖石裸露，坡陡流急，含沙量小；西支為濁漳，該區(qū)域為山丘盆地區(qū)，盆地內(nèi)黃土覆蓋較厚，植被較差，水土流失嚴重，洪水挾帶泥沙較多。

漳河上游流域內(nèi)以山區(qū)地貌為主，地形落差高達2 000多米，東鄰太行山脈，西鄰太岳山脈，中間是長治盆地。河流橫穿太行山脈，大地構(gòu)造處于太行山新華夏系第三隆起帶，由一系列北東向平緩復式群皺組成，巖性地層出露主要有震旦系石英砂巖、泥巖，寒武系紫紅色頁巖，奧陶系石灰?guī)r，河谷兩巖階地由第四系積物組成。

圖1 漳河上游流域位置

1.2 模型數(shù)據(jù)和參數(shù)

模型所需的空間數(shù)據(jù)主要包括數(shù)字高程模型(DEM)、土壤數(shù)據(jù)、土地利用類型數(shù)據(jù)和氣象數(shù)據(jù)。DEM使用國際科學數(shù)據(jù)服務網(wǎng)站提供的分辨率為30 m的數(shù)據(jù)，依照研究區(qū)經(jīng)緯度范圍下載，將多景的數(shù)據(jù)進行融合與拼接，使用數(shù)字化的研究區(qū)邊界裁切出邊界精準的漳河上游流域的DEM圖。土壤數(shù)據(jù)來源于世界糧農(nóng)組織網(wǎng)站的HWSD(Harmonized World Soil Database)數(shù)據(jù)。數(shù)據(jù)格式為柵格數(shù)據(jù)，分辨率為1 000 m，投影為WGS1984，土壤分類系統(tǒng)是FAO90。由于模型采用的是美國土壤分類系統(tǒng)，需要重新建立研究區(qū)的土壤屬性數(shù)據(jù)庫。土地利用類型數(shù)據(jù)由海河流域南部 1∶100 萬土地利用圖片數(shù)字化得到，并根據(jù)土地利用的一級分類，轉(zhuǎn)化為 SWAT 模型能識別的模型代碼，最終劃分為 5種土地利用類型：耕地、林地、草地、水域和其他用地。河網(wǎng)水系圖來源于國家基礎地理中心提供的分辨率為30 m的數(shù)據(jù)。降雨、溫度、相對濕度、風速及太陽輻射等氣象數(shù)據(jù)來源于SWAT 官方網(wǎng)站及中國氣象科學數(shù)據(jù)共享服務網(wǎng)。模型所需基礎數(shù)據(jù)的詳細描述見表1。利用Arcgis軟件對DEM、土壤類型圖和土地利用圖作Albers等積圓錐投影統(tǒng)一地理坐標系WGS1984，并統(tǒng)一生成相同的柵格大小。

表1 模型所需基礎數(shù)據(jù)

根據(jù)郝芳華[5]等的研究，由于現(xiàn)有的確定模型參數(shù)的最優(yōu)化方法對實測資料的依賴性很大，只能反映模擬值與實測值的擬合程度，而不能揭示參數(shù)的物理意義，因此本次直接采用模型默認的參數(shù)值進行模擬。

1.3 子流域劃分方案

本次研究基于最小集水面積( Critical Source Area，CSA)閾值劃分出不同數(shù)量的子流域。最小面積閾值(CSA)是形成河流的最小集水面積，也可以表示為其與流域總面積的百分比，CSA越大子流域劃分數(shù)量越少。

為深入分析流域離散化程度對SWAT模型模擬結(jié)果的影響，本次研究在前人成果的基礎上，結(jié)合研究區(qū)特點選擇6種子流域劃分方法，即CSA分別為0.4%、0.5%、0.7%、1%、2%和10%，并通過確定土地利用和土壤類型面積閾值分別為 5% 和10%來確定子流域中 HRU 的數(shù)量。表 2 列出了6種流域離散化方案的特征參數(shù)。

表2 子流域劃分及流域參數(shù)

根據(jù)CSA值對漳河上游流域進行子流域劃分,獲得6種不同的子流域劃分結(jié)果見圖2。

圖2 流域不同子流域劃分層次圖

定義相對偏差為：

Re=(Vi-Vmin)/Vmin×100%

式中：Vi為不同子流域劃分模擬得到的模擬值；Vmin為子流域數(shù)目最小時得到的模擬結(jié)果。

本次研究選擇的實測資料來自于課題組其他成員的分析成果。王喆[10]利用漳河上游近40年的氣象資料進行面降雨量的趨勢性分析，得出1979-1995年為枯水年，1996年為豐水年，1997-2013年為平水年。本次研究在時間尺度上選擇1995、1996年和1997年的氣象資料作為枯水年，豐水年和平水年的輸入數(shù)據(jù)進行模擬結(jié)果的分析。

子流域劃分方案下年均徑流、泥沙及營養(yǎng)物在平水年的相對偏差(Re)見表3。

表3 不同子流域劃分下年徑流、泥沙及營養(yǎng)物在平水年的相對偏差(Re)

注：Re1、Re2、Re3、Re4分別為年均徑流、泥沙、有機氮和有機磷相對偏差值。

2 結(jié)果與討論

2.1 子流域劃分方案對地表徑流的影響

由圖2可看出，不同子流域劃分結(jié)果最直接的表現(xiàn)是各子流域內(nèi)面積大小的變化及河網(wǎng)詳實程度的描述。由表2可看出，子流域劃分方案對徑流模擬結(jié)果有一定的影響，最大偏差值為11.23%，但和泥沙負荷相比而言，子流域劃分方案對徑流的影響并不大，尤其在子流域個數(shù)為58～111范圍內(nèi)，影響很小，這與前人的研究結(jié)論一致。這主要是因為本文所選用的產(chǎn)流方程為徑流曲線數(shù)法，其中決定產(chǎn)流量的敏感性參數(shù)為徑流曲線數(shù)CN，而不同的劃分方法對流域面積加權CN值的影響比較小。

圖3 不同水平年徑流量隨子流域數(shù)量變化曲線

圖3顯示了在不同水平年徑流量隨子流域數(shù)量變化的趨勢，總體表現(xiàn)為先增大，然后趨于穩(wěn)定，最后呈減小趨勢。具體表現(xiàn)為當子流域數(shù)量從5增加到25時，年徑流量隨之增大，豐水年增加的幅度最大，而枯水年呈略微減小趨勢，表明這兩者之間的子流域劃分并不是最佳的子流域劃分方案。當子流域劃分數(shù)量進一步增加，子流域個數(shù)大于25，小于111，CSA值為流域總面積的0.5%～2%時，年徑流量趨于穩(wěn)定，且不同水平年趨于穩(wěn)定的范圍一致。當子流域劃分數(shù)量超過111個時，徑流量明顯減少，尤其在豐水年減小更明顯，而枯水年呈略微的增加趨勢。這說明在不同水平年徑流量隨子流域劃分數(shù)量變化的趨勢并不完全一致。造成這種情況的原因可能是過多的子流域劃分，形成了較多的過于狹長的子流域，使得模型出現(xiàn)了對水文過程的不真實模擬[11]?？梢?，在漳河上游流域，合理的子流域劃分水平對徑流模擬的影響并不大，但存在一個使徑流模擬趨于穩(wěn)定的子流域劃分的上下限，當?shù)陀谙孪蓿哂谏舷迺r，模擬結(jié)果受劃分數(shù)量的影響較大。

2.2 子流域劃分方案對泥沙的影響

由表3可看出，泥沙負荷對子流域劃分方案的敏感性強于徑流量，最大偏差值高達30.31%。由圖4可看出，在不同水平年泥沙負荷隨子流域劃分數(shù)量的增加逐步呈增加趨勢，這和年均徑流量的變化趨勢并不完全一致。當子流域數(shù)量由5變到25時,泥沙負荷有增加趨勢，當子流域個數(shù)在25～111范圍時，泥沙負荷趨于穩(wěn)定，這和年徑流量趨于穩(wěn)定的范圍一致。但泥沙和年均徑流在不同水平年隨子流域數(shù)量變化的趨勢并不完全相同，當子流域個數(shù)超過模擬結(jié)果趨于穩(wěn)定的上限值111個時，年均徑流量明顯減少，而泥沙量明顯增加，尤其在豐水年表現(xiàn)的更明顯。這可能是因為子流域個數(shù)改變引起流域內(nèi)高程的空間變化，從而影響到流域的坡面坡長和坡度，而坡度和坡長參數(shù)(LS 因子)對于泥沙量是敏感因子，能顯著影響泥沙產(chǎn)量的預測。泥沙演算中的泥沙沉積和降解過程會影響產(chǎn)沙預測的結(jié)果，隨著子流域劃分數(shù)目的減小，流域描述簡化，河網(wǎng)密度減小，從而影響泥沙演算，并使預測的準確性降低[12]。

圖4 不同水平年泥沙隨子流域數(shù)量變化曲線

2.3 子流域劃分方案對營養(yǎng)物的影響

營養(yǎng)物氮、磷的變化原因比較復雜，除了受徑流、沉積負荷的影響外，還受到自身復雜計算算法的影響，而有機態(tài)養(yǎng)分的遷移主要是通過吸附在泥沙等固體上進行的。

隨著子流域劃分數(shù)目的變化，營養(yǎng)物有機氮和有機磷負荷發(fā)生明顯變化。由表3可看出，最粗略的子流域劃分和最精細劃分之間的相對誤差，最大偏差值分別為26.53%和45.30%，而且有機磷的相對偏差值總體都比有機氮的大。 由此可見，子流域劃分對營養(yǎng)物負荷的影響較大，這和前人的研究結(jié)果一致。

由圖5、6可看出，在不同水平年，隨著子流域劃分數(shù)量的增加，流域出口營養(yǎng)物有機氮、磷總體上呈現(xiàn)增加的趨勢，和泥沙負荷隨子流域劃分數(shù)量的變化趨勢類似，營養(yǎng)物負荷隨子流域數(shù)目變化趨于穩(wěn)定的范圍和泥沙負荷的范圍一致，都在25～111。

圖5 不同水平年有機氮隨子流域數(shù)量變化曲線

圖6 不同水平年有機磷隨子流域數(shù)量變化曲線

有機氮、磷負荷與泥沙負荷是成比例的，但是由于有機氮磷的河道傳輸與泥沙的河道傳輸算法的不同，所以出口有機氮磷與泥沙并不一定表現(xiàn)出共同的趨勢。徑流量和泥沙負荷都隨子流域劃分方案的不同而產(chǎn)生一定變化，且各自具有不同的變化趨勢和變化幅度(圖 3、圖 4)，因此營養(yǎng)物負荷量的變化原因很難歸結(jié)為某一簡單的原因。

2.4 最佳子流域劃分方案的確定

對于子流域劃分方案對SWAT模擬結(jié)果的影響，多數(shù)學者認為對徑流量的影響較小，對泥沙和營養(yǎng)物的影響較大，但是如何準確有效地確定子流域劃分數(shù)目至今并沒有統(tǒng)一的標準。一般來講，地形復雜的山區(qū)，由于地形破碎，空間差異大，需要更詳細的子流域劃分，以此充分體現(xiàn)該地區(qū)局部空間的差異。但CSA值過小，劃分過細時會產(chǎn)生較多過于狹長或過小的虛假子流域，導致不真實的水文模擬，而當亞流域劃分數(shù)量較少時，對流域描述不夠充分導致模擬輸出結(jié)果不穩(wěn)定，難以達到理想的預測精度。本次研究在不同子流域劃分方案下的模擬結(jié)果分析表明，在漳河上游流域，無論是徑流、泥沙、還是營養(yǎng)物有機氮磷，對子流域劃分方案趨于穩(wěn)定的范圍一致，都在25～111之間。當子流域的數(shù)量從5個增加到25個時，年均徑流、泥沙和有機氮磷相對來說變化幅度較大，且在不同水平年表現(xiàn)出的變化趨勢并不一致，這表明此區(qū)間并不是最佳的子流域劃分方案。隨著子流域劃分數(shù)量的進一步增加，年均徑流、泥沙和有機氮磷逐漸趨于穩(wěn)定，當超過某一子流域劃分方案時，劃分的HRU數(shù)量急劇增加，計算單元劃分過于精細，可能導致計算結(jié)果失真，同時花費更長的計算時間，降低計算效率。因此，在漳河上游流域，存在上下限閾值111、25，當子流域個數(shù)低于下限閾值，對空間的描述不夠，難以滿足流域規(guī)劃與管理的要求，高于上限閾值，模擬結(jié)果出現(xiàn)失真。通過綜合權衡下墊面特征、與其他數(shù)據(jù)的協(xié)調(diào)及模型運行能力等因素的基礎上，根據(jù)模型模擬結(jié)果，最后確定漳河上游流域最佳子流域劃分個數(shù)為79個，CSA取值為12 411 hm2，僅占流域總面積的0.7%，這說明漳河上游流域內(nèi)復雜的地理空間環(huán)境需要更加詳細的描述，這與實際需求相符。本次研究結(jié)果中子流域平均劃分面積占整個流域面積的比值與以往的研究結(jié)論并不一致，說明利用子流域面積閾值設定不同的子流域劃分方案存在地域差異，需要針對不同流域的實際地形、地貌、海拔等空間屬性對模型進行初步的徑流、泥沙和營養(yǎng)物質(zhì)的敏感性分析，合理確定CSA的取值及子流域的劃分方案，從而減少分布式水文模型模擬的不確定性，為進一步更準確地模擬水文過程奠定基礎。

3 結(jié) 語

本文以中國海河流域西南部的漳河上游流域為研究對象，建立流域SWAT年尺度產(chǎn)流產(chǎn)沙模型，進一步分析不同CSA取值下，子流域劃分數(shù)量對徑流、泥沙及營養(yǎng)物模擬的影響。主要結(jié)論如下：

(1)從總體來看, 子流域劃分層次對模擬結(jié)果的影響存在上限下限兩個閾值25、111，介于兩個閾值之間模擬結(jié)果相對穩(wěn)定，超出上限閾值111，模擬結(jié)果出現(xiàn)失真；低于下限閾值25，對空間的描述不夠，難以滿足流域規(guī)劃與管理的要求。子流域劃分數(shù)量對年均徑流量影響較小，對泥沙和營養(yǎng)物影響較大。

(2)在不同水平年，模擬結(jié)果對子流域劃分數(shù)量的敏感性不一致，在豐水年最敏感，在枯水年最不敏感，建議在徑流量較大的地區(qū)，要合理選擇子流域劃分個數(shù)。

(3)綜合考慮模擬精度和模擬效率，漳河上游流域 SWAT 模擬的最合理子流域劃分水平為79個，CSA值為12 411 hm2，占整個流域面積的0.7%。

□

[1] Manillapalli S. Effect of spatial variability on basin scale modeling[C]∥Third International Conference Workshop on Integrating GIS and Environmental Modeling. Proceedings, Santa Fe, New Mexico, 1996.

[2] Kim J G, Park Y, Yoo D, et al. Development of a SWAT patch for better estimation of sediment yield in steep sloping watersheds [J]. Journal of the American Water Resources Association，2009，45(4):963-972.

[3] FitzHugh TW,Mackay D S. Impacts of input parameter spatial aggregation on an agricultural nonpoint source pollution model[J]. Journal of Hydrology, 2000,(236):35- 53.

[4] 盧文喜，伊燕平，張 蕾，等． 不同亞流域劃分數(shù)量對SWAT 模型模擬結(jié)果的影響[J]．水電能源科學，2010，28(10):23-25.

[5] 郝芳華，張雪松，程紅光，等．分布式水文模型亞流域合理劃分水平芻議[J]．水土保持學報，2003，12(4)：75-77．

[6] 張雪松,郝芳華,程紅光,等.亞流域劃分對分布式水文模型模擬結(jié)果的影響 [J]. 水利學報,2004,(7):119-120.

[7] Tripathi M, Raghuwanshi N, Rao G. Effect of watershed subdivision on simulation of water balance components[J].Hydrological Processes,2005,20(5):1 137-1 151.

[8] Jha M, Gassman P W, Secchi S, et al. Effect of watershed subdivision on swat flow, sediment, and nutrient predictions[J].Journal of the American Water Resources Association,2007,40( 3):811-825.

[9] Soni M. Pradhanang, Russell D. B. Effects of critical source area on sediment yield and streamflow [J].Water and Environment Journal, 2014,28:222-232.

[10] 王 喆．氣候變化和人類活動對漳河上游地表水資源的影響研究[D]．河北邯鄲：河北工程大學城市建設學院，2015．

[11] 耿潤哲，李明濤，王曉燕，等. 基于SWAT模型的流域土地利用格局變化對面源污染的影響[J]．農(nóng)業(yè)工程學報，2015,31(16):241-250.

[12] 張召喜，羅春燕，張敬鎖，等．子流域劃分對農(nóng)業(yè)面源污染模擬結(jié)果的影響[J]．農(nóng)業(yè)環(huán)境科學學報，2012，31(10)：1 986-1 993．