import math
import matplotlib.pyplot as plt
%matplotlib inline

class TSP:
    
    def __init__(self):
        self.flat_mat = flat_mat
        self.n = 0
        self.melhor_dist = 1e11
        self.pontos = []
        self.melhores_pontos = []

    def busca_exaustiva(self, flat_mat, n, ite):
        
        if ite == n:
            dist = 0

            for j in range(1, n):
                dist += flat_mat[self.pontos[j-1] * n + self.pontos[j]]
            dist += flat_mat[self.pontos[n-1] * n + self.pontos[0]]

            if dist < self.melhor_dist:
                self.melhor_dist = dist
                self.melhores_pontos = self.pontos[:]   
            return
        
        for i in range(n):
            if self.pontos[i] == -1:   
                self.pontos[i] = ite
                self.busca_exaustiva(flat_mat, n, ite + 1)
                self.pontos[i] = -1
                
    def dist_mat(self):
        x = [] 
        y = []
        flat_mat = [] #matriz 1 dimensao contendo todas as distancias possiveis entre os pontos para facilitar cálculo.

        while True:
            try:    
                temp = input("Digite a coordenada x  y: ").split()
                x.append(float(temp[0]))
                y.append(float(temp[1]))
            except:
                break

        for i in range(len(x)):
            for j in range(len(y)):
                flat_mat.append(math.sqrt((x[i] - x[j])**2 + (y[i] - y[j])**2))
        return flat_mat, x, y

    def get_results(self):
        self.flat_mat, x, _ = self.dist_mat()
        self.n = len(x)
        self.pontos = [-1]*self.n
        self.pontos[0] = 0
        self.busca_exaustiva(self.flat_mat, self.n, 1)
        return self.melhor_dist, self.melhores_pontos

Tsp = TSP()
distancia, pontos = Tsp.get_results()
print("Melhor distancia encontrada: ", distancia)
print("Melhor caminho encontrado: ", pontos)

#plota gráfico
def connectpoints(x,y,p1,p2):
    x1, x2 = x[p1], x[p2]
    y1, y2 = y[p1], y[p2]
    plt.plot([x1,x2],[y1,y2],'ro-')

for i in range(1, len(pontos)):
    connectpoints(x,y,pontos[i-1],pontos[i])
connectpoints(x,y,pontos[len(x)-1],pontos[0])

plt.title("Percurso")
plt.show()

%%time
%%cmd
python TSP.py < in-1.txt
type out-1.txt
python TSP.py < in-2.txt
type out-2.txt
python TSP.py < in-3.txt
type out-3.txt
python TSP.py < in-4.txt
type out-4.txt

from qiskit import IBMQ
import numpy as np
#IBMQ.save_account('seu-tokenIBMQ-para-rodar-localmente')
IBMQ.load_account()

from qiskit import Aer
from qiskit.tools.visualization import plot_histogram
from qiskit.circuit.library import TwoLocal
from qiskit.optimization.applications.ising import max_cut, tsp
from qiskit.aqua.algorithms import VQE, NumPyMinimumEigensolver
from qiskit.aqua.components.optimizers import SPSA
from qiskit.aqua import aqua_globals
from qiskit.aqua import QuantumInstance
from qiskit.optimization.applications.ising.common import sample_most_likely
from qiskit.optimization.algorithms import MinimumEigenOptimizer
from qiskit.optimization.problems import QuadraticProgram

import logging
from qiskit.aqua import set_qiskit_aqua_logging

#Preparando os dados segundo os imputs do usuario para serem resolvidos pelo qiskit max 4 pontos por limitação de qubits
coord = []
flat_mat, x, y = TSP().dist_mat()
dist_mat = np.array(flat_mat).reshape(len(x),len(x))
for i, j in zip(x, y):
    coord.append([i,j])

ins = tsp.TspData('TSP_Q', dim=len(x), coord=np.array(coord), w=dist_mat)

qubitOp, offset = tsp.get_operator(ins)
print('Offset:', offset)
print('Ising Hamiltonian:')
print(qubitOp.print_details())

#Usando o numpyMinimumEigensolver como o solver do problema para resolver de forma quantica
ee = NumPyMinimumEigensolver(qubitOp)
result = ee.run()

print('energy:', result.eigenvalue.real)
print('tsp objective:', result.eigenvalue.real + offset)
x_Q = sample_most_likely(result.eigenstate)
print('feasible:', tsp.tsp_feasible(x_Q))
z = tsp.get_tsp_solution(x_Q)
print('solution:', z)
print('solution objective:', tsp.tsp_value(z, ins.w))

for i in range(1, len(z)):
    connectpoints(x,y,z[i-1],z[i])
connectpoints(x,y,z[len(x)-1],z[0])

plt.title("Percurso")
plt.show()

#instanciando o simulador ou o computador real importante lembrar que nao ira funcionar para mais de 4 pontos pelo numero de qubits disponibilizados pela IBM que sao apenas 16 para o simulador qasm e 15 para a maquina quantica
provider = IBMQ.get_provider(hub = 'ibm-q')
device = provider.get_backend('ibmq_16_melbourne')

aqua_globals.random_seed = np.random.default_rng(123)
seed = 10598
backend = Aer.get_backend('qasm_simulator')
#descomentar essa linha caso queira rodar na maquina real
#backend = device
quantum_instance = QuantumInstance(backend, seed_simulator=seed, seed_transpiler=seed)

#rodando no simulador quantico
spsa = SPSA(maxiter=10)
ry = TwoLocal(qubitOp.num_qubits, 'ry', 'cz', reps=5, entanglement='linear')
vqe = VQE(qubitOp, ry, spsa, quantum_instance=quantum_instance)

result = vqe.run(quantum_instance)

print('energy:', result.eigenvalue.real)
print('time:', result.optimizer_time)
x = sample_most_likely(result.eigenstate)
print('feasible:', tsp.tsp_feasible(x_Q))
z = tsp.get_tsp_solution(x_Q)
print('solution:', z)
print('solution objective:', tsp.tsp_value(z, ins.w))

