# A list of variables d_1, ..., d_6 for the coefficients
ds = [var("d%s"%i) for i in range(1,7)]

# A polynomial ring on all these generators
R = PolynomialRing(QQ, ds)


# Rewrite certain anticanonical coordinates in terms of the roots of W(E6)
Y123456=(d1-d2)*(d3-d4)*(d5-d6) 
Y142536=(d1-d4)*(d2-d5)*(d3-d6) 
Y162345=(d1-d6)*(d2-d3)*(d4-d5) 
Y123645=(d1-d2)*(d3-d6)*(d4-d5) 
Y162534=(d1-d6)*(d2-d5)*(d3-d4) 
Y142356=(d1-d4)*(d2-d3)*(d5-d6)
Y152634=(d1-d5)*(d2-d6)*(d3-d4)
Y132456=(d1-d3)*(d2-d4)*(d5-d6)
Y152436=(d1-d5)*(d2-d4)*(d3-d6) 
Y132645=(d1-d3)*(d2-d6)*(d4-d5)

# Two potential cubic equations
cubic = Y123456*Y142536*Y162345 - Y123645*Y162534*Y142356
cubic2 = Y152634*Y123645*Y132456 - Y123456*Y152436*Y132645

expand(R(cubic)) == 0 
# Returns True

expand(R(cubic2)) == 0 
# Returns True

load("../../General/settingsAndPermutations.sage")

cubic = S(Y123456*Y142536*Y162345 - Y123645*Y162534*Y142356)
cubic2 = S(Y152634*Y123645*Y132456 - Y123456*Y152436*Y132645)

ZZ= IntegerRing()

# Pick Random integers
# d = (ZZ.random_element(-1000,1000), ZZ.random_element(-1000,1000), ZZ.random_element(-1000,1000),ZZ.random_element(-1000,1000),ZZ.random_element(-1000,1000),ZZ.random_element(-1000,1000))
d = (-717, -662, 471, -136, -266, -182)

linears = load("../Output/270AnticanonicalLinearEquations.sobj")
linears0 = map(lambda f: f(d1=d[0], d2=d[1], d3=d[2], d4=d[3], d5=d[4], d6=d[5]), linears)

S0 = PolynomialRing(QQ, Xs+Ys)

I0 = S0.ideal(linears0)
S0(cubic) in I0 
#Returns False

S0(cubic2) in I0
#Returns False

load("../../General/orbit.sage")
load("../../General/removeDoubles.sage")

cubic = set([expand(S(cubic))])

allCubics = computeOrbit(cubic, homsStd)
len(allCubics)
#241

allCubicsNoDoubles = removeDoubles(allCubics)
len(allCubicsNoDoubles)
#121

#save(allCubicsNoDoubles,"../Output/121cubicEquations.sobj")
#save(allCubics, "../Output/241cubicEquations.sobj")

#Test that the elements in the W(E6)-orbit of the given cubic lies in the ideal generated by the specialized 270 linear equations and the input cubic
newIdeal = S0.ideal([S0(x) for x in linears0] + [S0(list(cubic)[0])])
all([S0(x) in newIdeal for x in list(allCubics)])
# Reuturns True

cubic2 = set([expand(S(cubic2))])

allCubics2 = computeOrbit(cubic2, homsStd)
len(allCubics2)
#241

allCubicsNoDoubles = removeDoubles(allCubics2)
len(allCubicsNoDoubles)
#121

list(cubic)[0] in allCubics
# Return True
list(cubic)[0] in allCubics2
# Return False

list(cubic2)[0] in allCubics2
# Return True
list(cubic2)[0] in allCubics
# Return False