Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
r""" Computation of Frobenius matrix on Monsky-Washnitzer cohomology
The most interesting functions to be exported here are :func:`matrix_of_frobenius` and :func:`adjusted_prec`.
Currently this code is limited to the case `p \geq 5` (no `GF(p^n)` for `n > 1`), and only handles the elliptic curve case (not more general hyperelliptic curves).
REFERENCES:
.. [Ked2001] Kedlaya, K., "Counting points on hyperelliptic curves using Monsky-Washnitzer cohomology", J. Ramanujan Math. Soc. 16 (2001) no 4, 323-338
.. [Edix] Edixhoven, B., "Point counting after Kedlaya", EIDMA-Stieltjes graduate course, Lieden (lecture notes?).
AUTHORS:
- David Harvey and Robert Bradshaw: initial code developed at the 2006 MSRI graduate workshop, working with Jennifer Balakrishnan and Liang Xiao
- David Harvey (2006-08): cleaned up, rewrote some chunks, lots more documentation, added Newton iteration method, added more complete 'trace trick', integrated better into Sage.
- David Harvey (2007-02): added algorithm with sqrt(p) complexity (removed in May 2007 due to better C++ implementation)
- Robert Bradshaw (2007-03): keep track of exact form in reduction algorithms
- Robert Bradshaw (2007-04): generalization to hyperelliptic curves
- Julian Rueth (2014-05-09): improved caching
"""
#***************************************************************************** # Copyright (C) 2006 William Stein <wstein@gmail.com> # 2006 Robert Bradshaw <robertwb@math.washington.edu> # 2006 David Harvey <dmharvey@math.harvard.edu> # 2014 Julian Rueth <julian.rueth@fsfe.org> # # Distributed under the terms of the GNU General Public License (GPL) # http://www.gnu.org/licenses/ #***************************************************************************** from __future__ import print_function
from sage.rings.all import Integers, Integer, PolynomialRing, PowerSeriesRing, Rationals, Rational, LaurentSeriesRing
from sage.rings.polynomial.polynomial_element import is_Polynomial
from sage.modules.module import Module from sage.structure.element import ModuleElement from sage.matrix.all import matrix from sage.modules.all import vector from sage.rings.ring import CommutativeAlgebra from sage.structure.element import CommutativeAlgebraElement from sage.structure.unique_representation import UniqueRepresentation from sage.structure.richcmp import richcmp from sage.misc.cachefunc import cached_method from sage.rings.infinity import Infinity
from sage.arith.all import binomial, integer_ceil as ceil from sage.functions.log import log from sage.misc.misc import newton_method_sizes
from sage.schemes.elliptic_curves.ell_generic import is_EllipticCurve from sage.schemes.elliptic_curves.constructor import EllipticCurve
class SpecialCubicQuotientRing(CommutativeAlgebra): r""" Specialised class for representing the quotient ring `R[x,T]/(T - x^3 - ax - b)`, where `R` is an arbitrary commutative base ring (in which 2 and 3 are invertible), `a` and `b` are elements of that ring.
Polynomials are represented internally in the form `p_0 + p_1 x + p_2 x^2` where the `p_i` are polynomials in `T`. Multiplication of polynomials always reduces high powers of `x` (i.e. beyond `x^2`) to powers of `T`.
Hopefully this ring is faster than a general quotient ring because it uses the special structure of this ring to speed multiplication (which is the dominant operation in the frobenius matrix calculation). I haven't actually tested this theory though...
.. TODO::
Eventually we will want to run this in characteristic 3, so we need to: (a) Allow `Q(x)` to contain an `x^2` term, and (b) Remove the requirement that 3 be invertible. Currently this is used in the Toom-Cook algorithm to speed multiplication.
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: R SpecialCubicQuotientRing over Ring of integers modulo 125 with polynomial T = x^3 + 124*x + 94
Get generators::
sage: x, T = R.gens() sage: x (0) + (1)*x + (0)*x^2 sage: T (T) + (0)*x + (0)*x^2
Coercions::
sage: R(7) (7) + (0)*x + (0)*x^2
Create elements directly from polynomials::
sage: A = R.poly_ring() sage: A Univariate Polynomial Ring in T over Ring of integers modulo 125 sage: z = A.gen() sage: R.create_element(z^2, z+1, 3) (T^2) + (T + 1)*x + (3)*x^2
Some arithmetic::
sage: x^3 (T + 31) + (1)*x + (0)*x^2 sage: 3 * x**15 * T**2 + x - T (3*T^7 + 90*T^6 + 110*T^5 + 20*T^4 + 58*T^3 + 26*T^2 + 124*T) + (15*T^6 + 110*T^5 + 35*T^4 + 63*T^2 + 1)*x + (30*T^5 + 40*T^4 + 8*T^3 + 38*T^2)*x^2
Retrieve coefficients (output is zero-padded)::
sage: x^10 (3*T^2 + 61*T + 8) + (T^3 + 93*T^2 + 12*T + 40)*x + (3*T^2 + 61*T + 9)*x^2 sage: (x^10).coeffs() [[8, 61, 3, 0], [40, 12, 93, 1], [9, 61, 3, 0]]
.. TODO::
write an example checking multiplication of these polynomials against Sage's ordinary quotient ring arithmetic. I can't seem to get the quotient ring stuff happening right now... """ def __init__(self, Q, laurent_series=False): """ Constructor.
INPUT:
- ``Q`` -- a polynomial of the form `Q(x) = x^3 + ax + b`, where `a`, `b` belong to a ring in which 2, 3 are invertible.
- ``laurent_series`` -- whether or not to allow negative powers of `T` (default=False)
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: R SpecialCubicQuotientRing over Ring of integers modulo 125 with polynomial T = x^3 + 124*x + 94
::
sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 + 2*t^2 - t + B(1/4)) Traceback (most recent call last): ... ValueError: Q (=t^3 + 2*t^2 + 124*t + 94) must be of the form x^3 + ax + b
::
sage: B.<t> = PolynomialRing(Integers(10)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + 1) Traceback (most recent call last): ... ArithmeticError: 2 and 3 must be invertible in the coefficient ring (=Ring of integers modulo 10) of Q """ raise TypeError("Q (=%s) must be a polynomial" % Q)
"coefficient ring (=%s) of Q" % base_ring)
# CommutativeAlgebra.__init__ tries to establish a coercion # from the base ring, by trac ticket #9138. The corresponding # hom set is cached. In order to use self as cache key, its # string representation is used. In otder to get the string # representation, we need to know the attributes _a and # _b. Hence, in #9138, we have to move CommutativeAlgebra.__init__ # further down: self._poly_ring = LaurentSeriesRing(base_ring, 'T') # R[T] else:
# Precompute a matrix that is used in the Toom-Cook multiplication. # This is where we need 2 and 3 invertible.
# (a good description of Toom-Cook is online at: # http://www.gnu.org/software/gmp/manual/html_node/Toom-Cook-3-Way-Multiplication.html)
1, 1, 1, 8, 4, 2])**(-1)).change_ring(base_ring).list()
def __repr__(self): """ String representation
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: print(R) SpecialCubicQuotientRing over Ring of integers modulo 125 with polynomial T = x^3 + 124*x + 94 """ (self.base_ring(), PolynomialRing(self.base_ring(), 'x')( [self._b, self._a, 0, 1]))
def poly_ring(self): """ Return the underlying polynomial ring in `T`.
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: R.poly_ring() Univariate Polynomial Ring in T over Ring of integers modulo 125 """
def gens(self): """ Return a list [x, T] where x and T are the generators of the ring (as element *of this ring*).
.. note::
I have no idea if this is compatible with the usual Sage 'gens' interface.
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: x, T = R.gens() sage: x (0) + (1)*x + (0)*x^2 sage: T (T) + (0)*x + (0)*x^2 """ self._poly_ring(1), self._poly_ring(0), check=False), SpecialCubicQuotientRingElement(self, self._poly_generator, self._poly_ring(0), self._poly_ring(0), check=False)]
def create_element(self, p0, p1, p2, check=True): """ Creates the element `p_0 + p_1*x + p_2*x^2`, where the `p_i` are polynomials in `T`.
INPUT:
- ``p0, p1, p2`` -- coefficients; must be coercible into poly_ring()
- ``check`` -- bool (default True): whether to carry out coercion
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: A, z = R.poly_ring().objgen() sage: R.create_element(z^2, z+1, 3) (T^2) + (T + 1)*x + (3)*x^2 """
def __call__(self, value): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: R(3) (3) + (0)*x + (0)*x^2 """
def _coerce_impl(self, value): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: R._coerce_impl(3) (3) + (0)*x + (0)*x^2 """ # coerce to underlying polynomial ring (possibly via base ring):
self._poly_ring(0), check=False)
class SpecialCubicQuotientRingElement(CommutativeAlgebraElement): """ An element of a SpecialCubicQuotientRing. """ def __init__(self, parent, p0, p1, p2, check=True): """ Constructs the element `p_0 + p_1*x + p_2*x^2`, where the `p_i` are polynomials in `T`.
INPUT:
- ``parent`` -- a SpecialCubicQuotientRing
- ``p0, p1, p2`` -- coefficients; must be coercible into parent.poly_ring()
- ``check`` -- bool (default True): whether to carry out coercion
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import SpecialCubicQuotientRingElement sage: SpecialCubicQuotientRingElement(R, 2, 3, 4) (2) + (3)*x + (4)*x^2 """ raise TypeError("parent (=%s) must be a SpecialCubicQuotientRing" % parent)
def coeffs(self): """ Returns list of three lists of coefficients, corresponding to the `x^0`, `x^1`, `x^2` coefficients. The lists are zero padded to the same length. The list entries belong to the base ring.
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: p = R.create_element(t, t^2 - 2, 3) sage: p.coeffs() [[0, 1, 0], [123, 0, 1], [3, 0, 0]] """
def __bool__(self): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: x, T = R.gens() sage: not x False sage: not T False sage: not R.create_element(0, 0, 0) True """
__nonzero__ = __bool__
def _richcmp_(self, other, op): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: x, t = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)).gens() sage: x == t False sage: x == x True sage: x == x + x - x True """
def _repr_(self): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: x, T = R.gens() sage: x + T*x - 2*T^2 (123*T^2) + (T + 1)*x + (0)*x^2 """
def _latex_(self): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: x, T = R.gens() sage: f = x + T*x - 2*T^2 sage: latex(f) (123 T^{2}) + (T + 1)x + (0)x^2 """ tuple([column._latex_() for column in self._triple])
def _add_(self, other): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: f = R.create_element(2, t, t^2 - 3) sage: g = R.create_element(3 + t, -t, t) sage: f + g (T + 5) + (0)*x + (T^2 + T + 122)*x^2 """ self._triple[0] + other._triple[0], self._triple[1] + other._triple[1], self._triple[2] + other._triple[2], check=False)
def _sub_(self, other): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: f = R.create_element(2, t, t^2 - 3) sage: g = R.create_element(3 + t, -t, t) sage: f - g (124*T + 124) + (2*T)*x + (T^2 + 124*T + 122)*x^2 """ self._triple[0] - other._triple[0], self._triple[1] - other._triple[1], self._triple[2] - other._triple[2], check=False)
def shift(self, n): """ Returns this element multiplied by `T^n`.
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: f = R.create_element(2, t, t^2 - 3) sage: f (2) + (T)*x + (T^2 + 122)*x^2 sage: f.shift(1) (2*T) + (T^2)*x + (T^3 + 122*T)*x^2 sage: f.shift(2) (2*T^2) + (T^3)*x + (T^4 + 122*T^2)*x^2 """ self._triple[0].shift(n), self._triple[1].shift(n), self._triple[2].shift(n), check=False)
def scalar_multiply(self, scalar): """ Multiplies this element by a scalar, i.e. just multiply each coefficient of `x^j` by the scalar.
INPUT:
- ``scalar`` -- either an element of base_ring, or an element of poly_ring.
EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: x, T = R.gens() sage: f = R.create_element(2, t, t^2 - 3) sage: f (2) + (T)*x + (T^2 + 122)*x^2 sage: f.scalar_multiply(2) (4) + (2*T)*x + (2*T^2 + 119)*x^2 sage: f.scalar_multiply(t) (2*T) + (T^2)*x + (T^3 + 122*T)*x^2 """ scalar * self._triple[0], scalar * self._triple[1], scalar * self._triple[2], check=False)
def square(self): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: x, T = R.gens()
::
sage: f = R.create_element(1 + 2*t + 3*t^2, 4 + 7*t + 9*t^2, 3 + 5*t + 11*t^2) sage: f.square() (73*T^5 + 16*T^4 + 38*T^3 + 39*T^2 + 70*T + 120) + (121*T^5 + 113*T^4 + 73*T^3 + 8*T^2 + 51*T + 61)*x + (18*T^4 + 60*T^3 + 22*T^2 + 108*T + 31)*x^2 """
def _mul_(self, other): """ EXAMPLES::
sage: B.<t> = PolynomialRing(Integers(125)) sage: R = monsky_washnitzer.SpecialCubicQuotientRing(t^3 - t + B(1/4)) sage: x, T = R.gens()
::
sage: f = R.create_element(1 + 2*t + 3*t^2, 4 + 7*t + 9*t^2, 3 + 5*t + 11*t^2) sage: g = R.create_element(4 + 3*t + 7*t^2, 2 + 3*t + t^2, 8 + 4*t + 6*t^2) sage: f * g (65*T^5 + 27*T^4 + 33*T^3 + 75*T^2 + 120*T + 57) + (66*T^5 + T^4 + 123*T^3 + 95*T^2 + 24*T + 50)*x + (45*T^4 + 75*T^3 + 37*T^2 + 2*T + 52)*x^2 """ return self.scalar_multiply(other)
# Here we do Toom-Cook three-way multiplication, which reduces the # naive 9 polynomial multiplications to only 5 polynomial multiplications.
# faster method if we're squaring
else:
# Now the product is c0 + c1 x + c2 x^2 + c3 x^3 + c4 x^4. # We need to reduce mod y = x^3 + ax + b and return result.
# todo: These lines are necessary to get binop stuff working # for certain base rings, e.g. when we compute b*c3 in the # final line. They shouldn't be necessary. Need to fix this # somewhere else in Sage.
-b*c3 + c0 + c3*T, -b*c4 - a*c3 + c1 + c4*T, -a*c4 + c2, check=False)
def transpose_list(input): """ INPUT:
- ``input`` -- a list of lists, each list of the same length
OUTPUT:
- ``output`` -- a list of lists such that output[i][j] = input[j][i]
EXAMPLES::
sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import transpose_list sage: L = [[1, 2], [3, 4], [5, 6]] sage: transpose_list(L) [[1, 3, 5], [2, 4, 6]] """
def helper_matrix(Q): """ Computes the (constant) matrix used to calculate the linear combinations of the `d(x^i y^j)` needed to eliminate the negative powers of `y` in the cohomology (i.e. in reduce_negative()).
INPUT:
- ``Q`` -- cubic polynomial
EXAMPLES::
sage: t = polygen(QQ,'t') sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import helper_matrix sage: helper_matrix(t**3-4*t-691) [ 64/12891731 -16584/12891731 4297329/12891731] [ 6219/12891731 -32/12891731 8292/12891731] [ -24/12891731 6219/12891731 -32/12891731] """
# Discriminant (should be invertible for a curve of good reduction)
# This is the inverse of the matrix # [ a, -3b, 0 ] # [ 0, -2a, -3b ] # [ 3, 0, -2a ]
[-9*b, -2*a**2, 3*b*a], [6*a, -9*b, -2*a**2]])
def lift(x): r""" Tries to call x.lift(), presumably from the `p`-adics to ZZ.
If this fails, it assumes the input is a power series, and tries to lift it to a power series over QQ.
This function is just a very kludgy solution to the problem of trying to make the reduction code (below) work over both Zp and Zp[[t]].
EXAMPLES::
sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import lift sage: l = lift(Qp(13)(131)); l 131 sage: l.parent() Integer Ring
sage: x=PowerSeriesRing(Qp(17),'x').gen() sage: l = lift(4+5*x+17*x**6); l 4 + 5*t + 17*t^6 sage: l.parent() Power Series Ring in t over Rational Field """
def reduce_negative(Q, p, coeffs, offset, exact_form=None): """ Applies cohomology relations to incorporate negative powers of `y` into the `y^0` term.
INPUT:
- ``p`` -- prime
- ``Q`` -- cubic polynomial
- ``coeffs`` -- list of length 3 lists. The `i^{th}` list [a, b, c] represents `y^{2(i - offset)} (a + bx + cx^2) dx/y`.
- ``offset`` -- nonnegative integer
OUTPUT: The reduction is performed in-place. The output is placed in coeffs[offset]. Note that coeffs[i] will be meaningless for i offset after this function is finished.
EXAMPLES::
sage: R.<x> = Integers(5^3)['x'] sage: Q = x^3 - x + R(1/4) sage: coeffs = [[10, 15, 20], [1, 2, 3], [4, 5, 6], [7, 8, 9]] sage: coeffs = [[R.base_ring()(a) for a in row] for row in coeffs] sage: monsky_washnitzer.reduce_negative(Q, 5, coeffs, 3) sage: coeffs[3] [28, 52, 9]
::
sage: R.<x> = Integers(7^3)['x'] sage: Q = x^3 - x + R(1/4) sage: coeffs = [[7, 14, 21], [1, 2, 3], [4, 5, 6], [7, 8, 9]] sage: coeffs = [[R.base_ring()(a) for a in row] for row in coeffs] sage: monsky_washnitzer.reduce_negative(Q, 7, coeffs, 3) sage: coeffs[3] [245, 332, 9] """
x = exact_form.parent().gen(0) y = exact_form.parent()(exact_form.parent().base_ring().gen(0))
# todo: the following divisions will sometimes involve # a division by (a power of) p. In all cases, we know (from # Kedlaya's estimates) that the answer should be p-integral. # However, since we're working over $Z/p^k Z$, we're not allowed # to "divide by p". So currently we lift to Q, divide, and coerce # back. Eventually, when pAdicInteger is implemented, and plays # nicely with pAdicField, we should reimplement this stuff # using pAdicInteger.
# need to lift here to perform the division else:
c0 = m[0]*a[0] + m[1]*a[1] + m[2]*a[2] exact_form += (c0 + c1*x + c2 * x**2) * y**(j+1)
except NotImplementedError: raise NotImplementedError("It looks like you've found a " "non-integral matrix of Frobenius! " "(Q=%s, p=%s)\nTime to write a paper." % (Q, p))
def reduce_positive(Q, p, coeffs, offset, exact_form=None): """ Applies cohomology relations to incorporate positive powers of `y` into the `y^0` term.
INPUT:
- ``Q`` -- cubic polynomial
- ``coeffs`` -- list of length 3 lists. The `i^{th}` list [a, b, c] represents `y^{2(i - offset)} (a + bx + cx^2) dx/y`.
- ``offset`` -- nonnegative integer
OUTPUT: The reduction is performed in-place. The output is placed in coeffs[offset]. Note that coeffs[i] will be meaningless for i offset after this function is finished.
EXAMPLES::
sage: R.<x> = Integers(5^3)['x'] sage: Q = x^3 - x + R(1/4)
::
sage: coeffs = [[1, 2, 3], [10, 15, 20]] sage: coeffs = [[R.base_ring()(a) for a in row] for row in coeffs] sage: monsky_washnitzer.reduce_positive(Q, 5, coeffs, 0) sage: coeffs[0] [16, 102, 88]
::
sage: coeffs = [[9, 8, 7], [10, 15, 20]] sage: coeffs = [[R.base_ring()(a) for a in row] for row in coeffs] sage: monsky_washnitzer.reduce_positive(Q, 5, coeffs, 0) sage: coeffs[0] [24, 108, 92] """
x = exact_form.parent().gen(0) y = exact_form.parent().base_ring().gen(0) # y = exact_form.parent()(exact_form.parent().base_ring().gen(0))
exact_form += Q.base_ring()(a[2].lift() / (3*j+9)) * y**(j+3)
# todo: see comments about pAdicInteger in reduceNegative()
# subtract off c1 of d(x y^j + 1), and else:
# subtract off c2 of d(x^2 y^j + 1) else:
exact_form += (c1*x + c2 * x**2) * y**(j+1)
def reduce_zero(Q, coeffs, offset, exact_form=None): """ Applies cohomology relation to incorporate `x^2 y^0` term into `x^0 y^0` and `x^1 y^0` terms.
INPUT:
- ``Q`` -- cubic polynomial
- ``coeffs`` -- list of length 3 lists. The `i^{th}` list [a, b, c] represents `y^{2(i - offset)} (a + bx + cx^2) dx/y`.
- ``offset`` -- nonnegative integer
OUTPUT: The reduction is performed in-place. The output is placed in coeffs[offset]. This method completely ignores coeffs[i] for i != offset.
EXAMPLES::
sage: R.<x> = Integers(5^3)['x'] sage: Q = x^3 - x + R(1/4) sage: coeffs = [[1, 2, 3], [4, 5, 6], [7, 8, 9]] sage: coeffs = [[R.base_ring()(a) for a in row] for row in coeffs] sage: monsky_washnitzer.reduce_zero(Q, coeffs, 1) sage: coeffs[1] [6, 5, 0] """
y = exact_form.parent()(exact_form.parent().base_ring().gen(0)) exact_form += Q.base_ring()(a[2] / 3) * y
def reduce_all(Q, p, coeffs, offset, compute_exact_form=False): """ Applies cohomology relations to reduce all terms to a linear combination of `dx/y` and `x dx/y`.
INPUT:
- ``Q`` -- cubic polynomial
- ``coeffs`` -- list of length 3 lists. The `i^{th}` list [a, b, c] represents `y^{2(i - offset)} (a + bx + cx^2) dx/y`.
- ``offset`` -- nonnegative integer
OUTPUT:
- ``A, B`` - pair such that the input differential is cohomologous to (A + Bx) dx/y.
.. note::
The algorithm operates in-place, so the data in coeffs is destroyed.
EXAMPLES::
sage: R.<x> = Integers(5^3)['x'] sage: Q = x^3 - x + R(1/4) sage: coeffs = [[1, 2, 3], [4, 5, 6], [7, 8, 9]] sage: coeffs = [[R.base_ring()(a) for a in row] for row in coeffs] sage: monsky_washnitzer.reduce_all(Q, 5, coeffs, 1) (21, 106) """
# exact_form = SpecialCubicQuotientRing(Q, laurent_series=True)(0) exact_form = PolynomialRing(LaurentSeriesRing(Q.base_ring(), 'y'), 'x')(0) # t = (Q.base_ring().order().factor())[0] # from sage.rings.padics.qp import pAdicField # exact_form = PolynomialRing(LaurentSeriesRing(pAdicField(p, t[1]), 'y'), 'x')(0) else:
else: return (coeffs[int(offset)][0], coeffs[int(offset)][1]), exact_form
def frobenius_expansion_by_newton(Q, p, M): r""" Computes the action of Frobenius on `dx/y` and on `x dx/y`, using Newton's method (as suggested in Kedlaya's paper [Ked2001]_).
(This function does *not* yet use the cohomology relations - that happens afterwards in the "reduction" step.)
More specifically, it finds `F_0` and `F_1` in the quotient ring `R[x, T]/(T - Q(x))`, such that
.. MATH::
F( dx/y) = T^{-r} F0 dx/y, \text{\ and\ } F(x dx/y) = T^{-r} F1 dx/y
where
.. MATH::
r = ( (2M-3)p - 1 )/2.
(Here `T` is `y^2 = z^{-2}`, and `R` is the coefficient ring of `Q`.)
`F_0` and `F_1` are computed in the SpecialCubicQuotientRing associated to `Q`, so all powers of `x^j` for `j \geq 3` are reduced to powers of `T`.
INPUT:
- ``Q`` -- cubic polynomial of the form `Q(x) = x^3 + ax + b`, whose coefficient ring is a `Z/(p^M)Z`-algebra
- ``p`` -- residue characteristic of the p-adic field
- ``M`` -- p-adic precision of the coefficient ring (this will be used to determine the number of Newton iterations)
OUTPUT:
- ``F0, F1`` - elements of SpecialCubicQuotientRing(Q), as described above
- ``r`` - non-negative integer, as described above
EXAMPLES::
sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import frobenius_expansion_by_newton sage: R.<x> = Integers(5^3)['x'] sage: Q = x^3 - x + R(1/4) sage: frobenius_expansion_by_newton(Q,5,3) ((25*T^5 + 75*T^3 + 100*T^2 + 100*T + 100) + (5*T^6 + 80*T^5 + 100*T^3 + 25*T + 50)*x + (55*T^5 + 50*T^4 + 75*T^3 + 25*T^2 + 25*T + 25)*x^2, (5*T^8 + 15*T^7 + 95*T^6 + 10*T^5 + 25*T^4 + 25*T^3 + 100*T^2 + 50) + (65*T^7 + 55*T^6 + 70*T^5 + 100*T^4 + 25*T^2 + 100*T)*x + (15*T^6 + 115*T^5 + 75*T^4 + 100*T^3 + 50*T^2 + 75*T + 75)*x^2, 7) """
# When we compute Frob(1/y) we actually only need precision M-1, since # we're going to multiply by p at the end anyway.
# Kedlaya sets s = Q(x^p)/T^p = 1 + p T^{-p} E, where # E = (Q(x^p) - Q(x)^p) / p (has integral coefficients). # Then he computes s^{-1/2} in S, using Newton's method to find # successive approximations. We follow this plan, but we normalise our # approximations so that we only ever need positive powers of T.
# Start by setting r = Q(x^p)/2 = 1/2 T^p s. # (The 1/2 is for convenience later on.)
# todo: this next loop would be clearer if it used the newton_method_sizes() # function
# We will start with a hard-coded initial approximation, which we provide # up to precision 3. First work out what precision is best to start with. # In this case there is no advantage to starting with precision three, # because we'll overshoot at the end. E.g. suppose the final precision # is 8. If we start with precision 2, we need two iterations to get us # to 8. If we start at precision 3, we will still need two iterations, # but we do more work along the way. So may as well start with only 2. else:
# Now compute the first approximation. In the main loop below, X is the # normalised approximation, and k is the precision. More specifically, # X = T^{p(k-1)} x_i, where x_i is an approximation to s^{-1/2}, and the # approximation is correct mod p^k. # approximation is 3/2 - 1/2 s # approximation is (15 - 10 s + 3 s^2) / 8 - (base_ring(20) * r).shift(p) + (base_ring(12) * r.square())) # The key to the following calculation is that the T^{-m} coefficient # of every x_i is divisible by p^(ceil(m/p)) (for m >= 0). Therefore if # we are only expecting an answer correct mod p^k, we can truncate # beyond the T^{-(k-1)p} term without any problems.
# todo: what would be really nice is to be able to work in a lower # precision *coefficient ring* when we start the iteration, and move up to # higher precision rings as the iteration proceeds. This would be feasible # over Integers(p**n), but quite complicated (maybe impossible) over a more # general base ring. This might give a decent constant factor speedup; # or it might not, depending on how much the last iteration dominates the # whole runtime. My guess is that it isn't worth the effort.
# Newton iteration loop # target_k = k' = precision we want our answer to be after this iteration
# This prevents us overshooting. For example if the current precision # is 3 and we want to get to 10, we're better off going up to 5 # instead of 6, because it is less work to get from 5 to 10 than it # is to get from 6 to 10.
# temp = T^{p(3k-2)} 1/2 s x_i^3
# We know that the final result is only going to be correct mod # p^(target_k), so we might as well truncate the extraneous terms now. # temp = T^{p(k'-1)} 1/2 s x_i^3
# X = T^{p(k'-1)} (3/2 x_i - 1/2 s x_i^3) # = T^{p(k'-1)} x_{i+1}
# Now k should equal M, since we're up to the correct precision
# We should have s^{-1/2} correct to precision M. # The following line can be uncommented to verify this. # (It is a slow verification though, can double the whole computation time.)
#assert (p * X.square() * r * base_ring(2)).coeffs() == \ # R(p).shift(p*(2*M - 1)).coeffs()
# Finally incorporate frobenius of dx and x dx, and choose offset that # compensates for our normalisations by powers of T.
def frobenius_expansion_by_series(Q, p, M): r""" Computes the action of Frobenius on `dx/y` and on `x dx/y`, using a series expansion.
(This function computes the same thing as frobenius_expansion_by_newton(), using a different method. Theoretically the Newton method should be asymptotically faster, when the precision gets large. However, in practice, this functions seems to be marginally faster for moderate precision, so I'm keeping it here until I figure out exactly why it is faster.)
(This function does *not* yet use the cohomology relations - that happens afterwards in the "reduction" step.)
More specifically, it finds F0 and F1 in the quotient ring `R[x, T]/(T - Q(x))`, such that `F( dx/y) = T^{-r} F0 dx/y`, and `F(x dx/y) = T^{-r} F1 dx/y` where `r = ( (2M-3)p - 1 )/2`. (Here `T` is `y^2 = z^{-2}`, and `R` is the coefficient ring of `Q`.)
`F_0` and `F_1` are computed in the SpecialCubicQuotientRing associated to `Q`, so all powers of `x^j` for `j \geq 3` are reduced to powers of `T`.
It uses the sum
.. MATH::
F0 = \sum_{k=0}^{M-2} \binom{-1/2}{k} p x^{p-1} E^k T^{(M-2-k)p}
and
.. MATH::
F1 = x^p F0,
where `E = Q(x^p) - Q(x)^p`.
INPUT:
- ``Q`` -- cubic polynomial of the form `Q(x) = x^3 + ax + b`, whose coefficient ring is a `\ZZ/(p^M)\ZZ` -algebra
- ``p`` -- residue characteristic of the `p`-adic field
- ``M`` -- `p`-adic precision of the coefficient ring (this will be used to determine the number of terms in the series)
OUTPUT:
- ``F0, F1`` - elements of SpecialCubicQuotientRing(Q), as described above
- ``r`` - non-negative integer, as described above
EXAMPLES::
sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import frobenius_expansion_by_series sage: R.<x> = Integers(5^3)['x'] sage: Q = x^3 - x + R(1/4) sage: frobenius_expansion_by_series(Q,5,3) ((25*T^5 + 75*T^3 + 100*T^2 + 100*T + 100) + (5*T^6 + 80*T^5 + 100*T^3 + 25*T + 50)*x + (55*T^5 + 50*T^4 + 75*T^3 + 25*T^2 + 25*T + 25)*x^2, (5*T^8 + 15*T^7 + 95*T^6 + 10*T^5 + 25*T^4 + 25*T^3 + 100*T^2 + 50) + (65*T^7 + 55*T^6 + 70*T^5 + 100*T^4 + 25*T^2 + 100*T)*x + (15*T^6 + 115*T^5 + 75*T^4 + 100*T^3 + 50*T^2 + 75*T + 75)*x^2, 7) """
# compute frobQ = Q(x^p) # anticipating the day when p = 3 is supported: # frobQ = x_to_p_cubed + Q[2]*x_to_p_squared + Q[1]*x_to_p + Q[0]*S(1)
# todo: Possible speedup idea, perhaps by a factor of 2, but # it requires a lot of work: # Note that p divides E, so p^k divides E^k. So when we are # working with high powers of E, we're doing a lot more work # in the multiplications than we need to. To take advantage of # this we would need some protocol for "lowering the precision" # of a SpecialCubicQuotientRing. This would be quite messy to # do properly over an arbitrary base ring. Perhaps it is # feasible to do for the most common case (i.e. Z/p^nZ). # (but it probably won't save much time unless p^n is very # large, because the machine word size is probably pretty # big anyway.)
def adjusted_prec(p, prec): r""" Computes how much precision is required in matrix_of_frobenius to get an answer correct to prec `p`-adic digits.
The issue is that the algorithm used in :func:`matrix_of_frobenius` sometimes performs divisions by `p`, so precision is lost during the algorithm.
The estimate returned by this function is based on Kedlaya's result (Lemmas 2 and 3 of [Ked2001]_), which implies that if we start with `M` `p`-adic digits, the total precision loss is at most `1 + \lfloor \log_p(2M-3) \rfloor` `p`-adic digits. (This estimate is somewhat less than the amount you would expect by naively counting the number of divisions by `p`.)
INPUT:
- ``p`` -- a prime = 5
- ``prec`` -- integer, desired output precision, = 1
OUTPUT: adjusted precision (usually slightly more than prec) """
# initial estimate: else:
# increase it until we have enough
def matrix_of_frobenius(Q, p, M, trace=None, compute_exact_forms=False): """ Computes the matrix of Frobenius on Monsky-Washnitzer cohomology, with respect to the basis `(dx/y, x dx/y)`.
INPUT:
- ``Q`` -- cubic polynomial `Q(x) = x^3 + ax + b` defining an elliptic curve `E` by `y^2 = Q(x)`. The coefficient ring of `Q` should be a `\ZZ/(p^M)\ZZ`-algebra in which the matrix of frobenius will be constructed.
- ``p`` -- prime = 5 for which E has good reduction
- ``M`` -- integer = 2; `p` -adic precision of the coefficient ring
- ``trace`` -- (optional) the trace of the matrix, if known in advance. This is easy to compute because it is just the `a_p` of the curve. If the trace is supplied, matrix_of_frobenius will use it to speed the computation (i.e. we know the determinant is `p`, so we have two conditions, so really only column of the matrix needs to be computed. it is actually a little more complicated than that, but that's the basic idea.) If trace=None, then both columns will be computed independently, and you can get a strong indication of correctness by verifying the trace afterwards.
.. warning::
THE RESULT WILL NOT NECESSARILY BE CORRECT TO M p-ADIC DIGITS. If you want prec digits of precision, you need to use the function adjusted_prec(), and then you need to reduce the answer mod `p^{\mathrm{prec}}` at the end.
OUTPUT:
2x2 matrix of frobenius on Monsky-Washnitzer cohomology, with entries in the coefficient ring of Q.
EXAMPLES:
A simple example::
sage: p = 5 sage: prec = 3 sage: M = monsky_washnitzer.adjusted_prec(p, prec) sage: M 5 sage: R.<x> = PolynomialRing(Integers(p**M)) sage: A = monsky_washnitzer.matrix_of_frobenius(x^3 - x + R(1/4), p, M) sage: A [3090 187] [2945 408]
But the result is only accurate to prec digits::
sage: B = A.change_ring(Integers(p**prec)) sage: B [90 62] [70 33]
Check trace (123 = -2 mod 125) and determinant::
sage: B.det() 5 sage: B.trace() 123 sage: EllipticCurve([-1, 1/4]).ap(5) -2
Try using the trace to speed up the calculation::
sage: A = monsky_washnitzer.matrix_of_frobenius(x^3 - x + R(1/4), ....: p, M, -2) sage: A [2715 187] [1445 408]
Hmmm... it looks different, but that's because the trace of our first answer was only -2 modulo `5^3`, not -2 modulo `5^5`. So the right answer is::
sage: A.change_ring(Integers(p**prec)) [90 62] [70 33]
Check it works with only one digit of precision::
sage: p = 5 sage: prec = 1 sage: M = monsky_washnitzer.adjusted_prec(p, prec) sage: R.<x> = PolynomialRing(Integers(p**M)) sage: A = monsky_washnitzer.matrix_of_frobenius(x^3 - x + R(1/4), p, M) sage: A.change_ring(Integers(p)) [0 2] [0 3]
Here is an example that is particularly badly conditioned for using the trace trick::
sage: p = 11 sage: prec = 3 sage: M = monsky_washnitzer.adjusted_prec(p, prec) sage: R.<x> = PolynomialRing(Integers(p**M)) sage: A = monsky_washnitzer.matrix_of_frobenius(x^3 + 7*x + 8, p, M) sage: A.change_ring(Integers(p**prec)) [1144 176] [ 847 185]
The problem here is that the top-right entry is divisible by 11, and the bottom-left entry is divisible by `11^2`. So when you apply the trace trick, neither `F(dx/y)` nor `F(x dx/y)` is enough to compute the whole matrix to the desired precision, even if you try increasing the target precision by one. Nevertheless, ``matrix_of_frobenius`` knows how to get the right answer by evaluating `F((x+1) dx/y)` instead::
sage: A = monsky_washnitzer.matrix_of_frobenius(x^3 + 7*x + 8, p, M, -2) sage: A.change_ring(Integers(p**prec)) [1144 176] [ 847 185]
The running time is about ``O(p*prec**2)`` (times some logarithmic factors), so it is feasible to run on fairly large primes, or precision (or both?!?!)::
sage: p = 10007 sage: prec = 2 sage: M = monsky_washnitzer.adjusted_prec(p, prec) sage: R.<x> = PolynomialRing(Integers(p**M)) sage: A = monsky_washnitzer.matrix_of_frobenius( # long time ....: x^3 - x + R(1/4), p, M) # long time sage: B = A.change_ring(Integers(p**prec)); B # long time [74311982 57996908] [95877067 25828133] sage: B.det() # long time 10007 sage: B.trace() # long time 66 sage: EllipticCurve([-1, 1/4]).ap(10007) # long time 66
::
sage: p = 5 sage: prec = 300 sage: M = monsky_washnitzer.adjusted_prec(p, prec) sage: R.<x> = PolynomialRing(Integers(p**M)) sage: A = monsky_washnitzer.matrix_of_frobenius( # long time ....: x^3 - x + R(1/4), p, M) # long time sage: B = A.change_ring(Integers(p**prec)) # long time sage: B.det() # long time 5 sage: -B.trace() # long time 2 sage: EllipticCurve([-1, 1/4]).ap(5) # long time -2
Let us check consistency of the results for a range of precisions::
sage: p = 5 sage: max_prec = 60 sage: M = monsky_washnitzer.adjusted_prec(p, max_prec) sage: R.<x> = PolynomialRing(Integers(p**M)) sage: A = monsky_washnitzer.matrix_of_frobenius(x^3 - x + R(1/4), p, M) # long time sage: A = A.change_ring(Integers(p**max_prec)) # long time sage: result = [] # long time sage: for prec in range(1, max_prec): # long time ....: M = monsky_washnitzer.adjusted_prec(p, prec) # long time ....: R.<x> = PolynomialRing(Integers(p^M),'x') # long time ....: B = monsky_washnitzer.matrix_of_frobenius( # long time ....: x^3 - x + R(1/4), p, M) # long time ....: B = B.change_ring(Integers(p**prec)) # long time ....: result.append(B == A.change_ring( # long time ....: Integers(p**prec))) # long time sage: result == [True] * (max_prec - 1) # long time True
The remaining examples discuss what happens when you take the coefficient ring to be a power series ring; i.e. in effect you're looking at a family of curves.
The code does in fact work...
::
sage: p = 11 sage: prec = 3 sage: M = monsky_washnitzer.adjusted_prec(p, prec) sage: S.<t> = PowerSeriesRing(Integers(p**M), default_prec=4) sage: a = 7 + t + 3*t^2 sage: b = 8 - 6*t + 17*t^2 sage: R.<x> = PolynomialRing(S) sage: Q = x**3 + a*x + b sage: A = monsky_washnitzer.matrix_of_frobenius(Q, p, M) # long time sage: B = A.change_ring(PowerSeriesRing(Integers(p**prec), 't', default_prec=4)) # long time sage: B # long time [1144 + 264*t + 841*t^2 + 1025*t^3 + O(t^4) 176 + 1052*t + 216*t^2 + 523*t^3 + O(t^4)] [ 847 + 668*t + 81*t^2 + 424*t^3 + O(t^4) 185 + 341*t + 171*t^2 + 642*t^3 + O(t^4)]
The trace trick should work for power series rings too, even in the badly- conditioned case. Unfortunately I don't know how to compute the trace in advance, so I'm not sure exactly how this would help. Also, I suspect the running time will be dominated by the expansion, so the trace trick won't really speed things up anyway. Another problem is that the determinant is not always p::
sage: B.det() # long time 11 + 484*t^2 + 451*t^3 + O(t^4)
However, it appears that the determinant always has the property that if you substitute t - 11t, you do get the constant series p (mod p\*\*prec). Similarly for the trace. And since the parameter only really makes sense when it is divisible by p anyway, perhaps this isn't a problem after all. """
raise ValueError("M (=%s) must be at least 2" % M)
# Expand out frobenius of dx/y and x dx/y. # (You can substitute frobenius_expansion_by_series here, that will work # as well. See its docstring for some performance notes.) #F0, F1, offset = frobenius_expansion_by_series(Q, p, M)
# we need to do all the work to get the exact expressions f such that F(x^i dx/y) = df + \sum a_i x^i dx/y F0_coeffs = transpose_list(F0.coeffs()) F0_reduced, f_0 = reduce_all(Q, p, F0_coeffs, offset, True)
F1_coeffs = transpose_list(F1.coeffs()) F1_reduced, f_1 = reduce_all(Q, p, F1_coeffs, offset, True)
# This implies that only one digit of precision is valid, so we only need # to reduce the second column. Also, the trace doesn't help at all.
# No trace provided, just reduce F(dx/y) and F(x dx/y) separately.
else: # Trace has been provided.
# In most cases this can be used to quickly compute F(dx/y) from # F(x dx/y). However, if we're unlucky, the (dx/y)-component of # F(x dx/y) (i.e. the top-right corner of the matrix) may be divisible # by p, in which case there isn't enough information to get the # (x dx/y)-component of F(dx/y) to the desired precision. When this # happens, it turns out that F((x+1) dx/y) always *does* give enough # information (together with the trace) to get both columns to the # desired precision.
# First however we need a quick way of telling whether the top-right # corner is divisible by p, i.e. we want to compute the second column # of the matrix mod p. We could do this by just running the entire # algorithm with M = 2 (which assures precision 1). Luckily, we've # already done most of the work by computing F1 to high precision; so # all we need to do is extract the coefficients that would correspond # to the first term of the series, and run the reduction on them.
# todo: actually we only need to do this reduction step mod p^2, not # mod p^M, which is what the code currently does. If the base ring # is Integers(p^M), then it is easy. Otherwise it is tricky to construct # the right ring, I don't know how to do it.
# make a copy, because reduce_all will destroy the coefficients:
# If the first entry is invertible mod p, then F(x dx/y) is sufficient # to get the whole matrix.
# using that the determinant is p: / F1_reduced[0]
else: # If the first entry is zero mod p, then F((x+1) dx/y) will be sufficient # to get the whole matrix. (Here we are using the fact that the second # entry *cannot* be zero mod p. This is guaranteed by some results in # section 3.2 of ``Computation of p-adic Heights and Log Convergence'' # by Mazur, Stein, Tate. But let's quickly check it anyway :-))
# Now G0_reduced expresses F((x+1) dx/y) in terms of dx/y and x dx/y. # Re-express this in terms of (x+1) dx/y and x dx/y.
# The thing we're about to divide by better be a unit.
# Figure out the second column using the trace... # ... and using that the determinant is p: / H0_reduced[1]
# Finally, change back to the usual basis (dx/y, x dx/y) H1_reduced[0] + H1_reduced[1]] H0_reduced[0] + H0_reduced[1] - F1_reduced[1]]
# One more sanity check: our final result should be congruent mod p # to the approximation we used earlier. (F1_reduced[0] - F1_modp_reduced[0]).is_unit() or (F1_reduced[1] - F1_modp_reduced[1]).is_unit() or F0_reduced[0].is_unit() or F0_reduced[1].is_unit()), msg
return matrix(base_ring, 2, 2, [F0_reduced[0], F1_reduced[0], F0_reduced[1], F1_reduced[1]]), f_0, f_1 else: F0_reduced[1], F1_reduced[1]])
#***************************************************************************** # This is a generalization of the above functionality for hyperelliptic curves. # # THIS IS A WORK IN PROGRESS. # # I tried to embed must stuff into the rings themselves rather than # just extract and manipulate lists of coefficients. Hence the implementations # below are much less optimized, so are much slower, but should hopefully be # easier to follow. (E.g. one can print/make sense of intermediate results.) # # AUTHOR: # -- Robert Bradshaw (2007-04) # #*****************************************************************************
import weakref
from sage.schemes.hyperelliptic_curves.constructor import HyperellipticCurve from sage.schemes.hyperelliptic_curves.hyperelliptic_generic import is_HyperellipticCurve from sage.rings.padics.all import pAdicField from sage.rings.all import QQ, IntegralDomain
from sage.rings.laurent_series_ring import is_LaurentSeriesRing
from sage.modules.free_module import FreeModule from sage.modules.free_module_element import is_FreeModuleElement
from sage.misc.profiler import Profiler from sage.misc.misc import repr_lincomb
def matrix_of_frobenius_hyperelliptic(Q, p=None, prec=None, M=None): """ Computes the matrix of Frobenius on Monsky-Washnitzer cohomology, with respect to the basis `(dx/2y, x dx/2y, ...x^{d-2} dx/2y)`, where `d` is the degree of `Q`.
INPUT:
- ``Q`` -- monic polynomial `Q(x)`
- ``p`` -- prime `\geq 5` for which `E` has good reduction
- ``prec`` -- (optional) `p`-adic precision of the coefficient ring
- ``M`` -- (optional) adjusted `p`-adic precision of the coefficient ring
OUTPUT:
`(d-1)` x `(d-1)` matrix `M` of Frobenius on Monsky-Washnitzer cohomology, and list of differentials \{f_i \} such that
.. MATH::
\phi^* (x^i dx/2y) = df_i + M[i]*vec(dx/2y, ..., x^{d-2} dx/2y)
EXAMPLES::
sage: p = 5 sage: prec = 3 sage: R.<x> = QQ['x'] sage: A,f = monsky_washnitzer.matrix_of_frobenius_hyperelliptic(x^5 - 2*x + 3, p, prec) sage: A [ 4*5 + O(5^3) 5 + 2*5^2 + O(5^3) 2 + 3*5 + 2*5^2 + O(5^3) 2 + 5 + 5^2 + O(5^3)] [ 3*5 + 5^2 + O(5^3) 3*5 + O(5^3) 4*5 + O(5^3) 2 + 5^2 + O(5^3)] [ 4*5 + 4*5^2 + O(5^3) 3*5 + 2*5^2 + O(5^3) 5 + 3*5^2 + O(5^3) 2*5 + 2*5^2 + O(5^3)] [ 5^2 + O(5^3) 5 + 4*5^2 + O(5^3) 4*5 + 3*5^2 + O(5^3) 2*5 + O(5^3)]
""" except AttributeError: raise ValueError("p and prec must be specified if Q is not " "defined over a p-adic ring") # extra_prec_ring = pAdicField(p, M) # SLOW!
# do reduction over Q in case we have non-integral entries (and it is so much faster than padics) # 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 coefficient size) # S._p = p # rational_S(F[0]).reduce_fast() # prof("reduce others")
# rational_S = S.change_ring(pAdicField(p, M))
# reduced = [F_i.reduce() for F_i in F]
# but the coeffs are WAY more precision than they need to be
# now take care of precision capping
class SpecialHyperellipticQuotientRing(UniqueRepresentation, CommutativeAlgebra): _p = None
def __init__(self, Q, R=None, invert_y=True): r""" Initialization.
TESTS:
Check that caching works::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import SpecialHyperellipticQuotientRing sage: SpecialHyperellipticQuotientRing(E) is SpecialHyperellipticQuotientRing(E) True
"""
# Trac ticket #9138: CommutativeAlgebra.__init__ must not be # done so early. It tries to register a coercion, but that # requires the hash being available. But the hash, in its # default implementation, relies on the string representation, # which is not available at this point. #CommutativeAlgebra.__init__(self, R) # moved to below.
raise NotImplementedError("Curve must be in Weierstrass " "normal form.")
raise NotImplementedError("Curve must be of form y^2 = Q(x).")
raise NotImplementedError("Polynomial must be monic.") else:
else: raise NotImplementedError("Must be an elliptic curve or polynomial " "Q for y^2 = Q(x)\n(Got element of %s)" % Q.parent())
# Trac ticket #9138: Initialise the commutative algebra here! # Below, we do self(self._poly_ring.gen(0)), which requires # the initialisation being finished.
def _repr_(self): """ String representation
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent() # indirect doctest SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = x^5 - 3*x + 1) over Rational Field """
def base_extend(self, R): """ Return the base extension of ``self`` to the ring ``R`` if possible.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().base_extend(UniversalCyclotomicField()) SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = x^5 - 3*x + 1) over Universal Cyclotomic Field sage: x.parent().base_extend(ZZ) Traceback (most recent call last): ... TypeError: no such base extension """ else:
def change_ring(self, R): """ Return the analog of ``self`` over the ring ``R``
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().change_ring(ZZ) SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = x^5 - 3*x + 1) over Integer Ring """
def __call__(self, val, offset=0, check=True): else:
def gens(self): """ Return the generators of ``self``
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().gens() (x, y*1) """
def x(self): r""" Return the generator `x` of ``self``
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().x() x """
def y(self): r""" Return the generator `y` of ``self``
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().y() y*1 """
def monomial(self, i, j, b=None): """ Returns `b y^j x^i`, computed quickly.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().monomial(4,5) y^5*x^4 """
else: else:
def monomial_diff_coeffs(self, i, j): r""" The key here is that the formula for `d(x^iy^j)` is messy in terms of `i`, but varies nicely with `j`.
.. MATH::
d(x^iy^j) = y^{j-1} (2ix^{i-1}y^2 + j (A_i(x) + B_i(x)y^2)) \frac{dx}{2y}
Where `A,B` have degree at most `n-1` for each `i`. Pre-compute `A_i, B_i` for each `i` the "hard" way, and the rest are easy. """ else: else: dg = self.monomial(i, j).diff() coeffs = [dg.extract_pow_y(j-1), dg.extract_pow_y(j+1)] self._monomial_diff_coeffs[i, j] = coeffs return coeffs
def monomial_diff_coeffs_matrices(self):
def _precompute_monomial_diffs(self):
def Q(self): """ Return the defining polynomial of the underlying hyperelliptic curve.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-2*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().Q() x^5 - 2*x + 1 """
def curve(self): """ Return the underlying hyperelliptic curve.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().curve() Hyperelliptic Curve over Rational Field defined by y^2 = x^5 - 3*x + 1 """
def degree(self): """ Return the degree of the underlying hyperelliptic curve.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().degree() 5 """
def prime(self): return self._p
def monsky_washnitzer(self):
def is_field(self, proof=True): """ Return False as ``self`` is not a field.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.parent().is_field() False """ SpecialHyperellipticQuotientRing_class = SpecialHyperellipticQuotientRing
class SpecialHyperellipticQuotientElement(CommutativeAlgebraElement):
def __init__(self, parent, val=0, offset=0, check=True): """ Elements in the Hyperelliptic quotient ring
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: MW = x.parent() sage: MW(x+x**2+y-77) # indirect doctest -(77-y)*1 + x + x^2 """ val, offset = val
def _richcmp_(self, other, op): """ Compare the elements.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x == x True sage: x > y True """
def change_ring(self, R): """ Return the same element after changing the base ring to R.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: MW = x.parent() sage: z = MW(x+x**2+y-77) sage: z.change_ring(AA).parent() SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = x^5 - 36*x + 1) over Algebraic Real Field """
def __call__(self, *x): """ Evaluate ``self`` at given arguments
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: MW = x.parent() sage: z = MW(x+x**2+y-77); z -(77-y)*1 + x + x^2 sage: z(66) 4345 + y sage: z(5,4) -43 """
def __invert__(self): """ Return the inverse of the element
The general element in our ring is not invertible, but `y` may be. We do not want to pass to the fraction field.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: MW = x.parent() sage: z = y**(-1) # indirect doctest sage: z.parent() SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = x^5 - 36*x + 1) over Rational Field
sage: z = (x+y)**(-1) # indirect doctest Traceback (most recent call last): ... ZeroDivisionError: Element not invertible """ else:
def __bool__(self): """ Return True iff ``self`` is not zero.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: bool(x) True """
__nonzero__ = __bool__
def __eq__(self, other): """ Return True iff ``self`` is equal to ``other``
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x == y # indirect doctest False """
def _add_(self, other): """ Return the sum of two elements
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x+y y*1 + x """
def _sub_(self, other): """ Return the difference of two elements
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: y-x y*1 - x """
def _mul_(self, other): """ Return the product of two elements
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-36*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: y*x y*x """ # over Laurent series, addition and subtraction can be # expensive, and the degree of this poly is small enough that # Karatsuba actually hurts significantly in some cases else:
def _rmul_(self, c):
def _lmul_(self, c):
def __lshift__(self, k):
def __rshift__(self, k):
def truncate_neg(self, n): """ Return ``self`` minus its terms of degree less than `n` wrt `y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: (x+3*y+7*x*2*y**4).truncate_neg(1) 3*y*1 + 14*y^4*x """
def _repr_(self): """ Return a string representation of ``self``.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: (x+3*y)._repr_() '3*y*1 + x' """
def _latex_(self): """ Return a LateX string for ``self``.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: (x+3*y)._latex_() '3y1 + x' """
def diff(self): """ Return the differential of ``self``
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: (x+3*y).diff() (-(9-2*y)*1 + 15*x^4) dx/2y """ # try: # return self._diff_x # except AttributeError: # pass
# d(self) = A dx + B dy # = (2y A + BQ') dx/2y # self._diff = self.parent()._monsky_washnitzer(two_y * A + dQ * B) # return self._diff
def extract_pow_y(self, k): r""" Return the coefficients of `y^k` in ``self`` as a list
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: (x+3*y+9*x*y).extract_pow_y(1) [3, 9, 0, 0, 0] """
def min_pow_y(self): """ Return the minimal degree of ``self`` w.r.t. y
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: (x+3*y).min_pow_y() 0 """
def max_pow_y(self): """ Return the maximal degree of ``self`` w.r.t. y
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: (x+3*y).max_pow_y() 1 """
def coeffs(self, R=None): """ Returns the raw coefficients of this element.
INPUT:
- ``R`` -- an (optional) base-ring in which to cast the coefficients
OUTPUT:
- ``coeffs`` -- a list of coefficients of powers of `x` for each power of `y`
- ``n`` -- an offset indicating the power of `y` of the first list element
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x,y = E.monsky_washnitzer_gens() sage: x.coeffs() ([(0, 1, 0, 0, 0)], 0) sage: y.coeffs() ([(0, 0, 0, 0, 0), (1, 0, 0, 0, 0)], 0)
sage: a = sum(n*x^n for n in range(5)); a x + 2*x^2 + 3*x^3 + 4*x^4 sage: a.coeffs() ([(0, 1, 2, 3, 4)], 0) sage: a.coeffs(Qp(7)) ([(0, 1 + O(7^20), 2 + O(7^20), 3 + O(7^20), 4 + O(7^20))], 0) sage: (a*y).coeffs() ([(0, 0, 0, 0, 0), (0, 1, 2, 3, 4)], 0) sage: (a*y^-2).coeffs() ([(0, 1, 2, 3, 4), (0, 0, 0, 0, 0), (0, 0, 0, 0, 0)], -2)
Note that the coefficient list is transposed compared to how they are stored and printed::
sage: a*y^-2 (y^-2)*x + (2*y^-2)*x^2 + (3*y^-2)*x^3 + (4*y^-2)*x^4
A more complicated example::
sage: a = x^20*y^-3 - x^11*y^2; a (y^-3-4*y^-1+6*y-4*y^3+y^5)*1 - (12*y^-3-36*y^-1+36*y+y^2-12*y^3-2*y^4+y^6)*x + (54*y^-3-108*y^-1+54*y+6*y^2-6*y^4)*x^2 - (108*y^-3-108*y^-1+9*y^2)*x^3 + (81*y^-3)*x^4 sage: raw, offset = a.coeffs() sage: a.min_pow_y() -3 sage: offset -3 sage: raw [(1, -12, 54, -108, 81), (0, 0, 0, 0, 0), (-4, 36, -108, 108, 0), (0, 0, 0, 0, 0), (6, -36, 54, 0, 0), (0, -1, 6, -9, 0), (-4, 12, 0, 0, 0), (0, 2, -6, 0, 0), (1, 0, 0, 0, 0), (0, -1, 0, 0, 0)] sage: sum(c * x^i * y^(j+offset) for j, L in enumerate(raw) for i, c in enumerate(L)) == a True
Can also be used to construct elements::
sage: a.parent()(raw, offset) == a True """
class MonskyWashnitzerDifferentialRing(UniqueRepresentation, Module): r""" A ring of Monsky--Washnitzer differentials over ``base_ring``. """ def __init__(self, base_ring): r""" Initialization.
TESTS:
Check that caching works::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: from sage.schemes.hyperelliptic_curves.monsky_washnitzer import SpecialHyperellipticQuotientRing, MonskyWashnitzerDifferentialRing sage: S = SpecialHyperellipticQuotientRing(E) sage: MonskyWashnitzerDifferentialRing(S) is MonskyWashnitzerDifferentialRing(S) True
"""
def invariant_differential(self): """ Returns `dx/2y` as an element of self.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.invariant_differential() 1 dx/2y """
def __call__(self, val, offset=0):
def base_extend(self, R): """ Return a new differential ring which is self base-extended to `R`
INPUT:
- ``R`` -- ring
OUTPUT:
Self, base-extended to `R`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.base_ring() SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = x^5 - 4*x + 4) over Rational Field sage: MW.base_extend(Qp(5,5)).base_ring() SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = (1 + O(5^5))*x^5 + (1 + 4*5 + 4*5^2 + 4*5^3 + 4*5^4 + O(5^5))*x + (4 + O(5^5))) over 5-adic Field with capped relative precision 5 """
def change_ring(self, R): """ Returns a new differential ring which is self with the coefficient ring changed to `R`.
INPUT:
- ``R`` -- ring of coefficients
OUTPUT:
Self, with the coefficient ring changed to `R`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.base_ring() SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = x^5 - 4*x + 4) over Rational Field sage: MW.change_ring(Qp(5,5)).base_ring() SpecialHyperellipticQuotientRing K[x,y,y^-1] / (y^2 = (1 + O(5^5))*x^5 + (1 + 4*5 + 4*5^2 + 4*5^3 + 4*5^4 + O(5^5))*x + (4 + O(5^5))) over 5-adic Field with capped relative precision 5 """
def degree(self): """ Returns the degree of `Q(x)`, where the model of the underlying hyperelliptic curve of self is given by `y^2 = Q(x)`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.Q() x^5 - 4*x + 4 sage: MW.degree() 5 """
def dimension(self): """ Returns the dimension of self.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: K = Qp(7,5) sage: CK = C.change_ring(K) sage: MW = CK.invariant_differential().parent() sage: MW.dimension() 4 """
def Q(self): """ Returns `Q(x)` where the model of the underlying hyperelliptic curve of self is given by `y^2 = Q(x)`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.Q() x^5 - 4*x + 4 """
@cached_method def x_to_p(self, p): """ Returns and caches `x^p`, reduced via the relations coming from the defining polynomial of the hyperelliptic curve.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.x_to_p(3) x^3 sage: MW.x_to_p(5) -(4-y^2)*1 + 4*x sage: MW.x_to_p(101) is MW.x_to_p(101) True """
@cached_method def frob_Q(self, p): """ Returns and caches `Q(x^p)`, which is used in computing the image of `y` under a `p`-power lift of Frobenius to `A^{\dagger}`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.frob_Q(3) -(60-48*y^2+12*y^4-y^6)*1 + (192-96*y^2+12*y^4)*x - (192-48*y^2)*x^2 + 60*x^3 sage: MW.Q()(MW.x_to_p(3)) -(60-48*y^2+12*y^4-y^6)*1 + (192-96*y^2+12*y^4)*x - (192-48*y^2)*x^2 + 60*x^3 sage: MW.frob_Q(11) is MW.frob_Q(11) True """
def frob_invariant_differential(self, prec, p): r""" Kedlaya's algorithm allows us to calculate the action of Frobenius on the Monsky-Washnitzer cohomology. First we lift `\phi` to `A^{\dagger}` by setting
.. MATH::
\phi(x) = x^p
\phi(y) = y^p \sqrt{1 + \frac{Q(x^p) - Q(x)^p}{Q(x)^p}}.
Pulling back the differential `dx/2y`, we get
.. MATH::
\phi^*(dx/2y) = px^{p-1} y(\phi(y))^{-1} dx/2y = px^{p-1} y^{1-p} \sqrt{1+ \frac{Q(x^p) - Q(x)^p}{Q(x)^p}} dx/2y
Use Newton's method to calculate the square root.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: prec = 2 sage: p = 7 sage: MW = C.invariant_differential().parent() sage: MW.frob_invariant_differential(prec,p) ((67894400*y^-20-81198880*y^-18+40140800*y^-16-10035200*y^-14+1254400*y^-12-62720*y^-10)*1 - (119503944*y^-20-116064242*y^-18+43753472*y^-16-7426048*y^-14+514304*y^-12-12544*y^-10+1568*y^-8-70*y^-6-7*y^-4)*x + (78905288*y^-20-61014016*y^-18+16859136*y^-16-2207744*y^-14+250880*y^-12-37632*y^-10+3136*y^-8-70*y^-6)*x^2 - (39452448*y^-20-26148752*y^-18+8085490*y^-16-2007040*y^-14+376320*y^-12-37632*y^-10+1568*y^-8)*x^3 + (21102144*y^-20-18120592*y^-18+8028160*y^-16-2007040*y^-14+250880*y^-12-12544*y^-10)*x^4) dx/2y """ # TODO, would it be useful to be able to take Frobenius of any element? Less efficient?
# cache for future use
# Q = self.base_ring()._Q # three_halves = Q.parent().base_ring()(Rational((3,2))) # one_half = Q.parent().base_ring()(Rational((1,2)))
# We are solving for t = a^{-1/2} = (F_pQ y^{-p})^{-1/2} # Newton's method converges because we know the root is in the same residue class as 1.
# t = self.base_ring()(1) # first iteration trivial, start with prec 2
# newton_method_sizes = [1, 2, ...] # binomial expansion is $\sum p^{k+1} y^{-(2k+1)p+1} f(x)$ # so if we are only correct mod p^prec, # can ignore y powers less than y_prec # t = (3/2) t - (1/2) a t^3
def frob_basis_elements(self, prec, p): """ Returns the action of a `p`-power lift of Frobenius on the basis
.. MATH::
\{ dx/2y, x dx/2y, ..., x^{d-2} dx/2y \}
where `d` is the degree of the underlying hyperelliptic curve.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: prec = 1 sage: p = 5 sage: MW = C.invariant_differential().parent() sage: MW.frob_basis_elements(prec,p) [((92000*y^-14-74200*y^-12+32000*y^-10-8000*y^-8+1000*y^-6-50*y^-4)*1 - (194400*y^-14-153600*y^-12+57600*y^-10-9600*y^-8+600*y^-6)*x + (204800*y^-14-153600*y^-12+38400*y^-10-3200*y^-8)*x^2 - (153600*y^-14-76800*y^-12+9600*y^-10)*x^3 + (63950*y^-14-18550*y^-12+1600*y^-10-400*y^-8+50*y^-6+5*y^-4)*x^4) dx/2y, (-(1391200*y^-14-941400*y^-12+302000*y^-10-76800*y^-8+14400*y^-6-1320*y^-4+30*y^-2)*1 + (2168800*y^-14-1402400*y^-12+537600*y^-10-134400*y^-8+16800*y^-6-720*y^-4)*x - (1596800*y^-14-1433600*y^-12+537600*y^-10-89600*y^-8+5600*y^-6)*x^2 + (1433600*y^-14-1075200*y^-12+268800*y^-10-22400*y^-8)*x^3 - (870200*y^-14-445350*y^-12+63350*y^-10-3200*y^-8+600*y^-6-30*y^-4-5*y^-2)*x^4) dx/2y, ((19488000*y^-14-15763200*y^-12+4944400*y^-10-913800*y^-8+156800*y^-6-22560*y^-4+1480*y^-2-10)*1 - (28163200*y^-14-18669600*y^-12+5774400*y^-10-1433600*y^-8+268800*y^-6-25440*y^-4+760*y^-2)*x + (15062400*y^-14-12940800*y^-12+5734400*y^-10-1433600*y^-8+179200*y^-6-8480*y^-4)*x^2 - (12121600*y^-14-11468800*y^-12+4300800*y^-10-716800*y^-8+44800*y^-6)*x^3 + (9215200*y^-14-6952400*y^-12+1773950*y^-10-165750*y^-8+5600*y^-6-720*y^-4+10*y^-2+5)*x^4) dx/2y, (-(225395200*y^-14-230640000*y^-12+91733600*y^-10-18347400*y^-8+2293600*y^-6-280960*y^-4+31520*y^-2-1480-10*y^2)*1 + (338048000*y^-14-277132800*y^-12+89928000*y^-10-17816000*y^-8+3225600*y^-6-472320*y^-4+34560*y^-2-720)*x - (172902400*y^-14-141504000*y^-12+58976000*y^-10-17203200*y^-8+3225600*y^-6-314880*y^-4+11520*y^-2)*x^2 + (108736000*y^-14-109760000*y^-12+51609600*y^-10-12902400*y^-8+1612800*y^-6-78720*y^-4)*x^3 - (85347200*y^-14-82900000*y^-12+31251400*y^-10-5304150*y^-8+367350*y^-6-8480*y^-4+760*y^-2+10-5*y^2)*x^4) dx/2y] """
def helper_matrix(self): """ We use this to solve for the linear combination of `x^i y^j` needed to clear all terms with `y^{j-1}`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: MW = C.invariant_differential().parent() sage: MW.helper_matrix() [ 256/2101 320/2101 400/2101 500/2101 625/2101] [-625/8404 -64/2101 -80/2101 -100/2101 -125/2101] [-125/2101 -625/8404 -64/2101 -80/2101 -100/2101] [-100/2101 -125/2101 -625/8404 -64/2101 -80/2101] [ -80/2101 -100/2101 -125/2101 -625/8404 -64/2101] """
# The smallest y term of (1/j) d(x^i y^j) is constant for all j. # must be using integer_mod or something to approximate self._helper_matrix = (~A.change_ring(QQ)).change_ring(A.base_ring()) else: MonskyWashnitzerDifferentialRing_class = MonskyWashnitzerDifferentialRing
class MonskyWashnitzerDifferential(ModuleElement):
def __init__(self, parent, val=0, offset=0): r""" Create an element of the Monsky-Washnitzer ring of differentials, of the form `F dx/2y`.
INPUT:
- ``parent`` -- Monsky-Washnitzer differential ring (instance of class :class:`~MonskyWashnitzerDifferentialRing`
- ``val`` -- element of the base ring, or list of coefficients
- ``offset`` -- if non-zero, shift val by `y^\text{offset}` (default 0)
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5 - 4*x + 4) sage: x,y = C.monsky_washnitzer_gens() sage: MW = C.invariant_differential().parent() sage: sage.schemes.hyperelliptic_curves.monsky_washnitzer.MonskyWashnitzerDifferential(MW, x) x dx/2y sage: sage.schemes.hyperelliptic_curves.monsky_washnitzer.MonskyWashnitzerDifferential(MW, y) y*1 dx/2y sage: sage.schemes.hyperelliptic_curves.monsky_washnitzer.MonskyWashnitzerDifferential(MW, x, 10) y^10*x dx/2y """
def _add_(left, right): """ Returns the sum of left and right, both elements of the Monsky-Washnitzer ring of differentials.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x + 4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: w + w 2*1 dx/2y sage: x*w + w (1 + x) dx/2y sage: x*w + y*w (y*1 + x) dx/2y """ left._coeff + right._coeff)
def _sub_(left, right): """ Returns the difference of left and right, both elements of the Monsky-Washnitzer ring of differentials.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: w-w 0 dx/2y sage: x*w-w (-1 + x) dx/2y sage: w - x*w - y*w ((1-y)*1 - x) dx/2y """ left._coeff - right._coeff)
def __neg__(self): """ Returns the additive inverse of self.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: -w -1 dx/2y sage: -((y-x)*w) (-y*1 + x) dx/2y """
def _lmul_(self, a): """ Returns `self * a`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: w*x x dx/2y sage: (w*x)*2 2*x dx/2y sage: w*y y*1 dx/2y sage: w*(x+y) (y*1 + x) dx/2y """
def _rmul_(self, a): """ Returns `a * self`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: x*w x dx/2y sage: 2*(x*w) 2*x dx/2y sage: y*w y*1 dx/2y sage: (x+y)*w (y*1 + x) dx/2y """
def coeff(self): r""" Returns `A`, where this element is `A dx/2y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: w 1 dx/2y sage: w.coeff() 1 sage: (x*y*w).coeff() y*x """
def __bool__(self): """ EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: not w False sage: not 0*w True sage: not x*y*w False """
__nonzero__ = __bool__
def _repr_(self): """ EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: w 1 dx/2y sage: (2*x+y)*w (y*1 + 2*x) dx/2y """
def _latex_(self): """ Returns the latex representation of self.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: latex(w) 1 \frac{dx}{2y} sage: latex(x*w) x \frac{dx}{2y} """ s = "\\left(%s\\right)" % s
def extract_pow_y(self, k): """ Returns the power of `y` in `A` where self is `A dx/2y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-3*x+1) sage: x,y = C.monsky_washnitzer_gens() sage: A = y^5 - x*y^3 sage: A.extract_pow_y(5) [1, 0, 0, 0, 0] sage: (A * C.invariant_differential()).extract_pow_y(5) [1, 0, 0, 0, 0] """
def min_pow_y(self): """ Returns the minimum power of `y` in `A` where self is `A dx/2y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-3*x+1) sage: x,y = C.monsky_washnitzer_gens() sage: w = y^5 * C.invariant_differential() sage: w.min_pow_y() 5 sage: w = (x^2*y^4 + y^5) * C.invariant_differential() sage: w.min_pow_y() 4 """
def max_pow_y(self): """ Returns the maximum power of `y` in `A` where self is `A dx/2y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-3*x+1) sage: x,y = C.monsky_washnitzer_gens() sage: w = y^5 * C.invariant_differential() sage: w.max_pow_y() 5 sage: w = (x^2*y^4 + y^5) * C.invariant_differential() sage: w.max_pow_y() 5 """
def reduce_neg_y(self): """ Use homology relations to eliminate negative powers of `y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-3*x+1) sage: x,y = C.monsky_washnitzer_gens() sage: (y^-1).diff().reduce_neg_y() ((y^-1)*1, 0 dx/2y) sage: (y^-5*x^2+y^-1*x).diff().reduce_neg_y() ((y^-1)*x + (y^-5)*x^2, 0 dx/2y) """ cs = [a/j for a in reduced.extract_pow_y(j-1)] else:
def reduce_neg_y_fast(self, even_degree_only=False): """ Use homology relations to eliminate negative powers of `y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^5-3*x+1) sage: x, y = E.monsky_washnitzer_gens() sage: (y^-1).diff().reduce_neg_y_fast() ((y^-1)*1, 0 dx/2y) sage: (y^-5*x^2+y^-1*x).diff().reduce_neg_y_fast() ((y^-1)*x + (y^-5)*x^2, 0 dx/2y)
It leaves non-negative powers of `y` alone::
sage: y.diff() (-3*1 + 5*x^4) dx/2y sage: y.diff().reduce_neg_y_fast() (0, (-3*1 + 5*x^4) dx/2y) """ # prof = Profiler() # prof("reduce setup")
# prof("extract coeffs")
# prof("loop %s"%self.min_pow_y()) else: # this is a total hack to deal with the fact that we're using # rational numbers to approximate fixed precision p-adics except AttributeError: pass # g = lin_comb[i] x^i y^j # self -= dg
# prof("recreate forms")
def reduce_neg_y_faster(self, even_degree_only=False): """ Use homology relations to eliminate negative powers of `y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-3*x+1) sage: x,y = C.monsky_washnitzer_gens() sage: (y^-1).diff().reduce_neg_y() ((y^-1)*1, 0 dx/2y) sage: (y^-5*x^2+y^-1*x).diff().reduce_neg_y_faster() ((y^-1)*x + (y^-5)*x^2, 0 dx/2y) """ # Timings indicate that this is not any faster after all...
return S(0), self
# See monomial_diff_coeffs # this is the B_i and x_to_i contributions respectively for all i
else: # this is a total hack to deal with the fact that we're using # rational numbers to approximate fixed precision p-adics try: v = coeffs[j-offset-1] for kk in range(len(v)): a = v[kk] ppow = S._p**max(-a.valuation(S._p), 0) v[kk] = ((a * ppow) % S._prec_cap) / ppow except AttributeError: pass
# reduced = self - f.diff()
def reduce_pos_y(self): """ Use homology relations to eliminate positive powers of `y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^3-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: (y^2).diff().reduce_pos_y() (y^2*1, 0 dx/2y) sage: (y^2*x).diff().reduce_pos_y() (y^2*x, 0 dx/2y) sage: (y^92*x).diff().reduce_pos_y() (y^92*x, 0 dx/2y) sage: w = (y^3 + x).diff() sage: w += w.parent()(x) sage: w.reduce_pos_y_fast() (y^3*1 + x, x dx/2y) """
def reduce_pos_y_fast(self, even_degree_only=False): """ Use homology relations to eliminate positive powers of `y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^3-4*x+4) sage: x, y = E.monsky_washnitzer_gens() sage: y.diff().reduce_pos_y_fast() (y*1, 0 dx/2y) sage: (y^2).diff().reduce_pos_y_fast() (y^2*1, 0 dx/2y) sage: (y^2*x).diff().reduce_pos_y_fast() (y^2*x, 0 dx/2y) sage: (y^92*x).diff().reduce_pos_y_fast() (y^92*x, 0 dx/2y) sage: w = (y^3 + x).diff() sage: w += w.parent()(x) sage: w.reduce_pos_y_fast() (y^3*1 + x, x dx/2y) """
# self -= c d(y^{j+1})
# the others are basis elements
# self -= c d(x^{i+1} y^{j-1})
def reduce(self): """ Use homology relations to find `a` and `f` such that this element is equal to `a + df`, where `a` is given in terms of the `x^i dx/2y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = (y*x).diff() sage: w.reduce() (y*x, 0 dx/2y)
sage: w = x^4 * C.invariant_differential() sage: w.reduce() (1/5*y*1, 4/5*1 dx/2y)
sage: w = sum(QQ.random_element() * x^i * y^j for i in [0..4] for j in [-3..3]) * C.invariant_differential() sage: f, a = w.reduce() sage: f.diff() + a - w 0 dx/2y """
def reduce_fast(self, even_degree_only=False): """ Use homology relations to find `a` and `f` such that this element is equal to `a + df`, where `a` is given in terms of the `x^i dx/2y`.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: E = HyperellipticCurve(x^3-4*x+4) sage: x, y = E.monsky_washnitzer_gens() sage: x.diff().reduce_fast() (x, (0, 0)) sage: y.diff().reduce_fast() (y*1, (0, 0)) sage: (y^-1).diff().reduce_fast() ((y^-1)*1, (0, 0)) sage: (y^-11).diff().reduce_fast() ((y^-11)*1, (0, 0)) sage: (x*y^2).diff().reduce_fast() (y^2*x, (0, 0)) """ # f1, reduced = self.reduce_neg_y() # f2, reduced = reduced.reduce_pos_y()
def coeffs(self, R=None): """ Used to obtain the raw coefficients of a differential, see :meth:`SpecialHyperellipticQuotientElement.coeffs`
INPUT:
- R -- An (optional) base ring in which to cast the coefficients
OUTPUT:
The raw coefficients of $A$ where self is $A dx/2y$.
EXAMPLES::
sage: R.<x> = QQ['x'] sage: C = HyperellipticCurve(x^5-4*x+4) sage: x,y = C.monsky_washnitzer_gens() sage: w = C.invariant_differential() sage: w.coeffs() ([(1, 0, 0, 0, 0)], 0) sage: (x*w).coeffs() ([(0, 1, 0, 0, 0)], 0) sage: (y*w).coeffs() ([(0, 0, 0, 0, 0), (1, 0, 0, 0, 0)], 0) sage: (y^-2*w).coeffs() ([(1, 0, 0, 0, 0), (0, 0, 0, 0, 0), (0, 0, 0, 0, 0)], -2) """
def coleman_integral(self, P, Q): r""" Computes the definite integral of self from $P$ to $Q$.
INPUT:
- P, Q -- Two points on the underlying curve
OUTPUT:
`\int_P^Q \text{self}`
EXAMPLES::
sage: K = pAdicField(5,7) sage: E = EllipticCurve(K,[-31/3,-2501/108]) #11a sage: P = E(K(14/3), K(11/2)) sage: w = E.invariant_differential() sage: w.coleman_integral(P,2*P) O(5^6)
sage: Q = E([3,58332]) sage: w.coleman_integral(P,Q) 2*5 + 4*5^2 + 3*5^3 + 4*5^4 + 3*5^5 + O(5^6) sage: w.coleman_integral(2*P,Q) 2*5 + 4*5^2 + 3*5^3 + 4*5^4 + 3*5^5 + O(5^6) sage: (2*w).coleman_integral(P, Q) == 2*(w.coleman_integral(P, Q)) True """
integrate = coleman_integral |