The diagonalization of matrices is done via routines specialized for the appropriate matrix type; this is exemplified by the Numerical Recipes approach, where the cases are classified in the following categories
- real, symmetric, tridiagonal - real, symmetric, banded - real, symmetric - real, nonsymmetric - complex, Hermitian - complex, non-Hermitian
We will need to diagonalize real, symmetric routines. We will use the Householder Method followed by the QR/QL algorithms. The former algorithm transforms the full matrix into a tridiagonal form, while the latter finds the eigenvalues & eigenvectors of the tridiagonal matrix; the eigenvalues & eigenvectors of the original matrix are calculated in a last step.
The routines TRED2 and TQLI are the routines we need.
The theory behind the Householder Method and QR/QL algorithms is explained in NR, Sections 11.0 - 11.3.
The routines TRED2 and TQLI, as well as utility and demonstration codes reside in a directory PHYS307_class_demo. This directory can be downloaded via this link: PHYS307_class_demo.tar.
What will be downloaded is a "tar" file, i.e., a concatenated form of the directory and all of the files it contains. You can expand back the directory via the command
tar -xvf PHYS307_class_demo.tar
This command will expand the directory within the location where your currently are. Be careful to be in the right directory to start with.
Help can be obtained on the Unix "tar" command via "man tar".
The directory you just built contains a "Makefile" file which will allow to build a demonstration code for matrix diagonalization using the NR TRED2 and TQLI routines.
Type:
make demo_diag
This command will build the demo_diag demonstration code.
The "make" Unix command performs pre-established tasks based on date stamps; for instance, in this case it will "compile" and "link" the demonstration main program demo_diag.c together with the NR diagonalization and utility routines. If you were to "make" the code again, nothing would happen. However, if you would "make" the code after modifying the source code, "make" would recognize the existence of more recent pieces of the code, recompile these pieces and link the executable again.
man make
will tell you more about the Unix "make" command.
The code sets a small arbitrary symmetric matrix and diagonalizes it.
The 1-D and 2-D arrays needed for the eigenvectors and symmetric matrix are built according to the NR convention. A 1-D array is first declared as a pointer, followed by a memory space allocation via a call to the NR "vector" routine. The latter sets an array with arbitrary offset (index of first position in array) and size. A 2-D array is first declared as a pointer to pointers (double pointers), followed by a call to the NR "matrix" routine to set aside memory space with arbitrary offset. This is well explained in Chapter 1 of NR.
Each time NR routines are used, the header file "nrutil.h" must be included. The routines "vector" and "matrix" are part of the NR "nrutil.c" library. The "tred2.c" and "tqli.c" are compiled separately (separate files).
The calls to "tred2.c" and "tqli.c" are well documented in NR Sections 11.0 - 11.3.
Table of Content