# ###################################################################
#
#       Program: tsp_recuit_simule.py
#        Author: Jean-Michel Richer
#  Organisation: Computer Science Department, 
#                University of Angers,
#	             France
#         Email: jean-michel.richer@univ-angers.fr
# Creation date: April, 2021
#  Modification: April, 2021
#
# ###################################################################
# 
# Aim:
#
#    This program is an implementation of a Simulated Annealing 
#    algorithm used to find a solution to the travelling salesman 
#    problem. 
#
# Objectif :
#
#    Ce programme repose sur l'implantation d'un algorithme de
#    Recuit Simulé afin de résoudre le problème du voyageur 
#    de commerce.
#
# ###################################################################
#
# License
#
#    This program is a free software you can modifiy it and 
#    redistribute it for non profitable use. 
#
# ###################################################################

#
# Sous cette forme l'algorithme peut parfois atteindre la solution
#
#[2, 3, 8, 0, 6, 5, 11, 7, 4, 9, 12, 1, 10], 18406
#=>  1 :  [2, 3, 7, 0, 6, 5, 11, 8, 4, 9, 12, 1, 10], 13120
#=>  2 :  [2, 3, 7, 0, 9, 5, 11, 8, 4, 6, 12, 1, 10], 10217
#=>  3 :  [2, 3, 7, 0, 9, 5, 11, 8, 1, 6, 12, 4, 10], 8906
#=>  4 :  [2, 3, 7, 0, 9, 5, 11, 1, 8, 6, 12, 4, 10], 8329
#=>  5 :  [2, 3, 7, 0, 10, 5, 11, 1, 8, 6, 12, 4, 9], 7791
#=>  6 :  [2, 9, 7, 0, 10, 5, 11, 1, 8, 6, 12, 4, 3], 7610
#=>  7 :  [2, 0, 7, 9, 10, 5, 11, 1, 8, 6, 12, 4, 3], 7320
#=>  8 :  [2, 7, 0, 9, 10, 5, 11, 1, 8, 6, 12, 4, 3], 7295
#=>  9 :  [2, 7, 0, 9, 5, 10, 11, 1, 8, 6, 12, 4, 3], 7293

#[8, 6, 10, 7, 4, 1, 2, 9, 12, 3, 5, 11, 0], 17395
#=>  1 :  [8, 6, 10, 7, 4, 0, 2, 9, 12, 3, 5, 11, 1], 13207
#=>  2 :  [8, 6, 10, 4, 7, 0, 2, 9, 12, 3, 5, 11, 1], 11063
#=>  3 :  [8, 6, 12, 4, 7, 0, 2, 9, 10, 3, 5, 11, 1], 8951
#=>  4 :  [8, 6, 12, 4, 7, 0, 2, 9, 3, 10, 5, 11, 1], 8101
#=>  5 :  [8, 6, 12, 4, 7, 0, 2, 3, 9, 10, 5, 11, 1], 7817
#=>  6 :  [8, 6, 12, 4, 2, 0, 7, 3, 9, 10, 5, 11, 1], 7736
#=>  7 :  [8, 6, 12, 4, 3, 0, 7, 2, 9, 10, 5, 11, 1], 7345
#=>  8 :  [8, 6, 12, 4, 3, 2, 7, 0, 9, 10, 5, 11, 1], 7295
#=>  9 :  [8, 6, 12, 4, 3, 2, 7, 0, 9, 5, 10, 11, 1], 7293

import random
import sys
import math
import numpy as np
import matplotlib.pyplot as plt


# ===================================================================
# Les données en entrée sont les suivantes :
# - une liste des villes sous forme de liste de chaines de caractères
# - les distances entre les villes sous forme d'une matrice d'entiers
#   (en l'occurrence il s'agit d'une liste de liste d'entiers)
# ===================================================================

villes = [ "1. New York, NY", 
    "2. Los Angeles, CA", 
    "3. Chicago, IL", 
    "4. Minneapolis, MN", 
    "5. Denver, CO", 
    "6. Dallas, TX", 
    "7. Seattle, WA", 
    "8. Boston, MA",  
    "9. San Francisco, CA", 
    "10. St. Louis, MI", 
    "11. Houston, TX", 
    "12. Phoenix, AZ", 
    "13. Salt Lake City, UT" ]
 

distances = [
      [0, 2451, 713, 1018, 1631, 1374, 2408, 213, 2571, 875, 1420, 2145, 1972],
      [2451, 0, 1745, 1524, 831, 1240, 959, 2596, 403, 1589, 1374, 357, 579],
      [713, 1745, 0, 355, 920, 803, 1737, 851, 1858, 262, 940, 1453, 1260],
      [1018, 1524, 355, 0, 700, 862, 1395, 1123, 1584, 466, 1056, 1280, 987],
      [1631, 831, 920, 700, 0, 663, 1021, 1769, 949, 796, 879, 586, 371],
      [1374, 1240, 803, 862, 663, 0, 1681, 1551, 1765, 547, 225, 887, 999],
      [2408, 959, 1737, 1395, 1021, 1681, 0, 2493, 678, 1724, 1891, 1114, 701],
      [213, 2596, 851, 1123, 1769, 1551, 2493, 0, 2699, 1038, 1605, 2300, 2099],
      [2571, 403, 1858, 1584, 949, 1765, 678, 2699, 0, 1744, 1645, 653, 600],
      [875, 1589, 262, 466, 796, 547, 1724, 1038, 1744, 0, 679, 1272, 1162],
      [1420, 1374, 940, 1056, 879, 225, 1891, 1605, 1645, 679, 0, 1017, 1200],
      [2145, 357, 1453, 1280, 586, 887, 1114, 2300, 653, 1272, 1017, 0, 504],
      [1972, 579, 1260, 987, 371, 999, 701, 2099, 600, 1162, 1200, 504, 0]
]


# ===================================================================
#	On définit un chemin (configuration) comme la liste des villes 
#	parcourues et le coût qui lui 	est associé
# ===================================================================
class Path(object):
    """
    a Path is the list of cities visited and the distance travelled
    """
    def __init__(self, chemin, distance):
        self.chemin = chemin
        self.distance = distance

    def __str__(self):
        return str(self.chemin) + ", " + str(self.distance)
    
    def copy(self):
        return type(self)(self.chemin.copy(), self.distance)
        
# ===================================================================
#	On définit une class TSP pour stocker et manipuler les données
#	et résoudre le problème
#	
# ===================================================================
class TSP(object):
	"""
	Class used to solve the problem that implements a Local Search
	algorithm
	"""
	
	# -----------------------------------------------------
	# Le constructeur modifie les données de manière à
	# utiliser des indices qui commencent à 1 et non à 0
	#
	# On modifie la matrice des distances en ajoutant une
	# ligne et une colonne de zéros, puis on convertit tous
	# les 0 en une distance égale à la valeur maximum trouvée
	# dans la matrice + 1, ce qui permettra de rechercher une
	# distance minimum sans obtenir de zéro
	# -----------------------------------------------------
	def __init__( self, villes, distances ):
		"""
		Constructor
		
		Parameters
		----------
		villes : list of string
		  list of cities names
		distances : list of list of int
		  matrix of distances between cities  
		"""
		self.nbr_villes = len(villes)
		self.villes = villes.copy()
		self.villes.insert(0, "-------")
		
		self.distances = np.array( distances )
		
		# add column of 0
		@@@...@@@
		@@@...@@@
		
		# add line of 0
		@@@...@@@
		@@@...@@@
		
		# find maximum distance 
		maximum = @@@...@@@
		
		# replace all 0 distances by maximum + 1
		@@@...@@@

		print("---- Matrice des distances modifiée ----")		
		print(self.distances)
	
		self.verbose = True

	# -----------------------------------------------------
	# Définition de la fonction objectif
	# -----------------------------------------------------
	def f_objectif( self, chemin ):
		distance_parcourue = 0
		
		@@@...@@@
		
		return distance_parcourue  
    
	# -----------------------------------------------------
	# Génération de voisins en échangeant des villes deux
	# à deux. Seules les configurations voisines amélio-
	# rantes sont gardées puis triées en ordre croissant.
	# En fonction du nombre d'itérations on élimine 
	# certaines configurations :
	# - si le nombre d'itérations est inférieu à 1000, on
	#   ne garde que les dix premières configurations 
	#   améliorantes
	# - au dela de 5000 itérations on ne garde que la 
	#   première configuration améliorante
	# -----------------------------------------------------	
	def voisinage( self, configuration, iteration ):

		# liste des voisins améliorants
		l = [ ]
		
		distance_a_ameliorer = int(configuration.distance * 1.1)
		
		for i in range( 1, self.nbr_villes-2):
		    for k in range( i+1, self.nbr_villes):
		        if i > self.nbr_villes or k > self.nbr_villes:
		            continue
		        
		        @@@...@@@
		        
		        if nouvelle_configuration.distance < distance_a_ameliorer:
		            l.append( nouvelle_configuration )
		        
		# tri des configuration en ordre croissant de la
		# distance parcourue        
		@@@...@@@
				
		return l
    
	# -----------------------------------------------------
	# Méthode de résolution basée sur le recuit simulé
	# -----------------------------------------------------
	def resolution_recuit_simule(self, t_i = 100.0, t_f = 0.001, alpha = 0.99):
	
		self.liste_evolution_courante = []
		self.liste_evolution_meilleur = []
		
		#
		# Génération de la configuration initiale
		#
		villes_dans_le_desordre = [ x for x in range(1, self.nbr_villes+1) ]
		random.shuffle( villes_dans_le_desordre )
		configuration_courante = Path( villes_dans_le_desordre, self. f_objectif( villes_dans_le_desordre ) )

		print("Configuration initiale=", configuration_courante)
		
		iteration = 1
		
		#
		# Paramètres de l'algorithme
		#
		
		# température initiale
		T_initiale = t_i
		# température finale
		T_finale   = t_f
		# facteur de refroidissement (alpha)
		facteur_de_refroidissement = alpha
		
		T = T_initiale
		meilleure_configuration = configuration_courante.copy()

		# Tant qu'on a pas atteint la température finale
		while T > T_finale:
			
			str_T = "{:.4f}".format(T)
			str_it = str(iteration)
			str_cc = str( configuration_courante.distance )
			str_mc = str( meilleure_configuration.distance)
			print("T=" + str_T + ",iter=" + str_it + " => " + str_cc + " ! " + str_mc, end='')
			
			# rechercher des voisins améliorants
			l_voisins_ameliorants = self.voisinage( configuration_courante, iteration )
			
			# si il n'y a qu'un seul voisin améliorant on le choisit
			# sinon on choisit au hasard
			if len(l_voisins_ameliorants) == 1:
				i = 0
			else:
				i = random.randint(0, len(l_voisins_ameliorants)-1 )
				
			configuration_voisine = l_voisins_ameliorants[ i ]
			
			str_cv = str(configuration_voisine.distance)
			str_lp = str(l_voisins_ameliorants[ 0 ].distance) 
			str_ld = str(l_voisins_ameliorants[ len(l_voisins_ameliorants)-1 ].distance) 
			print(",voisin=" + str_cv + " (" + str_lp + ", " + str_ld + ")", end='')    
			
			delta_f = configuration_voisine.distance - configuration_courante.distance

			str_u = ""
			str_exp = ""
			
			# si le voisin est meilleur que la configuration courante
			#    on le garde
			# sinon 
			#    on prend un voisin qui détériore la solution
			#    avec une certaine probabilité qui diminue à
			#    mesure que la température diminue
			if delta_f < 0.0:
				configuration_courante = configuration_voisine
			else:
				# générer une valeur entre 0 et 1.0
				u = @@@...@@@
				# calculer exp(-delta_f / T)
				e = @@@...@@@
				    
				str_u = "{:.4f}".format(u)
				str_exp = "{:.4f}".format(e)
				
				if (e > u):
				    configuration_courante = configuration_voisine
				    
			str_delta_f = "{:.4f}".format(delta_f)
			print(",delta_f="+str_delta_f+",u="+str_u+",exp="+str_exp)
				
			if  configuration_courante.distance < meilleure_configuration.distance:
				meilleure_configuration = configuration_courante.copy()
			
			# diminution de la température	
			T = facteur_de_refroidissement * T    
			iteration += 1

			self.liste_evolution_courante.append( configuration_courante.distance )
			self.liste_evolution_meilleur.append( meilleure_configuration.distance )
			
		return meilleure_configuration
				
	# -----------------------------------------------------
	# Affichage du chemin ainsi que des distances 
	# parcourues et de la distance totale
	# -----------------------------------------------------
	def print(self, chemin):
		print()
		print("chemin=",chemin)
		print()
		print("------------------------------------------")
		print("         ville            | dist. | totale")
		print("------------------------------------------")
		distance_totale = 0
		print( "{0:25s} | {1:5d} | {2:5d}".format(self.villes[ chemin[0] ], 0, 0))
		for ville in range(1, len(chemin)):	
			distance = self.distances[ chemin[ville-1] ][ chemin[ville] ]
			distance_totale += distance
			print( "{0:25s} | {1:5d} | {2:5d}".format(self.villes[ chemin[ville] ], distance, distance_totale))
		
		distance = self.distances[ chemin[ville] ][ chemin[0] ]
		distance_totale += distance
		print( "{0:25s} | {1:5d} | {2:5d}".format(self.villes[ chemin[0] ], distance, distance_totale))
		print()
	
	# -----------------------------------------------------
	# Affiche un graphique de l'évolution des scores 
	# (distances) de la configuration courante et de la
	# meilleure configuration
	# -----------------------------------------------------	
	def graphique(self):
		x = range( len( self.liste_evolution_courante ) )
		plt.grid()
		plt.plot( x, self.liste_evolution_courante, label='configuration courante' )
		plt.plot( x, self.liste_evolution_meilleur, label='meilleure configuration' )
		plt.xlabel('itérations')
		plt.legend()
		plt.title('Evolution du Recuit Simulé')
		plt.show()
			
# ===================================================================
# Programme principal
# ===================================================================

HELP_MESSAGE = """
Program arguments are:

-h or --help 
	to get this message

-q or --quiet 
	to remove verbose mode and avoid to display graphics
	
-i float or --t_init=float
	initial temperature (default is 100.0)
	
-f float or --t_final=float
	final temperature (default is 0.001)
	
-a float or --alpha=float		
	cooling factor (default is 0.99)
"""

# creer une instance du problème
tsp = TSP( villes, distances )

# température initiale
t_i = 100.0
# température finale
t_f = 0.001
#facteur_de_refroidissement
alpha = 0.99

#
# gestion des arguments du programme
#
if len(sys.argv) >= 1:
	try:
	
		i = 1
		while i < len(sys.argv):
			arg = sys.argv[i] 
			print(i, arg)
			if arg == "-h" or arg == "--help":
				print( HELP_MESSAGE )
				
			elif arg == "-q" or arg == "--quiet":
				tsp.verbose = False
				
			elif arg == "-i":
				i += 1
				t_i = float(sys.argv[i])
			elif arg.startswith("--t_init="):
				t_i = float(arg.split('=')[1])
				
			elif arg == "-f":
				i += 1
				t_f = float(sys.argv[i])
			elif arg.startswith("--t_final="):
				t_f = float(arg.split('=')[1])
				
			elif arg == "-a":
				i += 1
				alpha = float(sys.argv[i])
			elif arg.startswith("--alpha="):
				alpha = float(arg.split('=')[1])
				
			else:
				raise ValueError("Unexpected argument")
				
			i += 1
		
	except:
		print("="*40)
		print(sys.exc_info()[1])
		print("="*40)
		print( HELP_MESSAGE )
		sys.exit(1)

#
# recherche de la meilleure solution
#

configuration_finale = tsp.resolution_recuit_simule(t_i, t_f, alpha)

tsp.print(configuration_finale.chemin)
tsp.graphique()