在进入 DGEEV 时,参数编号 9 具有非法值
On entry to DGEEV parameter number 9 had an illegal value
我第一次尝试使用 C 中的 LAPACK 对角化矩阵,但我被卡住了。
我一直在尝试修改这个例子http://rcabreral.blogspot.co.uk/2010/05/eigenvalues-clapack.html from zgeev
to dgeev
. I have looked at the DGEEV
input parameters, http://www.netlib.org/lapack/explore-html/d9/d28/dgeev_8f.html,但似乎我还不够了解。
因此,下面的代码产生:
**** On entry to DGEEV parameter number 9 had an illegal value**
编辑: 错误发生在 dgeev
跨行 48 到(包括)53 的调用中。
编辑: 请注意,参数与此处的规范不同
http://www.netlib.org/lapack/explore-html/d9/d28/dgeev_8f.html
因为它们已被翻译成指针。在 C 中使用这些 Fortran 例程时,这是必要的,如下所述:
http://www.physics.orst.edu/~rubin/nacphy/lapack/cprogp.html
#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(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 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(i=0;i<N;i++){
for(j=0;j<N;j++) printf("%e", eigenVectors[i*N + j]);
printf("\n");
}
printf("\nEigenvalues \n");
for(i=0;i<N;i++) printf("%e", eigenValues[i] );
printf("\n------------------------------------------------------------\n");
return 0;
}
您不能直接从 zgeev
移植到 dgeev
。 zgeev
获取复数矩阵并计算复数特征值。而 dgeev
得到一个实数矩阵并计算复杂的特征值。为了保持一致 LAPACK
使用 WR
和 WI
用于每个特征值的实部和虚部。
所以请注意 dgeev
定义是
void dgeev_(char* JOBVL, char* JOBVR, int* N, double* A, int* LDA, double* WR, double* WI, double* VL, int* LDVL, double* VR, int* LDVR, double* WORK, int* LWORK, int* INFO);
我对你的例子的建议是删除:
#include <complex.h>
从双打矩阵中删除 I
:
double A[] = { 1. , 2. , 3 , 4. , 5. , 6. , 7., 8., 9.};
然后将特征值向量的大小加倍:
double eigenValues[2*N];
并调用 dgeev
使用 WR
和 WI
:
double *eigenvaluesWR = eigenvaluesW;
double *eigenvaluesWI = eigenvaluesW+N;
dgeev_(&JOBVL, &JOBVR, &N, AT, &N,
eigenvaluesWR, eigenvaluesWI,
VL, &LDVL,
eigenvectorsVR, &LDVR,
WORK, &LWORK, &INFO);
我第一次尝试使用 C 中的 LAPACK 对角化矩阵,但我被卡住了。
我一直在尝试修改这个例子http://rcabreral.blogspot.co.uk/2010/05/eigenvalues-clapack.html from zgeev
to dgeev
. I have looked at the DGEEV
input parameters, http://www.netlib.org/lapack/explore-html/d9/d28/dgeev_8f.html,但似乎我还不够了解。
因此,下面的代码产生:
**** On entry to DGEEV parameter number 9 had an illegal value**
编辑: 错误发生在 dgeev
跨行 48 到(包括)53 的调用中。
编辑: 请注意,参数与此处的规范不同 http://www.netlib.org/lapack/explore-html/d9/d28/dgeev_8f.html 因为它们已被翻译成指针。在 C 中使用这些 Fortran 例程时,这是必要的,如下所述: http://www.physics.orst.edu/~rubin/nacphy/lapack/cprogp.html
#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(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 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(i=0;i<N;i++){
for(j=0;j<N;j++) printf("%e", eigenVectors[i*N + j]);
printf("\n");
}
printf("\nEigenvalues \n");
for(i=0;i<N;i++) printf("%e", eigenValues[i] );
printf("\n------------------------------------------------------------\n");
return 0;
}
您不能直接从 zgeev
移植到 dgeev
。 zgeev
获取复数矩阵并计算复数特征值。而 dgeev
得到一个实数矩阵并计算复杂的特征值。为了保持一致 LAPACK
使用 WR
和 WI
用于每个特征值的实部和虚部。
所以请注意 dgeev
定义是
void dgeev_(char* JOBVL, char* JOBVR, int* N, double* A, int* LDA, double* WR, double* WI, double* VL, int* LDVL, double* VR, int* LDVR, double* WORK, int* LWORK, int* INFO);
我对你的例子的建议是删除:
#include <complex.h>
从双打矩阵中删除 I
:
double A[] = { 1. , 2. , 3 , 4. , 5. , 6. , 7., 8., 9.};
然后将特征值向量的大小加倍:
double eigenValues[2*N];
并调用 dgeev
使用 WR
和 WI
:
double *eigenvaluesWR = eigenvaluesW;
double *eigenvaluesWI = eigenvaluesW+N;
dgeev_(&JOBVL, &JOBVR, &N, AT, &N,
eigenvaluesWR, eigenvaluesWI,
VL, &LDVL,
eigenvectorsVR, &LDVR,
WORK, &LWORK, &INFO);