Hi,

I recently ported a simple toy program from F95 to C++/Eigen just to see what the buzz is all about. Generally I have to say I'm positively surprised. Eigen is quite easy to use, and somehow it even manages to spit out somewhat sensible error messages, in contrast to most heavily templated stuff I've seen.

Well, back to my toy program, which calculates the electronic binding energy of carbon structures using the tight binding approximation. The majority of the runtime in this program is spent in a symmetric eigenvalue solver. However, it seems that Eigen is not really competitive with LAPACK here. Timings for a 1600x1600 matrix, in all cases with -msse2 -DNDEBUG:

My laptop, g++ 4.4 (built from trunk 20090214), netlib reference BLAS:

Solving eigenvalues took 22 s 139000000 ns.
Solving eigenvalues with dsyevr took 7 s 970000000 ns.

Another machine, Intel C++ 10.1, reference BLAS:

Solving eigenvalues took 10 s 72000000 ns.
Solving eigenvalues with dsyevr took 3 s 174000000 ns.

Same machine, Intel C++ 10.1, MKL:

Solving eigenvalues took 10 s 540000000 ns.
Solving eigenvalues with dsyevr took 2 s 777000000 ns.

So you see that Eigen is, depending on the BLAS library, phase of the moon, etc. about a factor of 3-4 slower than LAPACK. Now this might not be a really fair comparison, as dsyevr AFAIK is the state of the art eigenvalue solver in LAPACK at the moment, but still.

The actual call to the eigenvalue solver in Eigen is:

VectorXd eigvals =
hamil.marked()
.eigenvalues();

As an aside, if I left the lower triangular part of the matrix unfilled, Eigen produced incorrect eigenvalues even though I have set UpperTrianguar above. However, I haven't looked into this that much myself yet, so unless this is a known bug I probably have missed something obvious.

yes, that's not much a surprise because 95% of the time is spent in the tridiagonalisation which is not super fast and highly limited by the memory (http://eigen.tuxfamily.org/index.php?title=Benchmark). It is still unclear whether we will provide an optional backend to LAPACK or if we'll manage to achieve similar performance ourself...

and, yes, currently the Eigen::UpperTriangular flag is not taken into account by the tridiagonalisation, and so does the eigensolver.

Well, I haven't much to say other than a deep rewrite of all the Householder/Givens operations is on the top of my TODO and will for sure be done before 2.1. I already saw a lot of places in our current code that are very far from optimal from both performance and numerical stability points of view, for instance the Givens rotation code. It's important to say that the impl of eigensolvers and SVD we have in 2.0 was done quickly to ensure we have something but isn't nearly as polished as the rest so if you want to "see what all the buzz is about" then it's definitely not about our eigensolvers.

jbl wrote:Hi,

I recently ported a simple toy program from F95 to C++/Eigen just to see what the buzz is all about. Generally I have to say I'm positively surprised. Eigen is quite easy to use, and somehow it even manages to spit out somewhat sensible error messages, in contrast to most heavily templated stuff I've seen.

Well, back to my toy program, which calculates the electronic binding energy of carbon structures using the tight binding approximation. The majority of the runtime in this program is spent in a symmetric eigenvalue solver. However, it seems that Eigen is not really competitive with LAPACK here. Timings for a 1600x1600 matrix, in all cases with -msse2 -DNDEBUG:

My laptop, g++ 4.4 (built from trunk 20090214), netlib reference BLAS:

Solving eigenvalues took 22 s 139000000 ns.
Solving eigenvalues with dsyevr took 7 s 970000000 ns.

Another machine, Intel C++ 10.1, reference BLAS:

Solving eigenvalues took 10 s 72000000 ns.
Solving eigenvalues with dsyevr took 3 s 174000000 ns.

Same machine, Intel C++ 10.1, MKL:

Solving eigenvalues took 10 s 540000000 ns.
Solving eigenvalues with dsyevr took 2 s 777000000 ns.

So you see that Eigen is, depending on the BLAS library, phase of the moon, etc. about a factor of 3-4 slower than LAPACK. Now this might not be a really fair comparison, as dsyevr AFAIK is the state of the art eigenvalue solver in LAPACK at the moment, but still.

The actual call to the eigenvalue solver in Eigen is:

VectorXd eigvals =
hamil.marked()
.eigenvalues();

As an aside, if I left the lower triangular part of the matrix unfilled, Eigen produced incorrect eigenvalues even though I have set UpperTrianguar above. However, I haven't looked into this that much myself yet, so unless this is a known bug I probably have missed something obvious.

ggael wrote:yes, that's not much a surprise because 95% of the time is spent in the tridiagonalisation which is not super fast and highly limited by the memory (http://eigen.tuxfamily.org/index.php?title=Benchmark). It is still unclear whether we will provide an optional backend to LAPACK or if we'll manage to achieve similar performance ourself...

Since I'm redoing the SVD anyway, I'll at least rewrite the Givens / Householder kernels, that'll at least bring some performance and numerical stability improvement.

I also plan to look closely at the tridiagonalization since that's what i'm doing in math currently anyway, but i make no guarantees that i can make the performance close to lapack...