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

【LKH算法体验】Python调用LKH算法求TSP问题,第1张

【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 ([email protected])
# 
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坐标
11.3042.312173.9182.179
23.6391.315184.0612.37
34.1772.244193.782.212
43.7121.399203.6762.578
53.4881.535214.0292.838
63.3261.556224.2632.931
73.2381.229233.4291.908
84.1961.044243.5072.376
94.3120.79253.3942.643
104.3860.57263.4393.201
113.0071.97272.9353.24
122.5621.756283.143.55
132.7881.491292.5452.357
142.3811.676302.7782.826
151.3320.695312.372.975
163.7151.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的算法,计算结果优于其他算法。计算时间效率很高!

欢迎分享,转载请注明来源:内存溢出

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

(0)
打赏 微信扫一扫 微信扫一扫 支付宝扫一扫 支付宝扫一扫
上一篇 2022-05-16
下一篇 2022-05-16

发表评论

登录后才能评论

评论列表(0条)

保存