## Sunday, May 2, 2010

### Eigenvalues: C/Lapack

Here is another example on the use of Lapack. This time the objective is to calculate the eigenvectors and eigenvalues of a complex matrix

`#include<stdio.h>#include<math.h>#include<complex.h>#include <stdlib.h>//.......................................................................................................void zgeTranspose( double complex *Transposed, double complex *M ,int n){int i,j;for(i=0;i<n;i++)  for(j=0;j<n;j++) Transposed[i+n*j] = M[i*n+j];}//......................................................................................................//  MatrixComplexEigensystem: computes the eigenvectors and eigenValues of input matrix A//  The eigenvectors are stored in columns//.....................................................................................................void MatrixComplexEigensystem( double complex *eigenvectorsVR, double complex *eigenvaluesW, double complex *A, int N){ int i;    double complex *AT = (double complex*) malloc( N*N*sizeof(double complex) );  zgeTranspose( AT, A , N);    char JOBVL ='N';   // Compute Right eigenvectors  char JOBVR ='V';   // Do not compute Left eigenvectors  double complex VL[1];  int LDVL = 1;   int LDVR = N; int LWORK = 4*N;  double complex *WORK =  (double complex*)malloc( LWORK*sizeof(double complex)); double complex *RWORK = (double complex*)malloc( 2*N*sizeof(double complex));  int INFO; zgeev_( &JOBVL, &JOBVR, &N, AT ,  &N , eigenvaluesW ,   VL, &LDVL,   eigenvectorsVR, &LDVR,    WORK,    &LWORK, RWORK, &INFO ); zgeTranspose( AT, eigenvectorsVR , N); for(i=0;i<N*N;i++) eigenvectorsVR[i]=AT[i];   free(WORK);  free(RWORK);  free(AT);}int main(){  int i,j;  const int N = 3;    double complex A[] = { 1.+I , 2. ,  3 , 4. , 5.+I , 6. , 7., 8., 9. + I};  double complex eigenVectors[N*N];  double complex eigenValues[N];    MatrixComplexEigensystem( eigenVectors, eigenValues, A, N);    printf("\nEigenvectors\n");  for(i=0;i<N;i++){    for(j=0;j<N;j++) printf(" (%f,%f) \t", eigenVectors[i*N + j]);    printf("\n");  }    printf("\nEigenvalues \n");  for(i=0;i<N;i++) printf("\n (%f, %f) \t",  eigenValues[i] );  printf("\n------------------------------------------------------------\n"); return 0;}`

1. A very useful example code that I have been searching for days. Could you also please add a few lines on how to compile this code. Also, are there any specialised functions to diagonalize real symmetric matrices using Lapack in C.
Thanks

2. The Makefile is something like:

#
CC = gcc

#edit LAPACK_PATH if necessary
LAPACK_PATH = /usr/lib64/atlas

a,out: main.c
\$(CC) main.c -L\$(LAPACK_PATH) -llapack -lblas -lgfortran -lm

3. Yes, Lapack contains methods for special kinds of matrices in order to optimize the code. You can find the function you want for example at http://www.netlib.org/lapack/double/.

My example uses zge functions.
z: double complex
ge: general matrices.

4. hi,
i am trying to compile this code but when i am doing
gcc main.c -llapack -lblas -lgfortran -lm
i am getting an error as
/usr/bin/ld: cannot find -llapack
Please if you can guide me in removing this error

Thanks
Rahi

5. Hi,

Your compiler cannot find the lapack library file
which can be
liblapack.so or liblapack.a

make sure you have those files.

You may need to edit
/etc/ld.so.conf
to include the path where the libraries are located
and execute ldconfig

6. Hi!! Thanks!! nice example!!

7. Thank you very much. This is quite useful.