Module: sage.schemes.elliptic_curves.ell_finite_field
Elliptic curves over finite fields
Author Log:
Class: EllipticCurve_finite_field
self, x, [y=None]) |
Special constructor for elliptic curves over a finite field
sage: EllipticCurve(GF(101),[2,3]) Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Finite Field of size 101
sage: F=GF(101^2, 'a') sage: EllipticCurve([F(2),F(3)]) Elliptic Curve defined by y^2 = x^3 + 2*x + 3 over Finite Field in a of size 101^2
Functions: abelian_group,
cardinality,
cardinality_exhaustive,
cardinality_from_group,
frobenius,
frobenius_order,
frobenius_polynomial,
gens,
order,
plot,
points,
random_element,
random_point,
trace_of_frobenius
self, [debug=False]) |
Returns the abelian group structure of the group of points on this elliptic curve.
WARNING: The algorithm is definitely *not* intended for use
with *large* finite fields! The factorization of the
orders of elements must be feasible. Also,
baby-step-giant-step methods are used which have space and
time requirements which are
.
Also, the algorithm uses random points on the curve and hence the generators are likely to differ from one run to another; but the group is cached so the generators will not change in any one run of Sage.
Note: This function applies to elliptic curves over arbitrary finite fields. The related function abelian_group_prime_field() uses the pari script, for prime fields only; it is now obsolete
Input:
Author: John Cremona
sage: E=EllipticCurve(GF(11),[2,5]) sage: E.abelian_group() (Multiplicative Abelian Group isomorphic to C10, ...
sage: E=EllipticCurve(GF(41),[2,5]) sage: E.abelian_group() (Multiplicative Abelian Group isomorphic to C22 x C2, ...
sage: F.<a>=GF(3^6,'a') sage: E=EllipticCurve([a^4 + a^3 + 2*a^2 + 2*a, 2*a^5 + 2*a^3 + 2*a^2 + 1]) sage: E.abelian_group() (Multiplicative Abelian Group isomorphic to C26 x C26, ...
sage: F.<a>=GF(101^3,'a') sage: E=EllipticCurve([2*a^2 + 48*a + 27, 89*a^2 + 76*a + 24]) sage: E.abelian_group() (Multiplicative Abelian Group isomorphic to C1031352, ...
The group can be trivial (but this is the only example of that up to isomorphism!)
sage: E=EllipticCurve(GF(2),[0,0,1,1,1]) sage: E.abelian_group() (Trivial Abelian Group, ())
Of course, there are plenty of points if we extend the field:
sage: E.cardinality(extension_degree=100) 1267650600228231653296516890625
This tests the patch for trac #3111, using 10 primes randomly selected:
sage: E = EllipticCurve('389a') sage: for p in [5927, 2297, 1571, 1709, 3851, 127, 3253, 5783, 3499, 4817]: ... G = E.change_ring(GF(p)).abelian_group() sage: for p in prime_range(10000): #long time (~20s) ... if p != 389: ... G=E.change_ring(GF(p)).abelian_group()
self, [algorithm=heuristic], [early_abort=False], [disable_warning=False], [extension_degree=1]) |
Return the number of points on this elliptic curve over an extension field (default: the base field).
Input:
The cardinality is cached.
Over prime fields, one of the above algorithms is used. Over non-prime fields, the serious point counting is done on a standard curve with the same j-invariant over the field GF(p)(j), then lifted to the base_field, and finally account is taken of twists.
For j=0 and j=1728 special formulas are used instead.
sage: EllipticCurve(GF(4,'a'),[1,2,3,4,5]).cardinality() 8 sage: k.<a> = GF(3^3) sage: l = [a^2 + 1, 2*a^2 + 2*a + 1, a^2 + a + 1, 2, 2*a] sage: EllipticCurve(k,l).cardinality() 29
sage: l = [1, 1, 0, 2, 0] sage: EllipticCurve(k,l).cardinality() 38
An even bigger extension (which we check against Magma):
sage: EllipticCurve(GF(3^100,'a'),[1,2,3,4,5]).cardinality() 515377520732011331036459693969645888996929981504 sage: magma.eval("Order(EllipticCurve([GF(3^100)|1,2,3,4,5]))") # optional -- requires magma '515377520732011331036459693969645888996929981504'
sage: EllipticCurve(GF(10007),[1,2,3,4,5]).cardinality() 10076 sage: EllipticCurve(GF(10007),[1,2,3,4,5]).cardinality(algorithm='sea') 10076 sage: EllipticCurve(GF(10007),[1,2,3,4,5]).cardinality(algorithm='bsgs') 10076 sage: EllipticCurve(GF(next_prime(10**20)),[1,2,3,4,5]).cardinality(algorithm='sea') 100000000011093199520
The cardinality is cached:
sage: E = EllipticCurve(GF(3^100,'a'),[1,2,3,4,5]) sage: E.cardinality() is E.cardinality() True sage: E=EllipticCurve(GF(11^2,'a'),[3,3]) sage: E.cardinality() 128 sage: EllipticCurve(GF(11^100,'a'),[3,3]).cardinality() 137806123398222701841183371720896367762643312000384671846835266941791510341 065565176497846502742959856128
self) |
Return the cardinality of self over the base field. Simply adds up the number of points with each x-coordinate: only used for small field sizes!
sage: p=next_prime(10^3) sage: E=EllipticCurve(GF(p),[3,4]) sage: E.cardinality_exhaustive() 1020 sage: E=EllipticCurve(GF(3^4,'a'),[1,1]) sage: E.cardinality_exhaustive() 64
self) |
Return the cardinality of self over the base field. Will be called by user function cardinality only when necessary, i.e. when the j_invariant is not in the prime field.
This function just calls abelian_group(), so results in the group structure and generators being cached as well as the group order.
sage: p=next_prime(10^3) sage: E=EllipticCurve(GF(p),[3,4]) sage: E.cardinality_from_group() 1020 sage: E=EllipticCurve(GF(3^4,'a'),[1,1]) sage: E.cardinality_from_group() 64
self) |
Return the frobenius of self as an element of a quadratic order
NOTES: This computes the curve cardinality, which may be time-consuming.
Frobenius is only determined up to conjugacy.
sage: E=EllipticCurve(GF(11),[3,3]) sage: E.frobenius() phi sage: E.frobenius().minpoly() x^2 - 4*x + 11
For some supersingular curves, Frobenius is in Z:
sage: E=EllipticCurve(GF(25,'a'),[0,0,0,0,1]) sage: E.frobenius() -5
self) |
Return the quadratic order Z[phi] where phi is the Frobenius endomorphism of the elliptic curve
NOTE: This computes the curve cardinality, which may be time-consuming.
sage: E=EllipticCurve(GF(11),[3,3]) sage: E.frobenius_order() Order in Number Field in phi with defining polynomial x^2 - 4*x + 11
For some supersingular curves, Frobenius is in Z and the Frobenius order is Z:
sage: E=EllipticCurve(GF(25,'a'),[0,0,0,0,1]) sage: R=E.frobenius_order() sage: R Order in Number Field in phi with defining polynomial x + 5 sage: R.degree() 1
self) |
Return the characteristic polynomial of Frobenius.
The Frobenius endomorphism of the elliptic curve has quadratic
characteristic polynomial. In most cases this is irreducible
and defines an imaginary quadratic order; for some
supersingular curves, Frobenius is an integer a and the
polynomial is
.
NOTE: This computes the curve cardinality, which may be time-consuming.
sage: E=EllipticCurve(GF(11),[3,3]) sage: E.frobenius_polynomial() x^2 - 4*x + 11
For some supersingular curves, Frobenius is in Z and the polynomial is a square:
sage: E=EllipticCurve(GF(25,'a'),[0,0,0,0,1]) sage: E.frobenius_polynomial().factor() (x + 5)^2
self) |
Returns a tuple of length up to 2 of points which generate the abelian group of points on this elliptic curve. See abelian_group() for limitations.
The algorithm uses random points on the curve, and hence the generators are likely to differ from one run to another; but they are cached so will be consistent in any one run of Sage.
Author: John Cremona
sage: E=EllipticCurve(GF(11),[2,5]) # random output sage: E.gens() ((0 : 4 : 1),) sage: EllipticCurve(GF(41),[2,5]).gens() # random output ((21 : 1 : 1), (8 : 0 : 1)) sage: F.<a>=GF(3^6,'a') sage: E=EllipticCurve([a,a+1]) sage: pts=E.gens() sage: len(pts) 1 sage: pts[0].order()==E.cardinality() True
self, [algorithm=heuristic], [early_abort=False], [disable_warning=False], [extension_degree=1]) |
Return the number of points on this elliptic curve over an extension field (default: the base field).
Input:
The cardinality is cached.
Over prime fields, one of the above algorithms is used. Over non-prime fields, the serious point counting is done on a standard curve with the same j-invariant over the field GF(p)(j), then lifted to the base_field, and finally account is taken of twists.
For j=0 and j=1728 special formulas are used instead.
sage: EllipticCurve(GF(4,'a'),[1,2,3,4,5]).cardinality() 8 sage: k.<a> = GF(3^3) sage: l = [a^2 + 1, 2*a^2 + 2*a + 1, a^2 + a + 1, 2, 2*a] sage: EllipticCurve(k,l).cardinality() 29
sage: l = [1, 1, 0, 2, 0] sage: EllipticCurve(k,l).cardinality() 38
An even bigger extension (which we check against Magma):
sage: EllipticCurve(GF(3^100,'a'),[1,2,3,4,5]).cardinality() 515377520732011331036459693969645888996929981504 sage: magma.eval("Order(EllipticCurve([GF(3^100)|1,2,3,4,5]))") # optional -- requires magma '515377520732011331036459693969645888996929981504'
sage: EllipticCurve(GF(10007),[1,2,3,4,5]).cardinality() 10076 sage: EllipticCurve(GF(10007),[1,2,3,4,5]).cardinality(algorithm='sea') 10076 sage: EllipticCurve(GF(10007),[1,2,3,4,5]).cardinality(algorithm='bsgs') 10076 sage: EllipticCurve(GF(next_prime(10**20)),[1,2,3,4,5]).cardinality(algorithm='sea') 100000000011093199520
The cardinality is cached:
sage: E = EllipticCurve(GF(3^100,'a'),[1,2,3,4,5]) sage: E.cardinality() is E.cardinality() True sage: E=EllipticCurve(GF(11^2,'a'),[3,3]) sage: E.cardinality() 128 sage: EllipticCurve(GF(11^100,'a'),[3,3]).cardinality() 137806123398222701841183371720896367762643312000384671846835266941791510341 065565176497846502742959856128
self) |
Draw a graph of this elliptic curve over a prime finite field.
Input:
sage: E = EllipticCurve(FiniteField(17), [0,1]) sage: P = plot(E, rgbcolor=(0,0,1))
self) |
All the points on this elliptic curve. The list of points is cached so subsequent calls are free.
sage: p = 5 sage: F = GF(p) sage: E = EllipticCurve(F, [1, 3]) sage: a_sub_p = E.change_ring(QQ).ap(p); a_sub_p 2
sage: len(E.points()) 4 sage: p + 1 - a_sub_p 4 sage: E.points() [(0 : 1 : 0), (1 : 0 : 1), (4 : 1 : 1), (4 : 4 : 1)]
sage: K = GF(p**2,'a') sage: E = E.change_ring(K) sage: len(E.points()) 32 sage: (p + 1)**2 - a_sub_p**2 32 sage: w = E.points(); w [(0 : 1 : 0), (0 : 2*a + 4 : 1), (0 : 3*a + 1 : 1), (1 : 0 : 1), (2 : 2*a + 4 : 1), (2 : 3*a + 1 : 1), (3 : 2*a + 4 : 1), (3 : 3*a + 1 : 1), (4 : 1 : 1), (4 : 4 : 1), (a : 1 : 1), (a : 4 : 1), (a + 2 : a + 1 : 1), (a + 2 : 4*a + 4 : 1), (a + 3 : a : 1), (a + 3 : 4*a : 1), (a + 4 : 0 : 1), (2*a : 2*a : 1), (2*a : 3*a : 1), (2*a + 4 : a + 1 : 1), (2*a + 4 : 4*a + 4 : 1), (3*a + 1 : a + 3 : 1), (3*a + 1 : 4*a + 2 : 1), (3*a + 2 : 2*a + 3 : 1), (3*a + 2 : 3*a + 2 : 1), (4*a : 0 : 1), (4*a + 1 : 1 : 1), (4*a + 1 : 4 : 1), (4*a + 3 : a + 3 : 1), (4*a + 3 : 4*a + 2 : 1), (4*a + 4 : a + 4 : 1), (4*a + 4 : 4*a + 1 : 1)]
Note that the returned list is an immutable sorted Sequence:
sage: w[0] = 9 Traceback (most recent call last): ... ValueError: object is immutable; please change a copy instead.
self) |
Returns a random point on this elliptic curve.
Returns the point at infinity with probability
where the base field has cardinality
.
sage: k = GF(next_prime(7^5)) sage: E = EllipticCurve(k,[2,4]) sage: P = E.random_element(); P (751 : 6230 : 1) sage: type(P) <class 'sage.schemes.elliptic_curves.ell_point.EllipticCurvePoint_finite_fi eld'> sage: P in E True
sage: k.<a> = GF(7^5) sage: E = EllipticCurve(k,[2,4]) sage: P = E.random_element(); P (a^4 + a + 5 : 6*a^4 + 3*a^3 + 2*a^2 + 4 : 1) sage: type(P) <class 'sage.schemes.elliptic_curves.ell_point.EllipticCurvePoint_finite_fi eld'> sage: P in E True
sage: k.<a> = GF(2^5) sage: E = EllipticCurve(k,[a^2,a,1,a+1,1]) sage: P = E.random_element(); P (a^4 : 0 : 1) sage: type(P) <class 'sage.schemes.elliptic_curves.ell_point.EllipticCurvePoint_finite_fi eld'> sage: P in E True
self) |
Returns a random point on this elliptic curve.
Returns the point at infinity with probability
where the base field has cardinality
.
sage: k = GF(next_prime(7^5)) sage: E = EllipticCurve(k,[2,4]) sage: P = E.random_element(); P (751 : 6230 : 1) sage: type(P) <class 'sage.schemes.elliptic_curves.ell_point.EllipticCurvePoint_finite_fi eld'> sage: P in E True
sage: k.<a> = GF(7^5) sage: E = EllipticCurve(k,[2,4]) sage: P = E.random_element(); P (a^4 + a + 5 : 6*a^4 + 3*a^3 + 2*a^2 + 4 : 1) sage: type(P) <class 'sage.schemes.elliptic_curves.ell_point.EllipticCurvePoint_finite_fi eld'> sage: P in E True
sage: k.<a> = GF(2^5) sage: E = EllipticCurve(k,[a^2,a,1,a+1,1]) sage: P = E.random_element(); P (a^4 : 0 : 1) sage: type(P) <class 'sage.schemes.elliptic_curves.ell_point.EllipticCurvePoint_finite_fi eld'> sage: P in E True
self) |
Return the trace of Frobenius acting on this elliptic curve.
NOTE: This computes the curve cardinality, which may be time-consuming.
sage: E=EllipticCurve(GF(101),[2,3]) sage: E.trace_of_frobenius() 6 sage: E=EllipticCurve(GF(11^5,'a'),[2,5]) sage: E.trace_of_frobenius() 802
The following shows that the issue from trac #2849 is fixed:
sage: E=EllipticCurve(GF(3^5,'a'),[-1,-1]) sage: E.trace_of_frobenius() -27
Special Functions: __getitem__,
__init__,
_cardinality_with_j_invariant_0,
_cardinality_with_j_invariant_1728,
_gp,
_magma_init_,
_pari_,
_points_via_group_structure
self, n) |
Return the n'th point in self's __points list. This enables users to iterate over the curve's point set.
sage: E=EllipticCurve(GF(97),[2,3]) sage: S=E.points() sage: E[10] (10 : 76 : 1) sage: E[15] (17 : 10 : 1) sage: for P in E: print P.order() 1 50 50 50 50 5 5 50 ...
self) |
Special function to compute cardinality when j=0.
An example with q=p=1 (mod 6)
sage: F=GF(1009) sage: [EllipticCurve(F,[0,0,0,0,11^i])._cardinality_with_j_invariant_0() for i in range(6)] [948, 967, 1029, 1072, 1053, 991]
An example with q=p=5 (mod 6)
sage: F=GF(1013) sage: [EllipticCurve(F,[0,0,0,0,3^i])._cardinality_with_j_invariant_0() for i in range(6)] [1014, 1014, 1014, 1014, 1014, 1014]
An example with
, p=5 (mod 6)
sage: F.<a>=GF(1013^2,'a') sage: [EllipticCurve(F,[0,0,0,0,a^i])._cardinality_with_j_invariant_0() for i in range(6)] [1028196, 1027183, 1025157, 1024144, 1025157, 1027183]
For examples in characteristic 2 and 3, see the function _cardinality_with_j_invariant_1728()
self) |
Special function to compute cardinality when j=1728.
An example with q=p=1 (mod 4)
sage: F=GF(10009) sage: [EllipticCurve(F,[0,0,0,11^i,0])._cardinality_with_j_invariant_1728() for i in range(4)] [10016, 10210, 10004, 9810]
An example with q=p=3 (mod 4)
sage: F=GF(10007) sage: [EllipticCurve(F,[0,0,0,5^i,0])._cardinality_with_j_invariant_1728() for i in range(4)] [10008, 10008, 10008, 10008]
An example with
, p=3 (mod 4)
sage: F.<a>=GF(10007^2,'a') sage: [EllipticCurve(F,[0,0,0,a^i,0])._cardinality_with_j_invariant_1728() for i in range(4)] [100160064, 100140050, 100120036, 100140050]
Examples with
, d odd (3 isomorphism classes):
sage: F.<a> = GF(2**15,'a') sage: ais = [[0,0,1,0,0],[0,0,1,1,0],[0,0,1,1,1]] sage: curves=[EllipticCurve(F,ai) for ai in ais] sage: all([all([e1==e2 or not e1.is_isomorphic(e2) for e1 in curves]) for e2 in curves]) True sage: [e._cardinality_with_j_invariant_1728() for e in curves] [32769, 33025, 32513]
Examples with
, d even (7 isomorphism classes):
sage: F.<a> = GF(2**16,'a') sage: b = a^11 # trace 1 sage: ais = [[0,0,1,0,0],[0,0,1,0,b],[0,0,1,b,0],[0,0,a,0,0],[0,0,a,0,a^2*b],[0,0,a^2,0,0],[0,0,a^2,0,a^4*b]] sage: curves=[EllipticCurve(F,ai) for ai in ais] sage: all([all([e1==e2 or not e1.is_isomorphic(e2) for e1 in curves]) for e2 in curves]) True sage: [e._cardinality_with_j_invariant_1728() for e in curves] [65025, 66049, 65537, 65793, 65281, 65793, 65281]
Examples with
, d odd (4 isomorphism classes):
sage: F.<a> = GF(3**15,'a') sage: b=a^7 # has trace 1 sage: ais=[[0,0,0,1,0],[0,0,0,-1,0],[0,0,0,-1,b],[0,0,0,-1,-b]] sage: curves=[EllipticCurve(F,ai) for ai in ais] sage: all([all([e1==e2 or not e1.is_isomorphic(e2) for e1 in curves]) for e2 in curves]) True sage: [e._cardinality_with_j_invariant_1728() for e in curves] [14348908, 14348908, 14342347, 14355469]
Examples with
, d even (6 isomorphism classes):
sage: F.<g>=GF(3^18,'g') sage: i=F(-1).sqrt() sage: a=g^8 # has trace 1 sage: ais= [[0,0,0,1,0],[0,0,0,1,i*a],[0,0,0,g,0],[0,0,0,g^3,0],[0,0,0,g^2,0], [0,0,0,g^2,i*a*g^3]] sage: curves=[EllipticCurve(F,ai) for ai in ais] sage: all([all([e1==e2 or not e1.is_isomorphic(e2) for e1 in curves]) for e2 in curves]) True sage: [E._cardinality_with_j_invariant_1728() for E in curves] [387459856, 387400807, 387420490, 387420490, 387381124, 387440173]
self) |
Return an elliptic curve in a GP/PARI interpreter with all Cremona's code for finite fields preloaded. This includes generators, which will vary from run to run.
The base field must have prime order.
sage: EllipticCurve(GF(41),[2,5])._gp() [Mod(0, 41), Mod(0, 41), Mod(0, 41), Mod(2, 41), Mod(5, 41), Mod(0, 41), Mod(4, 41), Mod(20, 41), Mod(37, 41), Mod(27, 41), Mod(26, 41), Mod(4, 41), Mod(11, 41), 44, [2, 2; 11, 1], [22, 2], ...
self) |
Return a Magma command that creates this curve.
sage: EllipticCurve(GF(41),[2,5])._magma_init_() # optional -- requires Magma 'EllipticCurve([_sage_[1]|GF(41)!0,GF(41)!0,GF(41)!0,GF(41)!2,GF(41)!5])' sage: magma(E) # optional -- requires Magma Elliptic Curve defined by y^2 = x^3 + 2*x + 5 over GF(41)
self) |
Return a GP/PARI elliptic curve
sage: EllipticCurve(GF(41),[2,5])._pari_() [Mod(0, 41), Mod(0, 41), Mod(0, 41), Mod(2, 41), Mod(5, 41), Mod(0, 41), Mod(4, 41), Mod(20, 41), Mod(37, 41), Mod(27, 41), Mod(26, 41), Mod(4, 41), Mod(11, 41), 0, 0, 0, 0, 0, 0]
self) |
Return a list of all the points on the curve, for prime fields only (see points() for the general case)
sage: S=EllipticCurve(GF(97),[2,3])._points_via_group_structure() sage: len(S) 100