Draft.

PROVING PRIMALITY IN ESSENTIALLY QUARTIC RANDOM TIME

DANIEL J. BERNSTEIN

Abstract. This paper presents an algorithm that, given a prime n, finds and verifies a proof of the primality of n in random time (lg n)^{4+o(1)}.

1. Introduction

This paper presents an algorithm that proves the primality of any prime n in random time (lg n)^{4+o(1)}:

• Section 3 defines certificates, and proves that n is a prime power if it has a certificate.

• Section 4 presents an algorithm that, given a prime n, finds a reasonably small certificate for n in random time (lg n)^{2+o(1)}.

• Section 6 presents an algorithm to verify a reasonably small certificate in

time (lg n)^{4+o(1)}. One can prove that n is not a perfect power in time (lg n)^{1+o(1)}, as explained in [6], so prime-power proving is tantamount to prime proving.

Section 7 discusses verification speed in more detail. Some of the complications in the certificate definition are irrelevant to the 4+o(1) result but produce speedups visible at the level of detail of Section 7. A simplified proof that n is a prime power under stronger assumptions, without these complications, appears in Section 2.

Genealogy. This algorithm uses an idea that one might call “proving primality with combinatorics.” This idea was introduced by Agrawal, Kayal, and Saxena in [4]. (Primitive forms of the idea were used by Fellows and Koblitz in [12], and by Konyagin and Pomerance in [16].) The Agrawal-Kayal-Saxena algorithm proves primality in polynomial time, using combinatorics in cyclotomic extensions of Z/n.

The algorithm in this paper replaces cyclotomic extensions with random Kummer extensions, so that it can twist x − 1 into ζx − 1, ζ^{2}x − 1, etc.; see the proof of Theorem 2.1. This idea was introduced by Berrizbeitia in [9], in the special case of Kummer extensions whose degrees are powers of 2. Berrizbeitia’s algorithm proves primality in random time (lg n)^{4+o(1)} for a sparse set of primes n, namely those for which n^{2} − 1 is divisible by some power of 2 near (lg n)^{2}.

Cheng in [11] adapted Berrizbeitia’s idea to prime degrees. Cheng’s algorithm proves primality in random time (lg n)^{4+o(1)} for a larger set of primes n, namely those for which n − 1 is divisible by a prime e ≈ (lg n)^{2}.

Date: 20030417.

2000 Mathematics Subject Classification. Primary 11Y16.

The author was supported by the National Science Foundation under grant DMS–0140542, and by the Alfred P. Sloan Foundation. He used the libraries at the Mathematical Sciences Research Institute, the University of California at Berkeley, and the American Institute of Mathematics.

This paper generalizes to arbitrary positive integers e ≈ (lg n)^{2} dividing n^{d} − 1 for any d ∈ n^{o(1)}. A standard result from analytic number theory implies that every prime n has a suitable pair (d, e); see Theorems 5.1 and 5.2. For practical purposes, the only interesting case is d = 1, as discussed in Section 7.

My generalization was independent of Cheng’s adaptation. I read Berrizbeitia’s paper on 26 January 2003 and promptly sent email to a few people saying how I expected it to generalize to any n. I was then told about Cheng’s paper, which had been published on 16 January 2003. I published a draft of this paper, with a detailed proof of Theorem 3.2, on 28 January 2003, and announced the result on the NMBRTHRY mailing list on 29 January 2003.

Competition. Another way to prove the primality of n is to exhibit a factor of the Jacobian group of a hyperelliptic curve over Z/n. Adleman and Huang in [2] proved that every prime n has a certificate of this type that can be found in random time (lg n)^{O(1)} and verified in time at most (lg n)^{3+o(1)}. The O(1) here is large.

A previous algorithm of Atkin, using small-discriminant complex-multiplication elliptic curves, is conjectured to find a certificate of the same type in time at most (lg n)^{5+o(1)}. An improved algorithm, pointed out by Shallit and reported in [17, page 711], is conjectured to find a certificate of the same type in time at most (lg n)^{4+o(1)}. As above, the certificates can be verified in time (lg n)^{3+o(1)}.

The algorithm in this paper is proven to find and verify certificates in random time at most (lg n)^{4+o(1)}.

For readers who want to actually prove the primality of various numbers n, rather than prove theorems about how quickly one can prove primality, the impact of the new algorithm is less clear. Is the (lg n)^{4+o(1)} time for the new algorithm smaller than the (lg n)^{4+o(1)} time to find elliptic-curve certificates? (Note also that, for small n, the slightly-superpolynomial-time method of [3], [10], etc. is faster than finding elliptic-curve certificates.) My current impression is that the answer is no, but that further results along the lines of [7] could change the answer.

2. The idea in a nutshell

Theorem 2.1. Let n, d, and e be positive integers such that 2^{e} − 1 ≥ n^{2d⌊\sqrt{e}⌋} and e divides n^{d} − 1. Let f be a monic polynomial in (Z/n)[y] of degree d. Define R as the ring (Z/n)[y]/f. Let r be an element of R such that r^{n^{d}−1} = 1 in R, r^{(n^{d}−1)/q} − 1 is a unit in R for each prime q dividing e, and r − 1 is a unit in R. If (x − 1)^{n^{d}} = r^{(n^{d}−1)/e}x − 1 in the ring R[x]/(x^{e} − r) then n is a power of a prime.

Theorem 2.1 improves in two ways upon the theorems of Berrizbeitia in [9] and Cheng in [11]:

• d is allowed to be any positive integer. Berrizbeitia considered only d ∈ {1, 2}, and Cheng considered only d = 1. Larger d’s are important for the (lg n)^{4+o(1)} result in this paper. On the other hand, as discussed in Section 7, the case d = 1 is the only important case in practice.

• e is allowed to be any positive divisor of n^{d} − 1. Berrizbeitia considered only powers of 2 (although with slightly more general moduli x^{2^{i}e} − r), and Cheng considered only primes e; the proofs relied on e having only one prime divisor. Arbitrary e’s are important for the (lg n)^{4+o(1)} result in this paper. Arbitrary e’s also save time in practice, because they allow many more n’s to be handled with d = 1.

Theorem 3.2 saves more time by allowing somewhat smaller e’s.

After I published the first draft of this paper, Mihailescu and Avanzi published a more complicated theorem that requires R to be a Galois extension of Z/n (see [18]) and that uses nth powers rather than n^{d}th powers; the same change can be made in Theorem 3.2. This change speeds up verification if d ≥ 2, for two reasons: first, checking nth powers is about d times faster than checking n^{d}th powers; second, the hypotheses of the theorem allow e to be chosen about d times smaller. (When d is 2 and e is a power of 2, this variant was already used by Berrizbeitia in [9].) In practice, however, one always has d = 1, so the extra complications do not save any time.

Proof. If n = 1 then n is a power of a prime, so assume that n > 1.

Step 1: Move to a field. R is a nonzero ring, so it maps onto a field k. Explicitly: find a prime p dividing n; find an irreducible polynomial g in F_{p}[y] dividing the image of f; then k = F_{p}[y]/g is a field.

Write N = #R = n^{d} and P = #k = p^{deg g}. Define ζ as the image of r^{(N−1)/e} in k. Then ζ has order e in k. (Indeed, r^{N−1} = 1 in R by hypothesis, so ζ^{e} = 1 in k. Furthermore, if q is a prime dividing e, then r^{(N−1)/q} − 1 is a unit in R by hypothesis, so its image ζ^{e/q} − 1 in k is a unit; hence ζ^{e/q} ≠ 1 in k.) Consequently e divides P − 1.

Step 2: Combinatorially enumerate many powers of x − 1. Define A as the ring k[x]/(x^{e} − r). By hypothesis (x − 1)^{N} = r^{(N−1)/e}x − 1 in R[x]/(x^{e} − r), so (x − 1)^{N} = ζx − 1 in A. Substitute ζ^{i}x for x: (ζ^{i}x − 1)^{N} = ζ^{i+1}x − 1 in the ring k[x]/((ζ^{i}x)^{e} − r) = A. Thus (x − 1)^{N^{i}} = ζ^{i}x − 1 in A for any integer i ≥ 0.

There are 2^{e} −1 vectors (a_{0}, a_{1}, …, a_{e−1}) ∈ {0, 1}^{e} such that ∑_{i} a_{i} ≤ e−1. Any product ∏_{i}(ζ^{i}x − 1)^{a_{i}} = (x − 1)^{∑_{i} N^{i}a_{i}} is a power of x − 1 in A. I claim that these products are distinct, so there are at least 2^{e} − 1 powers of x − 1 in A.

Indeed, say a, b are two such vectors with ∏_{i}(ζ^{i}x − 1)^{a_{i}} = ∏_{i}(ζ^{i}x − 1)^{b_{i}} in A. Then ∏_{i}(ζ^{i}x − 1)^{a_{i}} = ∏_{i}(ζ^{i}x − 1)^{b_{i}} in k[x]: distinct polynomials of degree at most e − 1 remain distinct when reduced modulo x^{e} − r. The polynomials x − 1, ζx − 1, …, ζ^{e−1}x − 1 are coprime in k[x], so a_{i} = b_{i} by unique factorization.

Step 3: Find colliding powers of x − 1. The nonzero element r^{(P −1)/e} of k has eth power r^{P −1} = 1; but the eth roots of 1 in k are exactly the powers of ζ. Thus r^{(P −1)/e} = ζ^{ℓ} in k for some integer ℓ. Now x^{P} = x^{P −1}x = r^{(P −1)/e}x = ζ^{ℓ}x in A, so x^{P^{j}} = ζ^{jℓ}x in A for any integer j ≥ 0. Thus (x − 1)^{N^{i}P^{j}} = ζ^{i+jℓ}x − 1 in A.

Consider the pairs (i, j) with 0 ≤ i ≤ ⌊\sqrt{e}⌋ and 0 ≤ j ≤ ⌊\sqrt{e}⌋. There are (⌊\sqrt{e}⌋ + 1)^{2} > e pairs (i, j), and only e possible powers ζ^{i+jℓ}, so there are distinct pairs (i, j), (i′, j′) with ζ^{i+jℓ} = ζ^{i′+j′ℓ}. Define u = N^{i}P^{j} and v = N^{i′}P^{j′}; then u and v are positive integers bounded by N^{2⌊\sqrt{e}⌋}, and (x − 1)^{u} = (x − 1)^{v}.

The remainder (x^{e} − r) mod (x − 1) = 1 − r is a unit in k, so x − 1 is a unit in A. Thus (x − 1)^{u−v} = 1 in the unit group A^{∗}. If u ≠ v then there are at most |u − v| powers of x − 1, but |u − v| < N^{2⌊\sqrt{e}⌋} ≤ 2^{e} − 1, contradiction.

Hence u = v; i.e., N^{i−i′} = P^{j′−j}. If i = i′ then j′ = j, contradiction. Thus a nontrivial power of n is a power of p; so n is a power of p. □

3. Certificates

Definition 3.1. Let n, d, and e be positive integers. Let c and c_{−} be integers. Let f be a monic polynomial in (Z/n)[y] of degree d. Define R as the ring (Z/n)[y]/f. Let r be an element of R. Let S be a subset of R. Assume that

• e divides n^{d} − 1;

• e > c ≥ c_{−} ≥ 0;

• r^{n^{d}−1} = 1 in R;

• r^{(n^{d}−1)/q} − 1 is a unit in R for each prime q dividing e;

• s is a unit in R for all s ∈ S;

• s^{e} − (s′)^{e} is a unit in R for all distinct s, s′ ∈ S;

• s^{e} − r is a unit in R for all s ∈ S;

• \binom{e#S}{c_{−}}\binom{c}{c_{−}}\binom{e#S−c_{−}+e−1−c}{e−1−c} ≥ n^{d⌈\sqrt{e/3}⌉}; and

• (x − s)^{n^{d}} = r^{(n^{d}−1)/e}x − s in the ring R[x]/(x^{e} − r) for all s ∈ S. Then (d, e, c, c_{−}, f, r, S) is a certificate for n.

For example, (1, 840, 419, 246, y, 17, {1}) is a certificate for

31415926535897932384626433832795028841;

(1, 2430, 1214, 928, y, 2, {1, 2}) is a certificate for

2718281828459045235360287471352662497757247093699959574966967627724076630353547594571;

and (1, 57449, 28724, 16826, y, 2, {1}) is a certificate for 2^{1024} + 643. Theorem 3.2. Let n, d, and e be positive integers. Let c and c_{−} be integers. Let f be a monic polynomial in (Z/n)[y] of degree d. Define R as the ring (Z/n)[y]/f. Let r be an element of R. Let S be a subset of R. Assume that (d, e, c, c_{−}, f, r, S) is a certificate for n. Then n is a power of a prime.

Readers who warmed up by reading the proof of Theorem 2.1 will recognize the overall structure, and many components, of the proof of Theorem 3.2. However, through the definition of a certificate, Theorem 3.2 includes three features not present in Theorem 2.1:

• It uses \sqrt{e/3} as suggested by Lenstra, instead of 2\sqrt{e}. This reduces the lower bound on e by a factor of about 12. (Intermediate results by various people are not discussed here.)

• It allows c > 0 as suggested by Voloch, with an optimization suggested by Vaaler; see the proof to understand the role of c. This reduces the lower bound on e by an extra factor of about 6.4673912111.

• It allows #S to vary. (Berrizbeitia and Cheng considered only #S = 1, or #S = 2^{i} for modulus x^{2^{i}e} − r.) This allows e to be chosen even smaller, at some cost; see Section 7.

Smaller e’s allow faster verification of the hypotheses; on the other hand, they make the following proof more complicated. See below for additional comments.

Proof. If n = 1 then n is a power of a prime, so assume that n > 1.

Step 1: Move to a field. R is a nonzero ring, so it maps onto a field k. Explicitly: find a prime p dividing n; find an irreducible polynomial g in F_{p}[y] dividing the image of f; then k = F_{p}[y]/g is a field.

[figure: k[x]/h; R[x]/(x^{e} − r); k[x]/(x^{e} − r); R = (Z/n)[y]/f; k = F_{p}[y]/g]

Write N = #R = n^{d} and P = #k = p^{deg g}. Note that P divides N. If N = P then n is a power of p, so assume that N > P. (Similarly, one can assume that N ≠ P^{2}; this allows a variation in the Minkowski argument below.)

Define ζ as the image of r^{(N−1)/e} in k. Then ζ has order e in k. (Indeed, r^{N−1} = 1 in R by hypothesis, so ζ^{e} = 1 in k. Furthermore, if q is a prime dividing e, then r^{(N−1)/q} − 1 is a unit in R by hypothesis, so its image ζ^{e/q} − 1 in k is a unit; hence ζ^{e/q} ≠ 1 in k.) Consequently e divides P − 1.

Step 2: Combinatorially enumerate many group elements. Find an irreducible polynomial h in k[x] dividing the image of x^{e} −r. Then k[x]/h is a field.

If s ∈ S then (x − s)^{N} = r^{(N−1)/e}x − s in R[x]/(x^{e} − r) by hypothesis, so (x − s)^{N} = ζx − s in k[x]/(x^{e} − r). Substitute ζ^{i}x for x: (ζ^{i}x − s)^{N} = ζ^{i+1}x − s in k[x]/((ζ^{i}x)^{e} − r) = k[x]/(x^{e} − r). Thus (x − s)^{N^{i}} = ζ^{i}x − s in k[x]/(x^{e} − r) for any integer i ≥ 0.

Note that ζ^{i}x−s is a unit in k[x]/h. (The remainder (x^{e}−r) mod (ζ^{i}x−s) = s^{e}−r is a unit in k, so ζ^{i}x − s is a unit in k[x]/(x^{e} − r), hence in k[x]/h.) Note also that ζ^{i}x − s and ζ^{i′}x − s′ are coprime in k[x] unless (ζ^{i}, s) = (ζ^{i′}, s′). (If they are not coprime then sζ^{i′} = s′ζ^{i} in k, so s^{e} = (s′)^{e} in k. If s ≠ s′ then s^{e} − (s′)^{e} is a unit in R by hypothesis, so it is a unit in k; contradiction. Thus s = s′; so sζ^{i′} = sζ^{i}; also s is a unit in R by hypothesis, so ζ^{i′} = ζ^{i}.)

Consider functions a : {0, 1, …, e − 1} × S → Z such that • #{(i, s) : a(i, s) < 0} = c_{−};

• ∑_{i,s} −a(i, s)[a(i, s) < 0] ≤ c; and

• ∑_{i,s} a(i, s)[a(i, s) ≥ 0] ≤ e − 1 − c.

There are \binom{e#S}{c_{−}}\binom{c}{c_{−}}\binom{e#S−c_{−}+e−1−c}{e−1−c} ≥ N⌈^{\sqrt{e/3}}⌉ ≥ N^{\sqrt{e/3}} such functions. I claim that the products ∏_{i,s}(ζ^{i}x − s)^{a(i,s)} are distinct in (k[x]/h)^{∗}; so there are at least

N^{\sqrt{e/3}} such products.

Indeed, assume that a, b are two such functions, and that ∏_{i,s}(ζ^{i}x − s)^{a(i,s)} = ∏_{i,s}(ζ^{i}x − s)^{b(i,s)} in (k[x]/h)^{∗}. Clear denominators to obtain polynomials

A = ∏_{i,s}(ζ^{i}x − s)^{a(i,s)[a(i,s)≥0]−b(i,s)[b(i,s)<0]} ∈ k[x],

B = ∏_{i,s}(ζ^{i}x − s)^{b(i,s)[b(i,s)≥0]−a(i,s)[a(i,s)<0]} ∈ k[x]

with A = B in k[x]/h. Now A(ζ^{j}x) = A^{N^{j}} = B^{N^{j}} = B(ζ^{j}x) in k[x]/h, since ζ^{i+j}x−s = (ζ^{i}x−s)^{N^{j}} in k[x]/(x^{e}−r). Thus A−B has roots x, ζx, ζ^{2}x, …, ζ^{e−1}x in k[x]/h; these roots are distinct, since x is invertible in k[x]/h; but A−B has degree at most c + e − 1 − c < e, so it cannot have e roots unless it is zero. Thus A = B in k[x]; so a(i, s)[a(i, s) ≥ 0]−b(i, s)[b(i, s) < 0] = b(i, s)[b(i, s) ≥ 0]−a(i, s)[a(i, s) < 0] by unique factorization into coprimes; so a(i, s) = b(i, s).

Step 3: Find colliding powers. The nonzero element r^{(P −1)/e} of k has eth power r^{P −1} = 1; but the eth roots of 1 in k are exactly the powers of ζ. Thus r^{(P −1)/e} = ζ^{ℓ} in k for some integer ℓ. Now x^{P} = x^{P −1}x = r^{(P −1)/e}x = ζ^{ℓ}x in k[x]/(x^{e} − r), so x^{P^{j}} = ζ^{jℓ}x in k[x]/(x^{e} − r) for any integer j ≥ 0. Thus (x − s)^{N^{i}P^{j}} = ζ^{i+jℓ}x − s in k[x]/(x^{e} − r).

Define L as the set of (α, β) ∈ Z × Z such that e divides α + (β − α)ℓ; then L is a lattice of determinant e. Define C as the set of (α, β) ∈ R × R such that max \{|α| lg(N/P), |β| lg P, |α lg(N/P) + β lg P |\} ≤ \sqrt{e/3} lg N; then C is a closed convex symmetric set of area 3(e/3)(lg N)^{2}/(lg P) lg(N/P) ≥ 4e. By Minkowski’s theorem, there is a nonzero point (α, β) ∈ L∩C. Assume without loss of generality that α ≥ 0.

(Variation: N ≠ P^{2} so the area of C is larger than 4e. Thus one can use a simpler form of Minkowski’s theorem, ignoring the fact that C is closed.)

If β ≥ 0, define u = (N/P)^{α}P^{β} and v = 1; then u and v are positive integers bounded by N^{\sqrt{e/3}}, and (x − s)^{uP^{α}} = (x − s)^{N^{α}P^{β}} = ζ^{α+βℓ}x − s = ζ^{αℓ}x − s = (x − s)^{P^{α}} = (x − s)^{vP^{α}} in k[x]/(x^{e} − r). If β < 0, define u = (N/P)^{α} and v = P^{−β}; then u and v are positive integers bounded by N^{\sqrt{e/3}}, and (x−s)^{uP^{α}} = (x−s)^{N^{α}} = ζ^{α}x − s = ζ^{(α−β)ℓ}x − s = (x − s)^{P^{α−β}} = (x − s)^{vP^{α}} in k[x]/(x^{e} − r).

P th powering is invertible on the powers of x − s, since (x − s)^{P^{e}} = x − s; so (x − s)^{u} = (x − s)^{v} in k[x]/(x^{e} − r). Consequently (x − s)^{u−v} = 1 in (k[x]/h)^{∗}. Take N^{i}th powers: (ζ^{i}x − s)^{u−v} = 1 in (k[x]/h)^{∗}. Thus each of the aforementioned products π = ∏_{i,s}(ζ^{i}x − s)^{a(i,s)} in (k[x]/h)^{∗} satisfies π^{u−v} = 1.

On the other hand, if u ≠ v, then the field k[x]/h has at most |u − v| < N^{\sqrt{e/3}} roots π of the polynomial π^{|u−v|} − 1; contradiction. Thus u = v; i.e., N^{α} = P^{α−β}. If α = 0 then P^{−β} = 0 so β = 0, but (α, β) was nonzero by construction; contradiction. Hence n is a power of p. □

Notes on the proof. Consider the subgroup G of (k[x]/(x^{e} − r))^{∗} generated by {x − s : s ∈ S}. The proof may be summarized as follows: G is large; G is cyclic; if n is not a power of p then G has small exponent, contradiction.

There are three different ways to prove that G is cyclic:

• Choose e and r so that k[x]/(x^{e} − r) is forced to be a field. This is the approach used by Berrizbeitia and Cheng; it is also the reason for their restrictions on e.

• Used in Theorem 2.1: Choose S = {1}, so that G is cyclic by definition. This is the simplest approach.

• Used in Theorem 3.2: The available equations for x − s imply that G is always isomorphic to its image in (k[x]/h)^{∗}. This idea was first published by Macaj in [19]; it was discovered independently by Agrawal.

The original approach of Agrawal, Kayal and Saxena in [4] was to work instead with the cyclic image of G in (k[x]/h)^{∗}, and force the degree of h to be fairly large. A subsequent improvement by Lenstra was to work with an isomorphic image of G in a product of copies of (k[x]/h)^{∗}. Both approaches are quantitatively worse than proving that G is cyclic.

Except in the simple cyclic-by-definition case, all of these approaches rely on the fact that a finite multiplicative subgroup of a field is cyclic. The proof of that fact starts from the observation that a large subgroup cannot have small exponent, and then does some extra work to construct a generator. The extra work is unnecessary for this application: the only reason to prove that G is cyclic is to prove that it does not have small exponent. This simplification was pointed out to me by Kiran Kedlaya (who was preparing to explain [4] to high-school students); it is used in the proof of Theorem 3.2.

The proof of Theorem 3.2 would still work if n^{d⌈\sqrt{e/3}⌉} were replaced by n^{d\sqrt{e/3}} in the definition of a certificate. However, this change would complicate certificate testing, and would have very little benefit.

4. Finding a certificate: the algorithm

Every prime n has a certificate of the form (d, e, 0, 0, f, r, {1}) with d ∈ (lg n)^{o(1)} and e ∈ (lg n)^{2+o(1)}. Furthermore, this certificate can be found in random time (lg n)^{2+o(1)}. This section discusses the construction of d and e, then the construction of f, and finally the construction of r.

As discussed in Section 6, one can then verify that this is a certificate for n in time (lg n)^{4+o(1)}. As discussed in Section 7, one can reduce the o(1) by choosing certificates more carefully.

Finding d and e. There is a positive integer d such that n^{d} −1 has a divisor e ≥ 6 between d^{2} ⌈lg n⌉^{2} and (d + 1)d^{2} ⌈lg n⌉^{2}, by Theorem 5.1. The smallest such d is in exp(O(lg lg lg n lg lg lg lg n)) by Theorem 5.2.

To compute the smallest d, one can try d = 1, then d = 2, etc.; success will occur within (lg n)^{o(1)} tries. For each d, there are (lg n)^{2+o(1)} possible divisors e between d^{2} ⌈lg n⌉^{2} and (d + 1)d^{2} ⌈lg n⌉^{2}, each e having (lg n)^{o(1)} bits. One can compute n^{d} − 1 modulo all these e’s simultaneously in time (lg n)^{2+o(1)}; see, e.g., [8, Section 18].

Finding f. For every prime number n and positive integer d, there is a monic irreducible polynomial f ∈ (Z/n)[y] of degree d.

One standard way to find f is to generate a uniform random monic polynomial f of degree d, see if it is irreducible, and try again if not. There are many choices of f that work: the expected number of trials is approximately d. If d is chosen as above then the expected number of trials is in (lg n)^{o(1)}.

Another standard way to find f is to search systematically through polynomials with small coefficients. This avoids randomness, and produces polynomials f that take very little space to write down. It has the disadvantage that the number of trials is no longer guaranteed to be small.

One way to check the irreducibility of a single f is to see whether f has factors in common with x^{n} − x, x^{n^{2}} − x, …, x^{n^{d−1}} − x. Each nth powering in (Z/n)[y]/f takes time (lg n)^{2+o(1)} if d ∈ (lg n)^{o(1)}, so the total time for an irreducibility test is (lg n)^{2+o(1)}.

There is much more to say about the construction of irreducible polynomials. I should cite some recent survey. I should say a little about the impact of GRH, and about the improvements available as d grows.

Finding r. For every prime number n, positive integer d, positive integer e dividing n^{d} − 1, and monic irreducible polynomial f ∈ (Z/n)[y] of degree d, there is an element r of the field R = (Z/n)[y]/f such that r^{(n^{d}−1)/e} has order e; for example, any generator r of R^{∗}. Furthermore, if e ≥ 6 and e ≥ d^{2} ⌈lg n⌉^{2} as in Theorems 5.1 and 5.2, then (d, e, 0, 0, f, r, {1}) is a certificate for n by Theorem 5.3.

Finding elements of specified order is analogous to (and in many ways tied to) finding irreducible polynomials of specified degree. One standard way to find r is to generate a uniform random element r of R − {0}, see if r^{(n^{d}−1)/e} has order e, and try again if not. There are many choices of r that work: the expected number of trials is the product of q/(q − 1) for primes q dividing e, which is in (lg n)^{o(1)} if e ∈ (lg n)^{2+o(1)}.

Another standard way to find r is to search systematically through elements of (Z/n)[y]/f with small coefficients. This avoids randomness, and produces r’s that take very little space to write down; the reader may have noticed that the examples of certificates in Section 3 are very short. It has the disadvantage that the number of trials is no longer guaranteed to be small.

One way to check whether r^{(n^{d}−1)/e} has order e is to check that r^{n^{d}−1} = 1 and that r^{(n^{d}−1)/q} ≠ 1 for each prime q dividing e. There are only (lg n)^{o(1)} such primes q if e ∈ (lg n)^{2+o(1)}, and all of them are easy to find since e is small; the main work is to compute r^{(n^{d}−1)/e} in the first place, which takes time (lg n)^{2+o(1)}.

I should, again, point to the literature: combining orders, using GRH, merging exponentiations, etc.

5. Finding a certificate: the theorems

Theorem 5.1. Let n be an integer with n ≥ 2. Then there exists a positive integer d such that n^{d} − 1 has a divisor e ≥ 6 with d^{2} ⌈lg n⌉^{2} ≤ e < (d + 1)d^{2} ⌈lg n⌉^{2}.

Proof. Observe that n^{2} −1 ≥ ⌈lg n⌉^{2} and n^{6} +n^{4} +n^{2} +1 ≥ 64. Thus e ≥ 64 ⌈lg n⌉^{2} where e = n^{8} − 1. Define d as the largest multiple of 8 with d^{2} ≤ e/ ⌈lg n⌉^{2}. Then d ≥ 8, and e divides n^{d}−1. Furthermore, d^{3} ≥ 16d+64, so d^{3}+d^{2} ≥ d^{2}+16d+64 = (d + 8)^{2} > e/ ⌈lg n⌉^{2}. □

Theorem 5.2. There are constants n_{0} and α such that, for every prime number n ≥ n_{0}, there is a positive integer d ≤ exp(α log(3 log lg n) log log(3 log lg n)) such that n^{d} − 1 has a divisor e ≥ 6 with d^{2} ⌈lg n⌉^{2} ≤ e < (d + 1)d^{2} ⌈lg n⌉^{2}.

This is a typical application of a well-known theorem of Odlyzko and Pomerance; see [3, Theorem 3]. The point is that the product of the small primes dividing n^{d}−1 grows, at a minimum, almost exponentially with d. Older theorems suffice for the bound d ∈ (lg n)^{o(1)}.

It is overkill to assume that n is prime; what matters is that n has no tiny prime divisors.

Proof. Choose α such that d below always exists, and choose n_{0} > 8 such that H below always exists.

Given n ≥ n_{0}, select a real number H > 16 such that H ≤ (lg n)^{3}, D + 1 < n, and H/D^{2} ≥ ⌈lg n⌉^{2}, where D = exp(α log log H log log log H). Asymptotically one can take H in (lg n)^{2+o(1)}, and thus D in (lg n)^{o(1)}, satisfying H/D^{2} ≥ ⌈lg n⌉^{2}, so the extra constraints H ≤ (lg n)^{3} and D + 1 < n are automatically satisfied for n large enough. Note that D ≤ exp(α log(3 log lg n) log log(3 log lg n)).

By [3, Theorem 3], there is a positive integer d ≤ D such that H ≤ π, where π is the product of the primes q with q − 1 dividing d.

Now d^{2} ⌈lg n⌉^{2} ≤ D^{2} ⌈lg n⌉^{2} ≤ H ≤ π. Find the smallest positive integer e ≥ d^{2} ⌈lg n⌉^{2} dividing π; note that e ≥ 6 since n > 8. Each prime q is at most d + 1, so e must be smaller than (d + 1)d^{2} ⌈lg n⌉^{2}.

Finally, e divides n^{d} − 1. Indeed, each prime q is at most d + 1 ≤ D + 1 < n, so q does not divide n, so q divides n^{q−1} − 1, hence n^{d} − 1. □

Theorem 5.3. Let n be a prime number. Let d be a positive integer. Let e ≥ 6 be a divisor of n^{d} −1 such that d^{2} ⌈lg n⌉^{2} ≤ e. Let f be a monic irreducible polynomial in (Z/n)[y] of degree d. Let r be an element of the ring (Z/n)[y]/f such that r^{(n^{d}−1)/e} has order e. Then (d, e, 0, 0, f, r, {1}) is a certificate for n.

Proof. Write R = (Z/n)[y]/f. Observe that R is a field.

By hypothesis, n, d, and e are positive integers; e divides n^{d} − 1; r^{n^{d}−1} = (r^{(n^{d}−1)/e})^{e} = 1; if q is a prime dividing e, then r^{(n^{d}−1)/q}−1 = (r^{(n^{d}−1)/e})^{e/q}−1 ≠ 0, so r^{(n^{d}−1)/q} − 1 is a unit; and e > 1, so r^{(n^{d}−1)/e} ≠ 1, so r ≠ 1, so 1^{e} − r = 1 − r is a unit.

Furthermore, e ≥ 6, so \binom{2e−1}{e−1} ≥ 2^{e} and e ≥ (\sqrt{e/3} + 1)^{2}. Thus (lg \binom{2e−1}{e−1})^{2} ≥ e^{2} ≥ (\sqrt{e/3} + 1)^{2}e ≥ (\sqrt{e/3} + 1)^{2}d^{2} ⌈lg n⌉^{2} ≥ ⌈\sqrt{e/3}⌉^{2} d^{2}(lg n)^{2}; i.e., \binom{2e−1}{e−1} ≥

n^{d⌈\sqrt{e/3}⌉}.

Finally, (x − 1)^{n^{d}} = x^{n^{d}} − 1 = x^{n^{d}−1}x − 1 = r^{(n^{d}−1)/e}x − 1 in R[x]/(x^{e} − r). □

6. Checking a certificate

This section presents an algorithm that decides whether (d, e, c, c_{−}, f, r, S) is a certificate for n, given positive integers n, d, e, integers c and c_{−}, a monic degree-d polynomial f ∈ (Z/n)[y], an element r of R = (Z/n)[y]/f, and a subset S of R.

This algorithm takes time (lg n)^{4+o(1)} for reasonably small inputs. “Reasonably small” means that d is in (lg n)^{o(1)}; #S is in (lg n)^{o(1)}; e is at most (lg n)^{2+o(1)}; and e > c ≥ c_{−} ≥ 0. Note that the certificates (d, e, 0, 0, f, r, {1}) constructed in Section 4, with d ∈ (lg n)^{o(1)} and e ∈ (lg n)^{2+o(1)}, are reasonably small.

The reader is assumed to be familiar with fast multiplication. See, e.g., [8].

The basic conditions. Computing n^{d} − 1, and checking that it is divisible by e, takes time (lg n)^{1+o(1)}. Checking that e > c ≥ c_{−} ≥ 0 takes time (lg n)^{o(1)}.

Multiplying in Z/n takes time (lg n)^{1+o(1)}. Thus multiplying in R takes time (lg n)^{1+o(1)}. Computing the n^{d} −1 power of r in R takes (lg n)^{1+o(1)} multiplications in R, hence time (lg n)^{2+o(1)}.

The units. There are (lg n)^{o(1)} primes q dividing e; finding them by trial division takes time (lg n)^{1+o(1)}. Computing the (n^{d} − 1)/q power of r in R takes time (lg n)^{2+o(1)}. Checking whether r^{(n^{d}−1)/q} − 1 is a unit in R takes time (lg n)^{1+o(1)}.

Computing s^{e} in R for each s ∈ S takes time (lg n)^{1+o(1)}. Checking all the remaining units takes time (lg n)^{1+o(1)}.

The binomial coefficients. Computing (e#S)! takes time at most (lg n)^{2+o(1)}, since e#S is at most (lg n)^{2+o(1)}. Similarly, computing c_{−}! and (e#S − c_{−})! takes time at most (lg n)^{2+o(1)}. Thus computing the binomial coefficient \binom{e#S}{c_{−}} takes time at most (lg n)^{2+o(1)}. Similar comments apply to the other binomial coefficients.

Computing n^{d⌈\sqrt{e/3}⌉} takes time (lg n)^{2+o(1)}. Checking whether n^{d⌈\sqrt{e/3}⌉} ≤ \binom{e#S}{c_{−}}\binom{c}{c_{−}}\binom{e#S−c_{−}+e−1−c}{e−1−c} takes time (lg n)^{2+o(1)}.

The big exponentiation. Multiplying in R[x]/(x^{e} − r) takes time (lg n)^{3+o(1)}. Computing each (x − s)^{n^{d}} in R[x]/(x^{e} − r) takes (lg n)^{1+o(1)} multiplications in R[x]/(x^{e} − r), hence time (lg n)^{4+o(1)}.

7. Optimizations and practical performance

This section looks at verification speed more closely in the important case d = 1.

Why d = 1 in practice. A substantial fraction of primes n have suitable divisors e of n − 1: divisors slightly above the lower bound in Theorem 3.2. What about the other primes n?

A single elliptic-curve-primality-proving step conjecturally takes time (lg n)^{3+o(1)} to reduce the problem of proving the primality of n to the problem of proving the primality of an auxiliary prime n′. Here n′ is slightly shorter than n and “looks random.” One can find several choices for n′ at similar speed.

Consequently one can parlay a fast algorithm for a substantial fraction of primes n into a fast algorithm for all primes n. This was suggested by Cheng in [11].

Consider, for example, a prime n that does not have any suitable divisors of n − 1, but that does have suitable divisors of n^{2} − 1. Here are two ways to prove that n is prime:

• Apply Theorem 3.2 to n with d = 2.

• Use an elliptic-curve-primality-proving step to locate an auxiliary n′ that has suitable divisors of n′ − 1. Apply Theorem 3.2 to n′ with d = 1.

Experiments confirm that the (lg n)^{4+o(1)} savings in moving from d = 2 to d = 1 is far above the (lg n)^{3+o(1)} cost of a single elliptic-curve-primality-proving step.

Of course, the same argument might mean that even the d = 1 case is of no practical interest: perhaps Theorem 3.2 will never be faster than a series of ellipticcurve-primality-proving steps. On the other hand, perhaps future improvements to Theorem 3.2 will make certificate verification so fast that the d = 2 case is worth considering again.

Choosing c and c_{−}. The choice c = c_{−} = 0 in Section 4 is far from optimal. If c ≈ αe and c_{−} ≈ βe then the product \binom{e#S}{c_{−}}\binom{c}{c_{−}}\binom{e#S−c_{−}+e−1−c}{e−1−c} is approximately exp(eγ) where γ = (#S − β + 1 − α) log(#S − β + 1 − α) + #S log #S + α log α − 2(#S − β) log(#S − β) − 2β log β − (α − β) log(α − β) − (1 − α) log(1 − α).

One can, either with a computer program or by hand, easily find α and β that maximize γ for any given #S. Any choice of c ≈ αe and c_{−} ≈ βe is reasonable; a small amount of additional searching will locate the optimal c and c_{−}.

It turns out that the optimal α and β have simple expressions: α = 1/2 and β = (#S + 1 − \sqrt{#S^{2} + 1})/2. For example, say #S = 1. The product of binomial coefficients is about 5.828427 …^{e} if one takes c ≈ e/2 and c_{−} ≈ (2 − \sqrt{2})e/2 = (0.2928932 …)e. For comparison: The product of binomial coefficients is about 4^{e} if one takes c = 0 and c_{−} = 0.

Choosing e and #S. Say there are many possibilities for (e, #S)—or, in the elliptic-curve context, many possibilities for an auxiliary (n, e, #S)—such that the maximized product of binomial coefficients exceeds n⌈^{\sqrt{e/3}}⌉. One should choose the possibility that minimizes verification time.

As a first approximation, this means minimizing e#S: verification time can be crudely modeled as (lg n)^{2}e#S. The following table shows e#S/(lg n)^{2} as a function of e/(lg n)^{2}, when #S is chosen as small as possible:

works for e/(lg n)^{2}

so e#S/(lg n)^{2}

#S between about and about

is between about

and about

1  0.051540 …  ∞  0.051540 …  ∞
2  0.027664 … 0.051540 …  0.055328 … 0.103081 …
3  0.020415 … 0.027664 …  0.061247 … 0.082992 …
4  0.016832 … 0.020415 …  0.067328 … 0.081663 …
5  0.014653 … 0.016832 …  0.073269 … 0.084160 …
6  0.013169 … 0.014653 …  0.079017 … 0.087923 …
7  0.012082 … 0.013169 …  0.084575 … 0.092187 …
8  0.011244 … 0.012082 …  0.089958 … 0.096658 …

If e drops substantially below 0.01(lg n)^{2} then e#S explodes: #S = 100 works for e/(lg n)^{2} down to about 0.004037 …; #S = 1000 works for e/(lg n)^{2} down to about 0.002164 …; #S = 10000 works for e/(lg n)^{2} down to about 0.001347 …; and so on.

A more precise model of verification time includes logarithmic factors that grow with e but not with #S. Reducing e at the expense of #S often saves time even if it increases e#S.

Multiplying quickly. One can square an element of (Z/n)[x]/(x^{e} − r) as follows:

• Lift to Z[x], obtaining polynomials of degree at most e − 1 with coefficients between −n/2 and n/2.

• Choose p so that 2^{p} > e(n/2)^{2}, and map to Z[x]/(x − 2^{p}) ≅ Z, obtaining an integer with approximately 2e lg n bits.

• Square in Z.

• Recover the product in Z[x].

• Reduce modulo x^{e} − r. This is particularly easy if r is small.

• Reduce each coefficient modulo n. The overall speed of certificate verification depends crucially on the details of these steps.

For example, one might square a polynomial with the following C code, using the GMP 4.1.2 library:

      mpz_set_si(t1,0);
      for (j = 0;j < e;++j)
        for (b = 0;b < nbits;++b)
          if (mpz_tstbit(poly[j],b))
            mpz_setbit(t1,j * padbits + b);
      mpz_mul(t1,t1,t1);

      for (j = 0;j < e;++j) {
        mpz_set_si(t2,0);
        for (b = 0;b < padbits;++b)
          if (mpz_tstbit(t1,(j + e) * padbits + b))
            mpz_setbit(t2,b);
        mpz_set_si(t3,r);
        mpz_mul(t2,t2,t3);
        mpz_set_si(t3,0);
        for (b = 0;b < padbits;++b)
          if (mpz_tstbit(t1,j * padbits + b))
            mpz_setbit(t3,b);
        mpz_add(t2,t3,t2);
        mpz_mod(poly[j],t2,n);
      }

This code uses approximately 2·10^{11} clock cycles on a Pentium III-800 to verify the aforementioned certificate (1, 2430, 1214, 928, y, 2, {1, 2}) for the prime ⌊10^{84} exp 1⌋. A large fraction of the time is spent testing and setting bits; GMP does not offer any good way to copy a stretch of bits from one number to another. Another large fraction of the time is spent in integer squaring, which can be sped up considerably. I would not be surprised to see an order of magnitude speed improvement.

Future improvements. It would be surprising for the group discussed in Section 3 to have size much smaller than p^{e}. Theorem 2.1 uses a lower bound near 2^{e}. Theorem 3.2 uses a lower bound near 5.828427 …^{e}. Can one do better?

An improvement from 2^{e} to 2^{γe} means that the lower bound on e drops by a factor of about γ^{2}. A lower bound close to p^{e}, with #S small, would reduce the lower bound on e to about (4 lg(n/p))/3 lg p < (4/3) lg n, saving a factor of (lg n)^{1+o(1)}. Mihailescu has pointed out that in this case one can further reduce the lower bound on e by using unit-group factors to quickly increase the lower bound on primes p dividing n.

Voloch in [22] suggested considering products in F_{p}[x] of degree somewhat larger than e, and applying the ABC theorem. I showed in [7] that, if #S = 1 and n does not have any tiny factors, then four distinct products of degree at most 1.1e cannot be congruent modulo x^{e} − r, so the group size is at least \frac{1}{3} \binom{⌊2.1e⌋}{e} ≈ 4.2768947738 …^{e}. Perhaps there are better results along these lines. (I suspect that any such results will work with multiple derivatives directly, rather than applying the ABC theorem.) If 1000 distinct products of degree at most 2e cannot be congruent modulo x^{e} − r then the group size is at least \frac{1}{999} \binom{3e}{e} ≈ 6.75^{e}.

References

[1] —, Proceedings of the 18th annual ACM symposium on theory of computing, Association for Computing Machinery, New York, 1986. ISBN 0–89791–193–8.

[2] Leonard M. Adleman, Ming-Deh A. Huang, Primality testing and abelian varieties over finite fields, Lecture Notes in Mathematics, 1512, Springer-Verlag, Berlin, 1992. ISBN 3– 540–55308–8. MR 93g:11128.

[3] Leonard M. Adleman, Carl Pomerance, Robert S. Rumely, On distinguishing prime numbers from composite numbers, Annals of Mathematics 117 (1983), 173–206. ISSN 0003–486X. MR 84e:10008.

[4] Manindra Agrawal, Neeraj Kayal, Nitin Saxena, PRIMES is in P (2002). Available from

    http://www.cse.iitk.ac.in/news/primality.html.

[5] A. O. L. Atkin, Francois Morain, Finding suitable curves for the elliptic curve method of factorization, Mathematics of Computation 60 (1993), 399–405. ISSN 0025–5718. MR 93k:11115.

[6] Daniel J. Bernstein, Detecting perfect powers in essentially linear time, Mathematics of Computation 67 (1998), 1253–1283. ISSN 0025–5718. MR 98j:11121. Available from http://cr.

yp.to/papers.html.

[7] Daniel J. Bernstein, Sharper ABC-based bounds for congruent polynomials. Available from

http://cr.yp.to/papers.html.

[8] Daniel J. Bernstein, Fast multiplication and its applications. Available from http://cr.yp.

to/papers.html.

[9] Pedro Berrizbeitia, Sharpening PRIMES is in P for a large family of numbers (2002). Avail-

able from http://arxiv.org/abs/math.NT/0211334.

[10] Wieb Bosma, Marc-Paul van der Hulst, Primality proving with cyclotomy, Ph.D. thesis, Universiteit van Amsterdam, 1990.

[11] Qi Cheng, Primality proving via one round in ECPP and one iteration in AKS (2003).

Available from http://www.cs.ou.edu/~qcheng/.

[12] Michael R. Fellows, Neal Koblitz, Self-witnessing polynomial-time complexity and prime factorization, Designs, Codes and Cryptography 2 (1992), 231–235. ISSN 0925–1022. MR 93e:68032.

[13] Shafi Goldwasser, Joe Kilian, Almost all primes can be quickly certified, in [1] (1986), 316–329; see also newer version in [14].

[14] Shafi Goldwasser, Joe Kilian, Primality testing using elliptic curves, Journal of the ACM 46 (1999), 450–472; see also older version in [13]. ISSN 0004–5411. MR 2002e:11182.

[15] Ronald L. Graham, Jaroslav Nešetřil (editors), The mathematics of Paul Erdős. I, Algorithms and Combinatorics, 13, Springer-Verlag, Berlin, 1997. ISBN 3–540–61032–4. MR 97f:00032.

[16] Sergei Konyagin, Carl Pomerance, On primes recognizable in deterministic polynomial time, in [15] (1997), 176–198. MR 98a:11184. Available from http://cr.yp.to/bib/entries.html#

1997/konyagin.

[17] Arjen K. Lenstra, Hendrik W. Lenstra, Jr., Algorithms in number theory, in [21] (1990), 673–715.

[18] Hendrik W. Lenstra, Jr., Galois theory and primality testing, in [20] (1985), 169–189. MR 87g:11171.

[19] Martin Macaj, Some remarks and questions about the AKS algorithm and related conjecture

(2002). Available from http://thales.doa.fmph.uniba.sk/macaj/aksremarks.pdf.

[20] I. Reiner, K. W. Roggenkamp (editors), Orders and their applications: proceedings of the conference held in Oberwolfach, June 3–9, 1984, Lecture Notes in Mathematics, 1142, Springer-Verlag, Berlin, 1985. ISBN 3–540–15674–7. MR 86g:16003.

[21] Jan van Leeuwen (editor), Handbook of theoretical computer science, volume A, Elsevier, Amsterdam, 1990. ISBN 0–444–88071–2. MR 92d:68001.

[22] José Felipe Voloch, On some subgroups of the multiplicative group of finite rings. Available

from http://www.ma.utexas.edu/users/voloch/preprint.html.

Department of Mathematics, Statistics, and Computer Science (M/C 249), The University of Illinois at Chicago, Chicago, IL 60607–7045

E-mail address: djb@cr.yp.to