各向异性含水层参数估计的差分进化算法

王晓晨，常安定，王静云，张 转(长安大学 理学院，陕西 西安 710064)

各向异性含水层参数估计的差分进化算法

王晓晨，常安定，王静云，张 转(长安大学 理学院，陕西 西安 710064)

【目的】 通过分析抽水试验数据进行各向异性含水层参数估计，为各向异性含水层参数估计提供新的方法。【方法】 以各向同性情况下的泰斯公式为基础,将差分进化算法应用于各向异性条件下含水层参数的估计，并将其计算结果与其他方法的参数估计结果进行比较，对不同初始参数范围下的参数估计值进行分析探讨。【结果】 差分进化算法估计的参数值与其他方法估计的参数值之间有较小的差异，相对误差不到4%，且目标函数值相对更小，达到0.003 93；对于不同的初始参数范围，利用该算法均能得到满意的参数估计结果。【结论】 基于抽水试验数据估计各向异性含水层参数的差分进化算法估计结果有效且可靠，算法稳定性高，寻优能力强。

各向异性；含水层参数；差分进化算法；参数估计；抽水试验

在研究与地下水运动有关的工程问题中，含水层的导水性能具有各向异性的特征，且含水层导水系数和储水系数作为最重要的参数，其准确有效地求解极为重要。为此许多学者在利用抽水试验确定各向异性含水层参数问题上做了很多的研究工作。Hantush[1]给出了有越流补给条件下，径向各向异性含水层中非完整井流三维问题的数学模型及相应的计算公式。Papadopulos[2]推导出了径向各向异性主值Kr与2个水平方向主渗透系数Kx、Ky大小和方向的关系式。之后，Way等[3]基于Hantush[1]和Papadopulos[2]的研究结果，提出了通过现场抽水试验确定径向各向异性渗透系数主值(Kr、Kz和Kx、Ky)的方法。Neuman等[4]利用张量分析的方法推导出了平面各向异性渗透系数主值和主方向的解析公式。Chen等[5]讨论了Hantush解中越流的影响及越流系数的计算方法。周志芳等[6]通过抽水试验，应用图解法和优化算法结合的半解析方法确定平面各向异性渗透系数主值和主方向，之后又推导出了有越流情况下各向异性含水层中地下水三维井流问题的计算公式[7]。刘燕等[8]以描述均质各向同性与各向异性情况下的非稳定流问题的解析解为基础，在满足直线图解法适用条件的情况下，分析了3组实际抽水试验数据，计算了各向同性与各向异性2种情况下的含水层参数。近年来，智能优化算法在经济和工程等各个领域有着广泛的应用。康瑞龙[9]于2013年成功地应用智能优化算法——人工根系算法求解了各向异性含水层参数问题。区别于传统的估计方法，智能优化算法求解速度快而且计算精度高。差分进化算法(Differential evolution,DE)是一种随机并行直接搜索的智能优化算法[10]，该算法容易实现，对参数的依赖程度低。鉴于此，本研究根据抽水试验数据，尝试将差分进化算法应用于确定各向异性条件下含水层参数的计算问题，并对算法和计算结果进行分析，以期为各向异性条件下含水层参数的准确估计提供方法上的支持。

1 各向异性条件下的含水层参数问题

1.1 抽水试验基本知识

在抽水试验中，假设在无限延伸含水层中，任一点处的水位降深随时间t变化用各向异性表达式表示：

(1)

式中：s为水位降深，m；Q为抽水井中的抽水流量，m3/s；Te可以定义为含水层的等效导水系数，m2/s；W(uxy)为与泰斯公式[11]形式相同的井函数。Te的表达式为：

(2)

式中：Txx和Txy分别为导水系数在当地坐标系x和y方向的张量分量，m2/s；函数W(uxy)中无量纲时间为：

(3)

式中：μ为含水层储水系数；Tyy为导水系数在当地坐标系下的张量分量，m2/s；x，y为当地坐标系的坐标分量，m。

在此，采用Srivastava[11]的泰斯井函数的近似表达式进行计算，具体表达式如下。W(u)=-lnu+a0+a1u+a2u2+a3u3+a4u4+a5u5，u≤1。

(4)

(5)

其中，a0=-0.577 22，a1=0.999 99，a2=-0.249 91，a3=0.055 19，a4=-0.009 76，a5=0.001 08；b0=0.267 77，b1=8.634 76，b2=18.059 02，b3=8.573 33；c0=3.958 50，c1=21.099 65，c2=25.632 96，c3=9.573 32。Srivastava[11]指出，以上的近似表达式与泰斯公式精确解相比，相对误差均小于 0.001%。

(6)

(7)

式中：X和Y为导水系数张量的全局坐标。

1.2 目标函数的建立

通常情况下，采用3个不在同一观测线上的观测孔中的降深-时间数据确定当地坐标系下的导水系数张量分量Txx、Tyy和Txy的值，然后再利用式(6)和(7)计算全局坐标系下的导水系数张量分量TX和TY值。因此，对于各向异性含水层参数的估计问题，可以将含水层储水系数μ及导水系数在当地坐标系下的张量分量Txx、Tyy和Txy设定为待估参数。

对于第i(i=1，2，3)个观测孔，其tj时刻的水位降深sj可由式(1)求得。为了使所求得的sj与真实水位降深相接近，所以构造目标函数：

(8)

对于各向异性条件下的含水层，u、Txx、Tyy、Txy为待估参数，其向量为θ=(u,Txx,Tyy,Txy)。为此，使得目标函数最小的参数值θ=(μ,Txx,Tyy,Txy)即为所求。

2 差分进化算法

2.1 DE算法的思想

DE算法是由Storn等[12]于1995年提出的,它是一种随机并行直接搜索算法，该算法的基本思想是应用当前种群个体的差异重组得到中间种群，然后应用子代个体与父代个体竞争获得新一代种群[13]。算法中采用由变异到交叉再到选择的操作顺序更新下一代种群。其中，变异操作是差异演化的关键步骤，交叉操作可以增加种群的多样性，选择决定着试验向量是否会成为下一代中的成员，确保下一代中的所有个体都比当前种群的对应个体更佳或者至少某一方面要好[14]。

2.2 DE算法的步骤及流程

步骤1：随机初始化种群，设置相关参数。

步骤2：计算各个体最优适应度fi，及其相应的位置gi，从中找出当代全局最优适应度f及其相应的位置g。

步骤4：若满足终止条件，则迭代停止；否则，k=k+1，转步骤3。

步骤5：输出全局最优适应度f、相应的位置g和迭代次数t。

差分进化算法流程图如图1所示。

3 数值试验

3.1 试验数据与试验条件

3.1.1 试验数据 试验数据引自文献[15]。假定在全局某正交各向异性含水层中，有一完整井以定流量Q=0.012 57 m3/s进行非稳定流抽水试验，有3个观测孔OW1、OW2和OW3，其在坐标系中的位置分别为(28.3，0)、(9.0，33.5)、(-19.3，-5.2)。当抽水井以定流量非稳定流抽水后，在3个观测孔中观测到的地下水降深随时间变化的试验数据如表1所示。

图1 差分进化算法流程图Fig.1 Flow diagram of DE

3.1.2 试验条件 初始种群规模为50；最大迭代次数2 000；达到精度eps=0.004；参数限定范围参考文献[9]：0≤μ≤1，0≤Txx≤600，0≤Tyy≤600，-600≤Txy≤0。对于差分进化算法中的参数CR的选取与文献[16]不同，其选取原则是使差分进化算法在解决各向异性条件下含水层参数问题具有较好的全局搜索能力，在此选取CR=0.8。

将差分进化算法程序反复执行100次，如果计算的目标函数值小于精度0.004，认为算法成功，即利用差分进化算法成功地估计了各向异性条件下的含水层参数；否则，认为算法不成功。记录算法每次的运行时间和算法成功的次数。

对于上述试验数据和试验条件，采用Matlab语言，根据算法流程实现DE算法。

表1 观测孔时间与降深数据Table 1 The time-drawdown data of the observation

3.2 试验结果与分析

3.2.1 试验结果 应用差分进化算法估计各向异性条件下的含水层参数，将其计算结果(迭代第80次)与相关文献中的结果进行对比，结果如表2所示。

表2 基于差分进化算法的各向异性含水层参数估计结果(迭代第80次)与相关文献结果的比较Table 2 The results (Running the 80th time) of the anisotropy aquifer parameters based on DE compare with the results of the related literature

3.2.2 算法有效性 由表2可以看出，在将差分进化算法应用于估计各向异性条件下含水层参数时，程序迭代到80次时达到精度要求，其目标函数值为0.003 93，计算出的相应参数与文献[15]中的参数相差较小，且计算出的相应参数间的相对误差均小于4%。因此，差分进化算法不论对待估参数μ、Txx、Tyy、Txy，还是通过式(6)、式(7)计算的TX和TY，均能得到可靠的计算结果。

图2是差分进化算法程序执行100次过程中，目标函数值的分布情况。由图2可以看出，在差分进化算法程序执行100次中，有93次成功，其寻优率达到93%，平均运行时间为1.928 35 s。由此可知，差分进化算法的寻优效率较高且运行时间较短，具有高效性。因此，运用差分进化算法能有效地估计各向异性条件下的含水层参数。

图2 差分进化算法执行100次过程中目标函数值的分布Fig.2 Distribution of the objective function value when DE is executed 100 times

3.2.3 算法稳定性 表3的试验结果是在算法参数设定不变，参数限定范围为文献[15]的不同倍数的情况下，应用差分进化算法得到的参数估计值。从表3中的时间数据可以看出，随着初始状态参数范围的不断增大，差分进化算法均能较好地收敛，得到较为满意的计算结果，即差分进化算法对初值选取的敏感性较低，寻优能力较强，具有较好的稳定性。在初始状态参数范围倍数增大的情况下，相应的所需运行时间有所增加，但是时间波动不大，其增加量不超过0.3 s，且运行时间较短，这说明差分进化算法能以较短的时间收敛于全局值。因此，差分进化算法解决各向异性条件下含水层参数问题是稳定的。

表3 参数初始范围为文献[15]参数值的不同倍数时DE算法对各向异性含水层参数的计算结果Table 3 The calculation results of DE under the initial range of parameters are setted to multiple the parameters of the literature [15]

3.2.4 算法收敛性 在不限制算法精度的条件下，让算法以迭代2 000次为终止条件，可以得到算法在执行189次时会收敛到0.003 402 81(图3)，由此可知差分进化算法有较快的收敛速度。

图3 差分进化算法执行2 000次时的收敛曲线Fig.3 Convergent curve of DE when DE is executed 2 000 times

4 结束语

鉴于各向异性含水层参数估计的复杂性，本文尝试应用差分进化算法对各参数进行估计。根据试验结果可以知道，用差分进化算法能够成功确定各向异性条件下的含水层参数，并且可以得到可靠的优化计算结果。本试验结果表明：(1)差分进化算法估计各向异性含水层参数时，收敛速度快，运行时间短，精度高而且具有较好的寻优率；(2)算法对不同的初值范围均能得到满意的参数计算结果，具有较高的稳定性和较强的寻优能力，且对初值的敏感程度较低。因此，差分进化算法可以成功、有效地用于求解各向异性含水层的参数问题。

[1] Hantush M S.A method for analyzing a drawdown test in anisotropic aquifers [J].Water Resources,1966,2(2):281-285.

[2] Papadopulos I S.Nonsteady flow to a well in an infinite anisotropic aquifer [C]//United Nations Educational Scientific and Cultural Organization.Symposium on Hydrology of Fractured Rocks,Paris:U.S.Geological Survey,1965:21-31.

[3] Way S C,Mckee C R.In-situ determination of three-dimensional aquifer permeabilities [J].Ground Water,1982,20(5):594-603.

[4] Neuman S P,Walter G R,Bentley W,et al.Determination of horizontal anisotropy with three wells [J].Ground Water,1984,22(1):66-72.

[5] Chen X,Ayers J F.Utilization of the Hantush solution for the simultaneous determination of aquifer paramenter [J].Ground Water,1997,35(5):751-756.

[6] 周志芳,朱学愚,李 艳.岩体渗透系数张量的半解析计算 [J].水利学报,1997(9):06-11. Zhou Z F,Zhu X Y,Li Y.Semi-analytic method for conductivity tensor in rock mass [J].Journal of Hydraulic Engineering,1997(9):6-11.(in Chinese)

[7] 周志芳.任意各向异性岩体渗透系数张量的半解析计算 [J].水利学报,1999(3):65-70. Zhou Z F.Semi-analytic method for conductivity tensor in random anisotropic rock mass [J].Journal of Hydraulic Engineering,1999(3):65-70.(in Chinese)

[8] 刘 燕,辛璐军,郭建青，等.抽水试验确定各向异性含水层参数的实例讨论 [J].勘察科技技术,2012(6):5-9. Liu Y,Xin L J,Guo J Q,et al.Discussion of determining anisotropic aquifer parameters by pumping tests data [J].Site Investigation Science and Technology,2012(6):5-9.(in Chinese)

[9] 康瑞龙.一种新型的智能优化算法：人工根系算法 [D].西安:长安大学,2013. Kang R L.A new intelligent optimization method:Artificial root optmization [D].Xi’an:Chang’an University,2013.(in Chinese)

[10] Stroll R,Price K.Differentifial evolution:A survey of the state of the art [J].IEEE Transaction on Evolutionary Computation,2011,15(1):4-31.

[11] Srivastava R.Implications of using approximation expressions for well function [J].Journal of Irrigation and Drainage Engineering,1995,121(6):459-462.

[12] Storn R,Price K.Differential evolution:A simple and efficient adaptive scheme for global optimization over continuous spaces [M].Berkeley:ICSI,1995:22-25.

[13] 万 东.差分进化算法研究及其应用 [J].科学技术与工程,2009,9(22):6672-6676. Wan D.Research and application based on differential evolution algorithm [J].Science Technology and Engineering,2009,9(22):6672-6676.(in Chinese)

[14] Storn R,Price K.Differential evolution:A simple and efficient heuristic for global optimization over continuous spaces [J].Journal of Global Optimization,1997,11(4):341-359.

[15] Kruseman G P,Ridder N A,Verweij J M.Analysis and evaluation of pumping test data [M].2ndEdition.Wageningen Netherlands:International Institute for Land Reclamation and Improvement,1994:133-145.

[16] 周艳平,顾幸生.差分进化算法研究进展 [J].化工自动化及仪表,2007,34(3):1-5. Zhou Y P,Gu X S.Development of differential evolution algorithm [J].Chemical Industry Automation and Instrumentation,2007,34(3):1-5.(in Chinese)

Differential evolution method for estimating the anisotropy aquifer parameters

WANG Xiao-chen,CHANG An-ding,WANG Jing-yun,ZHANG Zhuan

(CollegeofScience,ChanganUniversity,Xi’an,Shaanxi710064,China)

【Objective】 This study estimated the anisotropy aquifer parameters by the analysis of pumping test data. In order to provide a new method for estimating the parameter under the condition of anisotropy aquifer.【Method】 On the basis of the theis formula in the isotropic condition,differential evolution algorithm was applied to estimate the parameter under the condition of anisotropy aquifer.On the one hand,the calculation compare with other methods of parameter estimation and analyzes the results.On the other hand,discussing the estimated value of parameter that calculate under different initial scope.【Result】 The relative error that the values of the parameters are estimated between differential evolution algorithm and other method and the objective function value are relatively smaller.The relative error less than 4% and the objective function value is 0.003 93;For different initial parameters,differential evolution algorithm can get a satisfactory parameter estimation result.【Conclusion】 Based on pumping test data,the results of differential evolution method are more effective and reliable.It has higher stability and the stronger ability of optimization than other methods.

random anisotropic;aquifer parameters;differential evolution;parameter estimate;pumping test

2013-12-09

中央高校基本科研业务费专项资金(310829130225)

王晓晨(1988-),女,辽宁鞍山人,在读硕士,主要从事最优化理论与方法研究。E-mail:wangxiaochenbest@126.com

常安定(1964-),男,陕西大荔人,教授,硕士生导师,主要从事水文地质的数学分析方法研究。E-mail:chdanding@126.com

时间:2015-03-12 14:17

10.13207/j.cnki.jnwafu.2015.04.023网络出版地址:http://www.cnki.net/kcms/detail/61.1390.S.20150312.1417.023.html

1671-9387(2015)04-0223-06

TV211.1+2

A