上一篇我们讲了旅行商问题及遗传算法来解决此类问题：遗传算法(附代码)
今天介绍另外一种解决此类NP问题的方法是模拟退火算法(Simulated Annealing, SA)
模拟退火算法的思想借鉴于固体的退火原理，当固体的温度很高的时候，内能比较大，固体的内部粒子处于快速无序运动，当温度慢慢降低的过程中，固体的内能减小，粒子的慢慢趋于有序，最终，当固体处于常温时，内能达到最小，此时，粒子最为稳定。模拟退火算法便是基于这样的原理设计而成。
那么为什么在算法开始的时候要处于不稳定的状态呢？我们先看来一下爬山法吧
爬山法的问题
爬山法是一种贪婪算法，在当前位置的附近寻找一个更大的值，不断迭代这个过程直到处于稳定状态
如图中要寻找这个函数的最大值，采用爬山法的话，如果初始点选择D点，那么爬山法能够找到最大值B点，但是如果初始值选择C或者E点，那么爬山法找到的最大值就是局部最优A点，而不是全局最优点B点
这就是爬山法的局限性，如果选择的初始点不好，那么算法很有可能会陷入局部最优解，而模拟退火引入的温度，在温度越高的时候越有机会跳到一个相对不那么优的解，从而跳出局部最优。
模拟退火的过程
- 1.初始化一个温度T, 如T=100
- 2.老的解损失为: old_cost, 生成一个新的解,损失为: new_cost
- 3.计算当前是否采纳新的解概率: P = math.exp(-(new_cost-old_cost)/T)
- 4.如果损失new_cost<old_cost或者随机一个(0,1)之间数值<P, 则接受这个新解，否则还是使用老的解
- 5.当前的温度T下降一定的值(温度降低)
- 6.一直迭代2~5步，直到稳定状态
可以看到一开始温度T很高的时候P也很大，那么采纳新的解(损失比老的解更大)的概率也更高，从而跳出了局部最优
随着温度T的下降P也越来越小，会越来越多的采用损失大小来评价是否采纳新的解
代码实现
# -*- encoding: utf8 -*-
import random
import math
import matplotlib.pyplot as plt
class SA:
def __init__(self):
self.T = 10000.
self.cool = 0.95
def new_route(self, tour, data):
c1 = random.randint(0, len(tour['tour']) - 1)
c2 = random.randint(0, len(tour['tour']) - 1)
while c1 == c2:
c2 = random.randint(0, len(tour['tour']) - 1)
new_tour = []
for i in range(len(tour['tour'])):
if i == c1:
new_tour.append(tour['tour'][c2])
elif i == c2:
new_tour.append(tour['tour'][c1])
else:
new_tour.append(tour['tour'][i])
return self.get_tour_detail(new_tour, data)
def get_distance(self, c1, c2):
# 获取距离
xd = abs(c1['x'] - c2['x'])
yd = abs(c1['y'] - c2['y'])
distance = math.sqrt(xd * xd + yd * yd)
return distance
def get_cost(self, distance):
return distance
def get_tour(self, data):
# 随机获得一条路径
tour = []
for key, value in data.items():
tour.append(key)
random.shuffle(tour)
return self.get_tour_detail(tour, data)
def get_tour_detail(self, tour, data):
tmp = None
distance = 0
for item in tour:
if tmp is not None:
distance_tmp = self.get_distance(data[tmp], data[item])
distance += distance_tmp
tmp = item
return {'tour': tour, 'distance': distance, 'cost': self.get_cost(distance)}
def run(self, data):
route = self.get_tour(data)
print 'before route:'
print route
i = 0
while self.T > 0.1:
new_route = self.new_route(route, data)
old_cost, new_cost = route['cost'], new_route['cost']
if new_cost < old_cost or random.random() < math.exp(-(new_cost - old_cost) / self.T):
route = new_route
self.T = self.T * self.cool
i += 1
print 'total gen:', i
print route
return route
def init_data(num):
data = {}
def get_random_point():
# 随机生成坐标
return {
'x': random.randint(1, 99),
'y': random.randint(1, 99)
}
for i in range(num):
data[i + 1] = get_random_point()
return data
def show(citys, tour):
# 绘图
plt.bar(left=0, height=100, width=100, color=(0, 0, 0, 0), edgecolor=(0, 0, 0, 0))
plt.title(u'tsp')
plt.xlabel('total distance: %s m' % tour['distance'])
x = []
y = []
i = 0
for item in tour['tour']:
city = citys[item]
x.append(city['x'])
y.append(city['y'])
i += 1
if i == 1:
plt.plot(city['x'], city['y'], 'or')
else:
plt.plot(city['x'], city['y'], 'bo')
plt.plot(x, y, 'g')
plt.xlim(0.0, 100)
plt.ylim(0.0, 100)
plt.show()
if __name__ == '__main__':
scale, num, max_gen, pc, elite = 50, 10, 1000, 0.8, 0.2
data = init_data(num)
sa = SA()
new_fittest = sa.run(data)
show(data, new_fittest)
运行效果：
参考
[1] <<集体智慧编程>>
[2] https://www.cnblogs.com/heaad/archive/2010/12/20/1911614.html
[3] https://zhuanlan.zhihu.com/p/33184423
