import sys import copy import math import time from heapq import heappush, heappop from random import shuffle, sample, randint, SystemRandom import matplotlib.pyplot as plt from matplotlib.backends.backend_pdf import PdfPages # globals to contain reference solutions neos14 = [318, 324, 336, 319, 351, 311, 272, 361, 274, 322] neos15 = [313, 318, 281, 324, 378, 291, 348, 342, 353, 325] neos16 = [404, 353, 361, 349, 358, 344, 373, 355, 243, 330] neos36 = 464 neos = {14: neos14, 15: neos15, 16: neos16, 36: neos36} # returns a shuffled version of the cities list while maintaining the first & last elements def rand_order(cities): # shitty variable name suburbs = cities[1:-1] shuffle(suburbs) suburbs.insert(0, cities[0]) suburbs.append(cities[-1]) return suburbs # returns distance between 2 vertices def dist(city1, city2): return math.sqrt(int((city1[2]) - int(city2[2])) ** 2 + (int(city1[3]) - int(city2[3])) ** 2) # returns total cost/distance of a tour (list of cities) def length(cities): sum = 0 for i in range(len(cities) - 1): sum += dist(cities[i], cities[i+1]) return sum # returns the best ordered child state of a state (hill climbing, basic) & its cost def best_child(cities, cost): bestorder = cities bestcost = cost for i in range(1, len(cities) - 2): for j in range(i, len(cities) - 1): temporder = copy.deepcopy(cities) tempcity = temporder[i] temporder[i] = temporder[j] temporder[j] = tempcity templen = length(temporder) if (templen < cost): bestorder = temporder bestcost = templen return (bestorder, bestcost) # returns the best child using tabu # returns the best ordered child state of a state (hill climbing, basic) & its cost def bester_child(cities, cost, tabu): bestorder = cities bestcost = cost for i in range(1, len(cities) - 2): for j in range(i, len(cities) - 1): temporder = copy.deepcopy(cities) tempcity = temporder[i] temporder[i] = temporder[j] temporder[j] = tempcity if (get_hashable(temporder) in tabu): continue templen = length(temporder) if (templen < cost): bestorder = temporder bestcost = templen return (bestorder, bestcost) # returns a random child for simulated annealing search # returns the ordering and the cost of it def random_child(cities): #i = randbelow(len(cities) - 3) + 1 i = int(randint(1, len(cities) - 2)) j = i while (j == i): #j = randbelow(len(cities) - 3) + 1 j = int(randint(1, len(cities) - 2)) temporder = copy.deepcopy(cities) tempcity = temporder[i] temporder[i] = temporder[j] temporder[j] = tempcity return(temporder, length(temporder)) # returns the basic hill climbing solution for a given input def sir_edmund_hillary(cities): curcost = length(cities) curpath = copy.deepcopy(cities) steps = 0 # keep trying things til we can no longer try things while(1): steps += 1 child = best_child(curpath, curcost) if (child[1] < curcost): curcost = child[1] curpath = child[0] continue break return(curpath, curcost, steps) # retuns a tuple of the indices (ordered) in a list of cities def get_hashable(cities): newlist = [] for i in range(1, len(cities) - 1): newlist.append(cities[i][1]) return tuple(newlist) # returns the hill climbing solution w/tabu & sideways moves for a given input def tenzing_norgay(cities): curcost = length(cities) curpath = copy.deepcopy(cities) tabu = {get_hashable(cities)} laterals = 0 steps = 0 # keep trying things til we can no longer try things while(1): steps += 1 child = bester_child(curpath, curcost, tabu) if (child[1] < curcost): curcost = child[1] curpath = child[0] tabu.add(get_hashable(curpath)) continue if (child[1] == curcost and laterals < len(cities)): laterals += 1 curcost = child[1] curpath = child[0] tabu.add(get_hashable(curpath)) continue break return(curpath, curcost, steps) # returns the hill climbing solution w/ x random restarts for a given input # if you want x restarts, call with x + 1 def reinhold_messner(cities, x): start = time.process_time() curcost = length(cities) curpath = copy.deepcopy(cities) bestcost = length(cities) bestpath = copy.deepcopy(cities) # keep trying things til we can no longer try things while(x): child = best_child(curpath, curcost) if (child[1] < curcost): curcost = child[1] curpath = child[0] continue else: if (curcost < bestcost): bestcost = curcost bestpath = copy.deepcopy(curpath) curpath = rand_order(cities) curcost = length(curpath) x -= 1 end = time.process_time() return(bestpath, bestcost, end - start) # returns true or false randomly, based on the distribution e^(delta / temp) def bad_weather(delta, temp): return SystemRandom().random() < math.exp(delta / temp) # decrements temp by the option specified by opt: 1 is linear, 2 is logarithmic, 3 is exponential def colder(temp, opt, count): if (opt == 1): return temp - 15 if (opt == 2): return temp - (math.log(temp) + 2) else: return 500000 * math.exp(-0.0005 * count) - 1 # returns the hill climbing solution using simulated annealing w/ schedule opt def beck_weathers(cities, opt): start = time.process_time() temp = 500000 count = 1 curcost = length(cities) curpath = copy.deepcopy(cities) # keep trying things til we can no longer try things while(temp > 0): child = random_child(curpath) delta = curcost - child[1] if (delta > 0 or bad_weather(delta, temp)): curcost = child[1] curpath = child[0] temp = colder(temp, opt, count) count += 1 end = time.process_time() return(curpath, curcost, end - start) # solves the inputted TSP problem (filepath), using algorithm algo # 1 = hill climbing, 2 = tabu/sideways allowed, 3 = random restarts, 4 = annealing # x represents the optional parameters for 3 and 4 # specifically, the # of restarts and the annealing schedule, respectively def solver(problem, algo, x=None): # import map prob = problem with open(prob) as file: prob = file.read().splitlines() global n cities = prob[1:] counter = 0 # convert to list of list of ints # original format is: letter, coord, coord # output format is: iD#, letter, coord, coord for l in cities: temp = l.split() temp[1] = int(temp[1]) temp[2] = int(temp[2]) temp.insert(0, counter) counter += 1 cities[cities.index(l)] = temp # add in goal state cities.append(cities[0]) input = rand_order(cities) if (algo == 1): result = sir_edmund_hillary(cities) elif (algo == 2): result = tenzing_norgay(cities) elif (algo == 3): result = reinhold_messner(cities, x) else: result = beck_weathers(cities, x) return result # function to generate PDF results for basic hill climbing def sir_edmund_hillary_graph(): plt.ioff() plt.switch_backend('agg') solution_steps = [] solution_scores = [] good_sol_counts = [] for i in range(14, 17): steps = [] solution_steps.append(steps) scores = [] solution_scores.append(scores) good_sols = [] good_sol_counts.append(good_sols) for j in range(1, 11): filepath = "tsp_problems/" + str(i) + "/instance_" + str(j) + ".txt" prob_steps = 0 prob_scores = 0 good_solutions = 0 # run problem 100 times for k in range(0, 100): # path, cost, steps returned result = solver(filepath, 1) prob_steps += result[2] prob_scores += (result[1] / neos[i][j - 1]) if (result[1] <= neos[i][j- 1]): good_solutions += 1 steps.append(prob_steps / 100.0) scores.append(prob_scores / 100.0) good_sols.append(good_solutions) with PdfPages("hill_climbing.pdf") as pdf: figure, axes = plt.subplots(3, 1, True) figure.suptitle("Hill Climbing") for i in range(0, 3): axes[i].plot(range(1, 11), solution_steps[i], label = "Steps to Solution", color = 'b') axes[i].set_xlabel("Problem Instance w/ " + str(i + 14) + " cities") axes[i].set_ylabel("Steps Taken", color = 'b') axis_twin = axes[i].twinx() axis_twin.plot(range(1, 11), solution_scores[i], label = "Solution Quality", color = 'r') axis_twin.set_ylabel("Solution Quality", color = 'r') #axes[i].legend() #axis_twin.legend() plt.tight_layout() plt.subplots_adjust(top=0.92) pdf.savefig() plt.close() figure, axes = plt.subplots(3, 1, True) figure.suptitle("Hill Climbing Cont.") for i in range(0, 3): axes[i].plot(range(1, 11), good_sol_counts[i], color = 'b') axes[i].set_xlabel("Problem Instance w/ " + str(i + 14) + " cities") axes[i].set_ylabel("% of runs <= NEOS", color = 'b') plt.tight_layout() plt.subplots_adjust(top=0.92) pdf.savefig() plt.close() # now we can aggregate and replace the old values for i in range(0, 3): tempsteps = 0 tempqual = 0 tempgood = 0 for j in range(0, 10): tempsteps += solution_steps[i][j] tempqual += solution_scores[i][j] tempgood += good_sol_counts[i][j] solution_steps[i] = tempsteps / 10.0 solution_scores[i] = tempqual / 10.0 good_sol_counts[i] = tempgood / 10.0 figure, axes = plt.subplots(3, 1, True) figure.suptitle("Hill Climbing Aggegated (Averages)") axes[0].plot(range(14, 17), solution_steps, color = 'b') axes[0].set_xlabel("Problem Instances w/ x cities") axes[0].set_ylabel("Steps Used", color = 'b') axes[0].set_xticks([14, 15, 16]) axes[1].plot(range(14, 17), solution_scores, color = 'b') axes[1].set_xlabel("Problem Instances w/ x cities") axes[1].set_ylabel("Solution Quality", color = 'b') axes[1].set_xticks([14, 15, 16]) axes[2].plot(range(14, 17), good_sol_counts, color = 'b') axes[2].set_xlabel("Problem Instances w/ x cities") axes[2].set_ylabel("% of runs <= NEOS", color = 'b') axes[2].set_xticks([14, 15, 16]) plt.tight_layout() plt.subplots_adjust(top=0.92) pdf.savefig() plt.close() # function to generate PDF results for basic hill climbing w/ tabu and sideways moves def tenzing_norgay_graph(): plt.ioff() plt.switch_backend('agg') solution_steps = [] solution_scores = [] good_sol_counts = [] for i in range(14, 17): steps = [] solution_steps.append(steps) scores = [] solution_scores.append(scores) good_sols = [] good_sol_counts.append(good_sols) for j in range(1, 11): filepath = "tsp_problems/" + str(i) + "/instance_" + str(j) + ".txt" prob_steps = 0 prob_scores = 0 good_solutions = 0 # run problem 100 times for k in range(0, 100): # path, cost, steps returned result = solver(filepath, 2) prob_steps += result[2] prob_scores += (result[1] / neos[i][j - 1]) if (result[1] <= neos[i][j- 1]): good_solutions += 1 steps.append(prob_steps / 100.0) scores.append(prob_scores / 100.0) good_sols.append(good_solutions) with PdfPages("hill_climbing_tabu.pdf") as pdf: figure, axes = plt.subplots(3, 1, True) figure.suptitle("Hill Climbing - Tabu") for i in range(0, 3): axes[i].plot(range(1, 11), solution_steps[i], label = "Steps to Solution", color = 'b') axes[i].set_xlabel("Problem Instance w/ " + str(i + 14) + " cities") axes[i].set_ylabel("Steps Taken", color = 'b') axis_twin = axes[i].twinx() axis_twin.plot(range(1, 11), solution_scores[i], label = "Solution Quality", color = 'r') axis_twin.set_ylabel("Solution Quality", color = 'r') #axes[i].legend() #axis_twin.legend() plt.tight_layout() plt.subplots_adjust(top=0.92) pdf.savefig() plt.close() figure, axes = plt.subplots(3, 1, True) figure.suptitle("Hill Climbing - Tabu Cont.") for i in range(0, 3): axes[i].plot(range(1, 11), good_sol_counts[i], color = 'b') axes[i].set_xlabel("Problem Instance w/ " + str(i + 14) + " cities") axes[i].set_ylabel("% of runs <= NEOS", color = 'b') plt.tight_layout() plt.subplots_adjust(top=0.92) pdf.savefig() plt.close() # now we can aggregate and replace the old values for i in range(0, 3): tempsteps = 0 tempqual = 0 tempgood = 0 for j in range(0, 10): tempsteps += solution_steps[i][j] tempqual += solution_scores[i][j] tempgood += good_sol_counts[i][j] solution_steps[i] = tempsteps / 10.0 solution_scores[i] = tempqual / 10.0 good_sol_counts[i] = tempgood / 10.0 figure, axes = plt.subplots(3, 1, True) figure.suptitle("Hill Climbing - Tabu Aggegated (Averages)") axes[0].plot(range(14, 17), solution_steps, color = 'b') axes[0].set_xlabel("Problem Instances w/ x cities") axes[0].set_ylabel("Steps Used", color = 'b') axes[0].set_xticks([14, 15, 16]) axes[1].plot(range(14, 17), solution_scores, color = 'b') axes[1].set_xlabel("Problem Instances w/ x cities") axes[1].set_ylabel("Solution Quality", color = 'b') axes[1].set_xticks([14, 15, 16]) axes[2].plot(range(14, 17), good_sol_counts, color = 'b') axes[2].set_xlabel("Problem Instances w/ x cities") axes[2].set_ylabel("% of runs <= NEOS", color = 'b') axes[2].set_xticks([14, 15, 16]) plt.tight_layout() plt.subplots_adjust(top=0.92) pdf.savefig() plt.close() # function to generate PDF results for basic hill climbing w/ random restarts def reinhold_messner_graph(): plt.ioff() plt.switch_backend('agg') solution_times = [] solution_scores = [] # num of cities for i in range(14, 17): times = [] solution_times.append(times) scores = [] solution_scores.append(scores) # num of instances for j in range(1, 3): filepath = "tsp_problems/" + str(i) + "/instance_" + str(j) + ".txt" prob_times = [] times.append(prob_times) prob_scores = [] scores.append(prob_scores) # num of repeats for k in range(1, 16): ktimes = [] prob_times.append(ktimes) kscores = [] prob_scores.append(kscores) runtime = 0 qual = 0 # run problem 100 times for num in range(0, 100): # path, cost, runtime returned result = solver(filepath, 3, k) runtime += result[2] qual += (result[1] / neos[i][j - 1]) ktimes.append(runtime / 100.0) kscores.append(qual / 100.0) with PdfPages("hill_climbing_restarts.pdf") as pdf: figure, axes = plt.subplots(3, 2, True) figure.suptitle("Hill Climbing - Restarts") for i in range(0, 3): for j in range(0, 2): axes[i][j].plot(range(1, 16), solution_times[i][j], color = 'b') axes[i][j].set_xlabel("Problem Instance " + str(j + 1) + " w/ " + str(i + 14) + " cities") axes[i][j].set_ylabel("Time Taken (s)", color = 'b') axis_twin = axes[i][j].twinx() axis_twin.plot(range(1, 16), solution_scores[i][j], color = 'r') axis_twin.set_ylabel("Solution Quality", color = 'r') #axes[i].legend() #axis_twin.legend() plt.tight_layout() plt.subplots_adjust(top=0.92) pdf.savefig() plt.close() # function to generate PDF results for basic hill climbing w/ annealing def beck_weathers_graph(): plt.ioff() plt.switch_backend('agg') solution_times = [] solution_scores = [] # num of cities for i in range(14, 17): times = [] solution_times.append(times) scores = [] solution_scores.append(scores) # num of instances for j in range(1, 3): filepath = "tsp_problems/" + str(i) + "/instance_" + str(j) + ".txt" prob_times = [] times.append(prob_times) prob_scores = [] scores.append(prob_scores) # diff schedules for k in range(1, 4): ktimes = [] prob_times.append(ktimes) kscores = [] prob_scores.append(kscores) runtime = 0 qual = 0 # run problem 100 times for num in range(0, 100): # path, cost, runtime returned result = solver(filepath, 4, k) runtime += result[2] qual += (result[1] / neos[i][j - 1]) ktimes.append(runtime / 100.0) kscores.append(qual / 100.0) with PdfPages("hill_climbing_annealing.pdf") as pdf: figure, axes = plt.subplots(3, 2, True) figure.suptitle("Hill Climbing - Annealing") for i in range(0, 3): for j in range(0, 2): axes[i][j].plot(range(1, 4), solution_times[i][j], color = 'b') axes[i][j].set_title("Problem Instance " + str(j + 1) + " w/ " + str(i + 14) + " cities") axes[i][j].set_xlabel("Annealing Schedule") axes[i][j].set_ylabel("Time Taken (s)", color = 'b') axis_twin = axes[i][j].twinx() axis_twin.plot(range(1, 4), solution_scores[i][j], color = 'r') axis_twin.set_ylabel("Solution Quality", color = 'r') #axes[i].legend() #axis_twin.legend() plt.tight_layout() plt.subplots_adjust(top=0.85) pdf.savefig() plt.close() # main function that will call needed graphing function(s) def main(): #sir_edmund_hillary_graph() #tenzing_norgay_graph() #reinhold_messner_graph() beck_weathers_graph() main()