李鳳麗,徐嘉璐,張 游,韓 宇
(1.水發規劃設計有限公司,濟南 250013;2.山東省湖泊流域管理信息化工程技術研究中心,濟南 250013)
地下水庫起源于地下水人工補給,可以通過天然地下儲水空間攔蓄、調節和利用地下水[1, 2]。1972年,日本松尾氏采用灌漿法在長崎縣野母崎町樺島上建設了世界上第一座地下水庫,隨后又在福井、福岡、沖繩縣等地建設了十余座地下水庫,不僅提高了工、農業供水能力,還可以防止海水入侵[3, 4]。20世紀80年代以來,美國、瑞典、荷蘭和德國等地逐漸開展了人工儲水與回采工程,先將地表水回灌到地下水庫中,再通過排水渠和抽水井提取地下水,對水質進行適當處理后供城市用水[5, 6]。我國第一座地下水庫于1975年在河北省南宮市興建[7],之后北京西郊、山東龍口、大連旅順、廣西、貴州和福建等地[8-12]先后建設了地下水庫,為利用地下水提供了寶貴的經驗。在地下水庫水流數值模擬方面[13,14],藍盈盈[15]運用GMS軟件對洮兒河扇形地地下水庫水流進行了數值模擬,結果表明模型可以用于地下水位的預報;中科院利用回灌試驗和數值模擬研究了北京市西郊地下水庫永定河河道的入滲能力[16];為了研究秦嶺山前洪積扇地下水庫調蓄功能,康華等[17]在太平河洪積扇進行回灌試驗,并使用Visual Modflow進行數值模擬,結果表明適度加強開采可以激發地下水庫更大的調蓄潛力。
威海市區和文登區屬于資源性嚴重缺水的地區,隨著社會經濟的發展,水資源供需矛盾日益加劇,為了緩解這一矛盾,在母豬河下游擬建地下水庫,進行河道攔蓄,與位于其上游的米山水庫進行聯合調度,增加地表水和地下水的水資源利用率。本文結合母豬河下游的水文地質條件,利用Visual Modflow進行數值模擬,對地下水庫進行調蓄分析并選定出較為合適的開采方案。
威海市位于山東半島東端,地處北緯36°41′~ 37°35′、東經121°11′~122°42′。北、東、南三面瀕臨黃海,西與煙臺市接壤。威海市區、文登區水資源有年際變幅大、年內分配不均、地域分布不均、開發難度大、利用率低、人均水資源量少等特點,人均水資源占有量為858 m3,僅為全國人均占有量的32%,屬于資源性嚴重缺水地區。近年來隨著經濟社會的發展,各部門用水量急劇增加,水資源供需矛盾日益嚴重。研究開發新的供水工程,可增加威海市區、文登區供水能力,緩解水資源短缺狀況。
母豬河是威海市文登區第一大河,干流總長65 km,流域面積1 115.18 km2,占文登區總面積的63%,流域圖見圖1。母豬河中上游建有米山水庫,是一座以防洪、灌溉、城市及工業供水為主,兼顧發電和養殖等綜合效益的多年調節大(2)型水庫,控制流域面積440 km2,占母豬河流域面積的39.46%。母豬河地下水水庫擬建在母豬河下游,下游河段河谷寬闊,下游河床、河漫灘上伏6.0~10.0 m厚的第四系覆蓋層,以壤土及中粗砂為主構成,局部有中細砂、細砂,屬中等強透水層,地下水埋深一般為1.0~3.5 m,下伏基巖為二長花崗巖,地下截滲條件較好,水文地質剖面圖見圖2。母豬河地下水庫不僅可以發揮多年調節作用,還能攔蓄地表徑流用以補給地下水,從而增加地下水儲量和開采量。此外,文登區雖然有較多零星水系,但是大多獨立入海,無法攔蓄,使得水資源可利用率降低,且母豬河下游由于地下水長期超采,形成地下水位降落漏斗,在河口誘發海水入侵。興建地下水庫工程,通過“一截一蓄”,進行河道攔蓄,增加地表水的利用量;通過“以豐補歉”,在豐水期蓄積大量的水資源以備枯水期利用,還可以為米山水庫提供應急備用水源,緩解區域水資源供需矛盾;通過“蓄淡壓咸”,防止海水入侵,改善區域生態環境。

圖1 母豬河流域圖(單位:m)Fig.1 The basin map of Muzhu River

圖2 母豬河水文地質剖面圖Fig.2 The hydrogeological profile of Muzhu River
文登區母豬河地下水庫工程分為上級庫工程與下級庫工程,主要建筑物包括地下截滲墻、橡膠壩等。地下水庫庫容和水位參數見表1。上級地下截滲工程位于院東村東側,共長1.9 km,頂高程為9.0 m。下級地下截滲工程位于新建高速路南側,長共1.8 km,頂高程為6.0 m。該工程擬攔蓄母豬河右岸一級支流河水及兩岸的地下水。母豬河上下級庫示意圖見圖3。
考慮到降雨蒸發的時空分布對庫區的影響,為使模擬結果更加符合實際情況,以分水嶺為邊界對庫區的整體匯流區域進行模擬。圖4中紅線為母豬河地下水庫流域分水嶺,黃線為米山水庫流域分水嶺。由于區域②的匯流被米山水庫所截留,因此母豬河地下水庫匯流面積僅為①區域,選取①區域作為模擬范圍。

表1 母豬河地下水庫參數Tab.1 Parameters of underground reservoir of Muzhu River

圖3 母豬河上下級庫示意圖Fig.3 Schematic diagram of the upper and lower reservoirs of Muzhu River

圖4 流域水系圖Fig.4 The drainage map of basins
根據地下水埋藏條件,母豬河地下水庫模擬區含水層結構確定為潛水含水層。綜合考慮各地層的水力聯系程度、含水層滲透性大小以及巖性特征,對含水層結構自上而下進行概化。研究區地表以及下伏土壤多為中粗砂且滲透性相近,本次模擬將其概化為一層,該層滲透性強,富水性好,為母豬河地下水庫主要富水層,底層為強風化層花崗巖,透水率低,概化為不透水層。
通過對母豬河地下水庫地下水系統的分析,發現母豬河地下水庫為河道型地下水庫,其儲水地層均位于母豬河及支流河床上。其邊界條件概化如下:①底邊界:庫區基巖為下元古界荊山群變質巖,透水性較弱,富水性差,隔水性能良好,概化為隔水層,為零通量邊界。②上邊界:上邊界為潛水面,巖性主要為富水性和滲透性強的中粗砂,受大氣降雨、人為開采和蒸發排泄等因素的影響,為變通量邊界。③東部邊界和西部邊界:根據流域分水嶺的劃分東、西部邊界,均為隔水邊界。④流域上游邊界:上游邊界米山水庫以東為分水嶺設置為隔水邊界,米山水庫根據實際放流情況和地下水測流補給程度,設置為變通量邊界條件。⑤流域下游邊界:下游地勢較低為母豬河地下水庫的主要排泄通道,根據水位和地層高程,設置為自由排水邊界條件。
用于地下水流模型的水文地質參數主要分為兩類:第一類是用于計算地下水源匯項的經驗參數,如大氣降水入滲系數、灌溉入滲系數、蒸發系數等;第二類是含水層的水文地質參數,主要包括潛水含水層的滲透系數(K)、給水度(Sy)和有效孔隙度(η)等。在建模工作中,首先需要根據水文地質條件等因素確定參數的分區。經分析發現母豬河地下水庫流域的地貌特征主要分為丘陵和平原兩種,丘陵區巖性主要為下元古界荊山群變質巖,平原區巖性主要為中粗砂,二者參數差別較大。因此,模型分為兩個參數區域。圖5為母豬河地下水庫流域活動單元的三維圖,地勢較高的深色區域為丘陵區,地勢低平的淺色區域為平原區。根據實測鉆孔探測的底板高程和設置截滲墻的影響,確定母豬河地下水庫最大回水面積為圖中黑色線所圍區域,橘色線段為截滲墻所在位置。經過計算模擬范圍丘陵區的面積為138.2 km2,平原區的面積為80.7 km2。
按照參數的分區給定參數初值,通過參數率定得出各參數的分區值。參數分區初始值主要通過野外抽水試驗結合地區經驗值,經過綜合分析確定,具體見表2。
地下水系統的源匯項主要包括補給項與排泄項。補給項主要包括降雨入滲補給、河道側滲補給和山前側滲補給。排泄項主要包括地下水開采、潛水蒸發和下游排泄。

表2 主要分區參數初始值Tab.2 The initial parameters of areas
根據模擬區地下水系統水文地質概念模型,建立了地下水流數學模型。當不考慮弱透水層貯水能力時,淺層潛水系統地下水運動的數學模型概化為三維不穩定流動系統,其數學模型如下:
式中:K為邊界法線方向滲透系數,m/d;Sy為給水度;S1、S2為模擬區域第一類和第二類邊界;H0為含水層初始水頭,m;q為含水層第二類邊界單位面積過水斷面補給流量,m3/d;φ為源匯項強度(包括開采強度、入滲強度、蒸發強度等),1/d;D為滲流區域;n為滲流區邊界的單位外法線方向,矢量;W為源匯項,流進為正,流出為負;t為時間。
此次模擬采用矩形網格對滲流區進行離散化(剖分),將復雜的滲流問題處理成符合剖分單元內規則的簡單滲流問題。網格剖面見圖6。模型共剖分單元格100×100(行×列),單個單元格的尺寸是280 m×254 m,由于模擬區含水層結構為潛水含水層,此次網格剖分在垂直方向上分為1層。

圖6 網格剖分圖Fig.6 The grid sectional map
根據收集的研究區地下水系統資料,本次數值模型的模擬期為1979年1月到1979年12月,將整個模擬期劃分為12個應力期,每個應力期為一個自然月,時間步長為十天。在每個應力期中,所有外部源匯項的強度保持不變。初始條件:根據前期觀測孔的監測資料,插值計算潛水含水層的滲流場,并將其作為研究區的初始流場。模型潛水含水層初始流場見圖7。邊界條件:將所計算出來的各邊界流入流出量輸入到模型之中,通過邊界附近流場的擬合程度,適當調整邊界流入流出量。

圖7 初始流場圖Fig.7 The initial flow field map
選取葛家鎮小英村東21號井、澤頭鎮獸醫站院內25號井、葛家鎮糧食所院內26號井、葛家鎮鋪集固恒機械廠招待所院27號井和澤頭鎮道口村委院內東南角48號井共五眼井進行參數率定,位置分布見圖3。采用2016年逐日地下水位資料實測值與模型模擬值進行參數率定,結果表明實測值與模擬值擬合度較好,相關系數均大于95%。運用pest參數迭代調整參數,調整后的平原區中粗砂滲透系數K值為38.5 m/d,而丘陵區變質巖滲透系數K值為0.012 m/d 。平原區中粗砂的給水度為0.21,丘陵區變質巖給水度0.058。
將調整好的參數輸入到模型中運算,輸出結果與2017年的實測地下水位數據進行驗證,由于數據資料不全,只選取21號,26號和48號三眼井進行模型驗證,模擬值與實測值相關系數為98.7%,表明模型能夠較好的反應實際情況,精度較高。圖8和圖9分別是模型進行參數率定和驗證時48號井水位擬合圖。

圖8 2016年48號井水位實測值和模擬值擬合圖Fig.8 Fitting curve of measured and simulated water level of well 48 in 2016

圖9 2017年48號井水位實測值和模擬值擬合圖Fig.9 Fitting curve of measured and simulated water level of well 48 in 2017
根據建立的模型,結合水量調算,計算出母豬河地下水庫豐水年、平水年、枯水年和特枯年的最大調蓄水量,結果見表3。除利用模型計算最大調蓄水量,還利用水量平衡方程對地下水庫進行典型年調算,調算結果豐、平、枯和特枯年的分別是2 575.94、1 969.87、1 659.30和890.71 萬m3。通過對比可知,模型模擬結果和典型年調算結果相差較小,豐水年、平水年、枯水年結果相差均在10%之內,特枯年結果相差12.13%,再次驗證所建模型的合理性。從水量上分析,母豬河地下水庫建成后,有較大的開發潛力;從含水層水力條件及其富水性特征上分析,潛水含水層可調蓄空間適用性大;從長遠角度考慮,合理的調蓄、優化儲水用水方案是合理利用地下水資源的關鍵。
母豬河地下水庫建成后,擬將蓄積的地下水抽取并泵送至上游米山水庫,取代部分外調水源,優先利用當地優質的地下水資源。為合理開采地下水,結合上游米山水庫的放水情況和母豬河地下水庫的開采條件,初步擬定12種開采方案,見表4,并通過模型分析各方案水位變化情況和疏干情況,對比每種方案的優缺點。

表3 地下水庫調蓄計算表 萬m3

表4 不同開采方案分析結果Tab.4 The results of different scenarios of abstracting groundwater
對各方案水位變化進行分析,從模型的穩定性來看,流量為0.01 m3/s的方案比0.015 m3/s的方案更加穩定;但是在開采總量上,0.01 m3/s流量條件下的6種方案中只有55眼井全年開采的方案接近最大調蓄水量,同等條件下0.015 m3/s流量的開采方案能獲得較大的開采量;從工程經濟上分析,55眼井的方案工程量大、費用高,而27眼井又很難滿足開采需求。結合初步擬定的12種開采方案的分析結果,經過多次計算、調試和優化,最終確定36眼開采井(下級庫23眼、上級庫13眼),流量為0.015 m3/s且全年開采的方案,較為滿足工程實施和開采量的需求。開采井在模型中的分布如圖10所示,在開采井的位置上布設觀測井,觀測水位的實時變化。

圖10 井的位置分布Fig.10 The location of wells
開采井設置完畢,模擬枯水年條件下,開采地下水對地下水位的影響情況,模擬時段末模擬區地下水降深情況如圖11所示。模型運行過程中大多數井水位比較穩定,表明開采量和補給量達到平衡狀態,能形成穩定的降落漏斗,在當前開采量下可以持續抽取地下水。但一些離河床較遠的井(如5號井),水位有較大幅度的降低,如果持續抽水可能會導致單元格水位疏干。因此,在實際布井的過程中,應根據實際情況盡可能依河而建,且盡量避開村莊和農田。

圖11 模擬區地下水降深圖Fig.11 The drawdown of simulation zone
本文通過分析和計算,可以得出以下結論。
(1)模型參數率定時,實測值與模擬值相關系數均大于95%,模型驗證時,模擬值與實測值相關系數為98.7%,表明模型能夠較好地反映實際情況。
(2)通過所建模型對地下水庫進行調蓄分析,結果表明豐、平、枯和特枯年的最大調蓄水量分別是2 392.8、2 056.5、1 743.3和1 012.8 萬m3。母豬河地下水庫建成后,有較大的開發潛力,潛水含水層有較大的可調蓄空間。
(3)經過多次計算和調試,發現36眼開采井(下級庫23眼、上級庫13眼),流量為0.015 m3/s且全年開采的方案,能同時滿足工程要求和開采量的需求,為較為合適的開采方案。
本次研究成果可為威海市母豬河地下水庫建設和地下水資源開發利用提供參考,同時,也可以為我國地下水庫建設提供借鑒和思路。