See discussions, stats, and author profiles for this publication at: https://www.researchgate.

net/publication/221451662

Elliptic Curves with a Given Number of Points

Conference Paper  in  Lecture Notes in Computer Science · June 2004

DOI: 10.1007/978-3-540-24847-7_8 · Source: DBLP

CITATIONS READS
18 169

2 authors:

Reinier Broker Peter Stevenhagen

Brown University Leiden University
19 PUBLICATIONS   242 CITATIONS    54 PUBLICATIONS   466 CITATIONS   

SEE PROFILE SEE PROFILE

All content following this page was uploaded by Peter Stevenhagen on 28 May 2014.

The user has requested enhancement of the downloaded file.

 Elliptic curves with a given number of points

Reinier Bröker and Peter Stevenhagen

Mathematisch Instituut, Universiteit Leiden,

Postbus 9512, 2300 RA Leiden, The Netherlands
reinier@math.leidenuniv.nl psh@math.leidenuniv.nl

Abstract. We present a non-archimedean method to construct, given

an integer N ≥ 1, a finite field Fq and an elliptic curve E/Fq such that
E(Fq ) has order N .

1 Introduction
A classical theorem of Hasse from 1934 states that for an elliptic curve E defined
over the finite field Fq of q elements, the order of the group E(Fq ) of Fq -rational
points is an integer in the Hasse interval
√ √
Hq = [q + 1 − 2 q, q + 1 + 2 q]

around q. If E is given in some standard way, say by a Weierstrass equation

over Fq , there are several algorithms that compute the order of E(Fq ). The
1985 algorithm by Schoof [8, 9] runs in time polynomial in log q, and in small
characteristic p there are even faster p-adic algorithms due to Satoh [7] and
Kedlaya [5].
The situation is rather different in the case of the following problem, which
can be seen as an ‘inverse problem’ to the point counting problem.
Problem. Given an integer N ≥ 1, find a finite field Fq and an elliptic curve
E/Fq for which the number of Fq -rational points equals N .
As with other inverse problems, such as in Galois theory, this is mathematically
a natural question to ask. In this particular case, an efficient solution to the
problem would also be desirable in view of the need in current applications
to construct elliptic curves having point groups satisfying various smoothness
requirements with respect to their order. It is one of the reasons why we focus
on the order N , and do not specify the finite field Fq as being part of the input.
In addition, we will use the freedom with respect to the choice of a base field Fq
to our advantage.
A necessary condition for our problem to be solvable for given N is clearly
S
that N is contained in some Hasse interval Hq , so we would like the union q Hq
over all prime powers q to contain all positive integers. It is easy to see that the
contribution to the union coming from the ‘true’ prime powers q that are not
primes is negligible: it is contained in a zero density subset of Z≥1 . For this
reason, we may and will restrict in the sequel to the case where the base field Fq
is the prime field coming from a prime number q. In this particular case, all
integers in Hq actually do occur as the group order of E(Fq ) for some elliptic
curve E, so N ∈ Hq is sufficient to guarantee the existence of a solution. (For
arbitrary prime powers q there are often not enough supersingular curves to
realize all orders congruent toS1 modulo the characteristic.)
For the equality Z≥1 = q prime Hq we need to show that the primes are
not too far apart, i.e., that the gap between consecutive primes q and q 0 is
√
roughly bounded by 4 q for large q. This is more than what is currently known
to be true: even under assumption of the Riemann hypothesis the gap between
√
consecutive primes can only be shown to be of order O( q(log q)2 ). However,
from a practical, algorithmic point of view there are always lots of primes q
for which a large integer N is contained in Hq . Indeed, by the prime number
theorem, we expect 1 out of every log N integers around N to be prime,
√ so for
large N the set of primes q having N ∈ Hq is on average of size 4 N / log N ,
and finding such q is never a problem in practice.
Once we have found a prime q > 3 for which we have N ∈ Hq (we now
require N > 1), there is the following naive algorithm to find an elliptic curve
having exactly N rational points over Fq . Suppose that we are not in the easy
cases where we have N = q + 1 (then any supersingular curve over Fq will do)
or where one of the few curves with j-invariant 0 or 1728 has the right number
of points. Then we try
4a
Ea : y 2 = x3 + ax − a with j(Ea ) = 1728
4a + 27
for random a ∈ F∗q \ {− 27 4 } as the Weierstrass equation of the desired curve until
we find a curve having N points. More precisely, we write N = q + 1 − t and
check whether for our a the point (1, 1) ∈ Ea (Fq ) is annihilated by N = q + 1 − t
or q + 1 + t. If it is, we check whether the number of Fq -rational points is indeed
q+1±t. For order N = q+1−t we are done, for order q+1+t not Ea itself but its
quadratic twist has N points. Even though the distribution of the group orders
√
of elliptic curves over Fq is not quite uniform, we expect to examine O( q) =
√
O( N ) curves Ea before we hit a curve having exactly N points. As the amount
of time spent per a is usually very small, and certainly polynomial in log N , this
1
yields a probabilistic algorithm with expected running time O(N 2 +o(1) ). It is
quite practical for small values of N , but becomes unwieldy for N  1015 .
In the next section we briefly describe a classical deterministic algorithm
based on complex multiplication methods which, although not asymptotically
faster than the naive algorithm, can be improved in various ways. Our first
improvement is a p-adic approach to complex multiplication based on the re-
cent work of Couveignes and Henocq [3]. It is described in section 3, and illus-
trated by the explicit computation in section 4 of an ‘ANTS 6 curve’ having
2004061320040618 rational points. Our second improvement, in section 5, con-
sists of using ‘small’ modular functions in this p-adic context to push the limits
of what is feasible by p-adic methods. Although the resulting algorithm is still
far from polynomial, its power is illustrated in section 6 by the computation of
an elliptic curve having 1030 rational points.
2 Complex multiplication
A deterministic algorithm to produce an elliptic curve E over the prime field Fq
having N ∈ Hq points is provided by the theory of complex multiplication. One
writes N = q + 1 − t and observes that the Frobenius endomorphism Fq : E → E
on the desired curve E/Fq satisfies the quadratic relation Fq2 − tFq + q = 0 of
discriminant ∆ = t2 − 4q < 0 in End(E) = EndFq (E). Assume we are not in the
supersingular case t = 0. Then our observation gives rise to an embedding
√
End(E) −→ K = Q( √ ∆)
Fq 7−→ πq = (t + ∆)/2

that maps Fq to√a prime element πq of trace t and norm q in the quadratic order
O∆ = Z[(∆ + ∆)/2] ⊂ K of discriminant ∆. By the Deuring lifting theorem
[6, Chapter 13, Section 5], there exist a number field H ⊃ K and an elliptic
curve Ẽ/H such that
1. there exists φq ∈ End(Ẽ) satisfying φ2q − tφq + q = 0;
2. the prime q splits completely in H/Q, and for every prime q|q in H the
reduced curve Ẽ mod q is an elliptic curve having q + 1 − t points.
In fact, the reduction of the endomorphism φq ∈ End(Ẽ) above modulo a prime
q|q in H yields the Frobenius endomorphism of the curve Ẽ mod q over Fq .
The smallest field H ⊃ K over which a curve Ẽ satisfying 1 and 2 can
be defined is the Hilbert class field of K. More explicitly, from a list of reduced
binary quadratic forms [a, b, c] of discriminant D = disc(K), which can be viewed
as a list of elements of the class group Cl(D) of discriminant D, we can form the
class polynomial
!
Y  −b + √D 
FD = X −j ∈ Z[X]
2a
[a,b,c]∈Cl(D)

of discriminant D. Here j : H → C denotes the well-known elliptic modular

function on the complex upper half√plane. Using sufficiently accurate complex
D
approximations of the zeroes j( −b+ 2a ) of FD , one may exactly determine FD
as it has integral coefficients. Once we have FD , we are essentially done. Indeed,
any zero of the irreducible polynomial FD ∈ Z[X] generates the Hilbert class
field H of K over K, and modulo q the polynomial F̄D ∈ Fq [X] splits into linear
factors. (In fact, this property of FD is an excellent check for the correctness of
any algorithm to compute FD .) The zeroes of F̄D in Fq are the j-invariants of
the elliptic curves over Fq having endomorphism ring isomorphic to the ring of
integers OD of K. If j ∈ Fq is one of these zeroes, we write down a curve with
this j-invariant, and check whether it has q + 1 − t points. If it hasn’t, we have
found a curve with q + 1 + t points, and (for j 6= 0, 1728) its quadratic twist has
N = q + 1 − t points.
This deterministic algorithm, although relatively simple, is not much faster
than the naive algorithm, as for large D the class polynomial FD is of degree
 1 1
h(D) = O(D 2 +o(1) ) and has coefficients of size O(D 2 +o(1) ), making the total
running time of order O(D 1+o(1) ). However, if we have the freedom to pick q on
input N , it is often possible to find primes q for which the associated discriminant
1
∆ = ∆(q) = (q − 1 − N√)2 − 4N is not√only of size N 2 +o(1) (this is simply done by
picking q ∈ [N + 1 − 2 N , N + 1 + 2 N] close to the end points of the interval),
but moreover has a large square factor. This leads to a field discriminant D of K
which is quite a bit smaller than ∆, and makes the method feasible in cases
where the naive method would take too long. As√we currently √ cannot even prove
the existence of a single prime in q ∈ [N + 1 − 2 N , N + 1 + 2 N ], we certainly
cannot prove the existence of q for which ∆(q) has large square factors and D
1
is of order substantially smaller than O(N 2 +o(1) ).

3 A non-archimedean approach
The key feature of the complex multiplication method in the previous section
is the computation of the class polynomial FD ∈ Z[X] of the order OD for
suitable D. As the zeroes modulo q of FD are j-invariants of curves E/Fq for
which either E or its quadratic twist has exactly N points, this immediately
solves our problem. In this section we achieve the computation of FD in an
other way, using p-adic instead of complex approximations of the zeroes of F D .
Working in a non-archimedean setting has the advantage that we no longer have
to cope with the problem of rounding errors that arises in the complex case. It
does require a p-adic substitute for the complex analytic method to evaluate j
in CM-points of H using Fourier expansions, and this is provided by the recent
work of Couveignes and Henocq [3] explained in this section.
2
Let N = q + 1 − t and √ ∆ = t − 4q be as in the previous section, and D < −4
the discriminant of Q( ∆). We first construct an elliptic curve E over a finite
field Fp which has CM with OD . As we want p to be as small as possible, we let s
be the smallest positive integer of the same parity as D for which p = (s2 − D)/4
is prime. For D ≡ 1 mod 8 such p does not exist for parity reasons unless D = −7,
and we pick the smallest positive s for which p = (s2 − 4D)/4 is prime instead.
In practice we expect s to be small, at most a power of log |D|, so that p is of the
same order of magnitude as D. Unfortunately, even under GRH proven upper
bounds [3] for s are much weaker. √ √
As there is a prime element πp = (s + D)/2 (or πp = (s + 2 D)/2) of
norm p and trace s > 0 in the order OD = Z[πp ] (or O4D = Z[πp ]), there exists
an ordinary elliptic curve over Fp having CM by Z[πp ] and p + 1 ± s points
over Fp . We can find such a curve E/Fp by applying the naive algorithm, as we
saw in section 2 that |D| is much smaller than N . We have End(E) = OD for
D 6≡ 1 mod 8. For D ≡ 1 mod 8 the ring End(E) is either equal to Z[πp ] = O4D
or to OD ) Z[πp ]. We are in the second case if all 2-torsion of E is Fp -rational,
and in the first if it isn’t. As for D ≡ 1 mod 8 the class polynomials FD and F4D
have the same degree and both generate the Hilbert class field H of K, they are
both fine for our purposes. Alternatively, in case we have End(E) = O4D there
is a unique point P of order 2 in E(Fp ), and we can replace E by the 2-isogenous
curve E/hP i to achieve End(E) = OD . We assume that we have End(E) = OD
in the sequel. Note that D is by definition fundamental.
By the Deuring lifting theorem, there exists a curve Ẽ/H and a prime p|p
such that Ẽ reduces modulo p to E and such that we have End(Ẽ) ∼ = End(E).
As p splits completely in H/Q the curve Ẽ is actually defined over Qp . In fact,
Ẽ is the unique elliptic curve over Qp with reduction E and endomorphism
ring End(Ẽ) ∼ = End(E) ∼ = OD . It is this canonical lift Ẽ of E that we want to
compute, as its j-invariant is a zero of our class polynomial.
Let EllD (Fp ) be the set of Fp -isomorphism classes of elliptic curves over Fp
with endomorphism ring OD . The j-invariants of the elements in EllD (Fp ) are
the zeroes of (FD mod p), so EllD (Fp ) is finite of order h(D) = #Cl(D). It can
be identified with the similarly defined set EllD (Qp ) of Qp -isomorphism classes
of elliptic curves over Qp with endomorphism ring OD . The j-invariants of the
elements in EllD (Qp ) are the zeroes of FD . Let Cp be the completion of an
algebraic closure of Qp , and write XD (Cp ) for the set of isomorphism classes of
elliptic curves over Cp with the property that their reduction is an element of
EllD (Fp ). Then XD (Cp ) is a p-adic analytic space, as the j-invariant identifies
it with a subset of Cp . It consists of h(D) discs of p-adic radius 1. Each disc
contains exactly one element of EllD (Qp ), and this is the subset of Cp we want
to compute. It consists of the (j-invariants of the) isomorphism classes of elliptic
curves in XD (Cp ) having CM with OD .
Let I ⊂ OD be an OD -ideal prime to p and E/Fp an elliptic curve in
EllD (Fp ). Then there is a separable isogeny E → EI which has the subgroup
E[I] of I-torsion points of E as its kernel. In this way, we obtain a bijection
ρ̄I : EllD (Fp ) → EllD (Fp ) that sends the isomorphism class of E to that of EI .
We obtain an action of the group I(p) of OK -ideals prime to p on EllD (Fp ), and
since principal OK -ideals act trivially, this action factors via the quotient map
I(p)  Cl(D). This makes EllD (Fp ) into a principal homogeneous Cl(D)-space.
The fundamental idea in [3] is that the action of I(p) on EllD (Fp ) admits
a natural lift to an action of I(p) on XD (Cp ). More precisely, for I ∈ I(p) the
map ρI : XD (Cp ) → XD (Cp ) is a p-adic analytic map that lifts ρ̄I in the sense
that on EllD (Qp ) ⊂ XD (Cp ), the restriction of ρI is the standard Galois action
that factors via Cl(D). If I = (α) is principal with α ∈ OD \ Z, then ρα = ρI
stabilizes each disk around a CM-point in EllD (Qp ), and has the CM-point in
the disk as its unique fixed point.
It is shown in [3] that the derivative of the map ρα in a point of EllD (Qp )
equals α/ᾱ, and this can be used to compute the j-invariant of a CM-point in
EllD (Qp ) starting from an arbitrary point in the disk using a Newton iteration
process. If E1 is any lift of E to Cp , we put

j(ρα (Ek )) − j(Ek )

(1) j(Ek+1 ) = j(Ek ) − for k ∈ Z≥1 .
(α/ᾱ) − 1

If α is small with respect to p, then (α/ᾱ) − 1 is a unit in Zp and the sequence

(1) converges to the j-invariant of the canonical lift Ẽ. In each step, the p-adic
precision of the approximation is doubled.
XD (Cp ): e e e e

EllD (Fp ):

For the definition of ρI on the isomorphism class of E 0 in XD (Cp ), we note

that the subgroup E[I] of I-torsion points of the reduction E/Fp of E 0 lifts
canonically to a subgroup E 0 [I] of E 0 , and we put ρI (E 0 ) = EI0 , with E → EI0
the isogeny with kernel E 0 [I]. This provides a lift of E[I] to a group scheme over
the p-adic disk in XD (Cp ) lying over [E] ∈ EllD (Fp ). More algorithmically, if
E is given by a Weierstrass model, the subgroup E[I] can be described by a
separable polynomial f¯I ∈ Fp [X] having the x-coordinates of the affine points
in E[I] as its zeroes. If I has norm n, then p - n and f¯I divides the n-th division
polynomial of E in Fp [X]. Choosing a Weierstrass model for E 0 reducing to that
for E, we can lift f¯I uniquely by Hensel’s lemma to a factor fI of the n-th division
polynomial of E 0 , and this factor describes the subgroup E 0 [I]. Note that ρI is
the identity on XD (Cp ) if I is generated by an integer, and that, consequently,
ρI¯ is the inverse of ρI .
For the explicit computation of ρα (E 0 ), we take an element α = a + bπp in
OD \ Z that is sufficiently smooth, i.e., a product of OD -ideals L = (`, a + bπp )
of small prime norm ` 6= p. Such an element is found by sieving in the set

{a + bπp : a, b ∈ Z, b 6= 0, (a, b) = 1}.

Again, the smoothness properties are in practice much better than what can be
rigorously proved [3]. For the computation of ρL (E 0 ), we first compute the action
of ρ̄L on the reduction E = Ē 0 ∈ EllD (Fp ), i.e., the j-invariant of EL = ĒL0 .
The kernel E[L] of the isogeny E → EL is a cyclic subgroup of order ` of E[`]
that is an eigenspace of the Frobenius morphism with eigenvalue −b/a ∈ F` .
The corresponding polynomial f¯L ∈ Fp [X] can be computed by the techniques
used by Atkin and Elkies to improve Schoof’s original point counting algorithm,
see [9]. These techniques also yield a Weierstrass model for EL ; we only need
the j-invariant j(EL ) ∈ Fp of EL = ρ̄L (E). From the decomposition of ρα into
‘prime degree’ maps ρL , we obtain a cycle of isogenies
1 ρ̄L 2 ρ̄L t ρ̄L
(2) E −→ EL1 −→ EL1 L2 −→ . . . −→ E(α) = E.

To compute the action of L on the lift E 0 of E, we do not lift f¯L to some precision
to a divisor of the `-th division polynomial of E 0 , and then find a Weierstrass
model and the j-invariant for EL using Vélu’s formulas [11]. As the division
polynomial has degree (`2 −1)/2 for ` > 2, this would be rather time-consuming.
Instead, we exploit the `-th modular polynomial Φ` (X, Y ) ∈ Z[X, Y ], which is
of degree ` + 1 in each of the variables, and which we can precompute for a
number of small primes `. As there is an isogeny E 0 → EL0 of degree `, we know
that j(EL0 ) is a zero of Φ` (j(E 0 ), Y ) ∈ Zp [Y ]. Since we know the j-invariant of
the reduction of EL0 , we also know which root to approximate in Zp , and this
reduces the lifting process to a simple Hensel lift of a zero of a polynomial of
degree ` + 1.
For our Newton process in (1), we start from E ∈ EllD (Fp ), and compute
the cycle of Fp -isogenies in (2). We then lift E arbitrarily to a curve E1 over Qp ,
the ‘1 digit precision’ approximation of the canonical lift. Now we compute lifts
over Qp of our Fp -isogenies in 2 digit precision, using the modular polynomials,
and use the value of ρα (E1 ) obtained to update E1 as in (1) to a 2 digit preci-
sion approximation E2 of the canonical lift. We continue this process of making
(not really closed) cycles over Qp , doubling the precision of the computation at
each step, until we have the canonical lift with high enough accuracy (see [2,
Chapter 7, Section 6] for an estimate of the required accuracy).
If we know the canonical lift Ẽ, we compute its Galois conjugates again via
the modular polynomials. For this, we need small primes L that generate the
class group. Under GRH, this can be done using primes L not exceeding the
Bach bound 6 log2 (|D|), see [2, Chapter 5, Section 5]. In practice, this is never
a problem. If we have all the conjugates to the required precision, we find by
simple expansion of the product below the class polynomial
Y
FD = (X − j(ẼI )) ∈ Z[X].
[I]∈Cl(D)

4 An elliptic curve for ANTS 6

We illustrate the working of our p-adic method by computing a tailor made
elliptic curve for ANTS 6 having exactly N = 2004061320040618 points.
First we look for a small discriminant, so we write N = q + 1 − t for various
primes q and search for a large square dividing ∆ = t2 − 4q. In this example, the
choice
N = 2004061230508291 + 1 + 89532326
yields ∆ = −23 · 3 · 6192 · 22567, so we take D = −23 · 3 · 22567 = −541608 in
this section. The corresponding class group Cl(D) has order 132.
Our goal is to compute the class polynomial FD . We will do this p-adically, so
we first find an elliptic curve over some Fp which has CM with OD . The smallest
integer s > 0 for which (s2 − D)/4 is prime is s = 2, so we have D = 22 − 4p
with p = 135403. We fix this value of p for the rest of this section.
We now apply the naive method to find a curve Ea : y 2 = x3 + ax − a over
Fp with trace of Frobenius 2. We find that E1737 has trace 2 and take this as
our base curve E/Fp .
 Next we determine which element α ∈ OD \ Z we will use for the map ρα .
The ideal (α) = (22539 + 4πp ) of norm 510353281 = 192 · 292 · 412 factors as
(α) = L219 · L229 · L241 = (19, πp + 6)2 · (29, πp − 13)2 · (41, πp − 13)2 .
We compute the action for the prime ideal L19 . The eigenvalue of the action of
Frobenius on the 19-torsion is −6 ∈ F19 . If we evaluate the modular polynomial
Φ19 (X, Y ) in X = j(E) = 41556 ∈ Fp , we get a polynomial which has two roots
over Fp , namely 19533 and 54827. From this we deduce that L19 sends j(E) to
one of these two roots; we don’t know which one yet.
We just guess that the correct j-invariant is 54827 ∈ Fp . Following Elkies
[9], we now compute the eigenspace S of the 19-torsion corresponding to this
isogeny. We get the x-coordinates of the points on E in S as zeroes of
f¯L19 = X 9 +29873X 8 + 49874X 7 + 131130X 6 + 49222X 5 + 46538X 4+
111513X 3 + 68602X 2 + 126444X + 20947 ∈ Fp [X]
Since we know that the eigenvalue for L19 is −6, we can now just check whether
(X p , Y p ) = −6 · (X, Y )
holds for points in S, i.e., we compute both (X p , Y p ) and −6 · (X, Y ) in the ring
Fp [X, Y ]/(f¯L19 (X), Y 2 − X 3 + 1737X − 1737).
Note that the · means adding on the curve! In this example, it turns out that
(X p , Y p ) and −6 · (X, Y ) are not the same. It follows that the correct j-invariant
of the L19 -isogenous curve is the other value 19533 ∈ Fp .
The action of L19 on the curve with j-invariant 19533 is now easier to com-
pute: the modular polynomial has again two roots, but one of the roots has to
be j(E). This root corresponds to the action of L̄19 , so we pick the other root.
If we compute the entire cycle corresponding to L, we get:
L
19 19 L 29 L 29 L 41 L
41 L
41556 −→ 19533 −→ 100121 −→ 86491 −→ 40349 −→ 32517 −→ 41556.
We now lift E/Fp to E1 /Qp by lifting the coefficients of the Weierstrass equation
arbitrarily. The polynomial Φ19 (j(E1 ), Y ) ∈ Zp [Y ] has exactly two roots, one of
which reduces to 19533 modulo p. We compute this root, which is the value of
the L19 -action on j(E1 ), to two p-adic digits of precision. We continue to lift
the whole cycle to two p-adic digits of precision, and update j(E1 ) according to
formula (1) to obtain j(E2 ). Starting from j(E2 ), we now lift the cycle to four
p-adic digits of precision, compute j(E3 ) from this, and so on. We obtain
j(Ẽ) = 41556 + O(p)
= 41556 − 17953p + O(p2 )
= 41556 − 17953p − 51143p2 − 17793p3 + O(p4 )
= 41556 − 17953p − 51143p2 − 17793p3 + 45123p4 + 52596p5 + 18237p6
+42211p7 + O(p8 )
= 41556 − 17953p − 51143p2 − 17793p3 + 45123p4 + 52596p5 + 18237p6
+42211p7 + 45716p8 + 58788p9 + 18836p10 − 4101p11 − 60004p12
−24668p13 + 27527p14 − 58942p15 + O(p16 ).
From the estimate in [2, Chapter 7, Section 6], we know that we need approxi-
mately 3000 decimals of accuracy. As we have 103000 ≈ p585 , we compute j(Ẽ)
up to 585 p-adic digits. The class group Cl(D) ∼ = Z/2Z × Z/66Z is generated
by a prime of norm 2 and a prime of norm 29, so we find all 132 conjugates of
j(Ẽ) under Gal(H/K) up to 585 p-adic digits using the modular polynomials Φ2
and Φ29 . In the end, we expand the polynomial of degree 132 to find the class
polynomial Y
FD = (X − j(ẼI )) ∈ Z[X]
[I]∈Cl(D)

We now compute a root of FD ∈ Fq , and note that we can check our computa-
tions so far by testing whether FD splits completely in Fq [X]. One of the 132
roots is j = 5215470850369 ∈ Fq , and an elliptic curve with this j-invariant is
Ea : y 2 = x3 + ax − a with a = 27j/(4(1728 − j)) = 1460967812073632 ∈ F q .
We know that Ea has CM with OD and that its trace of Frobenius equals ±t.
To test whether it actually equals t, we look at the order of P = (1, 1) ∈ E a (Fq ).
We see that (q + 1 + t)P 6= O = (q + 1 − t)P , so Ea must have trace equal to t.
We conclude that the curve defined by

y 2 = x3 + 1460967812073632x + 543093418434659

over F2004061230508291 has exactly 2004061320040618 rational points.

5 Using class invariants

A serious drawback of the complex multiplication method is that it requires the

computation of the class polynomial FD , which grows rapidly in size with D,
and is already sizable for moderately small discriminants. Around 1900, Weber
[12] computed generating polynomials for Hilbert class fields using the values at
CM-points of modular functions of higher level instead of the j-function. These
techniques have become important again in an algorithmic context, and many
complications from Weber’s days are now well understood [4, 10]. The classical
theory of class invariants is firmly rooted in complex analytic arguments, but
much of it can be made to work in our non-archimedean setting. This section
gives an indication of what can be done, leaving a fuller treatment to [1].
The complex multiplication method to compute the class polynomial FD
is based on the fact that its zeroes are the j-invariants of the elliptic curves
in characteristic zero having complex multiplication with OD . In the complex
analytic setting, we simply list the complex lattices (up to scaling) giving rise
to such curves and compute their j-invariants to sufficient accuracy using the
q-expansion of j. In the p-adic setting, we first compute one such curve in the
finite set EllD (Fp ). Then we lift the action on EllD (Fp ) of the group I(p) of
OD -ideals prime to p to an action on the set X(Cp ) of their Cp -lifts in such
a way that the induced action on the subset EllD (Qp ) of canonical lifts is the
standard Galois action coming from Cl(D). This enables us to compute the finite
j
subset EllD (Qp ) ⊂ X(Cp ) −→ Cp in Cp consisting of the zeroes of FD . We use a
Newton type lift of our curve in EllD (Fp ) to EllD (Qp ) together with the explicit
Galois action of Cl(D) on EllD (Qp ), which we handle by exploiting modular
polynomials.
We now want to replace j by modular S functions of higher level. These are
elements of the modular function field F = n≥1 Fn as defined in [6, Chapter 6,
Section 3]. Here Fn denotes the field of modular functions of level n over Q. It can
be viewed as the function field of the modular curve X(n) over the cyclotomic
field Q(ζn ). Over F1 = Q(j), it is generated by the Fricke functions of level n,
which are normalized x-coordinates of n-torsion points on an elliptic curve with
j-invariant j.
The modular functions f we use are integral over Z[j], so they are given as
the zero of some irreducible polynomial Ψf (X, j) ∈ Z[j, X]. If we specialize j to
be the j-invariant of a curve Ẽ ∈ EllD (Qp ), then the roots of the polynomial
Ψf (X, j(Ẽ)) ∈ H[X], which has √ integral coefficients, lie in the ray class field
Hn of conductor n of K = Q( D), with n the level of f . It is known that for
many choices of ‘small’ f , one or more of these roots are class invariants that
actually lie in the Hilbert class field H = H1 ⊂ Hn of K. If we can determine
which roots end up in H, and compute the explicit Galois action of Cl(D) on
such a root, we can compute its irreducible polynomial over K or Q just like we
did this for j. In the complex analytic setting the tool to perform these tasks is
Shimura’s reciprocity law [4, 10]. It tells us in which points τ ∈ K ⊂ H a function
f should be evaluated to obtain a class invariant, and describes the conjugates
of f (τ ) over K as the values of conjugates of f over Q(j) in certain other points
τ 0 ∈ K. These values can be approximated in C using the q-expansions of f
and its conjugates. Once we have computed the irreducible polynomial of f (τ )
over K or Q, one can use the relation Ψf (f (τ ), j(τ )) = 0 to obtain information
on j(τ ) itself.
In a p-adic setting, we cannot deal with j and f as functions on the complex
upper half plane, and the expansion of modular functions as Fourier series in
q = e2πiτ has no non-archimedean analogue when dealing as we do with CM-
curves, which have integral j-invariants. What we do have is an action from class
field theory of the OD -ideals coprime to n on the roots of Ψf (X, j(Ẽ)) ∈ H[X]
in Hn associated to j(Ẽ). This action factors via the class group Cl(n2 D) of the
order of discriminant n2 D, which is a non-maximal order for n > 1.

Example 1. The modular function j : H → C has a holomorphic cube root

γ2 : H → C that is modular of level 3 and has Ψγ2 (X, j) = X 3 − j. It is the
unique root of X 3 − j having a rational q-expansion. If D is not divisible by 3
and we write OD = Z[τ ] with τ + τ̄ = 0 mod 3, then γ2 (τ ) is a class invariant,
and its ‘size’ is only one third of that of j(τ ).
If Ẽ : y 2 = x3 + ax + b is in EllD (Qp ) and c1 , . . . , c4 are the 4 roots of its
3-division polynomial, then

−48a
(3)
2a − 3(c1 c2 + c3 c4 )
is a cube root of j(Ẽ). As there are 3 ways to divide the roots c1 , . . . , c4 in two
sets of two roots each, formula (3) yields 3 distinct cube roots of j. Note that
there is no obvious way to single out an expression in (3) as being ‘the function’
corresponding to γ2 .
If I is an OD -ideal prime to 3p, the isogeny Ẽ → ẼI induces an isomorphism
∼
Ẽ[3] −→ ẼI [3] on the 3-torsion subgroups, and maps each possible cube root of
j(Ẽ) in (3) to some well-defined cube root of j(ẼI ). This is the Galois action
of the Artin symbol of I, which maps j(Ẽ) to j(ẼI ), on these cube roots. It
provides an extension of the map ρI : EllD (Qp ) → EllD (Qp ) to the set of
the cube roots of (the j-invariants of) the elements in EllD (Qp ). As all these
cube roots are p-adically integral, reduction of this map provides an extension of
ρ̄I : EllD (Fp ) → EllD (Fp ) to the set of cube roots in F̄p , provided we have p 6= 3.
In a similar way, we have an extension of the map ρI : XD (Cp ) → XD (Cp ) to
the set of cube roots.
Example 1 nicely illustrates that the cube roots of j are functions on the modular
curve X(3), the points of which can be viewed as isomorphism classes of elliptic
curves with complete 3-level structure: in order to have a well defined value not
only the j-invariant of the curve but also some ‘ordering’ of 3-torsion points is
required.
As the field Fn is generated over Q(j) by Fricke functions, a modular function
f of level n is always a rational expression in j and the roots of the n-th division
polynomial of a curve with j-invariant j. As in Example 1, the action ρI of
an OD -ideal I coprime to pn on XD (Cp ) naturally maps the roots of Ψf (X, j)
for j ∈ XD (Cp ) to the roots of Ψf (X, ρI (j)). This observation suffices to treat
modular functions providing class invariants by p-adic methods similar to that
in section 3.
Suppose f is an integral modular function of level n that is known to pro-
vide class invariants for OD at certain CM-points for OD in H. Then we know
that for every curve Ẽ ∈ EllD (Qp ), certain roots of Ψf (X, j(Ẽ)) ∈ Qp [X] are
class invariants, so they lie in Qp . If E ∈ EllD (Fp ) is the reduction of Ẽ, then
Ψf (X, j(E)) ∈ Fp [X] will have roots in Fp , and we want to know which roots
in Fp arise as the reduction of class invariants. Let β be a root, and assume for
simplicity that Ψf (X, j(E)) is separable. (This is usually the case in practice as
p is not too small, of size O(D 1+ε ).) The roots of Ψf (X, j(Ẽ)) ∈ Qp [X] all lie
in Hn , and they are class invariants exactly when they are fixed by the Galois
group
Gal(Hn /H) ∼ = (OD /nOD )∗ /OD ∗
= ker[Cl(n2 D) → Cl(D)].
It follows that β arises as the reduction of a class invariant if it is fixed by the
maps ρ̄x : EllD (Fp ) → EllD (Fp ) for a set of generators x of (OD /nOD )∗ /OD ∗
.
We can compute ρ̄x (j(E)) as before when x is sufficiently smooth, and we need
the extension of the action to the roots of Ψf (X, j(E)). In theory this can be
done by working with explicit Weierstrass models and an explicit description of
f in terms of n-torsion points as in Example 1. In practice we work with modular
polynomials for prime degree ` isogenies relating the roots of Ψf (X, j1 ) to that
of Ψf (X, j2 ) when j1 and j2 are j-invariants of `-isogenous curves. Such modular
polynomials are in many ways similar to the modular polynomials Φ` (X, Y ) aris-
ing for the j-function, and they may be found using complex analytic methods
involving q-expansions.
Suppose we find that a certain root β ∈ Fp of Ψf (X, j(E)) is the reduction
of a class invariant β̃. Then we lift the pair (E, β) to the pair (Ẽ, β̃) consisting
of the canonical lift Ẽ/Qp and the class invariant β̃. In theory this can be
done by lifting j(E) as in section 3, using cycles of smooth isogenies, and then
compute the Hensel lift of β to the root β̃. This method makes use of the modular
polynomials Φ` for j, which are big. In many cases the corresponding polynomials
Φf` for f are much more pleasant to work with, so it is better to lift isogeny
cycles not in terms of j-invariants but in terms of a root of Ψf (X, j), and find
the resulting j-invariant from its corresponding root.
All that remains is computing the conjugates of the pairs (Ẽ, β̃) under Cl(D).
This is done in exactly the same manner as before, using small primes that
generate Cl(D). We are only interested in the conjugates of β̃, so we use the
modular polynomials Φf` to compute these conjugates. In the end we expand the
polynomial Y
f
FD = (X − β̃ I ) ∈ Z[X] (or OD [X].)
[I]∈Cl(D)

This polynomial splits again completely in Fq [X], and from a root in Fq we

compute the corresponding j-value in Fq .
Example 2. We take for f the Weber function f, which is classically defined on H
−1
in terms of the Dedekind η-function as f(τ ) = ζ48 η( τ +1
2 )/η(τ ). It is a modular
function of level 48 of degree 72 over Q(j) with
Ψf (X, j) = (X 24 − 16)3 − jX 24 ∈ Z[j, X].
For discriminants D ≡ 1 mod 8 with 3 - D, it yields class invariants when eval-
uated at appropriate points τ ∈ H. In other cases small powers of f often have
the same property [4].
In principle we can express f in terms of x-coordinates of 48-torsion points,
but there is no need to do this. It suffices to find, for E/Fp a curve having CM
with OD , first a root β ∈ Fp of Ψf (X, j(E)) that is the reduction of a class
invariant β̃, then a good approximation of the root β̃ of Ψf (X, j(Ẽ)) in Qp ,
and finally good approximations in Qp of the conjugates of β̃ under Cl(D).
These questions ultimately reduce to computing the action of an ideal L of
prime norm ` - pn on pairs (E 0 , β 0 ), with E 0 in XD (Cp ) or EllD (Fp ) and β 0
a root of Ψf (X, j(E 0 )) in Qp or Fp . For E 0 we know how to compute j(ẼL0 ),
for (β 0 )L we use the fact that it is a zero of the modular polynomial Φf` (β 0 , X)
and of Ψf (X, j(EL0 )). Usually there are only two roots Φf` (β 0 , X) that we need
to consider, the correct one and the image under the action of L̄. It is of great
help that the modular polynomials Φf` are quite a bit smaller than the classical
modular polynomials for j. For small ` they are really small, like
Φf5 (X, Y ) = (X 5 − Y )(X − Y 5 ) + 5XY
Φf7 (X, Y ) = (X 7 − Y )(X − Y 7 ) + 7(XY − X 4 Y 4 ).
For ` = 13 it takes at least two of these pages to write down Φ` , but we have

Φf13 (X, Y ) = (X 13 − Y )(X − Y 13 ) + 5 · 13XY

+13(X 2Y 12 + X 12 Y 2 + 4X 10 Y 4 + 4X 10 Y 4 + 6X 6 Y 8 + 6X 8 Y 6 ).

6 An elliptic curve having 1030 points.

We illustrate the power of the method presented in the previous section by
constructing an elliptic curve having exactly N = 1030 rational points.
Just as in section 4, we first look for a suitable
√ discriminant. We write N =
q + 1 − t, and by looking at |t| slightly less than 2 N = 2 · 1015 , we find for trace
t = 1999999999167682 that the number q = 1030 + t − 1 is prime and that

∆ = t2 − 4q = −212 · 32 · 52 · 172 · 3672 · 3943 · 23537

has a large square factor leading to D = −92806391. The corresponding class

group Cl(D) has order 15610. Computing the Hilbert class polynomial for our D
would require an accuracy of 313618 decimals, which is clearly not practical.
Instead, we notice that we have D ≡ 1 mod 8 and D 6≡ 0 mod 3, so we can
use the classical Weber function f to compute a class polynomial for H(K).
We first compute an elliptic curve in EllD (Fp ) for a ‘small’ prime p. We have
D ≡ 1 mod 8, and the smallest s ∈ Z>0 with p = (s2 − 4D)/4 prime is s = 132,
leading to p = 92810747. The first curve Ea /Fp of trace 132 we encounter is
E1086 : y 2 = x3 + 1086x − 1086 of j-invariant 37202456. As E has all three of its
two-torsion points defined over Fp , its endomorphism ring is OD , not O4D .
We now have to determine which root β of the polynomial Ψf (X, j(E)) =
(X 24 − 16)3 − j(E)X 24 ∈ Fp [X] is the reduction of a class invariant β̃ ∈ Zp . We
are lucky since ±21677132 are the only two roots in Fp . Since −β̃ is also a class
invariant, it does not matter which root we pick. We take β = 21677132.
For the smooth ideal inducing ρα we pick (α) = (−420+πp), which factors as

(11, 8 + 2πp ) · (17, 4 + 2πp )2 · (23, 16 + 2πp )2 · (31, 13 + 2πp ) · (41, 30 + 2πp ).

Just as in section 4, we compute the cycle in Fp for the j-invariants:

11L 17 L 31 L 41 L
37202456 −→ 4967239 −→ · · · −→ 21402782 −→ 37202456.

Using this cycle, we can also compute the cycle for β. For instance, the modular
polynomial Φf11 (β, Y ) ∈ Fp [Y ] has two roots: 32604444 and 60476019. In order
to determine which root to take, we note that β is a root of Ψf (X, j(EL11 )) =
(X 24 − 16)3 − j(EL11 )X 24 ∈ Fp [X]. We find that 60476019 is the root we need.
Continuing like this, we get the following cycle for β in Fp :
11L 17 L 31 L 41 L
21677132 −→ 60476019 −→ · · · −→ 53004472 −→ 21677132.

Just as in section 4, we lift E/Fp to E1 /Qp by lifting the coefficients of its

Weierstrass equation. We could now compute the canonical lift Ẽ in two p-adic
digits just like we did in section 4, and use that information to compute β̃ ∈ Zp in
two p-adic digits accuracy. Indeed, β̃ is a root of Ψf (X, j(Ẽ)) ∈ Zp [X], and since
we have (β̃ mod p) = β, this root is just a Hensel lift of β. This approach has
the disadvantage that we have to use modular polynomials for the j-function.
Instead, we lift β to a root β1 of Ψf (X, j(E1 )) ∈ Zp [X]. Now we lift the cycle
that we had for β ∈ Fp to a cycle for β1 ∈ Zp by applying the small modular
(α) (α)
polynomials for f once more. Since we know that β1 is a root of Ψf (X, j(E1 )),
we can compute
(α)
(α) ((β1 )24 − 16)3
j(E1 ) = (α)
(β1 )24
and use this value to update j(E1 ) as in formula (1) to a value of the j-invariant of
the canonical lift Ẽ that is accurate to two p-adic digits. Knowing j(Ẽ) mod p2 ,
we can lift β ∈ Fp to a root of Ψf (X, j(Ẽ)) in Zp that is accurate to two p-adic
digits. We continue this process of doubling the precision until we have β̃ with
sufficient accuracy. The first four cycles yield:

β̃ = 21677132 + O(p)
= 21677132 + 28966941p + O(p2 )
= 21677132 + 28966941p + 7010373p2 + 31182954p3 + O(p4 )
= 21677132 + 28966941p + 7010373p2 + 31182954p3 − 33808617p4
+27519307p5 − 31601027p6 − 36195013p7 + O(p8 )
= 21677132 + 28966941p + 7010373p2 + 31182954p3 − 33808617p4
+27519307p5 − 31601027p6 − 36195013p7 − 8331811p8
−33957007p9 − 18191700p10 + 5895954p11 − 42670221p12
+23637278p13 − 40784695p14 + 7754196p15 + O(p16 ).

We expect to need d313618/72e = 4356 decimals of accuracy, so we compute β̃

upto 550 p-adic digits. The class group Cl(D), which is cyclic of order 15610, is
generated by a prime of norm 11. We can thus compute all the conjugates of β̃
under Gal(H/K) to 550 p-adic digits using the modular polynomial Φf11 . In the
end, we expand the polynomial of degree 15610 to find the class polynomial
f
Y
FD = (X − β̃ I ) ∈ Z[X],
[I]∈Cl(D)

which we reduce modulo q to compute a root γ ∈ Fq . From γ, we compute the

corresponding j-value and write down a curve Ē with that j-invariant. We then
know that Ē has CM with OD , but we still need to check whether the trace really
equals t. This turns out not to be the case, so we conclude that the quadratic
twist of Ē, given by

y 2 = x3 + 669397215131271955483581235905x + 363369366443977510319399421188

over Fq , with q = 1030 + 1999999999167681, has exactly 1030 rational points.

Checking that this curve indeed has the required number of points is an easy
matter for the current point counting algorithms.
 References
[1] Bröker, R.: Thesis, Universiteit Leiden, in preparation.
[2] Cohen, H.: A course in computational algebraic number theory, Graduate Texts in
Mathematics 138, (1993).
[3] Couveignes, J.-M. and Henocq, T.: Action of modular correspondences around CM
points in Algorithmic Number Theory Symposium V, Lecture Notes in Computer
Science 2369, (2002), 234–243.
[4] Gee, A.: Class invariants by Shimura’s reciprocity law , Journal de Théorie des
Nombres de Bordeaux 11, (1999), 45–72.
[5] Kedlaya, K.: Counting Points on Hyperelliptic Curves using Monsky-Washnitzer
Cohomology, Journal Ramanujan Mathematical Society 16, (2002), 323–338.
[6] Lang, S.: Elliptic functions, 2nd edition Graduate Texts in Mathematics 112,
(1987).
[7] Satoh, T.: The canonical lift of an ordinary elliptic curve over a finite field and its
point counting, Journal Ramanujan Mathematical Society 15, (2000), 247-270.
[8] Schoof, R.: Elliptic Curves over Finite Fields and the Computation of Square Roots
mod p, Mathematics of Computation 44, (1985), 483–494.
[9] Schoof, R.: Counting points on elliptic curves over finite fields, Journal de Théorie
des Nombres de Bordeaux 7, (1995), 219–254.
[10] Stevenhagen, P.: Hilbert’s 12th problem, complex multiplication and Shimura reci-
procity in Class field theory—its centenary and prospect (Tokyo, 1998) Adv. Stud.
Pure Math. 30, Math. Soc. Japan, (2001), 161–176.
[11] Vélu, J.: Isogénies entre courbes elliptiques, C.R. Math. Acad. Sc. Paris 273,
(1971), 238–241.
[12] Weber, H.: Lehrbuch der Algebra, Vol. III. Chelsea reprint, original edition, 1908.

View publication stats