山西省煤礦區地表水文過程的空間模擬研究

侯 帥，張吳平，賈若男

(山西農業大學資源環境學院，山西 太谷 030801)

0 引 言

長期以來，煤炭資源都是中國工業發展的重要支撐力量，特別是20世紀80年代以來，煤炭資源為我國的經濟建設和發展做出了巨大的貢獻，但長期的煤炭開采造成水循環途徑變化、土地和植被損毀、生態環境退化[1]。煤炭開采改變了地表形態，破壞了地表原有的水文過程，研究煤炭開采對地表水文過程的影響，可為煤礦區的環境保護、生態重建、科學開采等提供重要的理論依據。

自20世紀50年代以來，水文模型的研究取得了巨大的進展，其中分布式水文模型可更好地分析自然和人為因素對流域水環境的影響[2]，通常水文模型對水文過程的模擬主要包括植被截留、土壤入滲、徑流匯流等過程，需要大量的水文、氣象觀測資料，而無觀測站且受采煤擾動的長河流域難以滿足這些要求。

植被截留降雨是水文循環中極其重要的一部分，植被截留受到降雨強度、降雨歷時、植被類型、環境因素等諸多因素的影響，是一個復雜的混合過程[3]。國外的Horton、Rutter和Gash模型的應用較多，國內研究和應用較多的多為半經驗性理論模型，如李崇巍、李玉霞等以儀垂祥模型為基礎，建立了岷江上游林冠截持降雨遙感模型，結果較為理想[4-8]。但對農作物冠層截留的研究相對于林冠截留的研究較少且其結果差異較大[9]。SWAP模型中的作物截留計算采用了由Von Hoyningen-Hune和Braden提出的農田作物冠層截留量的一般方程，并且也被應用到其他研究水分運動的模型中去，如張圣微采用該方法建立了農田水分運動模型的降雨截留子模塊[10]。大多數水文模型中的植被截留過程模擬以經驗模型為主，且沒有區分天然植被和種植植被，對此過程的研究較為簡單。

土壤入滲指降雨或灌溉的水分從地表滲入土壤的運動過程，土壤入滲在水文循環過程中起著重要的紐帶作用。土壤入滲的研究模型主要有經驗模型，通過土壤入滲曲線得到模型參數如Kostiakov模型；半經驗模型，采用簡化的連續方程如Hortan、Holtan等模型；物理模型，基于Darcy定律和Richard方程，具有明確的物理意義，如Green-Ampt、Philip等模型,其中Green-Ampt模型廣泛應用于水土侵蝕、降雨入滲等領域的研究中，如Mein 等將其改進應用于恒定降雨條件下的土壤入滲研究，也被應用到SWAT模型中計算下滲量[11,12]。

降雨經植被截留、土壤入滲后，在地形作用下形成地表徑流。而數字高程模型(DEM)可反映地形特征，借助一定的算法可獲取水系特征，為徑流模擬提供了基礎數據，因此已被應用到SWAT、SHE、TOPMODEL等諸多水文模型中[13]。

筆者以30 m×30 m的柵格為計算單元，基于3S技術構建了長河流域地表水文過程模型，分析了煤炭開采下模型參數的空間分布特征，并模擬了不同降雨強度下的地表水文過程。

1 材料與方法

1.1 研究區概況

山西省長河流域位于晉城市澤州縣(圖1)，地理坐標北緯35°30′～35°38′，東經112°37′～112°46′，包括47個村莊，面積約為113 km2，海拔725～1 167 m，是典型的黃土丘陵區。氣候屬暖溫帶大陸性季風氣候，年平均氣溫11.2 ℃，年平均降雨量576 mm，長河于中部自北向南流過，長河以西分布著年產均在45 萬t以上的12座煤礦，長期的煤炭開采造成土地破壞、水資源污染、空氣污染等一系列問題，已在很大程度上影響了當地人民的生產生活和社會經濟的健康發展。

圖1 研究區地理位置及高程Fig.1 Geographical location and elevation of the study area

1.2 地表水文過程模型

研究區地處黃土高原，進行短時間徑流模擬可以不考慮蒸散發和壤中流的影響，因此本文重點模擬地表水文過程中的植被截留降雨、土壤入滲、徑流匯流過程。地表水文過程模型建立的條件為：以分水嶺劃定研究區邊界，研究區內居民以井水為主要水源，在情景模擬中無外界河流匯入研究區。

本文構建的地表水文過程模型包括植被截留模型、土壤入滲模型、徑流匯流模型，模型的總體結構及各模型間的聯系如圖2所示。由給定的降雨類型，首先根據各類型植被的截留降雨模型可得到截留的降雨量，落到地表的降雨再由土壤入滲模型得到入滲量，最后地表凈雨量經地表匯流過程模擬得到匯流累積量的空間分布。本文首次在長河流域構建地表水文過程模型，由于研究區內缺乏水文觀測資料，模型具體應用于下述降雨情景模擬中。

圖2 地表水文過程模型結構示意圖Fig.2 Schematic of the surface hydrological process model structure

1.2.1 植被截留模型

(1)天然植被截留模型。本部分采用儀垂祥提出的植被截留模型，對于給定的植被類型，植被冠層存在一個最大截留量，當降雨量P(mm)大于最大截留量時，最大截留降雨量PJ(mm)為[7]：

PJ=αVFCLAI

(1)

式中：α為葉表面最大持水深度。mm；VFC為植被覆蓋度；LAI為葉面積指數，m2/m2。

當降雨量小于最大截留量時，植被截留量由植被覆蓋度決定，其截留量PJ為：

PJ=VFCP

(2)

(2)種植植被截留模型。根據SWAP模型中計算作物截留的一般方程來模擬作物的截留量，隨著降雨量的增加，其截留量也逐漸達到最大值aLAI[14]，截留量可由下式計算：

(3)

式中：PJ為降雨截留量，mm；a為經驗系數；b為植被覆蓋度。

1.2.2 土壤入滲模型

采用以下由Mein和Larson將Green-Ampt模型推廣至降雨入滲情況的模型來對入滲過程進行模擬[15]。

假定有穩定雨強p，在降雨的開始階段，降雨全部滲入土壤，當p大于土壤的入滲能力時，地表開始形成積水，由模型可得出積水時的累積入滲量Fp(cm)：

Fp=Sf(θs-θi)/(p/Ks-1)

(4)

由此可得開始積水時間tp=Fp/p，只有當p>Ks時才可能積水，所以整個過程的入滲率f(cm/d)為：

(5)

式中：Ks為飽和導水率，cm/d；Sf為濕潤峰處的土壤水吸力值，cm；θs、θi分別為土壤的飽和含水率和初始含水率，cm3/cm3；F為開始積水以后的累積入滲量，cm。F的計算采用下面修正后的公式：

t>tp

(6)

式中：ts表示假設從t=0時開始積水，到入滲量為Fp時所需的時間，計算方法如下：

(7)

在實際應用中，首先根據降雨強度和降雨時間判斷地表是否可以形成積水，積水前的累積入滲量由式(4)得到，積水后的累積入滲量由式(6)、(7)得出。

1.2.3 徑流匯流模型

建立基于DEM的徑流匯流模型，DEM填洼處理后，采用單流向算法中應用較多的D8法確定水流方向，根據合適的閾值提取出河網，應用Strahler法對河網分級，并將整個流域劃分成若干子流域，從而建立研究區的數字水系[16]。

以植被截留、土壤入滲過程得到的地表凈雨量為賦值柵格，對生成的數字水系賦值可得到地表匯流量的空間分布及各子流域的匯流累積量。

1.3 數據來源與處理

1.3.1 植被截留模型參數

(1)遙感數據。來源于美國地質勘探局(http:∥www.usgs.gov/)的2000-2015年的Landsat影像，首先需要對原始影像進行預處理，主要是輻射校正和幾何校正。

(2)植被類型。采用研究區遙感影像結合當地土地利用規劃圖及《澤州縣林地保護利用規劃》對遙感影像進行解譯，并利用研究區1∶2 000正射影像對地物判別和修正得到。

(3)葉表面最大持水深度。根據各鄉鎮的統計資料及文獻[17]，并結合實地調查確定。

(4)植被覆蓋度。利用獲取的遙感影像來提取歸一化植被指數(NDVI)來估算[18]。

(5)葉面積指數。通過遙感影像結合植被冠層多次散射模型來獲取，具體步驟見文獻[19]。

1.3.2 土壤入滲模型參數

模型的關鍵參數土壤飽和導水率由土壤理化性質通過van Genuchten模型得到，濕潤峰處的土壤水吸力值根據張光輝等的研究得出[20,21]。所需的土壤的各項理化性質通過GPS確定1 km×1 km的采樣網格對土壤采樣后測定，共得到117組采樣點數據。

1.3.3 徑流匯流模型參數

DEM數據來源于澤州縣國土資源局，分辨率為10 m，在ArcGIS中經填洼、確定流向、提取河網、劃分子流域后可生成研究區的數字水系，其中含184個子流域，各徑流于研究區中部匯流，總體從東北流向西南，經流域西南出口流出。

2 結果與討論

2.1 參數的空間分布特征

2.1.1 植被覆蓋度的空間分布特征

植被覆蓋度是影響植被截留的重要因素，也反映了煤炭開采對植被的影響程度。2015年研究區的植被類型的空間分布如圖3，天然植被類型主要為暖溫帶落葉闊葉林和針闊混交林，農作物主要是玉米、小麥。草地、農作物的分布基本相等，灌木林地、闊葉林、針闊混交林在礦區的面積遠大于非礦區。

圖3 研究區植被類型分布圖Fig.3 Distribution map of vegetation types in the study area

將2000-2015年的植被覆蓋度分為5級：①無覆蓋，VFC≤0.1；②低覆蓋，0.1

圖4 2000-2015年研究區植被覆蓋度等級的空間分布圖Fig.4 Spatial distribution of vegetation coverage grades in the study area during 2000-2015

圖5 2000-2015年礦區各植被覆蓋度等級的面積Fig. 5 Area of the vegetation coverage grades in mining area during 2000-2015

表1 2000-2015年礦區植被覆蓋度等級轉移概率矩陣Tab.1 Transfer probability matrix of vegetation coveragegrades in mining area during 2000-2015

2.1.2 土壤飽和導水率的空間分布特征

土壤飽和導水率反映了土壤的滲透能力，是土壤水分運動的重要特征之一。表2是采樣點土壤飽和導水率Ks的基本統計特征，礦區Ks的最小值、最大值、均值均小于非礦區，且Ks屬中等變異，變異性較強，但礦區變異系數大于非礦區，說明受采煤影響，礦區Ks的變化情況大于非礦區。經Kolmogorov-Smirnov法檢驗，表明土壤飽和導水率符合正態分布(P=0.121>0.05)，通過地統計學方法對數據進行半變異函數的擬合，并通過交叉驗證檢驗其合理性和擬合精度，選取決定系數較大(0.568)、殘差較小(1.108×10-3)的最優擬合函數模型為指數模型，根據Kriging插值方法可得到土壤飽和導水率的空間分布(圖6)[23]。研究區土壤飽和導水率的范圍為4.337～139.321 cm/d，西部礦區的平均土壤飽和導水率34.375 cm/d小于東部非礦區的44.426 cm/d，土壤飽和導水率總體上呈現出西部小于東部的趨勢。說明煤炭開采改變了原有的土壤性質，土壤飽和導水率變小。

表2 礦區、非礦區土壤飽和導水率的基本統計特征Tab.2 Basic statistics of soil saturated hydraulic conductivityin mining area and non-mining area

圖6 研究區土壤飽和導水率的空間分布圖Fig.6 Spatial distribution of soil saturated hydraulic conductivity in the study area

2.2 基于多年降雨資料的降雨情景模擬

2.2.1 長河流域降雨趨勢分析

分析來源于澤州縣氣象中心的流域1986-2015年間的降雨資料可知，降雨主要集中在7、8月，兩月降雨量占年總降雨量比例的均值為41.524%。按照中國氣象局降雨量等級的劃分標準，降雨類型主要為小雨(0.1～9.9 mm)和中雨(10～24.9 mm)，分別占67.982%、20.175%，其余大雨(25～49.9 mm)占9.357%，暴雨(50～99.9 mm)、大暴雨(100～249.9 mm)共占2.486%。

2.2.2 情景模擬分析

在上述降雨資料分析的基礎上，以降雨強度分別為2.5、8、16 mm/h，降雨時間為1 h的降雨情景對研究區2015年7月的地表水文過程進行空間模擬。對于有植被覆蓋的柵格，首先根據植被截留模型判斷3種降雨情景下的降雨量是否達到各類型植被的最大截留降雨量：截留量未達最大截留量的柵格，降雨將不會落到地表；超過最大截留量的柵格根據土壤入滲模型中地表積水形成的時間，判斷落到地表的降雨在1 h內是否會形成地表積水，形成積水的柵格再由徑流匯流模型進行匯流累積量的計算。對于無植被覆蓋的柵格，降雨直接落到地表，同理根據土壤入滲模型判斷降雨在1 h內是否形成地表積水，是否進行徑流匯流的計算。計算結果表明：

(1)2.5 mm/h降雨強度下，降雨大部分被植被截留，較少降雨滲入土壤，如圖7(a)、(b)植被冠層平均截留量為0.160 mm，土壤平均入滲量分別為1.955 mm；

(2)8 mm/h降雨強度下，1 h內植被冠層截留量達到最大，其中農作物平均截留量為0.194 mm；闊葉林平均截留量為0.656 mm；針闊混交林平均截留量為0.442 mm；灌木林平均截留量為0.331 mm；草地平均截留量為0.197 mm。由于降雨強度小于土壤飽和入滲率，其余降雨經土壤入滲后難以形成地表徑流，如圖7(c)，(d)植被冠層平均截留量為0.198 mm，土壤平均入滲量分別為6.672 mm；

(3)16 mm/h降雨強度下，1 h內植被冠層截留量達到最大后，其余降雨一部分滲入土壤，一部分形成積水，植被截留量和土壤入滲量的空間分布如圖7(e)、(f)，土壤平均入滲量為13.457 mm。積水繼而形成徑流，經徑流匯流模擬可得到研究區地表徑流匯流量的空間分布(圖8)，匯流量由兩側向中部、由北向南逐漸增加，在流域西南即第184號子流域出口達到最大值。

圖7 不同降雨強度下研究區植被截留量、土壤入滲量的空間分布圖Fig.7 Spatial distribution of vegetation interception and soil infiltration in the study area under different rainfall intensities

圖8 研究區降雨匯流量的空間分布圖Fig.8 Spatial distribution of rainfall accumulation in the study area

3 結 語

(1)充分運用3S技術，通過構建植被截留模型、土壤入滲模型、徑流匯流模型，對長河流域的地表水文過程進行了空間模擬，模型具有明確的物理意義，其中充分考慮了不同植被類型對降雨的截留效應。

(2)2000-2015年礦區植被覆蓋度有所增加，可見隨著退耕還林、土地整治等措施的實施，植被有所恢復，植被的增加提高了植被的截留作用，減弱了地表土壤的侵蝕，提升了森林涵養水源的能力，引起一系列積極的生態效應。

(3)煤炭開采對土壤的破壞短時間內難以恢復，礦區的土壤飽和導水率小于非礦區，土壤滲透能力的改變降低了降雨和灌溉用水的利用效率，地表水流失加劇，生態環境和農業生產受到一定影響。

(4)長河流域在2.5 mm/h降雨強度下，降雨被植被截留；在8 mm/h降雨強度下，植被截留量達到最大后，其余降雨被土壤入滲，難以形成地表徑流；在16 mm/h降雨強度下，降雨會在短時間內形成積水并進行匯流，匯流量于流域西南出口達到最大值。

(5)長河流域礦區的生態治理工作需進一步加強，本研究可為流域的科學管理與生態恢復工作提供科學的參考依據。