This project has moved. For the latest updates, please go here.

EigenValue Decomposition - method failing

Feb 4, 2013 at 9:01 PM
Hi all,

I have taken to using the Math.NET Numerics library for some eigenvector work as part of a social network analysis project I am doing.

The Evd() method works fine for smaller cases however for some reason, when I try to run the method on a 650 dimensional matrix, the method never (I left it for 4 hours) seems to finish. A 649 dimensional matrix finished in ~6 seconds. In both cases, the matrices are square and defined as follows:
int dimension = 650;
Matrix<double> adjMatrix = new DenseMatrix(dimension, dimension, 0);
var eig = adjMatrix.Evd(); //fails with 650 dim matrix here
I then populate a "1" where a connection exists between nodes in my network.
Here are the matrices should they be of any use to you (loads of zeroes with some ones essentially):

649 (works) - http://txtup.co/ZhbPm

650 (fails) - http://txtup.co/RhHfY

Any help would be greatly appreciated.
Coordinator
Feb 4, 2013 at 9:42 PM
Hi, is your matrix by chance strictly triangular?
Feb 4, 2013 at 10:21 PM
Thanks for the response, my matrix is not triangular/ symmetrical.
Feb 6, 2013 at 9:34 AM
Could you post the positions of the "1" values for the 649 and 650 matrices? This would make it easier to test.

Thanks,
Marcus
Feb 6, 2013 at 3:04 PM
Hi Marcus,

As requested, here are the positions (as coordinates) of the "1" values in both matrices.

649 (works) - http://txtup.co/CZtcC

650 (fails) - http://txtup.co/kXzbR

As a bit of background, each dimension represents a person and if a "1" is in the location, it simply means they have communicated with that individual on that row/column.

I am adding a row into the 650 matrix at index 0 so you will notice a clear pattern as most of the coordinates are very similar (just incremented 1 in 650 to take into account the extra row).

In 650, there are just 2 extra "1" values located at (0, 129) and (129, 0) . Could the fact that this appears to be bidirectional cause problems? From a quick look through, bidirectional communications have been processed before in the working 649 version so it shouldn't matter.

Feel free to ask for anything else.
Feb 7, 2013 at 8:08 AM
Edited Feb 7, 2013 at 8:09 AM
Hi,

I tried to take a look, but the links have expired. Could email the data to me - marcus@cuda.org ?

Thanks,
Marcus
Feb 7, 2013 at 1:04 PM
Hi,

Sorry about that- how aggravating for it to have expired.

I'll email you Marcus but for anyone else's reference, here are the above matrices and positions uploaded to another file sharing site:

Matrices:

http://www.filedropper.com/649-completes

http://www.filedropper.com/650-fails

Positions:

http://www.filedropper.com/649positions

http://www.filedropper.com/650positions
Feb 7, 2013 at 1:42 PM
There is a bug in the EVD code where we are getting stuck in an infinite loop. Working on it.
Feb 7, 2013 at 4:17 PM
The code never converges to a solution. There is an iteration count, but is is never used. The FORTRAN HQR2 code (what JAMA and this code seem to be based on) errors out if the number of iterations is greater than 30*order. We should add something similar.

As to why this particular matrix doesn't converge, I don't know. I tinkered around with the convergence criteria it doesn't even seem to be close to converging.

Also, JAMA uses eps = Math.Pow(2,-52) but we use Math.Pow(2,-53) when determining if two numbers are equal. It doesn't make any difference in this case but I wanted to it out.
Feb 7, 2013 at 7:31 PM
Hi,

Thanks for the help. So you have no ideas as to why the matrix doesn't converge, is there something wrong with the matrix itself? I could really do with being able to display the eigenvalues for my matrix - ideas for alternative methods or anything you can do with Math.Net Numerics?

I had a go with another library from CodeProject ( http://www.codeproject.com/script/Articles/ViewDownloads.aspx?aid=5835 ) and found that I could get the eigenvalues for 650 dimensions fine but then the same issue occurred with 672 and it wouldn't converge. My knowledge of matrices is being pushed as it is so I don't understand any depths of why this wouldn't happen?
Feb 8, 2013 at 7:07 AM
Edited Feb 8, 2013 at 7:07 AM
found that I could get the eigenvalues for 650 dimensions fine but then the same issue occurred with 672
Interesting, since both sets of code seem to be a port of the JAMA code but with a different "eps". I'll poke are around a little more today and try to mimic the way the original FORTRAN code works for convergence (I couldn't find the Algol code that seems to the basis of all these implementations). I'll also see give MKL a spin and seem it converges.
Feb 8, 2013 at 12:55 PM
Edited Feb 8, 2013 at 12:56 PM
Playing around with the convergence criteria again didn't lead anywhere and I couldn't spot the difference between the two sets of code (ours and the the other JAMA port). Reading some comments on HQR2 suggest that there are just some matrices that it doesn't work with but I couldn't figure out what conditions cause it to fail.

MKL seems to handle the 650 matrix fine (and was quite quick). I'll add the LAPACK routine geev to the native wrapper and perhaps that will work for you.
Feb 8, 2013 at 3:10 PM
That would be great thanks, it might be enough to work for any other ones i will be throwing at it too.
Mar 13, 2013 at 9:04 AM
An update. I added EVD to the native wrapper in my fork, but the managed code needs to updated to work with a 1-D array instead of a jagged array. I will finish it when I get some free time.
Apr 24, 2013 at 6:18 AM
Native EVD is now up at https://github.com/mathnet/mathnet-numerics-native along with the updated native libraries.