Traveling Salesman Problem

Overview

The Traveling Salesman Problem is one of the most famous problems in computer science. In what follows, we'll describe the problem and show you how to find a solution.

The problem

Back in the days when salesmen traveled door-to-door hawking vacuums and encyclopedias, they had to plan their routes, from house to house or city to city. The shorter the route, the better. Finding the shortest route that visits a set of locations is an exponentially difficult problem: finding the shortest path for 20 cities is much more than twice as hard as 10 cities.

An exhaustive search of all possible paths would be guaranteed to find the shortest, but is computationally intractable for all but small sets of locations. For larger problems, optimization techniques are needed to intelligently search the solution space and find near-optimal solutions.

Mathematically, traveling salesman problems can be represented as a graph, where the locations are the nodes and the edges (or arcs) represent direct routes between the nodes. The weight of each edge is the distance between the nodes. The goal is to find the path with the shortest sum of weights. Below, we see a simple four-node graph and the shortest cycle that visits every node:

10 25 25 40 45 50 A B C D Sadly, your browser does not support inline SVG.

In addition to finding solutions to the classical Traveling Salesman Problem, OR-Tools also provides methods for more general types of TSPs, including the following:

  • Asymmetric cost problems—The traditional TSP is symmetric: the distance from point A to point B equals the distance from point B to point A. However, the cost of shipping items from point A to point B might not equal the cost of shipping them from point B to point A. OR-Tools can also handle problems that have asymmetric costs.
  • Prize-collecting TSPs, where benefits accrue from visiting nodes
  • TSP with time windows

Solving TSPs with the Google Directions API

Google also provides a way to solve simple TSPs of real-world locations without downloading OR-Tools. If you have a Google Directions API key, you can solve TSPs of real-world locations with the Directions API, providing the locations in a URL and getting the response back as JSON. You'll need your own free Directions API key for development, or an enterprise key for commercial use.

As an example, here's a URL that can be used to find a short tour of winemaking regions in South Australia, beginning in Adelaide. If you want to try this from your browser, replace API_KEY at the end of the URL with your key.

https://maps.googleapis.com/maps/api/directions/json?origin=Adelaide,SA&destination=Adelaide,SA&waypoints=optimize:true|Barossa+Valley,SA|Clare,SA|Connawarra,SA|McLaren+Vale,SA&key=API_KEY

The result will be a long JSON response detailing the solution, complete with Google Maps directions:

{
   "routes" : [
      {
         "bounds" : {
            "northeast" : {
               "lat" : -33.8347115,
               "lng" : 140.8547058
            },
            "southwest" : {
               "lat" : -37.3511758,
               "lng" : 138.4951576
            }
         },
         "copyrights" : "Map data ©2014 Google",
         "legs" : [
            {
               "distance" : {
                  "text" : "139 km",
                  "value" : 139119
               },
               "duration" : {
                  "text" : "1 hour 51 mins",
                  "value" : 6648
               },
               "end_address" : "Clare SA 5453, Australia",
               "end_location" : {
                  "lat" : -33.8333395,
                  "lng" : 138.6117283
               },
               "start_address" : "Adelaide SA, Australia",
               "start_location" : {
                  "lat" : -34.9285894,
                  "lng" : 138.5999429
               },
               "steps" : [
                  {
                     "distance" : {
                        "text" : "70 m",
                        "value" : 70
                     },
                     "duration" : {
                        "text" : "1 min",
                        "value" : 6
                     },
                     "end_location" : {
                        "lat" : -34.9285338,
                        "lng" : 138.6007031
                     },
                     "html_instructions" : "Head \u003cb\u003eeast\u003c/b\u003e on \u003cb\u003eReconciliation Plaza\u003c/b\u003e toward \u003cb\u003eVictoria Square\u003c/b\u003e",
...

Solving TSPs with OR-Tools

No solver can find the shortest paths for all problems. Lots of work has gone into techniques for quickly finding near-optimal solutions, and into proving bounds about how closely these techniques approach optimality.

You can solve TSPs using the OR-Tools vehicle routing library, a collection of algorithms designed especially for TSPs, and more general problems with multiple vehicles. The routing library is an added layer on top of the constraint programming solver. See RoutingModel for detailed information about the available methods for setting up and solving routing problems.

We next present an example that shows how to use OR-Tools to find the shortest route through the U.S. cities shown on the map below.

The following sections present a Python program that solves the TSP for these cities.

Create the data

Let's start with the data for the problem, which consists of the cities and the distances between them. The data is defined by an array, whose entry in row i and column j is the distance from city i to city j in miles, as shown below.

 [[   0, 2451,  713, 1018, 1631, 1374, 2408,  213, 2571,  875, 1420, 2145, 1972], # New York
  [2451,    0, 1745, 1524,  831, 1240,  959, 2596,  403, 1589, 1374,  357,  579], # Los Angeles
  [ 713, 1745,    0,  355,  920,  803, 1737,  851, 1858,  262,  940, 1453, 1260], # Chicago
  [1018, 1524,  355,    0,  700,  862, 1395, 1123, 1584,  466, 1056, 1280,  987], # Minneapolis
  [1631,  831,  920,  700,    0,  663, 1021, 1769,  949,  796,  879,  586,  371], # Denver
  [1374, 1240,  803,  862,  663,    0, 1681, 1551, 1765,  547,  225,  887,  999], # Dallas
  [2408,  959, 1737, 1395, 1021, 1681,    0, 2493,  678, 1724, 1891, 1114,  701], # Seattle
  [ 213, 2596,  851, 1123, 1769, 1551, 2493,    0, 2699, 1038, 1605, 2300, 2099], # Boston
  [2571,  403, 1858, 1584,  949, 1765,  678, 2699,    0, 1744, 1645,  653,  600], # San Francisco
  [ 875, 1589,  262,  466,  796,  547, 1724, 1038, 1744,    0,  679, 1272, 1162], # St. Louis
  [1420, 1374,  940, 1056,  879,  225, 1891, 1605, 1645,  679,    0, 1017, 1200], # Houston
  [2145,  357, 1453, 1280,  586,  887, 1114, 2300,  653, 1272, 1017,    0,  504], # Phoenix
  [1972,  579, 1260,  987,  371,  999,  701, 2099,  600, 1162,  1200,  504,   0]] # Salt Lake City
The order of the cities in the array is the following:
city_names = ["New York", "Los Angeles", "Chicago",
              "Minneapolis", "Denver", "Dallas", "Seattle", "Boston", "San Francisco",
              "St. Louis", "Houston", "Phoenix", "Salt Lake City"]

Declare the solver

The following code declares an instance of the solver and sets the default options for it.

  # Cities
  city_names = ["New York", "Los Angeles", "Chicago", "Minneapolis", "Denver", "Dallas", "Seattle",
                "Boston", "San Francisco", "St. Louis", "Houston", "Phoenix", "Salt Lake City"]
  tsp_size = len(city_names)
  num_routes = 1    # The number of routes, which is 1 in the TSP.
  # Nodes are indexed from 0 to tsp_size - 1. The depot is the starting node of the route.
  depot = 0

  # Create routing model
  if tsp_size > 0:
    routing = pywrapcp.RoutingModel(tsp_size, num_routes, depot)
    search_parameters = pywrapcp.DefaultRoutingSearchParameters()

    # Create the distance callback, which takes two arguments (the from and to node indices)
    # and returns the distance between these nodes.
    dist_between_nodes = CreateDistanceCallback()
    dist_callback = dist_between_nodes.Distance
    routing.SetArcCostEvaluatorOfAllVehicles(dist_callback)

The solver instance is declared by

 routing = pywrapcp.RoutingModel(tsp_size, 1)
The first argument (tsp_size) is the number of cities. The second is the number of routes allowed, which is 1 for the Traveling Salesman Problem.

The last three lines in the above code create the distance callback—a function that takes any pair of cities and returns the distance between them—and then passes the callback to the solver by the method SetArcCostEvaluatorOfAllVehicles.

Create the distance callback

The distance callback is created by the class CreateDistanceCallback shown below.

class CreateDistanceCallback(object):
  """Create callback to calculate distances between points."""
  def __init__(self):
    """Array of distances between points."""

    self.matrix = [
    [   0, 2451,  713, 1018, 1631, 1374, 2408,  213, 2571,  875, 1420, 2145, 1972], # New York
    [2451,    0, 1745, 1524,  831, 1240,  959, 2596,  403, 1589, 1374,  357,  579], # Los Angeles
    [ 713, 1745,    0,  355,  920,  803, 1737,  851, 1858,  262,  940, 1453, 1260], # Chicago
    [1018, 1524,  355,    0,  700,  862, 1395, 1123, 1584,  466, 1056, 1280,  987], # Minneapolis
    [1631,  831,  920,  700,    0,  663, 1021, 1769,  949,  796,  879,  586,  371], # Denver
    [1374, 1240,  803,  862,  663,    0, 1681, 1551, 1765,  547,  225,  887,  999], # Dallas
    [2408,  959, 1737, 1395, 1021, 1681,    0, 2493,  678, 1724, 1891, 1114,  701], # Seattle
    [ 213, 2596,  851, 1123, 1769, 1551, 2493,    0, 2699, 1038, 1605, 2300, 2099], # Boston
    [2571,  403, 1858, 1584,  949, 1765,  678, 2699,    0, 1744, 1645,  653,  600], # San Francisco
    [ 875, 1589,  262,  466,  796,  547, 1724, 1038, 1744,    0,  679, 1272, 1162], # St. Louis
    [1420, 1374,  940, 1056,  879,  225, 1891, 1605, 1645,  679,    0, 1017, 1200], # Houston
    [2145,  357, 1453, 1280,  586,  887, 1114, 2300,  653, 1272, 1017,    0,  504], # Phoenix
    [1972,  579, 1260,  987,  371,  999,  701, 2099,  600, 1162,  1200,  504,   0]] # Salt Lake City

  def Distance(self, from_node, to_node):
    return int(self.matrix[from_node][to_node])

The callback is the function Distance, which takes a pair of cities (from_node and to_node) and returns the corresponding entry of the distance array.

Note: Since the routing solver does all computations with integers, any function that creates a distance callback should convert the value returned by the distance callback to an integer. In this example, the conversion is done by the Python int() function in the line

return int(self.matrix[from_node][to_node])
While converting isn't essential in this case, because all distances are given as integers, it is a good practice to always convert the return value to an integer.

The solver section of the program sets the distance callback to be the Distance method in the CreateDistanceCallback class.

  dist_between_nodes = CreateDistanceCallback()
  dist_callback = dist_between_nodes.Distance

Finally, the program passes the callback to the solver:

routing.SetArcCostEvaluatorOfAllVehicles(dist_callback)
Notice that the method that specifies the distance callback is named SetArcCostEvaluatorOfAllVehicles. In this case, there is just one vehicle for the single salesman. But in more general routing problems, there can be multiple vehicles, and each node must be visited by one of them. In general, the function to be minimized is the total cost of all the routes. In this example, we are assuming that there is just one route, and the cost of each arc is proportional to its distance—so the minimum cost route is the same as the minimum distance route.

Invoke the solver and display the results

The following code invokes the solver (named routing) and displays the results.

    # Solve, returns a solution if any.
    assignment = routing.SolveWithParameters(search_parameters)
    if assignment:
      # Solution cost.
      print "Total distance: " + str(assignment.ObjectiveValue()) + " miles\n"
      # Inspect solution.
      # Only one route here; otherwise iterate from 0 to routing.vehicles() - 1
      route_number = 0
      index = routing.Start(route_number) # Index of the variable for the starting node.
      route = ''
      while not routing.IsEnd(index):
        # Convert variable indices to node indices in the displayed route.
        route += str(city_names[routing.IndexToNode(index)]) + ' -> '
        index = assignment.Value(routing.NextVar(index))
      route += str(city_names[routing.IndexToNode(index)])
      print "Route:\n\n" + route
    else:
      print 'No solution found.'
  else:
    print 'Specify an instance greater than 0.'

Displaying the route in the output

To display the route in the output, we first get the starting index of the route.

index = routing.Start(route_number)
index is the index of the variable that represents the starting node of the route, called the depot. The index of a node variable is not always the same as the index of the node itself (its row number in the locations array). So, when we add nodes to the route displayed in the output, we convert the indices of node variable to node indices, using the method routing.IndexToNode():
route += str(city_names[routing.IndexToNode(index)]) + ' -> '

The program iterates through the nodes of the route by

index = assignment.Value(routing.NextVar(index))
which takes an index for a node variable and returns the index for the next node variable in the route.

Here is the output of the program.

A depot must be specified, setting one at node 0
Total distance: 7293 miles

Route:

New York -> Boston -> Chicago -> Minneapolis -> Denver -> Salt Lake City -> Seattle ->
San Francisco -> Los Angeles -> Phoenix -> Houston -> Dallas -> St. Louis -> New York
The output gives the route and its total distance, 7293 miles. (The depot is the starting and ending point for the tour. Since we don't care where the depot is, OR-Tools uses the first location.)

The entire program

The entire program is shown below.

from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import routing_enums_pb2

# Distance callback
class CreateDistanceCallback(object):
  """Create callback to calculate distances between points."""
  def __init__(self):
    """Array of distances between points."""

    self.matrix = [
    [   0, 2451,  713, 1018, 1631, 1374, 2408,  213, 2571,  875, 1420, 2145, 1972], # New York
    [2451,    0, 1745, 1524,  831, 1240,  959, 2596,  403, 1589, 1374,  357,  579], # Los Angeles
    [ 713, 1745,    0,  355,  920,  803, 1737,  851, 1858,  262,  940, 1453, 1260], # Chicago
    [1018, 1524,  355,    0,  700,  862, 1395, 1123, 1584,  466, 1056, 1280,  987], # Minneapolis
    [1631,  831,  920,  700,    0,  663, 1021, 1769,  949,  796,  879,  586,  371], # Denver
    [1374, 1240,  803,  862,  663,    0, 1681, 1551, 1765,  547,  225,  887,  999], # Dallas
    [2408,  959, 1737, 1395, 1021, 1681,    0, 2493,  678, 1724, 1891, 1114,  701], # Seattle
    [ 213, 2596,  851, 1123, 1769, 1551, 2493,    0, 2699, 1038, 1605, 2300, 2099], # Boston
    [2571,  403, 1858, 1584,  949, 1765,  678, 2699,    0, 1744, 1645,  653,  600], # San Francisco
    [ 875, 1589,  262,  466,  796,  547, 1724, 1038, 1744,    0,  679, 1272, 1162], # St. Louis
    [1420, 1374,  940, 1056,  879,  225, 1891, 1605, 1645,  679,    0, 1017, 1200], # Houston
    [2145,  357, 1453, 1280,  586,  887, 1114, 2300,  653, 1272, 1017,    0,  504], # Phoenix
    [1972,  579, 1260,  987,  371,  999,  701, 2099,  600, 1162,  1200,  504,   0]] # Salt Lake City

  def Distance(self, from_node, to_node):
    return int(self.matrix[from_node][to_node])
def main():

  # Cities
  city_names = ["New York", "Los Angeles", "Chicago", "Minneapolis", "Denver", "Dallas", "Seattle",
                "Boston", "San Francisco", "St. Louis", "Houston", "Phoenix", "Salt Lake City"]
  tsp_size = len(city_names)
  num_routes = 1    # The number of routes, which is 1 in the TSP.
  # Nodes are indexed from 0 to tsp_size - 1. The depot is the starting node of the route.
  depot = 0

  # Create routing model
  if tsp_size > 0:
    routing = pywrapcp.RoutingModel(tsp_size, num_routes, depot)
    search_parameters = pywrapcp.DefaultRoutingSearchParameters()

    # Create the distance callback, which takes two arguments (the from and to node indices)
    # and returns the distance between these nodes.
    dist_between_nodes = CreateDistanceCallback()
    dist_callback = dist_between_nodes.Distance
    routing.SetArcCostEvaluatorOfAllVehicles(dist_callback)
    # Solve, returns a solution if any.
    assignment = routing.SolveWithParameters(search_parameters)
    if assignment:
      # Solution cost.
      print "Total distance: " + str(assignment.ObjectiveValue()) + " miles\n"
      # Inspect solution.
      # Only one route here; otherwise iterate from 0 to routing.vehicles() - 1
      route_number = 0
      index = routing.Start(route_number) # Index of the variable for the starting node.
      route = ''
      while not routing.IsEnd(index):
        # Convert variable indices to node indices in the displayed route.
        route += str(city_names[routing.IndexToNode(index)]) + ' -> '
        index = assignment.Value(routing.NextVar(index))
      route += str(city_names[routing.IndexToNode(index)])
      print "Route:\n\n" + route
    else:
      print 'No solution found.'
  else:
    print 'Specify an instance greater than 0.'

if __name__ == '__main__':
  main()

Add a distance function

In the previous example, the distances between cities are explicitly defined in an array. You can also define the distance callback by a function that calculates the distance between any two cities. The next example uses a function that calculates the distance between any two points on Earth from their latitudes and longitudes.

The advantage of this approach is that if you want to add a location to the network, you don't need to re-write the entire distance matrix—you simply add the latitude and longitude of the new location.

The distance function

The function that calculates the distance between locations is based on the Haversine formula, which navigators have used since the 1800s to calculate distances at sea. The function is shown below.

def distance(lat1, long1, lat2, long2):
    # Note: The formula used in this function is not exact, as it assumes
    # the Earth is a perfect sphere.

    # Mean radius of Earth in miles
    radius_earth = 3959

    # Convert latitude and longitude to
    # spherical coordinates in radians.
    degrees_to_radians = math.pi/180.0
    phi1 = lat1 * degrees_to_radians
    phi2 = lat2 * degrees_to_radians
    lambda1 = long1 * degrees_to_radians
    lambda2 = long2 * degrees_to_radians
    dphi = phi2 - phi1
    dlambda = lambda2 - lambda1

    a = haversine(dphi) + math.cos(phi1) * math.cos(phi2) * haversine(dlambda)
    c = 2 * math.atan2(math.sqrt(a), math.sqrt(1-a))
    d = radius_earth * c
    return d

def haversine(angle):
  h = math.sin(angle / 2) ** 2
  return h

Create the distance callback

As in the previous example, the distance callback is created in a class named CreateDistanceCallback, shown below.

# Distance callback

class CreateDistanceCallback(object):
  """Create callback to calculate distances between points."""

  def __init__(self):

    # Latitudes and longitudes of selected U.S. cities

    locations = [[40.71,  -74.01], # New York
                 [34.05, -118.24], # Los Angeles
                 [41.88,  -87.63], # Chicago
                 [44.98,  -93.27], # Minneapolis
                 [39.74, -104.99], # Denver
                 [32.78,  -96.89], # Dallas
                 [47.61, -122.33], # Seattle
                 [42.36,  -71.06], # Boston
                 [37.77, -122.42], # San Francisco
                 [38.63,  -90.20], # St. Louis
                 [29.76,  -95.37], # Houston
                 [33.45, -112.07], # Phoenix
                 [40.76, -111.89]] # Salt Lake City


    """Create the distance matrix."""
    size = len(locations)
    self.matrix = {}

    for from_node in xrange(size):
      self.matrix[from_node] = {}
      for to_node in xrange(size):
        if from_node == to_node:
          self.matrix[from_node][to_node] = 0
        else:
          x1 = locations[from_node][0]
          y1 = locations[from_node][1]
          x2 = locations[to_node][0]
          y2 = locations[to_node][1]
          self.matrix[from_node][to_node] = distance(x1, y1, x2, y2)

  def Distance(self, from_node, to_node):
    return int(self.matrix[from_node][to_node])

First, the program defines an array, locations, containing the latitudes and longitudes of the cities in the problem. Next, it creates the distance matrix, self.matrix. In this example, the program calculates the entries of the matrix by looping over each pair of cities and using the distance function.

Here's the output of the Python program:

Total distance: 7198 miles

Route:

New York -> Boston -> Chicago -> Minneapolis -> Denver -> Salt Lake City -> Seattle ->
San Francisco -> Los Angeles -> Phoenix -> Dallas -> Houston -> St. Louis -> New York

The distance function is not exact, because it assumes that the Earth is a perfect sphere, while in reality it is wider at the equator than at the poles. This accounts the slightly different values for the total distance of the route in the two examples. However, the route returned in both cases is the same.

To add a city to the problem in this example, you simply add its latitude and longitude to the array locations. If you hard coded the distances between cities, as in the previous example, you would need to add a row and column to the distance matrix, corresponding to the new city.

Example—drilling a circuit board

The next example is one of the problems in the TSPLIB collection: a280, a set of 280 abstract locations taken from a real-life problem of drilling holes in a circuit board with an automated drill. The problem is to find the shortest route for the drill to take on the board in order to drill all of the required holes .

Here's a Google Scatter Chart of the locations in a280:

The program that finds a solution for this TSP is very similar to the previous example. The main difference is the distance callback, which calculates distances between points on the circuit board using the standard Euclidean distance, rather than the distances between cities.

The program hard codes the data for the TSP as a set of points in the plane. (Ideally, there would be a parser for TSPLIB problems.) Under the hood, the OR-Tools routing library is written in C++, but it can be used from multiple languages:

Create the data

The data for the problem consists of 280 points in the plane, shown in the scatter chart above. The program creates the data in an array, the first few entries of which are shown below.

def create_data():
  locations = [[288, 149], [288, 129], [270, 133], [256, 141], [256, 157], [246, 157], ...

Create the distance callback

The class that creates the distance callback is the same as in the previous example, except that it calls the Euclidean distance function, instead of a function that calculates distances on Earth. The code is shown below.

def distance(x1, y1, x2, y2):

    dist = math.sqrt((x1 - x2)**2 + (y1 - y2)**2)

    return dist

# Distance callback

class CreateDistanceCallback(object):
  """Create callback to calculate distances between points."""

  def __init__(self):
    """Initialize distance array."""
    locations = create_data_array()
    size = len(locations)
    self.matrix = {}

    for from_node in xrange(size):
      self.matrix[from_node] = {}
      for to_node in xrange(size):
        if from_node == to_node:
          self.matrix[from_node][to_node] = 0
        else:
          x1 = locations[from_node][0]
          y1 = locations[from_node][1]
          x2 = locations[to_node][0]
          y2 = locations[to_node][1]
          self.matrix[from_node][to_node] = distance(x1, y1, x2, y2)

  def Distance(self, from_node, to_node):
    return int(self.matrix[from_node][to_node])

The function distance calculates the Euclidean distance between any two points (x1, y1) and (x2, y2) in the plane.

Output of the program

Here is the output of the program.

Total distance: 2790

Route:

0 -> 1 -> 279 -> 2 -> 278 -> 277 -> 247 -> 248 -> 249 -> 246 -> 244 -> 243 -> 242 -> 241 -> 240 ->
239 -> 238 -> 237 -> 236 -> 235 -> 234 -> 233 -> 232 -> 231 -> 230 -> 245 -> 250 -> 229 -> 228 ->
227 -> 226 -> 225 -> 224 -> 223 -> 222 -> 221 -> 220 -> 219 -> 218 -> 217 -> 216 -> 215 -> 214 ->
213 -> 212 -> 211 -> 210 -> 209 -> 208 -> 251 -> 254 -> 255 -> 257 -> 256 -> 253 -> 252 -> 207 ->
206 -> 205 -> 204 -> 203 -> 202 -> 142 -> 141 -> 146 -> 147 -> 140 -> 139 -> 265 -> 136 -> 137 ->
138 -> 148 -> 149 -> 177 -> 176 -> 175 -> 178 -> 179 -> 180 -> 181 -> 182 -> 183 -> 184 -> 186 ->
185 -> 192 -> 196 -> 197 -> 198 -> 144 -> 145 -> 143 -> 199 -> 201 -> 200 -> 195 -> 194 -> 193 ->
191 -> 190 -> 189 -> 188 -> 187 -> 163 -> 164 -> 165 -> 166 -> 167 -> 168 -> 169 -> 171 -> 170 ->
172 -> 105 -> 106 -> 104 -> 103 -> 107 -> 109 -> 110 -> 113 -> 114 -> 116 -> 117 -> 61 -> 62 ->
63 -> 65 -> 64 -> 84 -> 85 -> 115 -> 112 -> 86 -> 83 -> 82 -> 87 -> 111 -> 108 -> 89 -> 90 -> 91 ->
102 -> 101 -> 100 -> 99 -> 98 -> 97 -> 96 -> 95 -> 94 -> 93 -> 92 -> 79 -> 88 -> 81 -> 80 -> 78 ->
77 -> 76 -> 74 -> 75 -> 73 -> 72 -> 71 -> 70 -> 69 -> 66 -> 68 -> 67 -> 57 -> 56 -> 55 -> 54 ->
53 -> 52 -> 51 -> 50 -> 49 -> 48 -> 47 -> 46 -> 45 -> 44 -> 43 -> 58 -> 60 -> 59 -> 42 -> 41 ->
40 -> 39 -> 38 -> 37 -> 36 -> 35 -> 34 -> 33 -> 32 -> 31 -> 30 -> 29 -> 124 -> 123 -> 122 -> 121 ->
120 -> 119 -> 118 -> 156 -> 157 -> 158 -> 173 -> 162 -> 161 -> 160 -> 174 -> 159 -> 150 -> 151 ->
155 -> 152 -> 154 -> 153 -> 128 -> 129 -> 130 -> 131 -> 18 -> 19 -> 20 -> 127 -> 126 -> 125 -> 28 ->
27 -> 26 -> 25 -> 21 -> 24 -> 22 -> 23 -> 13 -> 12 -> 14 -> 11 -> 10 -> 9 -> 7 -> 8 -> 6 -> 5 ->
275 -> 274 -> 273 -> 272 -> 271 -> 270 -> 15 -> 16 -> 17 -> 132 -> 133 -> 269 -> 268 -> 134 ->
135 -> 267 -> 266 -> 264 -> 263 -> 262 -> 261 -> 260 -> 258 -> 259 -> 276 -> 3 -> 4 -> 0

Here is a graph of the resulting tour:

The OR-Tools library finds the above tour very quickly: in less than a second on a typical computer. The total length of the above tour is 2790.

Entire program

Here is the entire program for the circuit board example.

import math

from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import routing_enums_pb2
def distance(x1, y1, x2, y2):

    dist = math.sqrt((x1 - x2)**2 + (y1 - y2)**2)

    return dist

# Distance callback

class CreateDistanceCallback(object):
  """Create callback to calculate distances between points."""

  def __init__(self):
    """Initialize distance array."""
    locations = create_data_array()
    size = len(locations)
    self.matrix = {}

    for from_node in xrange(size):
      self.matrix[from_node] = {}
      for to_node in xrange(size):
        if from_node == to_node:
          self.matrix[from_node][to_node] = 0
        else:
          x1 = locations[from_node][0]
          y1 = locations[from_node][1]
          x2 = locations[to_node][0]
          y2 = locations[to_node][1]
          self.matrix[from_node][to_node] = distance(x1, y1, x2, y2)

  def Distance(self, from_node, to_node):
    return int(self.matrix[from_node][to_node])
def main():
  # Create the data.
  locations = create_data_array()
  tsp_size = len(locations)
  num_routes = 1      # The number of routes, which is 1 in the TSP.
  # Nodes are indexed from 0 to tsp_size - 1. The depot is the starting node of the route.
  depot = 0
  # Create routing model.
  if tsp_size > 0:
    routing = pywrapcp.RoutingModel(tsp_size, num_routes, depot)
    search_parameters = pywrapcp.DefaultRoutingSearchParameters()

    # Callback to the distance function. The callback takes two
    # arguments (the from and to node indices) and returns the distance between them.
    dist_between_locations = CreateDistanceCallback()
    dist_callback = dist_between_locations.Distance
    routing.SetArcCostEvaluatorOfAllVehicles(dist_callback)
    # Solve, returns a solution if any.
    assignment = routing.SolveWithParameters(search_parameters)
    if assignment:

      # Solution cost.
      print "Total distance: " + str(assignment.ObjectiveValue()) + "\n"

      # Inspect solution.
      # Only one route here; otherwise iterate from 0 to routing.vehicles() - 1.
      route_number = 0
      node = routing.Start(route_number)
      start_node = node
      route = ''

      while not routing.IsEnd(node):
        route += str(node) + ' -> '
        node = assignment.Value(routing.NextVar(node))
      route += '0'
      print "Route:\n\n" + route
    else:
      print 'No solution found.'
  else:
    print 'Specify an instance greater than 0.'
def create_data_array():
  locations = [[288, 149], [288, 129], [270, 133], [256, 141], [256, 157], [246, 157], [236, 169],
[228, 169], [228, 161], [220, 169], [212, 169], [204, 169], [196, 169], [188, 169], [196, 161],
[188, 145], [172, 145], [164, 145], [156, 145], [148, 145], [140, 145], [148, 169], [164, 169],
[172, 169], [156, 169], [140, 169], [132, 169], [124, 169], [116, 161], [104, 153], [104, 161],
[104, 169], [90, 165], [80, 157], [64, 157], [64, 165], [56, 169], [56, 161], [56, 153], [56, 145],
[56, 137], [56, 129], [56, 121], [40, 121], [40, 129], [40, 137], [40, 145], [40, 153], [40, 161],
[40, 169], [32, 169], [32, 161], [32, 153], [32, 145], [32, 137], [32, 129], [32, 121], [32, 113],
[40, 113], [56, 113], [56, 105], [48, 99], [40, 99], [32, 97], [32, 89], [24, 89], [16, 97],
[16, 109], [8, 109], [8, 97], [8, 89], [8, 81], [8, 73], [8, 65], [8, 57], [16, 57], [8, 49],
[8, 41], [24, 45], [32, 41], [32, 49], [32, 57], [32, 65], [32, 73], [32, 81], [40, 83], [40, 73],
[40, 63], [40, 51], [44, 43], [44, 35], [44, 27], [32, 25], [24, 25], [16, 25], [16, 17], [24, 17],
[32, 17], [44, 11], [56, 9], [56, 17], [56, 25], [56, 33], [56, 41], [64, 41], [72, 41], [72, 49],
[56, 49], [48, 51], [56, 57], [56, 65], [48, 63], [48, 73], [56, 73], [56, 81], [48, 83], [56, 89],
[56, 97], [104, 97], [104, 105], [104, 113], [104, 121], [104, 129], [104, 137], [104, 145],
[116, 145], [124, 145], [132, 145], [132, 137], [140, 137], [148, 137], [156, 137], [164, 137],
[172, 125], [172, 117], [172, 109], [172, 101], [172, 93], [172, 85], [180, 85], [180, 77],
[180, 69], [180, 61], [180, 53], [172, 53], [172, 61], [172, 69], [172, 77], [164, 81], [148, 85],
[124, 85], [124, 93], [124, 109], [124, 125], [124, 117], [124, 101], [104, 89], [104, 81],
[104, 73], [104, 65], [104, 49], [104, 41], [104, 33], [104, 25], [104, 17], [92, 9], [80, 9],
[72, 9], [64, 21], [72, 25], [80, 25], [80, 25], [80, 41], [88, 49], [104, 57], [124, 69],
[124, 77], [132, 81], [140, 65], [132, 61], [124, 61], [124, 53], [124, 45], [124, 37], [124, 29],
[132, 21], [124, 21], [120, 9], [128, 9], [136, 9], [148, 9], [162, 9], [156, 25], [172, 21],
[180, 21], [180, 29], [172, 29], [172, 37], [172, 45], [180, 45], [180, 37], [188, 41], [196, 49],
[204, 57], [212, 65], [220, 73], [228, 69], [228, 77], [236, 77], [236, 69], [236, 61], [228, 61],
[228, 53], [236, 53], [236, 45], [228, 45], [228, 37], [236, 37], [236, 29], [228, 29], [228, 21],
[236, 21], [252, 21], [260, 29], [260, 37], [260, 45], [260, 53], [260, 61], [260, 69], [260, 77],
[276, 77], [276, 69], [276, 61], [276, 53], [284, 53], [284, 61], [284, 69], [284, 77], [284, 85],
[284, 93], [284, 101], [288, 109], [280, 109], [276, 101], [276, 93], [276, 85], [268, 97],
[260, 109], [252, 101], [260, 93], [260, 85], [236, 85], [228, 85], [228, 93], [236, 93],
[236, 101], [228, 101], [228, 109], [228, 117], [228, 125], [220, 125], [212, 117], [204, 109],
[196, 101], [188, 93], [180, 93], [180, 101], [180, 109], [180, 117], [180, 125], [196, 145],
[204, 145], [212, 145], [220, 145], [228, 145], [236, 145], [246, 141], [252, 125], [260, 129],
[280, 133]]

  return locations

if __name__ == '__main__':
  main()

Changing the search strategy

Since TSPs are computationally intractable, the routing solver does not always find the optimal solution to a problem. For instance, the solution returned in the previous example is not the optimal route. To find a better solution, you can use a more advanced search strategy, called guided local search, which enables the solver to escape a local minimum—a solution that is shorter than all nearby routes, but which is not the global minumum. After moving away from the local minimum, the solver continues the search. Search strategies like this are called local search metaheuristics.

To set a guided local search, add the following option to search_parameters in the main method of the program:

search_parameters.local_search_metaheuristic = (
        routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH)
Note: When using GUIDED_LOCAL_SEARCH or other metaheuristics, you need to set a time limit for the solver—otherwise the solver will not terminate. The following code sets a time limit of 30000 milliseconds, or 30 seconds.

search_parameters.time_limit_ms = 30000

When you run the program after making these changes, you get the following solution, which is shorter than the previous solution.

Total distance: 2659

Route:

0 -> 3 -> 276 -> 4 -> 5 -> 6 -> 8 -> 7 -> 9 -> 10 -> 11 -> 14 -> 12 -> 13 -> 23 -> 22 -> 24 -> 21 ->
25 -> 26 -> 27 -> 28 -> 125 -> 126 -> 127 -> 20 -> 19 -> 130 -> 129 -> 128 -> 153 -> 154 -> 152 ->
155 -> 151 -> 150 -> 177 -> 176 -> 175 -> 180 -> 161 -> 160 -> 174 -> 159 -> 158 -> 157 -> 156 ->
118 -> 119 -> 120 -> 121 -> 122 -> 123 -> 124 -> 29 -> 30 -> 31 -> 32 -> 33 -> 34 -> 35 -> 36 ->
37 -> 38 -> 39 -> 40 -> 41 -> 42 -> 59 -> 60 -> 58 -> 43 -> 44 -> 45 -> 46 -> 47 -> 48 -> 49 ->
50 -> 51 -> 52 -> 53 -> 54 -> 55 -> 56 -> 57 -> 67 -> 68 -> 66 -> 69 -> 70 -> 71 -> 72 -> 73 ->
75 -> 74 -> 76 -> 77 -> 78 -> 80 -> 81 -> 88 -> 79 -> 92 -> 93 -> 94 -> 95 -> 96 -> 97 -> 98 ->
99 -> 100 -> 101 -> 102 -> 91 -> 90 -> 89 -> 108 -> 111 -> 87 -> 82 -> 83 -> 86 -> 112 -> 115 ->
85 -> 84 -> 64 -> 65 -> 63 -> 62 -> 61 -> 117 -> 116 -> 114 -> 113 -> 110 -> 109 -> 107 -> 103 ->
104 -> 105 -> 106 -> 173 -> 172 -> 171 -> 170 -> 169 -> 168 -> 167 -> 166 -> 165 -> 164 -> 163 ->
162 -> 187 -> 188 -> 189 -> 190 -> 191 -> 192 -> 185 -> 186 -> 184 -> 183 -> 182 -> 181 -> 179 ->
178 -> 149 -> 148 -> 138 -> 137 -> 136 -> 266 -> 267 -> 135 -> 134 -> 268 -> 269 -> 133 -> 132 ->
131 -> 18 -> 17 -> 16 -> 15 -> 270 -> 271 -> 272 -> 273 -> 274 -> 275 -> 259 -> 258 -> 260 -> 261 ->
262 -> 263 -> 264 -> 265 -> 139 -> 140 -> 147 -> 146 -> 141 -> 142 -> 145 -> 144 -> 198 -> 197 ->
196 -> 193 -> 194 -> 195 -> 200 -> 201 -> 199 -> 143 -> 202 -> 203 -> 204 -> 205 -> 206 -> 207 ->
252 -> 253 -> 256 -> 257 -> 255 -> 254 -> 251 -> 208 -> 209 -> 210 -> 211 -> 212 -> 213 -> 214 ->
215 -> 216 -> 217 -> 218 -> 219 -> 220 -> 221 -> 222 -> 223 -> 224 -> 225 -> 226 -> 227 -> 232 ->
233 -> 234 -> 235 -> 236 -> 237 -> 230 -> 231 -> 228 -> 229 -> 250 -> 245 -> 238 -> 239 -> 240 ->
241 -> 242 -> 243 -> 244 -> 246 -> 249 -> 248 -> 247 -> 277 -> 278 -> 2 -> 279 -> 1 -> 0

For more search options, see Routing Options.

from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import routing_enums_pb2

# Distance callback
class CreateDistanceCallback(object):
  """Create callback to calculate distances between points."""
  def __init__(self):
    """Array of distances between points."""

    self.matrix = [
    [   0, 2451,  713, 1018, 1631, 1374, 2408,  213, 2571,  875, 1420, 2145, 1972], # New York
    [2451,    0, 1745, 1524,  831, 1240,  959, 2596,  403, 1589, 1374,  357,  579], # Los Angeles
    [ 713, 1745,    0,  355,  920,  803, 1737,  851, 1858,  262,  940, 1453, 1260], # Chicago
    [1018, 1524,  355,    0,  700,  862, 1395, 1123, 1584,  466, 1056, 1280,  987], # Minneapolis
    [1631,  831,  920,  700,    0,  663, 1021, 1769,  949,  796,  879,  586,  371], # Denver
    [1374, 1240,  803,  862,  663,    0, 1681, 1551, 1765,  547,  225,  887,  999], # Dallas
    [2408,  959, 1737, 1395, 1021, 1681,    0, 2493,  678, 1724, 1891, 1114,  701], # Seattle
    [ 213, 2596,  851, 1123, 1769, 1551, 2493,    0, 2699, 1038, 1605, 2300, 2099], # Boston
    [2571,  403, 1858, 1584,  949, 1765,  678, 2699,    0, 1744, 1645,  653,  600], # San Francisco
    [ 875, 1589,  262,  466,  796,  547, 1724, 1038, 1744,    0,  679, 1272, 1162], # St. Louis
    [1420, 1374,  940, 1056,  879,  225, 1891, 1605, 1645,  679,    0, 1017, 1200], # Houston
    [2145,  357, 1453, 1280,  586,  887, 1114, 2300,  653, 1272, 1017,    0,  504], # Phoenix
    [1972,  579, 1260,  987,  371,  999,  701, 2099,  600, 1162,  1200,  504,   0]] # Salt Lake City

  def Distance(self, from_node, to_node):
    return int(self.matrix[from_node][to_node])
def main():

  # Cities
  city_names = ["New York", "Los Angeles", "Chicago", "Minneapolis", "Denver", "Dallas", "Seattle",
                "Boston", "San Francisco", "St. Louis", "Houston", "Phoenix", "Salt Lake City"]
  tsp_size = len(city_names)
  num_routes = 1    # The number of routes, which is 1 in the TSP.
  # Nodes are indexed from 0 to tsp_size - 1. The depot is the starting node of the route.
  depot = 0

  # Create routing model
  if tsp_size > 0:
    routing = pywrapcp.RoutingModel(tsp_size, num_routes, depot)
    search_parameters = pywrapcp.DefaultRoutingSearchParameters()

    # Create the distance callback, which takes two arguments (the from and to node indices)
    # and returns the distance between these nodes.
    dist_between_nodes = CreateDistanceCallback()
    dist_callback = dist_between_nodes.Distance
    routing.SetArcCostEvaluatorOfAllVehicles(dist_callback)
    # Solve, returns a solution if any.
    assignment = routing.SolveWithParameters(search_parameters)
    if assignment:
      # Solution cost.
      print "Total distance: " + str(assignment.ObjectiveValue()) + " miles\n"
      # Inspect solution.
      # Only one route here; otherwise iterate from 0 to routing.vehicles() - 1
      route_number = 0
      index = routing.Start(route_number) # Index of the variable for the starting node.
      route = ''
      while not routing.IsEnd(index):
        # Convert variable indices to node indices in the displayed route.
        route += str(city_names[routing.IndexToNode(index)]) + ' -> '
        index = assignment.Value(routing.NextVar(index))
      route += str(city_names[routing.IndexToNode(index)])
      print "Route:\n\n" + route
    else:
      print 'No solution found.'
  else:
    print 'Specify an instance greater than 0.'

if __name__ == '__main__':
  main()

The best algorithms can now routinely solve TSP instances with tens of thousands of nodes. (The record at the time of writing is the pla85900 instance in TSPLIB, a VLSI application with 85,900 nodes. For certain instances with millions of nodes, solutions have been found guaranteed to be within 1% of an optimal tour.)

Entire program

The entire program using GUIDED_LOCAL_SEARCH shown below.

import math

from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import routing_enums_pb2
def distance(x1, y1, x2, y2):

    dist = math.sqrt((x1 - x2)**2 + (y1 - y2)**2)

    return dist

# Distance callback

class CreateDistanceCallback(object):
  """Create callback to calculate distances between points."""

  def __init__(self):
    """Initialize distance array."""
    locations = create_data_array()
    size = len(locations)
    self.matrix = {}

    for from_node in xrange(size):
      self.matrix[from_node] = {}
      for to_node in xrange(size):
        if from_node == to_node:
          self.matrix[from_node][to_node] = 0
        else:
          x1 = locations[from_node][0]
          y1 = locations[from_node][1]
          x2 = locations[to_node][0]
          y2 = locations[to_node][1]
          self.matrix[from_node][to_node] = distance(x1, y1, x2, y2)

  def Distance(self, from_node, to_node):
    return int(self.matrix[from_node][to_node])

def main():
  # Create the data.
  locations = create_data_array()
  tsp_size = len(locations)
  num_routes = 1      # The number of routes, which is 1 in the TSP.
  # Nodes are indexed from 0 to tsp_size - 1. The depot is the starting node of the route.
  depot = 0
  # Create routing model.
  if tsp_size > 0:
    routing = pywrapcp.RoutingModel(tsp_size, num_routes, depot)
    search_parameters = pywrapcp.DefaultRoutingSearchParameters()
    search_parameters.local_search_metaheuristic = (
        routing_enums_pb2.LocalSearchMetaheuristic.GUIDED_LOCAL_SEARCH)
    # Set search time limit in milliseconds.
    search_parameters.time_limit_ms = 30000
    # Callback to the distance function. The callback takes two
    # arguments (the from and to node indices) and returns the distance between them.
    dist_between_locations = CreateDistanceCallback()
    dist_callback = dist_between_locations.Distance
    routing.SetArcCostEvaluatorOfAllVehicles(dist_callback)
    # Solve, returns a solution if any.
    assignment = routing.SolveWithParameters(search_parameters)
    if assignment:

      # Solution cost.
      print "Total distance: " + str(assignment.ObjectiveValue()) + "\n"

      # Inspect solution.
      # Only one route here; otherwise iterate from 0 to routing.vehicles() - 1.
      route_number = 0
      node = routing.Start(route_number)
      start_node = node
      route = ''

      while not routing.IsEnd(node):
        route += str(node) + ' -> '
        node = assignment.Value(routing.NextVar(node))
      route += '0'
      print "Route:\n\n" + route
    else:
      print 'No solution found.'
  else:
    print 'Specify an instance greater than 0.'
def create_data_array():
  locations = [[288, 149], [288, 129], [270, 133], [256, 141], [256, 157], [246, 157], [236, 169],
[228, 169], [228, 161], [220, 169], [212, 169], [204, 169], [196, 169], [188, 169], [196, 161],
[188, 145], [172, 145], [164, 145], [156, 145], [148, 145], [140, 145], [148, 169], [164, 169],
[172, 169], [156, 169], [140, 169], [132, 169], [124, 169], [116, 161], [104, 153], [104, 161],
[104, 169], [90, 165], [80, 157], [64, 157], [64, 165], [56, 169], [56, 161], [56, 153], [56, 145],
[56, 137], [56, 129], [56, 121], [40, 121], [40, 129], [40, 137], [40, 145], [40, 153], [40, 161],
[40, 169], [32, 169], [32, 161], [32, 153], [32, 145], [32, 137], [32, 129], [32, 121], [32, 113],
[40, 113], [56, 113], [56, 105], [48, 99], [40, 99], [32, 97], [32, 89], [24, 89], [16, 97],
[16, 109], [8, 109], [8, 97], [8, 89], [8, 81], [8, 73], [8, 65], [8, 57], [16, 57], [8, 49],
[8, 41], [24, 45], [32, 41], [32, 49], [32, 57], [32, 65], [32, 73], [32, 81], [40, 83], [40, 73],
[40, 63], [40, 51], [44, 43], [44, 35], [44, 27], [32, 25], [24, 25], [16, 25], [16, 17], [24, 17],
[32, 17], [44, 11], [56, 9], [56, 17], [56, 25], [56, 33], [56, 41], [64, 41], [72, 41], [72, 49],
[56, 49], [48, 51], [56, 57], [56, 65], [48, 63], [48, 73], [56, 73], [56, 81], [48, 83], [56, 89],
[56, 97], [104, 97], [104, 105], [104, 113], [104, 121], [104, 129], [104, 137], [104, 145],
[116, 145], [124, 145], [132, 145], [132, 137], [140, 137], [148, 137], [156, 137], [164, 137],
[172, 125], [172, 117], [172, 109], [172, 101], [172, 93], [172, 85], [180, 85], [180, 77],
[180, 69], [180, 61], [180, 53], [172, 53], [172, 61], [172, 69], [172, 77], [164, 81], [148, 85],
[124, 85], [124, 93], [124, 109], [124, 125], [124, 117], [124, 101], [104, 89], [104, 81],
[104, 73], [104, 65], [104, 49], [104, 41], [104, 33], [104, 25], [104, 17], [92, 9], [80, 9],
[72, 9], [64, 21], [72, 25], [80, 25], [80, 25], [80, 41], [88, 49], [104, 57], [124, 69],
[124, 77], [132, 81], [140, 65], [132, 61], [124, 61], [124, 53], [124, 45], [124, 37], [124, 29],
[132, 21], [124, 21], [120, 9], [128, 9], [136, 9], [148, 9], [162, 9], [156, 25], [172, 21],
[180, 21], [180, 29], [172, 29], [172, 37], [172, 45], [180, 45], [180, 37], [188, 41], [196, 49],
[204, 57], [212, 65], [220, 73], [228, 69], [228, 77], [236, 77], [236, 69], [236, 61], [228, 61],
[228, 53], [236, 53], [236, 45], [228, 45], [228, 37], [236, 37], [236, 29], [228, 29], [228, 21],
[236, 21], [252, 21], [260, 29], [260, 37], [260, 45], [260, 53], [260, 61], [260, 69], [260, 77],
[276, 77], [276, 69], [276, 61], [276, 53], [284, 53], [284, 61], [284, 69], [284, 77], [284, 85],
[284, 93], [284, 101], [288, 109], [280, 109], [276, 101], [276, 93], [276, 85], [268, 97],
[260, 109], [252, 101], [260, 93], [260, 85], [236, 85], [228, 85], [228, 93], [236, 93],
[236, 101], [228, 101], [228, 109], [228, 117], [228, 125], [220, 125], [212, 117], [204, 109],
[196, 101], [188, 93], [180, 93], [180, 101], [180, 109], [180, 117], [180, 125], [196, 145],
[204, 145], [212, 145], [220, 145], [228, 145], [236, 145], [246, 141], [252, 125], [260, 129],
[280, 133]]

  return locations

if __name__ == '__main__':
  main()

Send feedback about...

Optimization
Optimization