Module: sage.schemes.hyperelliptic_curves.jacobian_morphism
Jacobian `morphism' as a class in the Picard group.
This module implements the group operation in the Picard group of a hyperelliptic curve, represented as divisors in Mumford representation, using Cantor's algorithm.
A divisor on the hyperelliptic curve
is stored in
Mumford representation, that is, as two polynomials
and
such
that:
*
is monic,
*
divides
,
*
.
REFERENCES:
A readable introduction to divisors, the Picard group, Mumford representation, and Cantor's algorithm:
J. Scholten, F. Vercauteren. An Introduction to Elliptic and Hyperelliptic Curve Cryptography and the NTRU Cryptosystem. To appear in B. Preneel (Ed.) State of the Art in Applied Cryptography - COSIC '03, Lecture Notes in Computer Science, Springer 2004.
A standard reference in the field of cryptography:
R. Avanzi, H. Cohen, C. Doche, G. Frey, T. Lange, K. Nguyen, and F. Vercauteren, Handbook of Elliptic and Hyperelliptic Curve Cryptography. CRC Press, 2005.
The following curve is the reduction of a curve whose Jacobian has complex multiplication.
sage: x = GF(37)['x'].gen() sage: H = HyperellipticCurve(x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x); H Hyperelliptic Curve over Finite Field of size 37 defined by y^2 = x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x
At this time, Jacobians of hyperelliptic curves are handled differently than elliptic curves:
sage: J = H.jacobian(); J Jacobian of Hyperelliptic Curve over Finite Field of size 37 defined by y^2 = x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x sage: J = J(J.base_ring()); J Set of points of Jacobian of Hyperelliptic Curve over Finite Field of size 37 defined by y^2 = x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x defined over Finite Field of size 37
Points on the Jacobian are represented by Mumford's polynomials. First we find a couple of points on the curve:
sage: P1 = H.lift_x(2); P1 (2 : 11 : 1) sage: Q1 = H.lift_x(10); Q1 (10 : 18 : 1)
Observe that 2 and 10 are the roots of the polynomials in x, respectively:
sage: P = J(P1); P (x + 35, y + 26) sage: Q = J(Q1); Q (x + 27, y + 19)
sage: P + Q (x^2 + 25*x + 20, y + 13*x) sage: (x^2 + 25*x + 20).roots(multiplicities=False) [10, 2]
Frobenius satisfies
on the Jacobian of this reduction and the order of the Jacobian is
sage: 1904*P (1) sage: 34*P == 0 True sage: 35*P == P True sage: 33*P == -P True
sage: Q*1904 (1) sage: Q*238 == 0 True sage: Q*239 == Q True sage: Q*237 == -Q True
Module-level Functions
D1, D2, f, h, genus) |
sage: F.<a> = GF(7^2, 'a') sage: x = F['x'].gen() sage: f = x^7 + x^2 + a sage: H = HyperellipticCurve(f, 2*x); H Hyperelliptic Curve over Finite Field in a of size 7^2 defined by y^2 + 2*x*y = x^7 + x^2 + a sage: J = H.jacobian()(F); J Set of points of Jacobian of Hyperelliptic Curve over Finite Field in a of size 7^2 defined by y^2 + 2*x*y = x^7 + x^2 + a defined over Finite Field in a of size 7^2
sage: Q = J(H.lift_x(F(1))); Q (x + 6, y + 2*a + 2) sage: 10*Q # indirect doctest (x^3 + (3*a + 1)*x^2 + (2*a + 5)*x + a + 5, y + (4*a + 5)*x^2 + (a + 1)*x + 6*a + 3) sage: 7*8297*Q (1)
sage: Q = J(H.lift_x(F(a+1))); Q (x + 6*a + 6, y + 2) sage: 7*8297*Q # indirect doctest (1)
D1, D2, f, genus) |
Given
and
two reduced Mumford divisors on the Jacobian of the
curve
, computes a representative
.
WARNING: the representative computed is NOT reduced! Use
cantor_reduction_simple
to reduce it.
sage: x = QQ['x'].gen() sage: f = x^5 + x sage: H = HyperellipticCurve(f); H Hyperelliptic Curve over Rational Field defined by y^2 = x^5 + x
sage: F.<a> = NumberField(x^2 - 2, 'a') sage: J = H.jacobian()(F); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 + x defined over Number Field in a with defining polynomial x^2 - 2
sage: P = J(H.lift_x(F(1))); P (x - 1, y - a) sage: Q = J(H.lift_x(F(0))); Q (x, y) sage: 2*P + 2*Q # indirect doctest (x^2 + (-2)*x + 1, y + (-3/2*a)*x + 1/2*a) sage: 2*(P + Q) # indirect doctest (x^2 + (-2)*x + 1, y + (-3/2*a)*x + 1/2*a) sage: 3*P # indirect doctest (x^2 + (-25/32)*x + 49/32, y + (-45/256*a)*x - 315/256*a)
a, b, f, h, genus) |
Return the unique reduced divisor linearly equivalent to
on the
curve
.
See the docstring of sage.schemes.hyperelliptic_curves.jacobian_morphism for information about divisors, linear equivalence, and reduction.
sage: x = QQ['x'].gen() sage: f = x^5 - x sage: H = HyperellipticCurve(f, x); H Hyperelliptic Curve over Rational Field defined by y^2 + x*y = x^5 - x sage: J = H.jacobian()(QQ); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 + x*y = x^5 - x defined over Rational Field
The following point is 2-torsion:
sage: Q = J(H.lift_x(0)); Q (x, y) sage: 2*Q # indirect doctest (1)
The next point is not 2-torsion:
sage: P = J(H.lift_x(-1)); P (x + 1, y - 1) sage: 2 * J(H.lift_x(-1)) # indirect doctest (x^2 + 2*x + 1, y - 3*x - 4) sage: 3 * J(H.lift_x(-1)) # indirect doctest (x^2 - 487*x - 324, y - 10754*x - 7146)
a, b, f, genus) |
Return the unique reduced divisor linearly equivalent to
on the
curve
.
See the docstring of sage.schemes.hyperelliptic_curves.jacobian_morphism for information about divisors, linear equivalence, and reduction.
sage: x = QQ['x'].gen() sage: f = x^5 - x sage: H = HyperellipticCurve(f); H Hyperelliptic Curve over Rational Field defined by y^2 = x^5 - x sage: J = H.jacobian()(QQ); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 - x defined over Rational Field
The following point is 2-torsion:
sage: P = J(H.lift_x(-1)); P (x + 1, y) sage: 2 * P # indirect doctest (1)
Class: JacobianMorphism_divisor_class_field
self, parent, polys, [check=True]) |
Create a new Jacobian element in Mumford representation.
Input:
parent: the parent Homset
polys: Mumford's
and
polynomials
check (default: True): if True, ensure that polynomials define a
divisor on the appropriate curve and are reduced
WARNING: not for external use! Use J(K)([u, v])
instead.
sage: x = GF(37)['x'].gen() sage: H = HyperellipticCurve(x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x) sage: J = H.jacobian()(GF(37))
sage: P1 = J(H.lift_x(2)); P1 # indirect doctest (x + 35, y + 26) sage: P1.parent() Set of points of Jacobian of Hyperelliptic Curve over Finite Field of size 37 defined by y^2 = x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x defined over Finite Field of size 37 sage: type(P1) <class 'sage.schemes.hyperelliptic_curves.jacobian_morphism.JacobianMorphis m_divisor_class_field'>
Functions: scheme
self) |
Return the scheme this morphism maps to; or, where this divisor lives.
WARNING: although a pointset is defined over a specific field, the scheme returned may be over a different (usually smaller) field. The example below demonstrates this: the pointset is determined over a number field of absolute degree 2 but the scheme returned is defined over the rationals.
sage: x = QQ['x'].gen() sage: f = x^5 + x sage: H = HyperellipticCurve(f) sage: F.<a> = NumberField(x^2 - 2, 'a') sage: J = H.jacobian()(F); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 + x defined over Number Field in a with defining polynomial x^2 - 2
sage: P = J(H.lift_x(F(1))) sage: P.scheme() Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 + x
Special Functions: __cmp__,
__getitem__,
__init__,
__list__,
__neg__,
__tuple__,
_add_,
_latex_,
_printing_polys,
_repr_,
_sub_
self, other) |
Compare self and other.
TESTS:
sage: x = QQ['x'].gen() sage: f = x^5 - x sage: H = HyperellipticCurve(f); H Hyperelliptic Curve over Rational Field defined by y^2 = x^5 - x sage: J = H.jacobian()(QQ); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 - x defined over Rational Field
The following point is 2-torsion:
sage: P = J(H.lift_x(-1)); P (x + 1, y) sage: 0 == 2 * P # indirect doctest True sage: P == P True
sage: Q = J(H.lift_x(-1)) sage: Q == P True
sage: 2 == Q False sage: P == False False
Let's verify the same "points" on different schemes are not equal:
sage: x = QQ['x'].gen() sage: f = x^5 + x sage: H2 = HyperellipticCurve(f) sage: J2 = H2.jacobian()(QQ)
sage: P1 = J(H.lift_x(0)); P1 (x, y) sage: P2 = J2(H2.lift_x(0)); P2 (x, y) sage: P1 == P2 False
self, n) |
Return the
-th item of the pair
of polynomials
giving the Mumford representation of self.
TESTS:
sage: x = QQ['x'].gen() sage: f = x^5 + x sage: H = HyperellipticCurve(f) sage: F.<a> = NumberField(x^2 - 2, 'a') sage: J = H.jacobian()(F); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 + x defined over Number Field in a with defining polynomial x^2 - 2
sage: P = J(H.lift_x(F(1))) sage: P[0] # indirect doctest x - 1 sage: P[1] # indirect doctest a sage: P[-1] # indirect doctest a sage: P[:1] # indirect doctest [x - 1]
self) |
Return a list
of the polynomials giving the Mumford
representation of self.
TESTS:
sage: x = QQ['x'].gen() sage: f = x^5 + x sage: H = HyperellipticCurve(f) sage: F.<a> = NumberField(x^2 - 2, 'a') sage: J = H.jacobian()(F); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 + x defined over Number Field in a with defining polynomial x^2 - 2
sage: P = J(H.lift_x(F(1))) sage: list(P) # indirect doctest [x - 1, a]
self) |
Return the additive inverse of this divisor.
sage: x = GF(37)['x'].gen() sage: H = HyperellipticCurve(x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x) sage: J = H.jacobian()(GF(37)) sage: P1 = J(H.lift_x(2)); P1 (x + 35, y + 26) sage: - P1 # indirect doctest (x + 35, y + 11) sage: P1 - P1 # indirect doctest (1)
sage: H2 = HyperellipticCurve(x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x, x) sage: J2 = H2.jacobian()(GF(37)) sage: P2 = J2(H2.lift_x(2)); P2 (x + 35, y + 15) sage: - P2 # indirect doctest (x + 35, y + 24) sage: P2 - P2 # indirect doctest (1)
self) |
Return a tuple
of the polynomials giving the Mumford
representation of self.
TESTS:
sage: x = QQ['x'].gen() sage: f = x^5 + x sage: H = HyperellipticCurve(f) sage: F.<a> = NumberField(x^2 - 2, 'a') sage: J = H.jacobian()(F); J Set of points of Jacobian of Hyperelliptic Curve over Rational Field defined by y^2 = x^5 + x defined over Number Field in a with defining polynomial x^2 - 2
sage: P = J(H.lift_x(F(1))) sage: tuple(P) # indirect doctest (x - 1, a)
self, other) |
Return a Mumford representative of the divisor self + other.
sage: x = GF(37)['x'].gen() sage: H = HyperellipticCurve(x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x) sage: J = H.jacobian()(GF(37))
sage: P1 = J(H.lift_x(2)); P1 (x + 35, y + 26) sage: P1 + P1 # indirect doctest (x^2 + 33*x + 4, y + 13*x)
self) |
Return a LaTeXstring representing this Mumford divisor.
sage: F.<alpha> = GF(7^2) sage: x = F['x'].gen() sage: f = x^7 + x^2 + alpha sage: H = HyperellipticCurve(f, 2*x) sage: J = H.jacobian()(F)
sage: Q = J(0); print latex(Q) # indirect doctest \left(1\right) sage: Q = J(H.lift_x(F(1))); print latex(Q) # indirect doctest \left(x + 6, y + 2 \alpha + 2\right)
sage: print latex(Q + Q) \left(x^{2} + 5 x + 1, y + 3 \alpha x + 6 \alpha + 2\right)
self) |
Internal function formatting Mumford polynomials for printing.
TESTS:
sage: F.<a> = GF(7^2, 'a') sage: x = F['x'].gen() sage: f = x^7 + x^2 + a sage: H = HyperellipticCurve(f, 2*x) sage: J = H.jacobian()(F)
sage: Q = J(H.lift_x(F(1))); Q # indirect doctest (x + 6, y + 2*a + 2)
self) |
Return a string representation of this Mumford divisor.
sage: F.<a> = GF(7^2, 'a') sage: x = F['x'].gen() sage: f = x^7 + x^2 + a sage: H = HyperellipticCurve(f, 2*x) sage: J = H.jacobian()(F)
sage: Q = J(0); Q # indirect doctest (1) sage: Q = J(H.lift_x(F(1))); Q # indirect doctest (x + 6, y + 2*a + 2) sage: Q + Q # indirect doctest (x^2 + 5*x + 1, y + 3*a*x + 6*a + 2)
self, other) |
Return a Mumford representative of the divisor self - other.
sage: x = GF(37)['x'].gen() sage: H = HyperellipticCurve(x^5 + 12*x^4 + 13*x^3 + 15*x^2 + 33*x) sage: J = H.jacobian()(GF(37))
sage: P1 = J(H.lift_x(2)); P1 (x + 35, y + 26) sage: P1 - P1 # indirect doctest (1)
sage: P2 = J(H.lift_x(4)); P2 (x + 33, y + 34)
Observe that the
-coordinates are the same but the
-coordinates differ:
sage: P1 - P2 # indirect doctest (x^2 + 31*x + 8, y + 7*x + 12) sage: P1 + P2 # indirect doctest (x^2 + 31*x + 8, y + 4*x + 18) sage: (P1 - P2) - (P1 + P2) + 2*P2 # indirect doctest (1)
See About this document... for information on suggesting changes.