強降水條件下豫北平原地下水動態響應研究

朱洪生，王繼華，豆敬峰，朱瑞霞

（1.河南省地質環境監測院，河南鄭州450016；2.河南省地質環境保護重點實驗室，河南鄭州450016；3.平輿縣楊埠鎮農業服務中心，河南平輿463400）

近年來，隨著全球氣候變暖，水循環速度加快，汛期極端強降水事件增多，引起洪澇災害頻發。我國北方地區受地下水資源開采影響，地下水位多呈下降趨勢。降水作為地下水的主要補給來源，影響著地下水資源的動態變化特征，極端強降水的輸入導致地下水產生變化。目前，研究地下水動態與其影響因素的理論方法包括多元回歸分析法[1-2]、時間序列分析預測模型[3-4]、灰色系統模型[5-6]和地質統計學方法[7-8]等。平建華等[9]總結回顧了地下水動態預測模型，將各種模型分為確定性模型和隨機模型，并指出了各種模型的適用條件；陶虹等[10]研究了關中城市群50 a地下水動態變化及其影響因素，指出該地區地下水位整體呈下降趨勢，地下水開采是其主要影響因素；危潤初等[11]采用標準化降水指數研究了黑龍江建三江地區地下水位變化的突變點，指出地下水位下降受水稻種植面積以及降水影響；遲寶明等[12]利用基于遺傳算法的BP神經網絡模型預測了元寶山露天礦區地下水動態特征；張文化等[13]研究了氣候變化與人類活動對石羊河流域地下水動態變化的影響，指出人類活動在地下水動態變化中占主導因素。上述主要是對降水的周期性、頻率，以及降水的變化趨勢、空間分配等的研究，而對于典型區域強降水條件下地下水動態響應的研究較少，單次強降水引起的地下水動態響應和預測研究更少。

豫北平原位于河南省北部，地下水過量開采嚴重，造成地下水埋深很深，包氣帶厚度增大。近些年，該地區汛期強降水頻發，特別是2016年7月，華北地區出現5次明顯的降水過程，河北北部、北京、河南等地出現了50 mm以上的暴雨天氣。筆者利用豫北地區水文地質、地下水位統測資料和水質監測資料等，建立研究區地下水數值模型，對強降水條件下地下水動態響應進行研究。

1 研究區概況

研究區位于黃河以北的河南省北部。地形西高東低，高程42～1 725 m。研究區受地質構造控制明顯，地貌類型分為山高谷深、山勢陡峻雄偉的斷塊山地、堆積地貌，局部發育侵蝕、剝蝕型地貌。

研究地層主要有太古界、中太古界、古生界、中生界與新生界的古近系、新近系和第四系地層。研究區地下水類型可分為松散巖類孔隙水、碎屑巖類裂隙孔隙水和碳酸鹽巖類裂隙巖溶水，其中：松散巖類孔隙水在垂向上劃分為淺層地下水系統和深層地下水系統，在平面上分為漳衛河地下水系統（Ⅰ）、黃河地下水系統（Ⅱ）；碎屑巖類裂隙孔隙水可劃分為安陽丘陵區裂隙孔隙水系統（Ⅰ）和湯陰火龍崗裂隙孔隙水系統（Ⅱ）；碳酸鹽巖類裂隙巖溶水系統可分為焦作九里山泉巖溶水（Ⅰ）、輝縣百泉巖溶水系統（Ⅱ）、鶴壁許家溝泉巖溶洞水系統（Ⅲ）和安陽小南海—珍珠泉巖溶水系統（Ⅳ）。

2 數據來源和研究方法

2.1 數據來源

地下水位監測包括2016年1月—2017年12月豫北地區淺層地下水動態監測點34個，深層地下水監測點4個；水質監測點41個，2016年4月共采集41組水樣，2016年7月、9月10個加密點各采集1次共20組水樣。

2.2 地下水數值模擬

為了研究2016年7—9月強降水對地下水位的影響，建立地下水數值模擬模型，對地下水位變化趨勢進行研究。

2.2.1 水文地質概念模型

模擬區為整個豫北平原，總面積約21 460 km2。模擬計算的目的層為第四系全新統-下更新統孔隙含水層，含水層概化為2層，包括潛水含水層和承壓含水層。模型上部邊界為潛水面，為水量交換邊界，包括降水入滲、潛水蒸發、河渠入滲和灌溉水回滲等，下部多為第三系碎屑巖類，處理為隔水邊界。模型西部以山區與平原自然分界線為邊界，為地下水流入邊界。模型北側以行政區邊界為流量邊界，側向流入流出量根據達西定律計算得出。模型其他邊界為河流邊界，主要由黃河、金堤河和徒駭河組成。地下水接受大氣降水和地下水側向補給，整體流向為由西南流向東北，以蒸發、人工開采和地下水側向流出為主要排泄方式。

2.2.2 地下水流數學模型

模擬區第四系松散沉積物巖性在水平、垂直方向上都有較大變化，在同一點，水平方向滲透系數在各個方向大小一致，而垂直方向滲透系數小于水平方向滲透系數，因此將計算區域含水層概化為非均質各向異性介質。模擬區地下水分為潛水和深層承壓水，因此將研究區地下水流系統概化為非均質各向異性準三維非穩定地下水流系統，可用下面偏微分方程的定解問題來描述：

式中：K為滲透系數；S為自由面以下含水層的貯水率；W為源匯項；μ為潛水含水層在潛水面上的重力給水度；P 為潛水面上的降水入滲和蒸發量等；h（x，y，z，t）為水頭；h0（x，y，z）為初始壓力水頭；h1（x，y，z，t）為第一類邊界上的壓力水頭；f1（x，y，z，t）為第二類邊界上的水分通量；Kn為邊界法線方向的滲透系數；t為時間；Ω 為研究范圍；→n 為有效孔隙度；Γ0為潛水面；Γ1為第一類邊界；Γ2為第二類邊界。

2.2.3 模型的識別和驗證

模擬軟件采用Modflow軟件，初始（2016年6月）流場見圖1，第一含水層滲透系數、給水度分區和第二含水層滲透系數、貯水率分區見圖2。根據模擬區地下水位資料，將2016年6—12月作為模型識別驗證期。驗證期末刻（2016年12月）地下水流場擬合見圖3。驗證期模擬水位與觀測水位擬合較好，說明含水層結構、邊界條件的概化、水文地質參數的選取是合理的，所建立的數學模型能較真實地刻畫研究區地下水系統特征，識別后含水層參數見表1。

圖1 地下水初始流場（單位：m）

圖2 模擬區分區

圖3 2016年12月地下水流場擬合（單位：m）

表1 含水層參數

2.2.4 預測情景

2016年保持降水量、蒸發量、開采量等源匯項不變；其他年份保持降水量總量不變，7—9月增加降水量30%，其他月份降水量減少30%，其他源匯項保持不變。預測未來10、20 a地下水位變化情況，以及典型漏斗區地下水位變化情況。

3 強降水條件下地下水動態響應

3.1 地下水位的動態響應

利用6—8月、6—12月地下水位等值線得到地下水位變幅，將研究區淺層地下水劃分為基本平衡區（-0.5～0.5 m）、緩慢下降區（-0.5～-2.0 m）、急劇下降區（＜-2.0 m）、緩慢上升區（+0.5～2.0 m）和急劇上升區（＞2.0 m）5個水位變幅分區。

2016年6月淺層地下水平均水位為61.11 m，平均變幅0.75 m，最大降幅（-12.50 m）位于濮陽市胡村鄉北豆門村，最大升幅（13.81 m）位于鶴壁市浚縣淇門水文站。2016年6—8月豫北地區淺層地下水位變幅分區分布見圖4，其中：急劇上升區面積2 121.89 km2，占10.65%；緩慢上升區面積6 891.28 km2，占34.58%；基本平衡區面積9 505.54 km2，占47.69%，緩慢下降區面積682.53 km2，占2.90%；急劇下降區面積577.63 km2，占4.18%。與2016年6—12月豫北地區地下水位變幅分區（見圖5）相比較，急劇上升區面積減少200.6 km2，緩慢上升區面積增大1 373.97 km2，基本平衡區面積減少905.65 km2。2016年6—8月淺層地下水位急劇上升區面積與緩慢上升區面積之和比2016年6—12月淺層地下水位急劇上升區面積與緩慢上升區面積之和多1 173.37 km2，說明7—8月強降水過后水位埋深明顯減小。

3.2 地下水水質的動態響應

圖4 豫北地區淺層地下水位變幅分區（2016年6—8月）

圖5 豫北地區淺層地下水位變幅分區（2016年6—12月）

淺層地下水各單項組分變化情況見表2。根據2016年4—7月淺層地下水水質數據分析水化學類型和單個化學組分的變化情況。2016年4月HCO3-Ca·Mg·Na型水占80%；2016年7月HCO3-Ca·Mg·Na型水占70%，且與2016年4月相比井位有所改變。2016年4月、7月HCO3·SO4-Ca·Mg·Na型水都占20%，但2016年7月井位有所改變。SO4·Cl-Ca·Mg·Na型水僅在2016年7月測試中有1眼井為此類型。

單個化學組分變化情況：pH值為 7.3～8.3，與2016年4月相比，7月pH值增高井數為5眼，減少井數為2眼，不變井數為3眼。溶解性總固體含量7月低于4月的有9眼。與地下水質量標準Ⅲ類水相比，2016年4月和7月，氯化物含量均不超標，最高含量為231.49 mg/L；硫酸鹽超標率均為10%；硝酸鹽超標率均為50%，并且是相同井；氨氮均不超標；鐵含量超標率均為60%，最高含量為5.09 mg/L；耗氧量4月超標率為10%，7月不超標。與2016年4月相比，2016年7月指標含量減少的占比相對較多，其原因是，4月為相對枯水期，到7月受降水量增加影響，淺層地下水接受降水補給，地下水受到稀釋，指標含量減小。指標含量增加或不變的淺層地下水，可能是淺層地下水埋深大，降水入滲量小所致。

表2 淺層地下水各單項組分變化情況

3.3 地下水位動態響應預測

通過預測得出10 a和20 a后不同含水層地下水流場，見圖6。研究區整個地下水流場發生了一定程度的變化，整體地下水位有所下降，不同地區水位下降程度不同，但整體地下水流向未發生明顯變化，尤其是深層地下水。對于淺層地下水，在安陽市西北側，10 a后出現一定面積的疏干，但深層水未疏干；在安陽市東側和濮陽市北部，水位下降比較明顯，形成較大的地下水漏斗區；在模擬區南部，水位下降不明顯。研究區整體地下水位的下降說明地下水的補采不平衡，開采量大于補給量。對于當前存在的降落漏斗，如果按照當前開采量繼續開采，則20 a后漏斗范圍將持續擴大，地下水埋深將越來越大。

圖6 含水層地下水水位等值線（單位：m）

3.4 典型漏斗區地下水位的動態響應預測

為了詳細分析漏斗區地下水位變化情況，將安陽市東側和濮陽市大部分區域的漏斗區作為典型漏斗區進行分析。將典型漏斗區作為一個均衡區，同時在漏斗區的上、中、下3個地方安置水位觀察點，分別位于內黃縣、濮陽市和留固鎮，通過分析地下水位隨時間的變化了解漏斗區水位的變化。由圖7可知，不同觀察點淺層地下水位和深層地下水位均出現不同程度的下降，且下降速率不同。濮陽市淺層地下水位20 a下降了1.6 m，下降速率比較平穩；深層地下水位下降了3.3 m，前2 a下降速率較快，后期下降速率比較穩定。內黃縣淺層地下水位20 a變化不大，總體下降幅度不超過1.5 m；深層地下水位下降3.6 m，前期水位下降迅速，后期較為穩定。留固鎮淺層地下水位20 a下降了1.8 m，深層地下水位下降了3.7 m。按照當前開采速率，典型漏斗區淺層和深層地下水位將持續下降。

圖7 不同觀察點地下水位隨時間變化曲線

4 結 論

受強降水影響，研究區地下水位變化可劃分為基本平衡區、緩慢下降區、急劇下降區、緩慢上升區和急劇上升區5個水位變幅分區。受強降水影響，地下水化學類型和水化學組分均發生了變化，水質超標率減小。保持當前地下水開采量不變，降水量總量保持不變，7—9月增加降水量30%，其他月份減少降水量30%情景下，研究區地下水位將下降，不同地區水位下降程度不同，但地下水流向未發生明顯變化，尤其是深層地下水；在滲透性較好的地區，10 a后將出現一定面積的疏干，但深層地下水未疏干，在安陽市東側和濮陽市北部形成地下水漏斗區。典型漏斗區保持現狀開采，20 a后淺層地下水位和深層地下水位均將出現不同程度的下降，濮陽市淺層地下水位下降1.6 m，深層地下水位下降3.3 m；內黃縣淺層地下水位下降幅度不超過1.5 m，深層地下水位下降3.6 m；留固鎮淺層地下水位下降1.8 m，深層地下水位下降3.7 m。