內(nèi)陸河流域二維地下水運動模擬模型研究

內(nèi)陸河流域二維地下水運動模擬模型研究

在干旱和半干旱地區(qū)，水資源是最寶貴的自然資源。這不僅是沙漠綠洲形成、發(fā)展和穩(wěn)定的基礎(chǔ)，也是環(huán)境的必要和有機(jī)組成部分，也是生態(tài)環(huán)境的制約因素。而地下水作為水資源的重要組成部分和天然植被賴以生存的重要水源,有時甚至是唯一水源,在干旱區(qū)水資源開發(fā)利用中占據(jù)重要的地位。隨著干旱內(nèi)陸河流域水土資源矛盾的日益凸顯,再加之地表水資源在時間和空間上分布極不均勻,地下水資源已經(jīng)成了干旱內(nèi)陸河流域社會經(jīng)濟(jì)發(fā)展和生態(tài)環(huán)境建設(shè)必不可少的部分,是干旱區(qū)綠洲得以維持和存在的重要依靠,也是研究內(nèi)陸河流域水資源統(tǒng)一調(diào)度管理和優(yōu)化配置的關(guān)鍵要素。由此可見,地下水資源的可持續(xù)利用,緊密關(guān)系到國民經(jīng)濟(jì)的可持續(xù)發(fā)展和干旱區(qū)生態(tài)環(huán)境的重建。然而如何合理利用好地下水資源,是一個很值得深入探討的問題。由于內(nèi)陸河流域從河流源頭到下游地區(qū)可劃分為高山冰雪凍土帶、山區(qū)植被帶、平原綠洲帶和下游荒漠帶,而降水主要集中在高山冰雪凍土帶和山區(qū)植被帶,山前綠洲帶和下游荒漠帶降水稀少,基本不產(chǎn)流,因此,出山口徑流基本代表了流域的水資源量,且形成與出山口以上的高山冰雪帶和山區(qū)植被帶。山區(qū)徑流出山口后,在天然情況下部分洪水通過漫溢滲漏補給地下水,其余的通過河道滲漏補給河道兩側(cè)的地下水,最終消失于荒漠或尾閭湖中。在水資源大規(guī)模開發(fā)利用的情況下,人類活動增加了中下游的河道兩側(cè)徑流的側(cè)枝循環(huán),通過地表徑流灌溉和地下水開采,多次垂向交換,最終河水趨于消失。然而隨著計算機(jī)技術(shù)的高速發(fā)展及其在各個學(xué)科領(lǐng)域的廣泛應(yīng)用,地下水?dāng)?shù)值模擬技術(shù)也得到了長足進(jìn)展,國內(nèi)外已經(jīng)建立了許多城市和地區(qū)的地下水?dāng)?shù)值模擬模型和管理模型,例如MODFLOW、MT3D99、VizualMODFLOWF、VizualGroundwater和FEFLOW等,然而對干旱區(qū)內(nèi)陸河流域地下水系統(tǒng)進(jìn)行數(shù)值模擬相對較少。利用數(shù)值模型對地下水流和溶質(zhì)運移問題進(jìn)行模擬的方法以其有效性、靈活性和相對廉價性逐漸成為干旱區(qū)地下水研究領(lǐng)域的一種不可或缺的重要方法,在水資源評價問題、地下水污染問題、地下水管理與合理開發(fā)以及地下水-地表水聯(lián)合評價等方面受到越來越大的重視和廣泛的應(yīng)用。因此,本文將FEFLOW模型引用到干旱內(nèi)陸河流域,選取新疆臺蘭河流域作為研究區(qū),開展干旱內(nèi)陸流域地下水系統(tǒng)數(shù)值模擬研究。從而為內(nèi)陸河流域地下水的模擬提供一種可借鑒的方法。為了更準(zhǔn)確地評價研究區(qū)地下水資源,預(yù)測地下水系統(tǒng)狀態(tài)及其變化趨勢,本文在詳細(xì)分析研究區(qū)地質(zhì)和水文地質(zhì)條件的基礎(chǔ)上,建立研究區(qū)水文地質(zhì)概念模型和地下水流模擬模型,采用基于有限元方法的FEFLOW(Finiteelementsubsurfaceflowsystem)對所建立的二維地下水流數(shù)值方程進(jìn)行求解,完成模型的識別和驗證,從而有助于認(rèn)識臺蘭河流域在人類活動影響下地下水資源的發(fā)展趨勢,并預(yù)測給定條件下臺蘭河流域地下水流場變化情況,為干旱區(qū)地下水資源的合理開發(fā)利用提供科學(xué)依據(jù),旨在促進(jìn)干旱內(nèi)陸河流域地下水系統(tǒng)的研究。1地形地貌和地貌臺蘭河流域位于新疆阿克蘇地區(qū)溫宿縣境內(nèi),地理坐標(biāo)為80°21′44″～81°10′14″E,40°41′41″～42°15′13″N,地處北溫帶干旱區(qū),氣候干燥,日照充足,多風(fēng)沙,降水稀少,蒸發(fā)較大,晝夜溫差大,多年平均降水量為177mm,多年平均蒸發(fā)量為1800mm,屬典型干旱大陸性氣候。從地形地貌上看,流域分為丘陵區(qū)和沖洪積平原區(qū),總體地勢北高南低,可明顯分為3個帶,即上游沖洪積傾斜平原帶、中游沖洪積扇細(xì)土平原帶和下游沖洪積扇細(xì)土平原荒漠地帶。根據(jù)巖性結(jié)構(gòu)特征及相對關(guān)系和所處的地貌部位劃分,主要為上更新統(tǒng)沖洪積砂礫石層全新統(tǒng)沖洪積物。2水文干預(yù)模型的建立2.1地下水含水層分布根據(jù)臺蘭河流域水文地質(zhì)資料和不同時期地下水等水位線圖,確定以北部的沖洪積傾斜平原帶到南部的沖洪積扇的細(xì)粒土平原荒漠地帶為研究區(qū)。分析現(xiàn)有的探孔和生產(chǎn)井資料,臺蘭河流域?qū)贊撍?承壓水系統(tǒng)。自前山帶至溢出帶為潛水含水層。扇后緣靠山前帶為粗顆粒的砂卵石,至扇緣前為粗中顆粒相砂礫石,向平原過渡為細(xì)顆粒相沉積物。臺蘭河沖洪積扇第四系含水層均大于200m。但由于其上游補給源不同,各沖洪積扇的潛水含水層的富水程度差異也很大。潛水的埋藏深度,在前山帶大于100m,向扇緣方向逐漸變?yōu)?0～30m,30～0m。自扇緣帶以南到南部沙漠地帶為潛水-承壓含水層,在200m深度內(nèi),可見2～4個含水層。含水層巖性由砂礫石逐漸相變?yōu)楹[粗砂、粉細(xì)砂、亞砂、亞粘及粘土互層。南部沙漠邊緣一帶,含水層巖性為粉細(xì)砂,沖洪積成因的含水層厚度變薄。研究區(qū)內(nèi)廣泛分布第四系松散層孔隙潛水,含水層為砂碎石及中細(xì)砂層。由于受地質(zhì)構(gòu)造、地形地貌和沉積條件的制約,自北而南沉積物顆粒漸細(xì),地下水位埋深逐漸變淺,富水性漸弱,含水層由北部單一的潛水含水層,向南逐漸過渡為具有統(tǒng)一水力聯(lián)系的雙層或多層結(jié)構(gòu)的潛水-承壓含水層。由于沉積環(huán)境的變化交替,粘性土的分布在水平方向上都有較大的變化,隔水層呈斷續(xù)狀或透鏡體分布,很難在區(qū)域內(nèi)找到一層分布比較穩(wěn)定的隔水層,故承壓水和潛水在區(qū)域上有不可分割的水力聯(lián)系,構(gòu)成了統(tǒng)一的含水層系統(tǒng),可視為統(tǒng)一的多層介質(zhì)的潛水含水系統(tǒng)。本次模擬將研究區(qū)概化成非均質(zhì)各向同性二維非穩(wěn)定流潛水含水層。2.2流量邊界及邊界根據(jù)所確定的研究范圍及研究區(qū)地下水體系特征,模型的邊界條件(圖1)概化為:AB、MN、NA為第二類邊界(補給邊界),通量根據(jù)水文資料確定;HI為第二類邊界(排泄邊界),由于邊界處為平原荒漠地帶,可以認(rèn)為不受抽水井的影響,但由于缺少觀測井長期觀測資料,因此,很難確定荒漠中地下水位的變化,所以把該邊界定為定流量邊界(定流量值為0.0038m/d)。其余假設(shè)為第二類零流量邊界,概化后的計算面積為38.71×108m2。其北邊界以地下水埋深50m為界,南至荒漠邊緣為界,東西分別以鄰縣界為準(zhǔn)。研究區(qū)共自動剖分三角形單元1675個,結(jié)點897個。2.3基于fewelle-qued的調(diào)整根據(jù)收集地礦、新疆生產(chǎn)建設(shè)兵團(tuán)和水利部門鉆孔抽水試驗資料,確定研究區(qū)滲透系數(shù)K、給水度μ的初始值,結(jié)合研究區(qū)水文地質(zhì)圖和含水層剖面圖,直接由FEFLOW調(diào)用調(diào)試,調(diào)試時K和μ值的改變?nèi)匀徊捎肁rcView方便改動,并作出K和μ的分區(qū)圖(圖2)。由于研究資料缺少,無法準(zhǔn)確確定隔水層底板高程,但考慮到研究區(qū)域大,水平落差小(1/1000～1/1500),進(jìn)行概化認(rèn)定潛水層是100m的含水層,所以把模擬區(qū)的潛水含水層隔水底板高程全部給定為一個固定值,即區(qū)域地面最低處高程減去100m的值:110m。3數(shù)學(xué)模型的建立3.1h3水層隔水板內(nèi)h??x(Κh?Η?x)+??y(Κh?Η?y)+W(x,y,t)=μ?Η?tΗ(x,y,t)|t=0=Η0(x,y)(x,y)∈DΗ(x,y,t)|Γ1=Φ1(x,y,t)(x,y)∈Γ1?Η?n|Γ2=q2(x,y,t)(x,y)∈Γ2(1)式中,H為水位(或稱水頭)[L];H0:水頭初始條件[L];h:潛水含水層隔水底板至潛水自由面的距離[L];K:滲透系數(shù)[L/T];D:研究區(qū)域;μ:給水度;Г1:第一類邊界;Г2:第二類邊界;Φ1:第一類邊界條件,為給定邊界的水頭[L];q2:第二類邊界條件,為給定邊界上的單寬流量[L/T];n:為第二類邊界Г2的外法線方向;W(x,y,t):地下水源匯項,即垂直方向上地下水的抽出和入滲補給[L/T]。3.21補充排水項目研究區(qū)補給項主要包括大氣降水凝結(jié)水入滲、渠系水入滲、田間灌溉水入滲及地下水側(cè)向徑流;排泄項主要包括地下水開采、側(cè)向流出、蒸發(fā)蒸騰。3.2.1凝結(jié)水入滲量臺蘭河流域多年平均降水量177mm,并且降水次數(shù)少、有效的降水量低。因此,降水入滲對地下水的補給必須考慮地下水埋深因素,并采取分配在某些時段的方法。根據(jù)本區(qū)已有的研究,降水、凝結(jié)水入滲補給地下水限于地下水埋深小于5.0m的地區(qū)。降水入滲以全年每次降水大于10mm的降水量之和為基數(shù)(有效降水量),臺蘭河流域平原區(qū)多年平均為62mm,降水入滲系數(shù)取0.03～0.20(根據(jù)地下水埋深不同有所變化),分配在降水集中的7至9三個月。凝結(jié)水在地下水位埋深1～3m區(qū)全年為27mm,埋深3～5m區(qū)為35mm。凝結(jié)水主要發(fā)生在秋季溫差大、濕度大的降溫過程中,因此凝結(jié)水入滲的分配是:地下水埋深1～3m區(qū)從10-05至11-05的31d中每天入滲1mm,埋深3～5m區(qū)從10-05至11-10的36d中每天入滲1mm。3.2.2灌溉水入滲補給從理論上講,田間灌水入滲量的大小與灌溉定額、灌水量和灌溉地段的包氣帶巖性、厚度、濕度及氣象等多種因素有關(guān)。本文用入滲系數(shù)法計算入滲量。由于灌溉水入滲補給主要發(fā)生在耕地上,將2005年臺蘭河流域土地利用格局圖輸入ArcViewGIS中,將耕地圖和鄉(xiāng)級行政圖迭加,按照種植面積來確定灌溉制度,依各作物面積和灌溉制度算出每次灌水時流域某種作物的灌井水量和水庫水量,然后根據(jù)渠系利用系數(shù)算出應(yīng)抽井水量。同時根據(jù)不同地下水埋深耕地的面積,對迭加的圖進(jìn)行識別、剖分,并結(jié)合不同地下水埋深的入滲系數(shù),從而確定接受補給耕地單元面積,由此計算田間灌溉水的入滲量(表1)。3.2.3鄉(xiāng)單元的抽水量渠系滲漏補給量指渠系滲漏補給地下水的水量。將流域內(nèi)渠系分布圖輸入到ArcViewGIS中,劃分網(wǎng)格,確定渠系所在單元格,以備FEFLOW調(diào)用。在流域內(nèi)按每眼機(jī)井抽水量200m3/h,每天抽12h,由分布在作物面積上的機(jī)井?dāng)?shù)量算出鄉(xiāng)單元機(jī)井的日抽水強度,從而得到在鄉(xiāng)單元上每天的抽水量,進(jìn)一步根據(jù)式(2)計算渠系滲漏水入滲補給量:Q=qw×(1-N′)×0.9/SU(2)式中?Q為渠系滲漏補給量?qw為鄉(xiāng)單元某種作物灌水時機(jī)井每天抽水量(m3/d),Ν′為渠系有效利用系數(shù),0.9是渠系滲漏水的90%入滲補給地下水?SU為鄉(xiāng)單元耕地面積(m2)。3.2.4地下水側(cè)向徑流補給和排泄量計算根據(jù)文獻(xiàn)和《阿克蘇地區(qū)臺蘭河灌區(qū)規(guī)劃報告》中計算地下水側(cè)向流入(出)量?計算公式為∶Q=K×B×H×I×cosα(3)式中,K為地下水滲透系數(shù)(m/d),B為過水?dāng)嗝鎸挾?m),H為地下水含水層厚度(m),I為地下水水力坡度(‰),α為地下水流向與斷面線垂線夾角。3.2.5機(jī)井抽水量的測定臺蘭河流域機(jī)井?dāng)?shù)目多,是模擬的難點所在。按每眼機(jī)井抽水量200m3/h,每天抽12h,由分布在作物面積上的機(jī)井?dāng)?shù)量算出鄉(xiāng)單元機(jī)井的日抽水強度,從而得到在鄉(xiāng)單元上每天的抽水量。根據(jù)臺蘭河流域地下水開采強度的分區(qū),同時將流域地下水開采強度分區(qū)圖輸入到ArcViewGIS中,根據(jù)強度分區(qū)將地下水的開采量分布到相應(yīng)的網(wǎng)格單元上,以備FEFLOW調(diào)用,再結(jié)合鄉(xiāng)單元上每天的地下水抽取量,從而計算研究區(qū)的實際地下水開采量。3.2.6蒸發(fā)量的季節(jié)變化土面蒸發(fā)指的是沒有植被覆蓋或植被蓋度很小地區(qū)的蒸發(fā)。在穩(wěn)定蒸發(fā)條件下即為潛水蒸發(fā),可以按照不同地下水埋深計算潛水蒸發(fā)量。臺蘭河流域多年平均蒸發(fā)量最大在7月,其次在6、5月,然后在8、4月。據(jù)此把地下水埋深在1～3m區(qū)的陸面蒸發(fā)量按月分配,地下水埋深3～5m區(qū)的陸面蒸發(fā)量分配在7月。根據(jù)臺蘭河流域的氣候,平原區(qū)多年平均降雨量為62mm,植物生長旺盛季節(jié)在3—10月的8個月中,因此,此間的蒸散也最強烈,但不同季節(jié)的變化很明顯。因此林木草地蒸發(fā)引起的損失,在各時間段按比例分配到植物生長的幾個月中(表2)。4實測結(jié)果及分析以臺蘭河2005年度的土地利用類型圖及2004—2007年地下水位的觀測數(shù)據(jù)調(diào)試模型參數(shù)并校正模型,通過多次運行模型、調(diào)參、優(yōu)選,使得模型的計算結(jié)果與實測地下水?dāng)?shù)據(jù)變化趨勢盡可能一致,最終確定模型選用的參數(shù),利用10個長期觀測井的地下水位實測值與模擬值進(jìn)行對比。用2008年度地下水位的觀測數(shù)據(jù)來驗證模型。K和μ的調(diào)試結(jié)果見圖2。圖3是5個分布在流域上、中、下游的長期觀測井(井號1、2、5、8、10)2008年度的地下水位實測值與模擬值的對比?？梢钥闯?觀測井地下水位的模擬值與實測值的變化趨勢基本吻合,說明模型擬合的地下水位變化趨勢與研究區(qū)的水文地質(zhì)情況一致,模型的效果是較好的,可以用來模擬區(qū)域地下水的運動狀態(tài)。從模擬結(jié)果來看,臺蘭河流域?qū)崪y地下水位與模擬地下水位擬合過程比較理想,平均絕對誤差0.05m,最大誤差±0.5m。其它5個井的擬合情況見表3。根據(jù)各源匯項的計算結(jié)果,上述對地下水運動的數(shù)值模擬和模型校驗是建立在水均衡原理上的。對整個流域來說,流域進(jìn)水-流域消耗=流域排水。因此,流域消耗項即垂向源匯項的處理是地下水運動的模型校驗和數(shù)值模擬的關(guān)鍵。其中流域消耗項主要包括渠系入滲、田間入滲、降雨入滲、潛水蒸發(fā)、植物蒸騰及地下水開采。利用2004—2007年的資料計算的平均值與利用2008年的地下水位觀測資料驗證模型時的結(jié)果進(jìn)行比較,兩者的均衡量比較接近,均衡差率為-7.45%,絕對值小于10%,表明其結(jié)果符合規(guī)定的差率范圍(表4)。5建立水文模型利用ArcViewGIS與FEFLOW軟件,通過對內(nèi)陸河流域的潛水-微承壓水系統(tǒng)選用二維地下潛水流模型,基于水量均衡來處理地下水的源匯項,模擬了內(nèi)陸河流域地下水的演化特征。從地下水流數(shù)值模擬結(jié)果可以看出,所建立的水文地質(zhì)概念模型和數(shù)學(xué)模型是正確的,選取的水文地質(zhì)參數(shù)和計算的源匯項基本合理,符合臺蘭河流域地下水的實際情況,可用于該流域流場研究和水源地的開采規(guī)劃,在此基礎(chǔ)上可以為流域地