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