diff --git a/src/sage/rings/polynomial/multi_polynomial.pyx b/src/sage/rings/polynomial/multi_polynomial.pyx index ed665f1872d..4d157868272 100644 --- a/src/sage/rings/polynomial/multi_polynomial.pyx +++ b/src/sage/rings/polynomial/multi_polynomial.pyx @@ -31,6 +31,7 @@ from sage.modules.free_module_element import vector from sage.rings.rational_field import QQ from sage.arith.misc import gcd from sage.rings.complex_interval_field import ComplexIntervalField +from sage.rings.real_mpfr import RealField_class,RealField cdef class MPolynomial(CommutativeRingElement): @@ -1874,6 +1875,16 @@ cdef class MPolynomial(CommutativeRingElement): Traceback (most recent call last): ... ValueError: (=x^3*y*z + x^4 + y^2*z) must have two variables + + :: + + sage: R. = PolynomialRing(ZZ) + sage: F = - 8*x^6 - 3933*x^3*y - 725085*x^2*y^2 - 59411592*x*y^3 - 99*y^6 + sage: F.reduced_form(return_conjugation=False) + Traceback (most recent call last): + ... + ValueError: (=-8*x^6 - 99*y^6 - 3933*x^3*y - 725085*x^2*y^2 - + 59411592*x*y^3) must be homogenous """ from sage.matrix.constructor import matrix from sage.calculus.functions import jacobian @@ -1882,7 +1893,8 @@ cdef class MPolynomial(CommutativeRingElement): raise ValueError("(=%s) must have two variables"%self) if self.is_homogeneous() != True: raise ValueError("(=%s) must be homogenous"%self) - KK = ComplexIntervalField(prec=prec) + KK = ComplexIntervalField(prec=prec) # keeps trac of our precision error + JJ=RealField() R = self.parent() S = PolynomialRing(R.base_ring(),'z') phi = R.hom([S.gen(0),1],S)# dehomogenization @@ -1893,6 +1905,7 @@ cdef class MPolynomial(CommutativeRingElement): R=PolynomialRing(KK,'x,y') x,y = R.gens() Q =[] + # finds Stoll and Cremona's Q_0' for j in range(len(roots)): k = (1/(dF(roots[j]).abs()**(2/(n-2))))*((x-(roots[j]*y))*(x-(roots[j].conjugate()*y))) Q.append(k) @@ -1901,13 +1914,15 @@ cdef class MPolynomial(CommutativeRingElement): A = q.monomial_coefficient(x**2) B = q.monomial_coefficient(x*y) C = q.monomial_coefficient(y**2) + # need positive root try: z = (-B + ((B**2)-(4*A*C)).sqrt())/(2*A)# this is z_o except ValueError: raise ValueError("not enough precision") if z.imag()<0: z = (-B - ((B**2)-(4*A*C)).sqrt())/(2*A) - M = matrix(QQ, [[1,0], [0,1]])# this moves Z to fundamental domain + M = matrix(QQ, [[1,0], [0,1]]) # used to keep track of how our z is moved. + # this moves Z to our fundamental domain using the three steps laid out in the algorithim while z.real() < -0.5 or z.real() >= 0.5 or (z.real() <= 0 and z.abs() < 1)\ or (z.real() > 0 and z.abs() <= 1): if z.real() < -0.5: @@ -1923,7 +1938,7 @@ cdef class MPolynomial(CommutativeRingElement): elif (z.real() <= 0 and z.abs() < 1) or (z.real() > 0 and z.abs() <= 1): z = -1/z M = M * matrix(QQ, [[0,-1], [1,0]]) - q = q(tuple((M * vector([x, y])))) + q = q(tuple((M * vector([x, y])))) #our new Q A = q.monomial_coefficient(x**2) B = q.monomial_coefficient(x*y) C = q.monomial_coefficient(y**2) @@ -1944,7 +1959,7 @@ cdef class MPolynomial(CommutativeRingElement): d = (t-(L[j].real()))/((t-(L[j]))*(t-(L[j].conjugate()))+ u**2) c += d # return "A"broke here - #Newton's Method, error bound is while less than diameter of our z instaed of solve use newtons method + #Newton's Method, error bound is while less than diameter of our z err = z0.diameter() zz = z0.diameter() n = q.degree() @@ -1954,6 +1969,7 @@ cdef class MPolynomial(CommutativeRingElement): J = jacobian(G, [u,t]) v0 = vector([z0.imag(), z0.real()]) #z0 as starting point z = z0 + #finds our correct z while err <= zz: NJ = J.subs({u:v0[0], t:v0[1]}) NJinv = NJ.inverse() @@ -1967,8 +1983,9 @@ cdef class MPolynomial(CommutativeRingElement): w = z v0 = v0 - NJinv*G.subs({u:v0[0],t:v0[1]}) z = v0[1].coefficients()[0] + v0[0].coefficients()[0]*KK.gen(0) - err = z.diameter() # prec - zz = (w - z).abs() #diff in w and z + err = z.diameter() # precision + zz = (w - z).abs() #difference in w and z + # moves our z ro fundamental domain as before while z.real() < -0.5 or z.real() >= 0.5 or (z.real() <= 0 and z.abs() < 1)\ or (z.real() > 0 and z.abs() <= 1): if z.real() < -0.5: @@ -1987,7 +2004,7 @@ cdef class MPolynomial(CommutativeRingElement): err = z.diameter() zc = z.center() zr = zc.real() - if (zr.abs() - KK(0.5)).abs() <= err or (zc.abs() - 1).abs() <= err: + if (zr.abs() - JJ(0.5)).abs() <= err or (zc.abs() - 1).abs() <= err: raise ValueError("not enough precision") x,y = self.parent().gens() if return_conjugation: diff --git a/src/sage/schemes/projective/projective_morphism.py b/src/sage/schemes/projective/projective_morphism.py index 4931434541d..b20f8b1276c 100644 --- a/src/sage/schemes/projective/projective_morphism.py +++ b/src/sage/schemes/projective/projective_morphism.py @@ -3525,6 +3525,8 @@ def reduced_form(self, prec=100, return_conjugation=True): # add bool for matri the binary form reduction algorithm from Stoll and Cremona _[SC] to the periodic points (the dynatomic polynomial). + Implimented by Rebecca Lauren Miller as part of GSOC 2016. + INPUT:: - ``prec`` -- desired precision (default: 100) @@ -3559,16 +3561,15 @@ def reduced_form(self, prec=100, return_conjugation=True): # add bool for matri ) :: - - sage: PS. = ProjectiveSpace(QQ,1) + # failing why? + sage: PS. = ProjectiveSpace(ZZ,1) sage: H = End(PS) sage: f = H([x^2+ x*y,y^2]) #needs 3 periodic sage: m = matrix(QQ, 2, 2,[-221, -1, 1, 0]) sage: f = f.conjugate(m) sage: f.reduced_form(prec=200) ( - Scheme endomorphism of Projective Space of dimension 1 over Rational - Field + Scheme endomorphism of Projective Space of dimension 1 over Integer Ring Defn: Defined on coordinates by sending (x : y) to (-x^2 + x*y - y^2 : -y^2) , @@ -3592,22 +3593,34 @@ def reduced_form(self, prec=100, return_conjugation=True): # add bool for matri [-1 2] [ 2 -5] ) + + :: + # failing why? + sage: PS. = ProjectiveSpace(ZZ,1) + sage: H = End(PS) + sage: f = H([7*x^3 + 5*x*y^2 + y, 6*y^3]) + sage: m = matrix(QQ, 2, 2,[1 ,35, 0, 1])*matrix(QQ, 2, 2,[0, 1, -1, 0])*matrix(QQ,2 ,2, [1,0,131,1]) + sage: f = f.conjugate(m) + sage: f.reduced_form(prec=200) + Traceback (most recent call last): + ... + ValueError: polys (=[7*x^3 + 5*x*y^2 + y, 6*y^3]) must be homogeneous """ R = self.coordinate_ring() F = R(self.dynatomic_polynomial(1)) x,y = R.gens() - #F.parent() == R #needs ticket 21097, - F = R(F/gcd(F,F.derivative(x))) #remove multiple roots + F = R(F/gcd(F,F.derivative(x))) #removes multiple roots n = 2 - while F.degree() <= 2: # checks to have enough distinct roots need at least 3 - F = self.dynatomic_polynomial(n) # n period fix points - F = R(F/gcd(F,F.derivative(x)))#remove multiple roots + # Checks to make sure there are enough distinct, roots we need 3 + # if there are not if finds the nth period fix points until there are enough + while F.degree() <= 2: + F = self.dynatomic_polynomial(n) # finds n period fix points + F = R(F/gcd(F,F.derivative(x)))#removes multiple roots n += 1 G,m = F.reduced_form( prec=prec, return_conjugation=return_conjugation) if return_conjugation: return (self.conjugate(m), m) return self.conjugate(m) - # answers to examples have fractions?, 1st example now doenst have enough prec!, problem raising value error class SchemeMorphism_polynomial_projective_space_field(SchemeMorphism_polynomial_projective_space):