The traveling salesman problem is super easy to understand and (partly because of the problem’s simplicity) is a staple in most computer science algorithms courses.  The problem statement says:  Given the task of traveling through n cities and the restriction that each city must be visited exactly once, can we find a travel plan that minimizes the total distance traveled?

Below is a recently written algorithm that takes a natural number, n, for an argument and then solves an arbitrary traveling salesman problem with n nodes inside a unit square.  The algorithm then outputs the cities (ordered pairs) in the correct order and produces a visualization with a directed graph of the results.  The optimization mechanism for this result is a Monte Carlo Markov Chain with a fair amount of tuning (convergence speed vs absolute precision).

set.seed(321)

travel<-function(n){

#here we specify the number of points from the function argument
	nsim=n
	
#here we generate n random points in the unit square
	x_values<-c(rep(0,nsim))
	y_values<-c(rep(0,nsim))
	cities<-c(1:nsim)
	for(i in 1:nsim){x_values[i]=runif(1,0,1)}
	for(i in 1:nsim){y_values[i]=runif(1,0,1)}
	
#here we generate the distance matrix used to minimize the route
	distance<-matrix(0,nsim,nsim)
	candidate<-list()
	t=1
	s=1
	route_distance<-c()
	for(i in 1:nsim){
		for(j in 1:nsim){
			distance[i,j]=sqrt((x_values[i]-x_values[j])^2+(y_values[i]-
			y_values[j])^2)
		}
	}
	route<-sample(cities,nsim)


#here we repeatedly propose and accept/reject possible routes in order to arrive at a 
#global maximum.
	for(k in 1:100000){
		if(k>1000){s=10}
		if(k>2000){s=20}
		if(k>3000){s=50}
		if(k>4000){s=60}
		if(k>5000){s=70}
		if(k>6000){s=80}
		if(k>7000){s=90}
		if(k>8000){s=120}
#proposal step
	route_new<-route
	choose<-sample(route,2,replace=FALSE)
	route_new[choose[1]]=route[choose[2]]
	route_new[choose[2]]=route[choose[1]]
	
#acceptance/rejection step
	sumvector_TOP<-c(rep(0,nsim))
		for(i in 1:nsim-1){sumvector_TOP[i]=distance[route_new[i],route_new[i+1]]}
	sumvector_BOTTOM<-c(rep(0,nsim))
		for(i in 1:nsim-1){sumvector_BOTTOM[i]=distance[route[i],route[i+1]]}
		if(rbinom(1,1,min(1,(sum(sumvector_TOP)/sum(sumvector_BOTTOM))^(-s)))==1){
		route=route_new;
		candidate[[t]]<-route_new;
		route_distance[t]=sum(sumvector_TOP);
		t=t+1
		}
	}
#here we create a plot of the random points and connect them with directed edges
#the function requires the igraph package
best<-candidate[[which.min(route_distance)]]
output<-cbind(cities,x_values,y_values)
solution<-cbind(x_values,y_values)
best_odd<-c()
best_even<-c()
z<-c()

library(igraph)
layout<-solution
best_odd<-best
for(i in 2:nsim){
	best_even[i-1]<-best[i]
}
best_even[nsim]=best[1]
z=as.vector(rbind(best_odd,best_even))
G<-graph(z)
plot(G,layout=layout)

cat("the total distance for the shortest path is: ",min(route_distance),"\n")
cat("the route that generates this distance is: ",
candidate[[which.min(route_distance)]],"\n")
cat("the random points generated were: ","\n")
return(output)
}

#call the function
travel(20)

#the total distance for the shortest path is:  3.431226 
#the route that generates this distance is:  18 9 11 20 15 17 13 10 1 2 7 3 14 6 
#12 19 5 8 4 16 
#the random points generated were:  
#      cities   x_values   y_values
# [1,]      1 0.95589376 0.63596776
# [2,]      2 0.93728552 0.98987845
# [3,]      3 0.23822045 0.93117401
# [4,]      4 0.25507364 0.48550475
# [5,]      5 0.39051197 0.57468926
# [6,]      6 0.34117986 0.75123825
# [7,]      7 0.45238060 0.99272238
# [8,]      8 0.28993282 0.43086438
# [9,]      9 0.45067322 0.12436502
#[10,]     10 0.80659572 0.59426029
#[11,]     11 0.60572107 0.21050458
#[12,]     12 0.36284109 0.71025331
#[13,]     13 0.76633847 0.67582987
#[14,]     14 0.04508695 0.99290506
#[15,]     15 0.59217188 0.66287874
#[16,]     16 0.20163219 0.43870601
#[17,]     17 0.63273803 0.71822766
#[18,]     18 0.40358999 0.02558309
#[19,]     19 0.29050341 0.67232996
#[20,]     20 0.64045901 0.54095314

#   user  system elapsed 
#   2.248   0.036   2.961