Even faster integer multiplication

Polynomial multiplication can be done in subquadratic time using an evaluation-interpolation strategy: namely, to multiply two polynomials of degree less than n, we evaluate both polynomials at 2n points, compute 2n pointwise products, and recover the polynomial product by interpolation.

When one chooses roots of unity as the evaluation points, the evaluation and interpolation steps reduce to discrete Fourier transforms (forward and inverse, respectively), which can be done in O(n log n) arithmetic operations by means of the fast Fourier transform (the pointwise products require O(n) arithmetic operations).

Now, there are a couple of difficulties. How do you decompose the integers into polynomials? Is it better to use short polynomials with large coefficients, or long polynomials with small coefficients? Since the coefficients asymptotically grow with the length, you have to account for the bit complexity of the arithmetic operations in the FFTs and pointwise products (here you get recursive integer multiplications).

Which FFT algorithm should you use (there are many to choose from)?

Also, the ring of integers does not have enough roots of unity to do FFTs, so you have to lift the computation to a larger ring, such as the field of complex numbers (with numerical approximations), an appropriate finite field, or (as in the Schönhage-Strassen algorithm) a polynomial extension ring. Which choice gives the best complexity?

The improvements to integer multiplication published in the last few decades have basically been tightened analyses of such tradeoffs. I have not read this paper in detail yet, but it seems to provide a clear description in the introduction. Quote:

The idea of the new algorithm is remarkably simple. Given two n-bit integers, we split them into chunks of exponentially smaller size, say around log n bits, and thus reduce to the problem of multiplying integer polynomials of degree O(n/log n) with coeffcients of bit size O(log n). We multiply the polynomials using discrete Fourier transforms (DFTs) over C, with a working precision of O(log n) bits. To compute the DFTs, we decompose them into "short transforms" of exponentially smaller length, say length around log n, using the Cooley-Tukey method. We then use Bluestein's chirp transform to convert each short transform into a polynomial multiplication problem over C, and finally convert back to integer multiplication via Kronecker substitution. These much smaller integer multiplications are handled recursively.

You may also find the relevant section of http://arxiv.org/ useful as well. For example, http://arxiv.org/list/cs.SE/recent

work smart, not hard.

The Fürer algorithm change the loglog term of the complexity to 2^log* which doesn't change anything in practice for workable numbers but it also add a lot to the constant term which is not counted by the big O notation but is very important in practice.

    "We have implemented an unoptimised version of the algorithm
    from section 8 in the Mathemagix system [29] and found our
    implementation to be an order of magnitude slower than the
    Gmp library [23]. There is certainly room for improvement,
    but we doubt that even a highly optimised implementation of
    the new algorithm will be competitive in the near future."

And check out the table of multiplication algorithms in the paper's introduction, and the table in the GMP manual [0]. The FFT (Schönhage–Strassen) GMP uses is O(n log n log log n), and the algorithm in this paper is O(n log n 8^(log* n)). Those functions [1] diverge a bit slowly: they have the same value at n=2^2^2^21 (2^2^2,097,152).

[0] https://gmplib.org/manual/Multiplication-Algorithms.html#Mul...

[1] https://en.wikipedia.org/wiki/Iterated_logarithm