python蚁群算法解决TSP和CVRP

提示:文章写完后,目录可以自动生成,如何生成可参考右边的帮助文档


前言

蚁群算法,懂的自然懂,不懂的不会进来看


提示:以下是本篇文章正文内容,下面案例可供参考

一、TSP与CVRP的区别

TSP简单的说就是一个人以最短的路程走遍天下所有城市并回到出发城市不带重复的。
CVRP简单地说就是一个货车去各个城市接货,装满就拖车带回来,没装满就继续去写一个城市接货,直到所有城市的货都被他接收完,并且用的总路程最短。

二、使用步骤

1.引入库

import numpy as np
import matplotlib.pyplot as plt

2.直接上代码

代码如下(示例):

# -*- coding: utf-8 -*-

import os
import time
import math

os.getcwd()
#返回当前工作目录,用于保存文件到当前目录
import numpy as np
import matplotlib.pyplot as plt
#解决中文显示问题
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
""" 第一步:初始化参数-Start """



# 获取文件信息
def load_example():
    url = f'./data/xy.txt'
    file = open(url,'r')
    content = file.readlines()
    points =[]
    need=[]
    for data in content:
        n,x,y,l = data.split(" ")
        points.append([int(x),int(y)]) 
        need.append(int(l)) 
    distance = calculate_distance(points)
    points=np.array(points)
    return points,distance,need

def dis(point1,point2):
    return math.sqrt((point1[0]-point2[0]) ** 2 + (point1[1] - point2[1]) ** 2)

def calculate_distance(points):
    num = len(points)
    distMat = np.zeros([num,num])
    for i in range(len(points)):   
        for j in range(len(points)):
            distMat[i,j]= distMat[i][j] =np.linalg.norm(dis(points[i],points[j]))
    return distMat





points,distance,need = load_example()
distMat = distance  # 节点的距离矩阵
limit=20#配送车辆最大载重
coordinates= points
numAnt = int(len(coordinates)*1.5)  # 蚂蚁个数  节点数的1.5倍以内即可
points = len(coordinates)  #节点数
alpha = 1  # 信息素重要程度因子
beta = 1  # 启发函数重要程度因子
rho = 0.1  # 信息素的挥发速度
Q = 1  # 信息素常量
iter = 0  #迭代初始
iterMax = 150  #迭代总数
pheromoneTable = np.ones((points, points))  #信息素矩阵
pathTable = np.zeros((numAnt, points)).astype(int)  # 路径记录表,转化成整型
lengthAver = np.zeros(iterMax)  # 迭代数矩阵 ,存放每次迭代后,路径的平均长度,就是个长度为迭代数的数组
lengthBest = np.zeros(iterMax)  # 迭代数矩阵,存放每次迭代后,最佳路径长度
pathBest = np.zeros((iterMax, points))  # 迭代,节点数,存放每次迭代后,最佳路径城市的坐标
path = ''  #最优路径
pathText = np.zeros(iterMax)  #记录每次迭代产生的最短路径
thisLength = ''  #第t次迭代时当下产生的最短路径

""" 第一步:初始化参数-End """
""" 第二步:构建解空间-START """
""" 第三步:迭代-START """
time_start = time.perf_counter()




# 轮盘赌
def roulette(probability):
    probabilityTotal = np.zeros(len(probability))
    probabilityTmp = 0
    for i in range(len(probability)):
        probabilityTmp += probability[i]
        probabilityTotal[i] = probabilityTmp
    randomNumber = np.random.rand()
    result = 0
    for i in range(1, len(probabilityTotal)):
        if randomNumber < probabilityTotal[0]:
            result = 0
            # print("random number:",randomNumber,"<index 0:",probabilityTotal[0])
            break
        elif probabilityTotal[i - 1] < randomNumber <= probabilityTotal[i]:
            result = i
            # print("index ",i-1,":",probabilityTotal[i-1],"<random number:",randomNumber,"<index ",i,":",probabilityTotal[i])
    return result


while iter < iterMax:
       #m个蚂蚁随机放置于n个节点中,分次增加随机性
    pathTable[:points, 0] = np.random.permutation(range(points))[:]
    # 先放20个
    pathTable[points:, 0] = np.random.permutation(
        range(points))[:numAnt - points]  # 再把剩下的随机放完
    antcount = np.zeros(numAnt)
    length = np.zeros(numAnt)  # 1*蚂蚁个数的数组
    """ 计算k蚂蚁下一节点的概率-START """
    #本段程序算出每只/第i只蚂蚁转移到下一个城市的概率
    for k in range(numAnt):
        # i=0
        visiting = pathTable[k, 0]  # 当前所在的城市
        unVisited = set(range(points))  #未访问的城市集合
        unVisited.remove(visiting)  # 删除已经访问过的城市元素
        for j in range(1, points):  # 循环points-1次,访问剩余的所有points-1个城市
            # j=1
            # 每次用轮盘法选择下一个要访问的城市
            listUnisited = list(unVisited)  #未访问节点数,list
            probtrans = np.zeros(
                len(listUnisited))  #每次循环都初始化转移概率矩阵1*20,1*19,1*18,1*17....
            #以下是计算转移概率
            for p in range(len(listUnisited)):
                probtrans[p] = math.pow(pheromoneTable[visiting][listUnisited[p]], alpha) \
                               * math.pow(1/distMat[visiting][listUnisited[p]], beta ) #(ηij)eta-从城市i到城市j的启发因子 这是概率公式的分母   其中[visiting][listunvis[k]]是从本城市到k城市的信息素
            cumsumprobtrans=np.zeros(len(listUnisited))
            for l in range(len(listUnisited)):
                cumsumprobtrans[l] =  (probtrans[l] / sum(probtrans))

            # 以下使用轮盘赌,增加随机性
            p=listUnisited[roulette(cumsumprobtrans)]
            antcount+=antcount[k]+need[p]
            # if antcount[k]>limit:
                # print("超过最大载重限制")
            # 不干了,去他丫的CVRP


            pathTable[k, j] = p #采用禁忌表来记录蚂蚁i当前走过的第j城市的坐标,这里走了第j个城市.k是中间值
            unVisited.remove(p)  #将未访问城市列表中的K城市删去。
            length[k] += distMat[visiting][
                p]  #计算蚂蚁i从节点i到j的路径长度 d_ij并累加到之前的所有路径长度
            #计算本城市到K城市的距离
            visiting = p
        length[k] += distMat[visiting][pathTable[k, 0]]  #计算第i只蚂蚁所有路径长度和-L_i
    # print(pathTable)
    # print(f"L_k={length}")
    # 包含所有蚂蚁的一个迭代结束后,统计本次迭代的若干统计参数
    lengthAver[iter] = length.mean()
    thisLength = length.min()
    #本轮的平均路径
    """ 求出最佳路径 """
    # 对比本次迭代中的所有蚂蚁行走的路径长度,选择出最短的
    pathText[iter] = length.min()
    if iter == 0:
        lengthBest[iter] = length.min()
        pathBest[iter] = pathTable[length.argmin()].copy()
    #如果是第一轮路径,则选择本轮最短的路径,并返回索引值下标,并将其记录
    else:
        #后面几轮的情况,更新最佳路径
        # 如果此次迭代没有比以前更短的路径
        if length.min() > lengthBest[iter - 1]:
            lengthBest[iter] = lengthBest[iter - 1]
            pathBest[iter] = pathBest[iter - 1].copy()
        # 如果此次迭代有比以前更短的路径或者等于之前的最短路径,则选择本轮最短的路径,并返回索引值下标,并将其记录
        else:
            # print("禁忌表",pathTable)
            lengthBest[iter] = length.min()
            # print("输出信息素矩阵",pheromoneTable)
            pathBest[iter] = pathTable[length.argmin()].copy()
   

    #此部分是为了更新信息素-------∆τ_ij^k
    changepheromoneTable = np.zeros((points, points))
    for k in range(numAnt):  #更新所有的蚂蚁
        for j in range(points - 1):
            changepheromoneTable[pathTable[k, j]][pathTable[
                k, j + 1]] += Q / distMat[pathTable[k, j]][pathTable[k, j + 1]]
            #根据公式更新本只蚂蚁改变的节点间的信息素      Q/d_ij   其中d是从第j个城市到第j+1个城市的距离
        changepheromoneTable[pathTable[k, j + 1]][pathTable[
            k, 0]] += Q / distMat[pathTable[k, j + 1]][pathTable[k, 0]]
        #起点到终点 所有蚂蚁改变的信息素总和
        #∑_(k=1)^m *∆τ_ij^k
    #信息素更新公式p=(1-挥发速率)*现有信息素+改变的信息素
    pheromoneTable = (1 - rho) * pheromoneTable + changepheromoneTable
    # print('更新的信息素',pheromoneTable)





    #根据当前信息素更新下一时刻的信息素
    # τ_ij (t+1)=(1-ρ) τ_ij (t)+∆τ_ij
    print(f"迭代次数={iter+1}", f"本次迭代平均路径={lengthAver[iter]}",
          f"本次迭代最短的路径={thisLength}", f"历代迭代最短的路径={lengthBest[iter]}")
    iter += 1  # 迭代次数指示器+1
    # print('信息素',pheromoneTable)
    # print('禁忌表',pathTable)
#迭代完成----------------------------------------------------------------------------------------
time_end = time.perf_counter()
print('运行耗时', time_end - time_start, 's')
""" 第二步:构建解空间-end """
""" 第三步:迭代-end """
""" 输出最短路径 """
for i in pathBest[-1]:
    path += str(int(i)) + '>'
print("best_path=", path + str(int(pathBest[-1][0])),
      f'最短路径出现在第{list(pathText).index(pathText.min())+1}次迭代过程中',
      f'最短路径长度为={pathText.min()}')
#以下是做图部分
#做出平均路径长度和最优路径长度
fig, axes = plt.subplots(nrows=2, ncols=1, figsize=(12, 10))
axes[0].plot(lengthAver, 'k', marker='.')
axes[0].set_title('平均路径长度')
axes[0].set_xlabel(u'迭代次数')
#线条颜色black https://blog.csdn.net/ywjun0919/article/details/8692018
# axes[1].plot(pathText, 'k', marker='.')
axes[1].plot(lengthBest, 'k', marker='.')

axes[1].set_title('最短路径长度')
axes[1].set_xlabel(u'迭代次数')
fig.savefig('最短路径长度.png', dpi=500, bbox_inches='tight')
plt.show()
# 作出找到的最优路径图
bestpath = pathBest[-1]
plt.plot(coordinates[:, 0], coordinates[:, 1], 'ro')
plt.xlim([90, 200])
#x范围
plt.ylim([860, 950])
#y范围
for i in range(points - 1):
    #按坐标绘出最佳两两节点间路径
    m, n = int(bestpath[i]), int(bestpath[i + 1])
    plt.plot([coordinates[m][0], coordinates[n][0]],
             [coordinates[m][1], coordinates[n][1]], 'k')
plt.plot([
    coordinates[int(bestpath[0])][0], coordinates[int(
        bestpath[coordinates.shape[0] - 1])][0]
], [
    coordinates[int(bestpath[0])][1], coordinates[int(
        bestpath[coordinates.shape[0] - 1])][1]
], 'b')
ax = plt.gca()
ax.set_title("最优路径规划")
ax.set_xlabel('X-经度')
ax.set_ylabel('Y-纬度')
plt.savefig('最优路径规划.png', dpi=500, bbox_inches='tight')
plt.show()

总结

在写这个代码的时候,因为经常要测试,所以累了,现在发的这个是TSP的成品,CVRP没写完,但是对于TSP如果固定起点,就是说

while iter < iterMax:
       #m个蚂蚁随机放置于n个节点中,分次增加随机性
    pathTable[:, 0]=0
   使出发点永远是0节点的时候

就会出现老是产生局部最优解,可能是信息素问题没有优化好

  • wx:bybilibili
  • 如有疑问,问这…
  • 加油,同学

版权声明:本文为weixin_44700966原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接和本声明。