靖安水源地地下水流模擬與預(yù)測(cè)

范珊珊

(河北省地礦局秦皇島礦產(chǎn)水文工程地質(zhì)大隊(duì)，河北 秦皇島066001)

昌黎縣是河北省秦皇島市下屬的四縣之一，近年來，隨著社會(huì)經(jīng)濟(jì)的飛速發(fā)展，人口規(guī)模不斷擴(kuò)大，人民生活水平不斷提高，用水需求量逐年增加。為了滿足日益增長(zhǎng)的用水需求，擬在靖安鎮(zhèn)境內(nèi)建立新的地下水水源地。

靖安水源地位于灤河沖積扇東部的中上部，為傍河取水水源地，該區(qū)淺層含水層儲(chǔ)水空間大，調(diào)蓄能力強(qiáng)，遠(yuǎn)離隔水邊界，補(bǔ)給條件良好，附近多為非耕地，開采對(duì)農(nóng)田用水影響較小。研究區(qū)屬于暖溫帶半濕潤(rùn)季風(fēng)型大陸性氣候，四季分明，冬無嚴(yán)寒、夏無酷暑，日照充足，多年平均降水量為636.9 mm，年內(nèi)降水不均，7～9月降水量可占全年的 70%左右。利用GMS軟件對(duì)該水源地進(jìn)行地下水流數(shù)值模擬，按設(shè)計(jì)的開采方案進(jìn)行模擬預(yù)測(cè)，從而確定該水源地的允許開采量，為水源地的規(guī)劃建立提供依據(jù)。

1 地下水流模型

1.1 水文地質(zhì)概念模型

靖安水源地的研究區(qū)范圍為北至灤河出山口、安山一線，西南部以灤河為界，東以安山、新集一線為界，總面積305.76 km2。根據(jù)地層巖性、含水介質(zhì)特征及其間水力聯(lián)系等因素，將研究區(qū)的含水層劃分為四個(gè)含水組，即Ⅰ、Ⅱ、Ⅲ、Ⅳ含水組，它們的地層大體相當(dāng)于 Q4、Q3、Q2、Q1，其中以Ⅰ、Ⅱ含水組富水性最強(qiáng)且補(bǔ)給條件好，為研究的主要目的層。Ⅱ、Ⅲ含水組之間沉積有累計(jì)厚度達(dá)6～18 m的粘土、淤泥及亞粘土，使得Ⅱ、Ⅲ含水組之間的水力聯(lián)系微弱。灤河在Ⅰ含水組上部流過，屬未切穿含水層的非完整河。河床及漫灘上沉積有薄層亞粘土、粉細(xì)砂層，河水通過該層與Ⅰ含水組有較強(qiáng)的水力聯(lián)系。

圖1 研究區(qū)邊界條件示意圖

由于Ⅰ、Ⅱ含水組之間沒有連續(xù)的隔水層，將Ⅰ、Ⅱ含水組概化為淺層潛水含水層，所以潛水面為系統(tǒng)的上邊界，通過該邊界，系統(tǒng)與外界發(fā)生垂向水量交換，如大氣降水入滲補(bǔ)給以及蒸發(fā)排泄等;淺層含水層和下部含水組之間有較厚的粘土層，所以淺層含水層的底板概化為隔水邊界。研究區(qū)側(cè)向邊界包括流量邊界和水頭邊界(見圖1)，其中西南部以灤河為界，概化為水頭邊界;北部為側(cè)向補(bǔ)給邊界;東部安山鎮(zhèn)至新集為側(cè)向排泄邊界，其它為補(bǔ)給邊界。研究區(qū)地下水補(bǔ)給項(xiàng)包括大氣降水入滲補(bǔ)給、側(cè)向徑流補(bǔ)給和灌溉回歸補(bǔ)給等，排泄項(xiàng)包括人工開采、潛水蒸發(fā)和側(cè)向流出。

1.2 數(shù)學(xué)模型

由于含水層中的地下水流以水平運(yùn)動(dòng)為主，將研究區(qū)的地下水流概化為二維非均質(zhì)各向同性的非穩(wěn)定流地下水系統(tǒng)，建立如下的數(shù)學(xué)模型[1]來描述:

式中:H為含水層水位(m);H0(x，y)為初始水位(m);h(x，y，t)為第一類邊界 Γ1上的已知函數(shù)(m);q(x，y，t)為第二類邊界Γ2上的單寬流量(m2/d);K為淺層含水層(沿厚度)平均滲透系數(shù)(m/d);Hd為含水層底板標(biāo)高(m);w為單位時(shí)間在垂向上從單位水平面積含水層中流入或流出的水量(m/d);μ為給水度;Ω為模擬計(jì)算區(qū)域;Γ1為計(jì)算區(qū)域的第一類邊界(給定水位邊界);Γ2為計(jì)算區(qū)域的第二類邊界(給定流量邊界);n為第二類邊界Γ2的外法線方向。

1.3 數(shù)值模型

本文采用GMS軟件對(duì)上述數(shù)學(xué)模型進(jìn)行求解。首先對(duì)研究區(qū)進(jìn)行網(wǎng)格剖分，網(wǎng)格大小為237m×307m，共有有效單元格4 198個(gè)。由于2010年為平水年，本次模擬選擇2010年1月－2010年12月作為模擬期，以一個(gè)月為一個(gè)時(shí)間段，共有12個(gè)應(yīng)力期。含水層的頂板和底板高程分別根據(jù)地表高程和淺層含水層的底板標(biāo)高等值線確定的。本次模擬的初始流場(chǎng)為2010年1月的流場(chǎng)，是由實(shí)測(cè)水位資料通過插值方法得到的。

然后要對(duì)模型進(jìn)行識(shí)別和驗(yàn)證。通過調(diào)整水文地質(zhì)參數(shù)如滲透系數(shù)、給水度以及降水入滲系數(shù)等，使得觀測(cè)孔的計(jì)算水位和實(shí)測(cè)水位擬合程度不斷提高，從而使建立的模型能夠更為準(zhǔn)確地反映研究區(qū)的水文地質(zhì)條件[2]。

圖2為模型識(shí)別后典型觀測(cè)孔水位擬合曲線圖，從圖中可以看出觀測(cè)孔的計(jì)算水位和實(shí)測(cè)水位的變化趨勢(shì)相近，擬合程度較高，達(dá)到了模型調(diào)參識(shí)別的目的[3]。通過模型識(shí)別所得到的各分區(qū)參數(shù)值見表1。為進(jìn)一步驗(yàn)證所建立的地下水模型的可靠性，利用其它時(shí)期的資料對(duì)模型進(jìn)行驗(yàn)證。驗(yàn)證期為2013年1月－2013年12月，以一個(gè)月為一個(gè)時(shí)間段，共有12個(gè)應(yīng)力期。圖3為2013年12月計(jì)算流場(chǎng)和實(shí)際流場(chǎng)擬合圖，可以看出計(jì)算流場(chǎng)和實(shí)際流場(chǎng)擬合程度較高。由此可見，所建立的模型能夠較為真實(shí)地反映研究區(qū)的水文地質(zhì)條件，對(duì)于含水層結(jié)構(gòu)、邊界條件的概化以及水文地質(zhì)參數(shù)的選擇是合理的，可以用于地下水位預(yù)測(cè)[4]。

表1 水文地質(zhì)參數(shù)初始值和識(shí)別值一覽表

圖2 觀測(cè)孔水位擬合曲線圖

圖3 2013年12月計(jì)算流場(chǎng)和實(shí)際流場(chǎng)擬合圖

2 模型預(yù)測(cè)

2.1 預(yù)測(cè)方案

在降水保證率95%(特枯年)條件下，利用建立的地下水模型對(duì)不同的設(shè)計(jì)開采方案進(jìn)行模擬預(yù)測(cè)。本次模擬設(shè)計(jì)三種開采方案，分別為5萬 m3/d(方案一)、10萬 m3/d(方案二)和20萬 m3/d(方案三)。預(yù)測(cè)模型的邊界條件和水文地質(zhì)參數(shù)保持不變，補(bǔ)給項(xiàng)和排泄項(xiàng)按95%保證率重新計(jì)算，開采量根據(jù)不同方案的設(shè)計(jì)開采量確定。預(yù)測(cè)期為2013～2032年，以年為時(shí)間段，共20個(gè)應(yīng)力期，每個(gè)時(shí)間段內(nèi)包括不同的時(shí)間步長(zhǎng)，前10年的時(shí)間步長(zhǎng)為一天，后10年的時(shí)間步長(zhǎng)為一個(gè)月。初始流場(chǎng)選擇2012年12月末的流場(chǎng)。

2.2 預(yù)測(cè)結(jié)果

按方案一開采20年后的流場(chǎng)變化圖見圖4。從圖中可以看出，按方案一開采1年后，地下水位整體下降，開采3年、5年后，水位持續(xù)下降，直到10年后，地下水流場(chǎng)達(dá)到穩(wěn)定狀態(tài)，水位基本保持不變，與20年后相比降幅不超過0.01 m。

按方案二開采20年后的流場(chǎng)變化圖見圖5。從圖中可以看出，按方案二開采1年后，水源地位置水位降幅較大，開采3年、5年后，水位持續(xù)下降，直到10年后，地下水流場(chǎng)達(dá)到穩(wěn)定狀態(tài)，水位基本保持不變，與20年后相比降幅不超過0.02 m。

按方案三開采20年后的流場(chǎng)變化圖見圖6。從圖中可以看出，按方案三開采1年后，水源地中心地下水位大幅下降形成降落漏斗，開采3年、5年后，地下水位整體持續(xù)下降，漏斗面積不斷增大，直至10年后，地下水流場(chǎng)達(dá)到穩(wěn)定狀態(tài)，水位基本保持不變，與20年后相比降幅不超過0.02 m，漏斗面積不再增加。

圖4 方案一開采流場(chǎng)變化圖

圖7 為方案一條件下水源地及周邊地下水水位變化曲線，其中曾各莊、東蔡各莊和大周莊分別代表水源地上游、中心及下游位置。從圖中可以看出，水源地上游及中心的地下水位一開始都明顯下降，但中心年平均下降速率為1.2 m/a，上游下降速率為0.5 m/a。10年后水位基本達(dá)到穩(wěn)定狀態(tài)，20年后中心水位降深3.7 m。水源地下游的水位變化不明顯。

圖8為方案二條件下水位變化曲線。從圖中可以看出，水源地上游及中心的地下水位一開始都明顯下降，但中心年平均下降速率為 1.6 m/a，上游下降速率為 0.7 m/a。10年后水位基本達(dá)到穩(wěn)定狀態(tài)，20年后中心水位降深5.1 m。水源地下游的水位變化不明顯。

圖5 方案二開采流場(chǎng)變化圖

圖6 方案三開采流場(chǎng)變化圖

圖7 方案一開采水位變化曲線

圖8 方案二開采水位變化曲線

通過以上對(duì)預(yù)測(cè)結(jié)果的分析得到，水源地在三種方案持續(xù)開采20年后的中心水位降深分別為3.7 m、5.1 m及8.3 m，而水源地位置的含水層厚度為52 m，所以水位降深值遠(yuǎn)小于含水層厚度，不會(huì)造成含水層疏干，但方案三在水源地位置出現(xiàn)了面積較大的地下水降落漏斗，而方案一和方案二的地下水流場(chǎng)趨勢(shì)基本不變，只是水源地位置的水位降深與周邊相比較大。雖然三種方案均能使地下水系統(tǒng)達(dá)到穩(wěn)定狀態(tài)，補(bǔ)給量有保證，但是方案三的漏斗面積較大，可能會(huì)影響周邊的農(nóng)作物種植，與農(nóng)業(yè)生產(chǎn)爭(zhēng)水，對(duì)周圍的生態(tài)環(huán)境影響較大。所以為了盡可能的滿足昌黎縣的生產(chǎn)生活用水需求，最大程度地利用地下水資源，選擇方案二為合理的開采方案。

圖9 方案三開采水位變化曲線

3 結(jié)語

(1)應(yīng)用GMS軟件建立靖安水源地的地下水流數(shù)值模型，通過對(duì)模型進(jìn)行識(shí)別和驗(yàn)證，可以看出該模型能夠真實(shí)準(zhǔn)確地反映研究區(qū)的水文地質(zhì)概況，對(duì)邊界條件的分析以及參數(shù)的選擇是合理的，可以用于地下水位預(yù)測(cè)。

(2)利用建立的模型對(duì)水源地不同開采方案進(jìn)行模擬預(yù)測(cè)，預(yù)測(cè)結(jié)果表明，三種方案均能使地下水系統(tǒng)達(dá)到穩(wěn)定狀態(tài)，最大水位降深值遠(yuǎn)小于含水層厚度，不會(huì)造成含水層疏干，但方案三出現(xiàn)了面積較大的地下水降落漏斗，對(duì)周圍的生態(tài)環(huán)境影響較明顯。所以確定該水源地的允許開采量為10萬 m3/d。

(3)將GMS用于靖安水源地的地下水流模擬與預(yù)測(cè)，確定該水源地的允許開采量，可以為水源地的規(guī)劃建立提供科學(xué)依據(jù)，也為其它傍河型水源地的地下水資源評(píng)價(jià)提供了經(jīng)驗(yàn)和方法[5]。

[1]丁元芳，遲寶明，易樹平，等.Visual MODFLOW在李官堡水源地水流模擬中的應(yīng)用[J].水土保持研究.2006，13(5):99－105.

[2]陳冬琴.GMS軟件在杭嘉湖地下水資源評(píng)價(jià)中的應(yīng)用[J].軟件導(dǎo)刊.2007，(5):49－51.

[3]易立新，徐鶴.地下水?dāng)?shù)值模擬:GMS應(yīng)用基礎(chǔ)與實(shí)例[M].北京:化學(xué)工業(yè)出版社.2009.

[4]秦建甫，張海濱，萬偉鋒，等.哈頭才當(dāng)水源地地下水位模擬與預(yù)測(cè)[J].人民黃河.2012，34(5):52－54.

[5]劉記成，王現(xiàn)國(guó)，葛雁，等.Visual Modflow在鄭州沿黃水源地地下水資源評(píng)價(jià)中的應(yīng)用[J].地下水.2007，29(4):91－92.