芻議地下水取水影響模擬研究

費(fèi)國松許廣東李 丹閆 浩（.江蘇省水文水資源勘測局常州分局 常州 3000 .江蘇省水文水資源勘測局宿遷分局宿遷 3800）

芻議地下水取水影響模擬研究

費(fèi)國松1許廣東2李 丹1閆 浩1
（1.江蘇省水文水資源勘測局常州分局 常州 213000 2.江蘇省水文水資源勘測局宿遷分局宿遷 223800）

本文以取用東鵬煤礦礦坑水為例，采用數(shù)值模型分析地下水取水影響，根據(jù)數(shù)值模型的建立、求解、識別、驗(yàn)證，全方位闡述采用模型分析地下水取水影響的可靠性，為地下水取水影響探索新的分析模式。

地下水 數(shù)學(xué)模型 邊界概化 識別驗(yàn)證

1 項(xiàng)目取水基本情況

東鵬煤礦位于淮河流域，開采深度-350～-800m，礦井面積16.1km2。研究區(qū)內(nèi)孔隙地下水補(bǔ)給資源量為3.1 萬m3/d，可開采量為1.65萬m3/d，現(xiàn)狀條件下已開采資源量為1.22萬m3/d，仍有0.43萬m3/d的剩余可開采量。新建項(xiàng)目設(shè)計(jì)日取水量為0.2萬m3/d。

根據(jù)收集到的淮河流域東鵬煤礦區(qū)已有的水文地質(zhì)資料以及野外勘探等成果，通過概化出研究區(qū)水文地質(zhì)概念模型，結(jié)合地下水流在介質(zhì)中運(yùn)移的數(shù)學(xué)描述建立研究區(qū)地下水?dāng)?shù)學(xué)模型。

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

2.1 研究區(qū)邊界的概化

長巨斷層：落差2500m，東升西降，使區(qū)內(nèi)奧灰與區(qū)外前震旦紀(jì)變質(zhì)巖接觸，形成阻水邊界。

馬翔斷層：落差400～200m，西升東降，區(qū)內(nèi)奧灰與區(qū)外寒武系和前震旦紀(jì)變質(zhì)巖在深部接觸，大部地段為阻水邊界，但在北部地段為透水邊界。

東溝斷層：落差500m，南升北降，區(qū)內(nèi)奧灰與區(qū)外二疊系接觸，為阻水邊界。

鵬山斷層：落差2000m，北升南降，在劉鳧斷層以東寒武系及震旦系地層出露，在劉鳧斷層以西，區(qū)內(nèi)奧灰與區(qū)外上二疊統(tǒng)或侏羅系接觸，形成阻水邊界。

東鵬煤礦區(qū)的垂向邊界可以概化為：煤礦區(qū)的上邊界為地表，煤層下部為奧陶系灰?guī)r，因此下邊界（巖溶底板）可以定義為隔水邊界。

2.2 地下水的補(bǔ)、徑、排條件

地下水主要以3種形式得到補(bǔ)給：一是降水補(bǔ)給，鵬山北緣的灰?guī)r裸露區(qū)直接接受大氣降水的補(bǔ)給；二是側(cè)向補(bǔ)給，該水文單元東南的鵬山，有約300km2的寒武系灰?guī)r露頭區(qū)，地表溶溝、溶槽及節(jié)理較發(fā)育，較易得到補(bǔ)給，其所含地下水向北徑流側(cè)滲補(bǔ)給奧灰；三是上覆第四系水的下滲。

地下水的徑流方向基本與地表水徑流方向一致，由東北向西南，南部鵬山的寒武系灰?guī)r、奧灰裸露區(qū)接受的補(bǔ)給水，順巖層的傾斜方向向北運(yùn)動，與向南運(yùn)動的水在露頭附近匯合，形成強(qiáng)徑流帶。

地下水的排泄主要以泉的形式泄出地表，但隨著奧灰水開采量及礦井排水量的增加，區(qū)域范圍內(nèi)原有的地下水流場已發(fā)生較大的改變，人工排泄成為地下水的主要方式。

綜上所述，可將東鵬煤礦區(qū)概化為一個(gè)統(tǒng)一的非均質(zhì)、各向同性、三維非穩(wěn)定地下水流系統(tǒng)。

3 地下水?dāng)?shù)值模型的建立與求解

3.1 地下水?dāng)?shù)學(xué)模型的建立

地下水?dāng)?shù)學(xué)模型是用來描述地下水水頭、水溫及水質(zhì)等現(xiàn)象及其變化過程的數(shù)學(xué)表達(dá)公式，它是用數(shù)學(xué)方法表述經(jīng)過簡化和概化的地下水系統(tǒng)。其數(shù)學(xué)公式及定解條件如下：

式中：D—滲流區(qū)域；H—地下水水頭（m）；Kx，Ky，Kz—x，y，z方向上的滲透系數(shù)（m/d）；S—自由水面之下的含水層單位儲水系數(shù)（1/m）；W—承壓含水層源匯項(xiàng)（1/d）；p—潛水面的蒸發(fā)量和降水的補(bǔ)給量等（1/d）；μd—潛水含水層重力給水度；H0（x，y，z）—初始流場水頭分布值（m）；H1（x，y，z，t）—第一類邊界水頭分布值（m）；q（x，y，z，t）—第二類邊界單位面積流量（m3/m2·d）；Kn—邊界面外線方向的滲透系數(shù)（m/d）；Г1、Г2—滲流區(qū)第一、二類邊界。

3.2 地下水?dāng)?shù)學(xué)模型的求解軟件

為了求解上述數(shù)學(xué)模型，需要把模型中的二階非線性偏微分方程轉(zhuǎn)化為線性方程之后再求解，目前人們常采用數(shù)值法來求解，本文即是運(yùn)用基于有限差分法的Visual MODFLOW4.2軟件來求解上述地下水?dāng)?shù)學(xué)模型。

3.3 網(wǎng)格剖分

在整個(gè)模型中初始網(wǎng)格將具有均勻的網(wǎng)格單元大小，大網(wǎng)格剖分為50行50列5層。

根據(jù)淮河流域東鵬煤礦內(nèi)地層結(jié)構(gòu)特征、水力特征、流場特征、模擬精度及邊界條件等要求，將整個(gè)東鵬煤礦區(qū)剖分為5050的單元格，垂向?yàn)?層，其中第二層與第四層為隔水層。

3.4 模擬時(shí)間與初始水頭

根據(jù)前期收集到的東鵬煤礦區(qū)地下水水位動態(tài)觀測資料，選取2011年1月1日的地下水位作為識別階段的地下水初始流場，潛水含水層初始水頭見圖1。

3.5 水文地質(zhì)參數(shù)確定

此次模型模擬所選用的參數(shù)有滲透系數(shù)、給水度、蒸發(fā)強(qiáng)度系數(shù)以及降雨入滲系數(shù)等。上述這些參數(shù)是根據(jù)研究區(qū)內(nèi)多次抽水試驗(yàn)、多次勘探及觀測資料所確定。

通過水文地質(zhì)物探及鉆探或抽水試驗(yàn)得到的參數(shù)數(shù)據(jù)進(jìn)行匯總，就可以得到該區(qū)域含水層的水文地質(zhì)參數(shù)。根據(jù)東鵬煤礦區(qū)各初始參數(shù)的分布規(guī)律及滲流場的特征，將該研究區(qū)劃分為7個(gè)參數(shù)區(qū)域。

3.6 源匯項(xiàng)確定

東鵬煤礦區(qū)地下水系統(tǒng)的源匯項(xiàng)主要包括補(bǔ)給項(xiàng)與排泄項(xiàng)。研究區(qū)的地下水補(bǔ)給項(xiàng)包括降雨入滲補(bǔ)給、河流滲漏補(bǔ)給等；排泄項(xiàng)包括潛水蒸發(fā)排泄、人工開采排泄等。

3.6.1 補(bǔ)給項(xiàng)

a.降水入滲補(bǔ)給

降水入滲補(bǔ)給是東鵬煤礦區(qū)地下水的主要補(bǔ)給源。研究區(qū)降水入滲補(bǔ)給量的計(jì)算公式如下：

Q=α·F·P

式中：Q—大氣降水入滲量（m3/a）；α—降水入滲系數(shù)；F—降水入滲面積（m2）；P—年降水量（mm/a）。

東鵬煤礦區(qū)年內(nèi)降雨差異也較大，其中汛期6～9月份降雨量占全年的70%～76%，形成春秋兩旱夾一澇的自然規(guī)律。礦區(qū)多年平均降雨711.9mm，降水較為豐富但年際變化大且時(shí)空分布不均。

b.河流滲漏補(bǔ)給

地表河流的滲漏是研究區(qū)內(nèi)潛水的重要補(bǔ)給源。由于河床較高，地表水不能排入，河流與含水層的交換量可以用以下公式計(jì)算求得：

式中：QR—河流滲漏量（m3/a）；B—河床寬度（m）；L—河流長度（m）；KR—河底滲透系數(shù)；hR—河流水位（m）；H—地下水位（m）；S—河底積層厚度（m）；t—河流行水時(shí)間（d）。

c.側(cè)向徑流補(bǔ)給（Q側(cè)補(bǔ)）

圖1 東鵬煤礦區(qū)潛水含水層初始水頭等值線圖

圖2 潛水含水層觀測孔水位擬合曲線

圖3 第二含水層觀測孔水位擬合曲線

由孔隙地下水等水位線圖可知，研究區(qū)南部處于淮河流域城區(qū)降落漏斗范圍內(nèi)，僅于北部可接受側(cè)向徑流補(bǔ)給。經(jīng)長期監(jiān)測資料，該地段枯水期的水力坡度稍大于豐水期，受地下水降落漏斗影響，水力坡度較大，多年平均為I=1.5‰；地下水徑流斷面寬度為B=1230m；前面已確定導(dǎo)水系數(shù)T=880m2/d。因此，多年平均狀態(tài)下側(cè)向徑流補(bǔ)給量：

Q側(cè)補(bǔ)=T·B·I=1620 m3/d

d.灌溉回滲補(bǔ)給（Q灌）

根據(jù)調(diào)查統(tǒng)計(jì)，研究區(qū)灌溉面積7km2，采用面積比擬法估算研究區(qū)的灌溉回滲量為：Q灌=707m3/d。

3.6.2 排泄項(xiàng)

a.潛水蒸發(fā)排泄

研究區(qū)潛水含水層水位埋深較淺，所以潛水蒸發(fā)是該區(qū)排泄項(xiàng)之一。東鵬煤礦區(qū)多年平均蒸發(fā)量為1651.2mm。

b.人工開采排泄

研究區(qū)內(nèi)地下水開采主要包括農(nóng)業(yè)灌溉開采和居民生活用水開采兩部分。

經(jīng)調(diào)查統(tǒng)計(jì)，研究區(qū)內(nèi)耕地面積7km2，根據(jù)灌溉回滲補(bǔ)給量的計(jì)算結(jié)果，研究區(qū)農(nóng)業(yè)灌溉用水量為3533m3/d。區(qū)內(nèi)農(nóng)業(yè)灌溉用水約50%取自地下水，則農(nóng)業(yè)灌溉開采地下水量為1766m3/d。

據(jù)統(tǒng)計(jì)研究區(qū)內(nèi)人口約3.9萬人，人均日用水量按85L考慮，則研究區(qū)內(nèi)生活用水年開采量為3315m3/d。

c.側(cè)向流出量（Q側(cè)）

受水廠地下水集中開采影響，研究區(qū)西南部地下水側(cè)向流出至礦區(qū)。該流出邊界的徑流寬度B=3500m；平均水力梯度取I=2.0‰；導(dǎo)水系數(shù)T=880m2/d。由此計(jì)算側(cè)向流出量：Q側(cè)=6160m3/d。

d.越流排泄量

依據(jù)水均衡原理，Q補(bǔ)=Q排，即地下水的排泄量應(yīng)等于其補(bǔ)給量，由此計(jì)算得出，研究區(qū)的越流排泄量為10559m3/d。

4 模型的識別與驗(yàn)證

模型的識別與驗(yàn)證過程在整個(gè)研究區(qū)地下水模擬中具有極其重要的作用。一般來說，只有對源匯項(xiàng)和參數(shù)進(jìn)行多次的調(diào)整和修改，才能使整個(gè)模型得到比較理想的效果。此次數(shù)學(xué)模型的識別采用的是試估—校正法。

根據(jù)已有的東鵬煤礦區(qū)地下水動態(tài)觀測資料，選取2011年1月1日～6月30日為此次模型研究的識別階段，兩個(gè)月為一個(gè)時(shí)間段；選取2012年6月30日～2013年1 月1日為該模型的驗(yàn)證階段。但是為了能夠真實(shí)完整地反映東鵬煤礦區(qū)地下水的變化情況，在模擬和驗(yàn)證階段，應(yīng)滿足：模型輸出的地下水位計(jì)算值需要與實(shí)測水位擬合小于0.5m的絕對誤差值應(yīng)占已知觀測井點(diǎn)的70%以上；模型輸出的模擬地下水位等值線應(yīng)該與實(shí)際的等水位線在一定的誤差范圍內(nèi)吻合。如果模型計(jì)算的水位值與實(shí)際觀測水位值相差較大，則需要重新進(jìn)行調(diào)參，直到兩者之間有較好的擬合效果，最后得到的參數(shù)就是該模型校正后的最終參數(shù)。

煤礦區(qū)選擇了具有代表性的基巖觀測孔在群孔抽水試驗(yàn)過程當(dāng)中的地下水動態(tài)曲線進(jìn)行擬合。經(jīng)反復(fù)調(diào)整水文地質(zhì)參數(shù)等各種不確定因素，直到使各動態(tài)曲線擬合程度滿足要求為止。主要含水層觀測孔地下水動態(tài)擬合曲線，見圖2～3。

根據(jù)東鵬煤礦區(qū)地下水?dāng)?shù)值模擬的識別與驗(yàn)證階段的結(jié)果分析可知：該模型觀測井計(jì)算的水位值與實(shí)際觀測到的地下水位值的差值基本符合精度要求；該模型模擬的地下水流場與煤礦區(qū)內(nèi)實(shí)際地下水流場基本吻合。即此次的模型研究能夠真實(shí)反映東鵬煤礦區(qū)研究區(qū)地下水流規(guī)律和特征。

綜上所述可知，東鵬煤礦區(qū)這個(gè)非均質(zhì)、各向同性、三維非穩(wěn)定地下水?dāng)?shù)值模型符合煤礦區(qū)的實(shí)際水文地質(zhì)條件，因此可以利用該模型對東鵬煤礦區(qū)地下水流場變化趨勢進(jìn)行預(yù)測。

5 結(jié)論

本文以東鵬煤礦區(qū)為研究對象，針對地下水取水影響，構(gòu)建地下水?dāng)?shù)值模型，進(jìn)行模型識別與驗(yàn)證，研究數(shù)值模型對地下水流場變化趨勢進(jìn)行預(yù)測的可靠性，為取用地下水影響探索新的分析模式，可為地下水可持續(xù)發(fā)展提供基礎(chǔ)技術(shù)指導(dǎo)

（專欄編輯：顧 梅）