Test cases (#32)

* admonition and path rotation

* basic test cases

* docstrings

docs/src/index.md

@@ -17,11 +17,17 @@ A word of warning: the heuristics implemented are
 `TravelingSalesmanHeuristics` is a registered package, so you can install it with the REPL's `Pkg` mode: `]add TravelingSalesmanHeuristics`. Load it with `using TravelingSalesmanHeuristics`.
 
 ## Usage
-All problems are specified through a square distance matrix `D` where `D[i,j]` represents the cost of traveling from the `i`-th to the `j`-th city. Your distance matrix need not be symmetric and could probably even contain negative values, though I make no guarantee about behavior when using negative values.
+All problems are specified through a square distance matrix `D` where `D[i,j]` represents the cost of traveling from the `i`-th to the `j`-th city. Your distance matrix need not be symmetric or non-negative.
 
 !!! note
     Because problems are specified by dense distance matrices, this package is not well suited to problem instances with sparse distance matrix structure, i.e. problems with large numbers of cities where each city is connected to few other cities.
 
+!!! tip
+    TravelingSalesmanHeuristics supports distance matrices with arbitrary real element types: integers, floats, rationals, etc. However, mixing element types in your distance matrix is not recommended for performance reasons.
+
+!!! danger
+    TravelingSalesmanHeuristics does _not_ support distance matrices with arbitrary indexing; indices must be `1:n` in both dimensions for `n` cities.
+
 The simplest way to use this package is with the `solve_tsp` function:
 
 ```@docs

src/TravelingSalesmanHeuristics.jl

@@ -23,6 +23,9 @@ tend to lead to better solutions found at the cost of more computation time.
     typically the case. A `quality_factor` of 100 neither guarantees an optimal solution
     nor the best solution  that can be found via extensive use of the methods in this package.
 
+!!! danger
+    TravelingSalesmanHeuristics does not support distance matrices with arbitrary indexing; indices must be `1:n` in both dimensions for `n` cities.
+
 See also: individual heuristic methods such as [`farthest_insertion`](@ref)
 """
 function solve_tsp(distmat::AbstractMatrix{T}; quality_factor::Real = 40.0) where {T<:Real}

src/helpers.jl

@@ -11,12 +11,30 @@ function check_square(m, msg)
     return n
 end
 
+"""
+    rotate_circuit(circuit::AbstractArray{<:Integer}, desired_start::Integer)
+
+"Rotate" a closed circuit so that it starts (and ends) at `desired_start`. E.g.
+`rotate_circuit([1,2,3,1], 2) == [2,3,1,2]`. `circuit` must be a closed path.
+"""
+function rotate_circuit(circuit, desired_start)
+    if first(circuit) != last(circuit)
+        throw(DomainError(circuit, "Circuit passed to rotate_circuit is not closed."))
+    end
+    first(circuit) == desired_start && return copy(circuit)
+    start_ind = findfirst(circuit .== desired_start)
+    return vcat(
+        circuit[start_ind:end],
+        circuit[2:start_ind]
+    )
+end
+
 """
     legal_circuit(circuit::AbstractArray{<:Integer})
 
 Check that an array of integers is a valid circuit. A valid circuit over `n` locations has
-length `n+1`. The first `n` entries are a permutation of `1, ..., n`, and the `(n+1)`-st entry
-is equal to the first entry.
+length `n+1`. The first `n` entries are a permutation of `1, ..., n`, and the `(n+1)`-st
+entry is equal to the first entry.
 """
 function legal_circuit(circuit)
     n = length(circuit) - 1
@@ -134,7 +152,7 @@ to improvements with some minimum magnitude defined by the element type of the
 distance matrix."
 improvement_threshold(T::Type{<:Integer}) = one(T)
 improvement_threshold(T::Type{<:AbstractFloat}) = sqrt(eps(one(T)))
-improvement_threshold(T::Type{<:Real}) = sqrt(eps(1.0))
+improvement_threshold(::Type{<:Real}) = sqrt(eps(1.0))
 
 "Cost of inserting city `k` after index `after` in path `path` with costs `distmat`."
 function inscost(

test/performance_testing.jl

@@ -0,0 +1,98 @@
+#=
+This file contains a set of synthetic TSP test cases, plus some small helper functions
+for running the test cases. It's mostly intended for comparing algorithms in development
+and probably isn't of much use to end users. Note the dependence on DataFrames and
+Distributions as well, which aren't part of the TravelingSalesmanHeuristics environment.
+=#
+
+using Random
+using LinearAlgebra: norm
+
+using DataFrames
+using Distributions
+using TravelingSalesmanHeuristics
+TSP = TravelingSalesmanHeuristics
+
+struct TestCase
+    distmat
+    description::String
+end
+
+"""
+    runcase(alg, case)
+
+Test an algorithm on a TestCase. `alg` must be a function `distmat -> (path, cost)`.
+Return a NamedTuple describing the result. The `bound` and `excess` keys in the
+resulting tuple are respectively the vertwise bound and a simple scale-adjusted
+cost (lower is still better)."""
+function runcase(alg, case::TestCase)
+    time = @elapsed path, cost = alg(case.distmat)
+    path = TSP.rotate_circuit(path, 1)
+    bound = TSP.vertwise_bound(case.distmat)
+    excess = (cost - bound) / abs(bound)
+    return (;time, path, cost, bound, excess, description = case.description)
+end
+
+"""Run multiple test cases and return the result as a DataFrame.
+`alg` should be a function `distmat -> (path, cost)`."""
+function runcases(alg, cases)
+    return DataFrame(runcase.(alg, cases))
+end
+
+
+"Euclidean distance matrix from a matrix of points, one col per point, one row per dim"
+function _distmat_from_pts(pts)
+    n = size(pts, 2)
+    return [norm(pts[:,i] - pts[:,j]) for i in 1:n, j in 1:n]
+end
+
+"""Generate a list of Euclidean test cases. For each dimension and number of points,
+generate a test case in a unit square, a squashed rectangle, and a spherical shell."""
+function generate_cases_euclidean(;
+    seed = 47,
+    dimensions = [2, 5, 15],
+    ns = [5, 25, 100, 500],
+)
+    Random.seed!(seed)
+    cases = TestCase[]
+    for d in dimensions, n in ns
+        # hypercube
+        pts = rand(d, n)
+        dm = _distmat_from_pts(pts)
+        push!(cases, TestCase(dm, "euclidean_hypercube_$(d)d_$(n)n"))
+        # squashed hypercube
+        pts = rand(d, n)
+        pts[1:(d ÷ 2), :] ./= 10
+        dm = _distmat_from_pts(pts)
+        push!(cases, TestCase(dm, "euclidean_squashedhypercube_$(d)d_$(n)n"))
+        # hypersphere shell
+        pts = randn(d, n)
+        for i in 1:n
+            pts[:,i] ./= norm(pts[:,i])
+        end
+        dm = _distmat_from_pts(pts)
+        push!(cases, TestCase(dm, "euclidean_sphereshell_$(d)d_$(n)n"))
+    end
+    return cases
+end
+
+"""Generate non-Euclidean test cases by starting with the Eculidean cases
+and adding IID noise to each distance measurement, producing cases which are
+nonsymmetric and may have negative travel costs."""
+function generate_cases_euclidean_plus_noise(;
+    dists = [
+        Normal(),
+        Exponential(1),
+    ],
+    euclidean_args...
+)
+    euclidean_cases = generate_cases_euclidean(euclidean_args...)
+    return [
+        TestCase(
+            case.distmat .+ rand(dist, size(case.distmat)),
+            case.description * "_plus_" * string(dist)
+        )
+        for case in euclidean_cases
+        for dist in dists
+    ]
+end

test/runtests.jl

@@ -11,7 +11,7 @@ using LinearAlgebra
 # generate a Euclidean distance matrix for n points in the unit square
 function generate_planar_distmat(n)
     pts = rand(2, n)
-    dm = [norm(pts[:,i] - pts[:,j]) for i in 1:n, j in 1:n]
+    return [norm(pts[:,i] - pts[:,j]) for i in 1:n, j in 1:n]
 end
 
 # test that a path is acceptable:
@@ -147,6 +147,13 @@ function test_path_costs()
     end
 end
 
+function test_rotate()
+    # rotate a circuit
+    @test TSP.rotate_circuit([11,12,13,11], 12) == [12,13,11,12]
+    # not a circuit
+    @test_throws DomainError TSP.rotate_circuit([1,2,3,4], 2)
+end
+
 function test_solve_tsp()
     dm = rand(10,10)
     quality_factors = [-1, 0, 1.1, 35, 101]
@@ -195,6 +202,7 @@ test_farthest_insertion()
 test_simulated_annealing()
 test_bounds()
 test_path_costs()
+test_rotate()
 test_solve_tsp()
 test_two_opt()
 test_atypical_types()