王強(qiáng)民,靳德武,王文科,劉 基,3,楊 建,趙春虎,董興玲,尚宏波

(1.中煤科工集團(tuán)西安研究院有限公司，陜西 西安 710054; 2.陜西省煤礦水害防治技術(shù)重點(diǎn)實(shí)驗(yàn)室，陜西 西安 710054; 3.煤炭科學(xué)研究總院，北京 100013; 4.長安大學(xué) 環(huán)境科學(xué)與工程學(xué)院，陜西 西安 710054; 5.旱區(qū)地下水文與生態(tài)效應(yīng)教育部重點(diǎn)實(shí)驗(yàn)室，陜西 西安 710054)

榆神礦區(qū)是我國重要的煤炭生產(chǎn)基地和高強(qiáng)度開采區(qū)，同時(shí)也是典型的生態(tài)環(huán)境脆弱區(qū)[1]，區(qū)內(nèi)具有供水意義的2個(gè)含水層——薩拉烏蘇組含水層位于煤層之上，含水層底板距開采煤層3～50 m，受煤田開采引起的地面塌陷、地裂縫的影響，含水層結(jié)構(gòu)遭到破壞，地下水漏失，引起區(qū)域地下水位大幅度下降[2]。相對濕潤區(qū)，干旱礦區(qū)地下水和植被關(guān)系密切。當(dāng)前，全球氣候變化和人類強(qiáng)擾動礦區(qū)地下水環(huán)境變化對植被生態(tài)耗水的聯(lián)合影響已成為生態(tài)脆弱礦區(qū)亟待解決的關(guān)鍵科學(xué)問題之一[3]。

目前，很多學(xué)者在我國的干旱礦區(qū)開展了保水采煤的研究[4-6]，保水采煤是指在干旱半干旱地區(qū)煤層開采過程中，通過控制巖層移動維持具有供水意義和生態(tài)價(jià)值含水層水位變化在合理范圍內(nèi)[7]。具體其生態(tài)價(jià)值，是指潛水位埋較小時(shí)，土壤含水率較高，植被生長較好，覆蓋率高;當(dāng)潛水位埋深較大時(shí)，不能對土壤儲水量進(jìn)行有效補(bǔ)充，植被出現(xiàn)退化，導(dǎo)致生態(tài)環(huán)境惡化[8]。王文科、楊澤元等通過野外地下水和生態(tài)調(diào)查，建立了鄂爾多斯盆地地下埋深和植被蓋度的關(guān)系曲線，認(rèn)為植被蓋度和地下水埋深呈負(fù)相關(guān)[9]。同時(shí)，植被的生長也受到氣候變化影響，研究發(fā)現(xiàn)植被覆蓋度隨降雨量增加呈現(xiàn)非線性增加的趨勢[10]。鄂爾多斯盆地部分高強(qiáng)度開采礦區(qū)地下水位呈現(xiàn)出變化幅度大、范圍廣的現(xiàn)象[11]，同時(shí)氣候也呈現(xiàn)出較大變化[12]，2者對當(dāng)?shù)氐纳鷳B(tài)環(huán)境變化都有貢獻(xiàn)，但前人較少研究氣候變化和地下水位對植被生態(tài)耗水的聯(lián)合影響。

以往經(jīng)常采用野外調(diào)查和遙感技術(shù)研究植被和地下水關(guān)系，選取反映植被生長狀況的指標(biāo)大多為植被蓋度或植被覆蓋指數(shù)，選取反映地下水狀況的指標(biāo)為潛水埋深和包氣帶含水率等指標(biāo)，繼而建立植被蓋度(覆蓋指數(shù))與潛水埋深(土壤含水率)的關(guān)系曲線[9,13-14]。事實(shí)上，植被與地下水的關(guān)系不是簡單的線性關(guān)系，涉及降水、蒸發(fā)、土壤質(zhì)地等多個(gè)要素[15]，采用調(diào)查及遙感的方法，難以反映“地下水-土壤-植被-大氣連續(xù)體”之間的水動力過程，且花費(fèi)大量人力和物力，土壤水動力學(xué)模擬方法為研究地下水與植被的關(guān)系提供了一種可行途徑，目前已得到廣泛應(yīng)用[16]。

筆者以位于榆神礦區(qū)的典型植被為研究對象，利用野外調(diào)查、室內(nèi)測試和土壤水動力學(xué)模擬方法，研究不同干旱指數(shù)和地下水位埋深條件下典型植被的耗水狀況，并分析其聯(lián)合影響特征，引入單指數(shù)模型確定干旱礦區(qū)典型植被生態(tài)臨界水位，研究成果對豐富我國干旱礦區(qū)保水采煤和綠色開采技術(shù)內(nèi)涵具有重要意義。

1 研究區(qū)概況

榆神礦區(qū)位于陜西省神木縣、榆陽區(qū)境內(nèi)，其煤炭資源儲量豐富、煤質(zhì)優(yōu)良、地質(zhì)構(gòu)造簡單、開采技術(shù)條件優(yōu)越而引世人矚目，已成為我國重要的能源化工基地。主要地貌類型有風(fēng)沙地貌、黃土地貌、河谷地貌，總的地勢呈西北高而東南低，高程在1 200 m左右。其中，風(fēng)沙區(qū)地貌地下水埋深一般較小，基本在1～5 m。自然植被類型主要以苔草、沙蒿、沙柳等為主。本次研究以該區(qū)域的優(yōu)勢植被沙柳為研究對象，其周邊地下水位埋深在1 m左右，野外測定高度為2.50 m，葉面積指數(shù)為2.48。根系分布特征如圖1所示，沙柳主要根系分布在0～60 cm處，并隨著埋深的增大呈現(xiàn)指數(shù)減小的趨勢。

圖1 沙柳標(biāo)準(zhǔn)化根系分布函數(shù)Fig.1 Normalized root density distribution of Salix psammophila

由于該區(qū)域煤炭產(chǎn)量自20世紀(jì)90年代開始大幅度增長，因此本次研究收集并整理榆林市國家基本氣象站1991—2013年的基本氣象數(shù)據(jù)，包括:降水、蒸發(fā)量、最高氣溫、最低氣溫、相對濕度、平均風(fēng)速和日照時(shí)數(shù)，為模型的上邊界提供所需要的降水和蒸發(fā)數(shù)據(jù)。1991—2013年多年平均降水量為400 mm，多年平均潛在蒸散發(fā)量為1 153 mm，年內(nèi)降水量主要集中在5至10月，超過全年降水量的85%。為系統(tǒng)研究榆神礦區(qū)氣候變化和地下水位雙重因素對植被耗水的影響，筆者引入干旱指數(shù)[17]反映研究區(qū)氣候的整體變化，其定義為

AI=(ET0-P)/ET0

(1)

式中，AI為干旱指數(shù);ET0為潛在蒸散發(fā)量，mm;P為降水量，mm。

2 模型構(gòu)建

2.1 土壤水動力學(xué)模型及參數(shù)

有植被條件下，水文地質(zhì)概念模型如圖2所示，一維非飽和土壤水運(yùn)動的控制方程為

(2)

式中，C(h)為容水度，cm-1;h為土壤水壓力水頭，cm;t為時(shí)間，d;z為垂向坐標(biāo)，cm，向上為正;K(h)為非飽和滲透系數(shù)，cm/d;S(z,t)為植被根系吸水速率。

文中采用van-Genuchten-Mualem模型[18]描述土壤水分特征曲線和滲透系數(shù)曲線:

(3)

(4)

(5)

其中，θr為飽和含水量，cm3/cm3;θs為殘余含水量，cm3/cm3;Se為有效飽和度;Ks為飽和滲透系數(shù)，cm/d;α,m,n為土壤水分特征曲線形狀參數(shù)，其中m=1-n/1，n>1;l為彎曲度參數(shù)。根據(jù)室內(nèi)測試、物理試驗(yàn)及數(shù)值模擬相結(jié)合的方法，確定風(fēng)積沙土壤水分運(yùn)動參數(shù)見表1，該參數(shù)已經(jīng)得到模型校準(zhǔn)和驗(yàn)證，有較高準(zhǔn)確率。

圖2 水文地質(zhì)概念模型Fig.2 Hydrology geology concept model

θrθsα/cm-1Ks/(cm·h-1)nl0.0290.380.0523.8752.320.5

文中采用Feddes模型[20]描述植被根系系數(shù)模型:

S(z,t)=α(h)b(z)Tp

(6)

式中，α(h)為水分脅迫函數(shù);b(z)為根系分布函數(shù)，由野外調(diào)查獲取;Tp為植被潛在蒸騰量，cm/d。

2.2 初始及邊界條件

模擬的初始條件取野外原位試驗(yàn)場不同地下水位埋深條件下的土壤水分的多年觀測平均值，部分層位的土壤含水率通過插值獲取。模型的上邊界選取表層積水的“表面水庫”邊界條件，接受降雨補(bǔ)給和蒸發(fā)排泄，其降水和蒸發(fā)數(shù)據(jù)通過研究區(qū)氣象資料獲取。前人研究結(jié)果表明，當(dāng)?shù)叵滤裆畛^5 m時(shí)，植被和地下水的關(guān)系較小，因此選取初始地下水埋深5 m，并逐步減小(0.5 m/次)至地下水埋深1 m作為本次模擬計(jì)算的下邊界條件。因此，本次土壤水運(yùn)動模擬的邊界條件可總結(jié)為

(7)

式中，θ為初始含水量，cm3/cm3;q0(0,t)為凈入滲速率，即降雨量和蒸發(fā)量之差，cm/d。

綜上所述，本次研究共設(shè)置23種干旱指數(shù)情景(1991—2013年)，地下水埋深共設(shè)置1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5,5.0 m等9種情景。因此，需要搭建207個(gè)不同組合條件下的數(shù)值計(jì)算模型。

3 結(jié)果分析

3.1 干旱指數(shù)對植被耗水的影響

以干旱指數(shù)(AI)為橫坐標(biāo)，植被實(shí)際日蒸騰量(Ta)為縱坐標(biāo)，圖3為地下水埋深2.5 m時(shí)AI與Ta的散點(diǎn)圖?？傮w上，隨著干旱指數(shù)的增大，植被實(shí)際蒸騰量呈現(xiàn)減小趨勢，對上述散點(diǎn)圖進(jìn)行擬合，可以得到以下線性關(guān)系:y=-0.229 7x+0.310 8，該式的擬合相關(guān)系數(shù)達(dá)到0.608 2。

圖3 地下水埋深2.5 m時(shí)AI與Ta的散點(diǎn)圖Fig.3 Scatter diagram between AI and Ta when the water table depth is 2.5 m

同時(shí)，分別對地下水位埋深為1.0,1.5,2.0,2.5,3.0,3.5,4.0,4.5和5.0 m時(shí)的干旱指數(shù)和植被實(shí)際蒸騰量的散點(diǎn)圖進(jìn)行線性擬合，求得不同地下水埋深干旱指數(shù)AI與Ta線性式子的斜率，并畫出其斜率與地下水埋深的關(guān)系曲線圖，如圖4所示。由圖4可以看出，不同地下水埋深條件下，干旱指數(shù)對植被實(shí)際蒸騰量的影響差異較大，當(dāng)?shù)叵滤裆顬?.0 m時(shí)，其斜率為0，表明干旱指數(shù)對植被耗水基本不產(chǎn)生影響;當(dāng)?shù)叵滤宦裆钗挥?.0～2.5 m時(shí)，其斜率絕對值逐漸增大，表明干旱指數(shù)對植被耗水的影響隨著地下水埋深的增大呈現(xiàn)增大趨勢;當(dāng)?shù)叵滤裆畛^2.5 m時(shí)，其斜率絕對值基本不變，說明干旱指數(shù)對植被耗水的影響趨于穩(wěn)定，基本和地下水沒有關(guān)系。由此可以看出，植被耗水并不單一的受到氣候變化的影響，地下水也是影響植被耗水的關(guān)鍵要素。

圖4 不同地下水埋深條件下干旱指數(shù)與植被實(shí)際蒸騰量線性斜率變化特征Fig.4 Slope of the linear equation between AI and Ta at different water table depth

3.2 地下水埋深對植被耗水的影響

篩選出1990—2013年干旱指數(shù)的最小值、中間值和最大值，分別代表豐水年、平水年和枯水年，并繪制不同干旱指數(shù)條件下地下水埋深和植被實(shí)際耗水量的關(guān)系曲線，如圖5所示，圖5可以直觀反映地下水埋深對植被耗水的影響。不同干旱指數(shù)條件下，植被耗水量和地下水埋深呈現(xiàn)出典型的單指數(shù)分布規(guī)律，當(dāng)?shù)叵滤裆顬?.0～2.5 m時(shí)，隨著地下水埋深的增大，植被耗水量急劇減小，當(dāng)?shù)叵滤裆畲笥?.5 m時(shí)，植被耗水量基本趨于穩(wěn)定。另外，當(dāng)?shù)叵滤裆畲笥?.5 m時(shí)，圖5中的植被耗水量隨著干旱指數(shù)的增加呈現(xiàn)明顯的減小趨勢，這說明植被耗水量的多寡并不單純受到地下水埋深的控制，也受到氣候變化的影響。

圖5 不同地下水埋深條件下植被實(shí)際蒸騰量變化規(guī)律Fig.5 Variation trend of Ta under different water table depth

3.3 雙因素聯(lián)合影響植被耗水

圖6反映了干旱指數(shù)和地下水埋深聯(lián)合影響下的植被耗水特征，植被耗水量的高值點(diǎn)主要分布在地下水埋深為1.0～2.0 m處，即圖6中的A區(qū)域，此處植被耗水量的等值線圖基本與橫坐標(biāo)干旱指數(shù)平行，表明此區(qū)域植被耗水和干旱指數(shù)關(guān)系不大，植被生長主要受地下水控制。值得注意的是，植被耗水量的最大值并沒有出現(xiàn)在地下水埋深最小，干旱指數(shù)最小的區(qū)域，而是位于地下水埋深最小，干旱指數(shù)最大的區(qū)域，這是因?yàn)橹脖皇芨邼撍缓蛷?qiáng)降雨量的影響，植被根區(qū)的含水量處于近似飽和狀態(tài)，對植被的生長產(chǎn)生脅迫，因此出現(xiàn)局部區(qū)域植被耗水和干旱指數(shù)呈現(xiàn)正相關(guān)的關(guān)系;當(dāng)?shù)叵滤裆钐幱?.0～2.5 m時(shí)，即圖6中的B區(qū)域，此處植被耗水量的等值線圖即不平行于橫坐標(biāo)干旱指數(shù)，也不平行于縱坐標(biāo)地下水埋深，表明此區(qū)域植被耗水受地下水埋深和干旱指數(shù)的雙重影響;當(dāng)?shù)叵滤裆钐幱?.5～5.0 m時(shí)，即圖6中的C區(qū)域，此處植被耗水量的等值線圖近似平行于縱坐標(biāo)地下水埋深，表明此區(qū)域植被耗水和地下水埋深關(guān)系不大，主要受干旱指數(shù)影響。

圖6 干旱指數(shù)和地下水埋深聯(lián)合影響植被耗水分布特征Fig.6 Distribution of Ta in the coordinates composed of the aridity index and water table depth

3.4 礦區(qū)生態(tài)臨界水位估算

植被實(shí)際耗水量(Ta)和潛在耗水量(Tp)的比值，記為Ta/Tp，被認(rèn)為是衡量植被是否受到水分脅迫的關(guān)鍵指標(biāo)[21]。結(jié)合前人研究并考慮研究區(qū)的干旱氣候[16,22]，本文定義當(dāng)Ta/Tp<0.60時(shí)，植被生長受到水分限制;當(dāng)Ta/Tp<0.30時(shí)，植被出現(xiàn)退化現(xiàn)象。圖7反映了豐水年、平水年和枯水年不同地下水埋深條件下Ta/Tp比值，該圖可以直觀看出，Ta/Tp隨地下水埋深的變化趨勢和Ta的變化規(guī)律基本一致，與地下水埋深的關(guān)系可以用單指數(shù)模型來描述，即:

Ta/Tp=C+Aexp(-BH)

(8)

其中，Ta/Tp為植被實(shí)際蒸騰量與潛在蒸騰量比值;H為地下水埋深，m;A,B，C為擬合參數(shù)，參數(shù)值可以通過Matlab軟件進(jìn)行擬合，其相關(guān)數(shù)值見表2，可以看出擬合曲線相關(guān)系數(shù)R2高達(dá)0.99，因此，單指數(shù)模型可以較好地反映Ta/Tp與地下水埋深的關(guān)系。

圖7 不同地下水埋深條件下植被實(shí)際耗水量和潛在耗水量的比值Fig.7 Variation trend of Ta/Tp at different water table depths

干旱指數(shù)擬合參數(shù)ABCR2H1/mH2/m0.49(豐水年)2.8851.6010.2570.991.332.63 0.65(平水年)3.7301.7880.1780.991.221.91 0.80(枯水年)3.9941.6960.0550.991.171.65 平均值1.242.06

根據(jù)式(8)和擬合出的相關(guān)參數(shù)，可以利用式(9)推求不同干旱指數(shù)條件下的生態(tài)臨界水位:

(9)

其中，Hc為生態(tài)臨界水位，m;A,B，C為擬合參數(shù)，其參數(shù)值見表2;Ta/Tp分別取值0.60和0.30，當(dāng)取值為0.60時(shí)，計(jì)算出的Hc記為H1，定義地下水埋深大于H1時(shí)，植被受水分脅迫;當(dāng)取值為0.30時(shí)，計(jì)算出的Hc記為H2，定義地下水埋深大于H2時(shí)，植被出現(xiàn)退化現(xiàn)象。

計(jì)算出的H1和H2值見表2，程東會等[13]在該地區(qū)通過野外調(diào)查的手段研究了沙柳蓋度和地下水埋深的關(guān)系，認(rèn)為沙柳的適生地下水埋深為1.0～3.0 m，地下水埋深為1.0～1.5 m時(shí)沙柳蓋度最大，和本文研究結(jié)果基本一致，也再次印證了采用土壤水動力學(xué)模擬方法的準(zhǔn)確性。同時(shí)由表2可以看出，植被生態(tài)臨界水位并不是一個(gè)固定的數(shù)值，呈現(xiàn)出干旱指數(shù)越大，生態(tài)臨界水位埋深越小的變化規(guī)律。干旱指數(shù)較大條件下，降水量小，降水滿足不了植被的正常耗水需求，因此植被必須從地下水中吸收水分，從而表現(xiàn)出植被對地下水的依賴性較強(qiáng)，生態(tài)臨界水位埋深較小;相反，干旱指數(shù)較小條件下，降水量大，降水基本可以滿足植被的正常耗水需求，植被對地下水的依賴性較小，生態(tài)臨界水位埋深較大。

4 討 論

在干旱半干旱地區(qū)，地下水不僅是水資源的重要組成部分，同時(shí)也是維持當(dāng)?shù)厣鷳B(tài)安全的關(guān)鍵要素。榆神礦區(qū)位于陜北干旱半干旱區(qū)，該區(qū)域煤炭產(chǎn)量高，煤層埋深小，煤炭開采導(dǎo)致采空區(qū)上覆巖層直至地面，出現(xiàn)垮落帶、彎曲帶和裂縫帶[23]，改變了上覆含水層的結(jié)構(gòu)，影響了地下水的補(bǔ)徑排條件，使得地下水的運(yùn)動由采煤前的橫向運(yùn)動向垂向運(yùn)動過渡，表現(xiàn)為地下水采煤前的基流和潛流排泄變?yōu)榈V坑排水[24]，不可避免的破壞地下水資源，從而出現(xiàn)以采煤工作面為中心的地下水位降落漏斗。范立民等[6]分析了榆神府礦區(qū)高強(qiáng)度采煤對地下水的影響，得出高強(qiáng)度煤炭開采是礦區(qū)地下水位下降的主要驅(qū)動因素，認(rèn)為該礦區(qū)71.5%的地下水位明顯下降區(qū)(下降幅度大于8 m)為由高強(qiáng)度開采導(dǎo)致。但是，馬雄德等的研究結(jié)果認(rèn)為榆神府礦區(qū)植被尚未出現(xiàn)退化，甚至有轉(zhuǎn)好的趨勢[25]。結(jié)合本文研究成果，雖然榆神礦區(qū)地下水位出現(xiàn)大面積下降，但由于目前煤炭開采區(qū)采前地下水位埋深普遍大于2.5 m，影響礦區(qū)生態(tài)環(huán)境的主要控制因素是氣候變化(降水量)，考慮到近年來榆神礦區(qū)降水量有增大趨勢，因此出現(xiàn)“雖然地下水位明顯下降，但是生態(tài)環(huán)境局部轉(zhuǎn)好”的現(xiàn)象。

5 結(jié) 論

(1)干旱礦區(qū)植被生長受干旱指數(shù)和地下水埋深雙重影響，當(dāng)?shù)叵滤裆顬?.0～2.0 m處，植被生長主要受地下水控制;地下水埋深為2.0～2.5 m時(shí)，植被耗水受地下水和干旱指數(shù)的雙重影響;地下水埋深大于2.5 m時(shí)，植被耗水和地下水埋深關(guān)系不大，植被耗水主要受干旱指數(shù)影響。

(2)筆者引入單指數(shù)模型擬合地下水埋深和Ta/Tp的關(guān)系曲線，其相關(guān)系數(shù)高達(dá)0.99;同時(shí)利用單指數(shù)模型反求出枯水年、平水年和豐水年條件下的植被生態(tài)臨界地下水位，認(rèn)為當(dāng)?shù)叵滤宦裆畲笥?.24 m，植被受到水分脅迫，當(dāng)?shù)叵滤裆畲笥?.06 m，植被出現(xiàn)退化現(xiàn)象。

(3)采煤地下水位下降對植被生態(tài)的影響是有一定限度的，只有當(dāng)采前地下水位埋深為1.0～2.5 m時(shí)，地下水位下降才會引發(fā)植被生態(tài)退化。當(dāng)采前地下水位埋深大于2.5 m時(shí)，采煤地下水位下降基本對沙柳的生長不產(chǎn)生影響，此時(shí)植被生態(tài)退化主要受氣候變化影響。

本文以研究區(qū)優(yōu)勢植被沙柳為研究對象，因此估算出的生態(tài)臨界水位也僅針對沙柳或同類灌木，考慮到研究區(qū)有旱柳、小葉楊等根系深度較大的植被，其對地下水的依賴程度更高，生態(tài)臨界地下水埋深也更大。