diff --git a/src/sage/rings/multi_polynomial_element.py b/src/sage/rings/multi_polynomial_element.py index e880970e043..e62a4727d17 100644 --- a/src/sage/rings/multi_polynomial_element.py +++ b/src/sage/rings/multi_polynomial_element.py @@ -237,6 +237,9 @@ def __rpow__(self, n): def element(self): return self.__element + def change_ring(self, R): + return self.parent().change_ring(R)(self) + class MPolynomial_macaulay2_repr: """ diff --git a/src/sage/rings/multi_polynomial_ring.py b/src/sage/rings/multi_polynomial_ring.py index b3b2dd23fc7..fc63254098a 100644 --- a/src/sage/rings/multi_polynomial_ring.py +++ b/src/sage/rings/multi_polynomial_ring.py @@ -329,6 +329,10 @@ def latex_variable_names(self): self.__latex_variable_names = names return names + def change_ring(self, R): + from polynomial_ring_constructor import PolynomialRing + return PolynomialRing(R, self.variable_names(), order=self.term_order()) + class MPolynomialRing_polydict( MPolynomialRing_macaulay2_repr, MPolynomialRing_generic): """ Multivariable polynomial ring. diff --git a/src/sage/rings/power_series_ring_element.pyx b/src/sage/rings/power_series_ring_element.pyx index 057dd274238..ea7cf4eaf5d 100644 --- a/src/sage/rings/power_series_ring_element.pyx +++ b/src/sage/rings/power_series_ring_element.pyx @@ -837,6 +837,8 @@ cdef class PowerSeries(AlgebraElement): t^5 + O(t^25) """ denom = denom_r + if denom.is_zero(): + raise ZeroDivisionError, "Can't divide by something indistinguishable from 0" u = denom.valuation_zero_part() inv = ~u # inverse @@ -980,6 +982,8 @@ cdef class PowerSeries(AlgebraElement): 1 + t + t^2 + 7*t^3 + O(t^50) sage: sqrt(K(0)) 0 + sage: sqrt(t^2) + t AUTHOR: -- Robert Bradshaw @@ -991,8 +995,12 @@ cdef class PowerSeries(AlgebraElement): if val is not infinity and val % 2 == 1: raise ValueError, "Square root not defined for power series of odd valuation." + if self.degree() == 0: + x = self.valuation_zero_part()[0].sqrt() + return self.parent().base_extend(x.parent())([x], val/2) + prec = self.prec() - if prec is infinity and self.degree() > 0: + if prec == infinity: prec = self._parent.default_prec() a = self.valuation_zero_part() @@ -1045,6 +1053,8 @@ cdef class PowerSeries(AlgebraElement): raise ValueError, "Square root not defined for power series of odd valuation." elif not self[val].is_square(): raise ValueError, "Square root does not live in this ring." + elif is_Field(self.base_ring()): + return self.sqrt() else: try: return self.parent()(self.sqrt()) diff --git a/src/sage/rings/rational.pyx b/src/sage/rings/rational.pyx index b0170574d9e..655278a1650 100644 --- a/src/sage/rings/rational.pyx +++ b/src/sage/rings/rational.pyx @@ -163,6 +163,7 @@ cdef class Rational(sage.structure.element.FieldElement): def __set_value(self, x, unsigned int base): cdef int n + cdef Rational temp_rational if isinstance(x, Rational): set_from_Rational(self, x) @@ -236,6 +237,10 @@ cdef class Rational(sage.structure.element.FieldElement): if n or mpz_cmp_si(mpq_denref(self.value), 0) == 0: raise TypeError, "Unable to coerce %s (%s) to Rational"%(x,type(x)) + elif hasattr(x, 'rational_reconstruction'): + temp_rational = x.rational_reconstruction() + mpq_set(self.value, temp_rational.value) + else: raise TypeError, "Unable to coerce %s (%s) to Rational"%(x,type(x)) diff --git a/src/sage/schemes/elliptic_curves/ell_padic_field.py b/src/sage/schemes/elliptic_curves/ell_padic_field.py index 0fac5f57f9b..a3cf6f67617 100644 --- a/src/sage/schemes/elliptic_curves/ell_padic_field.py +++ b/src/sage/schemes/elliptic_curves/ell_padic_field.py @@ -25,10 +25,15 @@ from sage.rings.all import PowerSeriesRing, PolynomialRing, IntegerModRing, ZZ, QQ from sage.misc.functional import ceil, log +# Elliptic curves are very different than genus > 1 hyperelliptic curves, +# there is an "is a" relationship here, and common implementation with regard +# Coleman integration. +from sage.schemes.hyperelliptic_curves.hyperelliptic_padic_field import HyperellipticCurve_padic_field + import sage.databases.cremona -class EllipticCurve_padic_field(EllipticCurve_field): +class EllipticCurve_padic_field(EllipticCurve_field, HyperellipticCurve_padic_field): """ Elliptic curve over a padic field. """ @@ -40,7 +45,7 @@ def __init__(self, x, y=None): else: if isinstance(y, str): field = x - X = sage.databases.cremona.CremonaDatabase()[label] + X = sage.databases.cremona.CremonaDatabase()[y] ainvs = [field(a) for a in X.a_invariants()] else: field = x @@ -51,6 +56,7 @@ def __init__(self, x, y=None): EllipticCurve_field.__init__(self, [field(x) for x in ainvs]) self._point_class = ell_point.EllipticCurvePoint_field + self._genus = 1 def _pari_(self): try: @@ -62,63 +68,6 @@ def _pari_(self): return self.__pari -# The functions below were prototyped at the 2007 Arizona Winter School by -# Robert Bradshaw and Ralf Gerkmann, working with Miljan Brakovevic and -# Kiran Kedlaya -# All of the below is with respect to the Monsky Washnitzer cohomology. - - def local_analytic_interpolation(self, P, Q): - """ - Construct a linear interpolation from P to Q in a power series - in the local parameter t, with precision equal to the p-adic - precision of the underlying ring. - - Returns a point $X(t) = ( x(t) : y(t) : z(t) )$ such that - $X(0) = P and X(1) = Q$ - - This is implemented by linearly interpolating x, solving formally for y, - and letting z(t) = 1. - - For this to make sense, P and Q must be in the same residue series, neither equal to infinity. - - TODO: remove the last condition? - """ - prec = self.base_ring().precision_cap() - t = PowerSeriesRing(self.base_ring(), 't', prec).gen(0) - x = P[0]+t*(Q[0]-P[0]) - pts = self.lift_x(x) - if (pts[0][1] - P[1]).valuation() > 0: - return pts[0] - else: - return pts[1] - - def tiny_integrals(self, F, P, Q): - """ - Evaluate the integrals of $f_i dx/y$ from P to Q for each f_i in F - by formally integrating a power series in a local parameter $t$ - - P and Q MUST be in the same residue disk for this result to make sense. - """ - x, y, z = self.local_analytic_interpolation(P, Q) - dt = x.derivative() / y - integrals = [] - for f in F: - f = f(x,y) - I = (f*dt).integral() - integrals.append(I(1)) - return integrals - - def tiny_integrals_on_basis(self, P, Q): - """ - Evaluate the integrals of $dx/y$ and $x dx/y$ - by formally integrating a power series in a local parameter $t$ - - P and Q MUST be in the same residue disk for this result to make sense. - """ - R = PolynomialRing(self.base_ring(), ['x', 'y']) - return self.tiny_integrals([R(1), R.gen(0)], P, Q) - - def frobenius(self, P=None): try: @@ -152,35 +101,6 @@ def frob(P): else: return frob(P) - def teichmuller(self, P): - """ - Find a Teichm\:uller point in the same residue class of P. - - Because this lift of frobenius acts as $x \mapsto x^p$, - take the Teichmuler lift of $x$ and then find a matching y - from that. - - EXAMPLES: - sage: K = pAdicField(7, 5) - sage: E = EllipticCurve(K, [-31/3, -2501/108]) # 11a - sage: P = E(K(14/3), K(11/2)) - sage: E.frobenius(P) == P - False - sage: TP = E.teichmuller(P); TP - (0 : 2 + 3*7 + 3*7^2 + 3*7^4 + O(7^5) : 1 + O(7^5)) - sage: E.frobenius(TP) == TP - True - sage: (TP[0] - P[0]).valuation() > 0, (TP[1] - P[1]).valuation() > 0 - (True, True) - """ - x = padic_teichmuller(P[0]) - pts = self.lift_x(x) - if (pts[0][1] - P[1]).valuation() > 0: - return pts[0] - else: - return pts[1] - - def coleman_integrals_on_basis(self, P, Q): """ Return the coleman integral of dx/y and x dx/y from P to Q. @@ -254,16 +174,3 @@ def coleman_integrals_on_basis(self, P, Q): return P_to_TP + TP_to_TQ + TQ_to_Q - -# TODO: add this to padics (if it isn't there in the new version already). -def padic_teichmuller(a): - K = a.parent() - p = K.prime() - p_less_1_inverse = ~K(p-1) - one = K(1) - x = K(ZZ(a) % p) - if x == 0: - return x - for _ in range(ceil(log(K.precision_cap(),2))): - x = ((p-2)*x + x**(2-p)) * p_less_1_inverse - return x diff --git a/src/sage/schemes/elliptic_curves/ell_rational_field.py b/src/sage/schemes/elliptic_curves/ell_rational_field.py index d2e851227d8..79dd48c5cb9 100644 --- a/src/sage/schemes/elliptic_curves/ell_rational_field.py +++ b/src/sage/schemes/elliptic_curves/ell_rational_field.py @@ -4404,8 +4404,7 @@ def padic_sigma(self, p, N=20, E2=None, check=False, check_hypotheses=True): - def padic_E2(self, p, prec=20, check=False, check_hypotheses=True, - algorithm="auto"): + def padic_E2(self, p, prec=20, check=False, check_hypotheses=True, algorithm="auto"): r""" Returns the value of the $p$-adic modular form $E2$ for $(E, \omega)$ where $\omega$ is the usual invariant differential @@ -4531,6 +4530,38 @@ def padic_E2(self, p, prec=20, check=False, check_hypotheses=True, 1907 + 2819*3001 + 1124*3001^2 + O(3001^3) """ + frob_p = self.matrix_of_frobenius(p, prec, check, check_hypotheses, algorithm).change_ring(Integers(p**prec)) + + frob_p_n = frob_p**prec + + # todo: write a coercion operator from integers mod p^n to the + # p-adic field (it doesn't seem to currently exist) + # see trac #4 + + # todo: think about the sign of this. Is it correct? + + E2_of_X = output_ring( (-12 * frob_p_n[0,1] / frob_p_n[1,1]).lift() ) \ + + O(p**prec) + + # Take into account the coordinate change. + fudge_factor = (X.discriminant() / self.discriminant()).nth_root(6) + + # todo: here I should be able to write: + # return E2_of_X / fudge_factor + # However, there is a bug in SAGE (#51 on trac) which makes this + # crash sometimes when prec == 1. For example, + # EllipticCurve([1, 1, 1, 1, 1]).padic_E2(5, 1) + # makes it crash. I haven't figured out exactly what the bug + # is yet, but for now I use the following workaround: + fudge_factor_inverse = Qp(p, prec=(E2_of_X.precision_absolute() + 1))(1 / fudge_factor) + return output_ring(E2_of_X * fudge_factor_inverse) + + + def matrix_of_frobenius(self, p, prec=20, check=False, check_hypotheses=True, algorithm="auto"): + """ + See the parameters and documentation for padic_E2. + """ + # TODO, add lots of comments like the above if check_hypotheses: p = self.__check_padic_hypotheses(p) @@ -4617,29 +4648,7 @@ def padic_E2(self, p, prec=20, check=False, check_hypotheses=True, "Consistency check failed! (correct = %s, actual = %s)" % \ (correct_trace, trace_of_frobenius) - frob_p_n = frob_p**prec - - # todo: write a coercion operator from integers mod p^n to the - # p-adic field (it doesn't seem to currently exist) - # see trac #4 - - # todo: think about the sign of this. Is it correct? - - E2_of_X = output_ring( (-12 * frob_p_n[0,1] / frob_p_n[1,1]).lift() ) \ - + O(p**prec) - - # Take into account the coordinate change. - fudge_factor = (X.discriminant() / self.discriminant()).nth_root(6) - - # todo: here I should be able to write: - # return E2_of_X / fudge_factor - # However, there is a bug in SAGE (#51 on trac) which makes this - # crash sometimes when prec == 1. For example, - # EllipticCurve([1, 1, 1, 1, 1]).padic_E2(5, 1) - # makes it crash. I haven't figured out exactly what the bug - # is yet, but for now I use the following workaround: - fudge_factor_inverse = Qp(p, prec=(E2_of_X.precision_absolute() + 1))(1 / fudge_factor) - return output_ring(E2_of_X * fudge_factor_inverse) + return frob_p.change_ring(Zp(p, prec)) # This is the old version of padic_E2 that requires MAGMA: diff --git a/src/sage/schemes/elliptic_curves/monsky_washnitzer.py b/src/sage/schemes/elliptic_curves/monsky_washnitzer.py index 29040afebfb..15f02a27cd6 100644 --- a/src/sage/schemes/elliptic_curves/monsky_washnitzer.py +++ b/src/sage/schemes/elliptic_curves/monsky_washnitzer.py @@ -2021,9 +2021,17 @@ def matrix_of_frobenius_alternate(a, b, p, N): from sage.misc.profiler import Profiler from sage.misc.misc import repr_lincomb -def matrix_of_frobenius_hyperelliptic(Q, p, prec, M=None): +def matrix_of_frobenius_hyperelliptic(Q, p=None, prec=None, M=None): prof = Profiler() prof("setup") + if p is None: + try: + K = Q.base_ring() + p = K.prime() + prec = K.precision_cap() + #prec -= (adjusted_prec(p, prec) - prec) + except AttributeError: + raise ValueError, "p and prec must be specified if Q is not defined over a p-adic ring" if M is None: M = adjusted_prec(p, prec) extra_prec_ring = Integers(p**M) # pAdicField(p, M) # SLOW! @@ -2037,6 +2045,7 @@ def matrix_of_frobenius_hyperelliptic(Q, p, prec, M=None): # do reduction over Q in case we have non-integral entries (and it's so much faster than padics) rational_S = S.change_ring(QQ) # this is a hack until pAdics are fast + # (They are in the latest development bundle, but its not standard and I'd need to merge. # (it will periodically cast into this ring to reduce coefficent size) rational_S._prec_cap = p**M rational_S._p = p @@ -2060,7 +2069,7 @@ def matrix_of_frobenius_hyperelliptic(Q, p, prec, M=None): M += real_prec_ring(0).add_bigoh(prec) print prof # print len(S._monomials) - return M.transpose() + return M.transpose(), [f for f, a in reduced] @@ -2091,7 +2100,10 @@ def MonskyWashnitzerDifferentialRing(base_ring): class SpecialHyperellipticQuotientRing_class(CommutativeAlgebra): - def __init__(self, Q, R, invert_y=False): + def __init__(self, Q, R=None, invert_y=True): + if R is None: + R = Q.base_ring() + CommutativeAlgebra.__init__(self, R) x = PolynomialRing(R, 'x').gen(0) @@ -2099,13 +2111,13 @@ def __init__(self, Q, R, invert_y=False): E = Q if E.a1() != 0 or E.a2() != 0: raise NotImplementedError, "Curve must be in Weierstrass normal form." - Q = E.defining_polynomial()(x,0,1) + Q = (-E.defining_polynomial()).change_ring(R)(x,0,1) elif is_HyperellipticCurve(Q): C = Q if C.hyperelliptic_polynomials()[1] != 0: raise NotImplementedError, "Curve must be of form y^2 = Q(x)." - Q = E.hyperelliptic_polynomials()[0]()(x) + Q = C.hyperelliptic_polynomials()[0].change_ring(R) if is_Polynomial(Q): self._Q = Q.change_ring(R) @@ -2262,7 +2274,7 @@ def __init__(self, parent, val=0, offset=0): if isinstance(val, tuple): val, offset = val if isinstance(val, SpecialHyperellipticQuotientElement): - val = val.coeffs() + val, offset = val.coeffs() if isinstance(val, list) and len(val) > 0 and is_FreeModuleElement(val[0]): val = transpose_list(val) self._f = parent._poly_ring(val) @@ -2272,6 +2284,9 @@ def __init__(self, parent, val=0, offset=0): def change_ring(self, R): return self.parent().change_ring(R)(self.coeffs()) + def __call__(self, *x): + return self._f(*x) + def __invert__(self): """ The general element in our ring is not invertible, but y may be. diff --git a/src/sage/schemes/generic/homset.py b/src/sage/schemes/generic/homset.py index 029ae664aa7..141fe699623 100644 --- a/src/sage/schemes/generic/homset.py +++ b/src/sage/schemes/generic/homset.py @@ -3,6 +3,7 @@ """ import sage.rings.integer_ring +Z = sage.rings.integer_ring.ZZ # Some naive point enumeration routines for default. # AUTHOR: David R. Kohel diff --git a/src/sage/schemes/hyperelliptic_curves/constructor.py b/src/sage/schemes/hyperelliptic_curves/constructor.py index c5f9d8ff614..77154987689 100644 --- a/src/sage/schemes/hyperelliptic_curves/constructor.py +++ b/src/sage/schemes/hyperelliptic_curves/constructor.py @@ -13,11 +13,13 @@ from hyperelliptic_generic import HyperellipticCurve_generic from hyperelliptic_finite_field import HyperellipticCurve_finite_field from hyperelliptic_rational_field import HyperellipticCurve_rational_field +from hyperelliptic_padic_field import HyperellipticCurve_padic_field from hyperelliptic_g2_generic import HyperellipticCurve_g2_generic from hyperelliptic_g2_finite_field import HyperellipticCurve_g2_finite_field from hyperelliptic_g2_rational_field import HyperellipticCurve_g2_rational_field +from hyperelliptic_g2_padic_field import HyperellipticCurve_g2_padic_field -from sage.rings.all import is_FiniteField, is_RationalField, is_Polynomial +from sage.rings.all import is_FiniteField, is_RationalField, is_Polynomial, is_pAdicField def HyperellipticCurve(f,h=None,names=None,PP=None): r""" @@ -67,6 +69,11 @@ def HyperellipticCurve(f,h=None,names=None,PP=None): return HyperellipticCurve_g2_rational_field(PP, f, h, names=names, genus=g) else: return HyperellipticCurve_rational_field(PP, f, h, names=names, genus=g) + elif is_pAdicField(R): + if g == 2: + return HyperellipticCurve_g2_padic_field(PP, f, h, names=names, genus=g) + else: + return HyperellipticCurve_padic_field(PP, f, h, names=names, genus=g) else: if g == 2: return HyperellipticCurve_g2_generic(PP, f, h, names=names, genus=g) diff --git a/src/sage/schemes/hyperelliptic_curves/hyperelliptic_g2_padic_field.py b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_g2_padic_field.py new file mode 100644 index 00000000000..ff2abde395a --- /dev/null +++ b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_g2_padic_field.py @@ -0,0 +1,16 @@ +""" +Hyperelliptic curves over the rationals +""" + +#***************************************************************************** +# Copyright (C) 2007 Robert Bradshaw +# Distributed under the terms of the GNU General Public License (GPL) +# http://www.gnu.org/licenses/ +#***************************************************************************** + +import hyperelliptic_g2_generic, hyperelliptic_padic_field + +class HyperellipticCurve_g2_padic_field( + hyperelliptic_g2_generic.HyperellipticCurve_g2_generic, + hyperelliptic_padic_field.HyperellipticCurve_padic_field): + pass diff --git a/src/sage/schemes/hyperelliptic_curves/hyperelliptic_generic.py b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_generic.py index 832a080e8cb..62ae2dfcbf3 100644 --- a/src/sage/schemes/hyperelliptic_curves/hyperelliptic_generic.py +++ b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_generic.py @@ -55,6 +55,13 @@ def __init__(self, PP, f, h=None, names=None, genus=None): self._hyperelliptic_polynomials = (f,h) self._genus = genus + def chage_ring(self, R): + from constructor import HyperellipticCurve + f, h = self._hyperelliptic_polynomials + y = self._printing_ring.gen() + x = self._printing_ring.base_ring().gen() + return HyperellipticCurve(f.change_ring(R), h, "%s,%s"%(x,y)) + def _repr_(self): """ String representation hyperelliptic curves. @@ -85,6 +92,32 @@ def hyperelliptic_polynomials(self, K=None, var='x'): P = PolynomialRing(K, var) return (P(f),P(h)) + def lift_x(self, x): + f, h = self._hyperelliptic_polynomials + x += self.base_ring()(0) + one = x.parent()(1) + if h.is_zero(): + y2 = f(x) + if y2.is_zero(): + return self.point([x, y2, one], check=False) + elif y2.is_square(): + y = y2.square_root() + return [self.point([x, y, one], check=False), + self.point([x,-y, one], check=False)] + else: + return [] + else: + b = h(x) + D = b*b + 4*f(x) + if D.is_zero(): + return [self.point([x, -b/2, one], check=False)] + elif D.is_square(): + sqrtD = D.square_root() + return [self.point([x, (-b+sqrtD)/2, one], check=False), + self.point([x, (-b-sqrtD)/2, one], check=False)] + else: + return [] + def genus(self): return self._genus diff --git a/src/sage/schemes/hyperelliptic_curves/hyperelliptic_padic_field.py b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_padic_field.py new file mode 100644 index 00000000000..34126d8a910 --- /dev/null +++ b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_padic_field.py @@ -0,0 +1,200 @@ +""" +Hyperelliptic curves over a padic field. +""" + +#***************************************************************************** +# Copyright (C) 2007 Robert Bradshaw +# Distributed under the terms of the GNU General Public License (GPL) +# http://www.gnu.org/licenses/ +#***************************************************************************** + +import hyperelliptic_generic + +from sage.rings.all import PowerSeriesRing, PolynomialRing, ZZ, QQ +from sage.misc.functional import ceil, log + +class HyperellipticCurve_padic_field(hyperelliptic_generic.HyperellipticCurve_generic): + + +# The functions below were prototyped at the 2007 Arizona Winter School by +# Robert Bradshaw and Ralf Gerkmann, working with Miljan Brakovevic and +# Kiran Kedlaya +# All of the below is with respect to the Monsky Washnitzer cohomology. + + def local_analytic_interpolation(self, P, Q): + """ + Construct a linear interpolation from P to Q in a power series + in the local parameter t, with precision equal to the p-adic + precision of the underlying ring. + + Returns a point $X(t) = ( x(t) : y(t) : z(t) )$ such that + $X(0) = P and X(1) = Q$ + + This is implemented by linearly interpolating x, solving formally for y, + and letting z(t) = 1. + + For this to make sense, P and Q must be in the same residue series, neither equal to infinity. + + TODO: remove the last condition? + """ + prec = self.base_ring().precision_cap() + t = PowerSeriesRing(self.base_ring(), 't', prec).gen(0) + x = P[0]+t*(Q[0]-P[0]) + pts = self.lift_x(x) + if (pts[0][1] - P[1]).valuation() > 0: + return pts[0] + else: + return pts[1] + + def tiny_integrals(self, F, P, Q): + """ + Evaluate the integrals of $f_i dx/y$ from P to Q for each f_i in F + by formally integrating a power series in a local parameter $t$ + + P and Q MUST be in the same residue disk for this result to make sense. + """ + x, y, z = self.local_analytic_interpolation(P, Q) + dt = x.derivative() / y + integrals = [] + for f in F: + f = f(x,y) + I = (f*dt).integral() + integrals.append(I(1)) + return integrals + + def tiny_integrals_on_basis(self, P, Q): + """ + Evaluate the integrals of $dx/y$ and $x dx/y$ + by formally integrating a power series in a local parameter $t$ + + P and Q MUST be in the same residue disk for this result to make sense. + """ + R = PolynomialRing(self.base_ring(), ['x', 'y']) + x, y = R.gens() + return self.tiny_integrals([x**i for i in range(2*self.genus())], P, Q) + + + def teichmuller(self, P): + """ + Find a Teichm\:uller point in the same residue class of P. + + Because this lift of frobenius acts as $x \mapsto x^p$, + take the Teichmuler lift of $x$ and then find a matching y + from that. + + EXAMPLES: + sage: K = pAdicField(7, 5) + sage: E = EllipticCurve(K, [-31/3, -2501/108]) # 11a + sage: P = E(K(14/3), K(11/2)) + sage: E.frobenius(P) == P + False + sage: TP = E.teichmuller(P); TP + (0 : 2 + 3*7 + 3*7^2 + 3*7^4 + O(7^5) : 1 + O(7^5)) + sage: E.frobenius(TP) == TP + True + sage: (TP[0] - P[0]).valuation() > 0, (TP[1] - P[1]).valuation() > 0 + (True, True) + """ + x = padic_teichmuller(P[0]) + pts = self.lift_x(x) + if (pts[0][1] - P[1]).valuation() > 0: + return pts[0] + else: + return pts[1] + + def coleman_integrals_on_basis(self, P, Q): + import sage.schemes.elliptic_curves.monsky_washnitzer as monsky_washnitzer + from sage.misc.profiler import Profiler + prof = Profiler() + prof("setup") + K = self.base_ring() + p = K.prime() + from sage.modules.free_module import VectorSpace + from sage.matrix.constructor import matrix + dim = 2*self.genus() + V = VectorSpace(K, dim) + + prof("tiny integrals") + TP = self.teichmuller(P) +# print "TP", TP + if P == TP: + P_to_TP = V(0) + else: + P_to_TP = V(self.tiny_integrals_on_basis(P, TP)) +# print " P to TP:", P_to_TP[0] + + TQ = self.teichmuller(Q) +# print "TQ", TQ + if Q == TQ: + TQ_to_Q = V(0) + else: + TQ_to_Q = V(self.tiny_integrals_on_basis(TQ, Q)) +# print "TQ to Q:", TQ_to_Q[0] + + prof("mw calc") + try: + M_frob, forms = self._frob_calc + except AttributeError: + M_frob, forms = self._frob_calc = monsky_washnitzer.matrix_of_frobenius_hyperelliptic(self) + + prof("eval f") + R = forms[0].base_ring() + L = [f(R(TP[0]), R(TP[1])) - f(R(TQ[0]), R(TQ[1])) for f in forms] +# L = [f(TP[0], TP[1]) - f(TQ[0], TQ[1]) for f in forms] + b = 2*V(L) + + prof("lin alg") + M_sys = matrix(K, M_frob).transpose() - 1 + TP_to_TQ = M_sys**(-1) * b + +# print "TP to TQ: ", TP_to_TQ[0] +# print "\n" + prof("done") + print prof + return P_to_TP + TP_to_TQ + TQ_to_Q + + coleman_integrals_on_basis_hyperelliptic = coleman_integrals_on_basis + + + def monsky_washnitzer_gens(self): + import sage.schemes.elliptic_curves.monsky_washnitzer as monsky_washnitzer + S = monsky_washnitzer.SpecialHyperellipticQuotientRing(self) + return S.gens() + + def invariant_differential(self): + import sage.schemes.elliptic_curves.monsky_washnitzer as monsky_washnitzer + S = monsky_washnitzer.SpecialHyperellipticQuotientRing(self) + MW = monsky_washnitzer.MonskyWashnitzerDifferentialRing(S) + return MW.invariant_differential() + + def coleman_integral(self, w, P, Q): + import sage.schemes.elliptic_curves.monsky_washnitzer as monsky_washnitzer + S = monsky_washnitzer.SpecialHyperellipticQuotientRing(self, QQ) + MW = monsky_washnitzer.MonskyWashnitzerDifferentialRing(S) + print "w", w + w = MW(w) + f, rem = w.reduce_fast() + vec = rem.H1_vector() + basis_values = self.coleman_integrals_on_basis(P, Q) + dim = len(basis_values) + return f(P[0], P[1]) - f(Q[0], Q[1]) + sum([vec[i] * basis_values[i] for i in range(dim)]) # this is just a dot product... + + + +# TODO: add this to padics (if it isn't there in the new version already). +def padic_teichmuller(a): + K = a.parent() + p = K.prime() + p_less_1_inverse = ~K(p-1) + one = K(1) + x = K(ZZ(a) % p) + if x == 0: + return x + for _ in range(ceil(log(K.precision_cap(),2))): + x = ((p-2)*x + x**(2-p)) * p_less_1_inverse + return x + + + + + diff --git a/src/sage/schemes/hyperelliptic_curves/hyperelliptic_rational_field.py b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_rational_field.py index 7ea2b58c6ec..0fae2ff5783 100644 --- a/src/sage/schemes/hyperelliptic_curves/hyperelliptic_rational_field.py +++ b/src/sage/schemes/hyperelliptic_curves/hyperelliptic_rational_field.py @@ -9,6 +9,25 @@ #***************************************************************************** import hyperelliptic_generic +from sage.rings.all import is_pAdicField, is_pAdicRing, pAdicField class HyperellipticCurve_rational_field(hyperelliptic_generic.HyperellipticCurve_generic): - pass + + def matrix_of_frobenius(self, p, prec=20): + + # BUG: should get this method from HyperellipticCurve_generic + def my_chage_ring(self, R): + from constructor import HyperellipticCurve + f, h = self._hyperelliptic_polynomials + y = self._printing_ring.gen() + x = self._printing_ring.base_ring().gen() + return HyperellipticCurve(f.change_ring(R), h, "%s,%s"%(x,y)) + + import sage.schemes.elliptic_curves.monsky_washnitzer as monsky_washnitzer + if is_pAdicField(p) or is_pAdicRing(p): + K = p + else: + K = pAdicField(p, prec) + frob_p, forms = monsky_washnitzer.matrix_of_frobenius_hyperelliptic(my_chage_ring(self, K)) + return frob_p +