4 ms·

Eigen is still the best-in-class C++ library for linear algebra as far as I know. It's been around for a while, is still actively developed, and I don't think t

Eigen is still the best-in-class C++ library for linear algebra as far as I know. It's been around for a while, is still actively developed, and I don't think there's a risk of it becoming obsolete anytime soon. That being said, the learning curve for someone not familiar with C++ to write code using Eigen is much higher than writing the same code in Fortran. If you're not already a C++ programmer, you just want to do some math, and Numpy/Matlab are too slow or constrained, then Fortran is a solid choice.

enriquto 7y ago
> That being said, the learning curve for someone not familiar with C++ to write code using Eigen is much higher than writing the same code in Fortran. Moreover, the algorithms implemented inside eigen are hidden behind dozens of onion-like layers. Once you peel all these layers you find code such as this: https://eigen.tuxfamily.org/dox/JacobiSVD_8h_source.html https://eigen.tuxfamily.org/dox/JacobiSVD_8h_source.html I can read and write C++ and I teach numerical linear algebra for a living, but mother of god, this SVD implementation is horrendously unreadable to me. Most of the lines of code are stupid bureaucracy necessary just to add and multiply a few vectors. The same algorithm in linear algebra textbooks requires about fifteen lines of pseudo-code, and similarly for a fortran implementation.
- alexhutcheson 7y ago
  Yeah, libraries that make extensive use of templates are basically dark magic only comprehensible to C++ experts (see also: almost everything in Boost). The positive side of this is that it enables nice interfaces for the library, and makes a lot of the abstractions basically "free" (since the cost is paid at compile time rather than via pointer indirection at runtime). But it definitely makes those libraries inappropriate for students who want to look at an understandable implementation of the algorithm.
  - enriquto 7y ago
    > ...makes a lot of the abstractions basically "free" Maybe it is just me, but I do not really see the need for "abstractions" in a linear algebra library. Numbers are already an abstraction, you do not want to hide them! Ok, maybe you want to choose between "float" and "double", but this does not merit all that overcomplication.
    - vbarrielle 7y ago
      The abstractions inside Eigen are there to enable powerful compile time optimisations, such as compiling the addition of three vectors `Eigen::VectorXf x = a + b + c;` in a single loop instead of naively allocating a temporary to evaluate `Eigen::VectorXf tmp = b + c;` then evaluating `Eigen::VectorXf x = a + tmp;`. This logic is omnipresent in Eigen and is crucial to its performance. It is also used to have specialized versions of some decompositions when the matrix's size is known at compile time.
  - certik 7y ago
    The main advantage of this approach is that you can use any C++ compiler and the library (in this case Eigen) will work. The disadvantage of this approach is that the code is incomprehensible to average scientists, you must spend a lot of effort to become a C++ expert. The advantage of Fortran is that the code is simple and comprehensible to average scientist, and yet fast, because the effort goes into the compiler itself. I would also argue it is easier to implement an optimization pass in the compiler than implement optimizations on the template level in C++. The disadvantage is that you need a good Fortran compiler. If your compiler can't run your code on, say, GPU, then suddenly there is not a good path forward. While in C++, there is always a way forward via template metaprogramming: perhaps ugly, but at least it gets the job done.

enriquto 7y ago

> That being said, the learning curve for someone not familiar with C++ to write code using Eigen is much higher than writing the same code in Fortran. Moreover, the algorithms implemented inside eigen are hidden behind dozens of onion-like layers. Once you peel all these layers you find code such as this: https://eigen.tuxfamily.org/dox/JacobiSVD_8h_source.html https://eigen.tuxfamily.org/dox/JacobiSVD_8h_source.html I can read and write C++ and I teach numerical linear algebra for a living, but mother of god, this SVD implementation is horrendously unreadable to me. Most of the lines of code are stupid bureaucracy necessary just to add and multiply a few vectors. The same algorithm in linear algebra textbooks requires about fifteen lines of pseudo-code, and similarly for a fortran implementation.

alexhutcheson 7y ago
Yeah, libraries that make extensive use of templates are basically dark magic only comprehensible to C++ experts (see also: almost everything in Boost). The positive side of this is that it enables nice interfaces for the library, and makes a lot of the abstractions basically "free" (since the cost is paid at compile time rather than via pointer indirection at runtime). But it definitely makes those libraries inappropriate for students who want to look at an understandable implementation of the algorithm.
- enriquto 7y ago
  > ...makes a lot of the abstractions basically "free" Maybe it is just me, but I do not really see the need for "abstractions" in a linear algebra library. Numbers are already an abstraction, you do not want to hide them! Ok, maybe you want to choose between "float" and "double", but this does not merit all that overcomplication.
  - vbarrielle 7y ago
    The abstractions inside Eigen are there to enable powerful compile time optimisations, such as compiling the addition of three vectors `Eigen::VectorXf x = a + b + c;` in a single loop instead of naively allocating a temporary to evaluate `Eigen::VectorXf tmp = b + c;` then evaluating `Eigen::VectorXf x = a + tmp;`. This logic is omnipresent in Eigen and is crucial to its performance. It is also used to have specialized versions of some decompositions when the matrix's size is known at compile time.
- certik 7y ago
  The main advantage of this approach is that you can use any C++ compiler and the library (in this case Eigen) will work. The disadvantage of this approach is that the code is incomprehensible to average scientists, you must spend a lot of effort to become a C++ expert. The advantage of Fortran is that the code is simple and comprehensible to average scientist, and yet fast, because the effort goes into the compiler itself. I would also argue it is easier to implement an optimization pass in the compiler than implement optimizations on the template level in C++. The disadvantage is that you need a good Fortran compiler. If your compiler can't run your code on, say, GPU, then suddenly there is not a good path forward. While in C++, there is always a way forward via template metaprogramming: perhaps ugly, but at least it gets the job done.