模拟退火算法:TSP问题

针对中国31个省会城市的TSP问题:

通过模拟退火算法解决,坐标输入:

[[1304, 2312], [3639, 1315], [4177, 2244], [3712, 1399], [3488, 1535], [3326, 1556],
 [3238, 1229], [4196, 1004], [4312, 790], [4386, 570], [3007, 1970], [2562, 1756],
 [2788, 1491], [2381, 1676], [1332, 695], [3715, 1678], [3918, 2179], [4061, 2370],
 [3780, 2212], [3676, 2578], [4029, 2838], [4263, 2931], [3429, 1908], [3507, 2367],
 [3394, 2643], [3439, 3201], [2935, 3240], [3140, 3550], [2545, 2357], [2778, 2826],
 [2370, 2975]]

这里有3种扰动操作:1.两元素互换     2.任意长度片段进行逆序操作    3.选取任意元素进行插入操作

具体实现:

import sys
import math
import os
import random
from matplotlib import pyplot as plt

END_T = 1  # 结束温度
T = 100  # 温度
L = 100 * 31  # 马尔科夫链长
K = 0.9  # 降温系数

dist = []  # 邻接矩阵
city = [[1304, 2312], [3639, 1315], [4177, 2244], [3712, 1399], [3488, 1535], [3326, 1556],
        [3238, 1229], [4196, 1004], [4312, 790], [4386, 570], [3007, 1970], [2562, 1756],
        [2788, 1491], [2381, 1676], [1332, 695], [3715, 1678], [3918, 2179], [4061, 2370],
        [3780, 2212], [3676, 2578], [4029, 2838], [4263, 2931], [3429, 1908], [3507, 2367],
        [3394, 2643], [3439, 3201], [2935, 3240], [3140, 3550], [2545, 2357], [2778, 2826],
        [2370, 2975]]  # 各城市坐标
initial_path = list(range(0, 31))
losses = []


def __initial__():
    for i in city:
        dist.append([math.sqrt(math.pow(i[0] - j[0], 2) + math.pow(i[1] - j[1], 2)) for j in city])


def E(x: list):
    length = dist[x[0]][x[30]]
    i = 1
    while i < 31:
        length += dist[x[i]][x[i - 1]]
        i += 1
    return length


def accept_probability(xj, xi):
    return math.exp((E(xi) - E(xj)) / T)


def random_operation(path):
    """ 互换操作 """
    a, b = random.randint(1, 30), random.randint(1, 30)
    while a == b:  # 不能自己和自己交换,且从点0出发
        b = random.randint(1, 30)
    path[a], path[b] = path[b], path[a]
    return path


def random_operation2(path):
    """ 逆序操作 """
    a, b = random.randint(1, 30), random.randint(1, 30)
    if a > b:
        a, b = b, a
    while a < b:
        path[a], path[b] = path[b], path[a]
        a += 1
        b -= 1
    return path


def random_operation3(path):
    """ 插入操作 """
    a, b = random.randint(1, 30), random.randint(1, 30)
    if a == b:  # 没有操作
        return path
    if a < b:  # 后插
        i, temp = a, path[a]
        while i < b:
            path[i] = path[i + 1]
            i += 1
        path[b] = temp
    else:  # 前插
        i, temp = a, path[a]
        while i > b:
            path[i] = path[i - 1]
            i -= 1
        path[b] = temp
    return path


def 模拟退火(operator=None):
    if operator is None:
        operator = random_operation
    losses = []
    xi, xj = initial_path.copy(), None
    x_best = xi
    temperature = T
    while temperature > END_T:
        temperature *= K
        loop = L
        while loop > 0:
            loop -= 1
            xj = xi.copy()
            losses.append(E(x_best))
            operator(xj)  # 扰动操作
            if E(xj) < E(xi):
                xi = xj.copy()  # 距离更短则无条件接受
                if E(xj) < E(x_best):
                    x_best = xj.copy()  # 全局最小
                continue
            probablity = accept_probability(xj, xi)
            # 距离更长,按概率接受
            if random.random() < probablity:  # 小于概率,接受
                xi = xj.copy()
        # print('Dist:', E(xi), '\tTemperature=', temperature)
    return [x_best, losses]


def get_union_path_pic(ele1, ele2, ele3):
    eles = [ele1, ele2, ele3]
    titles = ['互换操作', '逆序操作', '插入操作']
    counter = 0
    plt.rcParams['font.sans-serif'] = ['SimHei']
    plt.rcParams['axes.unicode_minus'] = False
    while counter < 3:
        plt.subplot(3, 2, 2 * counter + 1)
        plt.title(titles[counter])
        i = 1
        while i < 31:
            plt.plot([city[eles[counter][0][i]][0], city[eles[counter][0][i - 1]][0]],
                     [city[eles[counter][0][i]][1], city[eles[counter][0][i - 1]][1]], color='blue')
            i += 1
        plt.plot([city[eles[counter][0][0]][0], city[eles[counter][0][30]][0]],
                 [city[eles[counter][0][0]][1], city[eles[counter][0][30]][1]], color='red')
        # line
        plt.subplot(3, 2, 2 * counter + 2)
        plt.plot(eles[counter][1])
        plt.title('Distance=' + str(E(eles[counter][0])))
        counter += 1
    plt.show()


def get_path_pic(path, losses):
    plt.subplot(1, 2, 1)
    i = 1
    while i < 31:
        plt.plot([city[path[i]][0], city[path[i - 1]][0]], [city[path[i]][1], city[path[i - 1]][1]], color='blue')
        i += 1
    plt.plot([city[path[0]][0], city[path[30]][0]], [city[path[0]][1], city[path[30]][1]], color='red')
    # line
    plt.subplot(1, 2, 2)
    plt.plot(losses)
    plt.title('Distance=' + str(E(path)))
    plt.show()


if __name__ == '__main__':
    __initial__()
    ele1 = 模拟退火()
    print('Dist:', E(ele1[0]), 'Path', ele1[0])
    ele2 = 模拟退火(operator=random_operation2)
    print('Dist:', E(ele2[0]), 'Path', ele2[0])
    ele3 = 模拟退火(operator=random_operation3)
    print('Dist:', E(ele3[0]), 'Path', ele3[0])
    get_union_path_pic(ele1, ele2, ele3)

运行效果:

3种扰动操作的运行效果

 

  • 0
    点赞
  • 9
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值