求解旅行商问题的联合算子模拟退火算法

李昌兴，黄 杉

(1.西安邮电大学 理学院,陕西 西安 710121；2.西安邮电大学 通信与信息工程学院,陕西 西安 710121)

模拟退火算法(Simulated Annealing Algorithm，SAA)不仅可以解决连续函数优化问题，也成功应用于求解组合优化问题[10]。对于组合优化问题，必须使用操作算子产生和接受新解。操作算子的策略非常重要，直接影响着模拟退火算法的优化性能[11]。在反序、移位和交换算子的基础上，文献[12]使用逆转操作算子，文献[13]使用强制随机邻域变换策略，这些操作算子改善了模拟退火算法的局部搜索能力。文献[14]提出融合最近、最远和随机插入操作的统一插入操作和最坏删除操作，通过自适应邻域提高了搜索性能。文献[15]提出自适应升温控制策略，平衡了全局搜索与局部搜索的关系，获得更好的全局优化能力。文献[16]建立了模拟退火算法的弛豫时间模型，给出退火过程Markov链长度的理论估计，提出了Markov链的动态调整策略，对于模拟退火算法的参数选择具有重要的理论指导意义。但是，以上研究并未充分利用不同操作算子的内在联系与相互作用，导致模拟退火算法精度不足，性能不高。

针对模拟退火算法中操作算子的生成策略问题，研究反序、移位和交换算子的特征与相互关系，发现交换操作与反序操作在机理上具有等价关系。利用这种等价关系，拟提出一种新的交换/反序联合算子(Joint Operator，JO)模拟退火算法，并利用不同规模的TSP问题对联合算子进行仿真实验。

1 模拟退火算法求解旅行商问题

1.1 模拟退火算法

模拟退火算法从较高的初始温度出发，在当前解的邻域不断通过随机扰动产生新的候选解，以一定的概率接受候选解。模拟退火算法包括内外两重循环：外循环是温度循环，控制退火温度降低的过程；内循环是在每一温度下迭代产生新解，循环次数也称Markov链长度。

模拟退火算法的基本步骤如下。

步骤1随机产生初始解X0。

步骤2在温度Tk下进行Lk次循环迭代。由当前解Xold产生新的候选解Xnew，并按Metropolis接受准则，以一定的概率接受候选解Xnew作为当前解。候选解与当前解目标函数值的差ΔE及候选解的接受概率P的表达式分别为

ΔE=E(Xnew)-E(Xold)

(1)

(2)

式中：k表示温度循环变量；E(Xold)表示当前解Xold的目标函数值；E(Xnew)表示新解Xnew的目标函数值。

Metropolis准则允许以一定的概率接受恶化解,使算法可以跳出局部最优解，是模拟退火算法能收敛于全局最优解的关键。

步骤3按照设定的退火控制策略缓慢降低退火温度，直到达到温度终值，完成退火过程，这时的当前解就是达到最低能量状态的最优解。

1.2 求解旅行商问题的操作算子

模拟退火算法要从当前解的邻域中产生新的候选解，解的表达形式和邻域结构对算法收敛非常重要。组合优化问题中，目标函数不仅依赖于解的取值，而且与解的结构、次序相关。旅行商问题的路径长度仅取决于解的排列次序，通常用城市编号序列表示问题的解。

为了满足组合优化问题的约束条件，需要通过操作算子对当前解序列进行变换操作，改变某些点的排列次序产生新解。基本的操作算子是交换算子、移位算子和反序算子[1，11]。

交换算子是将当前路径Snow中的第i个城市Ci与第j个城市Cj的位置交换，得到新路径Sswap，当前路径和新路径的表达式分别为

Snow=C1…Ci-1(Ci)Ci+1…Cj-1(Cj)Cj+1…Cn

Sswap=C1…Ci-1(Cj)Ci+1…Cj-1(Ci)Cj+1…Cn

移位算子是将当前路径Snow中的第i个城市Ci移动到第j个城市Cj之后的位置，得到新路径Smove，当前路径和新路径的表达式分别为

Snow=C1…Ci-1(Ci)Ci+1…Cj-1(Cj)Cj+1…Cn

Smove=C1…Ci-1Ci+1…Cj-1(Cj)(Ci)Cj+1…Cn

反序算子是将当前路径Snow中从第i个城市Ci到第j个城市Cj之间的城市排列顺序反向翻转，得到新路径Sinv，当前路径和新路径的表达式分别为

Snow=C1…Ci-1(CiCi+1…Cj-1Cj)Cj+1…Cn

Sinv=C1…Ci-1(CjCj-1…Ci+1Ci)Cj+1…Cn

1.3 多个新解的竞争机制

模拟退火算法本质上是基于随机搜索的优化方法，操作算子是模拟退火算法求解旅行商问题的核心机制。操作算子的运行策略不仅决定了新解的可行性，而且与获得优质解的概率密切相关，影响着模拟退火算法的优化性能。

交换算子、移位算子和反序算子所产生的新路径与当前路径相比，分别有8段、6段、4段路径的改变。随着当前路径的不断优化，通过这些算子获得更好路径的概率越来越小。增大循环次数即Mar-kov链长度可以提高优化性能，但这也带来了计算量的增大。

在每次循环产生的多个新解择优使用的竞争机制，可以提高获得更好路径的概率。文献[17]提出了多粒子模拟退火算法，每次扰动产生多个新解，选取最优者作为候选解。文献[18]在每次循环中通过断裂、倒置及重组操作产生6个新路径，选择最好路径作为候选解。类似地，文献[19]提出了多线程的并行算法，可以并行产生多个新解，从中找出最好的候选解。

引入产生多个新解的竞争机制，可以提高模拟退火算法的优化性能，但也增大一定的计算量。由于计算量的增长与新解个数的增加近似成正比，这与直接增大循环次数相比算法未必更为有效。只有在循环中产生多个新解，而又不增大计算量的前提下，才能真正提高模拟退火算法的总体性能。这就需要研究现有操作算子的运行特征与相互关系，构造同时产生多个新解的操作算子。

2 联合操作算子

2.1 交换、反序和移位操作的关系

模拟退火算法由新解与当前解的目标函数差ΔE计算新解的接受概率。在旅行商问题中，使用交换算子、反序算子或移位算子产生新路径时，并不需要直接计算新路径的总长度E(Xnew)，可以通过不同变换操作的特征直接计算ΔE，以减少计算量。

交换操作的路径差ΔEswap、反序操作的路径差ΔEinv与移位操作的路径差ΔEmove的表达式分别为

ΔEswap=[d(Ci-1,Cj)+d(Ci,Cj+1)]-
[d(Ci-1,Ci)+d(Cj,Cj+1)]+
[d(Cj,Ci+1)+d(Cj-1,Ci)]-
[d(Ci,Ci+1)+d(Cj-1,Cj)]

(3)

ΔEinv=[d(Ci-1,Cj)+d(Ci,Cj+1)]-
[d(Ci-1,Ci)+d(Cj,Cj+1)]

(4)

ΔEmove=[d(Ci-1,Ci+1)+d(Cj,Ci)+d(Ci,Cj+1)]-
[d(Ci-1,Ci)+d(Ci,Ci+1)+d(Cj,Cj+1)]

(5)

式中，d(Ci,Cj)表示城市Ci与Cj之间的距离。

对比这3种操作所产生的路径差，可以发现具有很多相同的距离项。特别是式(4)反序操作的路径差ΔEinv中的全部4段距离，都出现在式(3)交换操作的路径差ΔEswap中，并且正负号也相同。因此，在计算式(3)交换操作的路径差ΔEswap时，可以先计算式(4)中ΔEinv的4段距离，再计算剩下的4段距离。就能在计算ΔEswap的同时得到ΔEinv的值，并不增大计算量，同时还获得了2条新路径的ΔE。

(6)

于是，交换操作的路径差等于两个嵌套的反序操作的路径差之和，其表达式为

(7)

因此，可以将交换操作分解为两步：先将当前路径Snow中从城市Ci到城市Cj之间的城市排列顺序反向翻转得到新路径Sinv；再将路径Sinv中从城市Cj-1到城市Ci+1之间的城市排列顺序再次翻转，两次翻转后得到的路径就是通过交换算子所得到的新路径Sswap。第一步反向翻转路径Sinv和第二步反向翻转路径Sswap的表达式分别为

Sinv=C1…Ci-1Cj(Cj-1…Ci+1)CiCj+1…Cn
Sswap=C1…Ci-1Cj(Ci+1…Cj-1)CiCj+1…Cn

这说明交换操作等价于两个嵌套的反序操作的叠加复合。

2.2 交换-反序联合算子

3 实验结果及分析

3.1 实验方案

为检验提出的交换-反序联合算子的优化性能，从 TSPLib 案例库中选取 4 个不同规模的旅行商问题Eil51、Eil76、Eil101和Ch150作为测试案例，其城市数量分别为51、76、101、150，已知的最优路径长度分别为426、538、629和6 528。

所提算法实验环境为CPU：Intel(R)Core(TM) i5-4460S@2.90 GHz，内存为8 GB/1 600 MHz，操作系统为Window s7(64位/Service Pack 1)，采用 Python 3.8编程。

关于模拟退火算法的控制函数选择，控制温度按照Tk=αTk-1指数衰减，衰减系数取α∈(0.9,1.0)，按照式(2)Metropolis 准则接受新解。

考虑模拟退火算法的性能受到优化参数的影响，对每个测试案例分别采用不同的参数各进行两组实验。第1组参数为：初温100，终温1，衰减系数0.95，共90个温度状态，Markov链长度1 000；第2组参数为：初温100，终温1，衰减系数0.98，共228个温度状态，Markov链线性增大，平均链长1 000。由于模拟退火算法的优化结果有一定的随机性，在每组参数下各做10次重复实验，得到10次重复实验的最小值、平均值、最大值，具体对比结果分别如表1至表4所示。

表1 Eil51问题实验结果对比

表2 Eil76问题实验结果对比

表3 Eil101问题实验结果对比

表4 Ch150问题实验结果对比

3.2 实验结果的分析

对表1至表4所示的联合算子模拟退火算法的优化结果与使用现有的移位、交换、反序算子以及组合移位/交换/反序算子的模拟退火算法的优化结果进行比较分析，可以得到4个不同规模的旅行商问题在不同测试参数下的结果。使用交换-反序联合算子时最小值、平均值、最大值都是最好的，表明联合算子的性能优于移位、交换、反序算子及其组合方案，这得益于联合算子与其他算子相比的搜索数量更大。

使用联合算子和现有算子的模拟退火算法在各自的第2组的性能比第1组都有所提高，但计算时间也都更长，且时间增大与Markov链长度大致成正比。这也意味着如果通过增大循环迭代次数改善性能，将会直接付出计算量成比例增大的代价。各种方案在相同参数下运行时间随问题规模有所增大，但差异并不显著，这是由于算法中仅计算路径差，而不需要计算完整的路径长度。因此，计算量与问题规模没有太大关系。

在不同规模问题的测试中，联合算子模拟退火算法优化性能提高的幅度不同。较小规模的问题，现有算子也已经获得较好的结果，联合算子提高性能的幅度不大；较大规模的问题，复杂程度高，现有算子的优化结果较差，而联合算子提高性能的幅度较为显著。

联合算子模拟退火算法在Eil51问题找到了已知的最优路径，在Eil76问题很接近已知的最优路径，但对Eil101、Ch150问题的优化结果与已知的最优路径还有2%～3%的误差。这一方面说明大规模TSP问题的复杂程度更高，找到最优路径需要更多的搜索；另一方面通过调整联合算子的参数，还可以获得更好的优化路径。

在4个不同问题和两组不同参数条件下，使用联合算子的计算时间比使用其他算子时略有增大，但并不显著，这与前文的分析是一致的。而联合算子产生的新解数量是其他算子的3倍，远高于计算时间的增大。

通过以上分析，联合算子模拟退火算法不增加太多的计算量，而通过产生多个新解以增强搜索能力从而提高了算法的优化性能。

4 结语

分析了模拟退火算法操作算子的特征与相互关系，发现交换操作等价于两个嵌套的反序操作的叠加复合。提出了一种新的交换-反序联合算子，获得由交换操作和反序操作所产生的3条新路径，从中择优作为新解。联合算子既在循环中产生多个新解，又不增加过多计算量，提高了算法的优化性能。通过4个不同规模的TSP问题的实验，验证了联合算子的性能优于现有的移位、交换、反序算子，也优于这些算子的组合方案。交换-反序联合算子不仅可以用于模拟退火算法，也可以应用于其他进化算法，如遗传算法、粒子群算法中的变异操作以扩展搜索空间。此外，交换与反序算子存在等价关系，其他操作算子之间也具有不同程度的相关性。对于揭示和利用这些相关性提高算法效率，将是后续研究工作的方向。