Fast Fourier Transform

The fast Fourier transform is an algorithm for computing discrete Fourier transforms efficiently. In numerical quantum mechanics, FFTs are used to move between position and momentum grids, compute spectral derivatives, apply kinetic-energy phases, analyze wave packets, and implement split-operator time evolution.

The FFT is not a new Fourier convention. It is a fast way to compute a finite transform. The physics still depends on the chosen grid, normalization, sign convention, and mapping between array indices and physical momenta.

Discrete Fourier Transform

For $N$ complex samples $\psi_j$ , one common computational convention is

\widehat\psi_m = \sum_{j=0}^{N-1} \psi_j e^{-2\pi i jm/N}, \qquad m=0,\dots,N-1,

with inverse

\psi_j = \frac{1}{N} \sum_{m=0}^{N-1} \widehat\psi_m e^{2\pi i jm/N}.

Many libraries use this convention or a close variant. Some put the factor $1/N$ on the forward transform, and some use a unitary factor $1/\sqrt N$ in both directions. Always check the convention before interpreting amplitudes physically.

The continuous position-momentum convention used elsewhere is reviewed in Fourier Transform Conventions.

Quantum Fourier Transform applies a unitary discrete transform to amplitudes of a quantum register; unlike this classical FFT, it does not accept and return an explicitly readable array of all $N$ coefficients. The shared matrix does not make the two input-output contracts interchangeable.

Why the FFT Is Fast

A direct discrete Fourier transform computes $N$ output numbers, each as a sum over $N$ inputs, so it costs order $N^2$ operations.

The FFT exploits the periodic structure of the roots of unity. For many values of $N$ , especially powers of two or products of small primes, the transform can be decomposed into smaller transforms. This reduces the cost to order

N\log N.

That difference is decisive for wave packets, spectral methods, and multidimensional grids. The FFT is why Fourier-grid methods can apply kinetic energy in momentum space without forming a dense derivative matrix.

Position Grid and Momentum Grid

Take a periodic interval of length $L$ with

x_j = x_0+j\Delta x, \qquad \Delta x = \frac{L}{N}, \qquad j=0,\dots,N-1.

Do not include both endpoints $x_0$ and $x_0+L$ as distinct grid points. A periodic grid identifies them.

The corresponding wave numbers are spaced by

\Delta k = \frac{2\pi}{L}.

For the DFT ordering above, the physical wave number associated with index $m$ is commonly read as

k_m = \frac{2\pi}{L} \begin{cases} m, & 0\le m\lt N/2,\\ N/2, & m=N/2\ \text{for even }N,\\ m-N, & N/2\lt m\lt N. \end{cases}

The $m=N/2$ Nyquist mode for even $N$ needs convention care. Many libraries store it as a single real-frequency boundary mode for real-input transforms. For complex quantum wavefunctions, keep the full complex convention explicit.

Momentum values are

p_m = \hbar k_m.

This is the finite periodic-box version of the position-momentum Fourier relation.

Physical Normalization

Suppose $\psi_j\approx\psi(x_j)$ and the continuum state is normalized by

\int_0^L \lvert\psi(x)\rvert^2\,dx = 1.

The grid norm is

\sum_{j=0}^{N-1} \lvert\psi_j\rvert^2 \Delta x \approx 1.

If the periodic basis functions are $L^{-1/2}e^{ik_mx}$ , then the Fourier-series coefficient is approximated by

c_m \approx \frac{\Delta x}{\sqrt L} \sum_{j=0}^{N-1} \psi_j e^{-ik_mx_j}.

With the computational DFT, this is a scaled version of $\widehat\psi_m$ . The scale matters: raw FFT output is not automatically a normalized momentum-space wavefunction.

Discrete Parseval consistency should give

\sum_m \lvert c_m\rvert^2 \approx \sum_j \lvert\psi_j\rvert^2\Delta x.

If this check fails, the likely causes are a missing factor of $\Delta x$ , a misplaced factor of $N$ , or a sign convention mismatch.

Spectral Differentiation

On a periodic grid, differentiation is simple in Fourier space. If

\psi(x) \approx \sum_m c_m e^{ik_mx},

then

\frac{d\psi}{dx} \approx \sum_m ik_m c_m e^{ik_mx},

and

\frac{d^2\psi}{dx^2} \approx \sum_m (-k_m^2)c_m e^{ik_mx}.

Thus the kinetic-energy operator for a free particle is diagonal in the Fourier grid:

T_m = \frac{\hbar^2k_m^2}{2m_{\mathrm{phys}}}.

This is the computational version of the momentum-space simplification described in Momentum Representation and Spectral Methods.

Split-Operator Connection

For a Hamiltonian

H = T(p)+V(x),

the potential is diagonal in position space and the kinetic energy is diagonal in momentum space. A standard second-order split step is

\psi(t+\Delta t) \approx e^{-iV\Delta t/(2\hbar)} \mathcal F^{-1} e^{-iT\Delta t/\hbar} \mathcal F e^{-iV\Delta t/(2\hbar)} \psi(t).

Here $\mathcal F$ denotes the discrete Fourier transform with the chosen normalization. The two FFTs move the state to momentum space and back.

This formula is powerful because it avoids forming dense matrices. Its accuracy still depends on time-step size, operator noncommutativity, grid resolution, and boundary periodicity. The time-step issues are introduced in Time-Stepping Methods.

Periodicity and Wraparound

The DFT represents periodic data. A wave packet leaving the right edge of the interval reappears at the left edge unless the physical model, boundary treatment, or domain size prevents it.

This wraparound is not a small numerical detail. It is the boundary condition of the Fourier grid. For free-particle wave packets, choose a domain large enough that the packet does not hit the boundary during the simulated time, or use a method designed for open boundaries.

Endpoint mismatch also matters. If a nonperiodic function is placed on a Fourier grid, the periodic extension has a jump or kink. Fourier coefficients then decay slowly, spectral derivatives become contaminated, and Gibbs oscillations may appear.

Aliasing

Products in position space become convolutions in Fourier space. Multiplying $V(x)\psi(x)$ can create modes beyond the represented cutoff. Those modes can fold back into the resolved band as aliasing.

Aliasing is especially visible in nonlinear equations, but it also appears in linear quantum mechanics when the potential is rough or under-resolved. Diagnostics include increasing $N$ , smoothing only with physical justification, comparing with finite differences, and monitoring high-frequency Fourier coefficients.

Practical Checklist

Before trusting an FFT-based quantum calculation, record:

the physical interval length $L$ and grid spacing $\Delta x$ ;
whether endpoints are repeated;
the DFT sign and normalization convention;
the ordering of positive, negative, and Nyquist modes;
the mapping $p_m=\hbar k_m$ ;
the discrete norm check used in position and momentum space;
whether the wavefunction remains negligible near periodic boundaries;
the treatment of aliasing and high-frequency modes;
the time-step and split-operator error checks.

These details determine the physical meaning of the array.

Common Mistakes

Treating raw FFT output as a continuum momentum wavefunction without scaling.
Repeating both endpoints on a periodic grid.
Forgetting that FFT frequencies wrap from positive to negative values.
Mixing library DFT signs with the physics Fourier convention.
Ignoring the Nyquist mode for even $N$ .
Applying Fourier differentiation to a nonperiodic function without checking endpoint mismatch.
Letting a wave packet wrap around the periodic box and interpreting it as physical recurrence.
Assuming split-operator error is only a Fourier-grid issue rather than also a time-step issue.
Confusing wave number $k$ with physical momentum $p=\hbar k$ .

Cross-Links

References

J. W. Cooley and J. W. Tukey, “An algorithm for the machine calculation of complex Fourier series”, Mathematics of Computation 19, 297-301, 1965.
M. Frigo and S. G. Johnson, “The design and implementation of FFTW3”, Proceedings of the IEEE 93, 216-231, 2005.
L. N. Trefethen, Spectral Methods in MATLAB, SIAM, 2000.
J. P. Boyd, Chebyshev and Fourier Spectral Methods, 2nd ed., Dover, 2001.
J. M. Thijssen, Computational Physics, 2nd ed., Cambridge University Press, 2007.
D. J. Tannor, Introduction to Quantum Mechanics: A Time-Dependent Perspective, University Science Books, 2007.

Exercises

For a periodic interval of length $L$ with $N$ grid points, derive the wave-number spacing $\Delta k$ .

Solution

Periodic modes satisfy

e^{ik(x+L)} = e^{ikx}.

Thus $e^{ikL}=1$ , so $kL=2\pi n$ for integer $n$ . Adjacent allowed wave numbers differ by

\Delta k = \frac{2\pi}{L}.

Explain why the endpoint $x_0+L$ should not be stored as an additional point on a periodic FFT grid.

Solution

On a periodic interval, $x_0$ and $x_0+L$ represent the same physical point. Storing both duplicates one point, breaks the uniform $N$ -point periodic sampling assumed by the DFT, and changes the implied spacing and normalization.

A computational DFT uses no normalization on the forward transform and $1/N$ on the inverse. Why must a physical Fourier coefficient include factors involving $\Delta x$ or $L$ ?

Solution

The continuum coefficient is an integral against a normalized basis function. A grid approximation to the integral contains a quadrature factor $\Delta x$ , and the periodic basis contributes $1/\sqrt L$ . The raw DFT sum is only the unweighted finite sum, so it must be scaled before it has the physical units and normalization of a Fourier coefficient.

Why is the kinetic-energy step diagonal in the Fourier representation for a free particle?

Solution

Fourier modes satisfy

\frac{d^2}{dx^2}e^{ikx} = -k^2e^{ikx}.

Since the free-particle kinetic energy is $-\hbar^2(2m_{\mathrm{phys}})^{-1}d^2/dx^2$ , each Fourier mode is multiplied by

\frac{\hbar^2k^2}{2m_{\mathrm{phys}}}.

Therefore the kinetic operator is diagonal in the Fourier basis.