鄱陽湖典型洪泛區(qū)地下水數(shù)值模擬研究*

曹思佳，李云良，李寧寧，陳 靜，宋炎炎, 趙貴章，李志萍

(1：華北水利水電大學地球科學與工程學院，鄭州 450045) (2：中國科學院南京地理與湖泊研究所，南京 210008) (3：河海大學水文水資源與水利工程科學國家重點實驗室，南京 210009) (4：江西省地質局工程地質大隊，南昌 330000) (5：江西省水文監(jiān)測中心，南昌 330002)

洪泛濕地是河流與陸地之間典型的水文過渡帶，也是一種受洪水干擾頻繁的特殊下墊面類型，約占全球濕地面積15%[1]。從水循環(huán)角度而言，洪泛濕地是大氣降水、地表徑流以及地下水相互轉化的重要界面，因洪水過程的干擾，地表-地下水之間相互轉化通常具有高度的動態(tài)性和敏感性[2]。因此，洪泛濕地系統(tǒng)具有水文和生態(tài)環(huán)境研究的多重特色和意義，地表和地下水文過程的季節(jié)性變化對水流界面交互、營養(yǎng)物質交換以及沉積作用有重要影響作用，很大程度上促進了洪泛系統(tǒng)物質流、能量流和信息流的轉化和傳遞[3]。諸多研究表明，河流或者湖泊的水文情勢變化以及地表土地利用類型等對洪泛濕地地下水具有重要影響和貢獻，主要體現(xiàn)在地下水位的動態(tài)性和水動力場的變異性以及地下水對濕地的動態(tài)補給模式等[4]。地表和地下水文研究已成為全球洪泛濕地研究的重要問題，地表-地下水交互轉化是生態(tài)系統(tǒng)研究的基礎，是研究水環(huán)境防治的重要組成部分，同時也對洪泛濕地的管理和保護起至關重要的作用[5]。

在我國的長江中下游地區(qū)，以湖泊為核心的多類型復合濕地生態(tài)系統(tǒng)分布頗為廣泛。其中，鄱陽湖及其洪泛濕地系統(tǒng)具有水文和生態(tài)環(huán)境等諸多領域研究的區(qū)位優(yōu)勢，也是國際上公認的重要濕地。鄱陽湖洪泛濕地在維護系統(tǒng)平衡、植被群落分布等方面發(fā)揮著不可替代的作用和意義，也是我國重要的候鳥越冬棲息地[6]。然而，由于近年來人類活動的干預以及氣候條件變化的疊加影響，鄱陽湖水位呈現(xiàn)出洪季偏洪、枯季偏枯的總體態(tài)勢，尤其是鄱陽湖低枯水位的屢次出現(xiàn)，無疑干擾了湖泊洪泛的自然水文條件和生態(tài)狀況，導致洪泛濕地的淹沒動態(tài)、出露面積和地表形態(tài)等諸多因子發(fā)生顯著變化，將會使鄱陽湖洪泛濕地面臨外部復雜環(huán)境條件及水資源轉化的不確定性[7]。因而，以湖泊變化環(huán)境條件為背景，聚焦湖泊洪泛濕地地表-地下水轉化及其水量動態(tài)研究，對區(qū)域水資源優(yōu)化配置、濕地系統(tǒng)演變機制以及濕地生態(tài)系統(tǒng)健康等具有重要理論意義和實際應用價值。

近年來已有學者針對鄱陽湖洪泛濕地系統(tǒng)的地表和地下水方面開展了一些探索性研究工作，取得了一些重要進展和認識。先前研究發(fā)現(xiàn)，鄱陽湖水位變化對洪泛濕地系統(tǒng)的地下水動力場和湖水-地下水轉化關系具有不可忽視的影響作用，同時也對洪泛地表生態(tài)植被等方面具有重要影響[8]。Li等[9]基于野外觀測數(shù)據(jù)和水化學分析，發(fā)現(xiàn)了鄱陽湖與洪泛濕地的地下水之間存在密切轉化關系，但受資料和方法所限，還無法深入認識地下水對湖水變化的動力響應過程以及兩者之間的轉化通量。陳靜等[10]聯(lián)合應用野外觀測資料和能量守恒等計算方法，估算了鄱陽湖洪泛區(qū)典型碟形湖與地下水的交換通量，但由于原位監(jiān)測無法精確刻畫研究區(qū)的空間變異情況以及資料連續(xù)性相對較差，通量結果的估算誤差還無法確定。Zhang等[11]利用氫氧穩(wěn)定同位素示蹤技術分析了鄱陽湖洪泛濕地的土壤-植物-大氣連續(xù)體(SPAC)的水分傳輸，但研究結論僅是定性得出地下水對濕地水流遷移具有一定的促進作用，尚無法定量識別地下水與地表水之間的轉化與動態(tài)交互。Wang等[12]應用地下水數(shù)值模型開展了鄱陽湖贛撫平原區(qū)的地表-地下水相互作用研究，分析了地下水對河湖水體的貢獻作用及對氣候變化和人類活動的響應，但相關工作主要側重于流域下游的贛撫平原地區(qū)，對洪泛區(qū)地下水動力學過程和轉化機制的認識相對不足。上述分析可知，目前圍繞鄱陽湖洪泛系統(tǒng)地表-地下水轉化等方面已取得了一些共性認識，強調了地下水文的動態(tài)性及對濕地生態(tài)的貢獻，但研究方法以原位觀測和示蹤技術為主，開展洪泛區(qū)地下水動力場及其響應過程的數(shù)值模擬研究對深入理解和完整詮釋濕地生態(tài)環(huán)境效應內(nèi)涵具有重要意義。

在當前變化環(huán)境影響下，鄱陽湖洪泛濕地的水文和生態(tài)環(huán)境問題愈發(fā)突出。本文以鄱陽湖典型洪泛濕地為研究區(qū)，基于野外觀測數(shù)據(jù)和資料獲取，應用地下水數(shù)值模型，著重從地下水對湖水的貢獻作用的視角來分析鄱陽湖洪泛區(qū)的水資源轉化問題。本文主要研究目標為：(1)構建基于FEFLOW的洪泛濕地二維地下水動力學模型，分析地下水位的時空動態(tài)變化與特征；(2)定量識別湖泊不同水文時期地下水動力場的演變規(guī)律及對湖水位變化的響應過程；(3)基于地下水均衡分析，揭示不同水源組分對洪泛區(qū)地下水系統(tǒng)的相對貢獻。本研究可為后續(xù)水文地球化學循環(huán)、生態(tài)環(huán)境效應等相關研究提供重要數(shù)據(jù)資料和科學支撐。

1 材料與方法

1.1 研究區(qū)

鄱陽湖位于江西省北部及長江中下游南岸，是一個過水性、吞吐型湖泊，主要承接贛江、撫河、信江、饒河和修水五河來水，經(jīng)調蓄后由北部湖口注入長江干流(圖1a)。鄱陽湖流域南北長約173 km，東西最大寬度約74 km，北部入江水道最窄處寬度約2.8 km，湖岸線總長達1200 km(圖1c)。鄱陽湖流域面積占長江流域面積的9%，平均年徑流量為1450億m2，水量約占長江流域水量的5%[13]。鄱陽湖流域三面環(huán)山，周圍高中間低，地形高程變化為30～2200 m。鄱陽湖地處亞熱帶暖濕季風氣候區(qū)，年內(nèi)降水量分布不均，加上流域五河來水變化，導致其具有高度變異的湖泊水位變化特征[14]。湖泊洪、枯季節(jié)的水位差可達15 m，形成“洪水一片，枯水一線”的獨特洪泛特征。洪水期間湖區(qū)整體被淹沒，平均水深不足6 m，最大水深可達30 m，湖水面積可達3000 km2，而枯水期水域面積則萎縮至不足1000 km2[15]，大部分水體基本萎縮在主河道附近。在鄱陽湖洪泛濕地內(nèi)部，由于水文水動力、泥沙沖淤等原因形成眾多大大小小、形狀各異的洼地，共計約102個，這些洼地被稱為碟形湖[16]。在洪水期，這些碟形湖被洪水淹沒，與主湖融為一體，而在枯水期，湖水位快速下降，這些碟形湖逐漸形成相互獨立的子湖泊水體。碟形湖為洪泛濕地生態(tài)環(huán)境提供了重要水和物質保障。碟形湖附近的植被長勢良好，碟形湖草洲面積占整個濕地草洲面積的23.09%[17]。同時碟形湖中生存著大量的底棲動物，為候鳥越冬提供了物質基礎和適宜的生存條件。鄱陽湖受動態(tài)水情波動的影響，主湖區(qū)和洪泛濕地并沒有明確的界限，但普遍認為常年有水的地方為主湖區(qū)，受主湖區(qū)水位變化影響的周邊廣泛灘地則視為洪泛區(qū)。根據(jù)鄱陽湖洪枯季節(jié)的水面積變化范圍，本文典型洪泛區(qū)主要位于鄱陽湖主河道(星子-都昌-康山)和西側湖區(qū)邊界之間的洲灘濕地，面積約為1646 km2(圖1c)。

圖1 鄱陽湖流域主要水系分布(a)、湖區(qū)周邊流域地層巖性(b)和典型洪泛區(qū)觀測站點分布(c)Fig.1 Major rivers in the Lake Poyang catchment (a), formation lithology around the lake (b) and hydrological gauging stations distributed in the typical floodplain area of the lake (c)

1.2 基礎數(shù)據(jù)與獲取

鄱陽湖湖盆地形高程數(shù)據(jù)原始分辨率為5 m×5 m，主要用來刻畫洪泛區(qū)地形特點及構建地表-地下水數(shù)值模型。星子、都昌和康山水文站的日水位觀測數(shù)據(jù)來源于江西省水文局，主要用來描述湖泊水位的季節(jié)性變化及作為數(shù)值模型的邊界輸入條件。地下水觀測井(1#～5#)的日水位數(shù)據(jù)來源于加拿大Solinst傳感器的野外自動觀測，觀測點主要分布在吳城國家級自然保護區(qū)濕地以及南磯保護區(qū)濕地，該數(shù)據(jù)主要用來反映地下水位的波動狀況及驗證地下水模型。在修水和贛江的下游平原地區(qū)(未控區(qū)間)，鄱陽湖湖區(qū)主要巖性為第四紀松散巖類沉積物，周圍以變質巖為主，零星分布巖漿巖、碳酸巖及少量碎屑巖(圖1b)。未控區(qū)與湖區(qū)西側的交換通量根據(jù)未控區(qū)的地下水位資料和達西定律計算獲取，未控區(qū)地下水位數(shù)據(jù)來源于江西省水文局，通量數(shù)據(jù)用于作為模型的邊界輸入條件。降雨和蒸發(fā)的日觀測數(shù)據(jù)來源于中國科學院南京地理與湖泊研究所鄱陽湖湖泊濕地綜合研究站(星子站附近)，用于作為數(shù)值模型的大氣輸入條件。典型洪泛區(qū)碟形湖與地下水之間的日通量變化數(shù)據(jù)采用先前的估算結果[10]，作為數(shù)值模型的地表源匯條件。文中涉及的一些水文資料的具體觀測點位及分布可參考圖1c?？紤]到鄱陽湖水情變化特征的普適性以及所需不同數(shù)據(jù)資料的完整性，故選取2018年作為研究期來開展正常年份條件下的地表-地下水文過程分析。

1.3 數(shù)值模型介紹

FEFLOW(Finite Element subsurface FLOW system)最初是由德國水資源規(guī)劃與系統(tǒng)研究所WASY公司開發(fā)的地下水有限元數(shù)值模型，目前已成為丹麥水利研究所DHI模型系統(tǒng)的重要組成部分。模型采用有限單元離散技術，應用適應性更好的非結構性網(wǎng)絡，以此實現(xiàn)局部網(wǎng)格加密及其重點區(qū)域的精細化刻畫，進而來靈活應對各種地下水流場模擬與溶質運移模擬中可能出現(xiàn)的復雜的物理過程[18]。FEFLOW可以實現(xiàn)多孔介質達西流、非飽和流、潛水水流模擬和遷移、變密度流和裂隙等諸多實際問題。因此FEFLOW模型已被廣泛應用于地下水動態(tài)預測、地下水資源利用分配、水熱耦合運移、地下水污染運移、由于抽水引發(fā)的地面沉降以及海水入侵方面的研究，并取得了許多重要成果[19-22]。其地下水流運動數(shù)學模型為：

(1)

1.4 地下水模型構建

根據(jù)鄱陽湖洪泛區(qū)DEM地形高程特征，本文共計提取出41124個高程點，采用克里金插值方法，將插值結果作為數(shù)值模型的實際地表高程。結合湖區(qū)地形地貌分布格局，研究區(qū)共計剖分6329個三角形有限單元網(wǎng)格和8680個節(jié)點，三角形邊長變化范圍介于20～2000 m之間，很好地刻畫了洪泛區(qū)的復雜地形特點(圖2a)?？紤]到本文典型洪泛區(qū)實屬鄱陽湖區(qū)的一部分，相對大尺度流域而言，其地質類型和成因相對一致，加上洪泛區(qū)野外資料極其有限，其中一些水文地質參數(shù)未考慮分區(qū)特征(比如降雨入滲系數(shù)、滲透系數(shù)和給水度等)。根據(jù)研究區(qū)周邊的野外鉆孔資料(圖2b)，湖岸帶巖性主要為細砂和粉質黏土，粉質黏土中的粘粒含量較少，與細砂的巖性特征相差不大，且根據(jù)現(xiàn)場抽水試驗結果，細砂和粉質黏土層的滲透系數(shù)較為接近，在地表以下20 m處有一黏土層，透水性較差，其埋深情況相對均勻，可作為含水層底板。因此本文將地下含水層概化成一層，將潛水含水層底板埋深設為地表以下20 m處。研究區(qū)的東側邊界是鄱陽湖主河道，從北部星子延伸至南部康山，主要受湖泊水位變化控制，因此模型東側設置為給定水頭邊界條件，根據(jù)星子、都昌和康山水位觀測資料進行插值并分為3段來分別給定，以體現(xiàn)東部邊界的湖水位空間差異(圖2a)。洪泛區(qū)濕地的西側主要接受修水和贛江兩大區(qū)間的地下水補給作用，因此西側邊界根據(jù)修水和贛江的影響范圍，劃分為兩段并分別設置為給定流量邊界(圖2a)。

圖2 研究區(qū)模型概化示意圖(a)和都昌湖岸鉆孔柱狀圖(b)(鉆孔數(shù)據(jù)來源于江西省地質局)Fig.2 Conceptual map of groundwater flow model (a) and vertical borehole profile near the lakeshore of the Duchang gauging station (b) (data obtained from Jiangxi Geological Bureau)

研究區(qū)主要的源匯項為大氣降水、地表蒸發(fā)及碟形湖入滲補給或排泄。如前所述，本文選取2018平水年作為模擬期，降水蒸發(fā)等氣象條件與多年平均情況相當，具有很好的代表性。根據(jù)先前文獻資料，研究區(qū)的入滲系數(shù)取值范圍為0.06～0.1[23]，空間變化總體上較小，故本文將研究區(qū)的降雨入滲補給系數(shù)給定為0.1，即認為10%的有效降雨入滲補給潛水含水層。參照許秀麗等[24]的已有研究，鄱陽湖洪泛區(qū)濕地的土壤年蒸發(fā)量大概為100～200 mm，約為E601B型蒸發(fā)皿年蒸發(fā)總量的10%，考慮到植物蒸騰作用主要發(fā)生在每年3-5月且主要吸收土壤水，地下水的貢獻比重相對較小。因此，本文以E601B型蒸發(fā)皿數(shù)據(jù)為基礎，乘以折算系數(shù)0.1作為地下水的蒸發(fā)量。相對于鄱陽湖周邊的地層巖性分布而言，本研究所覆蓋的洪泛區(qū)巖性分布整體上較為均一(圖1b)，加上研究尺度相對較小，因此洪泛區(qū)的地下水數(shù)值模擬參數(shù)設置為均質各向同性。根據(jù)2018年鄱陽湖平水年遙感影像圖(Landsat)，對數(shù)值模型計算域內(nèi)的52個碟形湖進行定義和設置(圖2a)。在野外現(xiàn)場調研中，發(fā)現(xiàn)碟形湖底部存在明顯的淤泥質弱透水層，其厚度約為0.3～0.6 m，這直接約束了碟形湖對地下水的補給能力。根據(jù)先前研究結果，碟形湖以低枯水位季節(jié)滲漏補給地下水為主，而高洪水位季節(jié)兩者之間的交換通量較為微弱，其對地下水的最大補給通量約為0.1 m/d[16]。本文將52個碟形湖做同等處理，均將它們作為源匯項參與到模型計算，根據(jù)季節(jié)水文變化，設置0～0.1 m/d的動態(tài)補給通量?？紤]到鄱陽湖洪泛區(qū)實際情況，漲水期湖水主要通過地表漫灘(稱之為地表路徑)以及滲漏補給(稱之為地下路徑)兩種途徑來影響洪泛區(qū)地下水系統(tǒng)。因此，為進一步考慮溢出面對模型的影響，模型選取豐水期和枯水期兩幅遙感影像，通過水面積分布來概化不同時期洪泛區(qū)淹水對模型的影響，將其作為地表漫灘的水位邊界輸入模型(圖3a，b)。對于研究區(qū)的水文地質參數(shù)設定，根據(jù)現(xiàn)場鉆孔資料(圖2b)，盡量消除水文地質條件空間差異帶來的影響，模型將含水層介質主要作為細砂來考慮。本次地下水模型的時間步長設定為1 d，滲透系數(shù)和給水度是本次地下水模型的主要調整參數(shù)。為保證非穩(wěn)定流模擬結果的合理性，模型將多年平均的水文數(shù)據(jù)作為輸入條件，首先開展研究區(qū)的穩(wěn)定流模擬，將模型計算所得的穩(wěn)定流地下水位作為非穩(wěn)定流模型的初始條件。本文采用納什效率系數(shù)(Ens)、確定性系數(shù)(R2)以及均方根誤差(RMSE)對數(shù)值模型結果進行定量評價。

圖3 鄱陽湖洪泛區(qū)枯水期(a)和豐水期(b)水面分布(圖中淺綠色表示地表淹水)Fig.3 Water surface area during the dry (a) and flooding (b) seasons of Lake Poyang (the light green color in the figure indicates inundation distributions)

2 結果分析

2.1 模型驗證與可靠性評價

本研究中，結合地下水位觀測資料，采用手動試錯法對地下水模型的主要參數(shù)進行調整，調整后滲透系數(shù)取值為150 m/d，給水度為0.01。圖4繪制了基于野外5個地下水位觀測井日水位數(shù)據(jù)的模型驗證效果及定量評價結果，其中4個觀測井(1#～4#)位于蚌湖和沙湖濕地周邊灘地，1個觀測井(5#)位于南磯濕地。由圖4a～e可以看出，整個模擬期內(nèi)，地下水位模擬序列與觀測水位序列的年內(nèi)變化趨勢基本一致，尤其體現(xiàn)在地下水位相對較高的時間范圍內(nèi)。然而，模擬序列和觀測序列在冬季等一些低水位時期存在一定的偏差，水位偏差大約為0.3～0.4 m。這可能是因為野外地下水位觀測井受到降水或周邊其它水體匯入的影響，給觀測結果帶來一定的干擾，導致某些時刻地下水位急劇升高或者下降。其中南磯觀測井(5#)冬季模擬水位變化趨勢與其它觀測井存在一定差異，這可能是由于該觀測井位于鄱陽湖上游，其地質條件和水位變化情況不同于下游的幾個觀測井。通過5個觀測點的定量評價結果來看，水位擬合的Ens變化范圍介于0.71～0.92之間，確定性系數(shù)R2介于0.74～0.94之間，均方根誤差RMSE基本小于0.4 m，表明數(shù)值模擬結果較為理想，很好地再現(xiàn)了地下水位的時序動態(tài)特征以及變化幅度(圖4f)?？紤]到本文資料有限，且2018年的地下水模型是一個連續(xù)、完整的自然年模擬，因此認為當前驗證效果表明模型有能力反映地下水的響應變化?？偟膩碚f，本文所構建的地下水模型能夠很好地適應鄱陽湖洪泛區(qū)的水文變化特點，模擬結果的可靠性較高，可以用于洪泛區(qū)地下水位的季節(jié)性變化對外部環(huán)境的綜合響應。

圖4 洪泛區(qū)地下水位模擬效果驗證及評價Fig.4 Model validation and performance of simulated groundwater levels in the Lake Poyang floodplain

2.2 地下水動力場特征

2.2.1 地下水位時空分布 圖5選取了鄱陽湖“枯-漲-豐-退”4個典型水文時期，以此來分析湖泊洪泛區(qū)的地下水位的空間演變規(guī)律。考慮到多因素影響下地下水流場響應過程，本文每個典型期選擇兩個模擬時段加以綜合分析。在每年的枯水期，因鄱陽湖主湖區(qū)水位較低(約8～9 m)，加上研究區(qū)受到修水和贛江的地下水入流等補給，整個洪泛濕地通常保持著較高的地下水位(>11 m)，總體上地下水由洪泛濕地向東北部主湖區(qū)方向流動，即洪泛區(qū)地下水補給鄱陽湖主湖區(qū)(圖5a)。在漲水期，由于鄱陽湖主湖區(qū)水位快速提高，因湖水對周邊地下水的補給強度發(fā)生明顯變化，導致主湖區(qū)鄰近的大多數(shù)區(qū)域地下水位也隨之抬升(約13～15 m)，但同時可見洪泛區(qū)西北部的地下水位相對較低(<12 m)，可知湖水已逐漸向洪泛濕地進行水量傳輸和補給(圖5b)。對于鄱陽湖的高洪水位時期，湖水與整個洪泛濕地融為一體，洪泛區(qū)同樣受到主湖區(qū)水位的強烈影響作用，除了可見鄰近主湖區(qū)的大部分區(qū)域具有較高的地下水位，還可以發(fā)現(xiàn)主湖區(qū)周邊的地下水位要明顯高于其它區(qū)域，即豐水期仍表現(xiàn)為湖泊補給周邊地下水系統(tǒng)(圖5c)。退水時期，鄱陽湖主湖區(qū)水位逐漸降低(<10 m)，但大面積洪泛區(qū)仍保持一定的高地下水位(約11～12 m)，此時地下水由洪泛洲灘濕地迅速向主湖區(qū)排泄，即洪泛地下水補給湖泊(圖5d)。上述分析可得，鄱陽湖洪泛區(qū)的湖水和地下水位具有明顯的季節(jié)性轉化特征，湖水位的波動變化很大程度上決定了主湖區(qū)與周邊地下水之間的動態(tài)補排模式。總結可知，洪泛區(qū)地下水補給湖泊主要發(fā)生在鄱陽湖的低枯水位季節(jié)(例如枯水和退水階段)，而湖泊補給地下水主要發(fā)生在中高水位季節(jié)(漲水和洪水階段)。

圖5 鄱陽湖洪泛區(qū)枯水期(a)、漲水期(b)、豐水期(c)、退水期(d)的地下水位等值線圖(每個水文階段選擇2個代表日期的模擬結果)Fig.5 Contour maps of simulated groundwater levels for dry (a), rising (b), flooding (c), and recession (d) water periods of the Lake Poyang floodplain (simulated results of the two typical dates were selected for different water periods)

此外，通過圖6a、b可以清晰地發(fā)現(xiàn)，南北方向上，洪泛區(qū)地下水位總體上呈現(xiàn)出空間“南高北低”的分布格局，南北區(qū)域的地下水水位差可以達到2～4 m，這種空間格局主要與洪泛區(qū)地形地貌有關。而在東西方向上，地下水位受主湖區(qū)湖泊水位變化控制，漲水和豐水期呈現(xiàn)“東高西低”的變化格局，區(qū)域地下水位差異可以達到2～3 m，而退水和枯水期則呈現(xiàn)出“西高東低”的分布格局，區(qū)域地下水位差異基本小于2 m。通過圖6c洪泛區(qū)部分觀測點的地下水位變化狀況來看，整個洪泛區(qū)的地下水位年內(nèi)波動狀況幾乎一致，且與湖水位變化規(guī)律較為相似?？偟膩碚f，鄱陽湖洪泛區(qū)的地下水位年內(nèi)變幅基本介于2～5 m之間。另外，地下水位變幅較大的地方主要位于主湖區(qū)附近(例如3#和4#點位)，而地下水位變化幅度相對較小的地方則出現(xiàn)在遠離主湖區(qū)的西側灘地(例如1#、2#和5#點位)。上述分析表明，就鄱陽湖洪泛區(qū)濕地而言，地形地貌對整個洪泛區(qū)地下水位分布具有主導作用，但湖泊水位動態(tài)變化卻是一個關鍵的外部驅動要素，形成了地下水位時空響應的差異性。

圖6 鄱陽湖洪泛區(qū)南北和東西斷面的地下水位分布(a、b)、典型點位的地下水位變化(c)和剖面位置(d)Fig.6 Distribution of groundwater levels along the south-north and west-east transections (a-b), groundwater level variations for typical observation sites (c), and the corresponding transection within the Lake Poyang floodplain (d)

2.2.2 地下水流速場 為進一步揭示研究區(qū)的地下水動力場變化情況，本文選取了不同時期的地下水流速(達西流速)分布結果進行比較分析(圖7)。總的來說，整個研究區(qū)的地下水流速一般小于1 m/d, 部分區(qū)域的流速可達2 m/d，但局部地區(qū)因地形起伏變化較大，部分時段可達8～9 m/d。空間上，研究區(qū)北部的地下水流速要明顯大于南部地區(qū)的地下水流速，東部主湖區(qū)附近的地下水流速要明顯大于洪泛區(qū)(圖7a～d)。上述地下水流速的空間分布特征主要與地形地貌及其所導致的區(qū)域地下水力梯度有關，例如研究區(qū)北部的狹窄地形條件以及東部主湖水位對洪泛區(qū)地下含水層的激勵作用。對于鄱陽湖不同水文時期，地下水流速和運動方向仍存在著明顯差異(圖7a～d)?？菟竟?jié)，地下水主要向東北側的主湖區(qū)流動，主湖區(qū)附近的地下水流速較大，洪泛區(qū)濕地的地下水流速較小(<1 m/d)，在靠近西北側邊界附近的地下水流速則十分微弱(圖7a)。漲水季節(jié)，地下水流向發(fā)生明顯轉變，由主湖區(qū)向洪泛洲灘濕地流動，且由于洪泛區(qū)中部地形高程較低，在漲水初期，主湖區(qū)和修水、贛江的地下水均向洪泛區(qū)中部方向流動，整個空間的地下水最大流速可達1 m/d左右(圖7b漲水期選擇1時段)。但隨著主湖區(qū)水位的不斷上漲以及湖水-地下水的頻繁交換，此時主要表現(xiàn)為地下水整體上由主湖區(qū)向洪泛區(qū)西部流動(圖7b漲水期選擇2時段)。在豐水時期，地下水空間流速較漲水期有減小趨勢，整體上小于1 m/d，但地下水流動方向仍為由主湖區(qū)流向洪泛區(qū)，表明了主湖區(qū)水位波動的影響作用。因豐水期湖水淹沒導致整個研究區(qū)相對飽和，大部分區(qū)域的地下水流速低于0.1 m/d(圖7c)。退水時期，由于東側湖水位迅速下降，地下水流速的主要方向轉變?yōu)橛珊榉簠^(qū)流向主湖區(qū)，地下水流速的空間分布較為明顯(圖7d)。通過典型斷面的地下水流線可清晰發(fā)現(xiàn)，在鄱陽湖枯水期和退水期，由于主湖區(qū)水位總體上低于洪泛區(qū)地下水位，地下水主要由洪泛區(qū)流向主湖區(qū)，而漲水期和豐水期，地下水由主湖區(qū)流向洪泛區(qū)(圖8)。盡管如此，在研究區(qū)南部的部分區(qū)域，因地下水位常年低于周邊主湖區(qū)水位(康山站附近)，地下水基本以流向洪泛區(qū)為主。上述結果表明，鄱陽湖洪泛區(qū)地下水流速場在地形地貌和湖泊水位波動的疊加作用下，流速和流向均呈現(xiàn)明顯的季節(jié)性差異，這種差異主要體現(xiàn)在地下水在主湖區(qū)和洪泛區(qū)之間交互作用方式的動態(tài)轉變。

圖7 鄱陽湖洪泛區(qū)枯水期(a)、漲水期(b)、豐水期(c)、退水期(d)的地下水流速變化(每個水文階段選擇2個代表日期的模擬結果)Fig.7 Contour maps of simulated groundwater flow velocities for dry (a), rising (b), flooding (c), and recession (d) water periods of the Lake Poyang floodplain (simulated results of the two typical dates were selected for different water periods)

圖8 鄱陽湖洪泛區(qū)“枯-漲-豐-退”時期典型斷面的地下水流線變化(綠線所示的典型斷面處為流線起點)Fig.8 Changes of the groundwater streamline for selected transects in the Lake Poyang floodplain during the dry, rising, flooding and recession periods (green lines indicate the start points of the groundwater streamline)

2.3 地下水均衡分析

以整個洪泛區(qū)地下水系統(tǒng)為對象，本研究的收入項包括：東側及地表漫灘一類邊界水量補給、降水輸入、碟形湖補給地下水量以及西側二類通量邊界補給。主要支出項包括：東側及地表漫灘一類邊界的排泄水量、蒸發(fā)量以及西側二類邊界的排泄水量。因此，根據(jù)地下水入流量和地下水出流量來計算該地區(qū)地下水的儲量變化，構建基于月尺度的地下水均衡方程為：

ΔS=DBin-DBout+NBin-NBout+Qp+Ql-Qe

(2)

式中，ΔS為地下水儲量的變化量(m3/月)；DBin和DBout分別為一類邊界的流入量和流出量(m3/月)；NBin和NBout分別為二類邊界的流入量和流出量(m3/月)；Qp為降雨入滲補給量(m3/月)；Ql為碟形湖補給地下水量(m3/月)；Qe為潛水蒸發(fā)量(m3/月)。

根據(jù)地下水數(shù)值模擬結果，鄱陽湖洪泛區(qū)地下水均衡信息匯總如表1所示。由表1可得，2018年研究區(qū)的地下水總輸入水量約為1.91×108m3，地下水的總輸出水量約為1.95×108m3，水均衡計算誤差約2%，水量基本保持平衡。根據(jù)表1結果分析，降水和蒸發(fā)是研究區(qū)重要的水量平衡組分，年降水補給量為99.5×106m3，占總補給量的50%以上，地下水的年蒸發(fā)排泄量為140.5×106m3，約占總排泄量的70%以上。碟形湖是洪泛區(qū)地下水系統(tǒng)的一個穩(wěn)定補給水源組分，但碟形湖對地下水的年補給量為5.4×106m3，約占總補給量的3%，表明碟形湖對下伏地下水的貢獻比重相對較小。另外，第一類邊界和第二類邊界的年地下水補給量分別為74.8×106和11.2×106m3，占總補給的比重分別為40%和6%左右，同時可見第一類邊界的年地下水排泄量(24%)也要大于第二類邊界的地下水流出量(4%)，表明東側主湖區(qū)對地下水系統(tǒng)(主湖-洪泛區(qū))的貢獻作用要明顯強于西側的地下水交換(未控區(qū)-洪泛區(qū))?？偟膩碚f，鄱陽湖主湖區(qū)對研究區(qū)地下水平衡的影響要強于修水、贛江等平原地下水等輸入作用。

從表1月尺度水均衡變化結果可見，本研究區(qū)主要在春、夏季(3-7月)接受外部輸入和補給，例如3-5月明顯的降雨輸入(15.3×106～19.1×106m3/月)以及7月主湖區(qū)水量的貢獻作用(16.4×106m3/月)，而其它月份不同輸入組分的貢獻相對較小。地下水的蒸發(fā)排泄主要發(fā)生在夏、秋季(6-10月)，月蒸發(fā)量變化幅度約為13.6×106～20.4×106m3/月之間。地下水向主湖區(qū)的排泄則以秋、冬季為主(特別是12-1月)，排泄量總體上大于6.5×106m3/月。圖9清晰地呈現(xiàn)了洪泛區(qū)地下水系統(tǒng)的主要輸入-輸出條件和月平均水量動態(tài)變化?？傮w而言，降雨、蒸發(fā)和主湖區(qū)水量交換是研究區(qū)地下水均衡月尺度動態(tài)變化的主要組分(圖9a)，研究區(qū)地下水系統(tǒng)在春、夏季以補給狀態(tài)為主，秋、冬季以排泄狀態(tài)為主(圖9b)。

表1 鄱陽湖洪泛區(qū)地下水均衡分析(單位：106 m3/月)*Tab.1 Water budget analysis of the Lake Poyang floodplain groundwater system (unit: 106 m3/month)

圖9 鄱陽湖洪泛區(qū)地下水均衡組分相對貢獻(a)和月平均地下水補(正值)排(負值)量變化(b)Fig.9 Relative contribution of water budget components to the groundwater system (a) and monthly average recharge (positive values) and discharge (negative values) volume (b) for the Lake Poyang floodplain

3 討論

同其它水源類型相比，地下水是一種極其珍貴的隱藏資源，其水儲量要遠大于地表淡水河流和湖泊，但地下水的貢獻作用極易被忽視[25]。洪泛濕地是生物多樣性和生態(tài)系統(tǒng)研究的重要區(qū)域，具有攔蓄洪水、改善水質等多方面作用，但洪泛濕地的結構和功能對水文變化較為敏感，了解洪泛濕地的地表-地下水轉化和其水平衡信息對于探明濕地水文過程、維持生態(tài)平衡至關重要[26]。深入認識地下水的組成情況、水量來源及變化特征是研究洪泛濕地可持續(xù)發(fā)展以及水資源管理的基礎保障，更是水質分析、水污染溯源的重要前提[27]。關于洪泛區(qū)地下水-地表水交互作用及其重要意義，國外學者已針對亞馬遜Amazon湖泊洪泛區(qū)[28]、蘇格蘭Eddleston河流洪泛濕地[29]、澳大利亞Pike和Katarapko河流洪泛區(qū)[30]和Murray河流洪泛濕地[31]等開展了大量探索性工作。因此，洪泛濕地的地下水及其動力學機制研究是一個具有前瞻性和挑戰(zhàn)性的科學問題。

鄱陽湖復雜的水系統(tǒng)特征、多脈沖的外部輸入以及高度動態(tài)的水位生消過程，除了湖泊自身地表水動力過程研究，洪泛區(qū)地表-地下水轉化研究一直是目前眾多學者關注的重點。地表-地下水相互作用關系不僅僅是當前水文學、水文地球化學和生態(tài)水文學的研究熱點，更對其水資源管理、水資源可持續(xù)開發(fā)利用有著重要意義[3,11]。本文通過分析復雜河湖相互作用下地下水年內(nèi)動態(tài)特征，有效驗證了鄱陽湖地下水位與湖水位之間存在轉化機制。研究發(fā)現(xiàn)，鄱陽湖洪泛濕地的地下水位表現(xiàn)出明顯的季節(jié)性變化，模擬結果顯示，地下水位年內(nèi)變幅小于5 m，而地表水年內(nèi)變幅可達9 m以上。由此可見，研究區(qū)地下水年內(nèi)變幅小于地表水，地下水文過程的變化要相對穩(wěn)定。此外，湖水(低水位出現(xiàn)在第50～60天，高水位出現(xiàn)在第200～210天)與地下水(低水位出現(xiàn)在第60～65天，高水位出現(xiàn)在第205～215天)的高低水位出現(xiàn)時間比較接近，滯后期一般小于15天，說明研究區(qū)湖水-地下水文變化的同步性較好。通過進一步分析研究區(qū)地下水與湖水位的相互關系可知，越靠近主湖水體的區(qū)域，湖水位與地下水位的相關系數(shù)越大(可達0.95)，而遠離主湖水體的洪泛區(qū)，地下水與湖水位的相關系數(shù)較小(小于0.5)。相關性分析進一步表明，研究區(qū)南部的地下水位主要受康山段湖水位影響，而北部地區(qū)則主要受都昌段湖水位影響。另外，本文基于單因子分析并采用增加或減小10%、20%偏移量的方法分別對滲透系數(shù)K和給水度μ開展了敏感性分析。分析得知，給水度和滲透系數(shù)均可以引起較大的地下水位變化，說明模型對給水度和滲透系數(shù)的改變都存在著很敏感的響應(圖10)。

圖10 給水度敏感性分析(a)和滲透系數(shù)敏感性分析(b)Fig.10 Sensitivity analysis of the specific yield (a) and the hydraulic conductivity (b)

本文通過FEFLOW數(shù)值模擬分析了鄱陽湖洪泛區(qū)地下水動力場的演變及地下水的轉化通量，提升了對洪泛區(qū)地下水系統(tǒng)的理解，但在地下水模型概化過程中仍存在一定的不足：(1)模型的邊界條件是影響模型模擬精度的重要因素之一。模型邊界主要分為自然邊界和人為邊界兩類。研究區(qū)東側邊界為鄱陽湖主河道，可作為自然邊界考慮，因此將湖泊按照水位邊界進行概化，而研究區(qū)西側很難找到一個完整的自然邊界，且研究區(qū)西側主要接受修水、贛江流域地下水補給影響，但地下水流量較小，邊界位置對模型的模擬以及水均衡的影響相對有限。此外，本研究區(qū)西側邊界為鄱陽湖最大淹沒邊界，只有在最高湖水位才可能淹沒至此，大多情況下該邊界對模擬結果產(chǎn)生的影響有限；(2)碟形湖群是鄱陽湖濕地內(nèi)部的特色水體，模型在考慮碟形湖群對地下水的影響時，將碟形湖的滲漏量作為源匯項進行處理，盡管一定程度上體現(xiàn)了碟形湖群對地下水的影響，但后續(xù)工作應充分耦合地下水-碟形湖群的交互作用；(3)地下水數(shù)值模型十分依賴于研究區(qū)參數(shù)給定，同大尺度流域相比，本研究尺度相對有限且參數(shù)分區(qū)未做考慮。因資料獲取有限，模型將洪泛區(qū)含水層概化為細砂來考慮，但本文模型最終調整的滲透系數(shù)(150 m/d)和給水度(0.01)會與實際背景條件存在一定偏大或偏小的情況，可能很大程度上表征一種砂含礫石的介質特性。實際上，根據(jù)水文地質條件的基本特征，洪泛區(qū)的非均質性是必然存在的，本文空間均質假定會給地下水均衡組分的計算帶來一定的誤差。盡管本次模型驗證取得了理想的模擬效果，但下一步研究應在大量原位觀測和試驗基礎上，細化研究區(qū)含水層結構的空間異質性與水文地質參數(shù)分區(qū)。同時，針對鄱陽湖洪泛區(qū)當前所面臨的外部環(huán)境變化，下一步工作將圍繞湖泊-地下水轉化與生態(tài)水文的響應，應用氣候模式的未來預測結果，考慮鄱陽湖水位極端變化，結合鄱陽湖水利樞紐調度方案，模擬評估不同方案下洪泛區(qū)地下水動力場的演變及其對濕地的潛在影響，旨在為濕地生態(tài)系統(tǒng)健康的維護、濕地碳中和及未來可持續(xù)發(fā)展提供支撐。

4 結論

洪泛區(qū)地下水對地表河湖系統(tǒng)水資源分配以及生態(tài)環(huán)境功能等方面具有重要支撐和貢獻作用。鄱陽湖洪泛區(qū)濕地在長江中下游地區(qū)具有重要區(qū)位優(yōu)勢和研究特色，但當前變化環(huán)境已導致該系統(tǒng)水資源轉化和生態(tài)演化等朝著不確定性方向發(fā)展。本文構建了適應洪泛區(qū)特點的FEFLOW地下水流數(shù)值模擬，重點開展了洪泛區(qū)地表-地下水相互作用及動力學過程研究，綜合考慮了洪泛區(qū)地形地貌以及氣象水文等復雜條件，構建了地下水二維數(shù)值模型，模型驗證的納什效率系數(shù)變化范圍為0.71～0.92，確定性系數(shù)為0.74～0.94，均方根誤差為0.24～0.43 m，表明該模型能很好地模擬鄱陽湖洪泛濕地的地下水動力特征及對外部環(huán)境壓力的動態(tài)響應，并主要得出以下結論：

1)研究區(qū)湖水位的波動變化很大程度上決定了主湖區(qū)與周邊地下水之間的動態(tài)補排模式，即鄱陽湖洪泛區(qū)地下水補給湖泊主要發(fā)生在枯水和退水時期，地下水主要由洪泛區(qū)流向主湖區(qū)。而湖泊補給地下水主要發(fā)生在漲水和高洪水位時期，地下水由主湖區(qū)總體上流向洪泛區(qū)。但在研究區(qū)最南部，因地下水位常年低于周邊主湖區(qū)水位，地下水基本以流向洪泛區(qū)為主。

2)整個研究區(qū)的地下水流速一般小于1～2 m/d，北部地下水流速要明顯大于南部地區(qū)，東部主湖區(qū)附近的地下水流速要明顯大于洪泛區(qū)。研究區(qū)地下水流速場在地形地貌和湖泊水位波動的疊加作用下，呈明顯的季節(jié)性差異?？菟?、退水期和漲水期的地下水流速較大，基本介于1～2 m/d之間，而豐水期地下水流速較小，大部分區(qū)域低于0.1 m/d。

3)地下水均衡分析得出，降雨、蒸發(fā)和主湖區(qū)水量交換是洪泛區(qū)地下水均衡動態(tài)變化的主要組分，年降水占總補給量的50%以上，地下水年蒸發(fā)排泄量約占總排泄量的70%以上，湖水補給地下水約占總補給量的40%，而地下水補給湖水約占總排泄量的24%。此外，鄱陽湖洪泛區(qū)地下水系統(tǒng)主要在春、夏季接受外部降雨輸入和主湖區(qū)補給，而地下水的蒸發(fā)和向湖排泄主要發(fā)生在秋、冬季。

4)鄱陽湖主湖區(qū)對研究區(qū)地下水平衡的影響要強于修水、贛江等平原地下水的輸入作用。地形地貌對洪泛區(qū)地下水位分布以及流速場具有主導作用，但湖泊水位動態(tài)變化卻是一個關鍵的外部驅動要素，形成了地下水位時空響應的差異性以及地下水在主湖區(qū)和洪泛區(qū)之間交互方式的動態(tài)轉變。

致謝：感謝南京大學汪禎宸博士，河海大學劉波老師、王哲博士、殷曉然碩士、姜文瑜碩士在模型構建中給予的幫助。