How to multiply polynomials in Θ(n log n) time

http://agarri.ga/post/multiply-polynomials-n-log-n-time/

by rray • 11 years ago

80 23 11 years ago

glial• 11 years ago

This isn't limited to polynomials. As far as I know, most big integer libraries (like GMP) use FFT multiplication when the integers are big enough.

amelius• 11 years ago

Yes, but integers are a specialization of polynomials. For example, the number 263 can be written as 2x^2 + 6x + 3, where x=10.

(Note that for integers, a post-processing step is necessary to account for coefficients that are larger than the number base, but I guess this is relatively inexpensive).

fdej2• 11 years ago

Indeed, having a fast way to multiply in Z[x] gives you a fast way to multiply in Z and vice versa. In the reverse direction, you can encode the polynomial f = 2x^2 + 6x + 3 as 30602, and 30602^2 = 936482404 tells you that f^2 = 4x^4 + 24x^3 + 48x^2 + 36x + 9 (of course, on a computer you do this in binary as opposed to decimal).

This correspondence is often exploited in both directions. If you ask a sufficiently sophisticated computer algebra system (SageMath, for example) to multiply two large polynomials with integer or rational coefficients, chances are that it will reduce the problem to multiplying two large integers, then reduce that problem to multiplying two large polynomials (but different from the ones you started with), and then proceed using recursive polynomial and integer multiplications to compute that product... (i.e. going through a sequence of transformations Z[x] -> Z -> Z[x] -> ...)

jacobolus• 11 years ago

This blog post is doing precise integer arithmetic without rounding, and dealing with polynomials of tiny degree. But for more general uses...

Anyone trying to work with polynomials as continuous functions using floating point arithmetic (or trying to work efficiently with piecewise continuous functions in general) should check out Chebfun, http://www.chebfun.org

People should avoid working with polynomials expressed in the monomial basis if possible, unless dealing with a function defined on the complex unit circle. The Chebyshev basis is generally faster and more accurate for polynomials defined on some real interval, and makes it practical to compute with polynomials of very high degree. http://www.chebfun.org/examples/linalg/CondNos.html http://www.chebfun.org/examples/roots/RandomPolynomials.html

For a lot more see: http://www.chebfun.org/ATAP/atap-first6chapters.pdf http://www.chebfun.org/docs/guide/chebfun_guide.pdf

Here’s the guts of their multiplication routine, which multiplies two polynomials and returns a new polynomial which approximates the product to near machine precision: https://github.com/chebfun/chebfun/blob/development/%40chebt... (For the mathematical details see the (paywalled) paper http://epubs.siam.org/doi/abs/10.1137/120865458?journalCode=...). Multiplying values at Chebyshev points (similar to the blog post under discussion) would work too, but Chebfun natively stores an array of coefficients. In fact, they did multiply at values before April of last year, https://github.com/chebfun/chebfun/blob/05611dec557718b5b026...

Strilanc2• 11 years ago

This is technically not an O(n lg n) algorithm. Any O(n lg n) polynomial multiplication algorithm can be trivially turned into an O(n lg n) multiplication algorithm (just throw in some carries at the end), but O(n lg n) multiplication is an open problem.

In this case the issue is that you need more and more precision for the FFT as the inputs get larger, or else you'll get the wrong answer. I think it ends up being an O(n log^2 n) algorithm. Still pretty good; worth knowing about for sure. The constant factor on the precision might be bad enough that you can trigger the problem for reasonable inputs if you're using floats (e.g. a million coefficients, instead of an octillion coefficients).

On the other hand, many algorithms and data structures have hidden logarithmic factors (e.g. hash tables, sorting) that we ignore [1]. YMMV.

1: http://blog.computationalcomplexity.org/2009/05/shaving-logs...

pbsd• 11 years ago

The best complexity these days [1] is essentially n log n. More precisely, O( n log n 8^(log* n) ), where log* is the extremely slow-growing iterated logarithm. The original Schönhage-Strassen algorithm is O(n log n log log n).

[1] http://arxiv.org/abs/1407.3360

ot2• 11 years ago

IIRC if your polynomials have integer coefficients with a bounded number of digits you can use FFT on finite rings and get an actual O(n log n) algorithm.

fdej2• 11 years ago

This is not quite correct. A finite ring can only have finitely many primitive roots of unity, so the length of the FFTs you can compute in a ring of fixed size is bounded. To compute arbitrarily long convolutions over a finite ring, you need to extend the ring with more roots of unity (for example, using the Schönhage-Strassen trick), and this increases the complexity.

Strilanc2• 11 years ago

The number of digits needed to hold an output coefficient (and intermediate coefficients during the transform) increases with the number of input coefficients (because so many get added together). This forces the size of the ring to increase to accommodate the larger values, and pushes the runtime past O(n lg n).

abetusk• 11 years ago

Chapter 1 [1] of Chee Yap's "Fundamental Problems of Algorithmic Algebra" [2] gives integer multiplication at O(n lg(n) lg(lg(n))). To prove that bound a "discrete Fourier transform" is used using modular arithmetic over 2^L + 1 (for suitable L).

As you say, you can always convert to an integer multiplication problem so it stands at O(n lg(n) lg(lg(n))) (not O(n lg^2(n))).

[1] http://cs.nyu.edu/~yap/book/alge/ftpSite/l1.ps.gz

[2] https://www.cs.nyu.edu/yap/book/berlin/

tgb2• 11 years ago

Fun article, but it was missing a little bit that took me a minute to get. So if you're already familiar with fourier transforms, this might be the missing piece to this article:

Multiplication of polynomials when expressed in the domain of "array of coefficients" is a convolution operation: really you're summing (AB)_k = sum_{i=0}^k A_i B_{k-i} to find the kth coefficient of the result. The indices for B go backwards to those of A and they add up to the fixed quantity k, so we're just looking at a discrete convolution of these two arrays.

And what does Fourier Transform do? It turns convolutions into point-wise multiplication in the 'time' domain. The convolution started as an n^2 operation but the point-wise multiplication is n, so we're left with the FFT time of n log n. So this is another instance of the "formulate problem as a convolution, FFT, perform point-wise multiplication in the time domain, inverse FFT, done" class of algorithms.

(I got stuck a while trying to figure out how polynomial multiplication is convolution. Literally it's point-wise multiplication of two functions, not a convolution, so everything seemed backwards.)

rhaps0dy2• 11 years ago

This is mentioned in the article in the problem explanation section, but quickly and just passing by. I'll put it right at the beginning where it belongs. Thanks!

I had not understood properly on first read of your comment, but you did not mean what is on the post at all. I'll include this insight somehow :)

Don't you mean signal domain though?

tgb2• 11 years ago

Yeah I messed that domain up. Thanks.

nonsequitur2• 11 years ago

Can anybody guess in which way he used convolution for solving what is generally known as the change-making problem? Wikipedia [1] mentions a "probabilistic convolution tree", but that seems much more involved.

Edit: Solved. I've missed that the problem only deals with change amounts that can be reached with one or two coins. So a single convolution is sufficient in this specific case.

https://en.wikipedia.org/wiki/Change-making_problem

chaoxu• 11 years ago

probabilistic convolution tree is not much more involved. In some sense it's just multiple polynomial multiplication.

But anyway, coin change problem can be solved much faster. http://link.springer.com/article/10.1007%2Fs00453-007-0162-8

nonsequitur2• 11 years ago

Thanks! That paper is available on the author's website: http://gi.cebitec.uni-bielefeld.de/people/zsuzsa/papers/Algo...

charriu• 11 years ago

I realize that this not related to the topic of the post at all, but... why on earth would you (try to) make a pun and then add "(Haha, get it?)" at the end? At least for me, that takes all the fun out of getting a joke :(

vog• 11 years ago

I'm not a native English speaker and I had some trouble understanding that pun.

So how do we do it? We are going to have to talk a bit about the representation of polynomials in computer memory first. (Haha get it?)

Without that hint I would have totally missed it.

Also, this could be a joke on the joke. Probably the author didn't really think this pun is funny, and made fun of themselves by overstating how funny that pun is. Humor is a large field with lots of varieties. (And no, I will not spoil that this is a maths pun.)

acqq2• 11 years ago

Can please somebody explain what the pun is supposed to be for us who definitely don't see anything except the normal sentence? That he used "talk a bit" in the sentence that announces consideration of how computer memory works? I'm really disappointed then.

leni536• 11 years ago

This is the operation called discrete convolution, and it’s the equivalent of multiplying two polynomials in coefficient form.

As a physicist my first thought would be that convolution is simply multiplying in Fourier-space anyway, no need to think about polynomials. I didn't know the application of FFT for polynomials though so thanks for that.

How to multiply polynomials in Θ(n log n) time

Related Stories