基于1stOpt軟件計算Theis井函數的應用探討

梅 博，李 歡，楊云展，常林禎，張巖峰，劉普峰，刑化廬

(河北省地礦局第四水文工程地質大隊 地熱研究院，河北 滄州 061000)

進行地下水資源評價時，獲得含水層的水文地質參數不可或缺。非穩定流抽水試驗是預測地下水位變化，獲取水文地質參數的重要方法。利用非穩定流抽水試驗數據計算水文地質參數是建立在Theis井流模型基礎之上的。

定流量抽水試驗時，Theis井流模型假設條件如下：(1)含水層均質各向同性，等厚，側向無限延伸。產狀水平；(2)抽水前天然狀態下水利坡度為零；(3)完整井定流量抽水，井徑無限小；(4)含水層中水流符合Darcy定律；(5)水頭下降引起地下水從貯存量中的釋放是瞬間完成的[1]。根據上述假設條件，抽水后形成了以井軸為對稱軸的水位下降漏斗，如圖1所示。該抽水徑流歸納的數學模型如公式(1)所示。

(1)

該方程經Hankel變換，其公式變為公式(2)形式，

(2)

一般令W(u)代替式中的指數積分式，即

(3)

則(2)式改寫成經典的Theis公式：

(4)

前述所有公式中，s為抽水范圍內，任一點任一時刻的水位降深；Q為抽水井流量；T為導水系數；t為自抽水開始到計算時刻的時間；r為計算點到抽水井的距離；μ*含水層彈性釋水系數。

1 Theis公式求解方法概述

由于公式(3)沒有解析解，配線法和直線圖解法是在計算機技術沒有得到普遍應用時解決工程實際問題的經典方法，這兩種方法計算過程比較繁瑣，有很大的主觀隨意性，造成計算結果誤差很大，對于延續時間較長的抽水試驗，數據量多達幾百條，處理起來耗費大量時間和精力。多年來，國內外不少學者對其求解方式進行了分析研究，給出了井函數W(u)的近似解或者替代公式，其計算精度均控制在1%以內[2-9]。

圖1 承壓水完整井流

隨著計算及技術的發展和應用，采用智能化計算程序解決水文地質問題也得到推廣。對于常規配線法求解水文地質參數，Waterloo水文地質公司基于泰斯井函數理論、博爾頓模型理論和紐曼模型理論，開發出了專業求解軟件Aquifer Test，該軟件實現曲線自動擬合，還可進行人工可視化調參優化，兼顧自動與人工的優化效果，實現相應參數的求解[10]。基于該平臺軟件，張苗苗[11]，陶宗濤[12]，蔣輝[13]，唐輝等[14]均根據不同地方抽水試驗數據分析計算含水層水文地質參數，并與工配線法進行了對比表明AquiferTest軟件相比于人工配線法，具有快速、高效、準確以及受人為因素影響小的優點。

MATLAB是美國Math Works公司出品的可用于算法開發、數據可視化、數據分析以及數值計算商業數學軟件。廖梓龍等結合托克托縣地區抽水試驗，運用MATLAB軟件優化傳統的Theis配線法和Jacob直線圖解法在單孔非穩定流抽水試驗中確定了承壓含水層導水系數和滲透系數，計算結果符合當地水文地質條件，且比傳統方法精確可靠[15]。張勇等依托MATLAB軟件實現了人工智能算法完成水文地質參數的智能優化配線計算[16]。李林子使用MATLAB的數值積分工具和優化工具，實現了對水文地質參數的求解，所得的結果與應用配線法得到的結果基本一致，從而肯定了方法的應用價值[17]。

2 1stOpt軟件計算水文地質參數的應用效果

1stOpt是一款具有自主知識產權的的數學優化計算軟件平臺，在非線性曲線擬合領域和國際上常用的數學軟件相比具有其特點和優勢，程先云作為1stOpt軟件平臺的開發者，計算了地下水動力學中傳統降速配線法的經典案例，并且對比了其他平臺或者軟件計算結果表明1stOpt無論是精度和可操作性都有一定優勢[7]。但是將該計算軟件實際應用于水文地質參數求解相關報道幾近于無，本文作者長期從事水文地質、地熱地質野外一線工作，在地熱資源勘查等項目中多次將該軟件平臺應用于水文地質參數求解過程，對指導野外一線抽水試驗及計算求參工作取得了較好效果。

2.1 地熱井抽水試驗算例

在河北省滄州市東光縣開展的地熱井抽水試驗。地熱井利用熱儲層為新近系館陶組。抽水流量61.4 m3/h，觀測井為該抽水井配套的同層回灌井，觀測井距抽水井84.87 m。抽水試驗過程中，觀測井的時間-降深數據見表1，本次求參過程中方程采用Theis公式的簡化形式雅柯布方程。1stOpt代碼見表2。數據擬合圖見圖2。從圖2中可以出，經過計算，導水系數T=185.17 m2/d，壓力傳導系數a=187 731.91，數據擬合度較好，相關性系數R值為0.994，在工程上可以滿足需要。該方法需較少的工作量即可處理抽水延續時間較長的抽水試驗大量數據，節省了時間。

表1 地熱井抽水試驗基本數據表

表2 地熱井抽水1stOpt計算代碼

圖2 地熱井抽水試驗數據擬合圖

2.2 第四系承壓井抽水試驗算例

在河北省滄州市某處試驗場開展的一次多孔抽水試驗。抽水井深285 m，觀測井井深290 m，兩眼井利用段均為第Ⅲ含水組。抽水流量37 m3/h，觀測井距抽水井50 m。抽水試驗過程中，觀測井的時間-降深數據見表3，本次求參方程采用Theis公式。1stOpt代碼見表4，1stOpt軟件可以直接調用井函數進行計算。數據擬合圖見圖3。從圖3中可以出，相關性系數R值為0.989，也可滿足工程實際需要。

表3 第四系承壓井抽水試驗基本數據表

圖3 地熱井抽水試驗數據擬合圖

綜上，1stOpt軟件平臺無論在求解雅柯布方程和Theis方程都可以迅速得到計算結果，并且精度滿足要求，計算過程代碼比較簡單，對計算機語言的要求較低，優化算法可以解決直線圖解法和配線法的人為主觀因素影響，在實際的工程應用中可以替代常規的直線圖解法和配線法來計算水文地質參數。

表4 第四系承壓井抽水試驗1stOpt計算代碼

表5 坡度配線法算例抽水試驗數據

3 1stOpt軟件計算坡度配線數據的研究

坡度配線法又稱降深-距離配線法，其計算核心理論仍為Theis方程。一般應用于多個觀測孔參與的抽水試驗，取任意時刻的t為定值，利用多個觀測孔在時間t時刻的降深進行配線。坡度配線法精度比降速配線法更高。1stOpt軟件可以將不同時間點的多井降深數據擬合后進行全局綜合判斷，其最優解判定模型公式為各組數據擬合后的殘差之和的最小值。

本次計算過程中以地下水動力學經典案例的數據[18]進行計算：某承壓水進行了一次非穩定流抽水試驗，抽水流量為60 m3/h；觀測井為2#、15#、16#、10#、1#，距主孔距離分別為43 m、125 m、140 m、510 m、780 m。抽水試驗觀測資料如表4所示。

由于9#井數據間斷較多，計算時使用觀測時間120 min及之后的多個時間點的數據進行計算。計算過程代碼見表6。

表6 坡度配線數據試驗1stOpt計算代碼

該計算過程結果見圖4。在計算過程中，同時利用多組數據進行計算，擬合時曲線兼顧了每組數據中的每一個數據點，當每個數據點的計算數據與實測數據的殘差最小時，識別到最優解；1stOpt全局優化算法，不會舍棄抽水中的單個數據，且實現了較高的擬合度，加上可以直接調用的井函數，保證計算精度，結果比較穩定可靠。而傳統的坡度配線法只是對單個時間點的數據進行人工配線，計算繁瑣，在此基礎上若提高精度，獲得整體參數，只能對不同時間點配線后的結果再次進行平均值計算，工作量成倍增加，并且簡單的平均數計算獲得導水系數T，壓力傳導系數a，彈性釋水系數μ*建立不了地下水動力學中三個參數的關系式。Aquifer Test軟件雖然可以實現自動擬合數據，但是僅僅是針對單組數據；只是將人工配線方法實現了自動化，并未做到全局優化，多時間點和多井數據不能同時兼顧，若實現地層整體的水文參數也只能簡單求算平均數，獲取數據的精度仍不能有效得到保證。

圖4 坡度配線數據計算結果圖

4 結語

(1)1stOpt軟件是一款具有自主知識產權的的數學優化計算軟件平臺，在水文地質參數的計算求解中應用較少。針對作者開展的館陶組地熱井和第四系承壓井兩次抽水試驗的時間-降深數據進行計算，無論采用Theis公式還是其簡化形式的雅柯布方程，均能實現精確的求解，表明該軟件平臺可以得到廣泛應用，在實踐中替代傳統的降速配線法和直線圖解法。

(2)依據經典的地下水動力學案例數據，使用1stOpt軟件對人工坡度配線的公式進行計算表明1stOpt軟件全局優化算法可以充分利用每一個數據，實現較好的擬合度，加上可以直接調用的井函數，保證計算精度，結果比較穩定可靠。