【LKH算法体验】Python调用LKH算法求TSP问题_python

【LKH算法体验】Python调用LKH算法求TSP问题一、LKH算法简介

Keld Helsgaun 是丹麦罗斯基勒大学计算机科学专业的名誉副教授。他于 1973 年在哥本哈根大学获得DIKU 计算机科学硕士学位。他自 1975 年以来一直在罗斯基勒大学工作。他的研究兴趣包括人工智能（问题解决和启发式）和组合优化。

LKH 是Lin-Kernighan解决旅行商（TSP）问题启发式的有效实现。
计算实验表明，LKH 非常有效。尽管该算法是近似的，但以令人印象深刻的高效率产生最佳解决方案。LKH 已经为我们能够获得的所
有已解决问题提供了最佳解决方案；包括一个 109399 个城市的实例（最大的非平凡实例求解到最优）。此外，该算法改进了一系列具有未知最优值的大规模实例的已知解决方案，其中包括 1,904,711 个城市实例 ( World TSP )。

TSP（旅行售货商问题），又名为Hamilton环游问题。

对此问题感兴趣的可以阅读下面的网站：

http://www.seas.gwu.edu/~simhaweb/champalg/tsp/tsp.html

http://www.akira.ruc.dk/~keld/research/LKH/

K. Helsgaun 的 LKH (Lin-Kernighan-Helsgaun)方法是目前求解 TSP 问题的最有效方法。令人敬佩的是，为了促进研究，K. Helsgaun 在网站上发布了其算法的完整代码和他自己的研究论文。如果有时间和精力的话，好好的读读他的论文，将会很有收获的。感觉他的论文和在人脸识别中 P.Viola 的论文 “Robust Real-Time Face Detection” 一样的经典和让人兴奋。

K. Helsgaun 在其论文中，揭示了 TSP 最优解的若干特点：

最小生成树上的边出现在最优解中的概率非常高；

一个测试例子（532-city problem）发现，其最优解中有一条边e(from,to), 节点to 是 from 的第22个最近邻接点；

为了较好的度量一条边出现在最优解中的可能性，K. Helsgaun 引入了 1-tree 的概念并由此定义了边的 alpha-measure（测度）；

对于上面的 532-city 例子，若使用 alpha-measure ，则最优解上的边最多为第14个alpha-nearest。

进一步，K. Helsgaun 考虑了距离矩阵 C=( c i j c_{ij} cij) 的一种变换（每个节点 i 的关联边长增加一个常量 p i p_i pi）：

C = ( c i j ) = > D = ( d i j ) : d i j = c i j + p i + p j C=(c_{ij}) => D=(d_{ij}): d_{ij}= c_{ij} + p_i + p_j C=(cij)=>D=(dij):dij=cij+pi+pj
这种变换使得TSP最优解保持不变，但却会改变最小生成树和 1-tree。选取合适的这种变换，alpha-measure 可以表现的更好。例如上面的 532-city 例子，适当变换后其最优解上的边最多为第5个alpha-nearest。这正如矩阵计算中的条件数那样，适当的变换后条件数会好的多。K. Helsgaun 使用 subgradient optimization 方法来选取这种变换，并且在变换的过程中，有可能发现TSP问题的真实最有解。这是因为如果计算得到的最小的 1-tree 正好是个环游（tour）的话，那么这个tour就是最优解。

关于 1-tree 和 alpha-measure（测度）的高效计算也值得注意，展示了图论的灵活应用。

总之，对于研究TSP问题而言，LKH方法非常值得关注。

二、调用LKH的Python接口

如何用Python调用LKH这个比较牛X的TSP solver，我们可以在github上找到一位大神写的LKH的Pyhthon接口（网址链接：https://github.com/unr-arl/LKH_TSP），除了Python接口还有调用LKH的matlab接口。matlab接口，在上篇博文(https://blog.csdn.net/baidu/article/details/124672911)），已做介绍。

接下来废话不说，进入正题，在github上下载LKH的Python接口文件（文件名：InvokeLKH.py），原接口文件是在linux系统里运行的，我们在win10系统里运行存在一些问题，自己可以参照如下代码进行修改，或者简单复制一下。

#        __InvokeLKH__
#        Interface the TSP LKH Solver
# 
#         This script is a simple python interface to a compiled 
#        version of the LKH TSP Solver. It requires that the
#        solver is compiled at the given directories.
#
#
#        Example Syntax:
#        python InvokeLKH.py
#
#    This script is part of the "utils" section of the StructuralInspectionPlanner
#    Toolbox. A set of elementary components are released together with this 
#    path-planning toolbox in order to make further developments easier. 
#     
#    Authors: 
#    Kostas Alexis (kalexis@unr.edu)
# 
import os
import shutil
import math
import numpy as np

# Change with the Cost Matrix of your problem or 
# consider using it as an argument
CostMatrix = np.ones((10,10))*100

fname_tsp = "test"
user_comment = "a comment by the user"

# Change these directories based on where you have 
# a compiled executable of the LKH TSP Solver
lkh_dir = '/LKH/LKH-2.0.9/'
tsplib_dir = '/TSPLIB/'
lkh_cmd = 'LKH'
pwd=os.getcwd()


def writeTSPLIBfile_FE(fname_tsp,CostMatrix,user_comment):

    dims_tsp = len(CostMatrix)
    name_line = 'NAME : ' + fname_tsp + '\n'
    type_line = 'TYPE: TSP' + '\n'
    comment_line = 'COMMENT : ' + user_comment + '\n'
    tsp_line = 'TYPE : ' + 'TSP' + '\n'
    dimension_line = 'DIMENSION : ' + str(dims_tsp) + '\n'
    edge_weight_type_line = 'EDGE_WEIGHT_TYPE : ' + 'EXPLICIT' + '\n' # explicit only
    edge_weight_format_line = 'EDGE_WEIGHT_FORMAT: ' + 'FULL_MATRIX' + '\n'
    display_data_type_line ='DISPLAY_DATA_TYPE: ' + 'NO_DISPLAY' + '\n' # 'NO_DISPLAY'
    edge_weight_section_line = 'EDGE_WEIGHT_SECTION' + '\n'
    eof_line = 'EOF\n'
    Cost_Matrix_STRline = []
    for i in range(0,dims_tsp):
        cost_matrix_strline = ''
        for j in range(0,dims_tsp-1):
            cost_matrix_strline = cost_matrix_strline + str(int(CostMatrix[i][j])) + ' '

        j = dims_tsp-1
        cost_matrix_strline = cost_matrix_strline + str(int(CostMatrix[i][j]))
        cost_matrix_strline = cost_matrix_strline + '\n'
        Cost_Matrix_STRline.append(cost_matrix_strline)
    
    fileID = open((pwd + tsplib_dir + fname_tsp + '.tsp'), "w")
    print(name_line)
    fileID.write(name_line)
    fileID.write(comment_line)
    fileID.write(tsp_line)
    fileID.write(dimension_line)
    fileID.write(edge_weight_type_line)
    fileID.write(edge_weight_format_line)
    fileID.write(edge_weight_section_line)
    for i in range(0,len(Cost_Matrix_STRline)):
        fileID.write(Cost_Matrix_STRline[i])

    fileID.write(eof_line)
    fileID.close()

    fileID2 = open((pwd + tsplib_dir + fname_tsp + '.par'), "w")

    problem_file_line = 'PROBLEM_FILE = ' + pwd + tsplib_dir + fname_tsp + '.tsp' + '\n' # remove pwd + tsplib_dir
    optimum_line = 'OPTIMUM 378032' + '\n'
    move_type_line = 'MOVE_TYPE = 5' + '\n'
    patching_c_line = 'PATCHING_C = 3' + '\n'
    patching_a_line = 'PATCHING_A = 2' + '\n'
    runs_line = 'RUNS = 10' + '\n'
    tour_file_line = 'TOUR_FILE = ' + fname_tsp + '.txt' + '\n'

    fileID2.write(problem_file_line)
    fileID2.write(optimum_line)
    fileID2.write(move_type_line)
    fileID2.write(patching_c_line)
    fileID2.write(patching_a_line)
    fileID2.write(runs_line)
    fileID2.write(tour_file_line)
    fileID2.close()
    return fileID, fileID2

def copy_toTSPLIBdir_cmd(fname_basis):
    srcfile=pwd + '/' + fname_basis + '.txt'
    dstpath=pwd + tsplib_dir
    shutil.copy(srcfile, dstpath)

def run_LKHsolver_cmd(fname_basis):
    run_lkh_cmd =  pwd + lkh_dir  + lkh_cmd + ' ' + pwd + tsplib_dir + fname_basis + '.par'
    os.system(run_lkh_cmd)

def rm_solution_file_cmd(fname_basis):
    del_file=pwd + '/' + fname_basis + '.txt'
    os.remove(del_file)

def main(): 
    
    [fileID1,fileID2] = writeTSPLIBfile_FE(fname_tsp,CostMatrix,user_comment)
    run_LKHsolver_cmd(fname_tsp)
    copy_toTSPLIBdir_cmd(fname_tsp)
    rm_solution_file_cmd(fname_tsp)
    

if __name__ == "__main__":
    main()

解释一下：LKHdir就是在http://akira.ruc.dk/~keld/research/LKH/这个网站上下载的LKH.exe的存储位置。LKH.exe可以在画红色箭头的位置下载。

tsplib_dir 就是存放.par文件的路径，fname_tsp就是.tsp的文件名。

我们按接口的要求在当前文件夹下，建立LKH和TSPLIB文件夹，然后在LKH文件夹下建立LKH-2.0.9子文件夹。把LKH.exe文件，放在LKH-2.0.9目录下。

然后在命令行执行：python InvokeLKH.py

出现如下显示：

得到测试数据的最优值为Cost.min=1000，然后我们打开TSPLIB文件夹，发现里面有3个文件，如下图：

即：在TSPLIB文件夹中生成了3个文件：test.par（参数文件），test.tsp(距离矩阵数据文件)，test.txt（生成的最优路径文件）。

三、应用举例

使用LKH算法，求解31个城市的TSP问题：

假定31个城市的位置坐标下表所列。寻找出一条最短的遍历31个城市的路径。

城市编号	X坐标	Y坐标	城市编号	X坐标	Y坐标
1	1.304	2.312	17	3.918	2.179
2	3.639	1.315	18	4.061	2.37
3	4.177	2.244	19	3.78	2.212
4	3.712	1.399	20	3.676	2.578
5	3.488	1.535	21	4.029	2.838
6	3.326	1.556	22	4.263	2.931
7	3.238	1.229	23	3.429	1.908
8	4.196	1.044	24	3.507	2.376
9	4.312	0.79	25	3.394	2.643
10	4.386	0.57	26	3.439	3.201
11	3.007	1.97	27	2.935	3.24
12	2.562	1.756	28	3.14	3.55
13	2.788	1.491	29	2.545	2.357
14	2.381	1.676	30	2.778	2.826
15	1.332	0.695	31	2.37	2.975
16	3.715	1.678

第1步：参照"InvokeLKH.py"接口文件，编写如下代码文件，调用LKH算法进行处理：

#文件名：city031.py
#使用LKH算法求解31城市的TSP问题
import os
import shutil
import math
import numpy as np

def clac_distance(X, Y):
    """
    计算两个城市之间的欧氏距离，二范数
    :param X: 城市X的坐标.np.array数组
    :param Y: 城市Y的坐标.np.array数组
    :return:
    """
    distance_matrix = np.zeros((city_num, city_num))
    for i in range(city_num):
        for j in range(city_num):
            if i == j:
                continue

            distance = np.sqrt((X[i] - X[j]) ** 2 + (Y[i] - Y[j]) ** 2)
            distance_matrix[i][j] = distance

    return distance_matrix

#读取31座城市坐标
coord = []
with open("data.txt", "r") as lines:
    lines = lines.readlines()
for line in lines:
    xy = line.split()
    coord.append(xy)
coord = np.array(coord)
w, h = coord.shape
coordinates = np.zeros((w, h), float)
for i in range(w):
    for j in range(h):
        coordinates[i, j] = float(coord[i, j])
city_x=coordinates[:,0]
city_y=coordinates[:,1]
#城市数量
city_num = coordinates.shape[0]
CostMatrix = clac_distance(city_x, city_y)*1000    #将距离矩阵放大1000倍（LKH算法只能处理整数）

fname_tsp = "city031"
user_comment = "a comment by the user"

# Change these directories based on where you have 
# a compiled executable of the LKH TSP Solver
lkh_dir = '/LKH/LKH-2.0.9/'
tsplib_dir = '/TSPLIB/'
lkh_cmd = 'LKH'
pwd=os.getcwd()


def writeTSPLIBfile_FE(fname_tsp,CostMatrix,user_comment):

    dims_tsp = len(CostMatrix)
    name_line = 'NAME : ' + fname_tsp + '\n'
    type_line = 'TYPE: TSP' + '\n'
    comment_line = 'COMMENT : ' + user_comment + '\n'
    tsp_line = 'TYPE : ' + 'TSP' + '\n'
    dimension_line = 'DIMENSION : ' + str(dims_tsp) + '\n'
    edge_weight_type_line = 'EDGE_WEIGHT_TYPE : ' + 'EXPLICIT' + '\n' # explicit only
    edge_weight_format_line = 'EDGE_WEIGHT_FORMAT: ' + 'FULL_MATRIX' + '\n'
    display_data_type_line ='DISPLAY_DATA_TYPE: ' + 'NO_DISPLAY' + '\n' # 'NO_DISPLAY'
    edge_weight_section_line = 'EDGE_WEIGHT_SECTION' + '\n'
    eof_line = 'EOF\n'
    Cost_Matrix_STRline = []
    for i in range(0,dims_tsp):
        cost_matrix_strline = ''
        for j in range(0,dims_tsp-1):
            cost_matrix_strline = cost_matrix_strline + str(int(CostMatrix[i][j])) + ' '

        j = dims_tsp-1
        cost_matrix_strline = cost_matrix_strline + str(int(CostMatrix[i][j]))
        cost_matrix_strline = cost_matrix_strline + '\n'
        Cost_Matrix_STRline.append(cost_matrix_strline)
    
    fileID = open((pwd + tsplib_dir + fname_tsp + '.tsp'), "w")
    print(name_line)
    fileID.write(name_line)
    fileID.write(comment_line)
    fileID.write(tsp_line)
    fileID.write(dimension_line)
    fileID.write(edge_weight_type_line)
    fileID.write(edge_weight_format_line)
    fileID.write(edge_weight_section_line)
    for i in range(0,len(Cost_Matrix_STRline)):
        fileID.write(Cost_Matrix_STRline[i])

    fileID.write(eof_line)
    fileID.close()

    fileID2 = open((pwd + tsplib_dir + fname_tsp + '.par'), "w")

    problem_file_line = 'PROBLEM_FILE = ' + pwd + tsplib_dir + fname_tsp + '.tsp' + '\n' # remove pwd + tsplib_dir
    optimum_line = 'OPTIMUM 378032' + '\n'
    move_type_line = 'MOVE_TYPE = 5' + '\n'
    patching_c_line = 'PATCHING_C = 3' + '\n'
    patching_a_line = 'PATCHING_A = 2' + '\n'
    runs_line = 'RUNS = 10' + '\n'
    tour_file_line = 'TOUR_FILE = ' + fname_tsp + '.txt' + '\n'

    fileID2.write(problem_file_line)
    fileID2.write(optimum_line)
    fileID2.write(move_type_line)
    fileID2.write(patching_c_line)
    fileID2.write(patching_a_line)
    fileID2.write(runs_line)
    fileID2.write(tour_file_line)
    fileID2.close()
    return fileID, fileID2

def copy_toTSPLIBdir_cmd(fname_basis):
    srcfile=pwd + '/' + fname_basis + '.txt'
    dstpath=pwd + tsplib_dir
    shutil.copy(srcfile, dstpath)

def run_LKHsolver_cmd(fname_basis):
    run_lkh_cmd =  pwd + lkh_dir  + lkh_cmd + ' ' + pwd + tsplib_dir + fname_basis + '.par'
    os.system(run_lkh_cmd)

def rm_solution_file_cmd(fname_basis):
    del_file=pwd + '/' + fname_basis + '.txt'
    os.remove(del_file)

def main(): 
    
    [fileID1,fileID2] = writeTSPLIBfile_FE(fname_tsp,CostMatrix,user_comment)
    run_LKHsolver_cmd(fname_tsp)
    copy_toTSPLIBdir_cmd(fname_tsp)
    rm_solution_file_cmd(fname_tsp)
    

if __name__ == "__main__":
    main()

第2步：在命令终端运行：python city031.py，看到如下画面：

耗时约0.08s，计算效率很高。

在“TSPLIB”文件夹下，生成了如下3个文件。

第3步：编写代码，读取最优路径结果文件city031.txt，输出最优路径，还原出最优值（总距离）。

#读取最优路径，还原出最优值(总距离)
import os
import math
import numpy as np

def clac_distance(X, Y):
    """
    计算两个城市之间的欧氏距离，二范数
    :param X: 城市X的坐标.np.array数组
    :param Y: 城市Y的坐标.np.array数组
    :return:
    """
    distance_matrix = np.zeros((city_num, city_num))
    for i in range(city_num):
        for j in range(city_num):
            if i == j:
                continue

            distance = np.sqrt((X[i] - X[j]) ** 2 + (Y[i] - Y[j]) ** 2)
            distance_matrix[i][j] = distance

    return distance_matrix

def get_total_distance(x):
    dista = 0
    for i in range(len(x)):
        if i == len(x) - 1:
            dista += distance[x[i]][x[0]]
        else:
            dista += distance[x[i]][x[i + 1]]
    return dista

def print_path(best_line):
    result_cur_best=[]
    for i in best_line:
        result_cur_best+=[i]
    result_path = result_cur_best
    result_path.append(result_path[0])
    return result_path

if __name__ == "__main__":
    
    #读取31座城市坐标
    coord = []
    with open("data.txt", "r") as lines:
        lines = lines.readlines()
    for line in lines:
        xy = line.split()
        coord.append(xy)
    coord = np.array(coord)
    w, h = coord.shape
    coordinates = np.zeros((w, h), float)
    for i in range(w):
        for j in range(h):
            coordinates[i, j] = float(coord[i, j])
    city_x=coordinates[:,0]
    city_y=coordinates[:,1]
    #城市数量
    city_num = coordinates.shape[0]
    distance = clac_distance(city_x, city_y)

    #使用将距离矩阵放大1000倍后获得的最优路径，重新计算总距离
    x=[]
    txt=open('./TSPLIB/city031.txt','r').readlines()
    for i in range(len(txt)):        #读取路径数据
        if i>5 and i<=(city_num+5):
            x+=[eval(txt[i][:-1])]
    #输出结果
    print("最优路径：")
    print(print_path(x))
    x=[i-1 for i in x]   #将列表中每个元素减1
    grad=get_total_distance(x)
    print('总距离:',grad)

结果如下：

最优值为：15.3825

这应该是目前得到的最好结果，你如果得到比这个结果更优的值，不要忘记通知我一声。

LKH算法确实是目前比较牛X的算法，计算结果优于其他算法。计算时间效率很高！

欢迎分享，转载请注明来源：内存溢出

原文地址: http://outofmemory.cn/langs/918309.html

【LKH算法体验】Python调用LKH算法求TSP问题

发表评论

评论列表（0条）