import random
from IPython.core.display import SVG
import pyomo.environ as pyo
from pysat.solvers import Solver
from pysat.formula import CNF 
#import py_svg_combinatorics as psc
from ipywidgets import widgets, HBox
from collections import Counter
from pprint import pprint
from random import randint
import numpy as np
from IPython.display import IFrame
import IPython
from copy import copy
import os
from pathlib import Path
import itertools
nbname = ''
try:
    nbname = __vsc_ipynb_file__
except:
    if 'COCALC_JUPYTER_FILENAME' in os.environ:
        nbname = os.environ['COCALC_JUPYTER_FILENAME']    
title_ = Path(nbname).stem.replace('-', '_').title()    
IFrame(f'https://discopal.ispras.ru/index.php?title=Hardprob/{title_}&useskin=cleanmonobook', width=1280, height=300)

# не знаю как правильн генерировать граф с соблюдением неравенства треугольника поэтому генерирую рандомные точки
def dist(a, b):
    return int(np.ceil(np.sqrt((a[0] - b[0]) ** 2 + (a[1] - b[1]) ** 2)))

def out_of_verts(v, m, edges):
    for i in range(0, m):
        if i == v:
            continue
        if edges[v][i] == 1000000000:
            return False
    return True

def generate_random_graph(m, n):
    if n < m - 1 or n > m*(m-1) / 2:
        return [[]]
    vertexes = [(random.uniform(-100, 100), random.uniform(-100, 100)) for _ in range(m)]
    edges = [[1000000000] * m for _ in range(m)]
    edges[0][1] = edges[1][0] = dist(vertexes[0], vertexes[1])
    # обеспечение связности графа. поскольку вершины сами по себе рандомны, то можно идти вот так подряд - чтобы всегда был путь
    for i in range(2, m):
        edges[i][i - 1] = edges[i - 1][i] = dist(vertexes[i], vertexes[i - 1])
    #это чтобы всегда был хоть один замкнутый путь
    edges[0][m - 1] = edges[m - 1][0] = dist(vertexes[0], vertexes[m - 1])
    edges_left = n - m
    for _ in range(0, edges_left):
        v1 = randint(0, m - 1)
        while out_of_verts(v1, m, edges):
            v1 = randint(0, m - 1)
        v2 = randint(0, m - 1)
        while v2 == v1 or edges[v1][v2] < 1000000000:
            v2 = randint(0, m - 1)
        edges[v1][v2] = edges[v2][v1] = dist(vertexes[v2], vertexes[v1])
    return edges

def get_model(edges):
    m = pyo.ConcreteModel()
    N = len(edges)
    m.J = pyo.RangeSet(0, N-1)
    # m.J = edge
    m.N = N
    
    m.IJ = m.J * m.J
    m.x = pyo.Var(m.IJ, domain=pyo.Binary)
    # m.x = pyo.Var(m.J, m.J, domain=pyo.Binary)
    
    m.edges = edges
    
    m.obj = pyo.Objective(expr = sum(m.edges[i][j] * m.x[i, j] for i, j in m.IJ), sense=pyo.minimize)
    
    @m.Constraint(m.J)
    def не_должно_быть_селф_циклов(m, i):
        return m.x[i, i] == 0
    
    @m.Constraint(m.J)
    def только_один_вход(m, i):
        return sum(m.x[i, ...]) == 1
    
    @m.Constraint(m.J)
    def только_один_выход(m, i):
        return sum(m.x[..., i]) == 1
    
    m.I1 = pyo.RangeSet(1, N-1)
    m.IJ1 = pyo.Set(dimen=2, initialize=itertools.permutations(list(m.I1), 2))

    m.u = pyo.Var(m.I1, bounds=(0, N - 2))
    
    @m.Constraint(m.IJ1)
    def MillerTuckerZemlin(m, i: m.I1, j: m.I1):
        return m.u[i] - m.u[j] + m.N * m.x[i, j] <= m.N - 1
    
    return m
    
    
    

edges = generate_random_graph(5, 10)

model = get_model(edges)

model.pprint()

solver = pyo.SolverFactory('cbc')
solver.solve(model).write()
res = model.obj()
model.x.extract_values()

def solveByPyomo(edges):
    model = get_model(edges)
    solver = pyo.SolverFactory('cbc')
    solver.solve(model)
    res = model.obj()
    return model.x.extract_values()

def modelResToPerm(dictres, m):
    res = [0]*m
    current = 0
    for i in range(1, m):
        for item, v in dictres.items():
            if item[0] == current and v == 1.0:
                res[i] = item[1]
                current = item[1]
                break
    return res
        

perm = modelResToPerm(model.x.extract_values(), 5)

def sumFromPermutation(edges, perm):
    res = edges[0][len(perm) - 1]
    for i in range(len(perm) - 1):
        res += edges[i][i+1]
    return res

sumFromPermutation(edges, perm)

def primitiveSolving(edges):
    m = len(edges)
    res = 1000000000
    for perm in itertools.permutations(list(range(0, m)), m):
        res = min( sumFromPermutation(edges, perm) , res)
        
    return res

def primitiveTesting(n=100):
    for i in range(n):
        m = randint(5, 10)
        ed = randint(m, (m*(m-1)) // 2)
        edges = generate_random_graph(m, ed)
        if sumFromPermutation(edges, modelResToPerm(solveByPyomo(edges), m)) == primitiveSolving(edges):
            print(f'Test{i}: OK')
        else:
            print(f'Test{i}: Error | m = {m} | edg = {ed}')
            print(edges)
            print(solveByPyomo(edges))
            print(primitiveSolving(edges))