GMS在灤河三角洲淺層地下水?dāng)?shù)值模擬中的應(yīng)用

張佳寧

（中國(guó)市政工程?hào)|北設(shè)計(jì)研究總院有限公司，長(zhǎng)春 130021）

海水入侵是指由于自然因素或人為作用，破壞了沿海地區(qū)地下水中的淡水與海水中的咸水之間的平衡狀態(tài)，海水（含鹽量較大或礦化度較高）沿著與其具有水力聯(lián)系的含有淡水的地下水含水層向內(nèi)陸方向擴(kuò)散，使得原淡水含鹽量升高的過(guò)程[1]。導(dǎo)致海水入侵現(xiàn)象發(fā)生應(yīng)具備三個(gè)重要條件，前提條件是區(qū)域處于干旱或半干旱地區(qū)，降水較少，水資源量不足，存在潮汐與潮流；基礎(chǔ)條件是淡水含水層與海水之間存在水力聯(lián)系，決定性因素是人為導(dǎo)致的地下水超采或截?cái)嗪畬友a(bǔ)給源，三個(gè)因素綜合起來(lái)，使得濱海地區(qū)海水入侵范圍逐步擴(kuò)大[2]。

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

研究區(qū)位于燕山山前自北向南傾斜的沖洪積平原之上，地形起伏較小。地面絕對(duì)高程在10～50 m。環(huán)渤海海岸帶，灤河—洋河沖洪積平原及濱海沖海積平原組成了研究區(qū)的主要地貌，由海積和海積-沖積共同堆積共同作用形成了灤河—洋河沖洪積平原，地面高程在2～5 m。在河相沉積及海相沉積共同作用下，大多低洼地區(qū)排水不暢，加上海潮托頂，該區(qū)域鹽漬化較為普遍，澇災(zāi)頻發(fā)。海積平原絕對(duì)高程在1～3 m，地形起伏較小，坡比總體小于0.1‰。海積沖積平原地處上述兩種地貌單元之間，地面絕對(duì)高程3～5 m，主要不良地質(zhì)現(xiàn)象為土壤鹽漬化和沼澤化。根據(jù)不同巖性，海岸帶可以分為淤泥質(zhì)、砂質(zhì)及基巖海岸等類型。海岸帶地表出露巖性自東向西逐漸由細(xì)粒過(guò)度到粗粒。澗河口—大清河口多為淤泥質(zhì)海岸，地勢(shì)極低，岸線外潮間帶發(fā)育，并殘留一些砂島，鹽田、養(yǎng)蝦池遍布[3]。

本區(qū)基底為太古界和下元古界變質(zhì)巖系。第四紀(jì)沉積物主要為灤河水系建造，由一套陸相碎屑沉積物組成，濱海地區(qū)夾有2～6層海相層。Q1～Q4地層連續(xù)沉積，沒(méi)有明顯的沉積間斷。

評(píng)價(jià)區(qū)地下水主要由潛水及巖溶水組成，其中潛水主要賦存于第四系松散孔隙之中，巖溶水含水介質(zhì)為奧陶系灰?guī)r，由于研究區(qū)內(nèi)主要工業(yè)、農(nóng)業(yè)及生活用水來(lái)自于第四系松散地下水的開(kāi)采，因此第四系含水層為此次研究的目標(biāo)含水層。研究區(qū)主要由灤河補(bǔ)給的地下水系統(tǒng)，東至渤海西岸，西部主要通過(guò)黏性土組成的弱透水層與其下白河-薊運(yùn)河作為主要補(bǔ)給源的含水層在空間上進(jìn)行了相接。根據(jù)地下水含水層內(nèi)水質(zhì)不同，又可分為沖洪積淡水系統(tǒng)與沖海積咸水系統(tǒng)。

通常情況下，相對(duì)不透水層將上部孔隙潛水與下部巖溶水在空間上進(jìn)行了分割，二者之間幾乎不存在水力聯(lián)系。但是，隨著人類活動(dòng)的日益增多，研究區(qū)地下水得到廣泛利用（廣泛采用通天濾水管），致使天然隔水層作用減弱，發(fā)生兼之復(fù)雜的巖相變化，導(dǎo)致上下兩含水層水力聯(lián)系較為密切?？傮w上，大氣降水為研究區(qū)淺層地下水的主要補(bǔ)給來(lái)源，深層地下水主要為淺層地下水越流補(bǔ)給。

研究區(qū)地下水主要補(bǔ)給方式有大氣降水補(bǔ)給、河流側(cè)滲補(bǔ)給、地下水側(cè)向徑流補(bǔ)給及人工灌溉回滲補(bǔ)給。天然條件下地下水流場(chǎng)主要受研究區(qū)地形地貌影響，但研究區(qū)長(zhǎng)期大面積人工開(kāi)采地下水使得流場(chǎng)形態(tài)變化較大，第四系松散孔隙潛水主要滲流方向?yàn)樽员毕蚰?，從北部山區(qū)指向南部濱海平原，滲流方向于區(qū)內(nèi)灤河方向近似平行，水力梯度不大于2‰。在地形坡度與地下水賦存條件雙重影響下，南北兩側(cè)地下水水質(zhì)及其開(kāi)發(fā)利用程度有很大差異。模擬區(qū)北部屬低山丘陵區(qū)，含水介質(zhì)以圓角礫為主，顆粒較粗，地形坡度起伏較大，徑流條件較好，水質(zhì)較好，人為開(kāi)發(fā)利用程度相對(duì)較高；南部屬于濱海平原區(qū)，含水層介質(zhì)以粉細(xì)砂為主，地勢(shì)起伏較小，徑流滯緩條件極差，水質(zhì)較差，由于咸水原因，開(kāi)發(fā)利用程度較低。

根據(jù)地下水含鹽量，可將淺層地下水劃分為全淡水區(qū)和咸水區(qū)。全淡水區(qū)內(nèi)孔隙潛水為研究區(qū)內(nèi)生活飲用水主要來(lái)源，也是農(nóng)業(yè)灌溉的主要開(kāi)采層，開(kāi)采量較其他含水層大。隨著經(jīng)濟(jì)發(fā)展，深層巖溶水開(kāi)采程度逐漸加大，第四系松散孔隙潛水水位逐漸降低，開(kāi)采強(qiáng)度也日益增加，導(dǎo)致研究區(qū)內(nèi)局部出現(xiàn)小范圍的降落漏斗。在全咸水地區(qū)內(nèi)，咸水、鹵水在地下水中所占比例較大，開(kāi)發(fā)利用程度較差。受地形地貌及含水介質(zhì)滲透性的影響，南北部淺層地下水側(cè)向徑流也存在較大差異。

2 水文地質(zhì)模型

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

根據(jù)巖性及透水性不同，將研究區(qū)內(nèi)松散巖層劃分四層，自上而下依次為四個(gè)相對(duì)獨(dú)立的含水層。各含水層組之間存在相對(duì)不透水的黏土層，其分布較為分散，不連續(xù)，局部存在天窗，因此其無(wú)法構(gòu)成完整的隔水層。本次數(shù)值模擬采用了準(zhǔn)三維非穩(wěn)定地下水流模型。模型北邊界是根據(jù)項(xiàng)目工作區(qū)范圍劃定的邊界，垂直于區(qū)域多年平均等水位線，與邊界流線相重合，設(shè)定為二類零流量邊界。模型東邊界定在灤河西岸。模型南部邊界劃定為渤海北部海岸，海水與內(nèi)陸地下水含水層之間存在水量交換微弱的黏土層，將其概化為一類水頭邊界。根據(jù)此次研究的范圍，人為劃定了模型的西部邊界，由于該區(qū)域資料較為詳細(xì)，因此概化為二類零流量邊界。灤河位于模擬區(qū)的西北側(cè)，因此將其概化為一類水頭邊界。研究區(qū)下部為中風(fēng)化基巖，其透水性極差，構(gòu)成了此次模擬的下部邊界，將其視為隔水底板，模型向下流量為零。選取潛水頂面為上部邊界，其為開(kāi)放型邊界，參與研究區(qū)內(nèi)水文循環(huán)。

模擬區(qū)內(nèi)含水層的主要補(bǔ)給源主要有：以面狀形勢(shì)補(bǔ)給的大氣降水入滲、以線狀形式補(bǔ)給的河渠側(cè)向滲漏、以面狀形勢(shì)補(bǔ)給的農(nóng)田灌溉回滲等，由于模擬區(qū)內(nèi)內(nèi)河渠密布，縱橫交錯(cuò)，無(wú)法將其清晰區(qū)分，因此此次模擬將其視為面狀補(bǔ)給。

地下水流的數(shù)學(xué)模型如下：

式中：W為流入?yún)R或源項(xiàng)的水量；Ss為儲(chǔ)水率（承壓含水層），（給水度）潛水含水層；Kh為水平透系數(shù)；Kv為垂直透系數(shù)。

對(duì)于該模型來(lái)說(shuō)，式（2）為初始條件，式（3）為第一類邊界條件，式（4）為第二類邊界條件，式（5）為第三類邊界條件。

2.2 地下水流數(shù)值模型及識(shí)別

此次模擬利用GMS的Grid frame模塊采用六面體的矩形網(wǎng)格剖分法，將模擬區(qū)剖分為150行、150列、2層。根據(jù)評(píng)價(jià)區(qū)的主要巖性及地下水賦存條件，將淺層潛水含水層分為6個(gè)分區(qū)進(jìn)行水文地質(zhì)參數(shù)賦值，將深層潛水含水層分為3個(gè)分區(qū)進(jìn)行水文地質(zhì)參數(shù)賦值，這些經(jīng)驗(yàn)數(shù)據(jù)主要來(lái)自資料收集，基本可以反映研究區(qū)的實(shí)際情況。灤河三角洲數(shù)值模擬采用的時(shí)間單位是天，擬合時(shí)間為4年。采用這種時(shí)段劃分，模型計(jì)算水位是天平均值，與觀測(cè)水位吻合。根據(jù)搜集相關(guān)區(qū)域的入滲系數(shù)，將模擬區(qū)按4個(gè)分區(qū)進(jìn)行降水入滲補(bǔ)給系數(shù)的賦值。區(qū)內(nèi)小型河流及渠系成網(wǎng)狀密布，在區(qū)域性地下水流場(chǎng)中表現(xiàn)為面狀滲漏補(bǔ)給，只能以面狀入滲的方式進(jìn)行賦值。此次模擬將灌溉較為集中的4-6月的總?cè)霛B量平均分配到各個(gè)時(shí)段，將其作為模擬的初始流量，本次模擬是在保持深層開(kāi)采總量不變的條件下，對(duì)淺層及深層地下水開(kāi)采量的空間位置進(jìn)行了重新分配。

2.3 地下水流數(shù)值模型校正及水位預(yù)測(cè)

此次模擬選取多年水位觀測(cè)資料中的一個(gè)較為完整水文年的觀測(cè)數(shù)據(jù)，同時(shí)繪制時(shí)間-水位變化曲線圖，通過(guò)分析可以發(fā)現(xiàn)研究區(qū)內(nèi)部分水位觀測(cè)孔內(nèi)水位動(dòng)態(tài)變化曲線與區(qū)域地下水補(bǔ)給及排泄變化曲線相關(guān)性較好，整體趨于一致，較好地反映了模擬區(qū)地下水的年內(nèi)變化特征。但是，部分觀測(cè)井受人為地下水開(kāi)采影響，其地下水位動(dòng)態(tài)變化曲線與區(qū)域地下水補(bǔ)給及排泄變化曲線相關(guān)性較差，局部出現(xiàn)矛盾，不能用以反映地下水年內(nèi)變化特征，對(duì)于受附近開(kāi)采井影響的水位觀測(cè)孔，模型擬合時(shí)在觀測(cè)孔附近設(shè)置了抽水井。

結(jié)合研究區(qū)近期、中期、遠(yuǎn)期的地下水水位預(yù)測(cè)結(jié)果，筆者發(fā)現(xiàn)，研究區(qū)地下水流場(chǎng)基本穩(wěn)定，由于該區(qū)域地下水補(bǔ)排狀態(tài)呈負(fù)均衡，將導(dǎo)致降落漏斗的中心水位會(huì)進(jìn)一步加大降深，降落漏斗面積亦呈小幅度增長(zhǎng)，將會(huì)導(dǎo)致研究區(qū)由于區(qū)域降落漏斗引發(fā)的環(huán)境地質(zhì)問(wèn)題進(jìn)一步加劇，因此需針對(duì)這一問(wèn)題，在水資源利用方面做出合理分配。

3 控制海水入侵問(wèn)題的取水方案優(yōu)化

據(jù)調(diào)查，研究區(qū)內(nèi)主要河流等均發(fā)育在斷裂之上，各條河流直接入海，在其上游均修有水庫(kù)，從而減少了入海河水水量和泥砂量。另外，河流的下游入海河段，河床低而平緩，在漲潮時(shí)海水沿河道上溯，尤其是遇有天文大潮和風(fēng)暴潮時(shí)，其上溯距離更遠(yuǎn)（一般上溯3～4 km，最遠(yuǎn)出現(xiàn)在洋河可達(dá)10 km），充滿河床的海水與河水之混合水沿河床周邊滲漏補(bǔ)給地下水，使沿河兩岸地下水受到污染。因此，近海松散的砂礫質(zhì)含水層及構(gòu)造破碎帶與通海河流的砂礫質(zhì)河床是本區(qū)海水入侵的主要通道[4]。咸水分布區(qū)主要分布于研究區(qū)東南及南部，面積約3511 km2，從區(qū)域上看，自東北向西南，地下水的底板埋深同咸水體的厚度均呈加大的趨勢(shì)，而沿海地區(qū)分布最廣，其厚度自北向南、自西北向東南逐漸增大。將2015年年末水位與底面沉降速率等值線圖及地面沉降等值線圖進(jìn)行配準(zhǔn)，得出海水入侵范圍與地下水位降落漏斗分布特征。據(jù)分析，評(píng)價(jià)區(qū)東南部地區(qū)海水入侵情況相對(duì)較為嚴(yán)重，咸水體的分布位置與區(qū)內(nèi)的淺層地下水集中開(kāi)采中心基本吻合。據(jù)此可推測(cè)，研究區(qū)南部地下水降落漏斗為形成該區(qū)海水入侵的主要誘因，通過(guò)調(diào)整地下水開(kāi)采方案可以達(dá)到控制海水入侵范圍的目的。

本次地下水開(kāi)采方案的優(yōu)化以海水入侵范圍縮小或保持現(xiàn)有海水入侵范圍不再進(jìn)一步擴(kuò)大為目的，通過(guò)對(duì)海水入侵區(qū)域的空間演化進(jìn)行分析，以控制研究區(qū)南部淺層地下水降落漏斗降深。在減少現(xiàn)有地下水集中開(kāi)采區(qū)地下水開(kāi)采量的同時(shí)，在研究區(qū)中部及西部地下水未超采地區(qū)建立分散型地下水開(kāi)采水源地。

從以上優(yōu)化方案模擬的結(jié)果來(lái)看，2020年即優(yōu)化方案實(shí)施5年間原地下水位降落漏斗面積和中心水位降深變小，地下水流場(chǎng)開(kāi)始變化。至2030年原區(qū)域降落漏斗基本消失，形成新的面積及中心水位降深較小的多個(gè)小型降落漏斗，基本形成內(nèi)陸地區(qū)地下水位高于海水位，地下水流向由原來(lái)的海水向淺層降落漏斗流動(dòng)，基本改為內(nèi)陸地區(qū)向海洋排泄地下水。所以，該方案在改變現(xiàn)階段由于地下水超采造成的區(qū)域降落漏斗問(wèn)題上是可行的。

圖1 研究區(qū)2020年、2030年淺層地下水水位預(yù)測(cè)圖

4 結(jié)論

本文通過(guò)數(shù)值模擬方法，應(yīng)用GMS地下水流數(shù)值模擬軟件，針對(duì)冀東地區(qū)目前存在的與水資源開(kāi)發(fā)利用有關(guān)的海水入侵這一環(huán)境地質(zhì)問(wèn)題，通過(guò)模擬不同的水資源利用方案，得到了有利于改善當(dāng)前地下水超采問(wèn)題的水資源調(diào)度方案。但是，本次模擬仍然存在一些問(wèn)題。

本次模擬中監(jiān)測(cè)孔內(nèi)地下水水位變幅計(jì)算值整體小于地下水水位變化的實(shí)際值，未能準(zhǔn)確反映研究區(qū)黏性土釋水、儲(chǔ)水系數(shù)大的特點(diǎn)。此次模擬將模型西側(cè)概化為一類水頭邊界，屬于靜態(tài)邊界，不能很好地反映模擬區(qū)內(nèi)水位動(dòng)態(tài)變化特征；本次地下水?dāng)?shù)值模擬采用準(zhǔn)三維地下水流模型，計(jì)算值可以近似代替各含水層之間的越流補(bǔ)給量，但不能有效刻畫相對(duì)不透水層內(nèi)的水頭分布特征。準(zhǔn)三維模型雖然可以很好地反映該區(qū)域淺層含水層的水位變化，但是淺層對(duì)深層的越流補(bǔ)給、越流及壓縮釋水過(guò)程都不能得到較為準(zhǔn)確的刻畫；當(dāng)前收集的水文地質(zhì)資料，主要針對(duì)各主要含水層，對(duì)相對(duì)不透水的黏土層研究較少。由于目前資料有限，難以滿足高精度的模擬要求，在模型反演校正過(guò)程中，模型對(duì)水平滲透系數(shù)的識(shí)別精度不高。