【建模算法】基于蚁群算法（ACA）求解TSP问题（Python实现）
TSP (traveling salesman problem，旅行商问题)是典型的NP完全问题，即其最坏情况下的时间复杂度随着问题规模的增大按指数方式增长，到目前为止还未找到一个多项式时间的有效算法。本文探讨了基于蚁群算法求解TSP问题的Python实现。
一、问题描述
本案例以31个城市为例，假定31个城市的位置坐标如表1所列。寻找出一条最短的遍历31个城市的路径。
| 城市编号 | X坐标 | Y坐标 | 城市编号 | X坐标 | Y坐标 |
|---|---|---|---|---|---|
| 1 | 1.304 | 2.312 | 17 | 3.918 | 2.179 |
| 2 | 3.639 | 1.315 | 18 | 4.061 | 2.37 |
| 3 | 4.177 | 2.244 | 19 | 3.78 | 2.212 |
| 4 | 3.712 | 1.399 | 20 | 3.676 | 2.578 |
| 5 | 3.488 | 1.535 | 21 | 4.029 | 2.838 |
| 6 | 3.326 | 1.556 | 22 | 4.263 | 2.931 |
| 7 | 3.238 | 1.229 | 23 | 3.429 | 1.908 |
| 8 | 4.196 | 1.044 | 24 | 3.507 | 2.376 |
| 9 | 4.312 | 0.79 | 25 | 3.394 | 2.643 |
| 10 | 4.386 | 0.57 | 26 | 3.439 | 3.201 |
| 11 | 3.007 | 1.97 | 27 | 2.935 | 3.24 |
| 12 | 2.562 | 1.756 | 28 | 3.14 | 3.55 |
| 13 | 2.788 | 1.491 | 29 | 2.545 | 2.357 |
| 14 | 2.381 | 1.676 | 30 | 2.778 | 2.826 |
| 15 | 1.332 | 0.695 | 31 | 2.37 | 2.975 |
| 16 | 3.715 | 1.678 |
二、蚁群算法简介
2.1 蚁群算法基本原理
1、蚂蚁在行走过程中会依据信息素来选择道路，选择信息素较浓的路走，并且在行走的路径中会释放信息素，对于所有蚂蚁都没经过的路，则随机选择一条路走；
2、蚂蚁释放的信息素浓度与长度相关，通常与路径长度成反比；
3、信息素浓的路径会受到蚂蚁更大概率的选择，形成正向反馈，最短路径上的信息素浓度会越来越大，最终蚁群就都按这条最短路径走。
信息素计算公式、转移概率、信息素重要程度因子、启发函数重要程度因子、信息素挥发因子等详细介绍可参考"蚁群算法详细讲解"和“TSP解决之道——蚁群算法”
2.2 算法的两个关键步骤
1、选择：为蚂蚁选择下一个城市，信息素越多的路径被选中概率较大，可用轮盘赌算法实现；
2、信息素更新：一段时间后（蚂蚁走完一个城市或者走完整个路径后）重新计算信息素（计算方法：历史累计信息素-信息素挥发量+蚂蚁行走释放量），蚂蚁行走释放量的常见方法有三种：蚁周算法（ant-cycle，蚂蚁走完整个路径后，蚂蚁行走释放部分用Q/L计算，Q表示蚂蚁释放信息素的量，为常量，L表示路径总长度）、蚁密算法（ant-density，蚂蚁走完一个城市后，蚂蚁行走释放用Q表示）、蚁量算法(ant-quantity，蚂蚁走完一个城市后，蚂蚁行走释放用Q/dij表示，dij表示城市i和j之间的距离)。
三、蚁群算法设计
在本算法设计中，包含两层循环，外循环是迭代次数循环，内循环遍历每一只蚂蚁，其中信息素增量在每一只蚂蚁走完整体路径后对当前信息素更新，而信息素挥发只在每一代蚂蚁都走完后对当前信息素更新，流程如下：
四、求解结果
最优路线与最优值：
最优轨迹图：
五、Python源代码
#蚁群算法求解31座城市TSP问题完整代码：
import numpy as np
import pandas as pd
import math
import matplotlib.pyplot as plt
from time import perf_counter
class AntList(object):
def __init__(self,distfunc,getEtatable,numant=5,numcity=10,alpha=1,rho=0.1,Q=1):
""" 构造函数 """
self.numant = numant # 蚂蚁个数
self.numcity = numcity # 城市个数
self.alpha = alpha # 信息素重要程度因子
self.rho = rho # 信息素的挥发速度
self.Q = Q # 品质因子
self.distfunc=distfunc
self.getEtatable=getEtatable
self.bestantunit=None
self.population=[]
self.pathtable = np.zeros((self.numant,self.numcity)).astype(int) # 路径记录表
self.generation=0
def Init_eta_phe(self):
"""
函数名：Init_eta_phe(self)
函数功能： 对启发函数和信息素进行初始化
"""
self.etatable = self.getEtatable() # 启发函数矩阵，表示蚂蚁从城市i转移到矩阵j的期望程度
self.pheromonetable = np.ones((self.numcity,self.numcity)) # 信息素矩阵
def InitStartPosition(self):
"""
函数名：InitStartPosition(self)
函数功能： 初始化蚂蚁的起始位置
"""
# 随机产生各个蚂蚁的起点城市
if self.numant <= self.numcity: # 城市数比蚂蚁数多
self.pathtable[:,0] = np.random.permutation(range(0,self.numcity))[:self.numant]
else: # 蚂蚁数比城市数多，需要补足
self.pathtable[:self.numcity,0] = np.random.permutation(range(0,self.numcity))[:]
self.pathtable[self.numcity:,0] = np.random.permutation(range(0,self.numcity))[:self.numant-self.numcity]
def upDateInf(self):
"""
函数名：upDateInf(self)
函数功能： 对信息素进行更新
"""
changepheromonetable = np.zeros((self.numcity,self.numcity))
if self.population:
for antunit in self.population:
for i in range(self.numcity-1):
changepheromonetable[antunit.path[i]][antunit.path[i+1]] += self.Q/antunit.length
changepheromonetable[antunit.path[self.numcity-1]][antunit.path[0]] += self.Q/antunit.length
self.pheromonetable = (1-self.rho)*self.pheromonetable + changepheromonetable
else:
self.Init_eta_phe()
def getNextCity(self,unvisited,visiting):
"""
函数名：getNextCity(self,unvisited,visiting)
函数功能： 根据信息素和启发函数矩阵，通过轮盘赌法随机选下一个城市
输入 1 self：类自身
输入 2 unvisited：未走过的城市列表
输入 2 visited：已经走过的城市列表
输出 1 k：下一个城市的编号
其他说明：无
"""
listunvisited = list(unvisited)
probtrans = np.zeros(len(listunvisited))
for k in range(len(listunvisited)):
probtrans[k] = np.power(self.pheromonetable[visiting][listunvisited[k]],self.alpha)\
*np.power(self.etatable[visiting][listunvisited[k]],self.alpha)
cumsumprobtrans = (probtrans/sum(probtrans)).cumsum()
cumsumprobtrans -= np.random.rand()
k = listunvisited[np.where(cumsumprobtrans>0)[0][0]] # 下一个要访问的城市
return k
def GoOnePath(self,i):
"""
函数名：distance(self, path)
函数功能： 第i只蚂蚁从随机点出发找到一条路径
输入 1 self：类自身
输入 2 i：当代的第i只蚂蚁
输出 1 antunit：一个蚂蚁单元类
其他说明：无
"""
visiting = self.pathtable[i,0] # 当前所在的城市
unvisited = set(range(self.numcity))# 未访问的城市
unvisited.remove(visiting) # 删除元素
for j in range(1,self.numcity): # 循环numcity-1次，访问剩余的numcity-1个城市
# 每次用轮盘法选择下一个要访问的城市
k=self.getNextCity(unvisited,visiting)
self.pathtable[i,j] = k
unvisited.remove(k)
visiting = k
antunit=AntUnit(self.pathtable[i],self.distfunc(self.pathtable[i]))
if self.bestantunit:
if self.bestantunit.length>antunit.length:
self.bestantunit=antunit
else:
self.bestantunit=antunit
return antunit
def nextGeneration(self):
"""
函数名：nextGeneration(self)
函数功能： 产生下一代
"""
self.upDateInf()
newPopulation = [] # 新种群
for i in range(self.numant):
newPopulation.append(self.GoOnePath(i))
self.population = newPopulation
self.generation += 1
class AntUnit(object):
"""
类名：GAUnit
类说明： 遗传算法个体类
"""
def __init__(self, aPath = None,aLength = -1):
""" 构造函数 """
self.path = list(aPath) # 个体的基因序列
self.length = aLength # 初始化适配值
class Node:
"""
类名：Node
类说明： 城市节点类
"""
def __init__(self,CityNum):
"""
函数名：GetData()
函数功能： 从外界读取城市数据并处理
输入 无
输出 1 Position：各个城市的位置矩阵
2 CityNum：城市数量
3 Dist：城市间距离矩阵
其他说明：无
"""
self.visited=[False]*CityNum #记录城市是否走过
self.start=0 #起点城市
self.end=0 #目标城市
self.current=0 #当前所处城市
self.num=0 #走过的城市数量
self.pathsum=0 #走过的总路程
self.lb=0 #当前结点的下界
self.listc=[] #记录依次走过的城市
def GetData(datapath):
"""
函数名：GetData()
函数功能： 从外界读取城市数据并处理
输入 无
输出 1 Position：各个城市的位置矩阵
2 CityNum：城市数量
3 Dist：城市间距离矩阵
其他说明：无
"""
dataframe = pd.read_csv(datapath,sep=" ",header=None)
Cities = dataframe.iloc[:,1:3]
Position= np.array(Cities) #从城市A到B的距离矩阵
CityNum=Position.shape[0] #CityNum:代表城市数量
Dist = np.zeros((CityNum,CityNum)) #Dist(i,j)：城市i与城市j间的距离
#计算距离矩阵
for i in range(CityNum):
for j in range(CityNum):
if i==j:
Dist[i,j] = math.inf
else:
Dist[i,j] = math.sqrt(np.sum((Position[i,:]-Position[j,:])**2))
return Position,CityNum,Dist
def ResultShow(Min_Path,BestPath,CityNum,string):
"""
函数名：GetData()
函数功能： 从外界读取城市数据并处理
输入 无
输出 1 Position：各个城市的位置矩阵
2 CityNum：城市数量
3 Dist：城市间距离矩阵
其他说明：无
"""
print("基于"+string+"求得的最优路线：")
for m in range(CityNum):
print(str(BestPath[m])+"—>",end="")
print(BestPath[CityNum])
print("最优值："+str(Min_Path))
print()
def draw(BestPath,Position,title):
"""
函数名：draw(BestPath,Position,title)
函数功能： 通过最优路径将旅行商依次经过的城市在图表上绘制出来
输入 1 BestPath：最优路径
2 Position：各个城市的位置矩阵
3 title:图表的标题
输出 无
其他说明：无
"""
plt.rcParams['font.sans-serif'] = 'SimHei' # 设置中文显示
plt.rcParams['axes.unicode_minus'] = False
plt.title(title)
plt.plot(Position[BestPath, 0],Position[BestPath, 1], marker='>', mec='r', mfc='w',label=u'路线')
plt.legend() # 让图例生效
for i,city in enumerate(Position):
plt.text(city[0], city[1], str(i))
plt.xlabel('横坐标')
plt.ylabel('纵坐标')
plt.show()
def ant():
"""
函数名：ant()
函数功能：蚁群算法核心
"""
numant = 25 # 蚂蚁个数
numcity = CityNum # 城市个数
alpha = 1 # 信息素重要程度因子
rho = 0.1 # 信息素的挥发速度
Q = 1
iters = 0
itermax = 500
etatable = 1.0/(Dist+np.diag([1e10]*numcity)) # 启发函数矩阵，表示蚂蚁从城市i转移到矩阵j的期望程度
pheromonetable = np.ones((numcity,numcity)) # 信息素矩阵
pathtable = np.zeros((numant,numcity)).astype(int) # 路径记录表
lengthaver = np.zeros(itermax) # 各代路径的平均长度
lengthbest = np.zeros(itermax) # 各代及其之前遇到的最佳路径长度
pathbest = np.zeros((itermax,numcity)) # 各代及其之前遇到的最佳路径
while iters < itermax:
# 随机产生各个蚂蚁的起点城市
if numant <= numcity: # 城市数比蚂蚁数多
pathtable[:,0] = np.random.permutation(range(0,numcity))[:numant]
else: # 蚂蚁数比城市数多，需要补足
pathtable[:numcity,0] = np.random.permutation(range(0,numcity))[:]
pathtable[numcity:,0] = np.random.permutation(range(0,numcity))[:numant-numcity]
length = np.zeros(numant) # 计算各个蚂蚁的路径距离
for i in range(numant):
visiting = pathtable[i,0] # 当前所在的城市
unvisited = set(range(numcity)) # 未访问的城市
unvisited.remove(visiting) # 删除元素
for j in range(1,numcity): # 循环numcity-1次，访问剩余的numcity-1个城市
# 每次用轮盘法选择下一个要访问的城市
listunvisited = list(unvisited)
probtrans = np.zeros(len(listunvisited))
for k in range(len(listunvisited)):
probtrans[k] = np.power(pheromonetable[visiting][listunvisited[k]],alpha)\
*np.power(etatable[visiting][listunvisited[k]],alpha)
cumsumprobtrans = (probtrans/sum(probtrans)).cumsum()
cumsumprobtrans -= np.random.rand()
k = listunvisited[np.where(cumsumprobtrans>0)[0][0]] # 下一个要访问的城市
pathtable[i,j] = k
unvisited.remove(k)
length[i] += Dist[visiting][k]
visiting = k
length[i] += Dist[visiting][pathtable[i,0]] # 蚂蚁的路径距离包括最后一个城市和第一个城市的距离
# 包含所有蚂蚁的一个迭代结束后，统计本次迭代的若干统计参数
lengthaver[iters] = length.mean()
if iters == 0:
lengthbest[iters] = length.min()
pathbest[iters] = pathtable[length.argmin()].copy()
else:
if length.min() > lengthbest[iters-1]:
lengthbest[iters] = lengthbest[iters-1]
pathbest[iters] = pathbest[iters-1].copy()
else:
lengthbest[iters] = length.min()
pathbest[iters] = pathtable[length.argmin()].copy()
# 更新信息素
changepheromonetable = np.zeros((numcity,numcity))
for i in range(numant):
for j in range(numcity-1):
changepheromonetable[pathtable[i,j]][pathtable[i,j+1]] += Q/length[i]
changepheromonetable[pathtable[i,j+1]][pathtable[i,0]] += Q/length[i]
pheromonetable = (1-rho)*pheromonetable + changepheromonetable
# 迭代次数指示器+1
iters += 1
path_tmp=pathbest[-1]
BestPath=[]
for i in path_tmp:
BestPath.append(int(i))
BestPath.append(BestPath[0])
return BestPath,lengthbest[-1]
#########程序入口#########################################
if __name__ == "__main__":
Position,CityNum,Dist = GetData("data.txt")
start=perf_counter() #计时开始
BestPath,Min_Path = ant()
print("运行时间是: {:.5f}s".format(perf_counter()-start)) # 程序计时结束
print()
ResultShow(Min_Path,BestPath,CityNum,"蚁群算法")
draw(BestPath,Position,"轨迹图")
扫一扫
6182
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
