25 #include "matrixlapack.h" 39 int dgetrf_(
int *m,
int *n,
double *a,
int *lda,
int *ipiv,
int *info);
40 int dgetri_(
int *n,
double *a,
int *lda,
int *ipiv,
double *work,
int *lwork,
int *info);
43 int dsyevd_(
char *jobz,
char *uplo,
int *n,
double *a,
int *lda,
double *w,
double *work,
int *lwork,
int *iwork,
int *liwork,
int *info);
44 int M_DSYEV(
char *jobz,
char *uplo,
int *n,
double *a,
int *lda,
double *w,
double *work,
int *lwork,
int *iwork,
int *liwork,
int *info);
45 int zheevd_(
char *jobz,
char *uplo,
int *n,
double complex *a,
int *lda,
double *w,
double complex *work,
int *lwork,
double *rwork,
int *iwork,
int *liwork,
int *info);
47 int dsyev_(
char *jobz,
char *uplo,
int *n,
double *a,
int *lda,
double *w,
double *work,
int *lwork,
int *info);
48 int M_DSYEV(
char *jobz,
char *uplo,
int *n,
double *a,
int *lda,
double *w,
double *work,
int *lwork,
int *info);
50 int zheev_(
char *jobz,
char *uplo,
int *n,
double complex *a,
int *lda,
double *w,
double complex *work,
int *lwork,
double *rwork,
int *info);
53 int dsyevx_(
char *jobz,
char *range,
char *uplo,
int *n,
double *a,
int *lda,
double *vl,
double *vu,
54 int *il,
int *iu,
double *abstol,
int *m,
double *w,
double *z__,
int *ldz,
55 double *work,
int *lwork,
int *iwork,
int *ifail,
int *info);
98 int n, lda, lwork, info;
102 liwork = 5 * xNsize + 3;
103 iwork = (
int*)malloc(liwork *
sizeof(
double));
109 a = (
double*)malloc(xNsize*xNsize*
sizeof(
double));
110 w = (
double*)malloc(xNsize*
sizeof(
double));
111 work = (
double*)malloc(lwork*
sizeof(
double));
113 to_f(xNsize, xNsize, A, a);
119 M_DSYEV(&jobz, &uplo, &n, a, &lda, w, work, &lwork, &iwork, &liwork, &info);
122 M_DSYEV(&jobz, &uplo, &n, a, &lda, w, work, &lwork, &info);
132 for(k=0;k<xNsize;k++){
154 int DInv(
int xNsize,
double **xM,
double **xIM){
157 int m,n,lda,info,*piv,lwork;
161 m=n=lda=lwork=xNsize;
163 a = (
double*)malloc(xNsize*xNsize*
sizeof(
double));
164 work = (
double*)malloc(xNsize*
sizeof(
double));
165 piv = (
int*)malloc(xNsize*
sizeof(
int));
168 for(j=0;j<xNsize;j++){
169 for(i=0;i<xNsize;i++){
175 dgetrf_(&m, &n, a, &lda, piv, &info);
176 dgetri_(&n, a, &lda, piv, work, &lwork, &info);
185 for(k=0;k<xNsize*xNsize;k++){
186 xIM[k%xNsize][k/xNsize] = a[k];
214 int n, lda, lwork, info;
215 double *a, *w, *work;
218 liwork = 5 * xNsize + 3;
219 iwork = (
int*)malloc(liwork *
sizeof(
double));
220 lwork = 2 * xNsize * xNsize + 6 * xNsize + 1;
227 a = (
double*)malloc(xNsize*xNsize*
sizeof(
double));
228 w = (
double*)malloc(xNsize*
sizeof(
double));
229 work = (
double*)malloc(lwork*
sizeof(
double));
232 for(j=0;j<xNsize;j++){
233 for(i=0;i<xNsize;i++){
243 M_DSYEV(&jobz, &uplo, &n, a, &lda, w, work, &lwork, iwork, &liwork, &info);
246 dsyev_(&jobz, &uplo, &n, a, &lda, w, work, &lwork, &info);
250 for(i=0;i<xNsize;i++){
251 for(j=0;j<xNsize;j++){
264 for(k=0;k<xNsize;k++){
291 int DSEVXU(
int xNsize,
double **A,
double *r,
double **
X,
int xNev){
294 char jobz, range, uplo;
295 int n, lda, il, iu, m, ldz, lwork, *iwork, *ifail, info;
296 double *a, *w, *work, vl, vu, abstol, *z;
298 n = lda = ldz = xNsize;
301 a = (
double*)malloc(xNsize*xNsize*
sizeof(
double));
302 w = (
double*)malloc(xNsize*
sizeof(
double));
303 z = (
double*)malloc(xNsize*xNsize*
sizeof(
double));
304 work = (
double*)malloc(lwork*
sizeof(
double));
305 iwork = (
int*)malloc(5*xNsize*
sizeof(
int));
306 ifail = (
int*)malloc(xNsize*
sizeof(
int));
316 for(j=0;j<xNsize;j++){
317 for(i=0;i<xNsize;i++){
331 dsyevx_(&jobz, &range, &uplo, &n, a, &lda, &vl, &vu, &il, &iu, &abstol,
332 &m, w, z, &ldz, work, &lwork, iwork, ifail, &info);
345 r[k+xNsize-xNev] = w[k];
348 for(k=0;k<xNsize*xNev;k++){
349 X[k%xNsize][k/xNsize+xNsize-xNev] = z[k];
375 int ZHEEVall(
int xNsize,
double complex **A,
double complex *r,
double complex **
vec){
379 int n, lda, lwork, info, lrwork;
382 double complex *a, *work;
385 liwork = 5 * xNsize + 3;
386 iwork = (
int*)malloc(liwork *
sizeof(
double));
391 lwork = xNsize*xNsize + 2 * xNsize;
392 lrwork = 2 * xNsize*xNsize + 5*xNsize + 1;
398 a = (
double complex*)malloc(xNsize*xNsize*
sizeof(
double complex));
399 w = (
double*)malloc(xNsize*
sizeof(
double));
400 work = (
double complex*)malloc(lwork*
sizeof(
double complex));
401 rwork = (
double*)malloc(lrwork*
sizeof(
double));
404 for(j=0;j<xNsize;j++){
405 for(i=0;i<xNsize;i++){
415 int zheevd_(
char *jobz,
char *uplo,
int *n,
double complex *a,
int *lda,
double *w,
double complex *work,
int *lwork,
double *rwork,
int *iwork,
int *liwork,
int *info);
418 zheev_(&jobz, &uplo, &n, a, &lda, w, work, &lwork, rwork, &info);
430 for(i=0;i<xNsize;i++){
431 for(j=0;j<xNsize;j++){
438 for(k=0;k<xNsize;k++){
int dgetri_(int *n, double *a, int *lda, int *ipiv, double *work, int *lwork, int *info)
int M_DSYEV(char *jobz, char *uplo, int *n, double *a, int *lda, double *w, double *work, int *lwork, int *info)
int ZHEEVall(int xNsize, double complex **A, double complex *r, double complex **vec)
obtain eigenvalues and eigenvectors of Hermite matrix A
int DSEVvalue(int xNsize, double **A, double *r)
obtain eigenvalues of real symmetric A
double dlamch_(char *cmach)
void zheevd_(char *jobz, char *uplo, int *n, double complex *a, int *lda, double *w, double complex *work, int *lwork, double *rwork, int *lrwork, int *iwork, int *liwork, int *info)
int dsyevx_(char *jobz, char *range, char *uplo, int *n, double *a, int *lda, double *vl, double *vu, int *il, int *iu, double *abstol, int *m, double *w, double *z__, int *ldz, double *work, int *lwork, int *iwork, int *ifail, int *info)
int DSEVXU(int xNsize, double **A, double *r, double **X, int xNev)
obtain eigenvalues A
int DInv(int xNsize, double **xM, double **xIM)
obtain eigenvalues of inverse of real matrix xM
int zheev_(char *jobz, char *uplo, int *n, double complex *a, int *lda, double *w, double complex *work, int *lwork, double *rwork, int *info)
int DSEVvector(int xNsize, double **A, double *r, double **vec)
obtain eigenvalues and eigenvectors of real symmetric A
void to_f(int N, int M, double **A, double *a)
function for transforming Row-major matrix (C) to Column-major matrix (Fortran)
int dgetrf_(int *m, int *n, double *a, int *lda, int *ipiv, int *info)
int dsyev_(char *jobz, char *uplo, int *n, double *a, int *lda, double *w, double *work, int *lwork, int *info)