基于VBA和群井疊加原理的水文地質參數反演和降深預測

劉 軍

(水利部新疆維吾爾自治區水利水電勘測設計研究院，新疆 烏魯木齊 830000)

水文地質參數獲取是地下水資源評價中的重點之一，這些參數包括滲透系數、給水度和貯水系數等，而抽水試驗是獲取水文地質參數的最簡單方法之一。將抽水試驗按照抽水井數量可分為單井抽水和群井抽水試驗；單井定流量抽水試驗往往可采用穩定流的Dupuit公式或非穩定流的Theis公式計算含水層參數[1- 3]，在穩定流的單井抽水試驗中(觀測井小于1個)，需要引入影響半徑的經驗公式與Dupuit公式聯立方程組才可計算出滲透系數，在非穩定流的單井抽水試驗可通過Theis公式的Jacob直線圖解法求取含水層滲透參數；而群井抽水試驗通常采用數值模擬的方法進行參數反演及水位預測[4- 14]，但尚未較好地解決這一問題。

在基坑排水和水源地開采時往往是群井同時抽水，且各孔之間距離較近，這種情況下各井周圍含水層的水文地質條件差異較小，可視含水層為均質的，因此各井的水文地質參數也相同，在這種條件下，可以將含水層概化為均質無限含水層。本文基于烏魯木齊地鐵2號線群井抽水試驗結合疊加原理提出了一種參數反演的方法，并利用VBA編程預測群井附近含水層在任意時刻由于抽水導致的地下水位降深，再根據抽水前的初始流場即可得到抽水后任意時刻的流場。該方法的提出為水文地質參數獲取及基坑排水等問題提供了新思路。

1 材料及研究方法

本文主要利用單井抽水試驗和疊加原理提出了水文地質參數反演與水位預測的方法。具體內容如下：

1.1 單井抽水的Theis公式

在均質無限含水層的單井定流量抽水試驗中，抽水井周圍任意一點i的降深可用Theis公式[15]表示為：

(1)

式中，t—抽水試驗時間，min；si—任意點i在t時刻的降深，m；Q—抽水井的涌水量，m3/d；T—含水層導水系數，m2/d；ri—任意點距抽水井的距離，m；μ*—貯水系數。

1.2 疊加原理

由線性偏微分方程和線性定解條件組成的定解問題，可以采用疊加原理[15]，而承壓含水層的群井抽水試驗過程中，地下水位在承壓含水層頂板以上時，地下水運動的偏微分方程滿足疊加原理的條件。當承壓含水層中同時存在2口井同時定流量抽水時，如圖1所示，抽水任意一點i處產生的降深是兩口抽水井分別抽水時在i點處產生的降深疊加；將這一理推廣到計算N口抽水井同時定流量抽水時在任意一點i處的降深，可表達為：

圖1 群井疊加原理示意圖[15]

(2)

式中，N—抽水井個數；QN—第N口抽水井的涌水量，m3/d；rN—第i點與第N口抽水井之間的距離，m。

2 參數反演與水位預測

2.1 參數反演方法探討

單井抽水試驗可采用Theis公式的配線法將觀測孔數據與標準曲線擬合，從而計算含水層參數。根據這一特性，可將群孔抽水試驗所有孔的觀測降深與公式(2)的計算值擬合，從而反推導水系數T和貯水系數，這樣能夠很好的避免由于個別觀測孔數據異常而對結果產生較大影響。參數反演的具體方法如下：

(1)確定抽水孔個數N及各抽水孔涌水量QN，繪制各觀測孔的降深-時間散點圖。

(2)預先給定任意大于零的導水系數T和貯水系數μ*(根據經驗賦值能夠提高參數反演效率)，根據公式(2)計算各觀測孔(含有觀測數據的抽水孔)的疊加降深，繪制隨參數改變的動態降深-時間曲線。

(3)通過不斷調整導水系數T和貯水系數μ*，以達到所有觀測孔觀測到的降深-時間散點圖與計算的降深時間曲線達到最佳的擬合狀態，此時導水系數T和μ*就確定了。在Excel中通過公式和函數功能即可實現觀測孔曲線擬合，此處不再贅述。

采用這種方法擬合參數，個別觀測孔的數據可能無法完全擬合，這是因為在擬合過程中各個觀測孔的標準曲線是根據一個綜合的平均參數(滲透系數和貯水系數)計算出來的，因此擬合出的參數也是整個含水層的平均值，避免因個別井的含水層結構突變而計算出的參數不具有代表性。

2.2 水位降深預測方法

根據擬合后的參數及公式(2)，可計算群井抽水試驗中任意一點i(x，y)在t時刻的降深s(x，y，t)，因此只要給定一個范圍，即可得到該范圍內在任意時刻的降深或流場，其具體方法如下：

(1)確定一個所需要計算的流場范圍，需要給出橫坐標范圍(xmax,xmin)和縱坐標范圍(ymax,ymin)，設定一個水位計算點的間隔距離Δx和Δy，根據范圍和間隔采用VBA代碼提取出所有計算點的坐標值，計算點的間隔距越小插值出的降深場(流場)精度越高，但計算速度越慢。

(2)選定所需計算流場的時間點t，根據公示(2)采用VBA代碼實現計算所有計算點i(x，y)在各時間點t時的降深值。

(3)采用ArcGIS軟件可將各時刻的降深計算結果插值得到降深場，若需要繪制流場，將初始流場與降深場柵格化疊加即可得到各時刻的流場。

2.3 VBA程序實現

提取特定范圍內的計算點坐標i(x，y)可采用VBA代碼實現的具體代碼見表1。

表1 提取計算點坐標的代碼

提取出各計算點坐標后，給定一個抽水時間t，即可根據公式(2)在Excel中計算出該時刻各計算點的降深s(x，y，t)，通過插值軟件即可獲取任意時刻的降深場。

3 實例

3.1 含水層概況

本次研究選取了烏魯木齊市軌道交通2號線某試驗場地的承壓含水層，該含水層平面面積40000m2，含水層厚度在20～30m，壓力水頭高出地表2～3 m，含水層巖性單一；上覆隔水頂板厚度在20～25m左右；下層有穩定的隔水底板，含水層周圍無特殊構造，因此可將該含水層概化為均質、含水層水平無限延伸的承壓含水層，含水層結構示意圖如圖2所示。

圖2 含水層結構示意圖

3.2 群井抽水試驗概況

該含水層區域內有4口抽水井和4個觀測孔，抽水期間4口抽水井同時抽水，抽水孔均采取定流量抽水，各孔的涌水量見表2。本次抽水試驗時長為300min，4個觀測孔同時觀測水位并記錄觀測數據，各試驗井(孔)參數見表2。

表2 抽水試驗井(孔)參數表

3.3 參數反演

圖3 觀測孔降深擬合曲線圖

由擬合結果可看出，4個觀測孔的觀測降深與群井疊加原理的計算降深擬合結果較好，表明該方法能夠有效的反演水文地質參數，試算得到的導水系數和注水系數可視為該含水層的平均值，可用該參數進行流場預報。

3.4 水位降深預測

本次預報情景為4口抽水井同時開采，在每口涌水量均為2500m3/d的條件下，預報抽水30d內的水位變化情況，各觀測孔30d水位降深變化圖如圖4所示，抽水后1、10、30d降深場如圖5所示。由圖4觀測孔的降深變化曲線可以看出，隨著不斷地抽水，水位逐漸下降，水位下降速率逐漸變小。由圖5抽水試驗的不同時間點的降深場可以看出，抽水井群的降落漏斗連成一片，漏斗降深逐漸變大，等降深線逐漸向外擴散，這種特征滿足Theis公式的特點，也符合實際的抽水情況。

圖5 預測降深場(單位：m)

圖4 觀測孔降深變化曲線圖

4 結語

結合VBA程序和水文地質學中的群井疊加原理提出一種反演水文地質參數的方法，將該方法應用于烏魯木齊市地鐵2號線某站點的承壓含水沉群井抽水試驗中，得到以下結論：

(1)基于群井疊加原理反演水文地質參數的方法能夠有效利用所有觀測孔數據，同時可避免異常孔的數據。該方法計算出的水文地質參數為某一區域含水層的平均值，適用于小范圍內的基坑排水和水源地抽水條件下的求參。通過實例驗證了群井疊加原理反演水文地質參數的可靠性。

(2)基于VBA程序能夠有效提取設定范圍和設定密度的計算點坐標，結合群井疊加原理能夠預測不同抽水情景(定流量)下、不同時刻的降深場，通過對實例應用，預測結果較好。基于這一特征，還可將該方法用于解決“給定降深求解涌水量的問題”，為地下施工排水問題提供一種新的解決思路。

(3)群井疊加方法反演參數是基于Theis公式反演參數，因此該方法僅適用于水文地質結構較簡單的情況，含水層需要滿足Theis公式的假設條件。而在含水層結構復雜的條件下，應當采用數值法解決問題。

(4)與數值模擬相比，群井疊加的VBA方法預測流場操作簡便，結合其它插值軟件可繪制降深場或流場，在預測水位的過程中應當注意，當水位降至承壓頂板以后，地下水運動不再遵循承壓含水層的Theis公式，因此也不可使用于群井疊加方法預測水位降深。