经典地下水抽水试验水文地质参数解算方法

Calculation Method for Hydrogeological Parameters of Classic Groundwater Pumping Test

摘要: 地下水分潜水及承压水，抽水试验分稳定流及非稳定流抽水试验。随着我国建设工程的迅速发展，基坑工程规模越来越大，基坑深度越来越深，基坑降水对周边环境的影响越来越大，对查明场地水文地质条件及水文地质参数提出了更高的要求，通常需要进行抽水试验以查明含水层的水文地质参数。探索利用现代计算工具，对经典地下水稳定井流及Thies非稳定井流抽水试验进行水文地质参数解算具有一定意义。稳定流解算较为简单，可用Excel线性回归方法进行求解。非稳定流解算方法较为复杂，数值求解需采用非线性及最优化方法。Theis井函数计算，利用变量代换消除了被积函数奇点，将无限积分截断为有限积分，积分截断误差 < 10⁻¹⁵。Theis非稳定井流水文地质参数的反演求解，使用Mathcad数学软件求解局部最优解，利用1stOpt最优化专业软件求解全局最优解，两者最优化解算结果非常接近。

Abstract: The underground water consists of phreatic water and confined water, and pumping tests are divided into stable flow and unsteady flow pumping tests. With the rapid development of construction projects in our country, the scale of foundation pit engineering is getting larger and the depth of foundation pits is getting deeper. The impact of foundation pit dewatering on the surrounding environment is increasing, which puts higher requirements on identifying the hydrogeological conditions and parameters of the site. Usually, pumping tests are needed to determine the hydrogeological parameters of the aquifer. Exploring the use of modern computational tools for hydrogeological parameter calculation in classic groundwater stable well flow and Thies unsteady well flow pumping tests has certain significance. Stable flow calculation is relatively simple and can be solved using Excel linear regression method. The method for solving unsteady flow is relatively complex, and numerical solutions require the use of nonlinear and optimization methods. The well function calculation utilizes variable substitution to eliminate singular points in the integrand function, truncating infinite integrals into finite integrals with an integration truncation error of < 10⁻¹⁵. The inversion solution of the hydrogeological parameters of unstable well water flow is carried out by using Mathcad mathematical software to solve the local optimal solution and 1stOpt optimization professional software to solve the global optimal solution. The results of the two optimal solutions are very close.

1. 引言

随着我国经济建设的迅猛发展，基坑工程对周边环境的影响越来越大，对岩土勘察工作提出了更高的要求，需着重查明场地及其周边的水文地质条件及含水层的水文地质参数。由此，需要进行场地含水层的抽水试验。

地下水可分为潜水及承压水，抽水试验可分为稳定流试验及非稳定流试验。抽水试验水文地质参数的反演解算是岩土勘察工作者的一项重要技能。

19世纪80年代前，抽水试验反演解算水文地质参数，基本为手工方法。对稳定流，用两两孔，或s—logr人工直线图解法解算；对非稳定流，用s—logt人工直线图解法，或logs—logt手工配线法解算。两两孔解算，由于水位观测误差的影响，不同孔组解算的参数相差较大；人工直线图解法，也有一定的人为误差；手工配线法，人工随意性较大。

随着电脑的普及，抽水试验水文地质参数反演新技术应运而生。本文讨论了Thies井函数的计算，经典稳定井流及非稳定井流水文地质参数解算的基本原理与方法，利用现代计算工具Excel、Mathcad及1stOpt进行最优化参数求解，并用教科书经典实例予以解算验证，取得了良好效果。

2. 稳定流抽水试验水文地质参数解算方法

2.1. 稳定流抽水试验解算原理

在地下水运动中，当渗流为层流且符合达西定律时，其过水断面的渗流量为：

$Q=FV=F\cdot k\cdot I$ (1)

式中；$I=\frac{dh}{dl}$ 水力坡度

$h$ ——地下水水位高度

$l$ ——渗流长度

$F$ ——过水断面

$k$ ——渗透系数

对于抽水井流，当渗流达到稳定时，过任意径距$r$ 断面(详见图1)，其渗流量为[1] [2]：

(a) (b)

Figure 1. Groundwater stable flow pumping test diagram. (a) Phreatic aquifer; (b) Confined aquifer

图1. 地下水稳定流抽水试验图。(a) 潜水含水层；(b) 承压含水层

$\{\begin{array}{l}潜水：Q=2\pi rh\cdot k\cdot \frac{dh}{dr}=\pi k\cdot \frac{d{h}^2}{d\ln r}=\pi k\cdot \frac{\Delta h^2}{\Delta \ln r}\\ 承压水：Q=2\pi r\text{M}\cdot k\cdot \frac{dh}{dr}=2\pi k\text{M}\cdot \frac{dh}{d\ln r}=2\pi k\text{M}\cdot \frac{\Delta h}{\Delta \ln r}\end{array}$ (2)

从而有：

$\{\begin{array}{l}潜水：V=\Delta h^2=H_0^2-h^2=\frac{\text{Q}}{\pi k}\cdot \ln \frac{R}{r}\\ 承压水：s=\Delta h=H_0-h=\frac{Q}{2\pi k\text{M}}\cdot \ln \frac{R}{r}\end{array}$ (3)

式中：$H_0$ ——抽水前的初始水位

$h$ ——地下水水位高度

$s$ ——水位降深

$R$ ——地下水影响半径(圆柱岛状含水层半径，补给半径)

式(3)表现为$V—\ln r$ 或$s—\ln r$ 的线性方程。当抽水试验有2个以上(含2个)观测孔时，可用线性回归确定渗透系数$k$ 及影响半径$R$ 。回归方法可最大限度地消除由两两孔解算的参数误差。

潜水：$\{\begin{array}{l}V=c+m\cdot \ln r\\ m=-\frac{Q}{\pi k}，\Rightarrow k=-\frac{Q}{\pi \cdot m}\\ c=-m\cdot \ln R，\Rightarrow R=\exp \left(-\frac{c}{m}\right)\end{array}$ (4)

承压水：$\{\begin{array}{l}s=c+m\cdot \ln r\\ m=-\frac{Q}{2\pi kM}，\Rightarrow k=-\frac{Q}{2\pi M\cdot m}\\ c=-m\cdot \ln R，\Rightarrow R=\exp \left(-\frac{c}{m}\right)\end{array}$ (5)

式中：$m$ ——线性回归方程的斜率

$c$ ——线性回归方程的截距

其中，线性回归方程的斜率$m$ 及截距$c$ 可用Excel函数求解：

$m=\text{slope}\left(\text{y\_data, x\_data}\right)$

$c=\text{intercept}\left(\text{y\_data, x\_data}\right)$

2.2. 稳定井流解算实例

2.2.1. 潜水稳定井流解算

【实例1】选自：文献[3]，p. 40，例题6-3，表6-9；文献[4]，p. 307，例2，表15-2。

1961年山东某地小沽河谷冲积层中进行了一组多孔3降次稳定流抽水试验。河谷宽190 m，潜水含水层厚6.22 m，抽水井位于岸边，附近河流发生了弯曲，环绕抽水井呈半圆弧状河湾。抽水井，距河床20多m，井径0.250 m，观测孔平行于河流，1#及2#观测井距抽水井距离分别为5 m及50 m。抽水共进行了5个降次，最后3个降次的抽水试验数据详见表1。试求含水层渗透系数及影响半径。

由线性回归分析解算的水文地质参数详见表2及图2。解算的影响半径是河流补给距离的综合反映。

Table 1. Data of stable flow pumping test for phreatic water

表1. 潜水稳定流抽水试验数据

	孔号	0#	1#	2#	地下水位高度
井距	r/m	0.125	5	50	h
降次	Qm^3/d	s0/m	s1/m	s2/m	H0/m	h1/m	h2/m
1	3243	1.52	0.46	0.14	4.7	5.76	6.08
2	3517	1.78	0.54	0.19	4.44	5.68	6.03
3	4050	2.13	0.62	0.2	4.09	5.6	6.02

Table 2. Hydrogeological parameter calculation table for stable flow pumping test in phreatic water

表2. 潜水稳定流抽水试验水文地质参数解算表

	孔号	0#	1#	2#					文献解
降次	r/m	0.125	5	50	本文线性回归解				[3]文解	[4]文解
	Qm^3/d	V0	V1	V2	m	c	km/d	Rm	km/d	km/d
1	3243	16.5984	5.5108	1.7220	−2.5317	10.8485	407.74	72.61	341.34	345
2	3517	18.9748	6.4260	2.3275	−2.8367	12.4974	394.65	81.91	319.09	326
3	4050	21.9603	7.3284	2.4480	−3.3229	14.3914	387.96	76.02	321.64	336

表2中，文献解的渗透系数k，为观测孔解算值。

回归分析解是整体最小误差解，消除了两两孔解算的参数误差，解算值较为客观、优化。

Figure 2. Linear regression plot of v—lnr for stable flow pumping test in phreatic water

图2. 潜水稳定流抽水试验v—lnr线性回归图

2.2.2. 承压水稳定井流解算

【实例2】选自：文献[3]，p. 36，例题6-4，表6-2；文献[5]，p. 211，例2，表5-3-3。

1965年，在陕西咸阳渭河二级阶地冲洪积承压含水层中作了一次多降次多孔抽水试验。抽水试验场地位于渭南华家寨渭河北岸，含水层厚25 m，设置了A线及B线2组观测孔。典型剖面及抽水试验数据详见图3及表3。求含水层的渗透系数及影响半径。

Figure 3. Sectional view of observation hole for Huajiazhai pressurized water pumping test

图3. 华家寨承压水抽水试验观测孔剖面图

Table 3. Test data of pressurized water pumping in Huajiazhai

表3. 华家寨承压水稳定流抽水试验数据

观测孔	孔组		主井	A线组		B线组
	孔号		5	A2	A3	B1	B4	B5
	孔距r/m		0.20	17.34	60	5.3	60	300
阶段	降次	Qm^3/d	s5	sA2	sA3	sB1	sB4	sB5
抽水	1	5530	10.25	3.76	2.51	4.74	2.22	0.64
	2	4088	7.44	2.85	1.88	3.52	1.69	0.50
	3	1402	2.40		0.65	1.17	0.57	0.16
恢复	1	5530	10.25	3.76	2.51	4.72	2.22	0.64
	2	4088	7.40	2.80	1.82	3.47	1.61	0.49
	3	1402	2.35	0.99	0.65	1.15	0.52	0.17

由于5号主井存在较大水跃，本文仅选用观测井进行线性回归分析。线性回归分析解算的水文地质参数值详见表4及图4。

原文献渗透系数$k$ ，采用两两观测孔求算，然后计算平均值。文献[5]的影响半径$R$ ，3个降次的s—lnr直线收敛于$R$ 。

本文线性回归解算的参数为A线组及B线组的整体最小误差解，消除了原文献两两观测孔解算的参数误差，解算值较为客观、优化。

Table 4. Calculation of hydrogeological parameters of pressure water in Huajiazhai

表4. 华家寨承压水稳定流水文地质参数解算表

试验	抽水量		本文线性回归解				[3]文献A组解 B组解				[5]文献A组解 B组解
阶段	降次	Qm^3/d	m	c	km/d	Rm	kcpm/d	Rcpm	kcpm/d	Rcpm	kcpm/d	Rcpm	kcpm/d	Rcpm
抽水	1	5530	−1.0280	6.5598	42.81	590.53	33.73	594.00	34.65		42.19		42.19
	2	4088	−0.7606	4.8891	42.77	618.78	34.00	623.24	34.75		41.59	660	43.02	510
	3	1402	−0.2479	1.6019	45.00	639.98	35.91	592.71	35.49		41.89		42.09
恢复	1	5530	−1.0236	6.5396	42.99	595.08					42.19		43.84
	2	4088	−0.7536	4.8133	43.17	594.06					41.36		41.36
	3	1402	−0.2522	1.6247	44.24	627.91					42.09		40.75

(a) (b)

Figure 4. s—lnr linear regression chart of stable flow pumping test for pressurized water in Huajiazhai. (a) Pumping stage; (b) Recovery phase

图4. 华家寨承压水稳定流抽水试验s—lnr线性回归图。(a) 抽水阶段；(b) 恢复阶段

3. 非稳定井流抽水试验水文地质参数解算方法

3.1. 非稳定井流Theis模型

经典的非稳定Theis井流模型，含水层为无界承压含水层，假定含水层为水平无界、均质各向同性、等厚、无补给，抽水量为含水层存储水量的弹性释放。

1935年Theis给出了无界承压含水层的非稳定流抽水的降深解[1] [2] [6]：

$s=\frac{Q}{4\pi T}\cdot w\left(u\right)$ (6)

$w\left(u\right)={\displaystyle {\int }_u^{\infty }\frac{e^{-x}}{x}}dx$ (7)

$u=\frac{r^2}{4at}=\frac{S\cdot r^2}{4T\cdot t}$ (8)

$T=k\cdot M$ (9)

$a=\frac{T}{S}$ (10)

式中：$s$ ——降深/回升(m)

$Q$ ——井出水量(m³/d)

$M$ ——含水层厚度(m)

$r$ ——井距(m)

$t$ ——抽水延续时间/停抽恢复时间(min)

$k$ ——渗透系数(m/d)

$a$ ——导压系数(m²/d)

$T$ ——导水系数(m²/d)

$S$ ——弹性释水系数

当抽水时间足够长达到准稳定状态后，其随后的恢复可近似当作与抽水无关的独立恢复试验，与抽水试验等同。恢复阶段，因无其它干扰，测试的试验数据质量好于抽水阶段。

3.2. Theis井函数计算

Theis井函数(7)为指数积分。

1stOpt最优化软件自带指数积分函数$\mathrm{Expint}()$ ，数值计算可以直接调用该函数。

Mathcad也自带指数积分函数$\mathrm{Ei}()$ ，但该函数仅用于符号运算，数值计算不能直接调用。

Theis井函数(7)为无限积分，被积函数在x = 0为奇点，高精度的数值计算有不少困难。

Theis井函数的数值计算，不少学者进行过一些研究[7] [8]，或精度较低，或拟合式较复杂且有适用条件，使用不太方便。

本文，对Theis井函数的积分变量进行变量代换[9]：

$\{\begin{array}{l}y=\ln \left(x\right)\\ x=e^y\end{array}\Rightarrow \{\begin{array}{l}{w}_1\left(u\right)={\displaystyle {\int }_u^{\infty }\frac{\exp \left(-x\right)}{x}dx}\\ w_2\left(u\right)={\displaystyle {\int }_{\ln (u)}^{\infty }\exp \left(-e^y\right)dy}\end{array}\Rightarrow \{\begin{array}{l}{f}_1\left(x\right)=\frac{\exp \left(-x\right)}{x}\\ f_2\left(y\right)=\exp \left(-e^y\right)\end{array}$

变量代换后，消除了被积函数奇点，极大地改善了被积函数的强峰性态及数值积分条件，见图5。

(a) f₁(x) (b) f₂(y)

Figure 5. Curve plot of Theis well function as an integrand (double logarithmic coordinate)

图5. Theis井函数被积函数曲线图(双对数坐标)

在此基础上，将井函数的无限积分截断为有限积分，积分上限取为：

$\{\begin{array}{l}{y}^{+}=4\\ x^{+}=e^4\end{array}$

则井函数为：$w\left(u\right)={\displaystyle {\int }_{\ln u}^4\exp \left(e^y\right)}dy$ (7-1)

井函数(7-1)的积分截断误差(积分截断余项)$\epsilon <{10}^{-15}$ ，计算精度足以满足实际需求。

3.3. 非稳定Theis井流水文地质参数解算原理

Theis井流模型，降深$s$ 与水文地质参数$\left(T,S\right)$ 为非线性关系。传统方法是采用logs—logt手工配线法求解参数。当时间较长后，可使用s—logt直线法求解。直线法，有一定的适用条件，不能利用前期观测数据；人工配线法效率低、人为误差大。

用s—t观测数据进行最小二乘法拟合，可使用全部观测数据，消除人为误差，求解精度高。

由式(6)，降深为水文地质参数$\left(T,S\right)$ 的函数：

$s\left(T,S,t\right)=\frac{Q}{4\pi T}\cdot w\left(T,S,t\right)$ (11)

基于最小二乘法的参数求解，即为使实测降深$sg$ 与理论降深$s$ 的最小二乘目标函数E最小：

$\min E\left(T,S\right)={\displaystyle \sum _{i=1}^n{\left(s{g}_i-s{}_i\right)}^2}$ (12)

求解最小目标函数为非线性最优化问题，是目前的研究热点。

Mathcad，是一款所见即所得的数学软件，计算代码与数学表达式十分接近，兼及代数分析及数值计算。最优化问题可使用最小误差函数Minerr()求解局部最优解。实际求解时，需对求解参数赋多组初值，比较解算结果，最终确定符合水文地质条件的最优解。

1stOpt，是一款国内七维高科有限公司开发的最优化分析专业软件，使用简单，功能强劲。该软件的最大特色是不需要对求解参数赋初值即可大概率地找到全局最优解，性能与国外同类产品相比毫不逊色[10]。

3.3.1. Mathcad局部最优化关键代码

使用Mathcad的最小误差函数Minerr()，求解最优水文地质参数的局部最优解，其关键代码为：

1) 输入试验数据

$Q:=22.6*24$

$r:=117.85$

$tg:={\left(8\begin{array}{c}\end{array}10\begin{array}{c}\end{array}12\begin{array}{c}\end{array}15\begin{array}{c}\end{array}20\begin{array}{c}\begin{array}{c}\end{array}25\begin{array}{c}\end{array}30\begin{array}{c}\end{array}40\begin{array}{c}\end{array}\end{array}50\begin{array}{c}\end{array}60\begin{array}{c}\end{array}\cdots \right)}^{\text{T}}$

$sg:={\left(0.002\begin{array}{c}\end{array}0.005\begin{array}{c}\end{array}0.006\begin{array}{c}\end{array}0.007\begin{array}{c}\end{array}0.008\begin{array}{c}\begin{array}{c}\end{array}0.011\begin{array}{c}\end{array}0.020\begin{array}{c}\end{array}0.029\begin{array}{c}\end{array}\end{array}0.038\begin{array}{c}\end{array}0.050\begin{array}{c}\end{array}\cdots \right)}^{\text{T}}$

2) 定义函数

$u\left(T,S,t\right):=\frac{1440\cdot S\cdot r^2}{4\cdot T\cdot t}$

$w\left(T,S,t\right):={\displaystyle {\int }_{\ln \left(u\left(T,S,t\right)\right)}^4\exp \left(-e^y\right)}dy$

$s\left(T,S,t\right):=\frac{Q}{4\cdot \pi \cdot T}\cdot w\left(T,S,t\right)$

3) 赋参数初值

$\left(\begin{array}{l}T\\ S\end{array}\right):=\left(\begin{array}{l}{T}_0\\ S_0\end{array}\right)$ 可根据情况给出不同的初值$\left(\begin{array}{l}{T}_0\\ S_0\end{array}\right)$

4) 定义目标函数

$E\left(T,S\right):={\displaystyle \sum _{i=0}^{m-1}{\left(s_i-s{g}_i\right)}^2}$

5) 求解局部最优解

$\begin{array}{l}\text{Given}\\ \begin{array}{cc}& \end{array}E\left(T,S\right):=0\end{array}$

$\left(\begin{array}{l}T\\ S\end{array}\right):=\mathrm{Minerr}\left(T,S\right)$

6) 计算拟合度

$RMS:=\frac{\sqrt{E\left(T,S\right)}}{m}$ 拟合度越小，拟合效果越好

3.3.2. 1stOpt全局最优化关键代码

1stOpt最优化代码也很简洁，也由6部分构成：

1) 试验常量赋值

Constant Q=22.6*24 ;

Constant r=[117.850] ;

2) 定义变量

Variable t,s;

3) 字符串常量赋字

ConstStr u=1440*r^2*SS/(4*TT*t) ;

ConstStr w=Expint(u) ;

4) 参数赋初值与约束(参数赋初值非必需)

Parameters TT=80.00[0,] ;

Parameters SS=0.001[0,1] ;

5) 定义拟合函数(Theis降深函数)

Function s=Q/(4*pi*TT)*w ;

6) 给出抽水试验数据(本文仅给出抽水阶段数据)t—s

data;

// t

8,10,12,15,20,25,30,40,50,60,

80,100,120,150,180,210,240,300,360,480,

600,720,900,1200,1500,1800,2100,2400,3000,3600,

4200,4800,5820;

// s

0.002,0.005,0.006,0.007,0.008,0.011,0.020,0.029,0.038,0.050,

0.093,0.130,0.165,0.220,0.270,0.330,0.370,0.465,0.530,0.655,

0.755,0.880,1.000,1.150,1.220,1.320,1.390,1.450,1.510,1.670,

1.710,1.720,1.730;

注：因1stOpt不区分大小写字母，代码中TT为导水系数系数T，SS为弹性释水系数S。

3.4. 非稳定Theis井流水文地质参数解算实例

Table 5. Water pumping test data of Theis well flow in unbounded confined aquifer in Feng County, Jiangsu Province

表5. 江苏丰县无界承压含水层Theis井流抽水试验数据

日期	时间	延续时间	恢复时间	抽水流量	观1		延续时间	恢复时间	观1
					降深水位	回升水位			降深水位	回升水位
y-m-d	h:min	t/min	t'/min	Q/m^3/d	s/m	s'/m	t/min	t'/min	s/m	s'/m
1976/11/6	9:40	0		542.4	0.000	开抽	5825	5	1.729	0.001
		8		542.4	0.002		5826	6	1.728	0.002
		10		542.4	0.005		5828	8	1.727	0.003
		12		542.4	0.006		5830	10	1.725	0.005
		15		542.4	0.007		5832	12	1.721	0.009
		20		542.4	0.008		5835	15	1.719	0.011
		25		542.4	0.011		5840	20	1.715	0.015
		30		542.4	0.020		5845	25	1.710	0.020
		40		542.4	0.029		5850	30	1.706	0.024
		50		542.4	0.038		5855	35	1.702	0.028
		60		542.4	0.050		5860	40	1.698	0.032
		80		542.4	0.093		5865	45	1.692	0.038
		100		542.4	0.130		5870	50	1.685	0.045
		120		542.4	0.165		5880	60	1.670	0.060
		150		542.4	0.220		5900	80	1.640	0.090
		180		542.4	0.270		5920	100	1.605	0.125
		210		542.4	0.330		5940	120	1.570	0.160
		240		542.4	0.370		5970	150	1.515	0.215
		300		542.4	0.465		6000	180	1.465	0.265
		360		542.4	0.530		6030	210	1.400	0.330
		480		542.4	0.655		6060	240	1.360	0.370
		600		542.4	0.755		6120	300	1.273	0.457
		720		542.4	0.880		6180	360	1.200	0.530
		900		542.4	1.000		6300	480	1.080	0.650
		1200		542.4	1.150		6420	600	0.980	0.750
		1500		542.4	1.220		6540	720	0.895	0.835
		1800		542.4	1.320		6720	900	0.794	0.936
		2100		542.4	1.390		7020	1200	0.650	1.080
		2400		542.4	1.450		7320	1500	0.560	1.170
		3000		542.4	1.510		7620	1800	0.452	1.278
		3600		542.4	1.670		7920	2100	0.390	1.340
		4200		542.4	1.710		8520	2700	0.260	1.470
		4800		542.4	1.720		10,040	4220	0.090	1.640
		5820		542.4	1.730	停抽， 恢复

【实例3】选自：文献[1]，p.136，例题，抽水，表4-2；p.146，例题，恢复，表4-5。

1976年在江苏丰县进行了承压含水层的非稳定流抽水试验。抽水试验场地位于江苏丰县王沟公社曹楼大队冲洪积平原，承压含水层为粉细砂及亚砂土(粉土)，抽水井井深120 m，直径0.20 m，自1976年11月6日9时40分开始抽水，抽水后进行了恢复观测。抽水持续97 h，流量为22.60 m³/h。观测孔1距抽水井117.85 m。抽水试验数据详见表5。

注：原文献未给出含水层厚度，无法求取渗透系数。

抽水阶段及恢复阶段的解算结果详见表6及图6，Mathcad与1stOpt的最优化结果非常接近。

Table 6. Calculation of pumping test parameters for Theis well flow in unbounded confined aquifer

表6. 无界承压含水层Theis井流抽水试验参数解算表

解算方法	最优化方法	赋参数初值		抽水试验			恢复试验
		T	S	T	S	RMS	T	S	RMS
		m^2/d		m^2/d		m	m^2/d
Mathcad	局部最优解	10	0.1	84.92	1.452E−03	0.006	82.12	1.565E−03	0.003
		100	0.01	84.92	1.452E−03	0.006	84.51	1.521E−03	0.002
		1000	0.001	133.70	5.250E−04	0.030	84.53	1.520E−03	0.002
		10000	0.0001	3089.00	0.000E+00	0.105	3096.00	0.000E+00	0.084
		最优解		84.92	1.452E−03	0.006	84.51	1.521E−03	0.002
1stOpt	全局最优解			84.92	1.452E−03	0.006	84.36	1.531E−03	0.002
原文人工配线法				80.16	1.540E−03	0.007	88.56

(a) (b)

Figure 6. Fitting curve diagram of unsteady flow pumping test for Huajiazhai confined aquifer. (a) Pumping stage; (b) Recovery phase

图6. 华家寨承压含水层非稳定流抽水试验拟合曲线图。(a) 抽水阶段；(b) 恢复阶段

4. 几点讨论

非稳定井流试验，一般在抽水试验后会进行恢复水位观测。抽水过程中，流量及水位通常受到一定程度的干扰，波动较大，数据观测误差较大。恢复过程中，无设备干扰，无流量，水位没有波动，观测数据质量较高。

非稳定流参数解算，通常是分别对抽水阶段及恢复阶段进行单独解算。

抽水阶段，由于观测数据有一定程度的干扰，解算参数有一定程度的误差。

恢复阶段，由于没有抽水设备的干扰，观测数据误差较小，解算参数的精度相对较高。

实际抽水试验中，通常会在达到“稳定状态”一定时间(如4 h或8 h)就停止抽水转入恢复水位观测。

理论上，无界承压含水层的Theis井流，不会出现“稳定状态”，但抽水一定时间后可能会出现“假稳定状态”，水位变化很小，认为达到“稳定状态”。

实际上，“假稳定状态”下，恢复水位仍然受到抽水阶段的影响。单独用恢复试验资料进行最优化解算，恢复的前期水位并不能很好地拟合，见图6(b)拟合曲线图。

抽水试验，叠加恢复试验基本是常态。多降次、阶梯流量非稳定流抽水试验实际中也经常采用。

抽水迭加恢复，以及多降次、阶梯流量等非稳定井流抽水试验，其全时程全局最优化解算，目前尚未得到完全解决。利用现代计算工具，对非稳定井流抽水试验进行全时程全局最优化反演解答，仍需进一步深入研究探索。

5. 结语

稳定流的两两孔解算，以及s—logr人工直线图解法均有一定的解算误差。非稳定流logt—logs手工配线法，有较大的随意性，参数解算误差较大。

目前，抽水试验水文地质参数的反演有了新的工具，可利用现代计算工具进行最优化解算。

潜水及承压水稳定井流，参数反演求解比较简单。利用Excel回归方法解算的水文地质参数为整体最小误差解，可最大限度地消除由两两孔解算的参数误差。

Theis井函数计算，利用变量代换消除了被积函数奇点；将无限积分截断为有限积分，积分截断误差 < 10⁻¹⁵。

Theis非稳定井流，参数反演求解复杂。传统的s—logt直线图解法，只能利用部分观测数据；人工配线法，解算具有一定的随意性，误差较大。

使用Mathcad数学软件进行局部最优化解算，利用1stOpt最优化专业软件进行全局最优化解算，两者最优化解算结果非常接近。

潜水与承压水稳定井流以及Theis非稳定井流，均选择了教科书经典实例予以解算验证，参数解算结果客观、优化。

抽水迭加恢复，以及多降次、阶梯流量等非稳定井流抽水试验，其全时程全局最优化解算，目前尚未得到完全解决。利用现代计算工具，对非稳定井流抽水试验进行全时程全局最优化反演解答，仍需进一步深入研究探索。

NOTES

^*第一作者。

^#通讯作者。