/*    ATK IV - Numerical Programming (2006)      */
/*  Excercise 9.2 -  Eigenvalue Problems Cont'd  */

#include<stdio.h>
#include<math.h>
#define NRANSI
#include<nr.h>
#include<nrutil.h>

int main(void){
  float **a;
  float *d,*e;
  int n,nrot;
  float iprod=0.0;

  FILE *file;
  int i,j;

  file=fopen("matrix_sym.dat","rt");
  fscanf(file,"%d\n",&n);

  a=matrix(1,n,1,n);

  d=vector(1,n);
  e=vector(1,n);

  for(i=1;i<=n;i++){
    for(j=1;j<=n;j++){
      fscanf(file,"%f",&a[i][j]);
    }    
  }

  tred2(a,n,d,e);
  tqli(d,e,n,a);
  eigsrt(d,a,n);
  printf("\n");
  for(i=1;i<=n;i++){
    printf("Eigenvalue: %f\n",d[i]);
    printf("Corresponding eigenvector:  ");
    for(j=1;j<=n;j++){
      printf("%f ",a[j][i]);
    }
    printf("\n");
  }

  for(i=1;i<=n;i++){
    iprod += a[i][1]*a[i][2];
  }
  printf("\nInner product between 1st and 2nd eigenvector: %f\n",iprod);
  return 0;

}