最大最小蚁群算法求解TSP商旅问题

---

最大最小蚁群算法求解TSP商旅问题

---

算法要点

1. 每条路径上的信息素浓度都有一个最大最小值，最小信息素增加对于更优解探索的可能性，最大信息素浓度保证经验对于蚁群的启发性。

2. 对于蚂蚁的下一次选择使用轮盘赌的方式进行选择，每条路径的权重是根据启发公式进行计算，合适的alpha beta能够加速算法的收敛，这是经验性的参数

3. 蚂蚁的数量一般不宜过多，一般与城市数量差不多就可以了。经过多次迭代可以发现，当前的最优可行解就是信息素最大的那条路径。

4. 挥发系数的确定：在每轮迭代过程中，每一次迭代完毕后，下一次的路径的信息素的浓度会随着时间挥发。smell_new=smell_old*挥发系数+smell_added。在计算的同时要保证 smell_new在[smell_lowbound,smell_upbound]范围内，即信息素的浓度有边界（参见第一条）。

---

结论：以上条件能够使得算法具有很快的收敛速度，一般对于50个点的TSP问题，在200轮迭代后就能达到最优解，在图中可以看到，TSP的最优可行解是一个圈，不会 产生交叉。

---

代码：

'''蚁群算法求解商旅问题author: sunqinotice：使用前请先安装ffmpeg库algrithm: MMAS'''from matplotlib import pyplot as pltfrom matplotlib.animation import FuncAnimationimport numpy as npfrom matplotlib import collections as mcfrom matplotlib.ticker import  MultipleLocatorpoint_num = 30  # 商旅问题的点数rou = 0.05  # 蚁群算法中的蒸发系数alpha = 0.8  # 有关于信息素的权重系数beta = 2  # 有关于距离的权重系数MaxIterTimes = 10000;  # 允许的最大迭代次数Q = 15  # 单个蚂蚁所携带的信息素的总量，用于走完以后更新每条走过的边的信息素 delta=Q/L L:路径长度anti_num = 50;  # 蚂蚁的数量tao_min = 0.1/point_numtao_max = 1iter_count = 0  # 计数器计算迭代次数def generate_point(num):    '''    产生随机的点x,y坐标 个数num  坐标限制（0~100, 0~100）    '''    # 产生随机数，xy坐标    x, y = 100 * np.random.rand(2, num)    # ax1.scatter(x, y, marker='o')    # ax2.scatter(x, y)    return x, ydef point_matrix(x, y):    '''    处理产生的x,y坐标向量    返回:    point点对Array shape: pointnum,2    distance_matrix : distance from point_i to point_j    '''    distance_matrix = np.zeros((point_num, point_num))    point_mat = np.c_[x, y]    # print(point_mat.shape)    for i in range(point_num):        for j in range(i + 1, point_num):            distance_matrix[i, j] = np.linalg.norm(point_mat[i] - point_mat[j])            distance_matrix[j, i] = distance_matrix[i, j]    return point_mat, distance_matrixdef generate_smell_matrix():    '''    产生最开始的信息素矩阵    '''    smell_matrix = np.zeros((point_num, point_num))    # for i in range(point_num):    #     for j in range(point_num):    #         # 节点i to i之间是没有信息素的    #         smell_matrix1[i, j] = tao_max if i != j else 0    smell_matrix[::]=tao_max    for i in range(point_num):        smell_matrix[i,i]=0    return smell_matrixdef anti_choice(forbidList):    '''    蚂蚁在当前城市和已知禁忌表的情况下，下一步的去向    返回：nextPointIndex    '''    # 当前所在的城市    index = forbidList[-1]    # 获得可能的路径    possible_ends = []    for i in range(point_num):        if (i not in forbidList):            possible_ends.append(i)    # 计算概率    possibility = []    for j in possible_ends:        # 概率计算公式        p = smell_matrix[index, j] ** alpha /distance_matrix[index, j] ** beta        possibility.append(p)    # 计算归一化的概率    possibility = np.array(possibility) / sum(possibility)    # 以指定的概率进行选择路径    p = np.random.rand()    k = 0    for i in range(possibility.__len__()):        k += possibility[i]        if (p < k):            return possible_ends[i]    return possible_ends[-1]def calDistance(route1):    '''    计算蚂蚁走过的长度    '''    # 拷贝 不改变原来的参数    route = route1.copy()    route.append(route[0])    distance = 0    for i in range(len(route) - 1):        distance += distance_matrix[route[i], route[i + 1]]    return distance# 产生信息素矩阵标记每条路径 point_i to point_j的含有的信息素def update_smell_matrix(route):    '''    使得每条路径上的信息素浓度得到一定程度的衰减，同时增大最优路径上每条路径上的信息素    '''    global smell_matrix    smell_matrix = (1 - rou) * smell_matrix    newroute = route.copy()    newroute.append(newroute[0])    for i in range(len(newroute) - 1):        smell_matrix[newroute[i], newroute[i + 1]] += Q / distance_matrix[newroute[i],newroute[i + 1]]        #使对角线r_ij=r_ji 信息素含量相同        smell_matrix[newroute[i+1], newroute[i]]=smell_matrix[newroute[i], newroute[i + 1]]    for i in range(point_num):        for j in range(i+1,point_num):            if (i != j):                if (smell_matrix[i, j] < tao_min):                    smell_matrix[i, j] = tao_min                if (smell_matrix[i, j] > tao_max):                    smell_matrix[i, j] = tao_max            smell_matrix[j, i]=smell_matrix[i, j]    # 每条边的信息素的浓度都在tao_min~tao_max范围内    # smell_matrix[smell_matrix<tao_min]=tao_min    # smell_matrix[smell_matrix>tao_max]=tao_maxdef MMAS():    # 各个蚂蚁被随机分配的第一个城市    first_citys = np.random.randint(point_num, size=anti_num)    # 用于存放蚂蚁们的禁忌表,将第一个城市放入禁忌表中    anti_list = [[first_citys[i]] for i in range(anti_num)]    # 存储每一个蚂蚁走过的的长度    anti_distance_list = np.zeros((anti_num))    # 对每个蚂蚁进行迭代    for i in range(anti_num):        # 获取当前蚂蚁的禁忌表        forbid_list = anti_list[i]        # 每个蚂蚁都需要访问城市        for j in range(1,point_num):            # 每个蚂蚁在每个城市节点都面临选择            nextPoint = anti_choice(forbid_list)            forbid_list.append(nextPoint)  # 将选择的路径添加进禁忌表        anti_distance_list[i] = calDistance(forbid_list)    # 找出了当前最短的路径序号    minIndex = anti_distance_list.argmin()    # 获得全局最小长度    goodRoute=False   #如果此次的过程比历史最优解好就是一个好的选择    global minDistance,bestRoute,iter_count    if minDistance==None:        minDistance = anti_distance_list[minIndex]        bestRoute=anti_list[minIndex]    else:        dis=anti_distance_list[minIndex]        #只有该轮找到全局最小值才能对bestRoute进行更新        if dis<=minDistance:            minDistance=dis            bestRoute=anti_list[minIndex]    # 使用历史最优路径更新全局信息素的浓度    update_smell_matrix(bestRoute)    # 全局最优路径    iter_count+=1    #用一个列表记录time,distance 信息    time_list.append(iter_count)    distance_list.append(minDistance)    # print(bestRoute,minDistance,'\n')# 用于绘制第一个子图的路径曲线def routeIndex2mat(route1):    '''    :param route1:     :return: 返回坐标二维数组,将所有点围成一个回路    '''    route = route1.copy()    route.append(route[0])    mat = np.zeros((point_num + 1, 2))    k = 0    for i in route:        mat[k] = point_mat[i]        k += 1    # ax1.set_title("城市点数:{} 距离:{:.2f} 迭代次数:{}".format(point_num,distance_list[-1],iter_count))    ax1.set_title("Citys:{} Distance:{:.2f} Iter_times:{}".format(point_num,distance_list[-1],iter_count))    return mat# 绘制各条路径上的信息素浓度，用于第二张子图的绘制def LinesMatrix():    ax2.clear()    ax2.scatter(x,y)    for i in range(point_num):        for j in range(i+1,point_num):            value=smell_matrix[i,j]            if j!=None:                xdata=[point_mat[i,0],point_mat[j,0]]                ydata=[point_mat[i,1],point_mat[j,1]]                line=plt.Line2D(xdata,ydata,color='r',linewidth=value*4)                ax2.add_line(line)#更新ax3上的数据def update_Axes3():    limit=ax3.dataLim._get_bounds()    xlim=0    # if iter_count>1000:    #     x=iter_count-1000    #     ax3_line.set_data(time_list[-1000:],distance_list[-1000:])    # else:    ax3_line.set_data(time_list,distance_list)    ax3.set_xlim((0,iter_count*1.1))    ax3.set_ylim((0,distance_list[0]*1.1))    ax3.autoscale_view()def animate(i):    # 蚁群算法进行20次迭代 显示一次资源    for i in range(point_num):        MMAS()    #更新第一幅图内的最优路径    line=routeIndex2mat(bestRoute)    lc1.set_paths([line])    #更新第二幅图    LinesMatrix()    #更新第三个子图上的信息    update_Axes3()    # return ax1,ax2,ax3    return lc1,ax2,ax3_linedef init():    global lc1    # lc1 = mc.LineCollection([])    # ax1.add_collection(lc1)    # lc1.set_paths([[[10,10],[20,20]]])    lc1.set_paths([])    global ax3_line    ax3_line, = ax3.plot([], [])    ax3.autoscale_view()    # return ax1,ax2,ax3    return lc1,ax2,ax3_linesmell_matrix=generate_smell_matrix()minDistance = None;bestRoute = [];time_list=[]distance_list=[]# 展示蚁群算法性能的三个图# ax1:当前的最优路径 ax2:所有路径的信息素浓度图 ax3:最优路径下系统的消耗与迭代次数之间的关系plt.rcParams['font.family']='SimHei'x, y = generate_point(point_num)fig=plt.figure()ax1=plt.subplot(221)ax2=plt.subplot(222)ax3=plt.subplot(212)ax1.set_xlim(-10, 110)ax1.set_ylim(-10, 110)ax1.set_aspect('equal')#閉合曲线对象lc1 = mc.LineCollection([])ax1.add_collection(lc1)ax1.scatter(x, y)ax2.set_xlim(-10, 110)ax2.set_ylim(-10, 110)ax2.set_aspect('equal')lc2 = mc.LineCollection([])ax2.add_collection(lc2)ax2.scatter(x, y)ax3_line,=ax3.plot([],[])# ax3.set_xlabel('iter_times')ax3.set_xlabel('迭代次数')ax3.set_ylabel('代价')# ax3.set_ylabel('distance')ax3.xaxis.set_major_locator(MultipleLocator(50))point_mat, distance_matrix = point_matrix(x, y)# plt.rcParams["animation.convert_path"]=r"D:\MySoftWare\ImageMagick-7.0.5-9-portable-Q16-x64\magick.exe"anim1=FuncAnimation(fig,animate,init_func=init,frames=30)# anim1.save('demo.mp4',fps=2)    #放在前面会存储前面的帧数取消注释能够存储动态画面做成视频plt.show()# anim1.save('demo.gif',fps=2,dpi=80, writer='ffmpeg')

---

运行结果

---