## Wednesday, September 10, 2014

### a bit of fft: a short note on FFT

Mostly for my future references, and will be expanded later.

Fast Fourier Transform

For: Faster polynomial multiplication in $$O(N\lg{N})$$ time using a divide and conquer strategy, which makes us of a technique called polynomial interpolation, and a fact about the root of unity.

1. Polynomial Interpolation:
A polynomial $$P(x)$$ of degree $$n$$ can be represented by choosing $$n$$ distinct point on the curve $$P(x)$$. In other words, $$P(x)$$ can be constructed back solely by using $$n$$ points $$(x_0,y_0),(x_1,y_1),\ldots,(x_{n-1},y_{n-1})$$ where $$y_k = P(x_k)$$.
Proof that $$P(x)$$ constructed will be unique: Examine the Vandermonde matrix and its especially cool determinant to derive at the proof.

2. Root of unity
Given $$x^n = 1$$, there is a complex number $$w$$ the root of unity such that $$w = e^{\frac{2\pi k}{n}}$$ for $$k = 0,1,2,3,\ldots,n-1$$. Then further we have the Cancellation Lemma: $$w_{in}^{ik} = w_{n}^{k}$$ (direct substitution proof). In particular we have the Halving Lemma: $$(w_{n}^{2k+\frac{n}{2}})^2 = (w_{n}^{2k})^2$$ (Proof also by direct substitution and making use of the previous relationship).

3. FFT
A polynomial can be partitioned based on even-odd parity of its coefficient index:
$$P(x) = a_0+a_1x+a_2x^2+\ldots + a_{n-1}x^{n-1}$$
$$P(x) = A(x^2) + xB(x^2)$$
where
$$A(x) = a_0+a_2x+a_4x^2+\ldots$$ (even indexed coefficients)
$$B(x) = a_1+a_3x+a_5x^2+\ldots$$ (odd indexed coefficients)
Therefore we end up with a recursive relationship: to interpolate P(x), we can recursively interpolate A(x) and B(x) and combine its result in $$O(N)$$ thanks to the Cancellation Lemma and the Halving Lemma which gives rise to a convenient relationship between P and the results from A and B.

4. Pseudo-code:
rec_FFT( polynomial P ): (return interpolated P)
set n = deg(P)
if n == 1:
return P
set w_k = 1
set w = $$e^{\frac{2\pi}{n}}$$ (first root of unity)
set A = polynomial from even coefficient of P
set B = polynomial from odd coefficient of P
a[] = rec_FFT( A )
b[] = rec_FFT( B )
for i = 0 to floor(n/2):
p[i] = a[i] + w_k * b[i]
p[i + n/2] = a[i] - w_k * b[i]
w_k = w * w_k
return p