# Code from Chapter 9 of Machine Learning: An Algorithmic Perspective (2nd Edition) # by Stephen Marsland (http://stephenmonika.net) # You are free to use, change, or redistribute the code in any way you wish for # non-commercial purposes, but please maintain the name of the original author. # This code comes with no warranty of any kind. # Stephen Marsland, 2008, 2014 # A demonstration of four methods of solving the Travelling Salesman Problem import numpy as np def makeTSP(nCities): positions = 2*np.random.rand(nCities,2)-1; distances = np.zeros((nCities,nCities)) for i in range(nCities): for j in range(i+1,nCities): distances[i,j] = np.sqrt((positions[i,0] - positions[j,0])**2 + (positions[i,1] - positions[j,1])**2); distances[j,i] = distances[i,j]; return distances def exhaustive(distances): nCities = np.shape(distances)[0] cityOrder = np.arange(nCities) distanceTravelled = 0 for i in range(nCities-1): distanceTravelled += distances[cityOrder[i],cityOrder[i+1]] distanceTravelled += distances[cityOrder[nCities-1],0] for newOrder in permutation(range(nCities)): possibleDistanceTravelled = 0 for i in range(nCities-1): possibleDistanceTravelled += distances[newOrder[i],newOrder[i+1]] possibleDistanceTravelled += distances[newOrder[nCities-1],0] if possibleDistanceTravelled < distanceTravelled: distanceTravelled = possibleDistanceTravelled cityOrder = newOrder return cityOrder, distanceTravelled def permutation(order): order = tuple(order) if len(order)==1: yield order else: for i in range(len(order)): rest = order[:i] + order[i+1:] move = (order[i],) for smaller in permutation(rest): yield move + smaller def greedy(distances): nCities = np.shape(distances)[0] distanceTravelled = 0 # Need a version of the matrix we can trash dist = distances.copy() cityOrder = np.zeros(nCities) cityOrder[0] = np.random.randint(nCities) dist[:,cityOrder[0]] = np.Inf for i in range(nCities-1): cityOrder[i+1] = np.argmin(dist[cityOrder[i],:]) distanceTravelled += dist[cityOrder[i],cityOrder[i+1]] # Now exclude the chance of travelling to that city again dist[:,cityOrder[i+1]] = np.Inf # Now return to the original city distanceTravelled += distances[cityOrder[nCities-1],0] return cityOrder, distanceTravelled def hillClimbing(distances): nCities = np.shape(distances)[0] cityOrder = np.arange(nCities) np.random.shuffle(cityOrder) distanceTravelled = 0 for i in range(nCities-1): distanceTravelled += distances[cityOrder[i],cityOrder[i+1]] distanceTravelled += distances[cityOrder[nCities-1],0] for i in range(1000): # Choose cities to swap city1 = np.random.randint(nCities) city2 = np.random.randint(nCities) if city1 != city2: # Reorder the set of cities possibleCityOrder = cityOrder.copy() possibleCityOrder = np.where(possibleCityOrder==city1,-1,possibleCityOrder) possibleCityOrder = np.where(possibleCityOrder==city2,city1,possibleCityOrder) possibleCityOrder = np.where(possibleCityOrder==-1,city2,possibleCityOrder) # Work out the new distances # This can be done more efficiently newDistanceTravelled = 0 for j in range(nCities-1): newDistanceTravelled += distances[possibleCityOrder[j],possibleCityOrder[j+1]] distanceTravelled += distances[cityOrder[nCities-1],0] if newDistanceTravelled < distanceTravelled: distanceTravelled = newDistanceTravelled cityOrder = possibleCityOrder return cityOrder, distanceTravelled def simulatedAnnealing(distances): nCities = np.shape(distances)[0] cityOrder = np.arange(nCities) np.random.shuffle(cityOrder) distanceTravelled = 0 for i in range(nCities-1): distanceTravelled += distances[cityOrder[i],cityOrder[i+1]] distanceTravelled += distances[cityOrder[nCities-1],0] T = 500 c = 0.8 nTests = 10 while T>1: for i in range(nTests): # Choose cities to swap city1 = np.random.randint(nCities) city2 = np.random.randint(nCities) if city1 != city2: # Reorder the set of cities possibleCityOrder = cityOrder.copy() possibleCityOrder = np.where(possibleCityOrder==city1,-1,possibleCityOrder) possibleCityOrder = np.where(possibleCityOrder==city2,city1,possibleCityOrder) possibleCityOrder = np.where(possibleCityOrder==-1,city2,possibleCityOrder) # Work out the new distances # This can be done more efficiently newDistanceTravelled = 0 for j in range(nCities-1): newDistanceTravelled += distances[possibleCityOrder[j],possibleCityOrder[j+1]] distanceTravelled += distances[cityOrder[nCities-1],0] if newDistanceTravelled < distanceTravelled or (distanceTravelled - newDistanceTravelled) < T*np.log(np.random.rand()): distanceTravelled = newDistanceTravelled cityOrder = possibleCityOrder # Annealing schedule T = c*T return cityOrder, distanceTravelled def runAll(): import time nCities = 5 distances = makeTSP(nCities) print "Exhaustive search" start = time.time() print exhaustive(distances) finish = time.time() print finish-start print "Greedy search" start = time.time() print greedy(distances) finish = time.time() print finish-start print "Hill Climbing" start = time.time() print hillClimbing(distances) finish = time.time() print finish-start print "Simulated Annealing" start = time.time() print simulatedAnnealing(distances) finish = time.time() print finish-start runAll()