Jan 21. Polynomial multiplication using FFT

Polynomial multiplication

Given two degree n polynomials A(x) = a_0 + a_1 x + \cdots + a_n x^n and B(x) = b_0 + b_1 x + \cdots + b_n x^n, compute C(x) = A(x) B(x) = c_0 + c_1x + \cdots + c_{2n} x^{2n} where c_k = \sum_{i=0}^k a_i b_{k-i}.

Naive approach:  Computing c_k takes O(k) steps and there are 2n+1 steps, hence it would take \Theta(n^2) time. Using the divide and conquer approach from the integer multiplication algorithm would give us a running time of \Theta(n^{\log_2 3}). (Assume multiplying integers takes O(1) steps). Can we do better?

Fact: A degree n polynomial is uniquely characterized by its values at any distinct n+1 points.

Idea: Suppose we were given A and B evaluated at n+1 distinct points \{x_0,\cdots, x_n\} and wanted C evaluated at the same points, then we could compute it in time \Theta(n) by simply returning \{A(x_0)B(x_0), \cdots, A(x_n)B(x_n)\}. Thus an alternative approach would be:

  • Choose 2n+1 distinct points \{x_0, x_1, \cdots, x_{2n}\} and evaluate A and B on those points
  • Compute C(x_0), C(x_1), \cdots, C(x_{2n})
  • Compute c_0, c_1, \cdots, c_{2n} from C(x_0), \cdots, C(x_{2n}).

Notice that since A and B are of degree n, C has degree at most 2n and so by the above fact C(x_0), \cdots, C(x_{2n}) uniquely fixes C. Clearly the second step takes \Theta(n) time. Hence it remains to address the complexity of the first and third steps. The first step is called evaluation and the third step is called interpolation. Note also that we are free to choose any \{x_0, \cdots, x_{2n}\} as long as they are distinct.

Evaluation

Given a polynomial A(x) of degree n we wish to evaluate it at \{x_0, \cdots, x_n\}. Clearly one can do it in time \Theta(n^2) by evaluating each x_i serially. Can we do it faster?

For general \{x_0, \cdots, x_n\} we cannot; however a clever choice of x_i‘s can help. For example if our points were of the form \{\pm x_0, \pm x_1 ,\cdots, \pm x_{(n-1)/2}\} then  we could write A(x) as A(x) = A_e(x^2) + x A_o (x^2), where the degree of A_e(x), A_o(x) is at most \frac{n}{2}. Suppose we wish to evaluate A on x_i and -x_i. We have A(x_i) = A_e(x_i^2) + x_i A_o(x_i^2) and A(-x_i) = A_e(x_i^2) - x_i A(x_i^2). Thus in order to evaluate A on \{\pm x_0, \cdots, \pm x_{(n-1)/2}\} it suffices to evaluate A_e, A_o on \{x_0^2, x_1^2, \cdots, x_{(n-1)/2}^2\}. Note that we reduced our problem of evaluating a degree n polynomial on n points to the problem of evaluating two degree n/2 polynomials on \frac{n+1}{2} points. If we could ensure that \{x_0, \cdots , x_{(n-1)/2}\} are also of this special form then we could keep doing this trick further.

Choosing the points: Suppose we choose our points to be \{1, \omega, \omega^2, \cdots, \omega^{n}\} where \omega is a primitive n+1^{th} root of unity, i.e. \omega = e^{\frac{2\pi i}{n+1}}. If n is even then

  1. \omega^{\frac{n+1}{2} + j} = - \omega^j
  2. \omega^2 is a primitive \frac{n+1}{2}^{th} root of unity.

Notice that if n+1 is a power of 2 then the two properties ensure that the above trick works. Property 1 ensures that \{1, \omega, \omega^2, \cdots, \omega^n\} can be written as \{ \pm 1, \pm \omega, \cdots, \pm \omega^{\frac{n-1}{2}}\}. Property 2 ensures that we can keep doing this since \omega^2 is a primitive \frac{n+1}{2}^{th} root of unity.

This gives us the following algorithm FFT(A, \omega) for evaluating A on \{1, \omega, \cdots, \omega^n\}.

  • If \omega = 1 return A(1)
  • Let A(x) = A_e(x^2) + xA_o(x^2)
  • Recursively compute A_e(1), A_e(\omega^2), \cdots, A_e(\omega^{(n-1)/2}) by calling FFT(A_e, \omega^2).
  • Recursively compute A_o(1), A_o(\omega^2), \cdots, A_o(\omega^{(n-1)/2}) by calling FFT(A_o, \omega^2).
  • For j = 0 to n, compute A(\omega^j) as A_e(\omega^{2j}) + \omega^j A_o(\omega^{2j}).
  • Return A(1), A(\omega), \cdots, A(\omega^n).

If T(n) is the time to evaluate a degree n polynomial on \{1, \omega, \cdots, \omega^n\} then we get T(n) = 2 T(\frac{n}{2}) + O(n) implying T(n) = \Theta(n \log n).

Interpolation

Given A(1), A(\omega), \cdots, A(\omega^n) we wish to compute a_0, \cdots, a_n. Notice that for any x_0, \cdots, x_n we have

\begin{bmatrix} A(x_0) \\ A(x_1) \\ \vdots \\ A(x_n)\end{bmatrix} = \begin{bmatrix} 1 & x_0 & x_0^2 & \cdots & x_0^n \\ 1 & x_1 & x_1^2 & \cdots & x_1^n \\ \vdots \\ 1 & x_n & x_n^2 & \cdots & x_n^n\end{bmatrix} \begin{bmatrix} a_0 \\ a_1 \\ \vdots \\ a_n \end{bmatrix}

One can show that if the x_i‘s are distinct then this matrix is invertible. Let V(\omega) be this matrix where x_i = \omega^i. This matrix is called the Van der Monde matrix. It has the following nice property: V(\omega)^{-1} = \frac{1}{n} V(\omega^{-1}).

In other words to compute a_0, \cdots, a_n one needs to evaluate the polynomial \hat{A} at \{1, \omega^{-1}, \omega^{-2}, \cdots, \omega^{-n}\}, where \hat{A}(x) = A(1) + A(\omega) x + A(\omega^2) x^2 + \cdots + A(\omega^n) x^n. Observe that \omega^{-1} is also a primitive n+1^{th} root of unity. But we already know how to solve this problem! Thus interpolation is equivalent to evaluation and has time complexity \Theta(n\log{n}).

Combining the time complexity of the three steps gives us a \Theta(n\log{n}) algorithm to multiply two polynomials of degree n.

Advertisements

4 Comments

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s