darkfi/doc/src/crypto/fft.md

Discrete Fast Fourier Transform

Available code files:

• fft2.sage: implementation using vandermonde matrices illustrating the theory section below.
• fft3.sage: simple example with $n = 4$ showing 3 steps of the algorithm.
• fft4.sage: illustrates the full working algorithm.

Theory

$$ f(x)g(x) ∈ 𝔽{<2n}[x] $$ $$ fg = \sum{i + j < 2n - 2} a_i b_j x^{i + j} $$ Complexity: $O(n^2)$

Suppose $ω ∈ 𝔽$ is an nth root of unity.

Recall: if $𝔽 = 𝔽{p^k}$ then $∃N : 𝔽{p^N}$ contains all nth roots of unity.

$$ \DFT_ω : 𝔽^n → 𝔽^n $$ $$ \DFT_ω(f) = (f(ω^0), f(ω^1), …, f(ω^{n - 1})) $$

$$ V_ω = \begin{pmatrix} 1 & 1 & 1 & \cdots & 1 \ 1 & ω^1 & ω^2 & \cdots & ω^{n - 1} \ 1 & ω^2 & ω^4 & \cdots & ω^{2(n - 1)} \ \vdots \ 1 & ω^{n - 1} & ω^{2(n - 1)} & \cdots & ω^{(n - 1)^2} \ \end{pmatrix} $$

$$ \DFT_ω(f) = V_ω · f^T $$ since vandermonde multiplication is simply evaluation of a polynomial.

Lemma: $Vω^{-1} = \frac{1}{n} V{ω^{-1}}$

Use $1 + ω + ⋯ + ω^{n - 1}$ and compute $Vω V{ω^{-1}}$

Corollary: $\DFT_ω$ is invertible.

Definitions

1. Convolution $f * g = fg \mod (x^n - 1)$
2. Pointwise product $$(a0, …, a{n - 1})·(b0, …, b{n - 1}) = (a_0 b0, …, a{n - 1} b{n - 1}) ∈ 𝔽^n → 𝔽{<n}[x]$$

Theorem: $\DFT_ω(f*g) = \DFT_ω(f)·\DFT_ω(g)$

$$ fg = q'(x^n - 1) + f*g $$ $$ ⇒ f*g = fg + q(x^n - 1) $$ $$ \deg fg ≤ 2n - 2 $$

$$ \begin{align} (f*g)(ω^i) &= f(ω^i)g(ω^i) + q(ω^i)(ω^{in} - 1) \

       &= f(ω^i)g(ω^i)

\end{align} $$

Result

$$ f, g ∈ 𝔽_{<n/2}[x] $$ $$ fg = f*g $$ $$ \DFT_ω(f*g) = \DFT_ω(f)·\DFTω(g) $$ $$ fg = \frac{1}{n} \DFT{ω^{-1}} (\DFT_ω(f) · \DFT_ω(g)) $$

Finite Field Extension Containing Nth Roots of Unity

$$ μN = ⟨ω⟩, |𝔽{p^N}^×| = p^N - 1 $$ $$ \textrm{ord}(ω) = n | p^N - 1 $$ but $𝔽_{p^N}^×$ is cyclic.

For all $d | p^N - 1$, there exists $x ∈ 𝔽_{p^N}^×$ with ord$(x) = d$.

Finding $n | p^N - 1$ is sufficient for $ω ∈ 𝔽_{p^N}$ $$ n | p^N - 1 ⇔ \textrm{ord}(p) = (ℤ / nℤ)^× $$

FFT Algorithm Recursive Compute

We recurse to a depth of $\log n$. Since each recursion uses $ω^i$, then in the final step $ω^i = 1$, and we simply return $f^T$.

We only need to prove a single step of the algorithm produces the desired result, and then the correctness is inductively proven.

$$ \begin{align} f(X) &= a_0 + a_1 X + a2 X^2 + ⋯ + a{n - 1} X^{n - 1} \

    &= g(X) + X^{n/2} h(X)

\end{align} $$

Algorithm

Implementation of this algorithm is available in fft4.sage. Particularly the function called calc_dft().

function DFT($n = 2^d, f(X)$)
if $n = 1$ then
return $f(X)$
end
Write $f(X)$ as the sum of two polynomials with equal degree
$f(X) = g(X) + X^{n/2} h(X)$
Let $\mathbf{g}, \mathbf{h}$ be the vector representations of $g(X), h(X)$

$\mathbf{r} = \mathbf{g} + \mathbf{h}$
$\mathbf{s} = (\mathbf{g} - \mathbf{h})·(ω^0, …, ω^{n/2 - 1})$
Let $r(X), s(X)$ be the polynomials represented by the vectors $\mathbf{r}, \mathbf{s}$

Compute $(r(ω^0), …, r(ω^{n/2})) = \textrm{DFT}{ω^2}(n/2, r(X))$
Compute $(s(ω^0), …, s(ω^{n/2})) = \textrm{DFT}{ω^2}(n/2, s(X))$

return $(r(ω^0), s(ω^0), r(ω^2), s(ω^2), …, r(ω^{n/2}), s(ω^{n/2}))$
end

Sage code:

def calc_dft(ω_powers, f):
    m = len(f)
    if m == 1:
        return f
    g, h = vector(f[:m/2]), vector(f[m/2:])

    r = g + h
    s = dot(g - h, ω_powers)

    ω_powers = vector(ω_i for ω_i in ω_powers[::2])
    rT = calc_dft(ω_powers, r)
    sT = calc_dft(ω_powers, s)

    return list(alternate(rT, sT))

Even Values

$$ \begin{align} r(X) &= g(X) + h(X) \ \ f(ω^{2i}) &= g(ω^{2i}) + (ω^{2i})^{n/2} h(ω^{2i}) \

        &= g(ω^{2i}) +                h(ω^{2i}) \\
        &= (g + h)(ω^{2i}) \\

\end{align} $$

So then we can now compute $DFTω(f){k=2i} = DFT_{ω^2}(r)$ for the even powers of $f(ω^{2i})$.

Odd Values

For odd values $k = 2i + 1$

$$ \begin{align} s(X) &= (g(X) - h(X))·(ω^0, …, ω^{n/2 - 1}) \ \ f(X) &= a_0 + a_1 X + a2 X^2 + ⋯ + a{n - 1} X^{n - 1} \

          &= g(X) + X^{n/2} h(X) \\

f(ω^{2i + 1}) &= g(ω^{2i + 1}) + (ω^{2i + 1})^{n/2} h(ω^{2i + 1}) \ \end{align} $$ But observe that for any $n$th root of unity $ω^n = 1$ and $ω^{n/2} = -1$ $$ (ω^{2i + 1})^{n/2} = ω^{in} ω^{n/2} = ω^{n/2} = -1 $$ $$ \begin{align} ⇒ f(ω^{2i + 1}) &= g(ω^{2i + 1}) - h(ω^{2i + 1}) \

            &= (g - h)(ω^{2i + 1})

\end{align} $$

Let $\mathbf{s} = (\mathbf{g} - \mathbf{h})·(ω^0, …, ω^{n/2 - 1})$ be the representation for $s(X)$. Then we can see that $s(ω^{2i}) = (g - h)(ω^{2i + 1})$ as desired.

So then we can now compute $DFTω(f){k=2i + 1} = DFT_{ω^2}(s)$ for the odd powers of $f(ω^{2i + 1})$.

Example

Let $n = 8$ $$ \begin{align} f(X) &= (a_0 + a_1 X + a_2 X^2 + a_3 X^3) + (a_4 X^4 + a_5 X^5 + a_6 X^6 + a_7 X^7) \

    &= (a_0 + a_1 X + a_2 X^2 + a_3 X^3) + X^4 (a_4     + a_5 X   + a_6 X^2 + a_7 X^3) \\
    &= g(X) + X^{n/2} h(X) \\

g(X) &= a_0 + a_1 X + a_2 X^2 + a_3 X^3 \ h(X) &= a_4 + a_5 X + a_6 X^2 + a_7 X^3 \ \end{align} $$ Now vectorize $g(X), h(X)$ $$ \begin{align} \mathbf{g} &= (a_0, a_1, a_2, a_3) \ \mathbf{h} &= (a_4, a_5, a_6, a_7) \ \end{align} $$ Compute reduced polynomials in vector form $$ \begin{align} \mathbf{r} &= \mathbf{g} + \mathbf{h} \

       &= (a_0 + a_4, a_1 + a_5, a_2 + a_6, a_3 + a_7) \\

\mathbf{s} &= (\mathbf{g} - \mathbf{h})·(1, ω, ω^2, ω^3) \

       &= (a_0 - a_4, a_1 - a_5, a_2 - a_6, a_3 - a_7)·(1, ω, ω^2, ω^3) \\
       &= (a_0 - a_4, ω (a_1 - a_5), ω^2 (a_2 - a_6), ω^3 (a_3 - a_7)) \\

\end{align} $$ Convert them to polynomials from the vectors. We also expand them out below for completeness. $$ \begin{align} r(X) &= r_0 + r_1 X + r_2 X^2 + r_3 X^3 \

       &= (a_0 + a_4) + (a_1 + a_5) X + (a_2 + a_6) X^2 + (a_3 + a_7) X^3 \\

s(X) &= s_0 + s_1 X + s_2 X^2 + s_3 X^3 \

       &= (a_0 - a_4) + ω (a_1 - a_5) X + ω^2 (a_2 - a_6) X^2 + ω^3 (a_3 - a_7) X^3 \\

\end{align} $$ Compute $$ \textrm{DFT}{ω^2}(4, r(X)), \textrm{DFT}{ω^2}(4, s(X)) $$ The values returned will be $$

(r(1), s(1), r(ω^2), s(ω^2), r(ω^4), s(ω^4), r(ω^6), s(ω^6))

(f(1), f(ω), f(ω^2), f(ω^3), f(ω^4), f(ω^5), f(ω^6), f(ω^7)) $$ Which is the output we return.

Comparing Evaluations for $f(X)$ and $r(X), s(X)$

We can see the evaluations are correct by substituting in $ω^i$.

We expect that $s(X)$ on the domain $(1, ω^2, ω^4, ω^6)$ produces the values $(f(1), f(ω^2), f(ω^4), f(ω^6))$, while $r(X)$ on the same domain produces $(f(ω), f(ω^3), f(ω^5), f(ω^7))$.

Even Values

Let $k = 2i$, be an even number. Then note that $k$ is a multiple of 2, so $4k$ is a multiple of $n ⇒ ω^{4k} = 1$, $$ \begin{align} r(X) &= (a_0 + a_4) + (a_1 + a_5) X + (a_2 + a_6) X^2 + (a_3 + a_7) X^3 \ r(ω^{2i}) &= (a_0 + a_4) + (a_1 + a_5) ω^{2i} + (a_2 + a_6) ω^{4i} + (a_3 + a_7) ω^{6i} \ f(ω^k) &= (a_0 + a_1 ω^k + a_2 ω^{2k} + a_3 ω^{3k}) + ω^{4k} (a_4 + a_5 ω^k + a_6 ω^{2k} + a_7 ω^{3k}) \

       &= (a_0 + a_1 ω^k + a_2 ω^{2k} + a_3 ω^{3k}) +        (a_4     + a_5 ω^k   + a_6 ω^{2k} + a_7 ω^{3k}) \\
       &= (a_0 + a_4) + (a_1 + a_5) ω^k + (a_2 + a_6) ω^{2k} + (a_3 + a_7) ω^{3k} \\
       &= f(ω^{2i}) \\
       &= (a_0 + a_4) + (a_1 + a_5) ω^{2i} + (a_2 + a_6) ω^{4i} + (a_3 + a_7) ω^{6i} \\
       &= r(ω^{2i})

\end{align} $$

Odd Values

For $k = 2i + 1$ odd, we have a similar relation where $4k = 8i + 4$, so $ω^{4k} = ω^4$. But observe that $ω^4 = -1$. $$ \begin{align} s(X) &= (a_0 - a_4) + ω (a_1 - a_5) X + ω^2 (a_2 - a_6) X^2 + ω^3 (a_3 - a_7) X^3 \ s(ω^{2i}) &= (a_0 - a_4) + (a_1 - a_5) ω^{2i + 1} + (a_2 - a_6) ω^{4i + 2} + (a_3 - a_7) ω^{6i + 3} \ f(ω^k) &= (a_0 + a_1 ω^k + a_2 ω^{2k} + a_3 ω^{3k}) + ω^{4k} (a_4 + a_5 ω^k + a_6 ω^{2k} + a_7 ω^{3k}) \

       &= (a_0 + a_1 ω^k + a_2 ω^{2k} + a_3 ω^{3k}) -        (a_4     + a_5 ω^k   + a_6 ω^{2k} + a_7 ω^{3k}) \\
       &= f(ω^{2i + 1}) \\
       &= (a_0 + a_1 ω^{2i + 1} + a_2 ω^{4i + 2} + a_3 ω^{6i + 3}) -  (a_4     + a_5 ω^{2i + 1}   + a_6 ω^{4i + 2} + a_7 ω^{6i + 3}) \\
       &= 
       (a_0  - a_4)
       + (a_1 - a_5) ω^{2i + 1}
       + (a_2 - a_6) ω^{4i + 2}
       + (a_3 - a_7) ω^{6i + 3} \\
       &= s(ω^{2i})

\end{align} $$