def mynumerator(x): if parent(x) == R: return x return numerator(x) class fastfrac: def __init__(self,top,bot=1): if parent(top) == ZZ or parent(top) == R: self.top = R(top) self.bot = R(bot) elif top.__class__ == fastfrac: self.top = top.top self.bot = top.bot * bot else: self.top = R(numerator(top)) self.bot = R(denominator(top)) * bot def reduce(self): return fastfrac(self.top / self.bot) def sreduce(self): return fastfrac(I.reduce(self.top),I.reduce(self.bot)) def iszero(self): return self.top in I and not (self.bot in I) def isdoublingzero(self): return self.top in J and not (self.bot in J) def __add__(self,other): if parent(other) == ZZ: return fastfrac(self.top + self.bot * other,self.bot) if other.__class__ == fastfrac: return fastfrac(self.top * other.bot + self.bot * other.top,self.bot * other.bot) return NotImplemented def __sub__(self,other): if parent(other) == ZZ: return fastfrac(self.top - self.bot * other,self.bot) if other.__class__ == fastfrac: return fastfrac(self.top * other.bot - self.bot * other.top,self.bot * other.bot) return NotImplemented def __neg__(self): return fastfrac(-self.top,self.bot) def __mul__(self,other): if parent(other) == ZZ: return fastfrac(self.top * other,self.bot) if other.__class__ == fastfrac: return fastfrac(self.top * other.top,self.bot * other.bot) return NotImplemented def __rmul__(self,other): return self.__mul__(other) def __div__(self,other): if parent(other) == ZZ: return fastfrac(self.top,self.bot * other) if other.__class__ == fastfrac: return fastfrac(self.top * other.bot,self.bot * other.top) return NotImplemented def __pow__(self,other): if parent(other) == ZZ: return fastfrac(self.top ^ other,self.bot ^ other) return NotImplemented def isidentity(x): return x.iszero() def isdoublingidentity(x): return x.isdoublingzero() R. = PolynomialRing(QQ,23,order='invlex') I = R.ideal([ mynumerator((ux1^3+uy1^3+1)-(3*ud*ux1*uy1)) , mynumerator((ux1)-(uX1/uZ1)) , mynumerator((uy1)-(uY1/uZ1)) , mynumerator((uXX1)-(uX1*uX1)) , mynumerator((uYY1)-(uY1*uY1)) , mynumerator((uZZ1)-(uZ1*uZ1)) , mynumerator((uXY1)-(2*uX1*uY1)) , mynumerator((uXZ1)-(2*uX1*uZ1)) , mynumerator((uYZ1)-(2*uY1*uZ1)) , mynumerator((ux2^3+uy2^3+1)-(3*ud*ux2*uy2)) , mynumerator((ux2)-(uX2/uZ2)) , mynumerator((uy2)-(uY2/uZ2)) , mynumerator((uXX2)-(uX2*uX2)) , mynumerator((uYY2)-(uY2*uY2)) , mynumerator((uZZ2)-(uZ2*uZ2)) , mynumerator((uXY2)-(2*uX2*uY2)) , mynumerator((uXZ2)-(2*uX2*uZ2)) , mynumerator((uYZ2)-(2*uY2*uZ2)) , mynumerator((uZ2)-(1)) ]) J = I + R.ideal([0 , uXX1-uXX2 , uXY1-uXY2 , uXZ1-uXZ2 , uZZ1-uZZ2 , uYY1-uYY2 , uYZ1-uYZ2 , uX1-uX2 , uY1-uY2 , uZ1-uZ2 ]) ud = fastfrac(ud) ux2 = fastfrac(ux2) uy2 = fastfrac(uy2) ux1 = fastfrac(ux1) uy1 = fastfrac(uy1) uZZ1 = fastfrac(uZZ1) uYY1 = fastfrac(uYY1) uXX1 = fastfrac(uXX1) uY1 = fastfrac(uY1) uXZ1 = fastfrac(uXZ1) uYZ1 = fastfrac(uYZ1) uX1 = fastfrac(uX1) uXY1 = fastfrac(uXY1) uZ1 = fastfrac(uZ1) uZZ2 = fastfrac(uZZ2) uYY2 = fastfrac(uYY2) uXX2 = fastfrac(uXX2) uY2 = fastfrac(uY2) uXZ2 = fastfrac(uXZ2) uYZ2 = fastfrac(uYZ2) uX2 = fastfrac(uX2) uXY2 = fastfrac(uXY2) uZ2 = fastfrac(uZ2) uX3 = ((uYY1*uXZ2-uXZ1*uYY2)) uY3 = ((uXX1*uYZ2-uYZ1*uXX2)) uZ3 = ((uZZ1*uXY2-uXY1)) uXX3 = ((uX3^2)) uYY3 = ((uY3^2)) uZZ3 = ((uZ3^2)) uXY3 = (((uX3+uY3)^2-uXX3-uYY3)) uXZ3 = (((uX3+uZ3)^2-uXX3-uZZ3)) uYZ3 = (((uY3+uZ3)^2-uYY3-uZZ3)) ux3 = (((uy1^2*ux2-uy2^2*ux1)/(ux2*uy2-ux1*uy1))).reduce() uy3 = (((ux1^2*uy2-ux2^2*uy1)/(ux2*uy2-ux1*uy1))).reduce() print isidentity((ux3^3+uy3^3+fastfrac(1))-(fastfrac(3)*ud*ux3*uy3)) print isidentity((ux3)-(uX3/uZ3)) print isidentity((uy3)-(uY3/uZ3)) print isidentity((uXX3)-(uX3*uX3)) print isidentity((uYY3)-(uY3*uY3)) print isidentity((uZZ3)-(uZ3*uZ3)) print isidentity((uXY3)-(fastfrac(2)*uX3*uY3)) print isidentity((uXZ3)-(fastfrac(2)*uX3*uZ3)) print isidentity((uYZ3)-(fastfrac(2)*uY3*uZ3)) unified = True uXX4 = uXX3 uXY4 = uXY3 uXZ4 = uXZ3 uZZ4 = uZZ3 uYY4 = uYY3 uYZ4 = uYZ3 uX4 = uX3 uY4 = uY3 uZ4 = uZ3 ux4 = ((uy1*(fastfrac(1)-ux1^3)/(ux1^3-uy1^3))).reduce() uy4 = ((ux1*(uy1^3-fastfrac(1))/(ux1^3-uy1^3))).reduce() if unified: unified = isdoublingidentity((ux4^3+uy4^3+fastfrac(1))-(fastfrac(3)*ud*ux4*uy4)) if unified: unified = isdoublingidentity((ux4)-(uX4/uZ4)) if unified: unified = isdoublingidentity((uy4)-(uY4/uZ4)) if unified: unified = isdoublingidentity((uXX4)-(uX4*uX4)) if unified: unified = isdoublingidentity((uYY4)-(uY4*uY4)) if unified: unified = isdoublingidentity((uZZ4)-(uZ4*uZ4)) if unified: unified = isdoublingidentity((uXY4)-(fastfrac(2)*uX4*uY4)) if unified: unified = isdoublingidentity((uXZ4)-(fastfrac(2)*uX4*uZ4)) if unified: unified = isdoublingidentity((uYZ4)-(fastfrac(2)*uY4*uZ4)) if unified: print "Unified"