演算法思想

模仿自然界退火現象而得,利用了物理中固體物質的退火過程與一般優化問題的相似性從某一初始溫度開始,伴隨溫度的不斷下降,結合概率突跳特性在解空間中隨機尋找全域性最優解

具體應用

旅行商問題,即TSP問題(Travelling SalesmanProblem)又譯為旅行推銷員問題、貨郎擔問題,是數學領域中著名問題之一。
假設有一個旅行商人要拜訪n個城市,他必須選擇所要走的路徑,路徑的限制是每個城市只能拜訪一次,而且最後要回到原來出發的城市。
路徑的選擇目標是要求得的路徑路程為所有路徑之中的最小值。

吃今為止,這類問題中沒有一個找到有效演算法。傾向於接受NP完全問題(NP-Complet或NPC)和NP難題(NP-Hard或NPH)不存在有效演算法這一猜想,
認為這類問題的大型例項不能用精確演算法求解,必須尋求這類問題的有效的近似演算法。

TSP的應用領域

TSP “旅行商問題”的應用領域包括:如何規劃最合理高效的道路交通,以減少擁堵;如何更好地規劃物流,以減少運營成本;在網際網路環境中如何更好地設定節點,以更好地讓資訊流動等。

程式碼實現

import numpy as np
import matplotlib.pyplot as plt
import pdb

"旅行商問題 ( TSP , Traveling Salesman Problem )"
coordinates = np.array([[565.0, 575.0], [25.0, 185.0], [345.0, 750.0], [945.0, 685.0], [845.0, 655.0],
                        [880.0, 660.0], [25.0, 230.0], [525.0, 1000.0], [580.0, 1175.0], [650.0, 1130.0],
                        [1605.0, 620.0], [1220.0, 580.0], [1465.0, 200.0], [1530.0, 5.0], [845.0, 680.0],
                        [725.0, 370.0], [145.0, 665.0], [415.0, 635.0], [510.0, 875.0], [560.0, 365.0],
                        [300.0, 465.0], [520.0, 585.0], [480.0, 415.0], [835.0, 625.0], [975.0, 580.0],
                        [1215.0, 245.0], [1320.0, 315.0], [1250.0, 400.0], [660.0, 180.0], [410.0, 250.0],
                        [420.0, 555.0], [575.0, 665.0], [1150.0, 1160.0], [700.0, 580.0], [685.0, 595.0],
                        [685.0, 610.0], [770.0, 610.0], [795.0, 645.0], [720.0, 635.0], [760.0, 650.0],
                        [475.0, 960.0], [95.0, 260.0], [875.0, 920.0], [700.0, 500.0], [555.0, 815.0],
                        [830.0, 485.0], [1170.0, 65.0], [830.0, 610.0], [605.0, 625.0], [595.0, 360.0],
                        [1340.0, 725.0], [1740.0, 245.0]])


# 得到距離矩陣的函式
def getdistmat(coordinates):
    num = coordinates.shape[0]  # 52個座標點
    distmat = np.zeros((52, 52))  # 52X52距離矩陣
    for i in range(num):
        for j in range(i, num):
            distmat[i][j] = distmat[j][i] = np.linalg.norm(coordinates[i] - coordinates[j])
    return distmat


def initpara():
    alpha = 0.99
    t = (1, 100)
    markovlen = 1000
    return alpha, t, markovlen


num = coordinates.shape[0]
distmat = getdistmat(coordinates)  # 得到距離矩陣

solutionnew = np.arange(num)
# valuenew = np.max(num)

solutioncurrent = solutionnew.copy()
valuecurrent = 99000  # np.max這樣的原始碼可能同樣是因為版本問題被當做函式不能正確使用，應取一個較大值作為初始值
# print(valuecurrent)

solutionbest = solutionnew.copy()
valuebest = 99000  # np.max

alpha, t2, markovlen = initpara()
t = t2[1]

result = []  # 記錄迭代過程中的最優解
while t > t2[0]:
    for i in np.arange(markovlen):

        # 下面的兩交換和三角換是兩種擾動方式，用於產生新解
        if np.random.rand() > 0.5:  # 交換路徑中的這2個節點的順序
            # np.random.rand()產生[0, 1)區間的均勻隨機數
            while True:  # 產生兩個不同的隨機數
                loc1 = np.int(np.ceil(np.random.rand() * (num - 1)))
                loc2 = np.int(np.ceil(np.random.rand() * (num - 1)))
                ## print(loc1,loc2)
                if loc1 != loc2:
                    break
            solutionnew[loc1], solutionnew[loc2] = solutionnew[loc2], solutionnew[loc1]
        else:  # 三交換
            while True:
                loc1 = np.int(np.ceil(np.random.rand() * (num - 1)))
                loc2 = np.int(np.ceil(np.random.rand() * (num - 1)))
                loc3 = np.int(np.ceil(np.random.rand() * (num - 1)))

                if ((loc1 != loc2) & (loc2 != loc3) & (loc1 != loc3)):
                    break

            # 下面的三個判斷語句使得loc1<loc2<loc3
            if loc1 > loc2:
                loc1, loc2 = loc2, loc1
            if loc2 > loc3:
                loc2, loc3 = loc3, loc2
            if loc1 > loc2:
                loc1, loc2 = loc2, loc1

            # 下面的三行程式碼將[loc1,loc2)區間的資料插入到loc3之後
            tmplist = solutionnew[loc1:loc2].copy()
            solutionnew[loc1:loc3 - loc2 + 1 + loc1] = solutionnew[loc2:loc3 + 1].copy()
            solutionnew[loc3 - loc2 + 1 + loc1:loc3 + 1] = tmplist.copy()

        valuenew = 0
        for i in range(num - 1):
            valuenew += distmat[solutionnew[i]][solutionnew[i + 1]]
        valuenew += distmat[solutionnew[0]][solutionnew[51]]
        # print (valuenew)
        if valuenew < valuecurrent:  # 接受該解

            # 更新solutioncurrent 和solutionbest
            valuecurrent = valuenew
            solutioncurrent = solutionnew.copy()

            if valuenew < valuebest:
                valuebest = valuenew
                solutionbest = solutionnew.copy()
        else:  # 按一定的概率接受該解
            if np.random.rand() < np.exp(-(valuenew - valuecurrent) / t):
                valuecurrent = valuenew
                solutioncurrent = solutionnew.copy()
            else:
                solutionnew = solutioncurrent.copy()
    t = alpha * t
    result.append(valuebest)
    print(t)  # 程式執行時間較長，列印t來監視程式進展速度

# 用來顯示結果
plt.plot(np.array(result))
plt.ylabel("bestvalue")
plt.xlabel("t")
plt.show()