1 引言
地下水污染修复实际上是一个既考虑治理效果又考虑治理成本的复杂的多目标优化问题.当前,求解多目标问题最有效的方法是多目标智能优化算法,一方面,该算法可以并行地处理一组可能的解(群体),能在一次算法过程中找到Pareto最优集中的多个解;另一方面进化算法不局限于Pareto前端的形状和连续性,易于处理不连续的、凹形的多目标函数,能够有效地克服古典方法的局限性.如基于小生境的Pareto遗传算法(niched Pareto genetic algorithm,NPGA),快速非支配遗传算法(nondominated sorting genetic algorithm Ⅱ,NSGA Ⅱ)等.以上提出的多目标遗传算法具有良好的全局搜索能力,然而在搜索的中后期,算法收敛于真实最优解的效率大幅降低,而且遗传算法是以种群为单位的群体搜索方法,对于复杂的非线性规划问题不能保证收敛到最优解.鉴于以上不足,国内外研究者提出了多种对传统多目标智能优化算法的改进算法.本文对遗传算法初期搜索到的优化解进行局部搜索,形成了一种新的混合多目标优化模式.该混合算法即是利用遗传算法的全局搜索能力搜寻最优解,然后将进化的解作为初始解进行局部搜索以搜寻到真实的Pareto最优解.
本文将NSGAII与一种迭代式的局部搜索算法(Hill Climber with Step,HCS)相结合,开发了一种新的混合多目标遗传算法NSGAII-HCS.利用CONV1和ZDT6两个经典的多目标优化函数对NSGAII-HCS的性能进行测试.最后,将NSGAII-HCS与地下水流模拟软件MODFLOW和溶质运移模拟软件MT3DMS相耦合,并应用到一个理想的二维地下水污染修复管理模型中,结果分析表明该方法可为地下水污染治理提供多样的和收敛的Pareto管理策略,是一种稳定可靠的多目标优化方法.
2 混合多目标遗传算法
2.1 快速非支配遗传算法
NSGAII是一种可以快速进行Pareto排序的非支配多目标遗传算法,如图 1所示.该算法的优点主要包括以下3个方面：① 提出快速非支配排序法,降低算法的复杂度;② 提出拥挤度与拥挤度比较算子,代替了需要指定共享半径的适应度共享策略,并在同一级排序中将拥挤度作为个体的胜出标准,保持了种群的多样性;③ 引入精英策略,有利于保存种群进化过程中的优良个体,加速了种群的收敛.以上优点使NSGAII成为一种广泛应用的多目标智能优化算法,并作为与其它改进算法的比较基准.
图 1 NSGAII-HCS的计算流程图
2.2 局部搜索算法
HCS是一种基于最优搜寻方向的迭代式局部搜索方法.由于地下水污染修复管理模型是复杂的非线性规划问题,不能准确计算目标函数的梯度信息,因此本文采用无梯度计算的局部搜寻方法.
为了详细阐述HCS的具体搜索过程,可以将HCS作为一种对遗传算法进化的Pareto解进行局部搜索,然后输出局部最优解的独立的算法,而HCS与NSGAII的耦合过程如图 1所示.HCS的计算流程图如图 2所示.为了决定局部搜寻的方向,需要利用初始解的给定邻域半径内的邻域解.邻域解的定义即：
(1)
式中,X0=(x01,x02,…,x0m);Xn=(xn1,xn2,…,xnm).X0是初始解决策变量的向量形式;Xn 是邻域解决策变量的向量形式;m是决策变量的个数;x0i是初始解第i个决策变量;H是在寻优过程中选取邻域解的参考解集,以进化过程中父代与子代解集作为参考解集并从中选取符合要求的邻域解,这样可以减少因产生邻域解而需要的目标函数评价,增加算法的效率;r是邻域半径.若邻域解Pareto主导初始解并以Xn
(2)
式中,fi(X)表示第i个目标函数,k表示目标个数.
图 2 HCS的计算流程图
在流程图 2中,Nd是局部搜索的最大迭代次数;T是衡量决策变量变化的参数;Na是搜寻最优T值的最大试验次数;Xni表示第i次迭代选取的邻域解的决策变量;εy是计算最优变化量的参数,一般选取2~5;L0是近似表达目标函数梯度信息的变量;a是Nd次迭代搜寻过程中决策变量的平均变化量;其余符号如前所述.根据初始解与邻域解的Pareto主导关系可以将HCS分为Hill Climber与Step两个部分.当在Nd次迭代搜索内存在邻域解与初始解是Pareto主导关系,则进行最优方向的搜索.Lara et al提出以二次多项式的形式估计函数F(T),即F(T)=aT2+bT+c.当函数值F(T0),F(T1),F(T2)存在以下关系,即：
(3)
则认为存在最优值T=-b/(2a).如图 2所示,当经过Na次试验未满足式(3)的条件,则以最小变化量T0作为最优值.如果经过Nd次搜索,邻域解与初始解均处于非支配的关系,则认为初始解是局部最优Pareto解.为了保持解的多样性,Lara et al提出利用Nd个邻域解沿着Pareto解的锋面方向继续搜索,具体搜索过程如图 2所示.
HCS需要输入的主要参数是搜寻半径r与最大迭代搜寻次数N.N值越大,局部搜寻的解更可能达到最优解,但同时会增加目标函数的评价次数,因此,Lara et al.建议N值一般取5～10.由于是局部搜索,搜寻半径不宜过大,一般取值0.1～ 0.3.启发式搜寻最优T值的方法需要指定最大试验次数Na,Na值过大则降低算法的效率,本文根据Lara et al建议取值为3.
2.3 组合优化的方法
混合多目标算法的两个重要问题是如何选择个体作为局部搜寻的初始解与如何平衡全局搜索与局部搜索之间的关系.遗传算法的搜寻特点是在初期阶段的全局搜索能力可以快速找到接近最优解的个体,但是到后期阶段收敛到最优解的速度减慢.因此,在实施局部搜索时,应选取当前代进化得到最优个体作为初始解.全局搜索与局部搜索的平衡关系可以由当前选择个体作为初始解的概率来调节.如果在进化初期阶段选择概率较大,则搜寻的解易陷入局部最优的陷阱.
理想算例将局部搜寻设定在进化的10代以后,使遗传算法的全局搜索能力得到发挥.然后,每隔5代进行局部搜索,每一代搜索的概率服从以下函数分布：
(4)
式中,pl是局部搜索概率,rw,ra是函数分布参数,pmax是最大的局部搜索概率,m是进化过程中种群的局部搜索次数.在求解理想算例时,将pmax设定为0.3,m设定为20,函数图像如图 3所示.
图 3局部搜索概率分布曲线
从图 3中可以看到遗传算法进化的初始阶段局部搜寻概率小,有利于全局寻优,进化的中期增加局部搜索的能力改善解的精度,到后期可以减少局部搜索避免陷入局部最优的陷阱同时减少函数评价的次数,提高算法的效率.
2.4 标准测试函数的检验
为了测试NSGAII-HCS算法的性能,选用多目标问题的经典测试函数CONV1与ZDT6,其函数表达式分别是：
CONV1目标函数：
(5)
(6)
约束条件：
(7)
ZDT6目标函数：
(8)
(9)
约束条件：
(10)
(11)
标准测试函数的多目标遗传算法参数如表 1所示.如图 4a所示,对于凸形的多目标函数CONV1,NSGAII-HCS得到的解更趋近真实的Pareto解集,但是进化的总代数较小,因此不能使所有的解完全收敛到真实的Pareto锋面.如图 4b所示,对于凹形的多目标函数ZDT6,NSGAII-HCS在种群进化到第100代时,解完全收敛到真实的Pareto锋面,表明了该算法具有在保持解多样性的同时,可以使Pareto解达到局部最优性;而NSGAII在种群进化到第250代时,得到的Pareto锋面与真实的Pareto锋面的平均距离为0.16,仍与真实Pareto解存在较大的差距.标准函数的测试结果表明了NSGAII-HCS与NSGAII相比能搜寻到真实的Pareto解,具有明显的优势.
表 1 多目标遗传算法参数
图 4 NSGAII与NSGAII-HCS求解CONV1(a)与ZDT6(b)的数值计算结果
3 算例分析
3.1 算例概述
本理想算例来自文献Zheng and Wang,目的是利用NSGAII-HCS来设计一个抽取-处理(pump-and-treat,PAT)系统,以同时达到治理成本最小和含水层中剩余污染物最小这两个目标.算例场地信息可参见相关文献.该算例的地下水污染优化管理模型即：
目标函数：
(12)
(13)
约束条件：
(14)
(15)
其中,f1是治理成本,a1是安装抽水井的总费用(设为CNY150000),N是PAT系统控制井的数量,Nw是非零流量井的数量,a2是处理单位体积污水需要的费用(设为CNY 0.76· m-3),Qi是第i口井的流量(m3 · d-1),Δti是第i口井持续抽水的时间(d),f2是治理结束后剩余污染物质量与初始质量的百分比,massinitial是含水层初始污染物质量(kg),massend是剩余污染物质量(kg),C*是污染物浓度约束区域的最大允许浓度值(设为20 g · m-3),mc是污染物浓度约束区域节点的数量.
对于研究区的所有节点的浓度通过调用地下水水流模型MODFLOW与溶质运移模型MT3DMS计算,研究区的各点的污染物浓度取决于抽水量,抽水持续时间以及初始污染物浓度.函数表达式如下：
(16)
3.2 结果与讨论
本文将NSGAII-HCS和NSGAII分别与常用的地下水流模型MODFLOW和溶质运移模型MT3DMS相耦合,模拟优化理想条件下二维的地下水污染修复系统.两种算法在全局搜索的遗传算法上使用的参数相同,即种群大小为100,进化的总代数为100代,交叉概率为0.9,突变概率为0.05.对于局部搜索,最大的局部搜索概率为0.3,邻域搜索半径为0.2,每隔5代进行局部搜索,最大的迭代搜索次数为5次.本文引入多样性、收敛性和均匀性三个指标来评价NSGAII-HCS的寻优性能.
3.2.1 算法性能的衡量指标
多样性的衡量是通过将两种算法计算的解合并,并重新进行非支配解的排序得到最优解集,将每一种算法在最优解集中占有的Pareto解的个数作为多样性指标.收敛性指标是参照Deb提出的计算方法,即选取若干个真实的Pareto解作为参考点,然后计算多目标算法得到的权衡曲线与参考点构成的权衡曲线的距离.计算方程如下：
(17)
式中,N是种群大小,dij是第i个进化的个体到第j个参考点的距离,R是选取的参考点个数.Z值越大表示计算的解与真实解的差距越大,解的收敛性差,反之,收敛性好.由于地下水污染优化管理模型真实的Pareto解集未知,所以从两种算法合并后的最优解集中选取若干个均匀分布的Pareto解作为参考点.
均匀性指标的计算方程如下：
(18)
式中,N是Pareto解的个数,di是权衡曲线上相邻两点的距离, 是di的平均值.若解的分布比较均匀,则相邻点距离与平均距离相近,Δ较小.反之,Δ较大.
3.2.2 优化结果的对比与分析
将NSGAII与NSGAII-HCS计算得到的结果进行对比分析(图 5).在剩余污染物百分比大于15%时,NSGAII-HCS可以搜索到收敛性与多样性更优的解,同时保证所有的Pareto解在权衡曲线上均匀分布.因此,NSGAII-HCS能寻找到精度更高的Pareto解,并能提升解的多样性,具有更加明显的优势.
图 5 NSGAII与NSGAII-HCS运行到100代时的Pareto解集
为了比较在种群进化过程中两种不同算法的收敛性,可以计算NSGAII与NSGAII-HCS的每一代种群收敛性指标的差值ΔZ.ΔZ为正值表示NSGAII-HCS计算的解与参考点的距离较小.如图 6所示,在进化的初期由于遗传算法的全局搜索的能力,NSGAII算法发挥较大的作用,而局部搜索可能使解陷入局部最优值.在进化的中后期,遗传算法寻找到接近真实Pareto解的最优解而继续搜索真实解的效率逐渐降低.NSGAII-HCS具有的局部搜索能力可以很大程度上改善遗传算法在中后期搜索的Pareto解,使种群的收敛速度加快.如表 2所示,两种算法运行到100代时,NSGAII-HCS的Z值较小,表明解的精度更高.Pareto解集的多样性与均匀性的提升有助于决策者根据自己的要求选择最优的管理方案.从表 2中也可以看到,NSGAII-HCS在保证Pareto解局部最优性的同时,可以提供多样的解并且解的分布较均匀.具体参见 污水处理技术资料或污水技术资料更多相关技术文档。
图 6两种算法收敛性指标的差值与进化代数的关系
表 2 NSGAII与NSGAII-HCS运行到100代时Pareto解比较
4 结论与展望
1) 为了提高多目标遗传算法Pareto解的局部最优性,将NSGAII与一种迭代式的局部搜索算法HCS相结合,开发了一种新的混合多目标遗传算法NSGAII-HCS.HCS通过利用邻域解确定最优的搜寻方向,使种群朝着真实解的方向收敛,提高了收敛速度与解的精度.同时,在寻找到最优解后,HCS利用局部搜索的解沿着Pareto锋面的方向继续搜索,有助于改善解的多样性.通过选择两个经典的多目标测试函数检验算法的性能,结果表明NSGAII-HCS的解更加接近或完全收敛到真实的Pareto解.
2) 将NSGAII-HCS应用于理想的二维地下水污染修复管理模型中,从Pareto解的收敛性、多样性和均匀性方面分析,可以得到该算法与NSGAII相比能寻找到精度更高、跨度更大、沿权衡曲线分布更均匀的Pareto解.NSGAII-HCS同时具有遗传算法的全局搜索能力和HCS的局部搜索能力,而在进化过程中,NSGAII-HCS在中后期的局部搜索使收敛速度明显加快.在地下水污染修复优化管理模型中,精度更高的解可以节约污染物治理成本,同时有助于选择最高效的管理策略,使地下水污染治理的目标达到.NSGAII-HCS的优点则是能加速收敛到最优解,使决策者能选择最理想的治理方案.
3) 尽管本文算例是一理想的二维均质等厚各向同性承压含水层系统,但由于NSGAII-HCS耦合的水流和污染物运移模型适用于任何多孔介质的地下水系统,因此该优化模型同样可推广应用到实际含水层的地下水污染修复管理.当然,其优化效果如何,有待于后续应用实践的检验.另一方面,NSGAII-HCS每隔几代进行局部搜索的策略可能减小了算法的效率,如何进行自适应的局部搜索有待于一步研究.将NSGAII-HCS应用于实际场地条件下的地下水污染修复管理,充分发挥该算法的优点是未来的发展方向.