dist.m

function [ d ] = dist( x1, y1, x2, y2 )
% This function return the Euclidian distance of two points

d = sqrt((x2-x1)^2 + (y2-y1)^2);
end

output.txt

>> TravellingSalesman
The starting distance is: 22203.04
The ending distance is: 7044.30

simple_path_dist.m

function [d] = simple_path_dist(locations)
d = 0;
for i = 2:size(locations)
    x1 = locations(i-1, 1);
    y1 = locations(i-1, 2);
    x2 = locations(i, 1);
    y2 = locations(i, 2);
    d = d + dist(x1, y1, x2, y2);
end
end

TravellingSalesman.m

% Credit
% Understanding how to plot the United Stats map (no help was received from this page on the simulated annealing):
%   - http://www.mathworks.com/help/optim/ug/travelling-salesman-problem.html
%
% ---- Initial variable setup ----

clear all;

temp = 100;
cooling_rate = 0.01;

stops = 50;
locations = zeros(stops,2);
distances = [];

% ---- Setting up initial plot ----

figure;
load('usborder.mat','x','y','xx','yy');

% (Everything is scaled by 1000 to simulate something closer to real size
% of the US)

x = x * 1000;
y = y * 1000;
xx = xx * 1000;
yy = yy * 1000;

n = 1;
while (n <= stops)
    xp = rand * 1.5 * 1000;
    yp = rand * 1000;
    if inpolygon(xp,yp,xx,yy)
        locations(n,1) = xp;
        locations(n,2) = yp;
        n = n+1;
    end
end

title('Travelling Salesman')

subplot(1, 3, 1);
title('Randomized Solution')
xlabel('Longitude')
ylabel('Latitude')
hold on
plot(x,y,'Color','red'); 
plot(locations(:,1),locations(:,2),'*b')
for i = 2:length(locations)
    plot([locations(i, 1) locations(i-1, 1)], [locations(i, 2) locations(i-1, 2)], 'k')
end

hold off

% ---- Starting distance ----

distance = simple_path_dist(locations);
fprintf('The starting distance is: %0.2f\n', distance);

% ---- Simulated Annealing Algorithm ----

while(temp >= 0)
    % Choose first random city
    first_random_city = randi([1 length(locations)]);
    
    % Ensure that the second random city != the first random city
    % In the vast majority of cases, the while loop will execute once.
    second_random_city = first_random_city;
    while(second_random_city == first_random_city)
        second_random_city = randi([1 length(locations)]);
    end
    
    % Initialize current distance/locations matrix
    curr_locations = locations;
    curr_locations([first_random_city, second_random_city],:) = curr_locations([second_random_city, first_random_city],:);

    % Calculate the 
    curr_distance = simple_path_dist(curr_locations);
    
    
    if (distance > curr_distance)
        % If our new solution is better than previous
        % automatically accept the solution
        
        distance = curr_distance;
        locations = curr_locations;
    elseif (exp((distance - curr_distance) / temp) > rand())
        % Else, if a random probabilty is greater than exp((distance -
        % curr_distance) / temp), accept the worse local solution in order
        % to ensure the best global solution.
        
        distance = curr_distance;
        locations = curr_locations;
    end
    
    % Decrease temperature
    temp = temp - cooling_rate;
    
    % Keep up with distance
    distances = [distances distance];
end

% ---- Plot final graph ----

subplot(1, 3, 2);
title('Optimized Solution')
xlabel('Longitude')
ylabel('Latitude')
hold on;
plot(x,y,'Color','red'); 
plot(locations(:,1),locations(:,2),'*b')
for i = 2:length(locations)
    plot([locations(i, 1) locations(i-1, 1)], [locations(i, 2) locations(i-1, 2)], 'k')
end


% ---- Plot distances graph ----

subplot(1, 3, 3);
title('Minimization of distance over iterations')
xlabel('Iterations')
ylabel('Distance (mi)')
hold on;
plot(distances);
hold off

% ---- Final distance ----
fprintf('The ending distance is: %0.2f\n\n', distance);