譚廷靜，曹振東,3

(1.貴州省地質(zhì)礦產(chǎn)勘查開(kāi)發(fā)局111地質(zhì)大隊(duì)，貴州 貴陽(yáng) 550008；2.貴州地質(zhì)工程勘察設(shè)計(jì)研究院，貴州 貴陽(yáng) 550008；3.吉林大學(xué) 環(huán)境與資源學(xué)院，吉林 長(zhǎng)春 130021)

?

基于數(shù)值模擬技術(shù)的地下水流場(chǎng)變化特征分析

譚廷靜1,2，曹振東1,2,3

(1.貴州省地質(zhì)礦產(chǎn)勘查開(kāi)發(fā)局111地質(zhì)大隊(duì)，貴州 貴陽(yáng) 550008；2.貴州地質(zhì)工程勘察設(shè)計(jì)研究院，貴州 貴陽(yáng) 550008；3.吉林大學(xué) 環(huán)境與資源學(xué)院，吉林 長(zhǎng)春 130021)

礦山的開(kāi)采大多需要進(jìn)行疏干排水，而在疏干排水時(shí)持續(xù)的地下抽水勢(shì)必導(dǎo)致區(qū)內(nèi)含水層中地下水流場(chǎng)的變化，對(duì)周邊地下水的開(kāi)發(fā)利用造成影響。為了分析礦山開(kāi)采中地下水流場(chǎng)的變化特征，本文通過(guò)對(duì)興隆縣某金、鐵礦水文地質(zhì)條件分析建立研究區(qū)地下水系統(tǒng)的數(shù)值模擬模型，并通過(guò)對(duì)已知地下水動(dòng)態(tài)水位的擬合與檢驗(yàn)，確定模型的可靠性，然后利用識(shí)別后的模型預(yù)測(cè)興隆縣某金、鐵礦運(yùn)行投產(chǎn)含水層水位的變化特征，探討礦山開(kāi)采對(duì)地下水系統(tǒng)的影響。

數(shù)值模擬技術(shù)；礦山開(kāi)采；地下水流場(chǎng)

地球上蘊(yùn)藏著的種類(lèi)繁多的礦產(chǎn)資源，而礦產(chǎn)資源對(duì)人類(lèi)社會(huì)的發(fā)展與進(jìn)步有無(wú)可替代的作用，不斷促進(jìn)社會(huì)經(jīng)濟(jì)的發(fā)展，但礦產(chǎn)資源的開(kāi)發(fā)利用對(duì)地下水環(huán)境有著較大的影響作用。為了預(yù)防在礦山開(kāi)采中的突水事故，大多數(shù)礦山開(kāi)采時(shí)需要進(jìn)行疏干排水，疏干排水可能導(dǎo)致區(qū)域生態(tài)環(huán)境破壞、地下水位持續(xù)下降和水資源枯竭等[1]。隨著礦產(chǎn)資源開(kāi)采深度、開(kāi)采強(qiáng)度、開(kāi)采速度、開(kāi)采規(guī)模的增大，礦井突水事故頻繁發(fā)生，給人民生命財(cái)產(chǎn)造成重大損失[2]。礦山排水、礦井突水和水資源保護(hù)之間的矛盾突出，因而基于礦山排水下的地下水流場(chǎng)模擬預(yù)測(cè)較為重要。

國(guó)內(nèi)外很多研究者對(duì)地下水?dāng)?shù)值模擬方法進(jìn)行了大量探索，薛禹群等人建立露天礦排水疏干的二維水流模型，對(duì)疏干引起的地下水位變化情況進(jìn)行了分析[3]；張鳳娥等人利用地下水?dāng)?shù)值模擬技術(shù)對(duì)礦山開(kāi)采過(guò)程中地下水相關(guān)的問(wèn)題進(jìn)行研究[4]；周念清等人建立了地下水流三維數(shù)學(xué)模型，對(duì)地下水流場(chǎng)分布和礦坑涌水量進(jìn)行了預(yù)測(cè)[5]；楊彪等人在Visual Modflow中建立了地下水?dāng)?shù)值模擬模型，對(duì)金嶺南凡口鉛鋅礦鉛鋅礦地下水位動(dòng)態(tài)特征進(jìn)行了研究[6]；陳社明等人建立了地下水?dāng)?shù)值模擬模型，分析了河北某擴(kuò)建煤礦地下水流場(chǎng)變化特征[7]。日趨成熟的地下水模型軟件，如Visual Modflow、GMS等，也為這些模擬計(jì)算成果的準(zhǔn)確性提供了強(qiáng)有力的保證，使得地下水?dāng)?shù)值模擬已經(jīng)成為當(dāng)前地下水研究中最常用的手段之一。

本文采用數(shù)值模擬技術(shù)，在分析研究區(qū)水文地質(zhì)條件的基礎(chǔ)上，利用Visual Modflow軟件建立興隆縣某金、鐵礦的地下水?dāng)?shù)值模擬模型，并基于實(shí)測(cè)數(shù)據(jù)對(duì)模型進(jìn)行識(shí)別、驗(yàn)證，然后利用識(shí)別后的模型預(yù)報(bào)未來(lái)礦山開(kāi)采過(guò)程中含水層水位的變化特征，探討礦山開(kāi)采對(duì)地下水系統(tǒng)的影響。

1　研究區(qū)概況

興隆縣茅山鎮(zhèn)某金、鐵礦項(xiàng)目位于興隆縣青松嶺鎮(zhèn)螞蟻溝村及掛蘭峪鎮(zhèn)四拔子村一帶，礦區(qū)總面積為5.966 5 km2，中心坐標(biāo)為北緯40°15′2.25″，東經(jīng)117°32′51.79″。擬建工程礦區(qū)邊界西北緊鄰螞蟻溝村、距跑馬場(chǎng)村1 080 m，西側(cè)緊鄰道槽溝、距鐵爐子溝80 m，南側(cè)緊鄰?fù)跗綔衔?、壩尺溝，東側(cè)緊鄰四拔子村、大陰司溝、小陰司溝，距五拔子790 m，北側(cè)緊鄰黃砬棚、黃松峪。該礦山屬于“探轉(zhuǎn)采”新建礦山。

興隆縣茅山鎮(zhèn)螞蟻溝礦區(qū)位于六里坪山(1 442.8 m)南麓泃河與淋河之間的分水嶺上，處于北部花市花崗巖巖體(主峰1 007.2 m)與南部分布的長(zhǎng)城系分布區(qū)(高峰985.3 m)中間的片麻巖系中。礦區(qū)為中低山區(qū)，沿次級(jí)羽毛斷裂形成的次級(jí)溝谷很發(fā)育。礦區(qū)出露巖層有上太古界跑馬場(chǎng)組變質(zhì)巖系和第四系全新統(tǒng)松散堆積層。第四系坡洪積孔隙水巖組(Q4dl+pl)：主要分布在螞蟻溝內(nèi)的谷地及東部淋河的支溝中,其巖性為主溝中心為塊碎石夾砂，呈現(xiàn)棱角狀或半棱角狀，無(wú)分選性。第四系沖洪積孔隙含水層(Q4al+pl)：主要分布于泃河和淋河的主河道中,其巖性為礫卵石、漂石夾中粗砂,因?qū)偕絽^(qū)河谷,分選性不好,也無(wú)分帶性,厚一般7～10 m,局部可達(dá)13 m。透水不含水層：在山麓斜坡的緩坡處、山洼、山脊鞍部、溝谷上源都有坡殘積形成的亞粘土夾碎石，厚度1～3 m，最厚達(dá)5.0 m,其僅透水而不含水。

圖1　邊界條件概化

2　地下水?dāng)?shù)值模擬模型的建立

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

含水層結(jié)構(gòu)：根據(jù)水文地質(zhì)資料，本次模擬主要考慮淺層的裂隙水和孔隙水。對(duì)于裂隙水采用等效介質(zhì)原理，按孔隙水進(jìn)行模擬。模型在垂向上按一層處理，含水層類(lèi)型按照巖性進(jìn)行參數(shù)分區(qū)。

邊界條件：潛水含水層自由水面為系統(tǒng)上邊界，通過(guò)該邊界，潛水與系統(tǒng)外發(fā)生垂向水量交換(入滲補(bǔ)給、蒸發(fā)等)。模擬區(qū)的下部為基巖區(qū)，概化為隔水邊界。側(cè)向邊界的概化：根據(jù)模擬區(qū)地形條件，模擬區(qū)的西南側(cè)和東北側(cè)均為分水嶺，因此概化為隔水邊界；西北和東南側(cè)分布有泃河與淋河的河床，是區(qū)域內(nèi)地下水的流出通道，因此概化為流量邊界。

水力特征概化：評(píng)價(jià)區(qū)地下水系統(tǒng)符合質(zhì)量守恒定律；含水層分布廣，在常溫常壓下地下水運(yùn)動(dòng)符合達(dá)西定律。參數(shù)隨空間變化，體現(xiàn)了系統(tǒng)的非均質(zhì)性，但沒(méi)有明顯的方向性，所以參數(shù)概化成各向同性。由于地下水系統(tǒng)滲流運(yùn)動(dòng)要素隨時(shí)空變化，故地下水含水系統(tǒng)概化為非穩(wěn)定流。

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

對(duì)于非均質(zhì)、各向同性、空間三維結(jié)構(gòu)、非穩(wěn)定地下水流系統(tǒng)，可用如下微分方程的定解問(wèn)題來(lái)描述[8]：

式中：h(x，y，z，t)為含水層水頭(m)；h0(x，y，z)為含水層初始水頭(m)；K為含水層滲透系數(shù)(m/d)；μ為潛水含水層在潛水面上的重力給水度；S為自由面以下含水層的貯水率；W為源匯項(xiàng)(m/d)；P為潛水面上的降水入滲和蒸發(fā)等(m/d)；f1(x，y，z，t)為第二類(lèi)邊界上的水分通量(m/d)；Ω為計(jì)算區(qū)范圍；Γ1為流量邊界；Kn為邊界法線方向的滲透系數(shù)；n為邊界上的外法線。

圖2　網(wǎng)格平面和垂向剖分圖

2.3模型的離散及參數(shù)分區(qū)

2.3.1網(wǎng)格剖分

模擬面積為79.26 km2。根據(jù)Visual-MODFLOW的要求，在水平方向上對(duì)潛水含水層用相互垂直的平行線進(jìn)行網(wǎng)格剖分。將模擬區(qū)剖分成約500×500 m的單元格。模擬區(qū)單層總單元格為16 800個(gè)，其中有效單元格8 077個(gè)，剖分結(jié)果見(jiàn)圖2。剖分結(jié)果見(jiàn)圖2。

滲透系數(shù)及給水度通過(guò)收集以往的資料和本次水文地質(zhì)試驗(yàn)的成果獲取。考慮模擬區(qū)內(nèi)不同含水介質(zhì)類(lèi)型的，將模擬區(qū)淺層含水層根據(jù)滲透性分為4個(gè)區(qū)。其中：1和3分區(qū)為松散巖類(lèi)孔隙水；2分區(qū)為構(gòu)造和風(fēng)化裂隙水；4分區(qū)為透水而不含水區(qū)域。滲透系數(shù)及給水度分區(qū)及編號(hào)見(jiàn)圖3。

2.4模型的識(shí)別和驗(yàn)證

根據(jù)模擬區(qū)地下水統(tǒng)測(cè)資料，以2013年1月1日到2013年6月30日作為模型識(shí)別期，共計(jì)181 d。根據(jù)前期調(diào)查資料，識(shí)別期的補(bǔ)給量和排泄量根據(jù)實(shí)測(cè)值計(jì)算，以2013年1月1日水位為模型識(shí)別過(guò)程的初始水位。本次共選取5口觀測(cè)井作為識(shí)別和檢驗(yàn)的校準(zhǔn)數(shù)據(jù)點(diǎn)，觀察孔主要位于礦區(qū)附近，5口井的分布位置如圖4所示。

圖3　模擬區(qū)滲透系數(shù)及給水度分區(qū)圖

圖4　觀察點(diǎn)分布位置

圖5　識(shí)別期含水層地下水位擬合曲線

以觀測(cè)孔NO.1、NO.5的水位擬合情況為例(圖5、圖6)。觀測(cè)點(diǎn)在各個(gè)時(shí)段的模擬水位與觀測(cè)水位擬合程度較好，水位擬合差小于2 m。由于模擬區(qū)位于山區(qū)，地下水隨地形坡度變化幅度較大，因此該模擬精度能夠滿(mǎn)足生產(chǎn)需求，模型反應(yīng)了實(shí)際情況下水位隨時(shí)間的變化特征。

表1　含水層滲透系數(shù)識(shí)別結(jié)果

以2012年7月1日地下水位作為模型驗(yàn)證期初始水位，對(duì)模型進(jìn)行驗(yàn)證。圖5-7反映了觀測(cè)孔NO.1和NO.5淺層和深層含水層的水位擬合情況。通過(guò)模擬水位與觀察水位的對(duì)比看出，模擬效果也比較好，觀測(cè)井中最大誤差在允許范圍內(nèi)。擬合結(jié)果較好，能夠滿(mǎn)足預(yù)測(cè)要求。

圖6　驗(yàn)證期含水層地下水位擬合曲線

圖7　礦區(qū)開(kāi)采2 a、5 a、10 a后含水層地下水水位等值線

3　礦區(qū)開(kāi)采對(duì)地下水流場(chǎng)的影響分析

以識(shí)別和校正后的模型為基礎(chǔ)，源匯項(xiàng)的處理方法同識(shí)別期一致。其中降水量以多年平均降水量進(jìn)行計(jì)算。預(yù)測(cè)模型的初始水位取2014年1月1日統(tǒng)測(cè)水位。預(yù)測(cè)時(shí)段根據(jù)礦區(qū)服務(wù)年限進(jìn)行預(yù)報(bào)。依據(jù)某金、鐵礦項(xiàng)目設(shè)計(jì)，鐵礦的服務(wù)年限為8.39 a，金礦的服務(wù)年限為2.11年。本次模擬時(shí)段以鐵礦服務(wù)年限為準(zhǔn)，并考慮建設(shè)期等因素適當(dāng)延長(zhǎng)至10 a，即2014年到2024年。由于礦區(qū)13條鐵礦體零散賦存于堅(jiān)硬的片麻巖中，與其同步生成，礦體頂?shù)装逋暾?、穩(wěn)固，頂?shù)装辶严恫话l(fā)育，因此開(kāi)采過(guò)程中與淺層地下水受其影響有限。因此，在預(yù)報(bào)期僅對(duì)部分礦體的開(kāi)采增加新的地下水開(kāi)采。

通過(guò)2 a(2016年)、5 a(2019年)和10 a(2024年)后含水層地下水流場(chǎng)圖可以看出，2年 a(2016年)，按礦山建成后的開(kāi)采條件下，研究區(qū)淺層地下水流場(chǎng)與初始流場(chǎng)相比，在高海拔山區(qū)地下水位有所下降，在溝谷地區(qū)水位變化幅度不大，整個(gè)流場(chǎng)相比于初始流場(chǎng)更趨于平滑。在礦區(qū)所在位置，由于處于中間地下水分水嶺位置，一部分區(qū)域透水但不含水，其他區(qū)域則受礦井開(kāi)采影響地下水位出現(xiàn)一定幅度的下降，最大水位下降幅度在6～10 m左右。通過(guò)5年后(2019年)的流場(chǎng)圖可以看出，模擬區(qū)地下水流場(chǎng)與初始流場(chǎng)相比仍然未發(fā)生明顯變化。在礦區(qū)所在位置，最大水位下降幅度在13～18 m左右。通過(guò)10 a后(2024年)的流場(chǎng)圖可以看出，地下水等值線雖然與初始流場(chǎng)分布一致，但山區(qū)水位持續(xù)下降，水利坡度減小。但在溝谷地區(qū)，地下水位相比于初始流場(chǎng)仍未見(jiàn)明顯變化。在礦區(qū)所在位置，最大水位下降幅度在25～32 m左右。

4　結(jié)語(yǔ)

(1) 在充分了解模擬區(qū)水文地質(zhì)條件的基礎(chǔ)上，利用數(shù)值模擬技術(shù)分析礦山開(kāi)采過(guò)程中地下水流場(chǎng)的變化特征，具有精度高，效果好，簡(jiǎn)單可靠的特點(diǎn)。識(shí)別驗(yàn)證后的地下水?dāng)?shù)值模型能再現(xiàn)地下水的水位動(dòng)態(tài)、運(yùn)移規(guī)律，較好地反映了開(kāi)采條件下含水層地下水流場(chǎng)的狀況。

(2) 模擬區(qū)內(nèi)金礦、鐵礦開(kāi)采后，由于礦山排水，礦山附近地下水位出現(xiàn)下降，但地下水等值線雖然與初始流場(chǎng)分布一致。其中，山區(qū)水位持續(xù)下降，水利坡度減小；在溝谷地區(qū)，地下水位相比于初始流場(chǎng)仍未見(jiàn)明顯變化；在礦區(qū)所在位置，最大水位下降幅度在25～32 m左右。

(3) 建議在礦山期間中，注重礦山排水的再利用。由于裂隙含水層內(nèi)水質(zhì)較好，考慮將礦山排水達(dá)標(biāo)處理后，作為周邊工業(yè)和農(nóng)業(yè)的供水水源，從而減少地下水的開(kāi)采。

[1]王天平,解建倉(cāng),張建龍,等.基于突變理論的西王寨礦區(qū)礦井疏干水開(kāi)發(fā)利用風(fēng)險(xiǎn)評(píng)價(jià)[J].西安理工大學(xué)學(xué)報(bào).2010(4):417-423.

[2]張立新,李長(zhǎng)洪,趙宇.礦井突水預(yù)測(cè)研究現(xiàn)狀及發(fā)展趨勢(shì)[J].中國(guó)礦業(yè).2009(1):88-90.

[3]薛禹群,吳吉春,謝春紅,等.元寶山露天煤礦地下水疏干數(shù)值模擬[J].煤炭學(xué)報(bào).1996,21(3):255-260.

[4]張鳳娥,劉文生.煤礦開(kāi)采對(duì)地下水流場(chǎng)影響的數(shù)值模擬—以神府礦區(qū)大柳塔井田為例[J].安全與環(huán)境學(xué)報(bào).2002,2(4):30-33.

[5]周念清,徐小明,江思珉,等.武山銅礦地下水?dāng)?shù)值模擬與礦坑涌水量預(yù)測(cè)[J].勘察科學(xué)技術(shù).2012(2):34-38.

[6]楊彪,羅周全,王益?zhèn)?等.基于Visual Modflow的礦山地下水流場(chǎng)分析及預(yù)測(cè)[J].礦冶工程.2013,33(4):11-21.

[7]陳社明,王威,劉宏偉,等.河北省某礦建煤礦地下水流場(chǎng)變化數(shù)值模擬[J].人民黃河.2015(10):69-73.

[8]楊青春,盧文喜,馬洪云.Visual Modflow 在吉林省西部地下水?dāng)?shù)值模擬中的應(yīng)用[J].水文地質(zhì)工程地質(zhì).2005(3):67-69.

2016-08-14

譚廷靜(1987-)，女，貴州開(kāi)陽(yáng)人，助理工程師，主要從事水文地質(zhì)、環(huán)境地質(zhì)工作。

曹振東(1984-)，男，甘肅會(huì)寧人，工程師，主要從事水文地質(zhì)、地?zé)岬刭|(zhì)和淺層地溫能研究工作。

P641.2

B

1004-1184(2016)04-0214-03