- #1

Angelos K

- 48

- 0

Mod note: I revised the code below slightly, changing the loop control variable i to either j or k. The reason for this is that the browser mistakes the letter i in brackets for the BBCode italics tag, which causes some array expressions to partially disappear.

Hello,

I am trying for the first time to use LAPACK from C to diagonalize a matrix and I am stuck.

I try to modify this example http://rcabreral.blogspot.co.uk/2010/05/eigenvalues-clapack.html from zgeev to degeev. I have looked at the DGEEV input parameters, but it seems I don't understand the well enough.

Below is my code. Running it produces:

** On entry to DGEEV parameter number 9 had an illegal value

EDIT: The error occurs in the call of degeev spanning lines 48 to (including) 53.

Hello,

I am trying for the first time to use LAPACK from C to diagonalize a matrix and I am stuck.

I try to modify this example http://rcabreral.blogspot.co.uk/2010/05/eigenvalues-clapack.html from zgeev to degeev. I have looked at the DGEEV input parameters, but it seems I don't understand the well enough.

Below is my code. Running it produces:

** On entry to DGEEV parameter number 9 had an illegal value

EDIT: The error occurs in the call of degeev spanning lines 48 to (including) 53.

#include<stdio.h>

#include<math.h>

#include<complex.h>

#include <stdlib.h>

//...................

void dgeTranspose( double *Transposed, double *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 *eigenvectorsVR, double *eigenvaluesW, double *A, int N)

{

int i;

double *AT = (double *) malloc( N*N*sizeof(double ) );

dgeTranspose( AT, A , N);

char JOBVL ='N'; // Compute Right eigenvectors

char JOBVR ='V'; // Do not compute Left eigenvectors

double VL[1];

int LDVL = 1;

int LDVR = N;

int LWORK = 4*N;

double *WORK = (double *)malloc( LWORK*sizeof(double));

double *RWORK = (double *)malloc( 2*N*sizeof(double));

int INFO;

dgeev_( &JOBVL, &JOBVR, &N, AT , &N , eigenvaluesW ,

VL, &LDVL,

eigenvectorsVR, &LDVR,

WORK,

&LWORK, RWORK, &INFO );

dgeTranspose( AT, eigenvectorsVR , N);

for(j=0;j<N*N;j++) eigenvectorsVR[j]=AT[j];

free(WORK);

free(RWORK);

free(AT);

}int main()

{

int i,j;

const int N = 3; double A[] = { 1.+I , 2. , 3 , 4. , 5.+I , 6. , 7., 8., 9. + I};

double eigenVectors[N*N];

double eigenValues[N]; MatrixComplexEigensystem( eigenVectors, eigenValues, A, N);

printf("\nEigenvectors\n");

for(k=0;k<N;k++){

for(j=0;j<N;j++) printf("%e", eigenVectors[k*N + j]);

printf("\n");

}

printf("\nEigenvalues \n");

for(j=0;i<N;j++) printf("%e", eigenValues[j] );printf("\n------------------------------------------------------------\n");

return 0;}

Any help much appreciated!#include<math.h>

#include<complex.h>

#include <stdlib.h>

//...................

void dgeTranspose( double *Transposed, double *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 *eigenvectorsVR, double *eigenvaluesW, double *A, int N)

{

int i;

double *AT = (double *) malloc( N*N*sizeof(double ) );

dgeTranspose( AT, A , N);

char JOBVL ='N'; // Compute Right eigenvectors

char JOBVR ='V'; // Do not compute Left eigenvectors

double VL[1];

int LDVL = 1;

int LDVR = N;

int LWORK = 4*N;

double *WORK = (double *)malloc( LWORK*sizeof(double));

double *RWORK = (double *)malloc( 2*N*sizeof(double));

int INFO;

dgeev_( &JOBVL, &JOBVR, &N, AT , &N , eigenvaluesW ,

VL, &LDVL,

eigenvectorsVR, &LDVR,

WORK,

&LWORK, RWORK, &INFO );

dgeTranspose( AT, eigenvectorsVR , N);

for(j=0;j<N*N;j++) eigenvectorsVR[j]=AT[j];

free(WORK);

free(RWORK);

free(AT);

}int main()

{

int i,j;

const int N = 3; double A[] = { 1.+I , 2. , 3 , 4. , 5.+I , 6. , 7., 8., 9. + I};

double eigenVectors[N*N];

double eigenValues[N]; MatrixComplexEigensystem( eigenVectors, eigenValues, A, N);

printf("\nEigenvectors\n");

for(k=0;k<N;k++){

for(j=0;j<N;j++) printf("%e", eigenVectors[k*N + j]);

printf("\n");

}

printf("\nEigenvalues \n");

for(j=0;i<N;j++) printf("%e", eigenValues[j] );printf("\n------------------------------------------------------------\n");

return 0;}

Last edited by a moderator: