引言

选址问题是运筹学中的经典问题之一，广泛应用于生产、物流等领域，如利用一些特定的方法进行工厂、仓库、急救中心的选址等。现代网络选址理论最早由Hakimi^[1]提出，利用该理论可以在网络中找到一个或多个设施的最佳位置，并且在满足客户需求的前提下达到费用最小化、客户满意度最大化等目标。在众多选址目标中有一类设施，如垃圾填埋场^[2]、沼气工厂^[3]、生物柴油厂^[4]、水泥制造厂^[5]等是不受欢迎的，针对此类设施的选址问题称为厌恶设施选址问题。这类厌恶设施如选址不当，不仅会造成环境污染，还可能危害生命安全，因此，研究厌恶设施的选址问题有着非常重要的现实意义。

厌恶设施选址问题中主要有两大类模型：1-maximin模型和1-maxisum模型。二者的不同之处在于目标函数的设定。前者的选址目标是保证选址位置与所有客户权重距离中的最小权重距离达到最大；后者的目标是让选址位置与所有客户的权重距离之和达到最大。

目前针对1-maximin模型的研究较多。Melachrinoudis等^[6]证明了1-maximin模型的目标函数在每条边上是分段线性凹函数。Tamir^[7]用动态数据结构的方法，给出了通路和星图等特殊结构的求解算法。Colebrook等^[8]在文献[6]的工作基础上，基于去除无效边的思想提出了UnCenter算法，数值实验表明该算法找到最优解所需的时间有所减少。Nadirler等^[9]通过去除冗余的整数变量的方法简化了混合整数规划模型，并用改进的分支定界法来求解1-maximin模型。

20世纪70年代中期，Church等^[10]提出了1-maxisum模型，并基于有限的易于生成的可行解区域内至少有一个最优解的思想给出了1-maxisum模型的求解算法。Colebrook等^[11]利用新的上界方法求解了1-maxisum模型，并在数值实验部分给出了不同密度图下的计算结果。他们还对厌恶设施选址模型以及求解算法进行了总结^[12]。此外，Yamaguchi^[13]通过集体谈判的方式来选取厌恶设施的位置。Song等^[14]提出两种替代的数学优化模型目标是减少厌恶设施选址对人们情绪的总体影响程度。

以上文献中的厌恶设施选址位置都是用精确解的方法给出，但是其计算时间慢。近年来启发式算法(如蚁群算法、遗传算法和粒子群算法等等)逐渐成为求解优化问题的有效方法，其应用已深入到各个领域。粒子群算法是由Eberhart等^[15]提出的一种启发式算法，在求解传感器网络定位、常压塔稳态操作、车辆路径问题、交通流预测、数据网络等许多问题时表现出良好的性能^[16-20]。1-maximin模型的目标函数在每条边上是单峰函数^[6]，而黄金分割法是求解单峰函数极值点的很好选择。基于此，本文提出了粒子群与黄金分割法相结合的PSO-GS算法，用粒子群算法搜索最优解所在的边，再用黄金分割法在最优边上寻优，得到厌恶设施选址的最优位置，并将PSO-GS算法与UnCenter, Newalgorithm算法进行比较，以期为厌恶设施选址提供一种新的更高效的方法。

1 厌恶设施选址模型 1.1 距离的定义

令N=(V, E)是一个无向赋权简单联通图，其中V={v₁, v₂, …, v_n}是具有n个节点的节点集合，节点数表示客户数目；E={e_j=(v_s, v_t):v_s, v_t∈V, j=1, 2, …, m}是具有m条边的边集合，边数表示在道路系统中可以连通的路径数目。每个节点的权重w_i和每条边e_j的长度l_j分别由公式(1)、(2)计算

$ {v_i} \in V \to w\left( {{v_i}} \right) = {w_i} > 0 $

(1)

$ {e_j} = \left( {{v_{\rm{s}}},{v_{\rm{t}}}} \right) \in E \to l\left( {{e_j}} \right) = {l_j} > 0 $

(2)

图 1为5个客户8条边的简单联通图。图中每个点上括号外的数字是点的标号，括号内的数字是客户的权重；每条边上方括号中的数字表示该条边的长度。

从图 1中任意取两点v_s, v_t(v_s, v_t∈v₁…v₅)，定义其距离为两点之间的最短距离，用d(v_s, v_t)表示。边e=(v_s, v_t)∈E上的任意一点x到V上的点v_i的距离定义为

$ d\left( {x,{v_i}} \right) = \min \left\{ {x + d\left( {{v_{\rm{s}}},{v_i}} \right),{l_e} - x + d\left( {{v_{\rm{t}}},{v_i}} \right)} \right\} $

(3)

1.2 1-maximin和1-maxisum选址模型

对于图N中给定点x，定义

$ {f_{\min }}\left( x \right) = \mathop {\min }\limits_{{v_i} \in V} {w_i}d\left( {x,{v_i}} \right) $

(4)

$ {f_{{\rm{sum}}}}\left( x \right) = \sum\limits_{{v_i} \in V} {{w_i}d\left( {x,{v_i}} \right)} $

(5)

其中f_min(x)表示设施位置到所有客户的最小权重距离，f_sum(x)表示设施位置到所有客户的权重距离之和。

厌恶设施1-maximin模型为

$ \mathop {\max }\limits_{x \in N} \mathop {\min }\limits_{{v_i} \in V} {w_i}d\left( {x,{v_i}} \right) = \mathop {\max }\limits_{x \in N} {f_{\min }}\left( x \right) $

(6)

即找到设施所在的位置并使f_min(x)达到最大。

厌恶设施1-maxisum模型为

$ \mathop {\max }\limits_{x \in N} \sum\limits_{{v_i} \in V} {{w_i}d\left( {x,{v_i}} \right)} = \mathop {\max }\limits_{x \in N} {f_{{\rm{sum}}}}\left( x \right) $

(7)

即找到设施所在的位置并使f_sum(x)达到最大。

2 PSO-GS算法

考虑到1-maximin和1-maxisum选址模型的特点，本文采用启发式算法对其进行求解。本文提出的PSO-GS算法包含了3个环节：首先采用粒子群算法快速搜索最优解的大致位置，然后用加步探索法在第一阶段搜索到的位置上找到一个单峰函数区间，最后利用黄金分割法在加步探索法得到的区间内搜索，从而求得模型的最优解。

2.1 粒子群算法

标准的粒子群算法中速度、位置和惯性权重的更新公式分别如式(8)、(9)、(10)所示^[21]

$ v_i^{k + 1} = \omega v_i^k + {c_1}\xi \left( {p_i^k - x_i^k} \right) + {c_2}\eta \left( {p_g^k - x_i^k} \right) $

(8)

其中，v_i^k是第i个粒子第k次迭代的速度；x_i^k是第i个粒子第k次迭代的位置；p_i^k是第i个粒子第k次迭代时的个体最优值；p_g^k是第k次迭代中所有粒子经过的最优位置；c₁, c₂是学习因子，表示向个体最优和全局最优学习的能力；ξ, η是[0, 1]区间上均匀分布的随机数；ω是惯性权重。

$ x_i^{k + 1} = x_i^k + v_i^{k + 1} $

(9)

$ \omega = {\omega _{\max }} - \frac{{{\omega _{\max }} - {\omega _{\min }}}}{{{I_{\rm{p}}}}} \times c $

(10)

其中ω_max, ω_min分别是最大、最小权重；I_p为粒子群迭代的总次数；c为当前的迭代次数。从公式(10)中可以看出惯性权重随着迭代次数的增大而逐渐减小。

2.2 加步探索法

加步探索法是一种区间试探法，其基本思路是从一个初始点出发，按照一定的步长进行探索找到一个区间，使函数在此区间呈现低—高—低的形态。

2.3 黄金分割法

黄金分割法是求解单峰或单谷函数最优解的有效方法，其区间更新公式如式(11)、(12)

$ {x_{{\rm{l}}\left( {k + 1} \right)}} = {a_k} + 0.382 \times \left( {{b_k} - {a_k}} \right) $

(11)

$ {x_{{\rm{r}}\left( {k + 1} \right)}} = {a_k} + 0.618 \times \left( {{b_k} - {a_k}} \right) $

(12)

其中x_l(k+1)和x_r(k+1)分别是第k+1次迭代的左右分位点，a_k和b_k分别是第k次迭代的左右端点。通过比较函数值的大小，黄金分割法得到区间长度缩小为原来长度0.618倍的包含最优解的区间。

2.4 PSO-GS算法

本文基于1-maximin选址模型的特点以及粒子群算法和黄金分割法的优势，再结合加步探索法，提出的PSO-GS算法流程图如图 2。

3 数值实验 3.1 实验设置

将本文PSO-GS算法与UnCenter和Newalgorithm算法^{[8, 11-12]}进行了比较。用MATLAB2011编写程序，在联想M8400t-N000上运行，所用机器的主频3.4GHz，内存4G，操作系统为XP。算法中各参数设置如下：学习因子c₁=c₂=1.5，黄金分割法的精度ε=10^－6，步长增加的倍数α=2，加步探索的初始步长h=0.0001，最大最小权重分别为ω_max=0.9、ω_min=0.4，粒子群的迭代次数为I_p=5，总迭代次数为I_n=50。

3.2 PSO-GS与UnCenter，Newalgorithm算法性能比较

将PSO-GS算法与UnCenter和Newalgorithm算法的运行时间进行比较，统计的计算时间均不包括用弗洛伊德算法得到最短路的时间。

本文参照文献[8]、[11-12]生成数据给定客户数目后，根据不同的密度d(即点的数目n与边的数目m的比例关系)得到边的数目，其计算公式为m=d×n×(n－1)/2。本文取密度分别为1/2、1/4、1/8、1/16，不同密度下客户数目取值如表 1。计算时所选取的参数为：客户的权重范围1~10，边的长度范围1~50.

1-maxisum模型的计算结果如图 3所示，1-maximin模型的计算结果如图 4所示。

由图 3可知，不同密度的道路系统中，随着客户数目的增多，PSO-GS算法求解1-maxisum选址模型得到最优解的效率优于Newalgorithm算法。图 4表明，在低密度(1/16和1/8)系统中，客户数目在375~875时，PSO-GS算法的计算时间少于UnCenter算法，当客户数目增大到1000时两者的计算时间比较接近；在高密度(1/4和1/2)系统中，客户数目大于300时，PSO-GS算法的计算时间大约为UnCenter算法的一半，这表明PSO-GS算法更适用于求解高密度问题。在现实中，道路往往是错综复杂的，高密度系统较低密度系统更常见，因此，PSO-GS算法有广泛的应用。

4 结束语

本文提出了求解厌恶设施1-maximin和1-maxisum选址模型的混合启发式算法：PSO-GS算法。数值实验结果表明该方法的计算时间优于UnCenter和Newalgorithm算法，是求解1-maximin和1-maxisum选址模型的有效算法。