Replacing eigenvalue decomposition with Fourier Transform

Hello,

I am doing matrix inversion of circulant matrix using eigenvalue decomposition method. However, I found that the complexity is high which I know as O(n^3). Then, I found that the complexity can be reduced to O(n log n) using Fourier Transform (FFT) method.

Thus, I would be very happy if someone can explain to me the method how to replace Eigenvalue decomposition method with Fourier Transform (FFT) method in order to do matrix inversion.

Thank you.

Huda

[Christine Tobler]

There's an explanation of how the FFT method relates to circulant matrices here: https://en.wikipedia.org/wiki/Circulant_matrix#Other_properties (see the 4th bullet point).

[Nurulhuda Ismail]

Thank you

[Nurulhuda Ismail]

I tried to do as in the link that you shared. I got the eigenvalues of the circulant matrix, however the order of the diagonal entries is not the same as if we use eig(C) in MATLAB. For example:

Circulant matrix,

C =

16.7397 + 0.0000i -0.0426 + 0.6975i 2.7798 - 0.8504i -0.0426 - 0.6975i

-0.0426 - 0.6975i 16.7397 + 0.0000i -0.0426 + 0.6975i 2.7798 - 0.8504i

2.7798 + 0.8504i -0.0426 - 0.6975i 16.7397 + 0.0000i -0.0426 + 0.6975i

-0.0426 + 0.6975i 2.7798 + 0.8504i -0.0426 - 0.6975i 16.7397 + 0.0000i

My eigenvector of C based on DFT:

MyeigvalC = abs(diag(F*C(:,1))) =

19.4529 0 0 0

0 12.5936 0 0

0 0 19.6231 0

0 0 0 15.3785

Using built in function in MATLAB, eig(C)=

12.4623 0 0 0

0 15.1879 0 0

0 0 19.5652 0

0 0 0 19.7433

If we can see, the answer is in the different order. Is anybody here know how to solve this problem?

Thank you

[Matt J]

Why is it a problem that the eigenvalues are in a different order than as given by eig()? If you want them sorted, for example, why not just use sort()?

[Christine Tobler]

I agree with Matt, best to just call sort() on the eigenvalues computed using DFT. What order do you need the eigenvalues to be in?

I noticed something else: The eigenvalues should be much closer together than in your example above. The problem is that the numbers 2.7798 - 0.8504i and 2.7798 + 0.8504i are both used in the same "lane" of the matrix C, so this matrix is not circular. If I replace the values in that lane to all use "+", the eigenvalues are nearly identical

>> eig(C)
ans =
  12.5649 - 0.8504i
  15.3549 - 0.8504i
  19.4343 + 0.8504i
  19.6047 + 0.8504i
>> fft(C(:, 1))
ans =
  19.4343 + 0.8504i
  12.5649 - 0.8504i
  19.6047 + 0.8504i
  15.3549 - 0.8504i

Also, instead of multiplying with the matrix F, call the function fft instead. The behavior is the same, but it's much faster when the size of the matrix grows.

[Nurulhuda Ismail]

Hi Christine and Matt J,

Thank you for your comments. I did the corrections in my code.

Actually, I am interested to calculate the inversion of a circulant matrix using eigenvalue decomposition based on Fourier matrix. Circulant matrix can be decomposed as C= F^(-1) λF where F is the Fourier matrix, and λ is the eigenvalue of matrix C. So, what I did is:

1. Obtain F using 2 approaches; (i) F ==fft(F*C(:,1)), (ii) F= exp(-2*pi*j1*k1*i/K) for j1=0:K-1, k1=0:K-1
2. Obtain λ=diag(F*C(:,1))
3. Substitutes F and λ into Capp= F^(-1) λF
4. Calculate the inverse of C, (Capp)^(-1)= F^(-1) λ^(-1) F

As the results, the approximate inverse of C, (Capp)^(-1) gives exactly the same result as the actual inversion of C if I choose to use the first approach of F =fft(F*C(:,1)). And, there is not much differences with the second approach using F= exp(-2*pi*j1*k1*i/K) for j1=0:K-1, k1=0:K-1 , which I would say no difference. The different is only the sign of 0 element which happen in a few points. However, when I applied the approximate inversion to my system model using the second approach, for example bit error rate, unfortunately it gives me inaccurate result. You can see the details in the attachment.

What should I do?

[Nurulhuda Ismail]

Please someone help me for this.

[Christine Tobler]

In your PDF there are red circles around some numbers, but those numbers seem accurate up to a reasonable accuracy to me. Can you compute abs(C1 - C2)? I think the error will be very small.