################################################### # Checking genericity for the cubic and the Nodes # ################################################### # We aim to certify that the computations we have done for the anticanonical ideal and the nodes of each surface are valid as long as the surface is smooth and has no Eckardt points. # To this end we rewrite our basis for the P^3 with entries that are Laurent monomials in Yoshida and Cross functions and we check that one of the 4x4 minors is a Laurent monomial in the Yoshida and Cross functions. # In addition, we check that the solution to E1\cap F12 is always a multiple of the last vector of our chosen basis. # Finally, we analyze the computation of the cubic and show that this cubic does not vanish along P^3, so it cuts out the cubic surface in P^3. # We start by loading "OperatorsOnYoshidas.sage": load("OperatorsOnYoshidas.sage") # We load the action of W(E6) on all the variables. The action was computed with the Sage script "WE6_action.sage" allAction = load("Input/Symbolic_action_on_everything.sobj") Yoshidas = [var("Yoshida%s"%i) for i in range(0,40)] YRing = PolynomialRing(QQ, Yoshidas) YField = FractionField(YRing) YCRing = PolynomialRing(QQ, Yoshidas + cross_variables) YCField = FractionField(YCRing) Qd = PolynomialRing(QQ,ds) YM = load("../../Yoshida/Output/Yoshida_matrix.sobj") # The values of the Yoshida functions are stored in 'YOrdered' (loaded in "OperatorsOnYoshidas.sage") YOrdered = load("../../Yoshida/Scripts/Input/orderedYoshidaFunctionsDictionary.sobj") # We load the 4x45 matrix in Yoshidas corresponding to a basis of P^3 in P^44 MMInYoshidas = load("../../Yoshida/Output/YoshidaIntegerMM.sobj") # A = MatrixSpace(YField, 4, 45) # newMMInYoshidas = A.matrix(MMInYoshidas) ################################################################# # Rewriting the Solution matrix in Cross and Yoshida functions # ################################################################# ################################ # Using Roots in Factorization # ################################ # The following function transforms a monomial in Yoshidas into a vector of exponents in positive roots: def fromMonomialToRoots(f): if len(YRing(f).monomials())>1: print "The input is not a monomial" return None fl = SR(f).factor_list() expVector = sum([x[1]*vector(YM.rows()[Yoshidas.index(YRing(x[0]))]) for x in fl]) return expVector # The following will translate a polynomial in Yoshidas into the potentially non-monomial factor of its presentation in ds. def fromNonMonFactorsToDs(f): Monsf = YRing(f).monomials() # print len(Monsf) Coeffsf = YRing(f).coefficients() print len(Coeffsf) ExponentsMonsf = [fromMonomialToRoots(x) for x in Monsf] print ExponentsMonsf # we take the minimal of all exponents of positive roots to remove from each term. This will reduce the final degree of each term. minExps = vector({x: min([y[x] for y in ExponentsMonsf]) for x in range(0,36)}) print minExps # minExpsPoly = SR(VectorToYoshida(list(minExps), positive_roots)) # print minExpsPoly ScaledExponentsMonsf = [x-minExps for x in ExponentsMonsf] newTerms = [SR(VectorToYoshida(list(y), positive_roots)) for y in ScaledExponentsMonsf] return (sum([Coeffsf[i]*newTerms[i] for i in range(0, len(Coeffsf))]), minExps) ##################################### # Monomial and Non-monomial factors # ##################################### # The following function will be used to rewrite the factorization of an expression in Yoshidas into a Monomial in Yoshidas and Cross functions, and non-monomial factors. def fromBinomialsAndMonomialsToCrossesAndYoshidas(factorsList): nonMonomialFactors = [] newFactorsList = [] for x in factorsList: if len(YRing(x[0]).monomials())==2: newx = ([SR(y) for y in cross_to_yoshidas.keys() if SR(x[0]) in cross_to_yoshidas[y]][0], x[1]) newFactorsList = newFactorsList + [newx] elif len(YRing(x[0]).monomials())==1: newFactorsList = newFactorsList + [(SR(x[0]), x[1])] else: nonMonomialFactors = nonMonomialFactors + [(SR(x[0]), x[1])] return [newFactorsList, nonMonomialFactors] ################################################# # Rewriting using 16 Standard Yoshida functions # ################################################# # We analyze non-monomial factors separately. In order to find extra monomial factors, we rewrite each term of the non-monomial factor in terms of the 16 Yoshidas generating the row space of the Yoshida matrix. # Qs= polynomial ring on the parameters d1,..., d6 # YM = Yoshida matrix # f = irreducible polynomial in the Yoshida functions. def computingFixedExpressionsInYs(f,Qs,YM): newFactors = dict() Monsf = YRing(f).monomials() Coeffsf = YRing(f).coefficients() ExponentsMonsf = {x:fromMonomialToRoots(x) for x in Monsf} for x in Monsf: linear_vector = vector(ZZ,ExponentsMonsf[x]) if linear_vector in span(ZZ, YM.rows()): yoshida_vector = YM.solve_left(linear_vector) newFactors[Monsf.index(x)]= prod([Yoshidas[i]^yoshida_vector[i] for i in range(0, 40)]) else: print "we are in trouble" return None return sum([Coeffsf[k]*newFactors[k] for k in range(0, len(Monsf))]) ########################################################################### # Functions for writing ratios of Cross functions as monomial in Yoshidas # ########################################################################### # The next function expresses a product of positive roots as a Laurent monomial in the 16 exponent vectors of a Yoshida Q-basis, whenever possible. The result is correct up to sign, and we check the sign on a different function def constructMonomialFromProductFactored(YoshidaFactors, y_to_d, roots, basisOfYoshidas): v_exp = vector(36*[0]) for (f,e) in YoshidaFactors: print (f,e) v_exp = v_exp + e*exponent_vector(DField(SR(f).substitute(y_to_d)), roots) print v_exp soln = matrix(basisOfYoshidas).transpose() \ vector(v_exp) return soln # The next function finds the correct sign in the monomial expression of Yoshidas. We multiply the input vertex by 3 to avoid cube roots of Yoshidas. labels = [5]+[17..31] def findConstant(ratioInY): factorInYs = 1 tripleV = vector(ratioInY) for k in range(0,16): factorInYs = factorInYs*(Yoshidas[labels[k]].factor())^tripleV[k] return factorInYs os.chdir('../../ComputeLines/Scripts') DField = FractionField(PolynomialRing(QQ, ds)) YtoD = YField.hom([DField(y_to_d[y]) for y in y_to_d.keys()]) # The following is a basis of 16 vectors giving independent columns in the Yoshida matrix: uyvs_basis = load("Input/expVectorsZBasisIndex3.sobj") os.chdir('../../Yoshida/Scripts') ########################################################################### # Rewriting the matrix as Laurent monomial in Yoshida and Cross functions # ########################################################################### # We start by writing the factorization of the matrix in Yoshidas: MMInYCFactors = dict() for key in MMInYoshidas.keys(): print key test = MMInYoshidas[key] if test == 0: MMInYCFactors[key] = 0 else: MMInYCFactors[key] = fromBinomialsAndMonomialsToCrossesAndYoshidas(test.factor_list()) # Next, we confirm that no extraneous factor arise: [key for key in MMInYCFactors.keys() if MMInYCFactors[key] !=0 if MMInYCFactors[key][1] !=[]] [] # We use the above information to rewrite the MMInYoshidas as Laurent monomials in Yoshida and Cross functions: # MMInYC = dict() # for key in MMInYCFactors.keys(): # test = MMInYCFactors[key] # if test == 0: # MMInYC[key] = test # else: # MMInYC[key] = prod([z[0]^z[1] for z in MMInYCFactors[key][0]]) # We record the data: # save(MMInYC, "Input/MMInYC.sobj") # # We record the output as a plain text file: # f = file("../Output/SolutionMatrix_4x45_in_YC.txt",'w') # f.write('{') # for i in range(0,4): # for j in range(0,45): # if (i,j) == (3,44): # f.write(str((i,j)) + ':' + str(MMInYC[(i,j)]) + '}') # else: # f.write(str((i,j)) + ':' + str(MMInYC[(i,j)]) + ',') # f.close() # We load the matrix MMInYC = load("Input/MMInYC.sobj") ############################## # Checking the basis for P^3 # ############################## # In order to ensure that the linear equations are valid without any genericity assumptions, MatrixMMInYC = matrix(YCField, MMInYC) [S.gens()[x] for x in [0,1,5,6]] [X12, X13, X21, X23] minor_0_1_5_6 = MatrixMMInYC.matrix_from_columns([0,1,5,6]) nonSingDet = det(minor_0_1_5_6) nonSingDet Yoshida0^3*Yoshida1*Yoshida2^6*Yoshida4^14*Yoshida6^4*Yoshida8^8*Yoshida10^7*Yoshida15^7*Yoshida17^11/(Yoshida3^3*Yoshida5^7*Yoshida7^8*Yoshida9^12*Yoshida11^10*Yoshida12^14*Yoshida13^4*Cross11*Cross78*Cross106*Cross115) # This minor is a Laurent monomial in Yoshida and Cross functions. We rewrite it with standard Yoshida functions. # SR(nonSingDet).numerator().factor_list() # [(Yoshida0, 3), # (Yoshida1, 1), # (Yoshida10, 7), # (Yoshida15, 7), # (Yoshida17, 11), # (Yoshida2, 6), # (Yoshida4, 14), # (Yoshida6, 4), # (Yoshida8, 8)] # SR(nonSingDet).denominator().factor_list() # [(Cross106, 1), # (Cross11, 1), # (Cross115, 1), # (Cross78, 1), # (Yoshida11, 10), # (Yoshida12, 14), # (Yoshida13, 4), # (Yoshida3, 3), # (Yoshida5, 7), # (Yoshida7, 8), # (Yoshida9, 12)] allfactors = [x for x in SR(nonSingDet).numerator().factor_list()] + [(x[0],-x[1]) for x in SR(nonSingDet).denominator().factor_list()] print allfactors [(Yoshida0, 3), (Yoshida1, 1), (Yoshida10, 7), (Yoshida15, 7), (Yoshida17, 11), (Yoshida2, 6), (Yoshida4, 14), (Yoshida6, 4), (Yoshida8, 8), (Cross106, -1), (Cross11, -1), (Cross115, -1), (Cross78, -1), (Yoshida11, -10), (Yoshida12, -14), (Yoshida13, -4), (Yoshida3, -3), (Yoshida5, -7), (Yoshida7, -8), (Yoshida9, -12)] crossFactors = [x for x in allfactors if x[0] in cross_variables] yoshFactors = [x for x in allfactors if x[0] in Yoshidas] crossFactors [(Cross106, -1), (Cross11, -1), (Cross115, -1), (Cross78, -1)] crossPiece = prod([x[0]^x[1] for x in crossFactors]) crossPiece 1/(Cross106*Cross11*Cross115*Cross78) print yoshFactors [(Yoshida0, 3), (Yoshida1, 1), (Yoshida10, 7), (Yoshida15, 7), (Yoshida17, 11), (Yoshida2, 6), (Yoshida4, 14), (Yoshida6, 4), (Yoshida8, 8), (Yoshida11, -10), (Yoshida12, -14), (Yoshida13, -4), (Yoshida3, -3), (Yoshida5, -7), (Yoshida7, -8), (Yoshida9, -12)] YoshPieceNum = prod([SR(x[0])^x[1] for x in yoshFactors if x[1]>0]) YoshPieceDenom = prod([x[0]^x[1] for x in yoshFactors if x[1]<0]) testNum =SR(YoshPieceNum).substitute(y_to_d) testDenom =SR(YoshPieceDenom).substitute(y_to_d) newvalue = constructMonomialFromProductFactored(yoshFactors, y_to_d, roots, uyvs_basis) newvalue (0, -2, 0, -1, 0, 1, 4, -2, 0, 0, 0, -2, 1, 1, 2, 1) newExpression = findConstant(newvalue) newExpression Yoshida21*Yoshida22^4*Yoshida28*Yoshida29*Yoshida30^2*Yoshida31/(Yoshida17^2*Yoshida19*Yoshida23^2*Yoshida27^2) # Next, we need to test the sign: # We pick a generic value of ds: vald = {SR(d1): -717, SR(d2):-662, SR(d3):471,SR(d4): -136, SR(d5):-266, SR(d6):-182} numNewExpression = SR(newExpression.substitute(y_to_d)).substitute(vald) YoshPiece = prod([SR(x[0])^x[1] for x in yoshFactors]) numYoshPiece = SR(YoshPiece.substitute(y_to_d)).substitute(vald) numNewExpression/numYoshPiece 1 # There is no sign issue. # CONCLUSION: The matrix indeed yields a basis for the P^3. crossPiece*newExpression Yoshida21*Yoshida22^4*Yoshida28*Yoshida29*Yoshida30^2*Yoshida31/(Cross106*Cross11*Cross115*Cross78*Yoshida17^2*Yoshida19*Yoshida23^2*Yoshida27^2) ###################### # Checking the nodes # ###################### # We must certify that the computation of the formula for the nodes is valid for all surfaces, not the tropically stable ones. For this, we must check that the computation of a single pd yields a unique formula. By the action of W(E_6), it suffices to work with a single node. We choose E_1 \cap F_12. # The computation is done with the ds. We must check that the 4x9 system emanating from setting the 9 anticanonical coordinates vanishing for E1 \cap F12 to have a unique solution! This is done by finding a minor of this system whose value is written only in terms of quintics and roots. # We load the 4x45 matrix whose rows span the P^3 spanned by the del Pezzo cubic X inside P^44. ModIntegerMM = load("Input/ModIntegerMM.sobj") # We load the correspondence matrix: correspondence = load("Input/XYEFG_correspondence.sobj") # We load the 135 Nodes: all135Nodes = load("../../Yoshida/Scripts/Input/all135ClassicalNodes.sobj") # We use these functions as inputs for "CalculationsOfAllPds.sage" # Pij = Intersection of Ei and Fij # Gives a list of 0/Infinity entries where Infinity corresponds to anticanonical coordinates divisible by line1 or line2 def p(line1, line2): answer = dict() for a in correspondence.keys(): if T(line1).divides(correspondence[a]) or T(line2).divides(correspondence[a]): answer[a] = Infinity return answer.keys() def findMatrixForNodes(MatrixYC, node): colsCoords = p(node[0], node[1]) cols = sorted([S.gens().index(y) for y in colsCoords]) return MatrixYC.matrix_from_rows_and_columns([0,1,2,3],cols) # We compute the 4x9 matrices to determine the coordinates of each node all4x9NodeMatricesYC = dict() all4x9NodeMatricesds = dict() for node in all135Nodes.keys(): print node all4x9NodeMatricesYC[node] = findMatrixForNodes(MatrixMMInYC, node) all4x9NodeMatricesds[node] = findMatrixForNodes(ModIntegerMM, node) # We record the non-zero entries of each 4x9 matrix. nonZeroEntriesYC = dict() nonZeroEntriesds = dict() for node in all4x9NodeMatricesYC.keys(): print node nonZeroEntriesYC[node] = sorted([(i,j) for i in range(0,4) for j in range(0,9) if all4x9NodeMatricesYC[node][i][j] != 0]) nonZeroEntriesds[node] = sorted([(i,j) for i in range(0,4) for j in range(0,9) if all4x9NodeMatricesds[node][i][j] != 0]) # We check that the number of non-zero entries varies: set([len(nonZeroEntriesYC[node]) for node in nonZeroEntriesYC.keys()]) {15, 22, 24, 28, 30, 31, 32, 33} set([len(nonZeroEntriesds[node]) for node in nonZeroEntriesds.keys()]) {15, 22, 24, 28, 30, 31, 32, 33} # The following function print the shape of the submatrix def shape4x9(nonZeroEntriesPerNode): print " ", for j in range(0,9): print j, print "" for i in range(0,4): print i, for j in range(0,9): if (i,j) in nonZeroEntriesPerNode: print "*", else: print " ", print "" return None ############### # E1 \cap F12 # ############### node = (E1,F12) node in all135Nodes.keys() True print nonZeroEntriesYC[node] [(0, 0), (0, 2), (0, 3), (0, 4), (0, 6), (0, 7), (0, 8), (1, 1), (1, 2), (1, 3), (1, 4), (2, 5), (2, 6), (2, 7), (2, 8)] shape4x9(nonZeroEntriesYC[node]) 0 1 2 3 4 5 6 7 8 0 * * * * * * * 1 * * * * 2 * * * * 3 minor = all4x9NodeMatricesYC[node].matrix_from_rows_and_columns([0,1,2],[0,1,5]) det = minor.det() # det # Yoshida0^2*Yoshida1*Yoshida2^4*Yoshida4^10*Yoshida6^3*Yoshida8^6*Yoshida10^5*Yoshida15^5*Yoshida17^8/(Yoshida3^2*Yoshida5^5*Yoshida7^6*Yoshida9^9*Yoshida11^7*Yoshida12^10*Yoshida13^3*Cross11*Cross78*Cross115) minorFactorsDenom = SR(YCField(det).denominator()).factor_list() print minorFactorsDenom [(Cross11, 1), (Cross115, 1), (Cross78, 1), (Yoshida11, 7), (Yoshida12, 10), (Yoshida13, 3), (Yoshida3, 2), (Yoshida5, 5), (Yoshida7, 6), (Yoshida9, 9)] minorFactorsNum = SR(YCField(det).numerator()).factor_list() print minorFactorsNum [(Yoshida0, 2), (Yoshida1, 1), (Yoshida10, 5), (Yoshida15, 5), (Yoshida17, 8), (Yoshida2, 4), (Yoshida4, 10), (Yoshida6, 3), (Yoshida8, 6)] # This ensures that the unique solution to the system is e_3 as long as all Yoshidas and Cross functions do not vanish. ################################################# # Certifying the validity of the cubic equation # ################################################# # We start by writing four variables that will be used to write linear combinations of our basis for P^3. # Given the coordinates of P^3 and the binomial cubic equation of the universal del Pezzo cubic, we re-express the cubic equation from the 4 parameters defining the cubic. avars = [var("a0"), var("a1"), var("a2"), var("a3")] YRas = PolynomialRing(YRing, avars) YFas = PolynomialRing(YField,avars) YCRas = PolynomialRing(YCRing, avars) YCFas = PolynomialRing(YCField,avars) A = MatrixSpace(YCField, 4, 45) newMMInYC = A.matrix(MMInYC) # We create a generic vector for the P^3: AYCvector = vector(avars)*newMMInYC cubic = Cross107*Cross2*Cross9*X12*X23*X31 - Cross116*Cross80*Cross86*X13*X21*X32 ListXY = Xs + Ys # We create a dictionary to substitute the values of the anticanonical coordinates on the cubic. dictAYCVector = {x: AYCvector[S.gens().index(x)] for x in S.gens()} cubicInA = cubic.substitute(dictAYCVector) numCubic = YCFas(cubicInA).numerator() denomCubic = YCFas(cubicInA).denominator() denomCubic Yoshida1*Yoshida3^2*Yoshida5^6*Yoshida7^7*Yoshida9^10*Yoshida11^9*Yoshida12^12*Yoshida13^5*Cross0*Cross1*Cross11*Cross78*Cross84*Cross85*Cross106*Cross115 # Since the denom is a monomial in Yoshida and Cross functions, from now on we work with 'numCubic' numCubic Yoshida0*Yoshida1^2*Yoshida2^3*Yoshida3*Yoshida4^8*Yoshida5^2*Yoshida6^2*Yoshida7^2*Yoshida8^5*Yoshida9^2*Yoshida10^4*Yoshida11^4*Yoshida12^4*Yoshida13^2*Yoshida15^4*Yoshida17^7*Cross0*Cross1*Cross11*Cross80*Cross84*Cross85*Cross106^2*Cross116*a0*a1*a2 + Yoshida0^2*Yoshida1^2*Yoshida2^4*Yoshida3*Yoshida4^8*Yoshida5*Yoshida6^2*Yoshida7^2*Yoshida8^5*Yoshida9^2*Yoshida10^4*Yoshida11^4*Yoshida12^4*Yoshida13^3*Yoshida15^4*Yoshida17^7*Cross0*Cross1*Cross11*Cross80*Cross85*Cross86*Cross106*Cross116*a1^2*a2 + (-Yoshida1^2*Yoshida2^3*Yoshida3^2*Yoshida4^7*Yoshida5*Yoshida6^2*Yoshida7^3*Yoshida8^3*Yoshida9^3*Yoshida10^4*Yoshida11^7*Yoshida12^7*Yoshida13^4*Yoshida15^2*Yoshida17^5*Cross0*Cross1*Cross11*Cross80*Cross84*Cross86*Cross106*Cross116)*a1*a2^2 + (-Yoshida0^2*Yoshida1*Yoshida2^6*Yoshida4^13*Yoshida6^2*Yoshida8^7*Yoshida10^5*Yoshida15^7*Yoshida17^10*Cross0*Cross2*Cross9*Cross78*Cross84*Cross85*Cross107*Cross115)*a0^2*a3 + (-Yoshida0^2*Yoshida1*Yoshida2^5*Yoshida4^12*Yoshida5*Yoshida6^2*Yoshida8^7*Yoshida10^5*Yoshida13*Yoshida15^7*Yoshida17^10*Cross0*Cross1*Cross9*Cross78*Cross84*Cross85*Cross107*Cross115^2)*a0*a1*a3 + Yoshida0*Yoshida2^5*Yoshida3*Yoshida4^11*Yoshida6^2*Yoshida7*Yoshida8^6*Yoshida9*Yoshida10^4*Yoshida11^3*Yoshida12^3*Yoshida13^2*Yoshida15^5*Yoshida17^8*Cross0*Cross1*Cross9*Cross78^2*Cross84*Cross85*Cross107*Cross115*a0*a2*a3 + (-Yoshida0*Yoshida1^2*Yoshida2^3*Yoshida3*Yoshida4^7*Yoshida5^2*Yoshida6^2*Yoshida7^3*Yoshida8^4*Yoshida9^3*Yoshida10^4*Yoshida11^5*Yoshida12^5*Yoshida13^3*Yoshida15^3*Yoshida17^6*Cross0*Cross1*Cross11^2*Cross80*Cross84*Cross85*Cross106*Cross116)*a1*a2*a3 + Yoshida0^3*Yoshida1*Yoshida2^6*Yoshida4^12*Yoshida6^2*Yoshida8^6*Yoshida10^6*Yoshida13*Yoshida15^7*Yoshida17^10*Cross1*Cross2*Cross9*Cross78*Cross84*Cross85*Cross107*Cross115*a0*a3^2 numCubic.monomials() [a0*a1*a2, a1^2*a2, a1*a2^2, a0^2*a3, a0*a1*a3, a0*a2*a3, a1*a2*a3, a0*a3^2] numCubic.coefficients() [Yoshida0*Yoshida1^2*Yoshida2^3*Yoshida3*Yoshida4^8*Yoshida5^2*Yoshida6^2*Yoshida7^2*Yoshida8^5*Yoshida9^2*Yoshida10^4*Yoshida11^4*Yoshida12^4*Yoshida13^2*Yoshida15^4*Yoshida17^7*Cross0*Cross1*Cross11*Cross80*Cross84*Cross85*Cross106^2*Cross116, Yoshida0^2*Yoshida1^2*Yoshida2^4*Yoshida3*Yoshida4^8*Yoshida5*Yoshida6^2*Yoshida7^2*Yoshida8^5*Yoshida9^2*Yoshida10^4*Yoshida11^4*Yoshida12^4*Yoshida13^3*Yoshida15^4*Yoshida17^7*Cross0*Cross1*Cross11*Cross80*Cross85*Cross86*Cross106*Cross116, -Yoshida1^2*Yoshida2^3*Yoshida3^2*Yoshida4^7*Yoshida5*Yoshida6^2*Yoshida7^3*Yoshida8^3*Yoshida9^3*Yoshida10^4*Yoshida11^7*Yoshida12^7*Yoshida13^4*Yoshida15^2*Yoshida17^5*Cross0*Cross1*Cross11*Cross80*Cross84*Cross86*Cross106*Cross116, -Yoshida0^2*Yoshida1*Yoshida2^6*Yoshida4^13*Yoshida6^2*Yoshida8^7*Yoshida10^5*Yoshida15^7*Yoshida17^10*Cross0*Cross2*Cross9*Cross78*Cross84*Cross85*Cross107*Cross115, -Yoshida0^2*Yoshida1*Yoshida2^5*Yoshida4^12*Yoshida5*Yoshida6^2*Yoshida8^7*Yoshida10^5*Yoshida13*Yoshida15^7*Yoshida17^10*Cross0*Cross1*Cross9*Cross78*Cross84*Cross85*Cross107*Cross115^2, Yoshida0*Yoshida2^5*Yoshida3*Yoshida4^11*Yoshida6^2*Yoshida7*Yoshida8^6*Yoshida9*Yoshida10^4*Yoshida11^3*Yoshida12^3*Yoshida13^2*Yoshida15^5*Yoshida17^8*Cross0*Cross1*Cross9*Cross78^2*Cross84*Cross85*Cross107*Cross115, -Yoshida0*Yoshida1^2*Yoshida2^3*Yoshida3*Yoshida4^7*Yoshida5^2*Yoshida6^2*Yoshida7^3*Yoshida8^4*Yoshida9^3*Yoshida10^4*Yoshida11^5*Yoshida12^5*Yoshida13^3*Yoshida15^3*Yoshida17^6*Cross0*Cross1*Cross11^2*Cross80*Cross84*Cross85*Cross106*Cross116, Yoshida0^3*Yoshida1*Yoshida2^6*Yoshida4^12*Yoshida6^2*Yoshida8^6*Yoshida10^6*Yoshida13*Yoshida15^7*Yoshida17^10*Cross1*Cross2*Cross9*Cross78*Cross84*Cross85*Cross107*Cross115] # allFactorsCoeffs = {numCubic.monomials()[k]: YCRing(numCubic.coefficients()[k]).factor() for k in range(0,len(numCubic.monomials()))} # allFactorsCoeffs # {a0*a3^2: Cross115 * Cross107 * Cross85 * Cross84 * Cross78 * Cross9 * Cross2 * Cross1 * Yoshida13 * Yoshida1 * Yoshida6^2 * Yoshida0^3 * Yoshida10^6 * Yoshida8^6 * Yoshida2^6 * Yoshida15^7 * Yoshida17^10 * Yoshida4^12, # a1*a2*a3: (-1) * Cross116 * Cross106 * Cross85 * Cross84 * Cross80 * Cross1 * Cross0 * Yoshida3 * Yoshida0 * Cross11^2 * Yoshida6^2 * Yoshida5^2 * Yoshida1^2 * Yoshida15^3 * Yoshida13^3 * Yoshida9^3 * Yoshida7^3 * Yoshida2^3 * Yoshida10^4 * Yoshida8^4 * Yoshida12^5 * Yoshida11^5 * Yoshida17^6 * Yoshida4^7, # a0*a2*a3: Cross115 * Cross107 * Cross85 * Cross84 * Cross9 * Cross1 * Cross0 * Yoshida9 * Yoshida7 * Yoshida3 * Yoshida0 * Cross78^2 * Yoshida13^2 * Yoshida6^2 * Yoshida12^3 * Yoshida11^3 * Yoshida10^4 * Yoshida15^5 * Yoshida2^5 * Yoshida8^6 * Yoshida17^8 * Yoshida4^11, # a0*a1*a3: (-1) * Cross107 * Cross85 * Cross84 * Cross78 * Cross9 * Cross1 * Cross0 * Yoshida13 * Yoshida5 * Yoshida1 * Cross115^2 * Yoshida6^2 * Yoshida0^2 * Yoshida10^5 * Yoshida2^5 * Yoshida15^7 * Yoshida8^7 * Yoshida17^10 * Yoshida4^12, # a0^2*a3: (-1) * Cross115 * Cross107 * Cross85 * Cross84 * Cross78 * Cross9 * Cross2 * Cross0 * Yoshida1 * Yoshida6^2 * Yoshida0^2 * Yoshida10^5 * Yoshida2^6 * Yoshida15^7 * Yoshida8^7 * Yoshida17^10 * Yoshida4^13, # a1*a2^2: (-1) * Cross116 * Cross106 * Cross86 * Cross84 * Cross80 * Cross11 * Cross1 * Cross0 * Yoshida5 * Yoshida15^2 * Yoshida6^2 * Yoshida3^2 * Yoshida1^2 * Yoshida9^3 * Yoshida8^3 * Yoshida7^3 * Yoshida2^3 * Yoshida13^4 * Yoshida10^4 * Yoshida17^5 * Yoshida12^7 * Yoshida11^7 * Yoshida4^7, # a1^2*a2: Cross116 * Cross106 * Cross86 * Cross85 * Cross80 * Cross11 * Cross1 * Cross0 * Yoshida5 * Yoshida3 * Yoshida9^2 * Yoshida7^2 * Yoshida6^2 * Yoshida1^2 * Yoshida0^2 * Yoshida13^3 * Yoshida15^4 * Yoshida12^4 * Yoshida11^4 * Yoshida10^4 * Yoshida2^4 * Yoshida8^5 * Yoshida17^7 * Yoshida4^8, # a0*a1*a2: Cross116 * Cross85 * Cross84 * Cross80 * Cross11 * Cross1 * Cross0 * Yoshida3 * Yoshida0 * Cross106^2 * Yoshida13^2 * Yoshida9^2 * Yoshida7^2 * Yoshida6^2 * Yoshida5^2 * Yoshida1^2 * Yoshida2^3 * Yoshida15^4 * Yoshida12^4 * Yoshida11^4 * Yoshida10^4 * Yoshida8^5 * Yoshida17^7 * Yoshida4^8} # Conclusion: All coefficients are non-zero! ##################################### # Check cubic is always irreducible # ##################################### # In order to check that the cubic determines X, we must ensure that it is irreducible over P^3. # We use the shape of the cubic in the a variables to write down a potential factorization. The only option is: (var("W"), var("B"), var("C"), var("D"), var("E")) factoredCubic = ( a0 + B*a1 + C*a2 + D*a3)* (a0*a3 + E*a1*a2) factoredCubic.expand() B*E*a1^2*a2 + C*E*a1*a2^2 + D*E*a1*a2*a3 + E*a0*a1*a2 + B*a0*a1*a3 + C*a0*a2*a3 + D*a0*a3^2 + a0^2*a3 # This agrees with the original cubic up to a scalar, in particular, if we take four suitable ratios we should recover E. numCubic.monomials() [a0*a1*a2, a1^2*a2, a1*a2^2, a0^2*a3, a0*a1*a3, a0*a2*a3, a1*a2*a3, a0*a3^2] # This imposes four restrictions on the coefficients: coeffs = numCubic.coefficients() ########### # Ratio 1 # ########### # E = coeff(a1^2*a2)/coeff(a0*a1*a3) k = numCubic.monomials().index(a1^2*a2) l = numCubic.monomials().index(a0*a1*a3) e0 = coeffs[k]/coeffs[l] e0 Yoshida1*Yoshida3*Yoshida7^2*Yoshida9^2*Yoshida11^4*Yoshida12^4*Yoshida13^2*Cross11*Cross80*Cross86*Cross106*Cross116/(-Yoshida2*Yoshida4^4*Yoshida8^2*Yoshida10*Yoshida15^3*Yoshida17^3*Cross9*Cross78*Cross84*Cross107*Cross115^2) ########### # Ratio 2 # ########### # E = coeff(a1*a2^2)/coeff(a0*a2*a3) k = numCubic.monomials().index(a1*a2^2) l = numCubic.monomials().index(a0*a2*a3) e1 = coeffs[k]/coeffs[l] e1 (-Yoshida1^2*Yoshida3*Yoshida5*Yoshida7^2*Yoshida9^2*Yoshida11^4*Yoshida12^4*Yoshida13^2*Cross11*Cross80*Cross86*Cross106*Cross116)/(Yoshida0*Yoshida2^2*Yoshida4^4*Yoshida8^3*Yoshida15^3*Yoshida17^3*Cross9*Cross78^2*Cross85*Cross107*Cross115) ########### # Ratio 3 # ########### # E = coeff(a1*a2*a3)/coeff(a0*a3^2) k = numCubic.monomials().index(a1*a2*a3) l = numCubic.monomials().index(a0*a3^2) e2 = coeffs[k]/coeffs[l] e2 (-Yoshida1*Yoshida3*Yoshida5^2*Yoshida7^3*Yoshida9^3*Yoshida11^5*Yoshida12^5*Yoshida13^2*Cross0*Cross11^2*Cross80*Cross106*Cross116)/(Yoshida0^2*Yoshida2^3*Yoshida4^5*Yoshida8^2*Yoshida10^2*Yoshida15^4*Yoshida17^4*Cross2*Cross9*Cross78*Cross107*Cross115) ########### # Ratio 4 # ########### # E = coeff(a0*a1*a2)/coeff(a0^2*a3) k = numCubic.monomials().index(a0*a1*a2) l = numCubic.monomials().index(a0^2*a3) e3 = coeffs[k]/coeffs[l] e3 Yoshida1*Yoshida3*Yoshida5^2*Yoshida7^2*Yoshida9^2*Yoshida11^4*Yoshida12^4*Yoshida13^2*Cross1*Cross11*Cross80*Cross106^2*Cross116/(-Yoshida0*Yoshida2^3*Yoshida4^5*Yoshida8^2*Yoshida10*Yoshida15^3*Yoshida17^3*Cross2*Cross9*Cross78*Cross107*Cross115) # Equations: e0 = e1 = e2 = e3 # Check the denominators are Laurent monomials in Yoshidas and Cross functions, so we can ignore them: # (e0-e1).denominator() # -Yoshida0*Yoshida2^2*Yoshida4^4*Yoshida8^3*Yoshida10*Yoshida15^3*Yoshida17^3*Cross9*Cross78^2*Cross84*Cross85*Cross107*Cross115^2 # (e0-e2).denominator() # -Yoshida0^2*Yoshida2^3*Yoshida4^5*Yoshida8^2*Yoshida10^2*Yoshida15^4*Yoshida17^4*Cross2*Cross9*Cross78*Cross84*Cross107*Cross115^2 # (e0-e3).denominator() # Yoshida0*Yoshida2^3*Yoshida4^5*Yoshida8^2*Yoshida10*Yoshida15^3*Yoshida17^3*Cross2*Cross9*Cross78*Cross84*Cross107*Cross115^2 # We work with the numerators: e01 = (e0 - e1).numerator() e02 = (e0 - e2).numerator() e03 = (e0 -e3).numerator() e01Factors = SR(e01).factor_list() e02Factors = SR(e02).factor_list() e03Factors = SR(e03).factor_list() # We record the factors, since we only care about non-monomial ones. Only one appears on each list: # e01Factors # [(Cross115*Cross84*Yoshida1*Yoshida10*Yoshida5 - Cross78*Cross85*Yoshida0*Yoshida2*Yoshida8, # 1), # (Cross106, 1), # (Cross11, 1), # (Cross116, 1), # (Cross80, 1), # (Cross86, 1), # (Yoshida1, 1), # (Yoshida11, 4), # (Yoshida12, 4), # (Yoshida13, 2), # (Yoshida3, 1), # (Yoshida7, 2), # (Yoshida9, 2), # (-1, 1)] f01 = e01Factors[0][0] f01 Cross115*Cross84*Yoshida1*Yoshida10*Yoshida5 - Cross78*Cross85*Yoshida0*Yoshida2*Yoshida8 # e02Factors # [(Cross2*Cross86*Yoshida0^2*Yoshida10*Yoshida15*Yoshida17*Yoshida2^2*Yoshida4 - Cross0*Cross11*Cross115*Cross84*Yoshida11*Yoshida12*Yoshida5^2*Yoshida7*Yoshida9, # 1), # (Cross106, 1), # (Cross11, 1), # (Cross116, 1), # (Cross80, 1), # (Yoshida1, 1), # (Yoshida11, 4), # (Yoshida12, 4), # (Yoshida13, 2), # (Yoshida3, 1), # (Yoshida7, 2), # (Yoshida9, 2)] f02 = e02Factors[0][0] f02 Cross2*Cross86*Yoshida0^2*Yoshida10*Yoshida15*Yoshida17*Yoshida2^2*Yoshida4 - Cross0*Cross11*Cross115*Cross84*Yoshida11*Yoshida12*Yoshida5^2*Yoshida7*Yoshida9 # e03Factors # [(Cross2*Cross86*Yoshida0*Yoshida2^2*Yoshida4 - Cross1*Cross106*Cross115*Cross84*Yoshida5^2, # 1), # (Cross106, 1), # (Cross11, 1), # (Cross116, 1), # (Cross80, 1), # (Yoshida1, 1), # (Yoshida11, 4), # (Yoshida12, 4), # (Yoshida13, 2), # (Yoshida3, 1), # (Yoshida7, 2), # (Yoshida9, 2), # (-1, 1)] f03 = e03Factors[0][0] f03 Cross2*Cross86*Yoshida0*Yoshida2^2*Yoshida4 - Cross1*Cross106*Cross115*Cross84*Yoshida5^2 # Next, we factorize the expressions for f01, f02 and f03 in the ds. If for one of them the factors are all Eckardt quintics or roots, we know these will never vanish, so the cubic is always irreducible. We check that f01 has this property. f01ds = SR(SR(f01).substitute(cross_to_y)).substitute(y_to_d) f01dsFactors = SR(f01ds).factor_list() # We write down the factors: # f01dsFactors # [(d1^2*d3^3 - d2^2*d3^3 + d1*d3^4 - d2*d3^4 - d1^2*d2^2*d4 - d1^2*d2*d3*d4 - d1*d2^2*d3*d4 + d1^2*d3^2*d4 - d1*d2*d3^2*d4 + 2*d1*d3^3*d4 + d3^4*d4 - d1*d2^2*d4^2 - d1*d2*d3*d4^2 + d1*d3^2*d4^2 + d3^3*d4^2 - d1^2*d2^2*d5 - d1^2*d2*d3*d5 - d1*d2^2*d3*d5 + d1^2*d3^2*d5 - d1*d2*d3^2*d5 + 2*d1*d3^3*d5 + d3^4*d5 - 2*d1*d2^2*d4*d5 + d1^2*d3*d4*d5 - 2*d1*d2*d3*d4*d5 - d2^2*d3*d4*d5 + 3*d1*d3^2*d4*d5 - d2*d3^2*d4*d5 + 2*d3^3*d4*d5 - d2^2*d4^2*d5 + d1*d3*d4^2*d5 - d2*d3*d4^2*d5 + d3^2*d4^2*d5 - d1*d2^2*d5^2 - d1*d2*d3*d5^2 + d1*d3^2*d5^2 + d3^3*d5^2 - d2^2*d4*d5^2 + d1*d3*d4*d5^2 - d2*d3*d4*d5^2 + d3^2*d4*d5^2 - d1^2*d2^2*d6 - d1^2*d2*d3*d6 - d1*d2^2*d3*d6 + d1^2*d3^2*d6 - d1*d2*d3^2*d6 + 2*d1*d3^3*d6 + d3^4*d6 - 2*d1*d2^2*d4*d6 + d1^2*d3*d4*d6 - 2*d1*d2*d3*d4*d6 - d2^2*d3*d4*d6 + 3*d1*d3^2*d4*d6 - d2*d3^2*d4*d6 + 2*d3^3*d4*d6 - d2^2*d4^2*d6 + d1*d3*d4^2*d6 - d2*d3*d4^2*d6 + d3^2*d4^2*d6 - 2*d1*d2^2*d5*d6 + d1^2*d3*d5*d6 - 2*d1*d2*d3*d5*d6 - d2^2*d3*d5*d6 + 3*d1*d3^2*d5*d6 - d2*d3^2*d5*d6 + 2*d3^3*d5*d6 - d1^2*d4*d5*d6 - 2*d2^2*d4*d5*d6 + 2*d1*d3*d4*d5*d6 - 2*d2*d3*d4*d5*d6 + 3*d3^2*d4*d5*d6 - d1*d4^2*d5*d6 + d3*d4^2*d5*d6 - d2^2*d5^2*d6 + d1*d3*d5^2*d6 - d2*d3*d5^2*d6 + d3^2*d5^2*d6 - d1*d4*d5^2*d6 + d3*d4*d5^2*d6 - d1*d2^2*d6^2 - d1*d2*d3*d6^2 + d1*d3^2*d6^2 + d3^3*d6^2 - d2^2*d4*d6^2 + d1*d3*d4*d6^2 - d2*d3*d4*d6^2 + d3^2*d4*d6^2 - d2^2*d5*d6^2 + d1*d3*d5*d6^2 - d2*d3*d5*d6^2 + d3^2*d5*d6^2 - d1*d4*d5*d6^2 + d3*d4*d5*d6^2, # 2), # (d1 + d2 + d4, 2), # (d1 + d2 + d5, 2), # (d1 + d2 + d6, 2), # (d1 - d2, 1), # (d1 + d3 + d4, 1), # (d1 + d3 + d5, 2), # (d1 + d3 + d6, 1), # (d1 - d3, 1), # (d1 + d4 + d6, 1), # (d1 - d4, 1), # (d1 - d5, 1), # (d1 - d6, 1), # (d2 + d3 + d4, 1), # (d2 + d3 + d5, 1), # (d2 + d3 + d6, 2), # (d2 - d3, 1), # (d2 + d4 + d5, 1), # (d2 + d4 + d6, 1), # (d2 + d5 + d6, 1), # (d2 - d6, 1), # (d3 + d4 + d5, 1), # (d3 - d5, 1), # (d4 + d5 + d6, 3), # (d4 - d5, 2), # (d4 - d6, 2), # (d5 - d6, 1), # (-1, 1)] # noRoots = [x for x in f01dsFactors if x[0] not in roots] # noNegRoots = [x for x in noRoots if -1*x[0] not in roots] # len(noNegRoots) # 2 # quintics = [x for x in noNegRoots if x[0] not in quintics.values()] # len(noQuintics) # 2 # noNegQuintics = [x for x in noNegRoots if -1*x[0] not in quintics.values()] # noNegQuintics # [(-1,1)] # Conclusion: the factors are Eckardt quintics and roots only. # We write result as a Laurent monomial in Yoshidas and Cross functions, after multiplication by suitable roots. The functions for this are coded in "OperatorsOnYoshidas.sage" exp = SR(f01ds).factor_list() print exp [(d1^2*d3^3 - d2^2*d3^3 + d1*d3^4 - d2*d3^4 - d1^2*d2^2*d4 - d1^2*d2*d3*d4 - d1*d2^2*d3*d4 + d1^2*d3^2*d4 - d1*d2*d3^2*d4 + 2*d1*d3^3*d4 + d3^4*d4 - d1*d2^2*d4^2 - d1*d2*d3*d4^2 + d1*d3^2*d4^2 + d3^3*d4^2 - d1^2*d2^2*d5 - d1^2*d2*d3*d5 - d1*d2^2*d3*d5 + d1^2*d3^2*d5 - d1*d2*d3^2*d5 + 2*d1*d3^3*d5 + d3^4*d5 - 2*d1*d2^2*d4*d5 + d1^2*d3*d4*d5 - 2*d1*d2*d3*d4*d5 - d2^2*d3*d4*d5 + 3*d1*d3^2*d4*d5 - d2*d3^2*d4*d5 + 2*d3^3*d4*d5 - d2^2*d4^2*d5 + d1*d3*d4^2*d5 - d2*d3*d4^2*d5 + d3^2*d4^2*d5 - d1*d2^2*d5^2 - d1*d2*d3*d5^2 + d1*d3^2*d5^2 + d3^3*d5^2 - d2^2*d4*d5^2 + d1*d3*d4*d5^2 - d2*d3*d4*d5^2 + d3^2*d4*d5^2 - d1^2*d2^2*d6 - d1^2*d2*d3*d6 - d1*d2^2*d3*d6 + d1^2*d3^2*d6 - d1*d2*d3^2*d6 + 2*d1*d3^3*d6 + d3^4*d6 - 2*d1*d2^2*d4*d6 + d1^2*d3*d4*d6 - 2*d1*d2*d3*d4*d6 - d2^2*d3*d4*d6 + 3*d1*d3^2*d4*d6 - d2*d3^2*d4*d6 + 2*d3^3*d4*d6 - d2^2*d4^2*d6 + d1*d3*d4^2*d6 - d2*d3*d4^2*d6 + d3^2*d4^2*d6 - 2*d1*d2^2*d5*d6 + d1^2*d3*d5*d6 - 2*d1*d2*d3*d5*d6 - d2^2*d3*d5*d6 + 3*d1*d3^2*d5*d6 - d2*d3^2*d5*d6 + 2*d3^3*d5*d6 - d1^2*d4*d5*d6 - 2*d2^2*d4*d5*d6 + 2*d1*d3*d4*d5*d6 - 2*d2*d3*d4*d5*d6 + 3*d3^2*d4*d5*d6 - d1*d4^2*d5*d6 + d3*d4^2*d5*d6 - d2^2*d5^2*d6 + d1*d3*d5^2*d6 - d2*d3*d5^2*d6 + d3^2*d5^2*d6 - d1*d4*d5^2*d6 + d3*d4*d5^2*d6 - d1*d2^2*d6^2 - d1*d2*d3*d6^2 + d1*d3^2*d6^2 + d3^3*d6^2 - d2^2*d4*d6^2 + d1*d3*d4*d6^2 - d2*d3*d4*d6^2 + d3^2*d4*d6^2 - d2^2*d5*d6^2 + d1*d3*d5*d6^2 - d2*d3*d5*d6^2 + d3^2*d5*d6^2 - d1*d4*d5*d6^2 + d3*d4*d5*d6^2, 2), (d1 + d2 + d4, 2), (d1 + d2 + d5, 2), (d1 + d2 + d6, 2), (d1 - d2, 1), (d1 + d3 + d4, 1), (d1 + d3 + d5, 2), (d1 + d3 + d6, 1), (d1 - d3, 1), (d1 + d4 + d6, 1), (d1 - d4, 1), (d1 - d5, 1), (d1 - d6, 1), (d2 + d3 + d4, 1), (d2 + d3 + d5, 1), (d2 + d3 + d6, 2), (d2 - d3, 1), (d2 + d4 + d5, 1), (d2 + d4 + d6, 1), (d2 + d5 + d6, 1), (d2 - d6, 1), (d3 + d4 + d5, 1), (d3 - d5, 1), (d4 + d5 + d6, 3), (d4 - d5, 2), (d4 - d6, 2), (d5 - d6, 1), (-1, 1)] # The next function selects quintics and linear factors from a factorized expression in ds. def linearAndQuintic(exp): linear_factors = [] quintic_factors = [] for (f,e) in exp: if f.total_degree() == 1: linear_factors.append((f,e)) elif f.total_degree() == 5: quintic_factors.append((f,e)) return (linear_factors, quintic_factors) (linearF, quinticF) = linearAndQuintic([(Qd(x[0]),x[1]) for x in exp]) len(quinticF) [key for key in quintics.keys() if quintics[key] == Qd(-1*quinticF[0][0])] [X32] # ffdict[S(X32)] # [d1^2*d2*d4 + d1*d2^2*d4 + d1^2*d4^2 + d1*d2*d4^2 - d1^2*d2*d5 - d1*d2^2*d5 - d1^2*d5^2 - d1*d2*d5^2 - d2^2*d4*d6 - d2*d4^2*d6 + d2^2*d5*d6 + d2*d5^2*d6 - d2*d4*d6^2 - d4^2*d6^2 + d2*d5*d6^2 + d5^2*d6^2, # d1^2*d2*d4 + d1*d2^2*d4 + d1^2*d4^2 + d1*d2*d4^2 - d2^2*d4*d5 - d2*d4^2*d5 - d2*d4*d5^2 - d4^2*d5^2 - d1^2*d2*d6 - d1*d2^2*d6 + d2^2*d5*d6 + d2*d5^2*d6 - d1^2*d6^2 - d1*d2*d6^2 + d2*d5*d6^2 + d5^2*d6^2, # d1^2*d2*d5 + d1*d2^2*d5 - d2^2*d4*d5 - d2*d4^2*d5 + d1^2*d5^2 + d1*d2*d5^2 - d2*d4*d5^2 - d4^2*d5^2 - d1^2*d2*d6 - d1*d2^2*d6 + d2^2*d4*d6 + d2*d4^2*d6 - d1^2*d6^2 - d1*d2*d6^2 + d2*d4*d6^2 + d4^2*d6^2] test = ffdict[S(X32)] testFactored = [SR(x).factor_list() for x in test] sum([x[1] for x in linearF]) 35 # We check which of the four factors we must divide by to get a monomial in Yoshidas with integer exponents. There are 6 possible combinations since the quintic factor was squared. linearPart = dict() for i in range(0,3): for j in range(i,3): linearPart[(i,j)] = linearF + [(x[0], -1*x[1]) for x in testFactored[i]]+[ (x[0], -1*x[1]) for x in testFactored[j]] linearPart.keys() [(0, 1), (1, 2), (0, 0), (1, 1), (2, 2), (0, 2)] def constructMonomialFromProductRoots(rootFactors, roots, basisOfYoshidas): v_exp = vector(36*[0]) for (f,e) in rootFactors: print (f,e) v_exp = v_exp + e*exponent_vector(DField(SR(f)), roots) print v_exp try: soln = matrix(basisOfYoshidas).transpose() \ vector(v_exp) return soln except ValueError: print 'Expression is not monomial in cube roots of Yoshidas!' return False linearPartFactored =dict() for k in linearPart.keys(): print k linearPartFactored[k] = constructMonomialFromProductRoots(linearPart[k], roots,uyvs_basis) # linearPartFactored # {(0, 0): (5/3, -2, -7/3, 2, -2/3, 1/3, 5/3, 0, 0, 0, 1/3, -2/3, 1, 1, -1/3, 1), # (0, 1): (1, 0, -1, 1, -1, 1, 0, 1, -1, 1, 0, 0, -1, 1, 1, 0), # (0, 2): (4/3, -1, -5/3, 1, -4/3, 2/3, 1/3, 1, -1, 1, 2/3, -1/3, 0, 1, 1/3, 1), # (1, # 1): (1/3, 2, 1/3, 0, -4/3, 5/3, -5/3, 2, -2, 2, -1/3, 2/3, -3, 1, 7/3, -1), # (1, 2): (2/3, 1, -1/3, 0, -5/3, 4/3, -4/3, 2, -2, 2, 1/3, 1/3, -2, 1, 5/3, 0), # (2, 2): (1, 0, -1, 0, -2, 1, -1, 2, -2, 2, 1, 0, -1, 1, 1, 1)} # We pick one with integer exponents, namely (2,2): ratioInY = linearPartFactored[(2,2)] mon = findConstant(ratioInY) # mon # Yoshida21*Yoshida23^2*Yoshida25^2*Yoshida26*Yoshida29*Yoshida30*Yoshida31*Yoshida5/(Yoshida18*Yoshida20^2*Yoshida22*Yoshida24^2*Yoshida28) crossToFind = test[2]*quintics[S(X32)] cross_to_ds = load("Input/cross_to_ds.sobj") [cross for cross in cross_to_ds.keys() if cross_to_ds[cross] == crossToFind] [Cross86] finalMonomial = Cross86^2*mon finalMonomial Cross86^2*Yoshida21*Yoshida23^2*Yoshida25^2*Yoshida26*Yoshida29*Yoshida30*Yoshida31*Yoshida5/(Yoshida18*Yoshida20^2*Yoshida22*Yoshida24^2*Yoshida28) # The expression is correct up to sign. We must check if it is +1 or -1. For this we pick a generic value of d and compare the expressions after evaluation at this choice of d1,..., d6. finalMonNewExpression = SR(SR(finalMonomial.substitute(cross_to_y)).substitute(y_to_d)).substitute(vald) numf01ds = SR(f01ds).substitute(vald) # finalMonNewExpression/numf01ds # 1 # The sign is correct!