Friday, December 4, 2009

Computational Linear Algebra

How do we implement programs that require linear algebra?
For most of the cases there are no reasons to spent effort in developing our own libraries for linear algebra because we can use very good and freely available libraries.
  • The first layer is the Fortran BLAS library that implements matrix products. There are many variants of BLAS developed by certain companies and research institutions.
  • The second layer is the Fortran Lapack library, which implements more advanced routines such as matrix decompositions and computation of eigenvalues on the top of BLAS. However, there are no high level routines such as matrix inverse, determinants, etc. The reason of the absence of these high level routines is that there are many ways to implement them in terms of the Lapack routines but one has to choose a particular method according to the personal requirements for maximum efficiency. The implementation of programs in Lapack has high potential for efficiency optimization because one can specifies the type of matrix in the operation ie. real, complex, symmetric, etc and the decomposition that fits better for the final purpose.
  • The third layer is the Blas/Lapack wrapper, which implements the operations that are taught in a first curse of Linear algebra. The routines in Lapack are usually taught in a second course of linear algebra. The wrappers are easier to use but some of the flexibility is lost.
  • For those who work with C, GSL must be the best option. For C++, one has the choice to use the Boost library, which includes a lot more than linear algebra. Personally, the library I really feel as the most friendly is the Armadillo C++ library for linear algebra.

1 comment:

  1. I have found that Eigen 2 ( is a very good linear algebra library. It uses C++ templates, is extremely easy to use, and is very fast! Give it a spin when you have time :)