#include <numeric.h>
#define N 2

int main() {
array double a[2][2] = {3,      2,
                        1,      4};
     array double evalues[N];
     array double evectors[N][N];
     array double *v, vv[N][N];
     int i, status;

     printf("a[0][0] %f a[0][1] %f\n", a[0][0],a[0][1]);
     printf("a a =\n%e\n", a);
                        
     eigensystem(evalues, evectors, a);
     printf("eigennvalues of a =\n%e\n", evalues);
     printf("eigenvectors of a =\n%e\n", evectors);

     vv = transpose(evectors);
       for(i=0; i<N; i++) {
          v = (array double [:])(double [N])&vv[i][0];
          printf("V %e\n", v);
          printf("VV start %e\n", vv[i][0]);
          printf("A*v-evalues[i]*v = %f\n", a*v-evalues[i]*v);
       }
       

}