19
$\begingroup$

I have read a number of explanations of the steps involved in multiplying two polynomials using fast fourier transform and am not quite getting it in practice. I was wondering if I could get some help with a concrete example such as:

$$ p(x) = a_0 + a_2x^2 + a_4x^4 + a_6x^6 $$ $$ q(x) = b_0 + b_4x^4 + b_6x^6 + b_8x^8 $$

Since the $p(x)$'s degree is not a power of 2, do I pad it with zero coefficient terms and evaluate it as a degree 8 polynomial as well? I noticed that these polynomials don't have values for the odd terms; does this still multiply in O(nlog(n)) time complexity? I was unsure how to split it up.

When do I begin evaluating each term at the 8th roots of unity? i.e., $e^{\frac{i\pi k}{8}}$ for $k \in 2,4,6$.

Essentially, I seem to understand each component of component of the fft multiplication when I read it but I am yet to see a step by step concrete example of its process.

Can someone outline the steps for the multiplication of the above polynomials (or a similar simple multiplication) using fft? It would help me a lot.

$\endgroup$
  • 1
    $\begingroup$ It might interest you that "in practice", number theoretic transforms are used instead of "$e^{i\omega}$" transforms because you can't exactly represent a real number on an electronic computer. So integer multiplication uses number theory results for integers to replace the "$e^{i\omega}$". $\endgroup$ – DanielV Apr 22 '14 at 20:32
  • $\begingroup$ @DanielV except for too big polynomials or too big coefficients (like in the case of the Lucas-Lehmer test for Mersenne primes) ? $\endgroup$ – reuns Apr 2 '16 at 7:23
22
$\begingroup$

Let us take your concrete example:- $$p(x) = a_0 + a_2x^2 + a_4x^4 + a_6x^6$$ $$q(x) = b_0 + b_4x^4 + b_6x^6 + b_8x^8$$

We can express the above two polynomials in vector format:- $$\mathbf{p}=[a_0,0,a_2,0,a_4,0,a_6]$$ $$\mathbf{q}=[b_0,0,0,0,b_4,0,b_6,0,b_8]$$

where the $k$th entry (for $k\in\{0,1,2,..\}$) represents the coefficient of $x^k$.

Once we have the vectors, we need to pad $\mathbf{p}$ with $8$ zeros (the length of vector $\mathbf{q}$ minus $1$), and we need to pad $\mathbf{q}$ with $6$ zeros (the length of vector $\mathbf{p}$ minus $1$), so we have:-

$$\mathbf{p'}=[a_0,0,a_2,0,a_4,0,a_6,\color{red}{0,0,0,0,0,0,0,0}]$$ $$\mathbf{q'}=[b_0,0,0,0,b_4,0,b_6,0,b_8,\color{red}{0,0,0,0,0,0}]$$
The reason why we need zero padding is to ensure that when we carry out the FFT on the vectors, multiply the FFTs element wise and carry out an Inverse FFT (IFFT), the result will correspond to a linear convolution.

It is worth noting that the length of $\mathbf{p'}$ and $\mathbf{q'}$ is $15$. As the FFT is designed to optimally process an integer power of $2$ number of samples, we can pad $\mathbf{p'}$ and $\mathbf{q'}$ by an extra $0$ element, so that we have $16$ samples to process.

$$\mathbf{p''}=[a_0,0,a_2,0,a_4,0,a_6,\color{red}{0,0,0,0,0,0,0,0,0}]$$ $$\mathbf{q''}=[b_0,0,0,0,b_4,0,b_6,0,b_8,\color{red}{0,0,0,0,0,0,0}]$$

Next we carry out the FFT on vectors $\mathbf{p''}$ and $\mathbf{q''}$, and we denote the resulting vectors by $\mathbf{P}$ and $\mathbf{Q}$ respectively (each element will be a complex number):-

$$\mathbf{P}=[A_0,A_1,A_2,A_3,A_4,A_5,A_6,A_7,A_8,A_9,A_{10},A_{11},A_{12},A_{13},A_{14},A_{15}]$$ $$\mathbf{Q}=[B_0,B_1,B_2,B_3,B_4,B_5,B_6,B_7,B_8,B_9,B_{10},B_{11},B_{12},B_{13},B_{14},B_{15}]$$ We then multiply each element of $\mathbf{P}$ with the corresponding element of $\mathbf{Q}$, resulting in vector $$\mathbf{R}=[C_0,C_1,C_2,...,C_{15}]$$ where $C_k=A_k\times B_k$.

The final step is to calculate the IFFT of vector $\mathbf{R}$, resulting in $$\mathbf{r}=[r_0,r_1,r_2,...,r_{15}]$$ where element $r_k$ is the coefficient of $x^k$, for $\in{0,1,..,15}$. Note that $r_{15}$ should be $0$, as we zero padded both vectors by an extra $0$ term.

Assuming that the coefficients of the polynomials $p(x)$ and $q(x)$ were real, then the elements of $\mathbf{r}$ should be real (any deviations, in terms of a small imaginary component, could arise from the finite precision nature of how the FFT is implemented on a processor).

Given the nature of the polynomials $p(x)$ and $q(x)$, the elements $r_m$ for odd $m$ will be zero (or a very small number, also due to the finite precision effects).

$\endgroup$
  • $\begingroup$ How are the roots of unity for each polynomial used to reduce the problem size in this case? $\endgroup$ – alh Apr 23 '14 at 3:50
  • 1
    $\begingroup$ In terms of the $O(n\log n)$ time complexity of the FFT, the roots of unity for each polynomial are not used to reduce the problem size. The FFT complexity is reduced due to the fact that we have an integer power of $2$ number of samples. With this condition, the FFT recursively splits the processing into two smaller FFTs, each operating on half the number of samples. This is how we obtain the $O(n\log n)$ time complexity, rather than $O(n^2)$ time complexity. Note that the larger $n$ is, the greater the saving in computation. $\endgroup$ – Alijah Ahmed Apr 23 '14 at 6:09
2
$\begingroup$

You should read "Introduction to algorithms" [CLRS] it explains every process of the algorithm and proof the necesary lemma/theorems to derive it.

for example it answers your question:

"How are the roots of unity for each polynomial used to reduce the problem size?"

The basic idea is that it uses a divide and conquer approach (nlogn part) based on the halving and cancellation lemmas. The roots of unity are used only to get the polynomials in their point-value representation (any different points could be used to get such representation) but these allow to apply conveniently the theorems above and get the asymptotic complexity desired.

"do I pad it with zero coefficient terms and evaluate it as a degree 8 polynomial as well? ".

Taking 2 polynomials in their point-value representation and leting $x_k$ be the points you are going to evaluate both $ 0 \leq k \leq n$ will yield $n+1$ point-value pairs, but you need $2n + 1$ pairs to interpolate a unique polynomial of degree $2n$ (the degree of the product polynomial of both of degree $n$). So you must extend the point-value representations of the polynomials to a $2n + 1$ point-value pairs each.

"When do I begin evaluating each term at the 8th roots of unity?" http://imgur.com/NlZzXv0 This picture answers you "when" (line 4). and is not only evaluated the 8th roots of unity it's also evaluated the 4th roots of unity, 2nd roots of unity. (divide and conquer... O(nlgn)) the size of vector $\vec{a}$ is reduced by a half as you can notice.

"I noticed that these polynomials don't have values for the odd terms; does this still multiply in O(nlog(n)) time complexity? "

Keep in mind that this notation only states an upper bound for the algorithm, this means that if your polynomials have "some" terms the upper bound will be the same as if it has "all" the terms.

If you truly grasp the algorithm and theory behind this (fairly easy but long enough for a simple answer in here IMO) you would be able to "see" the process and answer your doubts such as "when" and "how", "do I ..." and "complexity".

Good luck and keep learning.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.