随着人们对地下水污染问题的高度重视，地下水污染修复治理工作成为地下水研究领域的热点^[1-3]。循环井技术作为一种原位修复技术，具有对含水层结构破坏小，易与其他技术结合使用，结构简单等优点，已在国外广泛应用^[4]。循环井技术通过抽注地下水，在井周围形成三维地下水环流场，增加地下水流动性以达到更佳的修复效果^[5]。循环井的运行效果受水文地质条件、循环井结构和运行条件等多种因素的综合影响^[6]。而循环井群运行时还存在井间干扰，这种干扰还受相邻循环井循环模式、井间距等因素影响^[7]。因此循环井间距及相邻循环井的循环方式也是循环井群设计中需要考虑的关键参数。因此，循环井技术的高效应用的关键是根据场地特定的水文地质条件对上述关键参数开展设计。

目前国内外多采用数值模拟的方法来评价循环井群在不同设计方案下的运行效果，且多以循环效率、污染物去除效率等作为衡量循环井修复效果的关键指标^[8]。如有学者利用有限差分法模拟循环井运行时地下水流特征，采用粒子回收率、影响半径等指标来评估循环效率的可靠性，分析各影响因素对循环效率的影响规律，并应用于循环井的初步设计^[8-10]。HUANG等^[11]利用粒子追踪方法，分析了相反抽注模式下循环井结构以及运行参数对双循环井运行效果的影响。许多学者根据修复过程中循环井对污染物的捕获率及对修复试剂的去除效率等开展试验研究^[12-13]。JOHN等^[14]分析了地下水的各向异性比对循环井群运行时污染物三氯乙烯 (TCE) 去除效率的影响。虽然国内外以污染物的去除对循环井运行效果开展了大量研究，受限于实际场地修复成本和环境等条件，如何在保证循环效果的同时使场地污染物去除效率最高，是应用双循环井优化设计时重点考虑的问题^[6,8]。上述研究均是采用数值模拟方法来预测循环井在不同因素影响下运行效果的变化特征，但是关于多种因素综合影响下最优方案设计的研究仍不成熟。因此将数值模拟方法和多目标优化方法相结合，构建模拟-优化耦合模型，可以在满足各种环境、技术等要求的前提下获得循环井优化设计方案。因此该技术作为一种解决优化设计方案问题的有效手段。

模拟-优化耦合模型需要将模拟模型作为目标函数编入优化模型，而在利用优化模型求解全局最优解时，通常需要多次调用数值模型，计算量过大，限制了模拟-优化耦合技术的使用^[15]。近年来，随着深度学习理念逐渐渗入到各研究领域，其优越性逐渐被人们所认可。深度学习方法在处理地下水相关研究中的高度非线性复杂问题时表现良好^[16]。LÄHIVAARA等^[17]使用卷积神经网络根据地震波数据获取地下水储量。HONG等^[18]提出一种3D-CNN模型，用于快速估计含水层的各向异性比，并将所设计的模型用以预测其他岩石的渗透率，预测结果的准确率高达0.91。通过深度学习模拟模型计算得到双循环井的输入-输出特征响应关系，建立替代模型，在优化模型中进行反复调用，节省模拟优化耦合模型的计算时间，减小计算负荷。由于替代模型求解容易，结构相对简单，大大减小了计算量，使模拟-优化耦合技术在地下水双循环井优化工作成为可能。

本研究基于深度学习方法形成模拟-优化耦合技术，以双循环井作为井群运行的典型单元，将双循环井运行效果的关键影响参数作为输入变量，利用FloPy调用MODFLOW、MODPATH和MT3DMS模块构建双循环井驱动下的地下水流、粒子追踪及溶质运移数值模型，将表征运行效果指标循环效率 (Pr) 和污染物去除效率 (η) 作为模型输出参数，构建数据集；利用卷积神经网络模型建立替代模型，利用NSGA-Ⅱ算法建立以双循环井运行效果的多目标优化模型；在此基础上对陕西省西安市某试验场地开展双循环井优化设计，以期为场地循环井群的结构及运行参数的确定提供理论支撑。

1.   方法

本研究的研究思路主要分为3个基本步骤，具体流程如图1所示。

图 1  优化过程流程图

Figure 1.  Flow chart of the optimization process

下载: 全尺寸图片 幻灯片

1.1   数据集的构建

1) 双循环井简介。双循环井具有多种结构形式，本研究以2个典型循环井组成的循环井对为研究单元。循环井包含2个由阻隔器分隔开的筛管 (抽水筛和注水筛) ，受污染的地下水通过抽水筛抽出，经处理后再通过注水筛以相同流量注入含水层。地下水在循环井周围形成环流场，而相邻两井的循环模式不同，引起的水流流向不同，且相邻井之间产生互相扰动。主要有两种不同的情况，即循环方式相同 (如图2(a)) 和循环方式相反 (如图2(b)) 。

图 2  双循环井水流流动

Figure 2.  Flow of groundwater around dual-circulation wells

下载: 全尺寸图片 幻灯片

2) 数据集构成。利用FloPy软件包和Python开发双循环井数值模型，利用该数值模型建立包含输入变量和模型输出参数的数据集。如图1所示，输入变量主要包括双循环井运行效果的影响因素 (水文地质参数、循环井结构参数以及运行参数^[19-20]) ，输出变量为模型输出双循环井运行180 d后^[10,21]的运行效果表征参数 (循环井群循环效率、污染物去除效率^{[10,19,21-22]}) 。其中循环井抽注方式 (相邻井抽注模式相同或相反) 的确定在本次研究中由循环井结构参数中抽注水筛的位置关系决定。循环井分为上抽下注和上注下抽2种抽注模式，即循环井抽水筛顶部与含水层顶部的距离小于注水筛顶部与含水层顶部的距离，则该循环井抽注模式为上抽下注模式，反之则称为上注下抽模式。

3) 数值模型求解。将输入变量输入到数值模型中，利用地下水流模拟程序 (MODFLOW) 、粒子追踪模拟程序 (MODPATH) 和溶质运移模拟程序 (MT3DMS) 分别进行双循环井地下水流、运动轨迹及污染物运移的数值模拟计算，根据MODPATH建立粒子追踪模型，以所有循环井的循环粒子数与投放总粒子数的比作为循环井群循环效率；此外根据MT3DMS模拟污染物在地下水中的运行过程，以运行时段末与运行前场地污染物的变化总量与初始污染物总量的比作为污染物去除效率^[23-25]。最终形成输入变量和模型输出参数共同组成的数据集。

1.2   基于卷积神经网络方法的替代模型构建

卷积神经网络 (convolutional neural networks，CNN) 是一种具有多个处理层的神经网络模型，主要由卷积层、池化层和全连接层这3类层组成，通过将这3类层按照不同数量及顺序叠加起来，即可构建出一个完整的卷积神经网络^[26]。本研究根据已获得的双循环井数据集，利用卷积神经网络模型学习双循环井影响因素与运行效果响应之间的映射关系。卷积神经网络模型的训练基于损失函数的计算结果，将误差进行反向传播到卷积神经网络的各个层，各层根据误差的反馈进行参数权值的更新。之后随着模型的训练、验证过程的进行，将卷积神经网络的结构及其他参数进行不断优化，从而使基于深度学习方法的替代模型与数值模拟模型达到最佳拟合。

根据前述样本空间的参数选取以及模型求解得到的双循环井性能指标，本研究共利用已有的800个数据集样本，随机抽取600组为训练集、100组为验证集和100组为测试集。其中训练集在卷积神经网络的训练过程中，通过梯度下降的方式减小训练误差，并得到模型内置参数，使模型获得预测双循环井运行效果表征参数响应的能力；验证集用于调整模型的网络结构和超参数，并对模型的能力进行初步评估；测试集用于最后模型精确度的评估。

在模型的构建上，本研究采用自适应动量优化 (adaptive moment，Adam) 算法优化卷积神经网络建立的替代模型，使用均方误差 (mean squared error，MSE) 作为损失函数，为了确定卷积神经网络模型中的网络结构和超参数，尝试了多种组合，并最终确定了误差最小的卷积神经网络模型结构，实验结果为建立的模型结构包括交替连接的3个卷积层和3个最大池化层、1个展平层和3个全连接层。另外该模型的核心参数配置：每批样本为50，迭代次数为3 000，激活函数均为线性整流函数 (rectified linear unit，ReLu) ，卷积神经网络内部参数如图3所示。

图 3  基于卷积神经网络的循环井运行效果预测模型

Figure 3.  The predictive model for GCW operation effect based on CNN

下载: 全尺寸图片 幻灯片

1.3   优化模型的构建

1) 目标函数。在本研究中，通过循环效率和污染物去除效率来衡量双循环井的运行效果，循环效率 (Pr) 是衡量双循环井抽水筛对地下水捕获能力的指标，可以体现水力环流的闭合性。较高的循环效率不仅可以保证修复试剂的充分利用，也可以保证污染组分再双循环井系统中得到降解，从而避免污染物从双循环井的环流场逃逸，造成二次污染。可以利用MODPATH计算，从而通过数值模型量化该指标。另外双循环群修复的最终目的是高效捕获和去除地下水中的污染物，所以将污染物去除效率 (η) 也作为指标。因此目标函数的数学表达式可以表示为式(1)和式(2)。

式中：${N}_{inj}$为注入的粒子数；${N}_{rec}$为随时间t变化的粒子回收数；$Conc$表示场地初始污染物浓度；$Conc'$表示双循环井运行180 d后污染物浓度。

2) 约束条件。循环井群主要受到以下3种类型的约束条件的影响，即循环井所在位置的约束条件，循环井筛管的约束条件和流量的约束条件。

①循环井位置的约束条件，本研究主要考虑双循环井的设计方案，假定以第一个井的位置为原点 (0，0) ，第二个井的位置则为 (0，d) ，则井间距的约束条件见式(3)。

式中：${d}_{{\mathrm{max}}}$根据待修复场地的具体大小得到的可布置的两个循环井之间的最大井间距，m，该参数由场地大小决定。

②筛管位置的约束条件，保证每个循环井的上下筛管位置均位于含水层中，是使用循环井技术的必要条件^[4]，则筛管位置的约束条件见式(4)。

式中：${d}_{{\mathrm{in}}}^{i}$表示第i口井注水筛顶部与含水层顶板的距离，m；${d}_{{\mathrm{out}}}^{i}$表示第i口井抽水筛顶部与含水层顶板的距离，m；${L}_{{\mathrm{max}}}$为筛管的最大长度，m；M表示含水层的厚度，m。

③循环井流量的约束条件，假设循环井每个井的抽注流量大小相等，则循环井流量的约束条件见式(5)。

式中：${Q}^{i}$表示为第i口循环井抽注流量，m³∙d⁻¹；${Q}_{{\mathrm{in}}}^{i}$表示第i口井的注水筛流量，m³∙d⁻¹；${Q}_{{\mathrm{out}}}^{i}$表示第i口井的抽水筛流量，m³∙d⁻¹；${Q}_{{\mathrm{max}}}$表示为抽注流量的最大值，该参数由水文地质参数决定^[27]，各参数的含义详见图2。

3) 基于NSGA-Ⅱ算法的双循环井参数多目标优化。双循环井优化问题实际上是一个具有多个变量以及多目标函数的复杂优化问题。目前用于解决多目标优化问题的方法很多，比如遗传算法，NSGA算法，粒子群算法等。NSGA-Ⅱ算法比其他优化算法适用性更加广泛，在复杂应用场合上效果更好^[28]。NSGA-Ⅱ算法的基本流程如下。

①首先随机初始化父代种群P₀，种群规模为P，对父代种群进行非支配排序后，经过选择、交叉和变异遗传算子生成第一个子代种群N₀。

②从第二个子代种群开始，将初始化父代种群P₀和第一个子代种群N₀合并为新种群K₀，对K₀种群进行非支配排序后，对每个非支配层的个体计算拥挤度，并根据非支配关系和拥挤度选取N个个体组成父代种群P₁，再采用遗传算子产生新一代种群N₁。

③以此类推，直到进化代数gen达到设定的最大遗传代数gen_max时停止计算，得到Pareto最优解集。

为了保证NSGA-Ⅱ算法的全局寻优能力，准确且快速的找到最优解集，模型参数设置如下：初始化种群个数100；迭代次数300；选择模拟二进制交叉算子 (SBX) ；交叉概率0.8；使用多项式变异 (PM) ；变异概率0.2。通过求解优化模型，即可获得约束条件下的双循环井污染物去除效率和循环效率。

2.   结果与讨论

2.1   数据集分析

本研究中通过数值模型建立了包含输入变量和模型输出参数共同组成的有效数据样本800个，数据集的主要指标见表1。该数据集是用于训练替代模型的训练样本，训练样本覆盖越广泛，样本数越多，越能保证训练效果，使替代模型能够更好的拟合数值模型中输入输出的响应关系。由图4可知，该数据集中各输入变量的分布基本覆盖了循环井适用的全部范围，表明通过数值模型获取的数据集具有代表性。

表 1  数值模型数据集样本参数一览表

Table 1.  List of parameters of numerical model datasets

数据集	类型	参数	符号	单位	循环井适用范围
输入变量	水文地质参数	含水层厚度	M	m	1.5~35
		水平渗透系数	KH	m·d−1	1~20
		水平渗透系数/垂向渗透系数	KH/KV	—	3~10
		给水度	μ	—	0.01~0.2
		弹性释水系数	ss	—	0.000 01~0.000 1
		孔隙度	n	—	0.01~0.5
	循环井结构参数	第i口井的注水筛顶部与含水层顶板的距离	diin	m	0~M
		第i口井的抽水筛顶部与含水层顶板的距离	diout	m	0~M
		第i口井的注水筛长度	Liin	m	0~M
		第i口井的抽水筛长度	Liout	m	0~M
	循环井运行参数	抽注水量	Q	m3·d−1	0~200
		井间距	d	m	场地大小
		循环井抽注方式	—	—	—
模型输出参数	表征指标	循环效率	Pr	%	0~100
		污染物去除效率	$\eta$	%	0~100

 | Show Table

DownLoad: CSV

图 4  数据集中输入变量数值分布

Figure 4.  Distribution of main parameters in the database

下载: 全尺寸图片 幻灯片

2.2   替代模型的构建及精确度评价

评估训练后替代模型的精确度，选取100组数据作为测试集，对训练后的替代模型进行测试，以评估模型泛化能力。此外，为证明本研究中卷积神经网络替代模型的性能，选取深度学习模型中的反向传播 (back propagation，BP) 神经网络模型和循环神经网络模型 (recurrent neural networks，RNN) 进行对比。

利用Tensorflow平台建立卷积神经网络、BP神经网络和循环神经网络3种深度学习模型，以决定系数 (r-square，R²) 、平均绝对百分比误差 (mean absolute percentage error，MAPE) 作为评价替代模型的精确度的指标，R²越接近于1，表示替代模型的精确度越高；MAPE越接近于0，表示替代模型的精确度越高。其计算式如式(6)~式(7)所示。

式中：n表示样本个数，i表示样本编号，${\hat{y}}_{i}$表示替代模型输出值，${y}_{i}$表示模拟模型输出值，${\stackrel{-}{y}}_{i}$表示模拟模型输出值的平均值。各模型的预测结果具体如图5、6所示，结果表明CNN模型效果预测效果越好。

图 5  替代模型与模拟模型预测循环效率结果对比

Figure 5.  Comparison of cycle efficiency predicted by alternative model and simulated model

下载: 全尺寸图片 幻灯片

图 6  替代模型与模拟模型计算预测污染物去除效率结果对比

Figure 6.  Comparison of predicted pollutant removal efficiency results form alternative and simulated model

下载: 全尺寸图片 幻灯片

几种预测模型对双循环井运行效果表征参数的精确度R²和MAPE如表2所示，可知3种预测模型对循环效率的预测效果均比污染物去除效率要好；另外，无论是预测循环效率还是污染物去除效率， CNN模型的精确度都要优于另外两种深度学习模型。CNN模型对于循环效率的精确度R²达到0.972，MAPE达到0.008；而对于污染物去除效率的精确度R²达到0.946，MAPE达到0.009。

表 2  替代模型精确度评价指标统计

Table 2.  Statistics of calculation accuracy of alternative model evaluation index

替代模型	循环效率			污染物去除效率
	R2	MAPE		R2	MAPE
CNN模型	0.972	0.008		0.946	0.009
BP模型	0.823	0.028		0.789	0.138
RNN模型	0.737	0.030		0.299	0.221

 | Show Table

DownLoad: CSV

综上所述，基于卷积神经网络的替代模型和数值模型具有较高的拟合精度并且具有强大的泛化能力，可以更好地反映数值模拟模型的输入输出变量响应关系，并且可以嵌入优化模型中，保证了双循环井优化方案的可靠性。

2.3   试验场地双循环井优化结果与分析

2.3.1   场地概况

本研究选取陕西省西安市某示范场地为研究区，针对场地实际水文地质条件开展双循环井优化设计。该场地位于关中盆地渭河河谷阶地中部，场地东西长180 m，南北宽60~70 m。试验目标含水层为第四系全新统、上更新统冲积孔隙潜水含水层。含水层以中砂为主，赋水性良好，含水层厚度 (M) 为13.5 m、水平渗透系数 (K_H) 为20 m∙d⁻¹、水平渗透系数/垂向渗透系数 (K_H/K_V) 为3、给水度 (μ) 为0.2、弹性释水系数 (ss) 为0.0001、孔隙度 (n) 为0.39。将该试验场地的水文地质参数输入到测试后的双循环井多目标优化设计模型中，以获取适合该场地的优化设计方案。

2.3.2   优化结果分析

将试验场地的水文地质参数输入双循环井运行效果预测模型中，并利用NSGA-Ⅱ算法求解适合该试验场地的双循环井结构及运行参数的最优组合，从而实现双循环井的优化设计。NSGA-Ⅱ优化算法求解结果可采用超体积指数指标 (Hypervolume) 指标来评价。如图7所示，在迭代70次后，优化算法收敛至较优水平，并逐渐寻找最优解，证明该算法可行。

图 7  Hypervolume随迭代次数变化

Figure 7.  Variation of Hypervolume by iteration times

下载: 全尺寸图片 幻灯片

在获得的Pareto最优解集中选择4个典型方案进行对比，优化方案如表3~表4所示。其中方案1和方案2的2个循环井抽注模式相同 (均为上抽下注模式) 。方案1循环井的循环效率低于方案2，而污染物去除效率高于方案2；2种方案的抽注水量相同，但方案2的井间距小于方案1。井间距越小，循环井间互相扰动越强，而多个循环井间的相互扰动对单井的循环效率起着促进作用，因而整体的循环效率变高。

表 3  优化模型计算最优设计方案参数统计表

Table 3.  Statistics of optimal design scheme parameters by optimal model

方案编号	第i口注水筛顶部与含水层顶板的距离/m			第i口井的抽水筛顶部与含水层顶板的距离/m			第i口井的注水筛长度/m			第i口井的抽水筛长度/m		抽注水量/(m3∙d−1)	井间距/m
	i=1	i=2		i=1	i=2		i=1	i=2		i=1	i=2
1	9.1	8.5		1.6	2.3		2.5	1.5		2.5	1.5	185	25
2	10.5	10.8		1.0	1.3		1.5	1.2		1.7	1.5	185	20
3	0.7	11.0		11.1	0.8		1.4	1.3		1.6	1.2	163	13
4	0.9	9.6		9.9	1.3		2.3	2.3		1.1	1.6	163	10

 | Show Table

DownLoad: CSV

表 4  模拟模型和优化模型计算结果对比表

Table 4.  Comparison of simulation and optimization model results %

方案编号	循环效率		循环效率相对误差	污染物去除效率		污染物去除效率相对误差
	模拟模型	优化模型		模拟模型	优化模型
1	87	90	3.4	29	30	3.4
2	98	99	1	26	27	3.8
3	92	95	2.5	33	34	1.6
4	98	97	1.1	30	31	2.9

 | Show Table

DownLoad: CSV

方案3和方案4的两个循环井抽注模式不同。方案3循环井循环效率低于方案4，而污染物去除效率高于方案4；2种方案的抽注水量相同，但方案4的井间距小于方案3。井间距越小，循环井间互相扰动越强，因此抽注模式相反方案中循环井间相互扰动同样对单井的循环效率有积极的影响。

多目标优化设计模型得到的优化设计方案均可以得到较高的循环效率和污染物去除效率，因此在双循环井的实际应用中，可进一步根据优化时的侧重点选择其中更符合实际需求的设计方案。例如当污染场地修复工期较长时，可侧重于保证修复试剂的充分利用，同时为防止污染物从双循环井的环流场逃逸而污染下游地下水，可选择循环效率最优的方案2作为最终优化设计方案；而当修复工期较短时，可侧重于保证场地污染物去除效率更高，则可选择污染物去除效率最高的方案3作为最终的优化设计方案。方案2和方案3的设计参数如图8所示。

图 8  优化方案示意图

Figure 8.  Schematic diagram for optimization scheme

下载: 全尺寸图片 幻灯片

将表3列出的优化设计参数代入到前述双循环井模拟模型中求解，将优化模型与模拟模型得到的结果进行相对误差计算，其计算式如式(8)所示。

式中：$\sigma$表示相对误差，${x}'$表示优化模型求解结果，x表示模拟模型求解结果。

计算对比结果如表4所示，循环效率的最大相对误差为3.4%，污染物去除效率的最大相对误差为3.8%，相对误差均低于±5%，证明该NSGA-Ⅱ优化模型的计算结果可靠。

3.   结论

1) 考虑到水文地质参数的不确定性，深度学习方法中基于卷积神经网络的替代模型与地下水数值模型有很大的适应性，可以稳健地反映地下水数值模型的输入输出关系。本研究选取确定性系数 (R²) 、平均绝对百分比误差 (MAPE) 2个指标去评价替代模型的精确度。结果表明替代模型一方面减少了计算负荷，同时在预测双循环井循环效率时的R²达到0.972，MAPE达到0.008；污染物去除效率的R²达到0.986，MAPE达到0.009，计算精度高。说明构建基于卷积神经网络的替代模型在双循环井的优化设计中是可行的。

2) 通过对双循环井优化设计方案的验证，发现得到的双循环井运行效果较好，验证了优化设计方案的可靠性，也表明了优化设计在双循环井实际使用工作中的必要性。本研究在获得的Pareto最优解集中选取4个典型方案进行对比分析，并最终确定了两个不同优化目标下的最优方案。本研究中的双循环井优化设计方法可以获取不同水文地质条件下双循环井结构及运行参数、井间距等指标，从而保证循环井的运行效果，达到该技术在地下水污染治理中的合理应用。