MODFLOW在滲渠取水工程中的應(yīng)用研究

張智印，朱 亮，羅奇斌

(1.西北大學(xué) 大陸動力學(xué)國家重點實驗室，陜西 西安710069;2.中國地質(zhì)科學(xué)院水文地質(zhì)環(huán)境地質(zhì)研究所，河北石家莊050061)

MODFLOW在滲渠取水工程中的應(yīng)用研究

張智印1，朱 亮2，羅奇斌1

(1.西北大學(xué) 大陸動力學(xué)國家重點實驗室，陜西 西安710069;2.中國地質(zhì)科學(xué)院水文地質(zhì)環(huán)境地質(zhì)研究所，河北石家莊050061)

安塞油田近幾年產(chǎn)量不斷攀升，給干旱少雨、地下水資源十分貧乏的安塞地區(qū)供水帶來了很大困難。為滿足采油廠的需水要求，需要對譚家營水源地的取水能力進(jìn)行評估。水源地的取水工程擬采用滲渠取水方式。根據(jù)延河河谷區(qū)概況建立延河水資源數(shù)學(xué)模擬模型，使用MODFLOW計算模擬滲渠在不同設(shè)計方案下的取水量。通過對滲渠取水地區(qū)地下水補給量、河流徑流量的分析，并和滲渠開采量進(jìn)行比較，確定滲渠取水量1 774 m3/d是有保證的。

MODFLOW;滲渠;數(shù)值模擬;取水量

隨著計算機技術(shù)和數(shù)值方法的發(fā)展，數(shù)值模擬已漸漸取代傳統(tǒng)的模擬技術(shù)成為研究地下水運動規(guī)律和定量評價地下水資源的主要手段而且其發(fā)展趨勢已遠(yuǎn)遠(yuǎn)超出作為一種計算手段的原有范疇，成為模擬一些水文地質(zhì)過程發(fā)生、發(fā)展的主要手段［1］。無論從理論上還是從實踐上，對于多孔介質(zhì)飽和帶中地下水運動機理的認(rèn)識是明確的，模擬飽和帶地下水流的模型是成熟的，應(yīng)用也是最廣泛的［2］。在眾多計算機軟件中，MODFLOW已成為應(yīng)用最為普及的地下水運動數(shù)值模擬計算程序［3］。MODFLOW是由美國地質(zhì)調(diào)查局于80年代開發(fā)的一套專門用于孔隙介質(zhì)中地下水三維有限差分法數(shù)值模擬軟件，具有程序結(jié)構(gòu)合理，離散方法簡單化、求解方法多樣化、易于理解和便于操作等特點［4］。對譚家營水源地的滲渠取水量，采用MODFLOW進(jìn)行評價，最終通過比較確定結(jié)果滿足要求。

1 研究區(qū)基本地質(zhì)概況

研究區(qū)地處西北內(nèi)陸，屬溫帶大陸性半干旱季風(fēng)氣候，春季干旱多風(fēng)，氣候多變;夏季炎熱維持時間短，干旱與雨澇交錯并存;秋季涼爽多雨;冬季嚴(yán)寒，少雨雪。多年平均降水量為505.3 mm，年最大降水量為793.9 mm，年最小降水量為289.6 mm，雨季集中在6—9月，占年總降水量的74.3%;多年平均蒸發(fā)量為1 054 mm，干旱指數(shù)為2.08;年平均氣溫為8.8℃，年極端最高氣溫為36.8℃，極端最低氣溫為零下23.6℃。

研究區(qū)位于鄂爾多斯盆地陜北斜坡中段東部，面積約2.24 km2。區(qū)域構(gòu)造輪廓為平緩西傾單斜。區(qū)域傾角一般小于半度，平均坡度8～10 m/km。局部構(gòu)造皆屬微弱起伏的鼻狀隆起，未見明顯的斷裂及背斜、向斜構(gòu)造。

區(qū)內(nèi)主要河流為延河，為黃河的一級支流，發(fā)源于靖邊縣天賜灣鄉(xiāng)的周山，由西北流向東南，流經(jīng)安塞縣、延安市、延長縣匯入黃河，在安塞縣境內(nèi)長約90 km流域面積2 649 km2，河道常流量為0.5～1.5 m3/s，多年平均徑流量約為1×108m3。依據(jù)安塞站1981－2008年28年的水文、降水資料，可以看出延河年徑流量多年呈豐、枯交替變化，并受降水量控制明顯，年徑流量與年降水量歷時均呈總體減少趨勢。

延河河谷區(qū)主要發(fā)育河漫灘與一級階地、二級階地、三級階地，階地分布較為連續(xù)，三級階地的部分地段被沖溝切割。河漫灘:連續(xù)分布于延河現(xiàn)代河床兩側(cè)，寬約10～100 m，高出河水位約1～2 m，地表出露多以粘性土或礫卵石為主;一級階地:呈帶狀連續(xù)分布于延河兩側(cè)，在河曲部位較寬，最寬處可達(dá)600 m左右。階面平坦，微向河谷傾斜，前緣以陡坎與漫灘相接，高出河床約5～6 m，出露地層具二元結(jié)構(gòu)，上部為粉質(zhì)粘土，下部為礫卵石層，總厚6～9 m;二級階地:分布于延河局部地段，多呈帶狀，寬約20～50 m，長約300～400 m，階面較平坦，向河谷微傾。地層具二元結(jié)構(gòu)，上部為粉質(zhì)粘土，下為礫卵石層，總厚約10～20 m;三級階地:位于河谷邊部，為堆積階地，階面以斜坡或臺階形式向河流方向傾斜。

區(qū)內(nèi)地層分為兩層:下部為侏羅系中統(tǒng)安定組(J2a)巖層，主要為砂泥巖，與下伏直羅組為連續(xù)沉積，厚度54～65 m;上部為第四系梁峁區(qū)黃土和河谷區(qū)的砂卵石層。

根據(jù)含水介質(zhì)的不同、水力性質(zhì)與埋藏條件的差異，區(qū)內(nèi)地下水分為第四系黃土裂隙孔洞潛水、第四系松散巖類孔隙—侏羅系碎屑巖類裂隙潛水:

(1)第四系黃土裂隙孔洞潛水:分布于黃土梁峁區(qū)。區(qū)內(nèi)黃土梁峁區(qū)地形破碎，沖溝發(fā)育，且多數(shù)沖溝已切入基巖，不利于地下水的儲存與富集。黃土潛水含水層分布不連續(xù)，含水性弱，透水性差，降水多以地表徑流形式排泄，補給作用微弱，地下水位埋藏深，泉水流量均小于9.5 m3/d，無實際集中供水意義。

(2)第四系松散巖類孔隙—侏羅系碎屑巖類裂隙潛水:分布于河谷區(qū)的漫灘與階地區(qū)。

第四系松散巖類孔隙潛水，含水層巖性為粉土、砂、砂礫卵石層，厚度變化較大，由河漫灘至階地后緣逐漸變厚，由幾米至十幾米不等。地下水位埋深一般0.4～5 m?；鶐r裂隙潛水，含水層主要為基巖的風(fēng)化裂隙帶和構(gòu)造裂隙帶，厚度一般僅幾米，地下水位埋深0.5～6 m。

第四系沖洪積砂、砂礫卵石層與侏羅系砂巖強風(fēng)化帶構(gòu)成了具有雙層介質(zhì)結(jié)構(gòu)的統(tǒng)一的潛水含水層(簡稱孔隙－裂隙潛水含水層)，據(jù)鉆孔抽水資料，水位埋0.45～1.87 m，含水層厚度 10.23～12.54 m，含水層滲透系數(shù) 0.05 ～8.03 m/d，單井涌水量3.37～62.04 m3/d，屬弱富水，為區(qū)內(nèi)主要供水目的層。

2 地下水系統(tǒng)模型

2.1 滲渠方案設(shè)計

譚家營水源地滲渠工程在延河現(xiàn)代河床主流線之下布設(shè)。滲渠結(jié)構(gòu):第四系松散層開挖深度(平均厚度)為2.5 m，開挖寬度為15 m;侏羅系砂泥巖平均開挖深度2.5 m，開挖寬度2.6～2.8 m。滲渠頂部之上為人工濾床，滲渠長度L為800 m，滲渠人工濾床寬度(開挖寬度)B為15 m，則滲渠的人工濾床水平面積F為1.2×104m2。人工濾床可以大量激發(fā)河水的入滲補給，同時可發(fā)揮濾床的過濾功能而濾清河水，達(dá)到既取水又凈沙的綜合目的?？紤]到河水泥沙含量大，人工濾床的滲透能力最終將取決于濾床的淤積與堵塞狀況;從供水安全考慮，人工濾床的淤塞系數(shù)β取濁水河流的經(jīng)驗值為0.3。確定滲渠人工濾床的單位面積垂向滲透水量Vz為0.22 m3/d·m2。人工濾床的河水入滲強度與河水入滲量估算值見表1。

表1 滲渠人工濾床的河水入滲強度與河水入滲量估算表

2.2 模型計算域及其邊界性質(zhì)

綜合考慮水源地水文、水文地質(zhì)特征、滲渠工程特點以及取水的影響范圍，確定模型計算域及其邊界性質(zhì)。

模型計算域:水平方向上，沿河谷的東西長度取1.5 km，垂直河谷的南北寬度取1.5 km(有效計算域為整個河谷區(qū))。垂直方向上，模型計算頂界面為地表面與河水面;水源地的含水層結(jié)構(gòu)為第四系松散層與侏羅系裂隙帶上下疊置的雙層結(jié)構(gòu)，根據(jù)區(qū)內(nèi)基巖風(fēng)化帶發(fā)育厚度及滲渠底設(shè)置深度等，模型計算深度取在河床下5 m處。

滲渠的上邊界:河床底部將會在滲渠運行期間有淤泥沉積，河床底部淤積層(泥膜)的滲透性能采用淤積層的垂向滲透系數(shù)與厚度的比值來表示。區(qū)內(nèi)河床淤積層(泥膜)厚1.5 cm，粉土夾泥膜的垂向滲透系數(shù)Kz取0.01 m/d(滲透試驗的最小值)，河水水深為30 cm，則滲渠頂部的河水入滲強度為 0.21 m3/d·m2。

模型區(qū)的水平周界與底界:模型計算域的東西邊界(河谷上下游邊界)距離滲渠取水工程位置較遠(yuǎn)，設(shè)定為一類定水頭邊界;模型有效計算區(qū)的南北兩側(cè)邊界設(shè)定在河谷邊界線，視為隔水邊界;模型底界(河床下5 m深度)為極弱風(fēng)化的砂泥巖，視為隔水邊界;模型頂部邊界為潛水面，且不考慮除滲渠抽水量外的其它垂向交換水量;滲渠頂部人工濾層的河水淹沒區(qū)，設(shè)置為河流補給邊界。模型計算區(qū)底板與潛水面均為順河谷方向由西向東傾斜，底板與水力坡度按地形坡度3.5‰設(shè)置。

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

延河常年流水不斷，河水與地下水之間具有密切的水力聯(lián)系。延河兩側(cè)的含水層是一個整體的滲流系統(tǒng)，設(shè)計在河床主流線之下的滲渠長期取水后，勢必會影響到滲渠兩側(cè)含水層中地下水的運動狀況，滲渠附近滲流的鉛直分速度不可忽視，故滲渠取水條件下的地下水運動采用三維流研究。建立滲流數(shù)學(xué)模型如下:

式中:H為地下水位標(biāo)高(m);K為滲透系數(shù)(m/d)，Kh水平滲透系數(shù)，Kv垂向滲透系數(shù)，Kr河床淤積層垂向滲透系數(shù);Ss為彈性釋水率(1/m);x，y，z，t為坐標(biāo)變量(m)，時間變量(d);f0為初始時刻的域內(nèi)地下水位標(biāo)高(m);f1為第一類(定水位)邊界水位標(biāo)高(m);Qr為滲渠開采量(m3/d);Hr為河水位標(biāo)高(m);A為河流淹沒面積(m2);Mr為河床淤積層厚度(m);Q為滲渠開采量(m3/d);F為滲渠進(jìn)水側(cè)面面積(m2);n為二類邊界或滲渠水平周界外法線方向;W為潛水面垂向交換水量(m3/d·m2);Γ1、Γ2、l0為一類邊界、二類邊界、滲渠水平周界;D為計算區(qū)范圍。

2.4 模型構(gòu)建與網(wǎng)格剖分

為準(zhǔn)確刻畫滲渠的復(fù)雜結(jié)構(gòu)，首先采用長方體單元對計算域進(jìn)行精細(xì)的剖分，模型區(qū)的剖分網(wǎng)格:模型區(qū)剖分為南北方向200行、東西方向200列，水平網(wǎng)格間距為7.5 m;模型區(qū)垂向上剖分為2層，第1層為第四系松散層與侏羅系基巖強風(fēng)化裂隙帶組成的雙層結(jié)構(gòu)含水層(層厚3 m)，第2層為侏羅系基巖裂隙含水層(層厚2 m)，滲渠部位的第2層為空心長方體，是匯集含水層地下水的蓄水與通水空間。

2.5 模型識別與驗證

模型的識別主要考慮含水層滲透系數(shù)和給水度。

模型調(diào)參的初始值根據(jù)鉆孔抽水試驗資料計算的結(jié)果，確定的模型區(qū)含水層水文地質(zhì)參數(shù):滲透系數(shù)K為0.07～0.59 m/d，給水度 μ 為 0.047～0.119。根據(jù)區(qū)域經(jīng)驗與試算結(jié)果，含水層彈性釋水率取1×10－5m－1。滲渠的濾層參數(shù)為:滲透系數(shù)0.68 m/d、給水度0.072。通過運行模型，反復(fù)調(diào)整有關(guān)參數(shù)，直至模型計算水位與抽水試驗的實際水位相擬合。最終確定水源地模型計算區(qū)的水文地質(zhì)參數(shù)(詳見表2)。

表2 水源地滲渠取水三維流數(shù)值模型參數(shù)表

3 模型預(yù)測

3.1 預(yù)測方案設(shè)計

譚家營水源地的滲渠取水?dāng)?shù)值模擬，考慮人工濾床為置換沙回填、原河沙回填兩種情況，以及平水期、偏枯水期、枯水期三種河流水文狀態(tài)(表1)，加上河流斷流的極端狀態(tài)，組合的每個水源地數(shù)值模型五種預(yù)測方案詳見表3。

第一預(yù)測方案:滲渠濾層頂部為置換沙回填，濾層的入滲能力為0.22 m3/d·m2;河流水文處于平水期狀態(tài)，即滲渠濾層全部被河水淹沒。

第二預(yù)測方案:滲渠濾層頂部為置換沙回填，濾層的入滲能力為0.22 m3/d·m2;河流水文處于偏枯水期狀態(tài)，即滲渠濾層被河水淹沒75%。

第三預(yù)測方案:滲渠濾層頂部為置換沙回填，濾層的入滲能力為0.22 m3/d·m2;河流水文處于枯水期狀態(tài)，即滲渠濾層被河水淹沒50%。

第四預(yù)測方案:滲渠濾層頂部為原河沙回填，濾層的入滲能力為0.06 m3/d·m2;河流水文處于平水期狀態(tài)，即滲渠濾層全部被河水淹沒。

第五預(yù)測方案:河流水文處于斷流狀態(tài)，即滲渠濾層無河水淹沒。

表3 水源地滲渠取水?dāng)?shù)值模型預(yù)測方案設(shè)置表

3.2 預(yù)測結(jié)果

水源地滲渠取水采用地下水流三維有限差分法求解，將各含水層有關(guān)參數(shù)和邊界條件等代入模型后，即可進(jìn)行計算。以2012年5月的水位為初始水位(詳見圖3)，依據(jù)所建立的數(shù)值模型，對各設(shè)計方案進(jìn)行滲渠工程取水預(yù)測，第一、二、三、四預(yù)測方案的預(yù)測時長為20 a，第五預(yù)測方案的預(yù)測時長為90 d。預(yù)測取水20 a后水位流場見圖4。預(yù)測取水20 a末典型網(wǎng)格單元的地下水位(頭)降深與降速詳見圖5。

對預(yù)測結(jié)果比較分析，方案二是比較合理的:

1)地下水位(頭)時空變化趨勢:地下水位(頭)下降趨勢及下降速度:滲渠內(nèi)的地下水位，在開采初期急速下降，并很快趨于穩(wěn)定。區(qū)域地下水位，在開采1年后降速開始變緩，開采10 a后變的更加平緩，河谷上下游地下水位已趨于穩(wěn)定，但河谷兩側(cè)在開采20 a未也未達(dá)到穩(wěn)定狀態(tài)。

2)地下水流場與降深場:水源地開采后并未影響地下水徑流的總方向，地下水仍然是由河谷上游向下游方向徑流，僅在滲渠周圍很小的范圍內(nèi)形成地下水局部封閉降落漏斗，封閉降落漏斗沿滲渠軸向延伸，垂直滲渠方向的封閉降落漏斗寬度在100～150 m之間，最大封閉降落漏斗面積約0.10 km2。水源地開采后的水位(頭)降落已擴展到河谷兩側(cè)邊界，但并未影響到模型計算區(qū)的河谷上下游邊界。

3.3 滲渠取水量評價

3.3.1 地下水補給量分析

區(qū)內(nèi)地下水主要補給源為大氣降水入滲補給、農(nóng)灌水入滲補給，以及滲渠開采激發(fā)河水的入滲補給;河谷上游斷面地下水徑流補給量小，忽略不計。

1)大氣降水入滲補給量:根據(jù)安塞縣氣象站1956—2007年間51 a的降水量資料，多年平均降水量487.1 mm;水源地河谷區(qū)面積631 843 m2，降水入滲系數(shù)取0.25，計算的多年平均降水入滲補給量為211 m3/d。

2)農(nóng)灌回歸補給量:水源地現(xiàn)有水澆地面積約1 200畝，農(nóng)灌定額為300 m3/a·畝，農(nóng)灌回歸補給系數(shù)為0.25，則計算區(qū)農(nóng)灌回歸補給量為9.0×104m3/a，平均到每天為247 m3/d。

3)滲渠開采激發(fā)河水的入滲補給:數(shù)值模型計算的第二開采方案河水入滲補給量1 776 m3/d。水源地第2開采方案的地下水總補給量為2 234 m3/d，滲渠取水量分別為1 774 m3/d，取水量分別占補給量的79.4%。水源地地下水補給量見表4。

由表知在斷流情況下，滲渠取水量為173 m3/d(斷流第10天的取水量)，而入滲補給量為458 m3/d，水源地的滲渠取水量小于地下水補給量。所以，只要地下水補給量有保障，則滲渠取水量就是有補給保證的。

3.3.2 河水徑流量分析

地下水補給量主要來源于滲渠取水過程中的激發(fā)河水入滲補給量，而河水徑流量是影響河水入滲量的重要因素之一。故，需要對河水徑流量進(jìn)行分析。

按照控制流域面積比法由延河安塞站推算水源地地區(qū)的延河徑流量，安塞站控制流域面積1 334 km2，譚家營水源地區(qū)控制流域面積808 km2。水源地第二開采方案的河水入滲量占的延河枯流量0.03 m3/s(即2 592 m3/d)的68.5%。顯然，滲渠取水量確定為第二開采方案的1 774 m3/d是有河水入滲補給保障的。

綜上，水源地的滲渠取水量確定為第二開采方案(置換沙回填，偏枯水年):滲渠取水量結(jié)論為1 774 m3/d。

表4 水源地地下水補給量與滲渠取水量計算成果表

4 結(jié)語

1)延河雖然徑流量較大，滿足需水要求，但含泥沙量大，故不宜直接取水利用。因此，滲渠取水成為水源地取水的有效方案。

2)考慮水源地水文地質(zhì)、滲渠工程特點，確定模型的邊界條件，建立滲渠取水三維流數(shù)值模型，預(yù)測滲渠取水量1 774 m3/d。分析可知滲渠取水量是有河水流量保證的，也有地下水補給保證的。

3)根據(jù)延河流量監(jiān)測資料，水源地河流無斷流發(fā)生，但連續(xù)枯流量時間可達(dá)5～10 d。建議為保證水源地運行安全，可以設(shè)定連續(xù)10 d的河流斷流時間，考慮河流斷流導(dǎo)致的滲渠取水量衰減而造成的供水量不足，水源地構(gòu)筑能夠調(diào)蓄10天供水的蓄水裝置(如調(diào)蓄水池等);為確保濾床有效滲水，每3～5 a清洗一次濾層，每次洪水過后都要清除濾床表層的泥皮。

［1］魏林宏，束金龍，郝振純.地下水流數(shù)值模擬的研究現(xiàn)狀和發(fā)展趨勢［J］.重慶大學(xué)學(xué)報(自然科學(xué)版).2000，(23):51.

［2］周仰效.國外地下水流及傳輸?shù)哪M［C］.現(xiàn)狀與模擬.1995，200－209.

［3］吳劍鋒，朱學(xué)愚.由MODFLOW淺談地下水流數(shù)值模擬軟件的發(fā)展趨勢［J］.工程勘察.2000，(2):12－15.

［4］李文躍，張博，洪梅等.Visual MODFLOW在大慶龍西地區(qū)地下水?dāng)?shù)值模擬中的應(yīng)用［J］.世界地質(zhì).2003，22(2):161 －165.

Application of MODFLOW in Channel Seepage Water Intake Works

ZHANG Zhi－ yin1，Zhu Liang2，Luo Qi－ bin1
(1.State Key Laboratory of Continental Dynamics，?Northwest University，Xi an，710069;2.The Institute of Hydrogeology and Environment Geology，Chinese Academy Geological Sciences，Shijiazhuang，Hebei，050061)

The rising production of Ansai oilfield has brought great difficulty for water supply in Ansai area，and it is short of rain and groundwater in the region.The paper makes an assessment on water intake capacity of Tan jiaying water source to evaluate the intake quantity which will meet requirement of production plant.In the evaluation test，water is gathered from seepage channels.Based on the mathematical model of Yanhe water resources，which is established with the consideration of general situation in Yanhe valley，infiltration gallery water withdrawal in different schemes are calculated by MODFLOW.The comparison between groundwater recharge and river runoff in tested area and seepage channel water withdrawal indicates that water supply is above 1 774 m3/d，which meets the requirement of production plant.

MODFLOW;seepage channel;numerical simulation and water withdrawal

TU991.3

A

1004－1184(2013)06－0026－04

2013－06－18

張智印(1986－)，男，山東菏澤人，在讀碩士研究生，研究方向為水文地質(zhì)。