The Order of the Group of Points

Point Counting

Magma contains an efficient implementation of the Schoof–Elkies–Atkin (SEA) algorithm, with Lercier’s extension to base fields of characteristic 2, for computing the number of points on an elliptic curve over a finite field. There are also implementations of the faster \(p\)-adic canonical lift methods and the \(p\)-adic deformation method in small characteristic. Calculations are performed in the smallest field over which the curve is defined, and the result is lifted to the original field.

Like all group functions on elliptic curves, these intrinsics really apply to a particular point set; the curve is identified with its base point set for the purposes of these functions. To aid exposition only the versions that take the curves are shown but an appropriate point set (over a finite field) may be used instead.

# H: SetPtEll -> RngIntElt
# E: CrvEll -> RngIntElt
Order(H): SetPtEll -> RngIntElt
Order(E): CrvEll -> RngIntElt: The order of the group of \(K\)-rational points of \(E\), where \(E\) is an elliptic curve defined over the finite field \(K\).

FactoredOrder(H): SetPtEll -> RngIntElt
FactoredOrder(E): CrvEll -> RngIntElt: This function returns the factorisation of the order of the group of \(K\)-rational points of \(E\), where \(E\) is an elliptic curve defined over the finite field \(K\). This factorisation is stored on \(E\) and is reused when computing the order of points on \(E\).

SEA(H : parameters): SetPtEll -> RngIntElt

SEA(E : parameters): CrvEll -> RngIntElt

Al        : MonStgElt                    Default: "Default"
MaxSmooth : RngIntElt                    Default: Infinity()
AbortLevel: RngIntElt                    Default: 0
UseSEA    : BoolElt                      Default: false

This is the internal algorithm used by the point counting routines, but it can be invoked directly on an elliptic curve \(E\) or a point set \(H\) if desired. The most obvious reason to do this is because of the parameters AbortLevel and MaxSmooth, which allow the computation to be stopped early if it can be shown not to satisfy certain conditions (such as the order of the group being prime).

Warning: Unlike the external functions which call it (such as Order and FactoredOrder), SEA does not store the computed group order back on the curve itself and so functions which need the group order will have to recompute it.

The parameter Al can be set to one of "Enumerate", "BSGS", or "Default". Setting Al to "Enumerate" causes the group order to be found by enumeration of all points on the curve and hence is only practical for small orders. Setting Al to "BSGS" is currently equivalent to setting Al to "Default", which uses point enumeration for suitably small finite fields, the full Schoof–Elkies–Atkin algorithm if the characteristic is not small, and canonical lift methods or deformation if it is. Currently, small characteristic means \(p < 1000\). For \(p < 37\) or \(p \in \{41,47,59,71\}\), canonical lift is used; deformation is the chosen method for other small \(p\).

The canonical lift methods use a range of techniques to lift a modular parameter of the curve to an unramified \(p\)-adic ring. If a Gaussian Normal Basis of type \(\le 2\) exists for that ring (see HasGNB) then the method of Lercier and Lubicz from [LL03] is used with some local adaptations by Michael Harrison to increase the speed for fields of moderate size. Otherwise, the MSST algorithm as described in [Gau02] is used for fields up to a certain size with Harley’s recursive version ([Hara]) taking over for larger fields.

For \(p < 10\) and \(p = 13\) the modular parameter used is a generator of the function field of a genus 0 modular curve \(X_0(p^r)\). The parameter gives a particularly nice modular polynomial \(\Phi\) for use in the above-mentioned iterative lifting processes (see [Koh03]). In these cases \(\Phi\) has a relatively simple factored form and evaluations are hand-coded for extra efficiency.

For \(p\) where \(X_0(p)\) is hyperelliptic (excepting the anomalous case \(p = 37\)) we use a new method, developed by Michael Harrison, which lifts a pair of modular parameters corresponding to the \(x\) and \(y\) coordinates in a simple Weierstrass model of \(X_0(p)\).

The deformation method uses code designed and implemented by Hendrik Hubrechts. It is based on ideas of Dwork, first used in the design of explicit computational algorithms by Alan Lauder, and works by computing the deformation over a family of curves of the Frobenius action on rigid cohomology by solving a certain differential equation \(p\)-adically. Details can be found in [Hub07]. The method is used for \(E\) over \({\mathbb{F}}_{p^n}\) with \(p < 1000\) for which the canonical lift method is not used and for \(n\) greater than a smallish bound that increases slowly with \(p\). For such \(p\) and smaller \(n\), either enumeration or small characteristic SEA is used.

To use Schoof–Elkies–Atkin in small characteristic (with Lercier’s improvements) instead of canonical lift or deformation, set UseSEA to true.

The parameter MaxSmooth can be used to specify a limit on the number of small primes that divide the group order. That is, we will consider the group order to be “smooth” if it is divisible by the product of small primes and this product is larger than the value of MaxSmooth. Then the possible values of AbortLevel have the following meanings:

•: 0 : Abort if the curve has smooth order.
•: 1 : Abort if both the curve and its twist have smooth order.
•: 2 : Abort if either the curve or its twist have smooth order.

If the early abort is triggered then the returned group order is \(0\).

One common use of these parameters is when searching for a curve with prime order. Setting MaxSmooth := 1 will cause the early termination of the algorithm if a small factor of the order is found. If, however, the algorithm did not abort then the group order is not necessarily prime — this just means that no small prime dividing the order was encountered during the computation.

Note that the canonical lift methods give no intermediate data on the order of \(E\), so MaxSmooth and AbortLevel have no effect when these methods are used.

Example: SEA Point Counting in Characteristic 2 (ex-05153b)

Run in calculator (opens in a new tab)

The first examples in characteristic 2 illustrate the increased speed when we work with fields that have GNBs available:

> //example with no Gaussian Normal Basis (GNB)
> K<w> := FiniteField(2, 160); // finite field of size 2^160
> f<x> := MinimalPolynomial(w); f;
x^160 + x^5 + x^3 + x^2 + 1
> E := EllipticCurve([K| 1, 0, 0, 0, w]);
> time #E;
1461501637330902918203686141511652467686942715904
Time: 0.050
> //example with a GNB of Type 1
> HasGNB(pAdicRing(2, 1), 162, 1);
true
> K<w> := FiniteField(2, 162); // finite field of size 2^162
> f<x> := MinimalPolynomial(w); f;
x^162 + x^7 + x^6 + x^5 + x^2 + x + 1
> E := EllipticCurve([K| 1, 0, 0, 0, w]);
> time #E;
5846006549323611672814742626685174784085558341632
Time: 0.010
> //example with a GNB of Type 2
> K<w> := FiniteField(2, 173); // finite field of size 2^173
> f<x> := MinimalPolynomial(w); f;
x^173 + x^8 + x^5 + x^2 + 1
> E := EllipticCurve([K| 1, 0, 0, 0, w]);
> time #E;
11972621413014756705924586057649971923911742202056704
Time: 0.030

and here is an example in characteristic 3 (with type 1 GNB)

> F := FiniteField(3, 100);
> j := Random(F);
> E := EllipticCurveWithjInvariant(j);
> time #E;
515377520732011331036461505619702343613256090042
Time: 0.010

Finally (for small characteristic), a case with \(p = 37\) where the deformation method is used.

> F := FiniteField(37, 30);
> j := Random(F);
> E := EllipticCurveWithjInvariant(j);
> time #E;
111186413610993811950186296693286649955782417767
Time: 0.700

Here we can see the early abort feature in action:

> p := NextPrime(10^9);
> p;
1000000007
> K := GF(p);
> E := EllipticCurve([K | -1, 1]);
> time SEA(E : MaxSmooth := 1);
0
Time: 0.010
> time #E;
1000009476
Time: 0.070

For curves this small the amount of time saved by an early abort is fairly small, but for curves over large fields the savings can be quite significant. As mentioned above, even when SEA does not abort there is no guarantee that the curve has prime order:

> E := EllipticCurve([K | 0, 1]);
> time SEA(E : MaxSmooth := 1);
1000000008
> IsSupersingular(E);
true
> E := EllipticCurve([K | -30, 1]);
> time SEA(E : MaxSmooth := 1);
1000035283
> Factorization($1);
[ <13, 1>, <709, 1>, <108499, 1> ]

In the first case the curve was supersingular and so the order was easily computed directly; thus no early abort was attempted. The latter case shows that small primes may not be looked at even when the curve is ordinary.

Next we find a curve with prime order:

> for k in [1..p] do
>   E := EllipticCurve([K | k, 1]);
>   n := SEA(E : MaxSmooth := 1);
>   if IsPrime(n) then
>     printf "Found curve of prime order %o for k = %o\\n", n, k;
>     break;
>   end if;
> end for;
Found curve of prime order 999998501 for k = 29
> E;
Elliptic Curve defined by y^2 = x^3 + 29*x + 1 over GF(1000000007)
> #E;
999998501

Run in calculator (opens in a new tab)

SetVerbose("SEA", v): MonStgElt, RngIntElt

This procedure changes the verbose printing level for the SEA algorithm which counts the number of points on an elliptic curve. Currently the legal values for \(v\) are false, true, or an integer between \(0\) and \(5\) (inclusive), where false is equivalent to \(0\) and true is equivalent to \(1\). A nonzero value gives information during the course of the algorithm, increasing in verbosity for higher values of \(v\).

N.B. The previous name for this verbose flag, "ECPointCount" is now deprecated and will be removed in a future release. It should continue to function in this release, but its verbose levels are different to those of "SEA".

Order(H, r): SetPtEll, RngIntElt -> RngIntElt
Order(E, r): CrvEll, RngIntElt -> RngIntElt: Returns the order of the elliptic curve \(E\) over the extension field of \(K\) of degree \(r\), where \(K\) is the base field of \(E\). The order is found by calculating the order of \(E\) over \(K\) and lifting.

Trace(H): SetPtEll -> RngIntElt
Trace(E): CrvEll -> RngIntElt
TraceOfFrobenius(E): CrvEll -> RngIntElt: Returns the trace of the Frobenius endomorphism on the elliptic curve \(E\), equal to \(q+1-n\) where \(q\) is the cardinality of the coefficient field and \(n\) is the order of the group of rational points of \(E\).

Trace(H, r): SetPtEll, RngIntElt -> RngIntElt
Trace(E, r): CrvEll, RngIntElt -> RngIntElt
TraceOfFrobenius(H, r): SetPtEll, RngIntElt -> RngIntElt
TraceOfFrobenius(E, r): CrvEll, RngIntElt -> RngIntElt: Returns the trace of the \(r\)-th power Frobenius endomorphism, where \(E\) is an elliptic curve over a finite field.

Example: Group Order and Trace of Frobenius (ex-e97387)

The computation of the order of a curve over a finite field invokes the SEA algorithm. This data determines the number of points over all finite extensions, and is equivalent to the data for the trace of Frobenius. The values of the trace of Frobenius and group order over all extensions are therefore a trivial computation.

> K<w> := GF(2^133);
> E := EllipticCurve([K| 1, 0, 0, 0, w]);
> time #E;
10889035741470030830869428720869868218880
Time: 8.830
> FactoredOrder(E);
[ <2, 9>, <3, 4>, <5, 1>, <29, 1>, <1103, 1>, <3793, 1>, <16855669, 1>,
  <25678053050475964927, 1> ]
> time TraceOfFrobenius(E);
-41441283053285452287
Time: 0.000
> time TraceOfFrobenius(E, 3);
1282596408769540337607822605365889523499049843333311465783809
Time: 0.000

Run in calculator (opens in a new tab)

Example: Twists and Trace of Frobenius (ex-94b05b)

Run in calculator (opens in a new tab)

The following code demonstrates the computation of twists of an elliptic curve over a finite field, and the relationship with the trace of Frobenius.

> E := EllipticCurveFromjInvariant(GF(101^2)!0);
> Twists(E);
[
   Elliptic Curve defined by y^2 = x^3 + 1 over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (61*w + 91) over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (98*w + 16) over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (9*w + 66) over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (59*w + 76) over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (87*w + 81) over GF(101^2)
]
> IsSupersingular(E);
true

Since \(E\) is supersingular, so are the twists of \(E\). Hence the traces are divisible by the prime:

> [ TraceOfFrobenius(F) : F in Twists(E) ];
[ -202, -101, 101, 202, 101, -101 ]
> [ TraceOfFrobenius(F) : F in QuadraticTwists(E) ];
[ -202, 202 ]

We see that the traces come in pairs; the trace of a curve is the negative of the trace of its quadratic twist. This is not just true of the supersingular curves:

> E := EllipticCurveFromjInvariant(GF(101^2)!12^3);
> Twists(E);
[
   Elliptic Curve defined by y^2 = x^3 + x over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (40*w + 53)*x over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (3*w + 99)*x over GF(101^2),
   Elliptic Curve defined by y^2 = x^3 + (92*w + 19)*x over GF(101^2)
]
> IsSupersingular(E);
false
> [ TraceOfFrobenius(F) : F in Twists(E) ];
[ -198, 40, 198, -40 ]
> [ TraceOfFrobenius(F) : F in QuadraticTwists(E) ];
[ -198, 198 ]

Run in calculator (opens in a new tab)

Zeta Functions

ZetaFunction(E): CrvEll -> FldFunRatUElt: Given an elliptic curve \(E\) over a finite field, the function returns the zeta function of \(E\) as a rational function in one variable.

Note that each of the function calls Order(E), Trace(E), and ZetaFunction(E) invoke the same call to the SEA algorithm and determine equivalent data.

Example: Zeta Function of an Elliptic Curve (ex-7a0b38)

The zeta function \(\zeta_E(t)\) of an elliptic curve \(E\) over a finite field \({\mathbb{F}}_q\) is a rational function whose logarithmic derivative has the power series expansion:

\[{{d\log \zeta_E(t)}\over{d\log t}} = \sum_{n=1}^{\infty} |E(F_{q^n})|\, t^n\]

The following example illustrates this property of \(\zeta_E(t)\).

> p := 11;
> E := EllipticCurve([GF(p) | 1, 2]);
> Z := ZetaFunction(E);
> Q<t> := LaurentSeriesRing(Rationals());
> t*Derivative(Log(Evaluate(Z, t))) + O(t^7);
16*t + 128*t^2 + 1264*t^3 + 14848*t^4 + 160976*t^5 + 1769600*t^6 + O(t^7)
> [ Order(E, n) : n in [1..6] ];
[ 16, 128, 1264, 14848, 160976, 1769600 ]

Run in calculator (opens in a new tab)

The Order of the Group of Points#

Point Counting#

Zeta Functions#

The Order of the Group of Points

Point Counting

Zeta Functions