[help]Data type in zheev.c example.

Shi__Yue · ‎02-09-2018

Hi,

In the example of

https://software.intel.com/sites/products/documentation/doclib/mkl_sa/11/mkl_lapack_examples/zheev_ex.c.htm

The to-be-diagonalized matrix is declared to be a dcomplex type

(struct _dcomplex { double re, im; };
typedef struct _dcomplex dcomplex;)

and initialized using:

dcomplex a[LDA*N] = {
           { 9.14,  0.00}, {-4.37,  9.22}, {-1.98,  1.72}, {-8.96,  9.50},
           { 0.00,  0.00}, {-3.35,  0.00}, { 2.25,  9.51}, { 2.57, -2.40},
           { 0.00,  0.00}, { 0.00,  0.00}, {-4.82,  0.00}, {-3.24, -2.04},
           { 0.00,  0.00}, { 0.00,  0.00}, { 0.00,  0.00}, { 8.44,  0.00}
        };

When I try to supply my own matrix elements to a[LDA*N] using

dcomplex a[LDA*N];
        for (int i=0;i<N;i++) {
        for (int j=0;j<LDA;j++) {
        int k=j+i*4*N;
        a = {std::real(hamiltonian),std::imag(hamiltonian)};
        }//j
        }//i

It complains:

error: expected an expression
a = {std::real(hamiltonian),std::imag(hamiltonian)};
^

when I tried to icc the code using mkl link.

I guess the question is: how to initialize/access members in a[LDA*N] in a more free way.

As you see I mostly use vectors to store my data. So I am not familiar with other data types.

Thanks a lot,

Yue

Shi__Yue · ‎02-09-2018

one update:

I have replace the initialization part with
a[0]={9.14,0.00};
a[1]={-4.37,9.22};
...
The code is appended.

The error persists for each declearation.
icc zheev_test00.cpp -mkl -o test
zheev_test00.cpp(93): error: expected an expression
a[0]={ 9.14, 0.00};
^

zheev_test00.cpp(94): error: expected an expression
a[1]={-4.37, 9.22};
^

zheev_test00.cpp(95): error: expected an expression
a[2]={-1.98, 1.72};
^

Thanks,
Yue

{\code}
#include <stdlib.h>
#include <stdio.h>

/* Complex datatype */
struct _dcomplex { double re, im; };
typedef struct _dcomplex dcomplex;

/* ZHEEV prototype */
extern "C" {
void zheev( char* jobz, char* uplo, int* n, dcomplex* a, int* lda,
double* w, dcomplex* work, int* lwork, double* rwork, int* info );
}
/* Auxiliary routines prototypes */
extern void print_matrix( char* desc, int m, int n, dcomplex* a, int lda );
extern void print_rmatrix( char* desc, int m, int n, double* a, int lda );

/* Parameters */
#define N 4
#define LDA N

/* Main program */
int main() {
/* Locals */
int n = N, lda = LDA, info, lwork;
dcomplex wkopt;
dcomplex* work;
/* Local arrays */
/* rwork dimension should be at least max(1,3*n-2) */
double w, rwork[3*N-2];
/*
dcomplex a[LDA*N] = {
{ 9.14, 0.00}, {-4.37, 9.22}, {-1.98, 1.72}, {-8.96, 9.50},
{ 0.00, 0.00}, {-3.35, 0.00}, { 2.25, 9.51}, { 2.57, -2.40},
{ 0.00, 0.00}, { 0.00, 0.00}, {-4.82, 0.00}, {-3.24, -2.04},
{ 0.00, 0.00}, { 0.00, 0.00}, { 0.00, 0.00}, { 8.44, 0.00}
};
*/
dcomplex a[LDA*N];
a[0]={ 9.14, 0.00};
a[1]={-4.37, 9.22};
a[2]={-1.98, 1.72};
a[3]={-8.96, 9.50};
a[4]={ 0.00, 0.00};
a[5]={-3.35, 0.00};
a[6]={ 2.25, 9.51};
a[7]={ 2.57, -2.40};
a[8]={ 0.00, 0.00};
a[9]={ 0.00, 0.00};
a[10]={-4.82, 0.00};
a[11]={-3.24, -2.04};
a[12]={ 0.00, 0.00};
a[13]={ 0.00, 0.00};
a[14]={ 0.00, 0.00};
a[15]={ 8.44, 0.00};
/* Executable statements */
printf( " ZHEEV Example Program Results\n" );
/* Query and allocate the optimal workspace */
lwork = -1;
zheev( "Vectors", "Lower", &n, a, &lda, w, &wkopt, &lwork, rwork, &info );
lwork = (int)wkopt.re;
work = (dcomplex*)malloc( lwork*sizeof(dcomplex) );
/* Solve eigenproblem */
zheev( "Vectors", "Lower", &n, a, &lda, w, work, &lwork, rwork, &info );
/* Check for convergence */
if( info > 0 ) {
printf( "The algorithm failed to compute eigenvalues.\n" );
exit( 1 );
}
/* Print eigenvalues */
print_rmatrix( "Eigenvalues", 1, n, w, 1 );
/* Print eigenvectors */
print_matrix( "Eigenvectors (stored columnwise)", n, n, a, lda );
/* Free workspace */
free( (void*)work );
exit( 0 );
} /* End of ZHEEV Example */

/* Auxiliary routine: printing a matrix */
void print_matrix( char* desc, int m, int n, dcomplex* a, int lda ) {
int i, j;
printf( "\n %s\n", desc );
for( i = 0; i < m; i++ ) {
for( j = 0; j < n; j++ )
printf( " (%6.2f,%6.2f)", a[i+j*lda].re, a[i+j*lda].im );
printf( "\n" );
}
}

/* Auxiliary routine: printing a real matrix */
void print_rmatrix( char* desc, int m, int n, double* a, int lda ) {
int i, j;
printf( "\n %s\n", desc );
for( i = 0; i < m; i++ ) {
for( j = 0; j < n; j++ ) printf( " %6.2f", a[i+j*lda] );
printf( "\n" );
}
}
{\code}

mecej4 · ‎02-09-2018

I am not sure what you are trying to do, and you have not stated which C++ standard you/your compiler is following. Perhaps, what you need is the "brace init" of C++11, as described in http://www.informit.com/articles/article.aspx?p=1852519 . Try the following code for the array initialization:

dcomplex a[LDA*N] { 
{ 9.14, 0.00}, {-4.37, 9.22}, {-1.98, 1.72}, {-8.96, 9.50},
{ 0.00, 0.00}, {-3.35, 0.00}, { 2.25, 9.51}, { 2.57, -2.40},
{ 0.00, 0.00}, { 0.00, 0.00}, {-4.82, 0.00}, {-3.24, -2.04},
{ 0.00, 0.00}, { 0.00, 0.00}, { 0.00, 0.00}, { 8.44, 0.00}
};

Note that there is no '=' between the variable and the initialization expression.