变邻域搜索算法

变邻域搜索算法(Variable Neighborhood Search,VNS)，于1997年由Hansen和Mladenovi首次提出，已经成为国内外的一个研究热点。作为一种经典的启发式算法，其在众多领域涌现了大量的研究成果。它的基本思想是在搜索过程中系统地改变邻域结构集来拓展搜索范围, 获得局部最优解, 再基于此局部最优解重新系统地改变邻域结构集拓展搜索范围找到另一个局部最优解的过程。
启发式方法一般用于生成NP-hard问题的初始解，好的启发式算法不仅能较快找出当前解周围的局部最优解，更能够跳出局部最优的桎梏，去更广阔的的天地寻找更满意的解，从而收敛于全局最优。VNS就是这样一种算法，它能够在当前邻域中找到最优解，又能跳出当前邻域寻找更好的解，因此，只要邻域结构集设置得较为合适，它有很大概率能够收敛于全局最优解。

VNS算法的步骤如下：
1.初始化 选择邻域结构集Nk (k=1, …, kmax) 和停止准则, 并给出初始解x。

2.重复如下步骤直到满足停止准则:

2.1 设置k=1;

2.2 直到k=kmax, 重复如下步骤:

2.2.1 随机搜索 在x的第k个邻域结构中随机产生x′ (x′∈Nk (x) ) ;

2.2.2 局部搜索 以x′为初始解, 应用一些局部搜索方法获得的局部最优解, 对应局部最优解为x*l。

2.2.3 更新 如果局部最优解优于当前最优解, 设置x=x*l, 继续在邻域结构N1内搜索;否则设置k=k+1。

对于经典的NP-hard问题–TSP问题，VNS算法是对其进行求解的一个有效方法。

#变邻域搜索解TSP问题
import numpy as np
import random as rd
import copy

def disCal(path):
    dis = 0
    for i in range(len(path)-1):
        dis += distmat[path[i]][path[i+1]]
    dis += distmat[path[0]][path[-1]]
    return dis

def shaking(solution):
    solutioni = []
    for i in range(0, 12, 3):
        lis = solution[i:i+3]
        solutioni.append(lis)
    sequence = [i for i in range(4)]
    rd.shuffle(sequence)
    sol = []
    for se in sequence:
        sol += solutioni[se]
    return sol

def variableNeighborhoodDescent(solution):
    i = 0
    dis , k = float("inf") , -1
    while i < 3:
        if i == 0:
            neiborSolution = neighborhoodOne(solution)
        elif i == 1:
            neiborSolution = neighborhoodTwo(solution)
        elif i == 2:
            neiborSolution = neighborhoodThree(solution)
        for j in range(len(neiborSolution)):
            if disCal(neiborSolution[j]) < dis:
                dis = disCal(neiborSolution[j])
                k = j
        if dis < disCal(solution):
            solution = neiborSolution[k]
            i = 0
        else:
            i += 1
    return disCal(solution),solution

def neighborhoodOne(sol):      #swap算子
    neighbor = []
    for i in range(len(sol)):
        for j in range(i+1,len(sol)):
            s = copy.deepcopy(sol)
            x = s[j]
            s[j] = s[i]
            s[i] = x
            neighbor.append(s)
    return neighbor

def neighborhoodTwo(sol):    #two_opt_swap算子
    neighbor = []
    for i in range(len(sol)):
        for j in range(i + 3, len(sol)):     #这里j从i+3开始是为了不产生跟swap算子重复的解
            s = copy.deepcopy(sol)
            s1 = s[i:j+1]
            s1.reverse()
            s = s[:i] + s1 + s[j+1:]
            neighbor.append(s)
    return neighbor

def neighborhoodThree(sol):  #two_h_opt_swap算子
    neighbor = []
    for i in range(len(sol)):
        for j in range(i+1,len(sol)):
            s = copy.deepcopy(sol)
            s = [s[i]] + [s[j]] + s[:i] + s[i+1:j] + s[j+1:]
            neighbor.append(s)
    return neighbor

distmat = np.array([[0,350,290,670,600,500,660,440,720,410,480,970],
                 [350,0,340,360,280,375,555,490,785,760,700,1100],
                 [290,340,0,580,410,630,795,680,1030,695,780,1300],
                 [670,360,580,0,260,380,610,805,870,1100,1000,1100],
                 [600,280,410,260,0,610,780,735,1030,1000,960,1300],
                 [500,375,630,380,610,0,160,645,500,950,815,950],
                 [660,555,795,610,780,160,0,495,345,820,680,830],
                 [440,490,680,805,735,645,495,0,350,435,300,625],
                 [720,785,1030,870,1030,500,345,350,0,475,320,485],
                 [410,760,695,1100,1000,950,820,435,475,0,265,745],
                 [480,700,780,1000,960,815,680,300,320,265,0,585],
                 [970,1100,1300,1100,1300,950,830,625,485,745,585,0]])

if __name__ == '__main__':
    currentSolution = [i for i in range(distmat.shape[0])]
    rd.shuffle(currentSolution)                #随机产生初始解
    shorterDisrance = disCal(currentSolution)
    iterx, iterxMax = 0, 5
    while iterx < iterxMax:
        currentSolution = shaking(currentSolution)
        currentDistance, currentSolution = variableNeighborhoodDescent(currentSolution)
        if currentDistance < shorterDisrance:
            shorterDisrance = currentDistance
            iterx = 0
        else:
            iterx += 1
    print(shorterDisrance)
    print(currentSolution)