基于GIS的地下水模型邊界地下水流估算方法

于國寶

(遼寧省朝陽水文局,遼寧 朝陽 122000)

1 概述

地下水和河流系統的生態完整性經常受到人類活動的威脅[1]。社會也面臨日益嚴重的水資源管理問題[2]。通常地下水建模是確定地下水管理策略的重要工具，然而，為了開發一個可靠的地下水模型，研究人員需要克服許多方法上的挑戰[3]。本研究論文集中于一種基于地理信息系統(GIS)的新方法，用于在地下水流模擬中估算研究區域邊界上的地下水流。

作為開發地下水流動模型的第一步，需要確定研究區域和模型邊界。模型邊界是研究區域和周圍環境之間的界面[4]。阿慧娟等[5]證明了分配正確邊界條件的重要性，并表明邊界位置的不準確性會導致校準過程中引入不必要的不均勻性。同時河流、湖泊和濕地與地下水的相互作用受水體相對于地下水流動系統的位置、河床的地質特征及其氣候環境的控制[6]。陳張建[7]全面總結了河流和地下水之間的相互作用。并認為在降雨量少的情況下，許多溪流中的基流大部分時間是來自于地下水的輸入。

GIS是為任何地下水流和污染物遷移建模工作開發地下水模型的重要工具。GIS提供數據管理和空間分析能力，可用于地下水建模。它提供數據插值、系統的模型參數分配、空間統計生成和模型結果的可視化顯示[8]。唐敏等[9]設計了一個GIS數據庫，研究了地下水脆弱性分析和水文地質建模。GIS在地下水流系統的概念化、特征化和數值模擬中也具有廣泛應用。龐園等[10]應用GIS技術對堅硬巖石地形的地下水補給進行定量建模。而且GIS還提供了許多工具，如創建剖面圖的工具、提取某一點的值、創建緩沖區的工具、插值工具等。因此，它可以和達西定律一起用來估算通過給定邊界的平均流量。

2 研究區

本研究選取遼寧省渾河的某支流作為研究對象。研究區的主要含水地層是沖積層，約占總面積的60%，其次是花崗巖和千枚巖。研究區的地下水位深度從平均海平面(MSL)的136～380m不等。MSL的地面高度從147～446m不等。該地區的平均降雨量為452mm。大部分降雨發生在6月至9月的季風季節。其余時間被認為是非季風季節。本文從水文、水文地質、地下水流方向和形態特征等方面對該區域進行了深入研究，并在此基礎上確定了地下水流的邊界條件。在研究區域的東側，是傾向NNE-SSW的山脈。可以將其視為一個重要的無流量邊界。在西側剩余的小部分以封閉研究區域的邊界。研究區域如圖1所示。選擇了研究區域內及其周圍的108個探井，如圖2所示，根據這些探井的歷史水位數據，用逆距離加權(IDW)插值法計算了1999—2003年季風和非季風季節的水位面。插值時選擇360m×360m的像元尺寸作為輸出圖像。圖2還顯示了2003年季風前地下水位的分布圖像。

圖1 研究區域

圖2 測點位置

根據獲得的28個探井記錄，繪制了研究區的水文地層圖。建造了3個沖積層、花崗巖和千枚巖的水文地層，并計算出所有三層的底部高程。ASTER高程值用于開發數字地形模型。研究區部分探井的詳細信息以及所有三層的井底高程見表1。利用鉆孔測井信息和IDW插值工具，為沖積層、花崗巖和千枚巖的水文地層底部高程準備表層。還計算了監測井所有季節的平均水位深度，并通過插值將其近似轉換為1個層。

表1 部分探井數據

3 研究方法與結果討論

本文采用了以下方法對邊界進行計算。

(1)在研究區的邊界上創建了一個5公里的緩沖區。邊界被分成22個橫截面，如圖3所示。西側的三個橫截面編號為F1、F2和F3。上側的部分從U1到U10順時針編號，底端的部分從L1到L9順時針編號。圖3顯示了其中一些部分的編號。在遇到邊界急劇彎曲的地方，應采用更多的橫截面。

圖3 橫截面示意圖

(2)利用三維分析工具，從1999—2003年的不同季風前和季風后水面剖面中提取了沿每個橫截面區域大約30—35個均勻分布點的地下水位。

(3)繪制這些點，并將四階多項式擬合到每條曲線，以確定精確邊界點的梯度，如圖4所示。觀察到四階多項式的適合度在所有情況下都非常好。2002年后季風的U4截面擬合相關性為0.998。從四階方程出發，通過簡單的微分，可以更容易地得到任意點的梯度(dy/dx)。因此，在每個這樣的輪廓的邊界點計算梯度。

圖4 2002年季風后沿U4斷面的水位變化

(4)盡管不同剖面之間的梯度值不同，但觀察到季風季節前后和不同年份之間給定剖面的梯度計算值沒有顯著變化。典型線段幾個橫截面不同年份的季風和非季風季節的梯度計算值見表2。流入研究區域邊界的梯度值保持為正，流出為負。可以明顯的發現，在季風前后各年份和各位置地下水位并沒有很明顯的變化，特別是對于L9而言，水位梯度變化值小于0.00001。

(5)因此，可以通過取所有季節所有梯度計算

表2 典型線段的梯度計算值

值的平均值來計算每個橫截面的梯度平均值。然后對所有橫截面獲得的邊界段的平均梯度值進行分析，并識別和合并具有幾乎相同梯度值的相鄰邊界段。通過這種方法，整個邊界被重新分組為8個邊界段，并計算梯度的平均值。對于這8個邊界段，基于單個線段的平均值。這些分段和梯度平均值如圖5所示。東部邊界由于靠近山脈，因此顯示值為“0”(無流量邊界)。

(6)在每個邊界段選擇5個點，并從相應的層中提取這些點的地面高度、水文地層底部高度和平均水深值。因此，在8個邊界段上選擇的40個邊界點的細節制成表格。然后用這些值計算三層含水層系統在每個點的等效橫向滲透系數，然后用達西定律計算各點邊界每米長度的凈流量，再計算每個區段的平均流量值，并乘以區段長度，計算每個區段的總凈流量，見表3，大部分的梯度值很低。梯度最大值(0.0024)發生在在4號段。用達西定律計算得到的流量結果表明單位長度下的凈流量比從-0.5040到0.9778之間變化，在第2段流入量最大，為35725m3/s。在第6段流出量最大，為29000m3/s。

表3 每個邊界段的凈流量

4 結論

本文在無法確定研究區域附近無流動邊界的情況下，構建了一種基于地理信息系統計算研究區域邊界的凈流量的方法。并認為地下水位梯度在不同季節和不同年份之間不會發生顯著變化，此假設也在地下水模型中得到驗證。在本研究中，研究區域邊界分為9個區段，并計算每個邊界區段的流量梯度值。通常，梯度值很低。梯度最大值(0.0024)發生在在4號段。然后使用達西定律來計算穿過這些邊界的流量。結果表明：在第2段流入量最大，為35725m3/s；在第6段流出量最大，為29000m3/s。本文提供了一種計算邊界平均凈流量的方法，以期對地下水污染研究和地下水模型建立提供一定的參考。

圖5 不同邊界段的最終水位梯度