#Version modifiée du programme de Peio
#Objet : Résoudre le TSP par le récuit simulé

from scipy import *
from math import *
from matplotlib.pyplot import *
import sys


# ###################################### Parametres du recuit ##################################
T0 = 150  # temperature initiale
Tmin = 1e-3 # temperature finale
tau = 1e4 # constante pour la décroissance de temperature
Alpha = 0.9 # constante pour la décroissance géométrique
Palier = 7 # nombre d'itérations sur un palier de température
IterMax = 15000 # nombre max d'itérations de l'algorithme
# ##############################################################################################

# fonction de fluctuation autour de l'etat "thermique" du systeme : echange de 2 points
# pre-condition :
#   - chemin : ordre de parcours des villes
#   - i,j indices des villes à permuter
# post-condition : nouvel ordre de parcours
def fluctuationDeux(chemin,i,j):
    nv = chemin[:]
    temp = nv[i]
    nv[i] = nv[j]
    nv[j] = temp
    return nv

# Factorisation des fonctions de calcul de la distance totale
# pre-condition : (x1,y1),(x2,y2) coordonnees des villes
# post-condition : distance euclidienne entre 2 villes
def distance(p1,p2):
    return sqrt((p1[0]-p2[0])**2 + (p1[1]-p2[1])**2)

# Fonction de calcul de l'energie du systeme,
# pre-condition : 
#   - coords : coordonnées des points du chemin
#   - chemin : ordre de parcours des villes
# post-condition : Pb du VC : la distance totale du trajet
def energieTotale(coords,chemin):
    energie = 0.0
    coord = coords[chemin]
    for i in range(-1,N-1): # on calcule la distance en fermant la boucle
        energie += distance(coord[i], coord[i+1])
    return energie

# fonction d'implementation de l'algorithme de Metropolis pour un trajet par rapport a son voisin
# pre-conditions :
#   - ch1 voisin, ch2 : trajet init,
#   - disti : distance de chaque trajet
#   - T : température actuelle du système
# post-condition : retourne la fluctuation retenue par le critère de Metropolis
def metropolis(ch1,dist1,ch2,dist2,T):
    global best_trajet, best_dist, x, y,cptx,cpty,cptz
    delta = dist1 - dist2 # calcul du differentiel
    if delta <= 0: # si ameliore,
        if dist1 <= best_dist: # comparaison au best, si meilleur, enregistrement et refresh de la figure
            best_dist = dist1
            best_trajet = ch1[:]
        cptx = cptx+1
        return (ch1, dist1) # la fluctuation est retenue, retourne le voisin
    else:
        if random.uniform() > exp(-delta/T): # la fluctuation n'est pas retenue selon la proba
            cpty = cpty+1
            return (ch2, dist2)              # trajet initial
        else:
            cptz = cptz+1
            return (ch1, dist1)              # la fluctuation est retenue, retourne le voisin

# Parsing du fichier de données
# pre-condition : nomfic : nom de fichier (doit exister)
# post-condition : (x,y) coordonnees des villes
def parse(nomfic):
    absc=[]
    ordo=[]
    with open(nomfic,'r') as inf:
        for line in inf:
            absc.append(float(line.split(' ')[1]))
            ordo.append(float(line.split(' ')[2]))
    return (array(absc),array(ordo))

# Instance du problème
#FIC="NUMERO.tsp"
if (len(sys.argv) > 1):
    FIC=sys.argv[1]
else:
    print("Aucun fichier spécifié...")
    sys.exit("USAGE : python RSTSP.py NUMERO_INSTANCE.tsp")



# initialisation des listes d'historique pour le graphique final
Henergie = []     # energie
Htemps = []       # temps
HT = []           # temperature
Hbest = []        # distance

# ##################################### INITIALISATION DE L'ALGORITHME ############################
# Construction des données depuis le fichier
(x,y) = parse(FIC) # x,y sont gardés en l'état pour l'affichage graphique
coords = array(list(zip(x,y)), dtype=float) # On contruit le tableau de coordonées (x,y)

# Paramètre du probleme
N = len(coords)    # nombre de villes

# defintion du trajet initial : ordre croissant des villes
trajet = [i for i in range(N)]
# calcul de l'energie initiale du systeme (la distance initiale a minimiser)
dist = energieTotale(coords,trajet)
# initialisation du meilleur trajet
best_trajet = trajet[:]
best_dist = dist

# boucle principale de l'algorithme du recuit
t = 0
T = T0
iterPalier = Palier
cptx = 0
cpty = 0
cptz = 0

print(best_trajet,best_dist)

# ##################################### BOUCLE PRINCIPALE DE L'ALGORITHME ############################

# Boucle de convergence sur critère de nb d'itération (pour tester les paramètres)
#for i in range(IterMax):
# Boucle de convergence sur critère de température
while T > Tmin:
    # palier de temperature TODO : ajouter un cpt d'amélioration
    while (iterPalier > 0): 
        # choix de deux villes au hasard
        i = np.random.randint(N)
        j = i
        while (i==j):
            j = np.random.randint(N)

        # creation de la fluctuation et mesure de l'energie
        voisin = fluctuationDeux(trajet,i,j)
        dist_voisin = energieTotale(coords,voisin)

        # application du critere de Metropolis pour determiner la fulctuation conservee
        (trajet, dist) = metropolis(voisin,dist_voisin,trajet,dist,T)

        iterPalier -= 1

    # règles de décroissances de températures
    #T = T0*exp(-t/tau)
    T = T*Alpha
    iterPalier = Palier

print(best_trajet,best_dist)