
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.



Hi, is your matrix by chance strictly triangular?



Thanks for the response, my matrix is not triangular/ symmetrical.



Could you post the positions of the "1" values for the 649 and 650 matrices? This would make it easier to test.
Thanks,
Marcus



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






There is a bug in the EVD code where we are getting stuck in an infinite loop. Working on it.



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.



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.



That would be great thanks, it might be enough to work for any other ones i will be throwing at it too.



An update. I added EVD to the native wrapper in my fork, but the managed code needs to updated to work with a 1D array instead of a jagged array. I will finish it when I get some free time.




