2

I want to find the eigenvalues and eigenvectors of a 3x3 matrix (mostly if not always symmetric!!). My numbers are stored in fixed-point format (16.16 to be exact).

Note that I don't mind much about the performance, but simply implementing an algorithm that does the job.

The code below, when you build it and run it (with libfixmath library) produces the right eigenvalues but not the correct eigenvectors.

If I understand correctly the algorithm the eigenvectors are the product of all the calculated Qs.

Does anyone know about what might be wrong? (Even correction about the code (writing style etc) anything you can think of, but ofcourse try to center your mind on the eigenvectors! :) :P

The actual loop is like 3 lines...and what it does is this:

eigenvectors = identity matrix

1) QR decomposition A = Q*R

2) Anew = R*Q (multiply the factors in the reverse order, and iterate)

3) eigenvectors = eigenvectors * Q

Thank you!! Oh and C newbie here....

the code:

artless noise
  • 21,212
  • 6
  • 68
  • 105
Franx
  • 87
  • 1
  • 10
  • C code link is broken. suggest using github gists – Jason S Mar 12 '15 at 20:34
  • @JasonS I see that the link is ok. You can not view but download the file. Nevertheless I will change it and start using github gists. Thank you! :) – Franx Mar 13 '15 at 09:17

1 Answers1

0

I have used for many years and with several benefits LAPACK library,

To compare results, you could use R (that uses the LINPACK routine DQRDC2 as a default)m even MATLAB.

And then you can use the qr() command in R with the LAPACK=TRUE option to use a LAPACK routine :

> QR <- qr(Mat,LAPACK=TRUE)
> QR

However, you should be aware that the function qr() in this case uses the LAPACK routine DGEQP3. Contrary to the DGEQRF routine you're using, DGEQP3 calculates the QR decomposition of a matrix with column pivoting.

If you get different results it could be that you're not using the same methods (do not know what method uses the code you posted).

You should remember that a QR decomposition is not a unique solution. To know whether your QR decomposition is correct, you can simply check if the Q and R matrices fulfill the requirements. eg in R :

> Q <- qr.Q(QR)
> round( t(Q) %*% Q , 10 )
     [,1] [,2] [,3] [,4]
[1,]    1    0    0    0
[2,]    0    1    0    0
[3,]    0    0    1    0
[4,]    0    0    0    1

> all.equal(Q %*% qr.R(QR),Mat)
[1] TRUE
edgarmtze
  • 24,683
  • 80
  • 235
  • 386
  • thank you for your answer. Firstly I already checked/compared the results with Octave and they are not the same. The libfixmath's QR algorithm uses gram schmit algorithm for the decomposition and it is correct as i mentioned (I have checked the Q and R). And isnt LAPACK for floating point calculations? – Franx Mar 02 '15 at 10:33