no_name
no_name

Reputation: 1375

Diagonalization of block circulant matrix with circulant blocks

Let A be a block circulant matrix with circulant blocks (i.e a BCCB matrix):

A = [1 2 3 4
     2 1 4 3
     3 4 1 2
     4 3 2 1]

that is:

 A = [C1 C2
      C2 C1]

where each block (C1, C2) is a circulant matrix. I've read (see here) that BCCB can be diagonalized by following the equation: A =F*·D·F where F is the 2-D discrete Fourier transform matrix, F* is the conjugate of F, and D is a diagonal matrix whose entries are the eigenvalues of A.

In MATLAB I use this code:

(conj(dftmtx(4))/16*(fft2(A))*dftmtx(4))

but the result is:

[1 4 3 2
 2 3 4 1
 3 2 1 4
 4 1 2 3]

Here the second and the fourth columns of A are switched. Where is the error?

Upvotes: 5

Views: 9345

Answers (1)

Eitan T
Eitan T

Reputation: 32930

Your source is a bit misleading. Diagonalizing a BCCB matrix with DFT is done as follows:

A = (FM⊗FN)*D(FM⊗FN)

where FN is the N-point DFT matrix, M is the number of Cj blocks and N is the size of each individual block (in your example M=2 and N=2). The "⊗" symbol denotes the tensor product.

Also note that F* = conj(F)T (F* is called the complex conjugate transpose matrix). In MATLAB it translates to F' instead of conj(F). Coincidentally, the DFT matrix F is symmetric, which means that F* = conj(F) is also true.

I'm not sure what you are trying to compute, but here's how the diagonalization of A is done in MATLAB:

M = 2; N = 2;
FF = kron(dftmtx(M), dftmtx(N)); %// Tensor product
D = FF' * A * FF / size(A, 1);   %// ' is the conjugate transpose operator

which yields:

D =
    10    0    0    0
     0   -2    0    0
     0    0   -4    0
     0    0    0    0

To diagonalize A using only 2-D FFT operations, you can do this instead:

c = reshape(A(:, 1), N, []);     %// First column of each block
X = fft2(c);
D = diag(X(:));  

or in a one-liner:

D = diag(reshape(fft2(reshape(A(:, 1), N, [])), [], 1));

All of these produce the same diagonal matrix D.

Hope this clarifies things for you!

Upvotes: 6

Related Questions