【建模算法】基于蚁群算法求解TSP问题（Python实现）

TSP (traveling salesman problem，旅行商问题)是典型的NP完全问题，即其最坏情况下的时间复杂度随着问题规模的增大按指数方式增长，到目前为止还未找到一个多项式时间的有效算法。本文探讨了基于蚁群算法求解TSP问题的Python实现。

一、问题描述

​ 本案例以31个城市为例，假定31个城市的位置坐标如表1所列。寻找出一条最短的遍历31个城市的路径。

二、蚁群算法简介

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,"轨迹图")