Skip to content

Instantly share code, notes, and snippets.

@kumanna
Created February 3, 2026 04:37
Show Gist options
  • Select an option

  • Save kumanna/c7c5118c204f9ff17a9e856b4d4330e0 to your computer and use it in GitHub Desktop.

Select an option

Save kumanna/c7c5118c204f9ff17a9e856b4d4330e0 to your computer and use it in GitHub Desktop.
import numpy as np
N = 8
h = np.array([1+1j, 0.5-0.4j, 0.2 + 0.3j])
x = np.array([1, 3, 2, 4, 7, -3, 4, 5])
assert(len(x) == N)
# Convolution
y = np.convolve(x, h)
print(y)
# Convolution matrix
H = np.zeros((N + len(h) - 1, N), dtype='complex')
H[0, 0] = h[0]
H[1, 0:2] = [h[1], h[0]]
H[2, 0:3] = h[::-1]
H[3, 1:4] = h[::-1]
H[4, 2:5] = h[::-1]
H[5, 3:6] = h[::-1]
H[6, 4:7] = h[::-1]
H[7, 5:8] = h[::-1]
H[8, 6:8] = [h[2], h[1]]
H[9, 7:8] = [h[2]]
assert(np.max(np.abs(H @ x - np.convolve(h, x)) < 1e-8))
# Construct the circulant matrix
H_circ = np.zeros((N, N), dtype='complex')
H_circ[0, 0] = h[0]
H_circ[0, -1] = h[1]
H_circ[0, -2] = h[2]
for i in range(1, N):
H_circ[i,:] = np.roll(H_circ[0,:], i)
# Construct Fourier matrix
FI = 1 / np.sqrt(N) * np.exp(1j * 2 * np.pi * np.outer(np.arange(N), np.arange(N)) / N)
F = FI.conj().T
# The circulant matrix is diagonalized by the DFT. That is:
# H_circ = IDFT matrix * L * DFT Matrix
# L = DFT Matrix * H_circ * IDFT Matrix
L = F @ H_circ @ FI
assert (np.max(np.abs(np.diag(L) == np.fft.fft(h, N))) < 1e-8)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment