I have coded two implementations of an algorithm to calculate all the eigenvalues and eigenvectors of a symmetric matrix. One implementation uses the REPA library
While the other implementation uses mutable Vectors and the ST monad
The algorithm is the Jacobi method which description can be found at http://en.wikipedia.org/wiki/Jacobi_eigenvalue_algorithm
I have tested the two implementations using 100 matrices of dimension 100 x 100, running the code sequentially and I have found the following times:
REPA Mutable Vectors
Total time(s) 6.7 28.5
GC (s) 0.2 1.2
The Jacobi Algorithm requires an iterative updating of some of the matrix entries, meaning that most of the matrix is left untouched between iterations. Therefore I would have guessed (erroneously) that the cost of copying a new matrix for each iteration in the REPA implementation will be bigger than the cost of mutating the matrix using the monad ST, because as far as I understand REPA doesn't mutate the array but makes a copy of it.
It is REPA fusing all the operation and avoiding to copy of a new array in each iteration? or it is something else?
Can someone please comment on this result?
U
and delayedD
arrays, and mixing these two properly seems to be important for good performance. I wonder which ones and where you used. -- For having fastST
code it's often important to use tight loops and prevent memory allocation inside (which I guess repa takes care of internally). And I guess copying data around a few times isn't that costly on modern machines, provided they're in one block to take advantage of memory caching. – Idealize