Advertisement

python模拟退火求解TSP问题

阅读量:

python模拟退火求解TSP问题

  • 问题描述

  • 思路

    • 模拟退火简介
    • TSP问题
  • 代码

  • 运行结果

问题描述

通过模拟退火算法解决30个城市组成的旅行商问题。通过直角坐标系中的两点距离公式计算两城市之间的距离,并使用以下具体的城市位置数据进行求解:
{41,94}; {37,84}; {54,67}; {25,62}; {7,64}; {2,99}; {68,58}; {71,44}; {54,62}; {83,69}
{64,60}; {18,54}; {22,60}; {83,46}; {91,38}; {25,38}; {24,42}; {58,69}; {71,71}; {74,78}
{87,76}; {18,40}; {13,40}; {82,7}; {62,32}; {58,35}; [此处应包含完整的坐标数据]

思路

模拟退火简介

  • (1)选取任意一个单元k,并赋予其一个随机的移动量q_k^*。
    • (2)计算由此产生的能量变化量ΔE_k = E(q_k^) - E(q_k),其中E表示系统的总能量。
      若ΔE_k ≤ 0,则此移动方案可行,并将调整后的新状态s' = q_k^
      + q_k作为下一步调整的基础状态;
      若ΔE_k > 0,则新状态s' = q_k^* + q_k被接受的概率为
运行结果

在式中,T代表温度.接着从服从均匀分布U(0,1)的概率空间中随机选取一个值RR.如果生成的概率值R小于Pk,则采用变换后的状态作为下一次迭代的基础;否则,则采用变换前的状态.

  • (3)转第(1)步继续执行,知道达到平衡状态为止。

TSP问题

模拟退火算法的核心机制在于,在特定的概率P下采用较差的状态作为转移状态的依据。这种情况下,“较好”的状态会影响转移的概率值。“较差”的状态对应较高的转移概率,“较好”的状态则对应较低的转移概率。

  1. 初始化,模拟退火中用温度来观察程序进程,这里模拟温度,设置t0=100,从100°开始降温,程序停止条件t0==1,温度下降的速率为alpha=0.9,在每个温度状态下,程序要尝试随机生成路径TIME次,TIME=1000
  2. 创建城市坐标矩阵pointspoint[i][0]表示第i个城市的x坐标,point[i][1]表示第i个城市的y坐标
  3. 先随机生成一条路径way,通过函数calWayDis(way)计算路径总长度waydis
  4. 再随机选择2个城市city1city2,交换way中这两个城市的位置,形成新的路径tempway
复制代码
    tempway[city1],tempway[city2]=tempway[city2],tempway[city1]
  1. 该函数调用tempway以计算其路径长度,并将其存储为tempdis
  2. 比较:
  3. 如果:
    新的路径较短,则取代当前路径。
    并将其与当前最短的路径进行比较;
    如果其较短,则判断是否小于现有最短路径。
复制代码
    if tempdis<waydis:
    	way = tempway.copy()
    	waydis = tempdis
    if tempdis<bestdis:
    	bestway = tempway.copy()
    	bestdis = tempdis
    	draw(bestway,bestdis)

2. 如果:
新的路径长度不小于 现有路径长度,则根据模拟退火的概率,去接受

复制代码
    if tempdis>waydis:
    if np.random.rand()<np.exp(-(tempdis-waydis)/(t0)):
    way = tempway.copy()
    waydis = tempdis

步骤4至6在初始状态t_0执行;经过TIME次循环,在每一次循环中将当前路径与临时路径进行比较;随后更新温度参数t_0 = t_0 \times \alpha;随后继续重复步骤4至6;当t_0降至终止温度时,则停止操作并输出目前最优的路径。

复制代码
    while t0>t[0]:
    	步骤4
    	步骤5
    	步骤6
    	t0 *= alpha
  1. 打印出得到的最优解,路径和路径长度

代码

复制代码
    import numpy as np
    import matplotlib.pyplot as plt 
    import pdb
    import imageio
    import shutil
    import os 
    
    # 生成距离矩阵
    def getDismat(points):
    	dismat = np.zeros((N,N))
    	for i in range(N):
    		for j in range(i,N):
    			dismat[i][j] = dismat[j][i] = np.linalg.norm(points[i]-points[j])
    	return dismat	
    
    def init():
    	alpha = 0.9
    	t = (1,100)
    	TIME = 1000
    	way = np.arange(N)
    	waydis = calWayDis(way)
    	return alpha,t,TIME,way,waydis
    
    # 计算路径长度
    def calWayDis(way0):
    	waydis = 0
    	for i in range(N-1):
    		waydis +=dismat[way0[i]][way0[i+1]] 
    	waydis += dismat[way0[N-1]][way0[0]]
    	return waydis
    
    def draw(way,dist):
    	global N,points ,TIMESIT, PNGFILE, PNGLIST
    	plt.cla()
    	plt.title('cross=%.4f' % dist)
    	xs = [points[i][0] for i in range(N)]
    	ys = [points[i][1] for i in range(N)]
    	plt.scatter(xs, ys, color='b')
    	xs = np.array(xs)
    	ys = np.array(ys)
    	# plt.plot(xs[[0, solutionpath[0]]], ys[[0, solutionpath[0]]], color='r')
    	# 连接路径
    	for i in range(N-1):
    		plt.plot(xs[[way[i], way[i + 1]]], ys[[way[i], way[i + 1]]], color='r')
    	# 将终点与起点连接
    	plt.plot(xs[[way[N - 1], 0]], ys[[way[N - 1], 0]], color='r')
    	plt.scatter(xs[0], ys[0], color='k', linewidth=10)
    	for i, p in enumerate(points):
    		plt.text(*p, '%d' % i)
    	plt.savefig('%s/%d.png' % (PNGFILE, TIMESIT))
    	PNGLIST.append('%s/%d.png' % (PNGFILE, TIMESIT))
    	TIMESIT += 1
    
    points = np.array([[41,94],[37,84],[54,67],[25,62],[7,64],[2,99],[68,58],[71,44],[54,62],[83,69],
    [64,60],[18,54],[22,60],[83,46],[91,38],[25,38],[24,42],[58,69],[71,71],[74,78],
    [87,76],[18,40],[13,40],[82,7],[62,32],[58,35],[45,21],[41,26],[44,35],[4,50]
    ])
    N = points.shape[0]
    dismat = getDismat(points)
    alpha,t,TIME,way,waydis=init()
    t0 = t[1]
    K=0.8
    
    # 画图初试化
    TIMESIT = 0
    PNGFILE = './png/'
    PNGLIST = []
    if not os.path.exists(PNGFILE):
    os.mkdir(PNGFILE)
    else:
    shutil.rmtree(PNGFILE)
    os.mkdir(PNGFILE)
    # 记录每次迭代的结果
    result = []
    tempway = way.copy()
    bestway = way.copy()
    bestdis = 10000
    
    while t0>t[0]:
    	for i in range(TIME):
    		if np.random.rand() > 0.5:
    		# 两点交换
    			while True:
    			# 随机生成不同2个点,
    				city1 = np.int(np.ceil(np.random.rand()*(N-1)))
    				city2 = np.int(np.ceil(np.random.rand()*(N-1)))
    				if city1!=city2:
    					break
    			# 交换
    	tempway[city1],tempway[city2]=tempway[city2],tempway[city1]
    		else:
    			# 3个点
    			while True:
    				city1 = np.int(np.ceil(np.random.rand()*(N-1)))
    				city2 = np.int(np.ceil(np.random.rand()*(N-1))) 
    				city3 = np.int(np.ceil(np.random.rand()*(N-1)))
    				if((city1 != city2)&(city2 != city3)&(city1 != city3)):
    					break
    			# 下面的三个判断语句使得city1<city2<city3
    			if city1 > city2:
    				city1,city2 = city2,city1
    			if city2 > city3:
    				city2,city3 = city3,city2
    			if city1 > city2:
    				city1,city2 = city2,city1
    			#下面的三行代码将[city1,city2)区间的数据插入到city3之后
    			temp = tempway[city1:city2].copy()
    			tempway[city1:city3-city2+1+city1] = tempway[city2:city3+1].copy()
    			tempway[city3-city2+1+city1:city3+1] = temp.copy()
    
    		tempdis = calWayDis(tempway)
    		if tempdis<waydis:
    			way = tempway.copy()
    			waydis = tempdis
    		if tempdis<bestdis:
    			bestway = tempway.copy()
    			bestdis = tempdis
    			draw(bestway,bestdis)
    		else:
    			if np.random.rand()<np.exp(-(tempdis-waydis)/(t0)):
    				way = tempway.copy()
    				waydis = tempdis
    				# 更新路径
    			else: tempway = way.copy()
    
    	t0 *= alpha
    	result.append(bestdis)
    
    print("最短路径:%s"%np.array(result[-1]))
    for i in bestway:
    	print(i,end="-->")
    print(bestway[0])
    
    #生成路径gif 
    generated_images = []
    for png_path in PNGLIST:
    generated_images.append(imageio.imread(png_path))
    # shutil.rmtree(PNGFILE)  # 可删掉
    generated_images = generated_images + [generated_images[-1]] 
    imageio.mimsave('Sovle_TSP01.gif', generated_images, 'GIF', duration=0.1)

运行结果

运行结果
生成路径图

全部评论 (0)

还没有任何评论哟~