Some readers will look at the title of this post and think “Ah yes, the FFT. I use it all the time. But what is this quadratic reciprocity?”

Others will look at the same title and think “Gauss called the quadratic reciprocity theorem the jewel in the crown of mathematics. But what is this FFT thing? I think I remember an engineer saying something about that.”

Gauss proved a theorem that relates quadratic reciprocity and the FFT:

left(frac{p}{q}right) left(frac{q}{p}right) = frac{text{Tr} {cal F}_{pq}}{ text{Tr} {cal F}_p, text{Tr} {cal F}_q}

I’ll spend the rest of this post unpacking this equation.

Legendre symbols

The expressions on the left are not fractions but rather Legendre symbols. The parentheses are not for grouping but are part of the symbol.

The Legendre symbol

left(frac{a}{r}right)

is defined to be 0 if a is a multiple of r, 1 if a has a square root mod r, and -1 otherwise.

Fourier transforms

The Discrete Fourier Transform (DFT) of a vector of length n multiplies the vector by the n by n Fourier matrix Fp whose j, k entry is equal to exp(2πijk/n). The Fast Fourier Transform (FFT) is a way to compute the DFT more quickly than directly multiplying by the Fourier matrix. Since the DFT is nearly always computed using the FFT algorithm, the DFT is commonly referred to as the FFT.

Matrix trace

The trace of a matrix is the sum of the elements along the main diagonal. So the trace of the Fourier matrix of size n is

text{Tr} {cal F}_n = sum_{j=1}^n exp(2pi ij^2/n)

Illustration

The quadratic reciprocity theorem, also due to Gauss, is usually stated as

left(frac{p}{q}right) left(frac{q}{p}right) = (-1)^{frac{p-1}{2}frac{q-1}{2}}

We won’t give a proof of the theorem at the top of the page but we will illustrate it numerically with the following Python code, using the quadratic reciprocity theorem to evaluate the product of the Legendre symbols.

from numpy import exp, pi

tr = lambda p: sum(exp(2j*pi*k**2/p) for k in range(1, p+1))
p, q = 59, 17
print( tr(p*q)/(tr(p)*tr(q)) )
print( (-1)**((p-1)*(q-1)/4) )

The first print statement produces (0.9999999999998136-1.4048176871018313e-13j) due to some loss of precision due to floating point calculations, but this is essentially 1, which is what the second print statement produces.

If we change q to 19, both statements print −1 (after rounding the first result).

Related posts

The post Connecting the FFT and quadratic reciprocity first appeared on John D. Cook.

Previous Post
Next Post