TAOCP 4.6.4 Exercise 57

Let x(u)=\sum_{j=0}^{n}x_j u^j,\qquad y(u)=\sum_{j=0}^{n}y_j u^j, and let

taocpmathematicsalgorithmsvolume-2math-medium

Section 4.6.4: Evaluation of Polynomials

Exercise 57. [M20] Show that fast Fourier transforms can be used to compute the coefficients of the product $x(u)y(u)$ of two given polynomials of degree $n$, using $O(n \log n)$ operations of (exact) addition and multiplication of complex numbers. [Hint: Consider the product of Fourier transforms of the coefficients.]

Verified: yes
Solve time: 4m31s

Solution

Let

$$ x(u)=\sum_{j=0}^{n}x_j u^j,\qquad y(u)=\sum_{j=0}^{n}y_j u^j, $$

and let

$$ z(u)=x(u)y(u)=\sum_{j=0}^{2n}z_j u^j. $$

The coefficients satisfy the convolution formula

$$ z_k=\sum_{j=0}^{k}x_jy_{k-j}, \qquad 0\le k\le 2n. $$

A direct evaluation of all coefficients requires $\Theta(n^2)$ arithmetic operations. The exercise asks for an $O(n\log n)$ method.

Choose an integer $N$ such that

$$ N>2n, $$

for example $N=2^m$ with $2^m>2n$. Let

$$ \omega=e^{2\pi i/N} $$

be a primitive $N$th root of unity.

Extend the coefficient sequences by zeros:

$$ x_j=y_j=0\qquad (j>n), $$

so that $x_j,y_j$ are defined for $0\le j<N$.

Define the discrete Fourier transforms

$$ X_r=\sum_{j=0}^{N-1}x_j\omega^{rj}, \qquad Y_r=\sum_{j=0}^{N-1}y_j\omega^{rj}, \qquad 0\le r<N. $$

The fast Fourier transform computes all $X_r$ from the coefficients $x_j$, and all $Y_r$ from the coefficients $y_j$, using $O(N\log N)$ additions and multiplications of complex numbers.

For each $r$, form

$$ Z_r=X_rY_r . $$

This requires $N$ complex multiplications.

Now compute the inverse discrete Fourier transform:

$$ z_k=\frac1N\sum_{r=0}^{N-1}Z_r\omega^{-rk}, \qquad 0\le k<N. $$

The inverse transform is again computable by FFT in $O(N\log N)$ operations.

It remains to prove that the resulting numbers $z_k$ are exactly the coefficients of $x(u)y(u)$.

Since

$$ Z_r = \left(\sum_{a=0}^{N-1}x_a\omega^{ra}\right) \left(\sum_{b=0}^{N-1}y_b\omega^{rb}\right) = \sum_{a,b}x_ay_b,\omega^{r(a+b)}, $$

the inverse transform gives

$$ \begin{aligned} z_k &= \frac1N \sum_{r=0}^{N-1} \sum_{a,b} x_ay_b,\omega^{r(a+b-k)} \ &= \sum_{a,b} x_ay_b \left( \frac1N \sum_{r=0}^{N-1} \omega^{r(a+b-k)} \right). \end{aligned} $$

The inner sum equals $1$ when

$$ a+b\equiv k \pmod N, $$

and equals $0$ otherwise, by the orthogonality of roots of unity. Hence

$$ z_k = \sum_{\substack{a,b\a+b\equiv k;(\bmod N)}} x_ay_b . $$

Because $a,b\le n$, we have $a+b\le 2n<N$. Therefore congruence modulo $N$ reduces to ordinary equality:

$$ a+b\equiv k\pmod N \iff a+b=k . $$

Thus

$$ z_k=\sum_{a+b=k}x_ay_b, $$

which is exactly the coefficient of $u^k$ in $x(u)y(u)$.

The computation consists of:

  1. One FFT of length $N$ for $(x_j)$;
  2. One FFT of length $N$ for $(y_j)$;
  3. $N$ pointwise multiplications;
  4. One inverse FFT of length $N$.

Each FFT requires $O(N\log N)$ complex additions and multiplications, and $N=O(n)$. Therefore the total number of exact complex arithmetic operations is

$$ O(N\log N)=O(n\log n). $$

Hence fast Fourier transforms compute all coefficients of the product $x(u)y(u)$ using $O(n\log n)$ exact additions and multiplications of complex numbers.

$$ \boxed{\text{The coefficients of }x(u)y(u)\text{ can be computed in }O(n\log n)\text{ complex operations.}} $$

This completes the proof.

∎