Background

It is not so hard to implement routines to find the eigenvalues and eigenvectors for real symmetric matrices, or to only calculate the eigenvalues for real non-symmetric matrices. You can find the solution in most of textbooks such as Numerical Recipes (NR). However, it is much harder to calculate eigenvectors for real dense non-symmetric matrices. C and Fortran programmers usually reply on ATLAS, LAPACK or GSL. The problem is frequently we just need a convenient routine to solve a small eigensystem, while we have to depend on a large library which is not always easy to compile and install. This is especially troublesome when we want others to use our source codes. I always prefer standalone single/double-file C libraries if they are efficient enough in a particular application.


A Failure

After quite some google searches, I found eigen.c in the PDL Perl module. Unfortunately, this piece of code does not work. On a small example, it segfaults and I do not know how to fix that. In case I was doing something stupid, I put my testing code at the end of this page.


A Success

Another google search directed me to feigen.c which is contained in the CD coming with the book Numerical Algorithms with C (NAC) by Frank Uhlig and his colleagues. This routine does work. Like always, I merged several source code files into one and redifined the API. My modified source code can be downloaded as eigeng.c. The example in the end of the file largely describes how to use the routine.


A Good Book

NAC is a truely amazing book. Although it is less descriptive and covers less topics than NR (no minimization for example), it provides some more advanced methods and the coding quality is also better in my view. NAC is a good complement to NR. Do not forget to read this book when NR does not meet your goal.

By the way, NAC gives an implementation of the order-4 Gear method for solving stiff ODEs. This method potentially outperforms other one-step methods for stiff equations (see also this page). I may test its performance when I have time. Currently LSODA looks quite satisfying to me.


Appendix: Testing Codes for PDL