嚴(yán)順君，俞學(xué)森，武蘭珍

(1.青海省水利水電勘測(cè)設(shè)計(jì)研究院，青海 西寧 810000；2.甘肅農(nóng)業(yè)大學(xué)水利水電工程學(xué)院，甘肅 蘭州 730000)

社會(huì)經(jīng)濟(jì)和城市化進(jìn)程的快速發(fā)展，加劇了水資源供需矛盾[1]，作為水資源之一的地下水，水質(zhì)較地表水不易受到污染，在一些地區(qū)被過(guò)度開采，使得資源環(huán)境問(wèn)題日益突出[2- 3]。地下水相關(guān)研究已日趨成熟，張宏仁[4- 7]對(duì)地下水滲流有限差分法差分格式及邊界條件進(jìn)行了研究，地下水資源承載力[2，8- 10]的概念也受到廣泛重視。且得益于計(jì)算機(jī)的快速發(fā)展，數(shù)值解法在地下水研究中得到廣泛應(yīng)用[11- 13]。格爾木河流域?qū)偾嗪Ｊ「駹柲臼?，該流域生態(tài)脆弱，近年來(lái)，維持該流域生態(tài)綠洲需水的泉水溢出量受到人工地下水開采影響，綠洲面積呈縮小趨勢(shì)，但相關(guān)研究不足?；贛odflow建立人工地下水開采對(duì)泉水溢出量影響模型，分析地下水開采量增加與泉水溢出量減少之間的相互關(guān)聯(lián)。

1 模型構(gòu)建

1.1 基本原理

當(dāng)不考慮地下水密度變化，假定滲透系數(shù)的主軸方向與坐標(biāo)軸的方向一致，則孔隙介質(zhì)中地下水三維滲流偏微分方程為：

(1)

式中，Kxx、Kyy、Kzz—滲透系數(shù)在x，y和z方向的分量；h—水頭；w—單位體積流量；Ss—孔隙介質(zhì)的貯水率；t—時(shí)間。

滲流方程的離散化多用有限元法和有限差分法，Modflow是基于后者的地下水?dāng)?shù)值模擬軟件，可用于孔隙介質(zhì)中三維地下水?dāng)?shù)值模擬，由美國(guó)地質(zhì)調(diào)查局開發(fā)。該軟件在地下水滲流分析中的準(zhǔn)確性已得到認(rèn)可，陳皓銳[14]應(yīng)用GMS-Modflow構(gòu)建了不同情景下地下水運(yùn)動(dòng)模型，借此分析研究區(qū)潛水位在氣候變化和人類活動(dòng)影響下的變化趨勢(shì)。雪力卡提·阿里木[15]利用Modflow軟件構(gòu)建了昌吉呼圖壁河流域地下水模型，模擬結(jié)果與實(shí)測(cè)值吻合較好。

1.2 模擬區(qū)概況

格爾木河發(fā)源于昆侖山北坡，最終注入達(dá)布遜湖等鹽湖。上游由雪水河、昆侖河兩大支流組成，匯流后由南向北至格爾木市區(qū)再次分成東西兩支。其地下徑流呈現(xiàn)山區(qū)地下水向河谷或山間河谷方向流，河谷和山間盆地地下水沿河流方向流動(dòng)的規(guī)律。流域地下水主要貯存在沖洪積扇，此處第四紀(jì)松散沉積物厚度較大，為地下水運(yùn)動(dòng)提供了空間?？紫稘撍话阌芍苓呄蚋叱梯^低的盆地中心滲流，其水平分帶具有狀分帶規(guī)律，與流域內(nèi)地貌、巖相規(guī)律一致。

格爾木河流域沖洪積平原可劃分為三部分：山前至沖洪積扇中部、沖洪積扇中部及下部和沖洪積扇前緣。其中，沖洪積扇前緣含水層巖性以砂質(zhì)為主，山前至沖洪積扇中部、沖洪積扇中部及下部主要為砂礫石、卵礫石含水層。地下水埋深從山前至沖洪積扇前緣依次減小，山前至中部地下水埋深不小于50m，中部至沖洪積扇前緣地下水埋深一般為15～35m，沖洪積扇前緣地下水埋深較淺，且由于砂質(zhì)含水層透水性較好，地下水開采井大多建于此部位。格爾木流域含水層厚度較大，一般大于200m。表1為含水層特性表。

表1 含水層巖性特性

1.3 模型建立

流域泉水補(bǔ)給綠洲帶植被生長(zhǎng)需水及尾閭湖泊部分需水，格爾木河出山后，在綠洲帶、鹽化草甸帶上邊界形成兩個(gè)明顯的泉水溢出帶，因此選取南起昆侖山前第四系基巖交界處，北至達(dá)布遜鹽湖區(qū)上邊界，東西兩側(cè)以格爾木河沖洪積扇邊緣為界區(qū)域作為模擬范圍，并根據(jù)特性將其劃分為戈壁帶、綠洲帶和鹽化草甸區(qū)三部分，整個(gè)研究區(qū)概化為一個(gè)均質(zhì)含水層中夾有一個(gè)薄隔水層。模擬區(qū)總面積5400km2，在計(jì)算核心區(qū)，網(wǎng)格剖分為800m×1000m，其它區(qū)域網(wǎng)格剖分為1600m×1000m。模擬范圍內(nèi)戈壁區(qū)、綠洲區(qū)、鹽化草甸區(qū)水文地質(zhì)參數(shù)取值見(jiàn)表2。

表2 含水層參數(shù)分區(qū)與取值

由位于戈壁帶和綠洲區(qū)上游1986年7個(gè)實(shí)測(cè)長(zhǎng)觀測(cè)孔資料，綠洲區(qū)下游以及整個(gè)鹽化草甸區(qū)20世紀(jì)80年代末完成的水文地質(zhì)圖31個(gè)地下水埋深點(diǎn)數(shù)值，建立模擬區(qū)初始水位，如圖1所示。綠洲區(qū)和鹽化草甸區(qū)非地帶性植被的蒸騰蒸發(fā)作為研究區(qū)主要的排泄項(xiàng)，其空間蒸發(fā)量分布如圖2所示。

圖1 模擬區(qū)初始水位

圖2 模擬區(qū)蒸發(fā)排泄分布位置圖

研究區(qū)多年平均降雨量約80～120mm，綠洲區(qū)降雨入滲補(bǔ)給量為0.0271億m3。河流入滲主要考慮垂直向，不同分區(qū)河流入滲補(bǔ)給，經(jīng)巖性、模型調(diào)試分析，荒漠戈壁區(qū)滲漏系數(shù)0.19，綠洲區(qū)0.15，鹽化草甸區(qū)0.1。由遙感影像并考慮其季節(jié)性變化，確定河流寬度荒漠區(qū)為125m，河流進(jìn)入綠洲區(qū)后，分東西兩支，寬度分別為80、90m，鹽化草甸區(qū)入流處為60m，出流處為30m。

1.4 地下水開采

研究區(qū)主要用水對(duì)象為格爾木市的工農(nóng)業(yè)發(fā)展與生活用水。地下水開采主要取自潛水層，現(xiàn)狀共133口開采井，設(shè)計(jì)最大供水能力28.41萬(wàn)m3/d。主要水源地有位于格爾木河西側(cè)的格爾木市自來(lái)水公司二期水源地和位于格爾木河?xùn)|側(cè)的鹽湖一期水源。概化現(xiàn)狀地下水取水井抽水量自西向東分別為14.04萬(wàn)、15.75萬(wàn)、14.04萬(wàn)m3/d，開采井概化位置分布圖如圖3所示。為進(jìn)一步研究地下水開采量與泉水溢出量關(guān)聯(lián)性，假設(shè)在現(xiàn)狀開采井井位不變的前提下，增大地下水開采量，分別取地下水開采量為1.6億、2.0億、2.5億m3模擬計(jì)算泉水溢出量。

圖3 模擬區(qū)人工開采井概化位置分布圖

2 結(jié)果分析

現(xiàn)狀地下水開采情況下，模型720d流場(chǎng)圖如圖4所示，從圖4中可以得出：

(1)由于鹽化草甸區(qū)蒸發(fā)強(qiáng)烈，而且植被分布呈現(xiàn)塊狀，地下水的流場(chǎng)并非似初始流場(chǎng)概化出完全平行的直線，局部地下水流場(chǎng)呈閉合狀。通過(guò)模擬計(jì)算，模擬區(qū)儲(chǔ)水量的變化量為-339萬(wàn)m3，結(jié)果為負(fù)均衡。

(2)山前戈壁區(qū)總補(bǔ)給量為3.32億m3，其中，河道和灌渠測(cè)滲補(bǔ)給量為2.32億m3；綠洲區(qū)的補(bǔ)給量2.04億m3，其中，河道滲漏補(bǔ)給1.63億m3，農(nóng)田灌溉補(bǔ)給量0.41億m3；鹽化草甸區(qū)河流滲漏補(bǔ)給量0.492億m3；綠洲區(qū)和鹽化草甸區(qū)植被蒸發(fā)量為3.69億m3；綠洲區(qū)和鹽化草甸區(qū)泉水溢出量為0.443億m3。

(3)地下水開采量為1.6億、2.0億、2.5億m3時(shí)，洪積扇前緣泉水溢出量分別為0.443億、0.224億、0.068億m3，地下水開采量與泉水溢出量呈現(xiàn)負(fù)相關(guān)，為滿足流域生態(tài)綠洲需水，減小地下水開采量應(yīng)是首要考慮因素。

圖4 模擬720d流場(chǎng)圖

3 結(jié)論

綠洲作為格爾木流域主要的植被生態(tài)系統(tǒng)，其生態(tài)脆弱，對(duì)水資源的依賴程度高。泉水溢出量減少，將直接威脅綠洲生態(tài)?；贛odflow地下水模型顯示，現(xiàn)狀地下水開采使得該區(qū)域儲(chǔ)水量呈負(fù)均衡，且地下水開采量的增加會(huì)進(jìn)一步減少研究區(qū)沖積扇前緣泉水溢出量，加劇綠洲退化趨勢(shì)，此后，應(yīng)采取挖掘節(jié)水潛力，嚴(yán)格地下水開采管理，利用其它水源替代地下水等措施，逐步減小地下水開采量。流域生態(tài)需水、耗水是一整體系統(tǒng)，地下水開采其響應(yīng)因素應(yīng)是泉水溢出量、徑流量、蒸發(fā)量等，不僅僅是引起泉水溢出量變化，這是本文的不足之處，尚需進(jìn)一步研究分析。