从理论到实践：用Python实现模拟退火算法求解经典TSP

1. 模拟退火算法与旅行商问题简介第一次听说模拟退火算法时我脑海中浮现的竟然是打铁匠淬火的场景。这种将物理现象转化为数学模型的思路正是算法设计的精妙之处。模拟退火算法(Simulated Annealing)最早由Kirkpatrick等人在1983年提出灵感来源于固体退火过程——金属加热到高温后缓慢冷却原子最终会排列成能量最低的稳定状态。旅行商问题(TSP)则是组合优化领域的经典难题给定一系列城市和它们之间的距离找到一条最短的路径让旅行商访问每个城市一次并最终回到起点。别看问题描述简单当城市数量增加到30个时可能的路径组合就超过了天文数字(30!/2≈4.4×10^30)。传统穷举法在这里完全失效而模拟退火这类启发式算法却能给出不错的近似解。为什么模拟退火适合解决TSP我总结三点关键优势跳出局部最优相比贪心算法容易陷入局部最优模拟退火通过概率性接受劣质解的机制有机会跳出局部陷阱参数可控通过调整温度参数可以灵活平衡探索(全局搜索)与开发(局部优化)的关系实现简单核心算法用几十行Python就能实现不需要复杂的数学工具2. 算法核心原理拆解2.1 物理退火与算法对应关系理解模拟退火最好的方式就是从物理原型出发。在金属退火过程中温度(T)算法中的核心参数控制接受劣解的概率能量(E)对应目标函数值TSP中就是路径总长度状态即当前解TSP中是一个城市访问序列算法运行时温度会从初始高温T0逐渐降低到终止温度Tf。高温时系统活跃容易接受能量更高的状态更长的路径随着温度降低系统趋于稳定只接受更好的或轻微变差的状态。2.2 Metropolis准则的实现细节这个看似简单的概率公式藏着算法的精髓deltaE new_energy - current_energy if deltaE 0: accept True # 新解更优直接接受 else: p math.exp(-deltaE / T) # 计算接受概率 accept random.random() p # 概率性接受劣解实际编码时我踩过一个坑温度T必须大于0否则math.exp会溢出。建议设置Tf1e-8这样的小数而非绝对0。2.3 邻域搜索算子设计TSP问题中如何从当前解产生新解尤为关键。常见有三种算子交换算子(Swap)def swap(tour): i, j random.sample(range(len(tour)), 2) new_tour tour.copy() new_tour[i], new_tour[j] new_tour[j], new_tour[i] return new_tour反序算子(2-opt)def reverse(tour): i, j sorted(random.sample(range(len(tour)), 2)) return tour[:i] tour[i:j1][::-1] tour[j1:]移位算子(Move)def move(tour): i, j random.sample(range(len(tour)), 2) city tour.pop(i) tour.insert(j, city) return tour实测发现对于50个城市以内的问题反序算子效果最好。因为它能同时改变多个边搜索效率更高。3. Python完整实现3.1 算法参数设置参数调优是门艺术经过多次实验我总结出这些经验值def init_parameters(): params { T0: 1000, # 初始温度(与问题规模相关) Tf: 1e-8, # 终止温度 alpha: 0.98, # 降温系数(0.95-0.99) Lk: 1000, # 每个温度的迭代次数 max_stay: 20 # 最优解持续未改进次数 } return params特别说明alpha的选择取值0.95时降温快适合简单问题快速求解取值0.99时降温慢适合复杂问题精细搜索我习惯从0.95开始逐步调高观察效果3.2 核心算法框架完整算法流程如下def simulated_annealing(cities): # 初始化 current_tour random_permutation(cities) current_energy calc_distance(current_tour) best_tour, best_energy current_tour.copy(), current_energy params init_parameters() T, stay_counter params[T0], 0 while T params[Tf] and stay_counter params[max_stay]: for _ in range(params[Lk]): # 产生新解 new_tour get_neighbor(current_tour) new_energy calc_distance(new_tour) # Metropolis准则 if accept_new_solution(current_energy, new_energy, T): current_tour, current_energy new_tour, new_energy # 更新最优解 if new_energy best_energy: best_tour, best_energy new_tour.copy(), new_energy stay_counter 0 # 降温 T * params[alpha] stay_counter 1 return best_tour, best_energy3.3 可视化与性能分析用matplotlib绘制优化过程能直观理解算法行为def plot_optimization(record_best, record_now): plt.figure(figsize(10,5)) plt.plot(record_best, r-, labelBest Solution) plt.plot(record_now, b--, labelCurrent Solution) plt.xlabel(Iterations) plt.ylabel(Tour Length) plt.legend() plt.title(SA Optimization Process) plt.show()典型优化曲线会呈现三个阶段高温震荡期解的质量波动剧烈算法广泛探索中温过渡期逐渐收敛偶尔接受劣解低温稳定期基本只接受优化最终稳定4. 实战技巧与性能优化4.1 初始解生成策略完全随机初始解虽然简单但可能让算法需要更长时间收敛。我常用的改进方法def greedy_init(cities): tour [random.choice(range(len(cities)))] unvisited set(range(len(cities))) - {tour[0]} while unvisited: last tour[-1] next_city min(unvisited, keylambda x: distance(cities[last], cities[x])) tour.append(next_city) unvisited.remove(next_city) return tour这种贪心初始化能让算法起步更快平均减少20%的迭代时间。4.2 自适应降温策略固定降温系数有时效率不高我实现了动态调整方案def adaptive_cooling(T, improvement_ratio): if improvement_ratio 0.1: # 改进明显加快搜索 return T * 0.95 elif improvement_ratio 0.01: # 改进缓慢深入挖掘 return T * 0.99 else: # 正常降温 return T * 0.974.3 并行化加速技巧对于大规模TSP问题可以结合多进程加速from multiprocessing import Pool def parallel_SA(args): # 每个进程独立运行SA return simulated_annealing(*args) if __name__ __main__: with Pool(4) as p: results p.map(parallel_SA, [(cities,)]*4) best min(results, keylambda x: x[1])这种实现能在4核CPU上获得接近线性的加速比特别适合需要多次运行取最优解的场景。