Comparison of cost vs quality for 50 random cities

Eric P. Hanson

We will use TravelingSalesmanExact to compute the exact cost and compare to the estimated best costs found by TravelingSalesmanHeuristics with various settings of quality. First we load the packages:

using TravelingSalesmanExact, GLPK
using TravelingSalesmanHeuristics
using TravelingSalesmanExact: ATT, euclidean_distance
using TravelingSalesmanBenchmarks
using Plots
using Printf, Random
gr(fmt=:svg)
Plots.GRBackend()

For plotting, we will override TravelingSalesmanExact's plot_cities function to use Plots, since we load that package anyway. We also add a mutating version.

function get_plt_coords(cities)
    n = length(cities)
    inc(a) = a == n ? one(a) : a + 1
    return [cities[inc(j)][1] for j = 0:n], [cities[inc(j)][2] for j = 0:n]
end
TravelingSalesmanExact.plot_cities(cities; kwargs...) = plot(get_plt_coords(cities)...; kwargs...)
plot_cities!(cities; kwargs...) = plot!(get_plt_coords(cities)...; kwargs...)

function plot_tours(cities, pairs; kwargs...)
    plts = []
    colors = sequential_palette(0, length(pairs)+1)[2:end]
    for (index, (tour, label)) in enumerate(pairs)
        plt = plot_cities(cities[tour], label = label, linewidth = 2,
                            color = colors[index])
        plot_cities!(cities; title="Comparison of tours", st=:scatter, label="City locations", markersize = 5, kwargs...)
        push!(plts, plt)
    end
    return plts
end
plot_tours (generic function with 1 method)

Then we generate a cost matrix by choosing cities at random, although with a fixed seed for reproducibility:

Random.seed!(155);
N = 50
cities = [100*rand(2) for _ = 1:N]
cost = [ euclidean_distance(cities[i], cities[j]) for i=1:N, j=1:N ]
50×50 Array{Float64,2}:
  0.0     44.8088  19.9692   43.5758  …   64.6586    54.5324   52.8702
 44.8088   0.0     58.114    71.448       86.0398    43.4927   26.5045
 19.9692  58.114    0.0      24.8275      46.5257    74.4112   71.2316
 43.5758  71.448   24.8275    0.0         21.7496    96.4314   90.0815
 45.1792  44.3448  40.5373   36.709       44.1871    81.7696   69.0097
 54.7367  85.363   34.8922   13.9695  …   18.9868   108.721   103.388 
 17.971   27.3493  31.5108   49.1805      67.5469    48.5724   40.9607
 23.3844  68.1842  20.3087   41.8838      63.189     72.111    74.3188
 39.2596  13.9264  56.1032   73.6303      90.7237    30.5753   16.8647
 55.6061  77.5743  37.8824   13.5463       9.18445  106.725    98.5898
  ⋮                                   ⋱                               
 29.4587  17.4684  45.9701   63.9585      81.7009    36.0023   26.1233
 44.3408  81.3876  24.813    18.2641      34.1703    98.8478   96.0092
 47.7985  92.5511  41.9794   56.4328      74.7942    90.9141   96.2662
 47.2316  30.7208  49.3069   51.5112      60.4885    71.934    56.8522
 70.4713  67.1173  61.5265   47.1347  …   41.136    107.333    93.1449
 21.0135  56.5979   3.46318  22.7644      44.2661    74.9837   70.8998
 64.6586  86.0398  46.5257   21.7496       0.0      115.892   107.542 
 54.5324  43.4927  74.4112   96.4314     115.892      0.0      19.3796
 52.8702  26.5045  71.2316   90.0815     107.542     19.3796    0.0

Now we will compute an optimal tour and cost, and plot these versus those found by the tsp_solve function of TravelingSalesmanHeuristics.

t_exact, c_exact = get_optimal_tour(cost, with_optimizer(GLPK.Optimizer))

c(q) = solve_tsp(cost; quality_factor = q)[2]

qs = range(10, stop = 100, step = 10)

plot(qs, c, xlabel="quality", ylabel="Cost", label="solve_tsp 1", title="random cities")
for j = 2:5
    plot!(qs, c, label="solve_tsp $j")
end
hline!([c_exact], label="Exact cost")

We've run tsp_solve five times for each quality, since the cost will vary from run to run due to the randomness of the heuristics.

We can compare plots of the tours:

t_heuristic_100, c_heuristic_100 = solve_tsp(cost; quality_factor = 100)
t_heuristic_40, c_heuristic_40 = solve_tsp(cost; quality_factor = 40)

plts = plot_tours(cities, [
    (t_heuristic_40, @sprintf("Heuristic; quality 40, cost=%.2f",   c_heuristic_40)),
    (t_heuristic_100, @sprintf("Heuristic; quality 100, cost=%.2f", c_heuristic_100)),
    (t_exact, @sprintf("Optimal tour; cost=%.2f", c_exact))
    ])
for plt in plts
    display(plt)
end

We can repeat this several times to compare various choices of random cities.

function compare_cities(N)
    cities = [100*rand(2) for _ = 1:N]
    cost = [ euclidean_distance(cities[i], cities[j]) for i=1:N, j=1:N ]
    t_exact, c_exact = get_optimal_tour(cost, with_optimizer(GLPK.Optimizer))

    c(q) = solve_tsp(cost; quality_factor = q)[2]

    qs = range(10, stop = 100, step = 10)

    line_plt = plot(qs, c, xlabel="quality", ylabel="Cost", label="solve_tsp 1", title="random cities")
    for j = 2:5
        plot!(qs, c, label="solve_tsp $j")
    end
    hline!([c_exact], label="Exact cost")

    t_heuristic_100, c_heuristic_100 = solve_tsp(cost; quality_factor = 100)
    t_heuristic_40, c_heuristic_40 = solve_tsp(cost; quality_factor = 40)
    city_plts = plot_tours(cities, [
    (t_heuristic_40, @sprintf("Heuristic; quality 40, cost=%.2f",   c_heuristic_40)),
    (t_heuristic_100, @sprintf("Heuristic; quality 100, cost=%.2f", c_heuristic_100)),
    (t_exact, @sprintf("Optimal tour; cost=%.2f", c_exact))
    ])
    return line_plt, city_plts...
end

for j = 1:5
    plts = compare_cities(50)
    for plt in plts
        display(plt)
    end
end

Appendix

These benchmarks are a part of the TravelingSalesmanBenchmarks.jl repository, found at: https://github.com/ericphanson/TravelingSalesmanBenchmarks.jl, based on the weave branch of DiffEqBenchmarks.

To locally run this tutorial, do the following commands:

using TravelingSalesmanBenchmarks
TravelingSalesmanBenchmarks.weave_file("random_50_cities_cost_vs_quality.jmd")

Computer Information:

Julia Version 1.1.0
Commit 80516ca202 (2019-01-21 21:24 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin14.5.0)
  CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, broadwell)
Environment:
  JULIA_EDITOR = "/Applications/Visual Studio Code.app/Contents/Frameworks/Code Helper.app/Contents/MacOS/Code Helper"

Package Information:

Status: `/Users/eh540/Dropbox (Personal)/LinkedFolders/Julia/dev/TravelingSalesmanBenchmarks/Project.toml`
[a93c6f00-e57d-5684-b7b6-d8193f3e46c0] DataFrames 0.17.1
[60bf3e95-4087-53dc-ae20-288a0d20c6a6] GLPK 0.9.1
[2e9cd046-0924-5485-92f1-d5272153d98b] Gurobi 0.6.0
[7073ff75-c697-5162-941a-fcdaad2a7d2a] IJulia 1.18.1
[91a5bcdd-55d7-5caf-9e0b-520d859cae80] Plots 0.24.0
[737fac7d-4440-55ef-927e-002196e95561] TravelingSalesmanExact 0.3.2
[8c8f4381-2cdd-507c-846c-be2bcff6f45f] TravelingSalesmanHeuristics 0.3.1
[44d3d7a6-8a23-5bf8-98c5-b353f8df5ec9] Weave 0.9.0
[ddb6d928-2868-570f-bddf-ab3f9cf99eb6] YAML 0.3.2
[b77e0a4c-d291-57a0-90e8-8db25a27a240] InteractiveUtils 
[d6f4376e-aef5-505a-96c1-9c027394607a] Markdown 
[44cfe95a-1eb2-52ea-b672-e2afdf69b78f] Pkg 
[de0858da-6303-5e67-8744-51eddeeeb8d7] Printf 
[9a3f8284-a2c9-5f02-9a11-845980a1fd5c] Random 
[10745b16-79ce-11e8-11f9-7d13ad32a3b2] Statistics