数学建模大赛第一日

模拟退火算法！

from random import*
import numpy as np
from math import*
from matplotlib import pyplot as plt

citys = np.loadtxt("data.csv",dtype=int,delimiter=",")

number_of_citys = 25
distance = np.zeros((number_of_citys,number_of_citys))
for i in range(number_of_citys):
    for j in range(number_of_citys):
        distance[i][j] = sqrt((citys[i][0]-citys[j][0])**2+(citys[i][1]-citys[j][1])**2)

print(distance)
#初始化参数
iteration1 = 2000                #外循环迭代次数
T0 = 100000                      #初始温度，取大些
Tf = 1                           #截止温度，可以不用
alpha = 0.95                     #温度更新因子
iteration2 = 10                  #内循环迭代次数
fbest = 0                        #最佳距离

#初始化初解
x = []
for i in range(100):
    x.append(i)
np.random.shuffle(x)
x = np.array(x)
for j in range(len(x) - 1):
    fbest = fbest + distance[x[j]][x[j + 1]]
fbest = fbest + distance[x[-1]][x[0]]
xbest = x.copy()
f_now = fbest
x_now = xbest.copy()

for i in range(iteration1):
    for k in range(iteration2):
        #生成新解
        x1 = [0 for q in range(number_of_citys)]
        n1,n2 = randint(0,number_of_citys-1),randint(0,number_of_citys-1)
        n = [n1,n2]
        n.sort()
        n1,n2 = n
        #n1为0单独写
        if n1 > 0:
            x1[0:n1] = x_now[0:n1]
            x1[n1:n2+1] = x_now[n2:n1-1:-1]
            x1[n2+1:number_of_citys] = x_now[n2+1:number_of_citys]
        else:
            x1[0:n1] = x_now[0:n1]
            x1[n1:n2+1] = x_now[n2::-1]
            x1[n2+1:number_of_citys] = x_now[n2+1:number_of_citys]
        s = 0;
        for j in range(len(x1) - 1):
            s = s + distance[x1[j]][x1[j + 1]]
        s = s + distance[x1[-1]][x1[0]]
        #判断是否更新解
        if s <= f_now:
            f_now = s
            x_now = x1.copy()
        if s > f_now:
            deltaf = s - f_now
            if random() < exp(-deltaf/T0):
                f_now = s
                x_now = x1.copy()
        if s < fbest:
            fbest = s
            xbest = x1.copy()

    T0 = alpha * T0                #更新温度

    # if T0 < Tf:                  #停止准则为最低温度时可以取消注释
    #     break

#打印最佳路线和最佳距离
print(xbest)
print(fbest)

#绘制结果
plt.title('SA_TSP')
plt.xlabel('x')
plt.ylabel('y')
plt.plot(citys[...,0],citys[...,1],'ob',ms = 3)
plt.plot(citys[xbest,0],citys[xbest,1])
plt.plot([citys[xbest[-1],0],citys[xbest[0],0]],[citys[xbest[-1],1],citys[xbest[0],1]],ms = 2)
plt.show()