David Harvey
email: d.harvey@unsw.edu.au

Joris Van Der Hoeven

Grégoire Lecerf
email: lecerf@lix.polytechnique.fr

Faster polynomial multiplication over finite fields

Keywords: Polynomial multiplication, finite field, algorithm, complexity bound, FFT A.C.M. subject classification: G.1.0 Computer-arithmetic, F.2.1 Computations in finite fields A.M.S. subject classification: 68W30, 68Q17, 68W40

Let p be a prime, and let M p (n) denote the bit complexity of multiplying two polynomials in F p [X] of degree less than n. For n large compared to p, we establish the bound M p (n) = O(n log n 8 log * n log p), where log * is the iterated logarithm. This is the first known Fürertype complexity bound for F p [X], and improves on the previously best known bound M p (n) = O(n log n log log n log p).

Introduction

Given a ring R, a basic problem in complexity theory is to find an upper bound for the cost of multiplying two polynomials in R[X] of degree less than n. Several complexity models may be considered. In algebraic complexity models, such as the straight-line program model [START_REF] Bürgisser | Algebraic complexity theory[END_REF]Chapter 4], we count the number of ring operations in R, denoted M R (n). In this model, the best currently known bound M R (n) = O(n log n log log n) was obtained by Cantor and Kaltofen [START_REF] Cantor | On fast multiplication of polynomials over arbitrary algebras[END_REF]. More precisely, their algorithm performs O(n log n log log n) additions and subtractions and O(n log n) multiplications in R. Their bound may be viewed as an algebraic analogue of the Schönhage-Strassen result for integer multiplication [START_REF] Schönhage | Schnelle Multiplikation großer Zahlen[END_REF] and generalises previous work by Schönhage [START_REF] Schönhage | Schnelle Multiplikation von Polynomen über Körpern der Charakteristik 2[END_REF]. These algorithms all rely on suitable incarnations of the fast Fourier transform (FFT) [START_REF] Cooley | An algorithm for the machine calculation of complex Fourier series[END_REF]. For details we refer the reader to [START_REF] Zur Gathen | Modern Computer Algebra[END_REF]Chapter 8].

In this paper we are mainly interested in the case that R is the finite field F p = Z / p Z for some prime p. The standard bit complexity model based on deterministic multitape Turing machines [START_REF] Papadimitriou | Computational Complexity[END_REF] is more realistic in this setting, as it takes into account the dependence on p. We write M p (n) for the bit complexity of multiplying two polynomials in F p [X] of degree less than n. Two basic approaches are known for obtaining good asymptotic bounds for M p (n).

If n is not too large compared to p, then one may use Kronecker substitution, which converts the problem to integer multiplication by packing the coefficients into large integers. Let I(n) denote the bit complexity of n-bit integer multiplication. Throughout the paper, we make the customary assumption that I(n)/n is increasing, and for convenience we define I(x) := I(⌈x⌉) for x ∈ R, x > 0. According to our recent sharpening [START_REF] Harvey | Even faster integer multiplication[END_REF] of Fürer's bound [START_REF] Fürer | Faster integer multiplication[END_REF], we may take

I(n) = O(n log n 8 log * n ), (1.1)
where log * x denotes the iterated logarithm for x ∈ R, i.e., log * x := min {k ∈ N:

log •k x 1}, (1.2) log •k := log • ••• k× • log.
If log n = O(log p), then Kronecker substitution yields M p (n) = O(I(n log p)) (see section 2.6), so we have M p (n) = O n log p log (n log p) 8 log * (n log p) .

However, Kronecker substitution leads to inferior bounds when n is large compared to p, due to coefficient growth in the lifted product in Z[X]. In this situation the best existing method is the algebraic Schönhage-Strassen algorithm (i.e., the Cantor-Kaltofen algorithm). In the Turing model this yields the bound M p (n) = O(n log n log log n log p + n log n I(log p)).

The first term dominates for large enough n, say for log log p 8 log * p = O(log log n), and we get simply M p (n) = O(n log n log log n log p).

(1.3) Since Fürer's breakthrough result, there has been a gap between what is known for integer multiplication and for multiplication in F p [X]. Namely, it has remained an open question whether the factor log log n appearing in (1.3) can be improved to K log * n for some K > 1. It seems that Fürer's construction does not work in F p [X].

Our contributions

The main result of this paper is a new algorithm that completely closes the gap mentioned above. holds uniformly for all n 2 and all primes p.

The bound follows immediately via Kronecker substitution if log n = O(log p). The bulk of the paper concerns the reverse case where n is large compared to p. For example, if log log p = O(log n), the bound becomes simply M p (n) = O(n log n 8 log * n log p).

As promised, this replaces the factor log log n in (1.3) by 8 log * n .

The basic idea of the new algorithm is as follows. We first construct a special extension F p κ of F p , where κ is exponentially smaller than n, and for which p κ -1 has many small divisors. These divisors are themselves exponentially smaller than n, yet they are so numerous that their product is comparable to n. We now convert the given multiplication problem in F p [X] to multiplication in F p κ[Y ], by cutting up the polynomials into small pieces, and we then multiply in F p κ[Y ] by using FFTs over F p κ. Applying the Cooley-Tukey method to the small divisors of p κ -1, we decompose the FFTs into transforms of exponentially shorter lengths (not necessarily powers of two). We use Bluestein's method to convert each short transform to a convolution over F p κ, and then use Kronecker substitution (from bivariate to univariate polynomials) to convert this to multiplication in F p [X]. These latter multiplications are then handled recursively. We continue the recursion until n is comparable to p, at which point we switch to using Kronecker substitution (from polynomials to integers).

In many respects this approach is very similar to the algorithm for integer multiplication introduced in [START_REF] Harvey | Even faster integer multiplication[END_REF], but there are various additional technical issues to address. We recommend [START_REF] Harvey | Even faster integer multiplication[END_REF] as a gentle introduction to the main ideas.

We also prove the following conditional result. holds uniformly for all n 2 and all primes p.

Of these conjectures, Conjecture 8.3 seems to be in reach for specialists in analytic number theory; Conjecture 8.5 is at least as hard as Artin's conjecture on primitive roots, but may be tractable under the generalised Riemann hypothesis (GRH); and Conjecture 8.1 is almost as strong as the Lenstra-Pomerance-Wagstaff conjecture on the distribution of Mersenne primes, and seems to be inaccessible at present.

Outline of the paper

The paper is structured as follows. In section 2, we start with a survey of relevant basic techniques: discrete Fourier transforms (DFTs) and FFTs, Bluestein's chirp transform, and Kronecker substitution. Most of this material is repeated from [START_REF] Harvey | Even faster integer multiplication[END_REF]Section 2] for the convenience of the reader; however, sections 2.5 and 2.6 differ substantially.

In section 3, we recall basic complexity results for arithmetic in finite fields. In particular, we consider the construction of irreducible polynomials in F p [X], algorithms for finding roots of unity in F p κ, and the cost of arithmetic in F p κ and F p κ[Y ].

In section 4, we show how to construct special extensions of F p whose multiplicative group has a large subgroup of highly smooth order, i.e., is divisible by many small integers. The main tool is [START_REF] Adleman | On distinguishing prime numbers from composite numbers[END_REF]Theorem 3].

Section 5 gives complexity bounds for functions that satisfy recurrence inequalities involving postcompositions with "logarithmically slow" functions. The prototype of such an inequality is T (n) K T (log n) + L, where K and L are constants. The definitions and theorems are duplicated from [START_REF] Harvey | Even faster integer multiplication[END_REF]Section 5]; for the proofs, see [START_REF] Harvey | Even faster integer multiplication[END_REF].

To minimise the constant K in the bound M p (n) = O(n log n K log * n log p) (for n large relative to p), we need one more tool: in section 6, we present a polynomial analogue of the Crandall-Fagin convolution algorithm [START_REF] Crandall | Discrete weighted transforms and large-integer arithmetic[END_REF]. This allows us to convert a cyclic convolution over F p of length n, where n is arbitrary, into a cyclic convolution over F p κ of somewhat smaller length N , where N is prescribed and where κ ≈ 2 n/ N . (We can still obtain K = 16 without using this Crandall-Fagin analogue, but we do not know how to reach K = 8 without it.) Section 7 is devoted to the proof of Theorem 1.1. Section 8 gives the proof of Theorem 1.2, and discusses the three number-theoretic hypotheses of that theorem.

The last section offers some final notes and suggested directions for generalisation. We first quickly dispense with the bit complexity of multiplication in F q [X] for prime powers q, and in (Z / m Z)[X] for arbitrary integers m 1. We then sketch some algebraic complexity bounds for polynomial multiplication over F p -algebras and (Z / m Z)-algebras, especially in the straightline program model. Our techniques also give rise to new strategies for polynomial evaluationinterpolation over F q . This may for instance be applied to the efficient multiplication of polynomial matrices over F q . Although we have not implemented any of our algorithms yet, we conclude by a few remarks on their expected efficiency in practice. (1) (see [START_REF] Zur Gathen | Modern Computer Algebra[END_REF]Chapter 25,Section 7] for details). The symbol R denotes the set of non-negative real numbers, and N denotes {0, 1, 2, ...}. For a ring R and n ∈ N, we write R[X] n := {P ∈ R[X]: deg P < n}.

Notations. We use Hardy's notations

f ≺ g for f = o(g), and f ≍ g for f = O(g) and g = O(f ). The "soft-Oh" notation f (n) = O ˜(g(n)) means that f (n) = g(n) (log(3 + g(n))) O
We will write lg n := ⌈log n/log 2⌉. In expressions like log log p or lg lg p, we tacitly assume that the function is adjusted so as to take positive values for small primes such as p = 2.

What is significant from a complexity point of view is that occasionally we must switch representations, to access an array (say 2-dimensional) by "rows" or by "columns". In the Turing model, we may transpose an n 1 × n 2 matrix of b-bit elements in time O(b n 1 n 2 lg min (n 1 , n 2 )), using the algorithm of [START_REF] Bostan | Linear recurrences with polynomial coefficients and application to integer factorization and Cartier-Manin operator[END_REF]Appendix]. Briefly, the idea is to split the matrix into two halves along the "short" dimension, and transpose each half recursively.

We will also require more complex rearrangements of data, for which we resort to sorting. Suppose that X is a totally ordered set, whose elements are represented by bit strings of length b, and suppose that we can compare elements of X in time O(b). Then an array of n elements of X may be sorted in time O(b n lg n) using merge sort [START_REF] Knuth | The art of computer programming, volume 3: Sorting and Searching[END_REF], which can be implemented efficiently on a Turing machine.

Discrete Fourier transforms

Let R be a commutative ring with identity and let n 1. An element ω ∈ R is said to be a principal n-th root of unity if ω n = 1 and k=0 n-1

(ω i ) k = 0 (2.1)
for all i ∈ {1, ..., n -1}. In this case, we define the discrete Fourier transform (or DFT) of an n-tuple a = (a 0 , ..., a n-1 ) ∈ R n with respect to ω to be DFT ω (a) = a ˆ= (a ˆ0, ..., a ˆn-1 ) ∈ R n where

a ˆi := a 0 + a 1 ω i + ••• + a n-1 ω (n-1)i .
That is, a ˆi is the evaluation of the polynomial A(X)

:= a 0 + a 1 X + ••• + a n-1 X n-1 at ω i .
If ω is a principal n-th root of unity, then so is its inverse ω -1 = ω n-1 , and we have

DFT ω -1(DFT ω (a)) = n a.
Indeed, writing b := DFT ω -1(DFT ω (a)), the relation (2.1) implies that

b i = j =0 n-1 a ˆj ω -ji = j =0 n-1 k=0 n-1 a k ω j (k-i) = k=0 n-1 a k j =0 n-1 ω j (k-i) = k=0 n-1 a k (n δ i,k ) = n a i ,
where δ i,k = 1 if i = k and δ i,k = 0 otherwise.

Remark 2.1. In all of the new algorithms introduced in this paper, we actually work over a field, whose characteristic does not divide n. In this setting, the concept of principal root of unity coincides with the more familiar primitive root of unity.

The Cooley-Tukey FFT

Let ω be a principal n-th root of unity and let n = n 1 n 2 where 1 < n 1 < n. Then ω n1 is a principal n 2 -th root of unity and ω n2 is a principal n 1 -th root of unity. Moreover, for any i 1 ∈ {0, ..., n 1 -1} and i 2 ∈ {0, ..., n 2 -1}, we have

a ˆi1n2+i2 = k1=0 n1-1 k2=0 n2-1 a k2n1+k1 ω (k2n1+k1)(i1n2+i2) = k 1 =0 n1-1 ω k1i2 k 2 =0 n2-1 a k2n1+k1 (ω n1 ) k2i2 (ω n2 ) k1i1 . ( 2

.2)

If A 1 and A 2 are algorithms for computing DFTs of length n 1 and n 2 , we may use (2.2) to construct an algorithm A 1 ⊙ A 2 for computing DFTs of length n as follows.

For each k 1 ∈ {0, ..., n 1 -1}, the sum inside the brackets corresponds to the i 2 -th coefficient of a DFT of the n 2 -tuple (a 0n 1 +k 1 , ..., a (n 2 -1)n 1 +k 1 ) ∈ R n2 with respect to ω n1 . Evaluating these inner DFTs requires n 1 calls to A 2 . Next, we multiply by the twiddle factors ω k1i2 , at a cost of n operations in R. (Actually, fewer than n multiplications are required, as some of the twiddle factors are equal to 1. This optimisation, while important in practice, has no asymptotic effect on the algorithms discussed in this paper.) Finally, for each i 2 ∈ {0, ..., n 2 -1}, the outer sum corresponds to the i 1 -th coefficient of a DFT of an n 1 -tuple in R n1 with respect to ω n2 . These outer DFTs require n 2 calls to A 1 .

Denoting by F R (n) the number of ring operations needed to compute a DFT of length n, and assuming that we have available a precomputed table of twiddle factors, we obtain

F R (n 1 n 2 ) n 1 F R (n 2 ) + n 2 F R (n 1 ) + n. For a factorisation n = n 1 ••• n d , this yields recursively F R (n) i=1 d n n i F R (n i ) + (d -1) n. (2.3)
The corresponding algorithm is denoted

A 1 ⊙ ••• ⊙ A d .
The ⊙ operation is neither commutative nor associative; the above expression will always be taken to mean

(•••((A 1 ⊙ A 2 ) ⊙ A 3 ) ⊙ •••) ⊙ A d .
The above discussion requires several modifications in the Turing model. Assume that elements of R are represented by b bits.

First, for A 1 ⊙ A 2 , we must add a rearrangement cost of O(b n lg min (n 1 , n 2 )) to efficiently access the rows and columns for the recursive subtransforms (see section 2.1). For the general case

A 1 ⊙ ••• ⊙ A d , the total rearrangement cost is bounded by O( i b n lg n i ) = O(b n lg n).
Second, we will sometimes use non-algebraic algorithms to compute the subtransforms, so it may not make sense to express their cost in terms of F R . The relation (2.3) therefore becomes

F(n) i=1 d n n i F(n i ) + (d -1) n m R + O(b n lg n), (2.4)
where F(n) is the (Turing) cost of a transform of length n over R, and where m R is the cost of a single multiplication in R.

Finally, we point out that A 1 ⊙ A 2 requires access to a table of twiddle factors ω i1i2 , ordered lexicographically by (i 1 , i 2 ), for 0 i 1 < n 1 , 0 i 2 < n 2 . Assuming that we are given as input a precomputed table of the form 1, ω, ..., ω n-1 , we must show how to extract the required twiddle factor table in the correct order. We first construct a list of triples (i The corresponding cost for A 1 ⊙ ••• ⊙ A d is determined as follows. Assuming that the table 1, ω, ..., ω n-1 is given as input, we first extract the subtables of (n

1 , i 2 , i 1 i 2 ), ordered by (i 1 , i 2 ), in time O(n lg n); then sort by i 1 i 2 in time O(n lg 2 n) (
1 ••• n i )-th roots of unity for i = d -1, ..., 2 in time O((n 1 ••• n d + ••• + n 1 n 2 ) (b + lg n)) = O(n (b + lg n)). Extracting the twiddle factor table for the decomposition (n 1 ••• n i-1 ) × n i then costs O(n 1 ••• n i lg n (b + lg n)); the total over all i is again O(n lg n (b + lg n)).
Remark 2.2. An alternative approach is to compute the twiddle factors directly in the correct order. When working over C, as in [START_REF] Harvey | Even faster integer multiplication[END_REF]Section 3], this requires a slight increase in the working precision. Similar comments apply to the root tables used in Bluestein's algorithm in section 2.5.

Fast Fourier multiplication

Let ω be a principal n-th root of unity in R and assume that n is invertible in R. Consider two polynomials

A = a 0 + ••• + a n-1 X n-1 and B = b 0 + ••• + b n-1 X n-1 in R[X]. Let C = c 0 + ••• + c n-1 X n-1 be the polynomial defined by c := 1 n DFT ω -1(DFT ω (a) DFT ω (b)),
where the product of the DFTs is taken pointwise. By construction, we have c ˆ= a ˆb ˆ, which means that C(ω i ) = A(ω i ) B(ω i ) for all i ∈ {0, ..., n -1}. The product

S = s 0 + ••• + s n-1 X n-1 of A and B modulo X n -1 also satisfies S(ω i ) = A(ω i ) B(ω i ) for all i. Consequently, s ˆ= a ˆb ˆ, s = DFT ω -1(s ˆ)/n = c, whence C = S.
For polynomials A, B ∈ R[X] with deg A < n and deg B < n, we thus obtain an algorithm for the computation of A B modulo X n -1 using at most 3 F R (n) + O(n) operations in R. Modular products of this type are also called cyclic convolutions. If deg (A B) < n, then we may recover the product A B from its reduction modulo X n -1. This multiplication method is called FFT multiplication.

If one of the arguments (say B) is fixed and we want to compute many products A B (or cyclic convolutions) for different A, then we may precompute DFT ω (b), after which each new product A B can be computed using only 2 F R (n) + O(n) operations in R.

Bluestein's chirp transform

We have shown above how to multiply polynomials using DFTs. Inversely, it is possible to reduce the computation of DFTs -of arbitrary length, not necessarily a power of two -to polynomial multiplication [START_REF] Bluestein | A linear filtering approach to the computation of discrete Fourier transform[END_REF], as follows.

Let ω be a principal n-th root of unity and consider an n-tuple a ∈ R n . First consider the case that n is odd. Let

f i := ω (i 2 -i)/2 , f i ′ := ω (i 2 +i)/2 , g i := ω (-i 2 -i)/2 .
Note that (i2 i)/2 and (i 2 + i)/2 are integers, and that

g i+n = ω (-(i+n) 2 -(i+n))/2 = ω (-i 2 -i)/2 ω n(-i-(n+1)/2) = g i . Let F := f 0 a 0 + ••• + f n-1 a n-1 X n-1 , G := g 0 + ••• + g n-1 X n-1 and C := c 0 + ••• + c n-1 X n-1 ≡ F G modulo X n -1. Then f i ′ c i = j =0 n-1 f i ′ (f j a j ) g i-j = j =0 n-1 ω (i 2 +i)/2 ω (j 2 -j)/2 ω (-(i-j) 2 -(i-j))/2 a j = j =0 n-1 ω ij a j = a ˆi.
In other words, the computation of a DFT of odd length n reduces to a cyclic convolution product of the same length n, together with O(n) additional operations in R. Notice that the polynomial G is fixed and independent of a in this product. Now suppose that n is even. In this case we require that 2 be invertible in R, and that ω n/2 = -1. Let σ := (-1) n/2 , and put

f i := ω i 2 , f i ′ := ω i 2 +i , g i := ω -i 2 + ω -i 2 -i .
Then

g i+ n 2 = σ ω -i 2 -ω -i 2 -i , g i+n = g i , and
1 2 g i + σ g i+ n 2 = ω -i 2 , 1 2 g i -σ g i+ n 2 = ω -i 2 -i .
As above, let

F := f 0 a 0 + ••• + f n-1 a n-1 X n-1 , G := g 0 + ••• + g n-1 X n-1 and C := c 0 + ••• + c n-1 X n-1 ≡ F G modulo X n -1. Then for i ∈ 0, ..., n 2 -1 we have 1 2 f i c i + σ c i+ n 2 = 1 2 f i j =0 n-1 f j a j g i-j + σ g i-j + n 2 = j =0 n-1 ω i 2 ω j 2 ω -(i-j) 2 a j = j =0 n-1 ω 2ij a j = a ˆ2i ,
and similarly

1 2 f i ′ c i -σ c i+ n 2 = j=0 n-1 ω i 2 +i ω j 2 ω -(i-j) 2 -(i-j) a j = j =0 n-1 ω (2i+1) j a j = a ˆ2i+1 .
Again, the DFT reduces to a convolution of length n, together with O(n) additional operations in R.

The only complication in the Turing model is the cost of extracting the f i , f i ′ and g i in the correct order. For example, consider the f i in the case that n is odd. Given as input a precomputed table 1, ω, ω 2 , ..., ω n-1 , we may extract the f i in time O(n lg n (b + lg n)) by applying the strategy from section 2.3 to the pairs (i, Remark 2.3. Bluestein's original formulation used the weights f i := ω i 2 /2 , g i := ω -i 2 /2 . This has two drawbacks in our setting. First, it requires the presence of 2 n-th roots of unity, which may not exist in R. Second, if n is odd, it leads to negacyclic convolutions, rather than cyclic convolutions. The variants used here avoid both of these problems.

Remark 2.4. It is also possible to give variants of the new multiplication algorithms in which Bluestein's transform is replaced by a different method for converting DFTs to convolutions, such as Rader's algorithm [START_REF] Rader | Discrete Fourier transforms when the number of data samples is prime[END_REF].

Kronecker substitution

Multiplication in F p [X] may be reduced to multiplication in Z using the classical technique of Kronecker substitution [START_REF] Zur Gathen | Modern Computer Algebra[END_REF]Corollary 8.27]. More precisely, let n > 0 and suppose that we are given two polynomials A, B ∈ F p [X] of degree less than n. Let A ˜, B ˜∈ Z[X] be lifts of these polynomials, with coefficients A i and B i satisfying 0 A i , B i < p. Then for the product C ˜= A ˜B ˜we have 0 C i < n p 2 , and the coefficients of C ˜may be read off the integer product

C ˜(2 N ) = A ˜(2 N ) B ˜(2 N )
where N := 2 lg p + lg n. We deduce the coefficients of C := A B by dividing each C i by p. This shows that M p (n) = O(I(n (lg p + lg n)) + n I(lg p)). Using the assumption that I(n) /n is increasing, we obtain

M p (n) = O(I(n (lg p + lg n))).
If we also assume that lg n = O(lg p), i.e., that the degree is not too large, then this becomes simply

M p (n) = O(I(n lg p)). (2.5)
A second type of Kronecker substitution reduces bivariate polynomial multiplication to the univariate case. Indeed, let n 1 and κ 1, and suppose that

A, B ∈ F p [X , Z], where deg X A, deg X B < n and deg Z A, deg Z B < κ. We may then recover C := A B from the univariate product C(Y 2κ , Y ) = A(Y 2κ , Y ) B(Y 2κ , Y ) in F p [Y ]. Note that A(Y 2κ , Y ) and B(Y 2κ , Y ) have degree less than 2 n κ, so the cost of the bivariate product is M p (2 n κ) + O(n κ lg p).
The same method works for computing cyclic convolutions: to multiply A, B ∈ F p [X , Z] / (X n -1), the same substitution leads to a product in

F p [Y ] / (Y 2nκ -1). The cost is thus M p ′ (2 n κ) + O(n κ lg p), where M p ′ (d) denotes the cost of a multiplication in F p [X]/(X d -1).

Arithmetic in finite fields

Let p be a prime and let κ 1. In this section we review basic results concerning arithmetic in F p κ and

F p κ[Y ].
We assume throughout that F p κ is represented as F p [Z]/ P for some irreducible monic polynomial P ∈ F p [Z] of degree κ. Thus an element of F p κ is represented uniquely by its residue modulo P , i.e., by a polynomial F ∈ F p [Z] of degree less than κ. Lemma 3.1. Let p be a prime and let κ 1. We may compute a monic irreducible polynomial

P ∈ F p [Z] of degree κ in time O ˜(κ 4 p 1/2 ). Proof. See [33, Theorem 3.2].
The above complexity bound is very pessimistic in practice. Better complexity bounds are known if we allow randomised algorithms, or unproved hypotheses. For instance, assuming GRH, the bound reduces to (κ log p) O(1) [START_REF] Adleman | Finding irreducible polynomials over finite fields[END_REF].

We now consider the cost of arithmetic in F p κ, assuming that P is given. Addition and subtraction in F p κ have cost O(κ lg p). For multiplication we will always use the Schönhage-Strassen algorithm. Denote by S p (κ) the cost of multiplying polynomials in F p [Z] κ by this method, i.e., using the polynomial variant [START_REF] Schönhage | Schnelle Multiplikation von Polynomen über Körpern der Charakteristik 2[END_REF] for the polynomials themselves, followed by the integer version [START_REF] Schönhage | Schnelle Multiplikation großer Zahlen[END_REF] to handle the coefficient multiplications. Then we have S p (κ) = O(κ lg κ lg lg κ lg p lg lg p lg lg lg p).

Of course, we could use the new multiplication algorithm recursively for these products, but it turns out that Schönhage-Strassen is fast enough, and leads to a simpler complexity analysis in section 7. Let D p (κ) denote the cost of dividing a polynomial in F p [Z] 2κ by P , returning the quotient and remainder. Using Newton's method [START_REF] Zur Gathen | Modern Computer Algebra[END_REF]Chapter 9] we have D p (κ) = O(S p (κ)). Thus elements of F p κ may be multiplied in time O(S p (κ)).

Let N 1 be a divisor of p κ -1. To compute a DFT over F p κ of length N , we must first have access to a primitive N -th root of unity in F p κ. In general it is very difficult to find a primitive root of F p κ, as it requires knowledge of the factorisation of p κ -1. However, to find a primitive N -th root, Lemma 3.3 below shows that it is enough to know the factorisation of N . The construction relies on the following existence result. Lemma 3.2. Let ℓ be a positive integer such that p ℓ > c 1 r 4 (log r + 1) 4 κ 2 , where r is the number of distinct prime divisors of p κ -1, and where c 1 > 0 is a certain absolute constant. Then there exists a monic irreducible polynomial Θ ∈ F p [Z] of degree ℓ such that Θ modulo P is a primitive root of unity for

F p κ = F p [Z]/P. Proof. This is [34, Theorem 1.1].
Lemma 3.3. Assume that N divides p κ -1 and that the factorisation of N is given. Then we may compute a primitive N-th root of unity in F p [Z]/P in time O ˜(κ 9 p).

Proof. Testing whether a given α ∈ F p κ is a primitive N -th root of unity reduces to checking that α N = 1 and α N /s = / 1 for every prime divisor s of N . According to Lemma 3.2, it suffices to apply this test to α := (Θ mod P ) (p κ -1)/N for Θ ∈ F p [Z] running over all monic polynomials of degree ℓ := ⌈log(c 1 r 4 (log r + 1) 4 κ 2 )/log p⌉, where r = O(κ lg p) is as in the lemma. The number of candidates is at most p ℓ c r 4 (log r + 1) 4 κ 2 p = O ˜(κ 6 p), and each candidate can be tested in time O ˜(r lg N (κ lg p)) = O ˜(κ 3 lg 3 p) using binary exponentiation.

Of course, we can do much better if randomised algorithms are allowed, since ζ (p κ -1)/N is reasonably likely to be a primitive N -th root for randomly selected ζ.

Finally, we consider polynomial multiplication over F p κ. Let n 1, and let A, B ∈ F p κ[X ] n . Let A ˜, B ˜∈ F p [X , Z] be their natural lifts, i.e., of degree less than κ with respect to Z. The bivariate product C ˜:= A ˜B ˜may be computed using Kronecker substitution (section 2.6) in time

M p (2 n κ) + O(n κ lg p). Writing C ˜= i=0 2n-2 C ˜i(Z) X i ,
M p κ(n) = M p (2 n κ) + O(n D p (κ)) = M p (2 n κ) + O(n S p (κ)).
As in section 2.6, the same method may be used for cyclic convolutions. If M p κ ′ (n) denotes the cost of multiplication in F p κ[X ]/(X n -1), then we get

M p κ ′ (n) = M p ′ (2 n κ) + O(n S p (κ)).
Remark 3.4. In the other direction, we can also reduce multiplication in

F p [X] to multiplication in F p κ[Y ]
by splitting the inputs into chunks of size ⌊κ / 2⌋. This leads to the bound M p (n) M p κ(⌈n / ⌊κ / 2⌋⌉) + O(n lg p) for n κ. A variant of this procedure is developed in detail in section 6.

Preparing for DFTs of large smooth orders

The aim of this section is to prove the following theorem, which allows us to construct "small" extensions F p λ of F p containing "many" roots of unity of "low" order. Recall that a positive integer is said to be y-smooth if all of its prime divisors are less than or equal to y.

Theorem 4.1. There exist computable absolute constants c 3 > c 2 > 0 and n 0 ∈ N with the following properties. Let p be a prime and let n n 0 . Then there exists an integer λ in the interval

(lg n) c2lg lg lg n < λ < (lg n) c3lg lg lg n ,
and a (λ + 1)-smooth integer M n, such that M | p λ -1. Moreover, given p and n, we may compute λ and the prime factorisation of M in time O((lg n) lg lg n ).

The proof depends on a series of preparatory lemmas. For λ 2, define

H λ := q prime, q -1|λ q.
For example,

H 36 = 1919190 = 2 • 3 • 5 • 7 • 13 • 19 • 37 and H 37 = 2.
Note that H λ is always squarefree and (λ + 1)-smooth. We are interested in finding λ for which H λ is unusually large; that is, for which λ has many divisors d such that d + 1 happens to be prime. Most of the heavy lifting is done by the following remarkable result.

Lemma 4.2. Define λ 0 (k) := min λ ∈ N: H λ k √
. There exist computable absolute constants c 5 > c 4 > 0 such that for any integer k > 100 we have

(log k) c4log log log k < λ 0 (k) < (log k) c5log log log k .
Proof. This is part of [2, Theorem 3]. (The threshold k √ could of course be replaced by any fixed power of k. It is stated this way in [START_REF] Adleman | On distinguishing prime numbers from composite numbers[END_REF] because that paper is concerned with primality testing.)

The link between H λ and F p λ is as follows.

Lemma 4.3. Let p be a prime and let λ 2. Then there exists a

(λ + 1)-smooth integer M H λ /(λ + 1) such that M | p λ -1.
Proof. We take M := H λ / p if p divides H λ , and otherwise M := H λ . In the former case we must have p λ + 1, so in both cases M H λ / (λ + 1). To see that M | p λ -1, consider any prime divisor q = / p of H λ . Then q -1 | λ, so q | p λ -1 by Fermat's little theorem.

Remark 4.4. The integer M constructed in Lemma 4.3 only takes into account the structure of H λ , and ignores p λ -1 itself. In practice, p λ -1 will often have small factors other than those in H λ , possibly including repeated factors (which are never detected by H λ ). For example,

H 6 = 2 • 3 • 7, but 19 6 -1 = 2 3 • 3 3 • 5 • 7 3 • 127.
In an implementation, one should always incorporate these highly valuable "accidental" factors into M . We will ignore them in our theoretical discussion.

Next we give a simple sieving algorithm that tabulates approximations of log H λ for all λ up to a prescribed bound.

|ℓ λ -log H λ | 1, (1 λ m). (4.1)
Proof. Initialise a table of ℓ λ with ℓ λ := 0 for λ = 1, ..., m. Since log H λ q λ+1 log q = O(λ), it suffices to set aside O(lg m) bits for each ℓ λ .

For each integer q = 2, ..., m, perform the following steps. First test whether q is prime, discarding it if not. Using trial division, the cost is q 1/2+o (1) , so m 3/2+o (1) = O(m 2 ) overall. Now assume that q is found to be prime. Using a fast algorithm for computing logarithms [START_REF] Brent | Fast multiple-precision evaluation of elementary functions[END_REF], compute an integer L q such that |L q -2 r log q| 1, where r := 1 + lg m, in time O (lg m) 1+o (1) . Scan through the table, replacing ℓ λ by ℓ λ + L q for those λ divisible by q -1, i.e., every (q -1)-th entry, in time O(m lg m). The total cost for each prime q is O(m lg m), and there are O(m/lg m) primes, so the overall cost is O(m 2 ). At the end, we divide each ℓ λ by 2 r , yielding the required approximations satisfying (4.1).

Proof of Theorem 4.1. We are given p and n as input. Applying Lemma 4.2 with k := n 3 , we find that for large enough n there exists some λ with (log (n 3 )) c4log log log (n 3 ) < λ < (log (n 3 )) c5log log log (n 3 ) and such that H λ n 3/2 . Choose any positive c 2 < c 4 log 2 and any c 3 > c 5 log 2. Since log log(n 3 ) log log log(n 3 ) ∼ (log 2) log lg n lg lg lg n, we have (lg n) c2lg lg lg n < (log(n 3 )) c4log log log(n 3 ) < λ < (log(n 3 )) c5log log log(n 3 ) < (lg n) c3lg lg lg n for large n. Using Lemma 4.5, we may find one such λ in time O(m 2 ), where m := (lg n) c3lg lg lg n .

Let M be as in Lemma 4.3. Then M divides p λ -1, and

M H λ λ + 1 n 3/2 m + 1 n
for large n. We may compute the prime factorisation of M by simply enumerating the primes q λ + 1, q = / p, and checking whether q -1 divides λ for each q. This can be done in time

O(λ 2 ) = O(m 2 ).
The main multiplication algorithm depends on a reduction of a "long" DFT to many "short" DFTs. It is essential to have some control over the long and short transform lengths. The following result packages together the prime divisors of the abovementioned M , to obtain a long transform length N near a given target L, and short transform lengths N i near a given target S.

Theorem 4.6. Let p, n, λ, M be as in Theorem 4.1. Let L and S be positive integers such that λ < S < L < M. Then there exist (λ + 1)-smooth integers N 1 , ..., N d with the following properties:

a) N := N 1 ••• N d divides M (and hence divides p λ -1). b) L N (λ + 1) L. c) S N i S 3 for all i.
Given λ, S, L, and the prime factorisation of M, we may compute such N 1 , ..., N d (and their factorisations) in time O ˜(λ 3 ).

Proof. Let M = N 1 ••• N s be the prime decomposition of M . Taking d minimal with N 1 ••• N d L,
we ensure that (a) and (b) are satisfied. At this stage we have N i λ + 1 S for all i. As long as the tuple (N 1 , ..., N d ) contains an entry N i with N i < S, we pick the two smallest entries N i and N j and replace them by a single entry N i N j . Clearly, this does not alter the product N of all entries, so (a) and (b) continue to hold. Furthermore, as long as there exist two entries N i and N j with N i < S and N j < S, new entries N i N j will always be smaller than S 2 . Only at the very last step of the loop, the second smallest entry N j might be larger than S. In that case, the product of the entry N i with N i < S with N j is still bounded by S • S 2 = S 3 . This shows that condition (c) is satisfied at the end of the loop.

Determining d requires at most s = O(λ) multiplications of integers less than M . There are at most s iterations of the replacement loop. Each iteration must scan through at most s integers of bit size O(log M ) = O(λ) and perform one multiplication on such integers.

Logarithmically slow recurrence inequalities

Let Φ: (x 0 , ∞) → R be a smooth increasing function, for some x 0 ∈ R. We say that Φ

* : (x 0 , ∞) → R is an iterator of Φ if Φ * is increasing and if Φ * (x) = Φ * (Φ(x)) + 1 (5.1)
for all sufficiently large x.

For instance, the standard iterated logarithm log * defined in (1.2) is an iterator of log. An analogous iterator may be defined for any smooth increasing function Φ: (x 0 , ∞) → R for which there exists some σ x 0 such that Φ(x) x -1 for all x > σ. Indeed, in that case,

Φ * (x) := min {k ∈ N: Φ •k (x) σ }
is well-defined and satisfies (5.1) for all x > σ. It will sometimes be convenient to increase x 0 so that Φ(x) x -1 is satisfied on the whole domain of Φ.

We say that Φ is logarithmically slow if there exists an ℓ ∈ N such that

(log •ℓ • Φ • exp •ℓ )(x) = log x + O(1) (5.2)
for x→ ∞. For example, the functions log (2 x), 2 log x, (log x) 2 and (log x) log log x are logarithmically slow, with ℓ = 0, 1, 2, 3 respectively.

Lemma 5.1. Let Φ: (x 0 , ∞) → R be a logarithmically slow function. Then there exists σ x 0 such that Φ(x) x -1 for all x > σ. Consequently all logarithmically slow functions admit iterators.

In this paper, the main role played by logarithmically slow functions is to measure size reduction in multiplication algorithms. In other words, multiplication of objects of size n will be reduced to multiplication of objects of size n ′ , where n ′ Φ(n) for some logarithmically slow function Φ(x).

The following result asserts that, from the point of view of iterators, such functions are more or less interchangeable with log x. Lemma 5.2. For any iterator Φ * of a logarithmically slow function Φ, we have

Φ * (x) = log * x + O(1).
The next result is our main tool for converting recurrence inequalities into actual asymptotic bounds for solutions.

Proposition 5.3. Let K > 1, B 0 and ℓ ∈ N. Let x 0 exp •ℓ (1)
, and let Φ: (x 0 , ∞) → R be a logarithmically slow function such that Φ(x) x -1 for all x > x 0 . Then there exists a positive constant C (depending on x 0 , Φ, K, B and ℓ) with the following property.

Let σ x 0 and L > 0. Let S ⊆ R, and let T : S → R be any function satisfying the following recurrence. First, T (y) L for all y ∈ S, y σ. Second, for all y ∈ S, y > σ, there exist y 1 , ..., y d ∈ S with y i Φ(y), and weights γ 1 , ..., γ d 0 with i γ i = 1, such that

T (y) K 1 + B log •ℓ y i=1 d γ i T (y i ) + L.
Then we have T (y) C L K log * y -log * σ for all y ∈ S, y > σ.

The Crandall-Fagin algorithm

Consider the problem of computing the product of two n-bit integers modulo m := 2 n -1. The naive approach is to compute their ordinary 2 n-bit product and then reduce modulo m. The reduction cost is negligible because of the special form of m. If n is divisible by a high power of two, one can save a factor of two by using the fact that FFTs naturally compute cyclic convolutions.

An ingenious algorithm of Crandall and Fagin [START_REF] Crandall | Discrete weighted transforms and large-integer arithmetic[END_REF] allows for the gain of this precious factor of two for arbitrary n. Their algorithm is routinely used in the extreme case where n is prime, in the large-scale GIMPS search for Mersenne primes (see http://www.mersenne.org/).

A variant of the Crandall-Fagin algorithm was a key ingredient of our algorithm for integer multiplication that conjecturally achieves the bound Section 9]. In this section we present yet another variant, for computing products in F p [X]/(X n -1).

I(n) = O(n log n 4 log * n ) [18,

The Crandall-Fagin reduction for polynomials

Let n, N and κ be positive integers with 1 N n and κ 2 ⌈n / N ⌉. Our aim is to reduce multiplication in F p [X] / (X n -1) to multiplication in F p κ[Y ] / (Y N -1). In the applications in subsequent sections, N will be a divisor of p κ -1 with many small factors, so that we can multiply efficiently in F p κ[Y ]/(Y N -1) using FFTs over F p κ. Suppose that F p κ is represented as F p [Z]/P , for some irreducible P ∈ F p [Z] of degree κ. The reduction relies on the existence of an element θ ∈ F p κ such that θ N = Z. (This is the analogue of the real N -th root of 2 appearing in the usual Crandall-Fagin algorithm, which was originally formulated over C.) It is easy to see that such θ may not always exist; for example, there is no cube root of Z in F 16 = F 2 [Z]/(Z 4 + Z + 1). Nevertheless, the next result shows that if κ is large enough, then we can always find some modulus P for which a suitable θ exists. Proposition 6.1. Suppose that p κ/2 > N. Then we may compute an irreducible polynomial P ∈ F p [Z] of degree κ, and an element θ ∈ F p κ = F p [Z]/P such that θ N = Z, in time O ˜(κ 9 p).

Proof. We first observe that if ζ ∈ F p κ is a primitive root, then ζ N cannot lie in a proper subfield of F p κ. (This property is independent of P .) Indeed, if ζ N ∈ F p ℓ for some proper divisor ℓ of κ, then every N -th power in F p κ lies in F p ℓ. This contradicts the fact that the number of N -th powers in F p κ is at least 1 + (p κ -1)/N p κ /N > p κ/2 p ℓ . Now we give an algorithm for computing P and θ. We start by using Lemma 3.1 to compute an irreducible P ˜∈ F p [U ] of degree κ in time O ˜ κ 4 p 1/2 , and we temporarily represent F p κ as F p [U ]/P ˜. Using Lemma 3.2, we may construct a list of O ˜(κ 6 p) candidates for primitive roots of F p κ (see the proof of Lemma 3.3). For each candidate ζ, we compute ζ N and its powers 1, ζ N , ..., ζ (κ-1)N . If these are linearly dependent over F p then ζ N belongs to a proper subfield (and so, as shown above, ζ is not actually a primitive root). Thus we must eventually encounter some ζ for which they are linearly independent. We take P to be the minimal polynomial of this ζ N , which may be computed by finding a linear relation among 1, ζ N , ..., ζ κN . Let ϕ: F p [Z]/P → F p [U ]/ P ˜be the isomorphism that sends Z to ζ N . The matrix of ϕ with respect to the standard bases 1, Z , ..., Z κ-1 and 1, U , ..., U κ-1 is given by the coefficients of 1, ζ N , ..., ζ (κ-1)N . The inverse of this matrix yields the matrix of ϕ -1 :

F p [U ] / P ˜→ F p [Z] / P . We then set θ = ϕ -1 (ζ), so that θ N = ϕ -1 (ζ N ) = Z.
For each candidate ζ, the cost of computing the necessary powers of ζ is O ˜((κ + lg N ) κ lg p) = O ˜(κ 2 lg 2 p), and the various linear algebra steps require time O ˜(κ 3 lg p) using classical matrix arithmetic.

In the remainder of this section, we fix some P and θ as in Proposition 6.1, and assume that F p κ is represented as F p [Z] / P . Suppose that we wish to compute the product of u, v ∈ F p [X] / (X n -1). The presentation here closely follows that of [18, Section 9.2]. We decompose u and v as

u = i=0 N -1 u i X ei , v = i=0 N -1 v i X ei , (6.1)
where

e i := ⌈n i/N ⌉, u i , v i ∈ F p [X] e i+1 -e i .
Notice that e i+1e i takes only two possible values: ⌊n/N ⌋ or ⌈n/N ⌉.

For 0 i < N , let

c i := N e i -n i, (6.2)
so that 0 c i < N . For any 0 i 1 , i 2 < N , define δ i 1 ,i 2 ∈ Z as follows. Choose σ ∈ {0, 1} so that i := i 1 + i 2σ N lies in the interval 0 i < N , and put

δ i 1 ,i 2 := e i 1 + e i 2 -e i -σ n.
From (6.2), we have

c i1 + c i2 -c i = N (e i1 + e i2 -e i ) -n (i 1 + i 2 -i) = N δ i1,i2 .
Since the left hand side lies in the interval (-N , 2 N ), this shows that δ i1,i2 ∈ {0, 1}. Now, since e i1 + e i2 = e i + δ i1,i2 (mod n), we have

u v = i1=0 N -1 i2=0 N -1 u i 1 v i 2 X ei 1 +ei 2 = i=0 N -1 w i X ei (mod X n -1),
where

w i := i1+i2=i (mod N ) X δi 1 ,i 2 u i1 v i2 .
Since u i1 ∈ F p [X] ⌈n/N ⌉ and similarly for v i2 , we have w i ∈ F p [X] 2⌈n/N ⌉ . Note that we may recover u v from w 0 , ..., w N -1 in time O(n lg p), by a standard overlap-add procedure (provided that N = O(n/lg n)). Now, regarding u i and v i as elements of F p κ by sending X to Z, define polynomials U , V ∈ F p κ[Y ] / (Y N -1) by U i := θ ci u i and V i := θ ci v i for 0 i < N , and let

W = W 0 + ••• + W N -1 Y N -1 := U V be their (cyclic) product. Then w ˜i := θ -ci W i = i 1 +i 2 =i (mod N ) θ -ci U i1 V i2 = θ ci 1 +ci 2 -ci u i1 v i2 = Z δi 1 ,i 2 u i1 v i2
coincides with the reinterpretation of w i as an element of F p κ. Moreover, we may recover w i unambiguously from w ˜i, as κ 2 ⌈n / N ⌉ and w i ∈ F p [X] 2⌈n/N ⌉ . Altogether, this shows how to reduce multiplication in

F p [X]/(X n -1) to multiplication in F p κ[Y ]/(Y N -1).
Remark 6.2. The pair (e i+1 , c i+1 ) can be computed from (e i , c i ) in O(lg n) bit operations, so we may compute the sequences e 0 , ..., e N -1 and c 0 , ..., c N -1 in time O(N lg n). Moreover, since c i+1c i takes on only two possible values, we may compute the sequence θ c0 , ..., θ cN -1 using O(N ) multiplications in F p κ.

The main algorithm

Consider the problem of computing t 1 products u 1 v, ..., u t v with u 1 , ..., u t , v ∈ F p [X]/(X n -1), i.e., t products with one fixed operand. Denote the cost of this operation by C p,t (n). Our algorithm for this problem will perform t + 1 forward DFTs and t inverse DFTs, so it is convenient to introduce the normalisation

C p (n) := sup t 1 C p,t (n) 2 t + 1 .
This is well-defined since clearly C p,t (n) t C p,1 (n). Roughly speaking, C p (n) may be thought of as the notional cost of a single DFT. The problem of multiplying two polynomials in F p [X] of degree less than k may be reduced to the above problem by using zero-padding, i.e., taking n := 2 k and t := 1. Since C p,1 (n) 3 C p (n), we obtain M p (k) 3 C(2 k) + O(k lg p). Thus it suffices to obtain a good bound for C p (n).

The next theorem gives the core of the new algorithm for the case that n is large relative to p.

Theorem 7.1. There exist x 0 2 and a logarithmically slow function Φ: (x 0 , ∞) → R with the following property. For all integers n > x 0 , there exist integers n 1 , ..., n d Φ(n), and weights γ 1 , ..., γ d 0 with i γ i = 1, such that

C p (n) n lg p lg (n lg p) 8 + O 1 lg lg n i=1 d γ i C p (n i ) n i lg p lg(n i lg p) + O(1), (7.1)
uniformly for n > max (x 0 , p 2 ).

Proof. We wish to compute t 1 products u 1 v, ..., u t v with u 1 , ..., u t , v ∈ F p [X] / (X n -1), for some sufficiently large n. We assume throughout that p 2 < n.

Choose parameters. Using Theorem 4.1, we obtain integers λ and M with

(lg n) c2lg lg lg n < λ < (lg n) c3lg lg lg n ,
and so that M n, M | p λ -1, and M is (λ + 1)-smooth. We choose long and short target transform lengths L := ⌈n/λ 3 ⌉ and S := λ (lg lg n) 2 . For large enough n we then have λ < S < L < n M , so we may apply Theorem 4.6. This yields (λ + 1)-smooth integers N 1 , ..., N d , with known factorisations, such that N := N 1 ••• N d divides p λ -1 and lies in the range L N (λ + 1) L, and such that S N i S 3 for all i. Finally we set κ := ⌈2 ⌈n/ N ⌉/λ⌉ λ, so that N also divides p κ -1, and we put n i := 2 N i κ. All of these parameters may be computed in time O((lg n) lg lg n ) = O(n).

For large n, we observe that the following estimates hold. First, we have

n λ 3 N 2 n λ 2 .
Since log λ ≍ lg lg n lg lg lg n, it follows that

log 2 N = 1 + O 1 lg lg n lg n,
and also that log κ ≍ lg lg n lg lg lg n. We have

N κ = 2 n + O(N λ) = 2 n + O(n/λ), so N κ = 2 + O 1 lg n n.
To estimate d, note that N S d and log S ≍ (lg lg n) 3 lg lg lg n, so

d = O lg n (lg lg n) 3 lg lg lg n .
Since p n we have Long transforms. The factorisation of N is known, and N divides p κ -1, so by Lemma 3.3 we may compute a primitive N -th root of unity ω ∈ F p κ in time O ˜(κ 9 p) = O(n).

lg (n lg p) = 1 + O 1 lg lg n lg n, ( 7
We will multiply in F p κ[Y ]/ (Y N -1) by using DFTs with respect to ω. The table of roots 1, ω, ..., ω N -1 may be computed in time O(N S p (κ)). In a moment, we will describe an algorithm A i for computing a "short" DFT of length N i with respect to ω i := ω N /Ni ; we then use the algorithm Let D denote the complexity of A, and for t ′ 1 let D i,t ′ denote the cost of performing t ′ independent DFTs of length N i using A i . Then by (2.4) we have

A := A 1 ⊙ ••• ⊙ A d for
D i=1 d D i,N /N i + O(d N S p (κ)) + O((κ lg p) N lg N ).
The last term is simply O(n lg n lg p).

Bluestein conversion. We now begin constructing A i , assuming that t ′ 1 independent transforms are sought. We first use Bluestein's algorithm (section 2.5) to convert each DFT of length N i to a multiplication in F p κ[X ]/(X Ni -1). We must check that 2 is invertible in

F p κ if N i is even; indeed, if N i is even, then so is p κ -1, so p = / 2. The Bluestein conversion contributes O(t ′ N i S p (κ)) to the cost of A i .
We must also compute a suitable table of roots, once at the top level. We first extract the table 1, ω i , ..., ω i Ni-1 from the top level table in time O(N κ lg p) = O(n lg p), and then sort them into the correct order (and perform any necessary additions) in time

O(N i lg N i (κ lg p + lg N )) = O(S 3 lg N i κ lg p) = O(lg N i (n lg p)). Over all i the cost is O(lg N (n lg p)) = O(n lg n lg p).
Kronecker substitution. We finally convert each multiplication in F p κ[X ]/(X Ni -1) to a multiplication in F p [X] / (X ni -1) using Kronecker substitution (see section 2.6). The latter multiplications have cost C p,t ′(n i ), since one argument is fixed. After the multiplications, we must also perform t ′ N i divisions by P to recover the results in

F p κ[X ]/(X Ni -1), at cost O(t ′ N i S p (κ)).
Consolidating the estimates for the Bluestein conversion and Kronecker substitution, we have

D i,t ′ C p,t ′(n i ) + O(t ′ N i S p (κ)). (2 t ′ + 1) C p (n i ) + O(t ′ N i S p (κ)).
For t ′ = N /N i N /S 3 ≻ lg lg n, this becomes

D i,N /Ni (2 + O(1/lg lg n)) (N /N i ) C p (n i ) + O(N S p (κ)).

Conclusion. Summing over i yields

D 2 + O 1 lg lg n i=1 d N N i C p (n i ) + O(d N S p (κ)) + O(n lg n lg p). Since d N S p (κ) = O lg n (lg lg n) 3 lg lg lg n N κ lg κ lg lg κ lg p lg lg p lg lg lg p = O n lg n lg p (lg lg lg n) 2 lg lg n = O(n lg n lg p) and N N i = 2 N κ 2 N i κ = 4 + O 1 lg n n n i , this becomes D 8 + O 1 lg lg n i=1 d n n i C p (n i ) + O(n lg n lg p).
Let γ i := log N i /log N , so that i γ i = 1. Since

lg n i = log 2 N i + O(log κ) = 1 + O log λ log S log 2 N i = 1 + O 1 lg lg n log 2 N i and lg n = (1 + O(1/lg lg n)) log 2 N , we get D 8 + O 1 lg lg n i=1 d γ i n lg n n i lg n i C p (n i ) + O(n lg n lg p).
To compute the desired t products, we must execute t + 1 forward transforms, and t inverse transforms. Each product also requires O(N ) pointwise multiplications in F p κ and O(N ) multiplications by 1/N . These have cost O(N S p (κ)), which is absorbed by the O(n lg n lg p) term. Thus we obtain

C p,t (n) (2 t + 1) D + O(t n lg n lg p).
Dividing by (2 t + 1) n lg n lg p, taking suprema over t 1, and using (7.2) and (7.3), yields the bound (7.1). Finally, since

n i = O(S 3 λ 3 ) = O (lg n) 3c3lg lg lg n((lg lg n) 2 +1) = O (lg n) (lg lg n) 3 ,
we have log n i = O((log log n) 4 ) and hence log log log n i log log log log n + C for some constant C and all large n. We may then take Φ(x) := exp •3 (log •4 x + C).

Now we may prove the main theorem announced in the introduction.

Proof of Theorem 1.1. For n 2 define

T p (n) := C p (n) n lg p lg(n lg p) .
It suffices to prove that T p (n) = O(8 log * (nlg p) ), uniformly in n and p. Let x 0 and Φ(x) be as in Theorem 7.1. Increasing x 0 if necessary, we may assume that Φ(x) x -1 for x > x 0 , and that x 0 exp (exp (1)).

Let σ p := max (x 0 , p 2 ) for each p. We will consider the regions n σ p and n > σ p separately. First consider the case n σ p . There are only finitely many primes p < (x 0 ) 1/2 , so we may assume that p 2 x 0 and that n σ p = p 2 . In this region we use Kronecker substitution: by (2.5) and (1.1) we have

M p (n) = O(I(n lg p)) = O n lg p lg (n lg p) 8 log * (nlg p)
, and since

C p (n) = sup t 1 C p,t (n) 2 t + 1 sup t 1 t C p,1 (n) 2 t + 1 = O(C p,1 (n)) = O(M p (n)),
we get T p (n) = O(8 log * (n lg p) ) for n σ p . In fact, this even shows that T p (n) = O(8 log * p ) uniformly for n σ p . Now consider the case n > σ p . Here we invoke Theorem 7.1 to obtain absolute constants B, L > 0 such that for every n > σ p , there exist n 1 , ..., n d Φ(n) and γ 1 , ..., γ d such that

T p (n) 8 1 + B log log n i=1 d γ i T p (n i ) + L. Set L p := max (L, max 2 n σp T p (n)). Applying Proposition 5.3 with ℓ := 2, K := 8, S := {2, 3, ...}, we find that T p (n) = O(L p 8 log * n-log * σp ) uniformly for n > σ p . Since log * σ p = log * p + O(1) and L p = O(8 log * p ), we conclude that T p (n) = O(8 log * n ) = O(8 log * (nlg p) ) uniformly for n > σ p .

Conjecturally faster multiplication

Recall that in [START_REF] Harvey | Even faster integer multiplication[END_REF] we established the bound I(n) = O(n lg n K log * n ) with K = 8 for the complexity of integer multiplication. We also proved that this can be improved to K = 4, assuming the following slight weakening of the Lenstra-Pomerance-Wagstaff conjecture on the distribution of Mersenne primes [START_REF] Pomerance | Recent developments in primality testing[END_REF][START_REF] Wagstaff | Divisors of Mersenne numbers[END_REF] (see [START_REF] Harvey | Even faster integer multiplication[END_REF]Section 9] for further discussion).

Conjecture 8.1. Let π m (x) be the number of Mersenne primes less than x. Then there exist constants 0 < a < b such that for all x > 3,

a log log x < π m (x) < b log log x.
The source of the (conjectural) speedup is as follows. The K = 8 algorithm of [START_REF] Harvey | Even faster integer multiplication[END_REF]Section 6] computes DFTs over C, and so we encounter the "short product obstruction". Namely, to compute the product of two real numbers with p significant bits using FFT algorithms, we are forced to compute a full product of two p-bit integers, and then truncate the result to p bits. To achieve K = 4, we replaced the coefficient ring C by the finite field F p [i], where p = 2 q -1 is a Mersenne prime. Multiplication modulo 2 q -1 is a "cyclic product" rather than a short product, and this saves a factor of two at each recursion level.

The aim of this section is to outline a credible strategy for achieving the same improvement, from K = 8 to K = 4, in the context of multiplication in F p [X], at least under certain plausible number-theoretic hypotheses.

Consider the problem of multiplying polynomials in F p [X] of degree less than n. In the Kronecker substitution region, i.e., for lg n = O(lg p), we can of course achieve K = 4 if we assume Conjecture 8.1. Now consider the region where n is much larger than p. In the algorithm of section 7, we reduced the multiplication problem to computing products in F p κ[Y ]/ (Y Ni -1), where F p κ = F p [Z]/P for some monic irreducible P ∈ F p [Z] of degree κ, and where N i is a "short" transform length. These multiplications were in turn converted to multiplications in F p [X] / (X 2κNi -1) via Kronecker substitution. It is exactly this factor of two in zero-padding that we wish to avoid.

Taking our cue from the integer case, we observe that if the modulus P is of a particularly special form, then this factor of two can be eliminated. Indeed, suppose that

P (Z) = Z κ+1 -1 Z -1 = Z κ + ••• + Z + 1. Then F p κ[Y ] / (Y Ni -1) is isomorphic to F p [Y , Z] / (Y Ni -1, Z κ + ••• + 1), which is a quotient of F p [Y , Z] / (Y Ni -1, Z κ+1 -1)
. Assuming further that gcd(N i , κ + 1) = 1, the latter ring is isomorphic to F p [X]/ (X (κ+1)Ni -1). In other words, we can reduce a cyclic convolution over F p κ of length N i to a cyclic convolution over F p of length (κ + 1) N i ∼ κ N i , rather than length 2 κ N i . This saves a factor of two at this recursion level; if we can manage something similar at every recursion level, we reduce K from 8 to 4.

The essential question is therefore how to choose κ so that Proof. If ℓ is a nontrivial factor of κ + 1, then Z κ+1 -1 is divisible by Z ℓ -1, so P has the nontrivial factor Z ℓ-1 + ••• + Z + 1. On the other hand, suppose that κ + 1 is prime. If κ + 1 = p, then P = (Z p -1) / (Z -1) = (Z -1) p-1 is not irreducible; otherwise, the number of irreducible factors of P over F p is exactly κ / m, where m is the order of p in (Z/ (κ + 1) Z) * [36, Theorem 2.13 and Proposition 2.14] (this last statement is sometimes called the "cyclotomic reciprocity law").

P = Z κ + ••• + Z + 1
A first attempt to reach K = 4 might run as follows. First choose λ as in the proof of Theorem 7.1, so that F p λ has many roots of unity of smooth order. Then use the same multiplication algorithm as before, but now working over

F p κ = F p [Z] / (Z κ + ••• + 1)
, where κ is the smallest positive multiple of λ such that q := κ + 1 is prime and such that p is a primitive root modulo q.

Unfortunately this plan does not quite work, for reasons that we now explain. First, if p | λ, there may be a trivial obstruction to the existence of suitable κ. For example, if p = 5 and 5 | λ, then for any q = 1 (mod λ) we have (5/ q) = (q /5) = 1 by quadratic reciprocity, so 5 cannot be a primitive root modulo q. The only way to avoid this seems to be to insist that λ not be divisible by p in the first place. Therefore we propose the following strengthening of Lemma 4.2.

Conjecture 8.3. Define λ 0 (k, p) := min λ ∈ N: H λ k √ , p ∤ λ p odd , min λ ∈ N: H λ k √ , 8 ∤ λ p = 2.
There exist computable absolute constants c 5 ′ > c 4 ′ > 0 such that for any integer k > 100 and any prime p we have (log k) c 4 ′ log log log k < λ 0 (k, p) < (log k) c 5 ′ log log log k .

It seems likely that this conjecture should be accessible to specialists in analytic number theory. Experimentally, the "missing" factor p does not seem to have much effect on the propensity of H λ to have many divisors d such that d + 1 is prime. We could not find a quick way to prove the conjecture directly from Lemma 4.2.

The second problem is more serious. A special case of Artin's conjecture on primitive roots (see [START_REF] Moree | Artin's primitive root conjecture-a survey[END_REF] for a survey) asserts that for any prime p, there are infinitely many primes q for which p is a primitive root modulo q. Unfortunately, there is not even a single prime p for which Artin's conjecture is known to hold! Therefore, in general we cannot prove existence of a suitable κ.

However, we do have the following result of Lenstra, which states that the only obstruction to the existence of κ is the trivial one noted above, provided we are willing to accept GRH. Lemma 8.4. Assume GRH. Assume that p ∤ λ if p is odd, or that 8 ∤ λ if p = 2. Then there exist infinitely many primes q = 1 (mod λ) such that p is a primitive root modulo q.

Proof. This is a special case of [START_REF] Lenstra | On Artin's conjecture and Euclid's algorithm in global fields[END_REF]Theorem 8.3]. (For a description of exactly which number fields are supposed to satisfy GRH for this result to hold, see [START_REF] Lenstra | On Artin's conjecture and Euclid's algorithm in global fields[END_REF]Theorem 3.1].)

The third problem is that existence of κ is not enough; we also need a bound on its size, as a function of λ and p. The literature provides little guidance on this question. Under the hypotheses of Lemma 8.4 (including GRH), it is known that the number of primes q < x such that q = 1 (mod λ) and for which p is a primitive root modulo q grows as

C x log x + O x log log x log 2 x (8.1)
for some positive C (see for example [START_REF] Moree | On primes in arithmetic progression having a prescribed primitive root[END_REF]Theorem 2]). Both C and the implied big-O constant depend on λ and p. While the value of C is reasonably well understood, we could not find in the literature a similarly precise description of the implied constant, so we have been unable to derive a bound for κ.

In the interests of getting our conjectural algorithm off the ground, we make the following guess. Define q 0 (λ, p) to be the smallest prime q such that q = 1 (mod λ) and such that p is a primitive root modulo q (assuming that such q exists).

Conjecture 8.5. There exists a constant β ∈ N and a logarithmically slow function Ψ(x) with the following property. Suppose that p ∤ λ if p is odd, or that 8 ∤ λ if p = 2, and suppose that log * λ β + log * p. Then q 0 (λ, p) exp Ψ(λ).

We now apply the Crandall-Fagin algorithm to reduce multiplication in F p [X] / (X n -1) to multiplication in F p κ[Y ]/(Y µN -1). The prerequisite κ 2 ⌈n/ µ N ⌉ is satisfied (since κ is even). Instead of using Proposition 6.1 to construct θ ∈ F p κ satisfying θ µN = Z, we may simply take θ := Z (µN ) -1 mod (κ+1) . The required modular inverse exists, as µ and κ + 1 are distinct primes, and gcd(N , κ + 1) = 1 since N is (λ + 1)-smooth.

We now take advantage of the isomorphism

F p κ[Y ] / (Y µN -1) = ∼ R[Y ] / (Y N -1)
, where R := F p κ[U ]/(U µ -1). We multiply in R[Y ]/(Y N -1) by using DFTs over R: one DFT of length N over R reduces to µ DFTs of length N over F p κ. The latter are handled by decomposing into short transforms of length N i , which are in turn converted to multiplications in F p κ[Y ]/(Y Ni -1) via Bluestein's algorithm. Finally -and this is where all the hard work pays off -each such multiplication reduces to a multiplication in F p [X]/(X (κ+1)Ni -1), since gcd(N i , κ + 1) = 1. The multiplications in R itself are handled using the algebraic Schönhage-Strassen algorithm.

We omit the rest of the complexity argument, which is essentially the same as that of Theorem 7.1 and Theorem 1.1. We mention only that S was chosen sufficiently large that the multiplications in R and in F p κ make a negligible contribution overall.

Notes and generalisations

In this section we outline some directions along which the results in this paper can be extended. We also provide some hints concerning the practical usefulness of the new ideas. Our treatment is more sketchy and we plan to provide more details in a forthcoming paper.

Multiplication of polynomials over F p κ

Recall that M p κ(n) denotes the cost of multiplying polynomials in F p κ[X ] of degree less than n, where we assume that some model F p [Z]/P for F p κ has been fixed in advance. Indeed, we saw in section 3 that

M p κ(n) = M p (2 n κ) + O(n D p (κ)),
where D p (κ) denotes the cost of dividing a polynomial of degree less than 2 κ by P . Having established the bound M p (n) = O(n lg p lg(n lg p) 8 log * (nlg p) ), it is now permissible to assume that M p (n)/n is increasing, so the usual argument for Newton iteration shows that D p (κ) = O(M p (κ)). Using again that M p (n)/n is increasing, we obtain n D p (κ) = O(M p (n κ)). For lg n = O(α lg p) we may simply use Kronecker substitution. For α lg p = O(lg n) we must modify the algorithm of section 7. Recall that in that algorithm we reduced multiplication in

Multiplication of polynomials over

F p [X] / (X n -1) to multiplication F p κ[Y ] / (Y N -1)
, where N is a transform length dividing p κ -1. Our task is to define a ring R, analogous to F p κ, so that multiplication in (Z/ p α Z)[X]/(X n -1) can be reduced to multiplication in R[Y ]/(Y N -1).

Let π be the natural projection Z / p α Z → F p , and write π also for the corresponding map

(Z / p α Z)[Z] → F p [Z]. Let P ∈ F p [Z]
be any monic irreducible polynomial of degree κ, and let P ˜∈ (Z/ p α Z)[Z] be an arbitrary lift via π, also monic of degree κ. Then P ˜is irreducible, and we will take R := (Z/ p α Z)[Z]/P ˜. Any primitive N -th root of unity ω ∈ F p κ has a unique lift to a principal N -th root of unity ω ˜∈ R. (This can be seen, for example, by observing that ω ˜must be the image in R of a Teichmüller lift of ω in the p-adic field whose residue field is F p κ.) Given ω, we may compute ω ˜efficiently using fast Newton lifting [START_REF] Cohen | Handbook of elliptic and hyperelliptic curve cryptography[END_REF]Section 12.3].

Moreover, if we choose P and θ ∈ F p κ = F p [Z] / P as in Proposition 6.1, so that θ N = Z, then fast Newton lifting can also be used to obtain θ ˜∈ R such that θ ˜N = Z in R.

Multiplication of polynomials over Z/ m Z

Complexity bounds for straight-line programs

Until now we have considered only complexity bounds in the Turing model. The new techniques also lead to improved bounds in algebraic complexity models. In what follows, A is always a commutative ring with identity.

For the simplest case, first assume that A is an F p -algebra for some prime p. In the straightline program model [9, Chapter 4], we count the number of additions and subtractions in A, the number of scalar multiplications (multiplications by elements of F p ), and the number of nonscalar multiplications (multiplications in A). As pointed out in the introduction, the best known bound for the total complexity was previously O(n lg n lg lg n), which was achieved by the Cantor-Kaltofen algorithm.

Theorem 9.4. Let A be an F p -algebra. We may multiply two polynomials in A[X] of degree less than n using O(n lg n 8 log * n ) additions, subtractions, and scalar multiplications, and O(n 4 log * n ) nonscalar multiplications. These bounds are uniform over all primes p and all F p -algebras A.

The idea of the proof is to use the same algorithm as in section 7, but instead of switching to Kronecker substitution when we reach n ≈ p, we simply recurse all the way down to n = 1. The role of the extension F p κ is played by A ⊗ Fp F p κ = ∼ A[Z] / P , where P ∈ F p [Z] is monic and irreducible of degree κ (here we have used the fact that F p [X] may be viewed as a subring of A[X], since A contains an identity element and hence a copy of F p ). Thus we reduce multiplication in A[X] / (X n -1) to multiplication in (A ⊗ Fp F p κ)[Y ]/(Y N -1). The latter multiplication may be handled by DFTs over A ⊗ F p F p κ, since any primitive N -th root of unity in F p κ corresponds naturally to a principal N -th root of unity in A ⊗ Fp F p κ.

The O(n lg n 8 log * n ) bound covers the cost of the DFTs, which are accomplished entirely using additions, subtractions and scalar multiplications. Nonscalar multiplications are needed only for the pointwise multiplication step. To explain the O(n 4 log * n ) bound, we observe that at each recursion level the total "data size" grows by a factor of four: one factor of two arises from the Crandall-Fagin splitting, and another factor of two from the Kronecker substitution. Faster polynomial multiplication over finite fields Note that in the straight-line program model, we give a different algorithm for each n. Thus all precomputed objects, such as the defining polynomial P and the required roots of unity, are obtained free of cost (they are built directly into the structure of the algorithm). The uniformity in p follows from the fact that the bounds for λ and M in Theorem 4.1 are independent of p, although of course the algorithm will be different for each p. For each fixed p, we may use the same straight-line program for all F p -algebras A.

Assuming Conjectures 8.3 and 8.5, the bounds may be improved to respectively O(n lg n 4 log * n ) and O(n 2 log * n ). However, in this case we lose uniformity in p, due to the issue raised in Remark 8.6. Theorem 9.4 can be generalised to (Z/ p α Z)-algebras, along the same lines as section 9.2. It is also possible to handle (Z/m Z)-algebras for any integer m 1, i.e., rings of finite characteristic, but we cannot proceed exactly as in section 9.3, because the straight-line model has no provision for a "reduction modulo p α " operation. Instead, we use a device introduced in [START_REF] Cantor | On fast multiplication of polynomials over arbitrary algebras[END_REF]. Suppose that m = p 1 α1 ••• p l αl and that A is now a (Z/m Z)-algebra. For each i we may construct a straight-line program M i that takes as input polynomials f , g ∈ B[X] of degree less than n, where B is any (Z / p i αi Z)-algebra, and computes f g. By replacing every constant in Z / p i αi Z by a compatible constant in Z/ m Z, we obtain a straight-line program M i ′ that takes as input f , g ∈ A[X] of degree less than n, and computes h i ∈ A[X] such that h if g ∈ p i αi A[X]. By the Chinese remainder theorem we may choose e i ∈ Z/m Z such that e i = 0 (mod m/ p i αi ) and e i = 1 (mod p i αi

) for each i. The linear combination i e i h i is then equal to f g in A[X]. We conclude that we may multiply polynomials in A[X] using O(n lg n 8 log * n ) additions, subtractions and scalar multiplications (by elements of Z/m Z), and O(n 4 log * n ) nonscalar multiplications (i.e., multiplications in A). These bounds are not uniform in m, but for each m they are uniform over all (Z/m Z)-algebras.

Other algebraic models

If we wish to take into account the cost of precomputations, and give a single algorithm that works uniformly for all p and n, we may use a more refined complexity model such as the BSS model [START_REF] Blum | Complexity and real computation[END_REF]. In this model, a machine over A is, roughly speaking, a Turing machine in which the tape cells hold elements of A. Actually, we need a multi-tape version of the model described in [START_REF] Blum | Complexity and real computation[END_REF]. The machine can perform arithmetic operations on elements of A in unit time, but can also manipulate data such as index variables in the same way as a Turing machine, and must deal with data layout in the same way as a Turing machine. In this model we obtain similar bounds to Theorem 9.4, but without uniformity in p. Alternatively, we could obtain bounds uniform in p if we added extra terms to account for the precomputations.

Another point of view in Theorem 9.4 is that we have described a new evaluation-interpolation strategy for polynomials over an F p -algebra A. We refer the reader to [START_REF] Van Der Hoeven | Newton's method and FFT trading[END_REF] for classical examples of evaluation-interpolation schemes, and [START_REF] Gao | Additive fast Fourier transforms over finite fields[END_REF] for algorithms specific to finite fields. Such schemes are characterised by two quantities: the evaluation/interpolation complexity E(n) and the number N(n) of evaluation points. The new algorithms yield the bounds E(n) = O(n lg n 8 log * n ) and N(n) = O(n 4 log * n ), and these can be used to prove complexity bounds for problems more general than polynomial multiplication. For example, we can multiply r × r matrices with entries in A[X] n using O(r 2 n lg n 8 log * n + r ω n 4 log * n ) ring operations, where ω is an exponent of matrix multiplication. One of the main advantages of our algorithms is that N(n) is almost linear, contrary to synthetic FFT methods [START_REF] Cantor | On fast multiplication of polynomials over arbitrary algebras[END_REF][START_REF] Schönhage | Schnelle Multiplikation von Polynomen über Körpern der Charakteristik 2[END_REF] derived from Schönhage-Strassen multiplication [START_REF] Schönhage | Schnelle Multiplikation großer Zahlen[END_REF], which achieve only N(n) = O(n lg n).

In the setting of bilinear complexity [START_REF] Bürgisser | Algebraic complexity theory[END_REF]Chapter 14], the new algorithms do not improve asymptotically on the best known bounds. For example, it is known [START_REF] Pieltant | New uniform and asymptotic upper bounds on the tensor rank of multiplication in extensions of finite fields[END_REF] that for any n there exist F 2 -linear maps a i , b i : F 2 [X] n → F 2 and polynomials c i ∈ F 2 [X] 2n for i ∈ {1, ..., k(n)}, with k(n) = (189 / 22 + o(1)) n, such that u v = i=1 k(n) a i (u) b i (v) c i for all u, v ∈ F 2 [X] n . The new method yields the asymptotically inferior bound k(n) = O(n 4 log * n ). Nevertheless, these bounds are asymptotic and do not take into account that our new algorithm would rather work over (say) F 2 60 instead of F 2 . In practice, the new method might therefore outperform the bilinear algorithms from [START_REF] Pieltant | New uniform and asymptotic upper bounds on the tensor rank of multiplication in extensions of finite fields[END_REF].

Theorem 1 . 1 .

 11 The boundM p (n) = O n log p log(n log p) 8 log * (nlog p)

Theorem 1 . 2 .

 12 Assume Conjectures 8.1, 8.3 and 8.5. Then the bound M p (n) = O n log p log(n log p) 4 log * (nlog p)

  ) in time O(n lg n (b + lg n)). The total cost of the extraction is thus O(n lg n (b + lg n)).

  to recover A B we must divide each C ˜i by P . Denoting by M p κ(n) the complexity of the original multiplication problem, we obtain

Lemma 4 . 5 .

 45 Let m 2. In time O(m 2 ), we may compute integers ℓ 1 , ..., ℓ m with

. 2 )

 2 and as n i N i S (lg n) lg lg n , we similarly obtain lg (n i lg p) = 1 + O 1 lg lg n lg n i . (7.3) Crandall-Fagin reduction. Let us check that the hypotheses of section 6.1 are satisfied, to enable the reduction to multiplication in F p κ[Y ] / (Y N -1). We certainly have 1 N n and κ 2 ⌈n / N ⌉. For Proposition 6.1, observe that κ λ 2 ≻ lg n ∼ lg N , so p κ/2 2 κ/2 > N for large n. Thus we obtain an irreducible P ∈ F p [Z] of degree κ, and an element θ ∈ F p κ = F p [Z]/P such that θ N = Z, in time O ˜(κ 9 p) = O ˜(λ 27 n 1/2 ) = O(n). Each multiplication in F p κ has cost O(S p (κ)) (see section 3). Computing the sequences e i and c i costs O(N lg n) = O(n lg n), and computing the sequence θ ci costs O(N S p (κ)). The initial splitting and final overlap-add phases require time O(t n lg p), and the multiplications by θ ci and θ -ci cost O(t N S p (κ)).

  the main transform of length N (see section 2.3). The corresponding twiddle factor tables may be extracted in time O(N lg N (κ lg p + lg N )) = O(n lg n lg p).

Theorem 9 . 1 .

 91 We have M p κ(n) = O n κ lg p lg (n κ lg p) 8 log * (nκlg p) , uniformly for all n, κ 1 and all primes p. Assuming Conjectures 8.1, 8.3 and 8.5, we may replace K = 8 by K = 4.

  Z/ p α Z For any prime p and any integer α 1, denote by M p,α (n) the bit complexity of multiplying polynomials in (Z/ p α Z)[X] n . Theorem 9.2. We have M p,α (n) = O n α lg p lg (n α lg p) 8 log * (nαlg p) , uniformly for all n, α 1 and all primes p. Assuming Conjectures 8.1, 8.3 and 8.5, we may replace K = 8 by K = 4.

For any m 1 ,

 1 denote by M m (n) the bit complexity of multiplying polynomials in (Z/m Z)[X] n . Theorem 9.3. We have M m (n) = O n lg m lg (n lg m) 8 log * (nlg m) , uniformly for all n, m 1. Assuming Conjectures 8.1, 8.3 and 8.5, we may replace K = 8 by K = 4. We use the isomorphism Z / m Z = ∼ i (Z / p i αi Z), where m = p 1 α1 ••• p l αl is the prime decomposition of m. The cost of converting between the Z / m Z and i (Z / p i αi Z) representations is O(I(lg m) lg l) = O(I(lg m) lg lg m) [17, Section 10.3]. By Theorem 9.2 we get M m (n) = O n lg m lg (n lg m) 8 log * (nlg m) + n I(lg m) lg lg m . The first term dominates if n m. If n < m we may simply use Kronecker substitution.

  see section 2.1); then merge with the given root table to obtain a table (i 1 , i 2 , ω i1i2 ), ordered by i 1 i 2 , in time O(n (b + lg n)); and finally sort by (i 1 , i 2

  is irreducible over F p . The following lemma gives a simple number-theoretic characterisation of such κ. Lemma 8.2. Let κ 2. The polynomial P = Z κ + ••• + Z + 1 is irreducible over F p if and only if κ + 1 is prime and p is a primitive root modulo κ + 1.

David Harvey, Joris van der Hoeven, Grégoire Lecerf

(i 2i) mod n) for 0 i < n. The other sequences are handled similarly. For the g i in the case that n is even, we also need to perform O(n) additions in R; assuming that additions in R have cost O(b), this is already covered by the above bound.

Faster polynomial multiplication over finite fields

Acknowledgments. We would like to thank Karim Belabas for directing us to the paper [2], and Igor Shparlinski for answering some analytic number theory questions.

We have deliberately given an absurdly weak formulation of this conjecture, so as to provide the largest possible target for analytic number theorists. For example, it would be enough to prove that q 0 (λ, p) exp ((log λ) 100 ) for all λ exp (exp (exp (exp p))). We suspect that the conjecture is accessible under GRH; presumably a proof would require analysing the implied constant in (8.1).

Remark 8.6. The reason we need λ to be large compared to p is that it is easy to construct, using the Chinese remainder theorem, a prime p which fails to be a primitive root modulo q for all primes q up to a prescribed bound.

Proof of Theorem 1.2. The proof is similar to that of Theorem 1.1. We will only give a sketch, highlighting the main differences.

In the region lg n = O(lg p), we use Kronecker substitution together with [18, Theorem 2], which states that I(n) = O(n lg n 4 log * n ) under Conjecture 8.1. Now assume that lg p = O(lg n). We will pursue a strategy similar to Theorem 7.1, recursing from n down to an exponentially smaller Φ(n). Assume that we wish to compute t 1 products in F p [X]/(X n -1) with one fixed operand.

Using an appropriate modification of Theorem 4.1, in which the use of Lemma 4.2 is replaced by Conjecture 8.3, we find integers λ and M such that

Now we wish to apply Conjecture 8.5 to construct suitable κ. However, that conjecture requires that log * λ β + log * p. If we have instead log * λ < β + log * p, we simply use the algorithm of Theorem 7.1. This yields the expansion factor K = 8, but only for O(1) recursion levels (depending on β), since log * λ = log * n + O(1). So it is permissible to assume that log * λ β + log * p, losing only a constant factor in the main complexity bound for the last few recursion levels.

Applying Conjecture 8.5 and Lemma 8.2, we obtain a positive multiple κ of λ, such that κ + 1 is prime, κ < exp Ψ(λ), and

Notice that κ < Ψ ′ (n) for some logarithmically slow function Ψ ′ (x). Indeed, there exists ℓ ∈ N and C , C ′ > 0 such that Ψ(λ) exp •ℓ (log •(ℓ+1) (λ) + C) and log log log λ log log log log n + C ′ for large n. Increasing ℓ if necessary, we get κ < exp Ψ(λ) exp •(ℓ+1) (log •(ℓ+2) (n) + C ′′ ), the latter being a logarithmically slow function of n.

In the proof of Theorem 7.1, we selected the transform length N first, and then chose κ to fine-tune the total bit size. Here we have less control over κ, so we must use a different strategy. Choose target long and short transform lengths L := ⌈n/(λ 2 κ)⌉ and S := κ (lg lg n) 2 lg lg κ . Applying Theorem 4.6, we obtain (λ + 1)-smooth integers N 1 , ..., N d such that

for each i, and such that n / (λ 2 κ) N 2 n/(λ κ) for large n. Let µ ′ := ⌈2 n/ (N κ)⌉ λ, and let µ be the smallest prime greater than µ ′ + 1 and different to κ + 1. By [START_REF] Baker | The difference between consecutive primes[END_REF] we have

Implementations

The new complexity bounds have a quite theoretic flavour: in practice, slowly growing functions such as log log n and 8 log * n really behave as constants. Furthermore, for concrete applications, it would be more relevant to redefine log * as log * n := min {k ∈ N:

With this definition, we observe that log * n 1 for all feasible values of n in the foreseeable future. Nevertheless, we think that variants of the ideas in this paper could still be useful in practice.

Let us briefly sketch how an efficient multiplication algorithm over F 2 could be implemented. Various modern processors offer a special instruction for the multiplication of two polynomials in F 2 [X] 64 . The practical translation of the results from section 4 is to work over the field F 2 60 = F 2 [Z] / P where P := Z 60 + ••• + Z + 1. This has two advantages. First, if we represent elements of F 2 60 by their residues modulo Z 61 -1, i.e., using a redundant representation, then each multiplication in F 2 60 requires a single hardware multiply instruction, followed by simple shift and XOR instructions. Second, there exist primitive roots of unity of high smooth order

Transforms of small prime lengths p | N over F 2 60 be implemented using specialised codelets, to those implemented in FFTW3 [START_REF] Frigo | The design and implementation of FFTW3[END_REF]. Roughly speaking, the constant factor of a radix-p DFT is p / log p. In our case, (13 / log 13) / (2 / log 2) ≈ 1.8, so specialised codelets should remain reasonably efficient for p 13. Consequently, through the mere use of such small radices, we obtain efficient algorithms for computing DFTs of lengths dividing 3 2 • 5 2 • 7 • 11 • 13 = 225225. To handle larger primes, we may use any algorithm for converting DFTs to convolutions. Rather than use Bluestein's algorithm as in the theoretical portion of the paper, it is more economical in this case to use Rader's algorithm, which reduces a DFT of length p to a cyclic convolution of length p -1, together with O(p) additions in F 2 60 (see also Remark 2.4).

More precisely, consider the remaining transform lengths p = 31, 41, 61, 151, 331 and 1321. We have 30

so in each case, the DFT reduces to a two-dimensional convolution of size 2 i × d, where d is a divisor of 225225 and 1 i 3. Such a convolution in turn reduces to 2 i+1 DFTs of length d together with d convolutions of length 2 i . The latter cannot be handled using DFTs over F 2 60, since there are no roots of unity of order 2, 4 or 8, but in any case they are all very short convolutions and could be performed directly, or by using Karatsuba's algorithm.

Transforms of other lengths dividing N may be reduced to these cases via the Cooley-Tukey algorithm. Finally, to handle arbitrary input lengths efficiently, we may use a mixed-radix generalisation of the truncated Fourier transform [START_REF] Van Der Hoeven | The truncated Fourier transform and applications[END_REF][START_REF] Van Der Hoeven | Notes on the Truncated Fourier Transform[END_REF].

We also reemphasise the fact that DFTs of this kind can be used as an evaluation-interpolation technique. This approach is for instance attractive for the multiplication of polynomial matrices over F 2 .

We have not yet implemented any of these algorithms. We expect that much fine-tuning will be necessary to make them most effective in practice. We intend to report on this issue in a future work.

Bibliography