姚邦杰 劉 琦② 任 標(biāo) 王 涵
( ①同濟(jì)大學(xué)土木工程學(xué)院地下建筑與工程系 上海 200092)
( ②巖土及地下工程教育部重點(diǎn)實(shí)驗(yàn)室( 同濟(jì)大學(xué)) 上海 200092)
巖溶石漠化是指在巖溶環(huán)境下,受到人類不合理的社會(huì)活動(dòng)和經(jīng)濟(jì)活動(dòng)等的干擾和破壞造成的土壤嚴(yán)重侵蝕、土壤肥力和生產(chǎn)力大幅下降,并導(dǎo)致基巖大面積裸露且出現(xiàn)類似于荒漠化景觀的生態(tài)環(huán)境問題,它是土地荒漠化的主要類型之一(王世杰等,2003)。
巖溶地區(qū)因?yàn)樗畬?duì)可溶性巖石的侵蝕以及兩者之間的能量物質(zhì)交換形成了該地區(qū)特有的地形地貌形態(tài),這種地形地貌形態(tài)可以分為地下和地表兩個(gè)子系統(tǒng)。其中管道、溶洞、地下暗河等組成了地下子系統(tǒng),而峰叢、溶溝、溶槽、溶隙、巖溶漏斗、天窗、落水洞等是地表子系統(tǒng)的組成部分( Ford,1989; 盧耀如,1999) 。巖溶水徑流在地表發(fā)生并分配,然后在地下子系統(tǒng)中輸送和調(diào)蓄。巖溶系統(tǒng)存在孔隙、裂隙和管道3 種介質(zhì)空間,孔隙流、裂隙流和管道流并存。因此,巖溶地下水在結(jié)構(gòu)上表現(xiàn)為功能上的耗散結(jié)構(gòu)、三維的空間地域結(jié)構(gòu)、二元流場(chǎng)形態(tài)結(jié)構(gòu)還有不均一的雙重含水介質(zhì)結(jié)構(gòu)( White,1969) 。巖溶地下水一般不符合達(dá)西定律,表現(xiàn)為紊流或者混合流( Ford,1989; 石朋等,2012) ,正是由于巖溶水系統(tǒng)循環(huán)演化的特殊性,水不僅會(huì)攜帶土壤發(fā)生地表流失,而且會(huì)發(fā)生地下漏失,導(dǎo)致巖溶石漠化現(xiàn)象突出。巖溶水系統(tǒng)循環(huán)演化受到自然和人為因素控制。自然因素包括降雨量、蒸發(fā)量、年徑流量、地形地貌、地質(zhì)構(gòu)造、巖性、水文地質(zhì)條件等; 人為因素包括土地利用形式、水利水電工程、向河流排放污染水源等( Stringfield et al.,1979; Rosenberry et al.,2008) 。國內(nèi)外學(xué)者進(jìn)行了大量的研究,其研究方法主要包括:室內(nèi)實(shí)驗(yàn)法、野外調(diào)查試驗(yàn)法、理論分析法和數(shù)值模擬法( 史運(yùn)良等,1992; 盧衛(wèi)中,1995; 袁道先等,1996; 王臘春等,2000; 韓寶平等,2004; 章程等,2004; 劉延惠等,2005; 曾科,2012) 。但目前少有對(duì)典型巖溶石漠化地區(qū)的巖溶水系統(tǒng)演化問題進(jìn)行系統(tǒng)性研究,而這一循環(huán)過程恰恰是導(dǎo)致石漠化地區(qū)土壤地表流失和地下漏失方向和趨勢(shì)的前提條件和關(guān)鍵性問題。
本文以典型巖溶石漠化地區(qū)貴州貞豐—關(guān)嶺花江高原峽谷區(qū)為研究區(qū),對(duì)該區(qū)的大氣降水、河水、泉水進(jìn)行水化學(xué)、氫氧同位素分析,研究巖溶水的補(bǔ)、徑、排特征,水化學(xué)特征的變化; 利用Modflow對(duì)研究區(qū)進(jìn)行地下水?dāng)?shù)值模擬,用PEST 法( 參數(shù)評(píng)估) 進(jìn)行模型識(shí)別與驗(yàn)證,通過構(gòu)建研究區(qū)的地下水流場(chǎng),計(jì)算研究區(qū)的降雨入滲量、地下暗河的排泄量、地下水泄流量、水均衡情況等,以定量研究區(qū)內(nèi)的巖溶水系統(tǒng)循環(huán)演化過程,從而為進(jìn)一步研究石漠化地區(qū)土壤地表流失和地下漏失的方向和趨勢(shì)提供相關(guān)理論依據(jù)。
研究區(qū)位于貴州省安順市關(guān)嶺縣的北盤江北岸與打邦河西岸所夾的三角地帶,東經(jīng)105°18'45″~105°52'30″、北緯25°33'45″~25°50'00″,研究區(qū)內(nèi)巖溶地貌發(fā)育,地形高差極大,地形復(fù)雜。研究區(qū)年均降雨量約為1100 mm,年均溫約18 ℃,海拔376 ~1828 m,研究區(qū)河流兩岸高差由于河流的下切作用非常大,是典型的巖溶高原峽谷地貌。該區(qū)石漠化發(fā)育程度強(qiáng)烈,石漠化面積約占研究區(qū)總面積的49%,且主要是中度和強(qiáng)度石漠化,已經(jīng)對(duì)人類的生存發(fā)展造成威脅。
該區(qū)主要出露三疊系、白堊系及第四系地層,其中三疊系地層出露面積占95%以上,碳酸鹽巖出露面積約占85%,根據(jù)巖性、巖層的水理性質(zhì)及地下水賦存狀態(tài),將研究區(qū)地下水劃分為3 類,包括碳酸鹽巖裂隙溶洞水、碎屑巖碳酸鹽巖溶洞裂隙水和層狀巖類裂隙水。研究區(qū)內(nèi)共有主要地下河兩條,巖溶管道總長度約39 km,巖溶大泉20 個(gè)。研究區(qū)主要發(fā)育有4 組裂隙,并且主要發(fā)育1 組北西向的裂隙,巖溶發(fā)育的主控方向即主要受北西向裂隙控制,幾條主要地下暗河走向均為北西向,主要分布于研究區(qū)打邦河西岸,各含水類型及地下河的分布見圖1。
2018 年9 月,在研究區(qū)進(jìn)行野外調(diào)查并采集雨水樣4 組,在采集雨水樣2 d 后采集泉水樣20 組,河水樣2 組( 采樣點(diǎn)見圖1) ,基本覆蓋研究區(qū)泉水出露點(diǎn),采樣后進(jìn)行氫氧同位素和水化學(xué)成分測(cè)定,結(jié)果見表1。
圖1 研究區(qū)巖溶水系統(tǒng)水文地質(zhì)圖Fig. 1 Hydrogeological map of karst groundwater system in the study area
因隨著海拔的增加,氣溫相應(yīng)降低,此時(shí)由于蒸發(fā)分餾作用的影響,大氣降水中氫氧同位素組成會(huì)發(fā)生線型變化,因此研究區(qū)大氣降水中氫氧同位素的組成也具有顯著的高程效應(yīng)( 陳中笑等,2010) ,據(jù)此可推算出地下水的補(bǔ)給高程,它是判斷地下水的補(bǔ)給位置以及補(bǔ)給范圍的基礎(chǔ),巖溶地區(qū)地下水循環(huán)更替較快,該方法更加有效( 黃荷等,2016) 。通過研究區(qū)水樣分析計(jì)算結(jié)果表明,花江雨水中δD值與海拔成正相關(guān)關(guān)系: H=-28.976 δD-998.25( 圖2) 。每個(gè)主要巖溶泉的氫氧同位素值代入研究區(qū)雨水δD 值與降雨高程的關(guān)系公式,計(jì)算出主要巖溶泉的大致補(bǔ)給高程,并根據(jù)結(jié)果標(biāo)注到地圖上( 圖3) 。大部分泉水的補(bǔ)給范圍與巖溶發(fā)育主控方向( 北西—南東) 一致,泉水的補(bǔ)給區(qū)域主要受斷層、地下管道控制,其中16 號(hào)泉水( 地下暗河) 補(bǔ)給范圍最大,從西北邊的上關(guān)鎮(zhèn)到東北部的壩坎周圍均為其匯水范圍。根據(jù)泉水補(bǔ)給范圍推測(cè)出研究區(qū)地下水流向主要為北西—南東向,靠近南部北盤江時(shí)地下水流向轉(zhuǎn)為近北—南向,局部地區(qū)如花江鎮(zhèn)東南部地下水流向因地形高差原因會(huì)有一定變化。
為了深入了解研究區(qū)巖溶水資源轉(zhuǎn)化量的關(guān)系,本文使用MODFLOW 軟件進(jìn)行地下水?dāng)?shù)值模擬,根據(jù)室內(nèi)外資料構(gòu)建地下水流場(chǎng),計(jì)算研究區(qū)的水均衡情況、降雨入滲量、地下河管道的排泄量、地下水泄流量等,以定量研究區(qū)內(nèi)巖溶水轉(zhuǎn)化關(guān)系。
表1 研究區(qū)泉水、地表水、雨水同位素及水化學(xué)分析數(shù)據(jù)Table 1 Spring,surface water,rainwater isotope and water chemical analysis data of the study area
圖2 研究區(qū)雨水δD 值的高程效應(yīng)Fig. 2 Elevation effect of rainfall δD value in the study area
按照研究區(qū)的水文地質(zhì)條件圈定研究區(qū)的邊界,并將研究區(qū)的所有匯水面積圈定為模型范圍,最終模型橫向長31.7 km,縱向長38.3 km。其中大面積發(fā)育的碳酸鹽巖為模型的主要含水介質(zhì)。研究區(qū)的邊界主要分為河流和斷層兩類,北部由兩條大型壓性斷層構(gòu)成邊界,南部由北盤江和打邦河兩條大型河流構(gòu)成研究區(qū)邊界。模型的邊界條件主要分為4 類,第1 類為大型河流構(gòu)成的河流邊界,分為兩段,模型南面為北盤江,東面為北盤江支流打邦河,大部分地下水補(bǔ)給入這兩條河流,因此將這兩條河流概化為模型的河流邊界和流出邊界; 第2 類為地下暗河構(gòu)成的排水溝邊界,模型邊界內(nèi)主要發(fā)育有許凹地下河,運(yùn)用MODFLOW 中的“排水溝”子程序來模擬該條地下河,為流出邊界; 第3 類為模型西南邊海拔最高處的定流量邊界; 第4 類為補(bǔ)給邊界。在穩(wěn)定流計(jì)算時(shí)該地區(qū)的降雨補(bǔ)給量設(shè)定為400 mm·a-1。非穩(wěn)定流計(jì)算時(shí),降雨補(bǔ)給量按氣象站實(shí)測(cè)值的0.33 倍給出。模擬范圍、研究區(qū)邊界及模型邊界條件如圖4 所示。
圖3 研究區(qū)主要巖溶泉推測(cè)補(bǔ)給范圍Fig. 3 Speculation range of major karst springs in the study area
表2 水文地質(zhì)參數(shù)初值Table 2 Initial values of hydrogeological parameters
依據(jù)上述水文地質(zhì)概念模型,可將研究區(qū)的裂隙巖溶地下水流問題概化為非均質(zhì)各向異性介質(zhì)中的三維非穩(wěn)定流問題,研究區(qū)在垂向上分為6 層,含水層均為潛水含水層,水文地質(zhì)參數(shù)初值如表2 所示。對(duì)模擬區(qū)進(jìn)行網(wǎng)格剖分。模型在垂向上除了原有的分區(qū)以外不進(jìn)行更細(xì)致的劃分。模型的網(wǎng)格剖分總體上采用154×159 m2的形式。剖分后活動(dòng)單元格數(shù)約為24 500 個(gè),模擬面積1214.1 km2。以研究區(qū)內(nèi)主要巖溶泉的高程值作為該點(diǎn)地下水位的擬合值,采用MODFLOW 中的PEST 算法進(jìn)行參數(shù)反演。
3.2.1 地下水位變化
研究區(qū)地下水流場(chǎng)模擬結(jié)果如圖5 所示,地下水流向大致呈北西—南東向,分別匯入南部北盤江和東部打邦河,地下水流向與該地區(qū)巖溶發(fā)育主控方向一致,較好反映了研究區(qū)的地下水流場(chǎng)。OW-4觀測(cè)點(diǎn)的水位變化情況( 圖6) ,該區(qū)地下水位年變化幅度較大,季節(jié)性強(qiáng),3 ~5 月為地下水枯水期,8~10 月為地下水豐水期,地下水位變化情況相對(duì)降雨量變化情況有一定的滯后性。
3.2.2 地下暗河流量
圖4 模型范圍、研究區(qū)邊界及模型邊界條件Fig. 4 Model range,study area boundaries,and model boundary conditions
模擬得到地下暗河流量變化情況( 圖7) ,地下暗河流量變化與降雨量變化情況存在一定的滯后性,地下暗河流量隨雨季來臨而上升速度加快,模擬地下暗河流量的回落速度在降雨量減小時(shí)較慢。結(jié)合模擬地下河出口的總排泄流量與實(shí)際觀測(cè)流量的擬合效果分析,模擬雨季最大流量小于實(shí)測(cè)雨季最大流量,而模擬枯季最小流量大于實(shí)測(cè)枯季最小流量,可見模擬流量變化情況比實(shí)際情況稍緩,但誤差不大,模擬過程和實(shí)際條件下地下水系統(tǒng)排泄量的變化特征基本吻合。
3.2.3 水均衡
模擬計(jì)算出的研究區(qū)年均地下水均衡情況( 表3) ,本次模擬期中降雨補(bǔ)給量是主要的補(bǔ)給來源; 側(cè)向補(bǔ)給量是次要補(bǔ)給來源。河流側(cè)向排泄量占排泄量的60.21%,雖然占比較大,但唯一一條地下河的管道排泄量占總排泄量的24.92%,符合巖溶地區(qū)的實(shí)際情況。
將2018 年每個(gè)月研究區(qū)的水均衡情況歸納如圖8 所示。在雨季,大量增加的地下水無法迅速通過地下暗河或河流排泄,其中的大部分被地下含水層儲(chǔ)存,而地下水向地下暗河及河流的排泄量迅速增大。在旱季,因該地區(qū)幾乎沒有大氣降水的補(bǔ)給,地下含水層中儲(chǔ)存的地下水開始發(fā)揮補(bǔ)給功能,而地下水向地下暗河及河流的排泄量迅速減小。一年之中,側(cè)向排泄量和側(cè)向補(bǔ)給量變化較小,地下暗河的排泄量變化幅度顯著大于河流排泄量的變化幅度。
根據(jù)模擬計(jì)算結(jié)果看出,研究區(qū)33%的大氣降水通過表層巖溶裂隙帶的下滲作用轉(zhuǎn)化為地下水,其中24.9%的地下水通過地下暗河排泄、60.2%的地下水通過河流排泄作用進(jìn)入北盤江和打邦河成為地表水,其余的14.9%的地下水側(cè)向排泄成為區(qū)域外的地下水。
通過野外調(diào)查、水樣分析及數(shù)值模擬等方法,對(duì)貴州貞豐—關(guān)嶺花江高原峽谷石漠化治理示范區(qū)巖溶水系統(tǒng)的補(bǔ)、徑、排特征,水化學(xué)特征及水資源量等進(jìn)行分析,得到如下結(jié)論:
圖5 OW-4 觀測(cè)點(diǎn)地下水位模擬值Fig. 5 Groundwater level simulation value of OW-4 observation point
圖6 模擬地下水流場(chǎng)和地下水位擬合點(diǎn)Fig. 6 Simulated groundwater flow field and groundwater level fitting point
表3 年均地下水均衡表Table 3 Annual average groundwater balance table
(1) 研究區(qū)內(nèi)大氣降水的氫氧同位素組成具有顯著的高程效應(yīng),據(jù)此計(jì)算出了研究區(qū)主要巖溶泉的補(bǔ)給高程。綜合水化學(xué)分析推斷出的補(bǔ)給區(qū)地層,以及氫氧同位素分析計(jì)算出的補(bǔ)給高程,得到了研究區(qū)主要巖溶泉的補(bǔ)給范圍。得出研究區(qū)河水上下游的水化學(xué)、氫氧同位素變化與河流兩岸的海拔、巖性的變化高度一致,可見地下水與河水之間聯(lián)系密切。
圖7 地下暗河流量變化模擬情況Fig. 7 Simulation of underground dark river flow changes
圖8 研究區(qū)2018 年水均衡變化情況Fig. 8 Water balance changes in the study area in 2018
( 2) 建立了關(guān)嶺—貞豐花江地區(qū)MODFLOW 數(shù)值模型,根據(jù)實(shí)際情況將含水巖層分為6 層,模擬過程和實(shí)際條件下地下水系統(tǒng)排泄量的變化特征基本吻合,根據(jù)模擬結(jié)果得到,研究區(qū)33%的大氣降水通過表層巖溶裂隙帶的下滲作用轉(zhuǎn)化為地下水,其中24.9%的地下水通過地下暗河排泄、60.2%的地下水通過河流排泄作用進(jìn)入北盤江和打邦河成為地表水,其余的14.9%的地下水側(cè)向排泄成為區(qū)域外的地下水。其結(jié)果符合巖溶地區(qū)的實(shí)際情況,但由于巖溶地區(qū)地下水流情況十分復(fù)雜,且地層資料精度較差、地下水位觀測(cè)數(shù)據(jù)缺乏,地下水的模擬仍有一定誤差。
分析結(jié)果基本反映了實(shí)際巖溶水系統(tǒng)排泄量變化及地下水流向特征,初步闡明了該區(qū)域巖溶水系統(tǒng)循環(huán)演化特征,為巖溶石漠化地區(qū)的土壤地表流失和地下漏失的方向和趨勢(shì)研究及防治提供理論依據(jù)。