#include "header.h"

int main(void)
{
	gsl_matrix *A,*B;
	gsl_vector *e,*ans;
	float tol;
	int maxite,n,m,flag=0,i;

	scanf("%d",&n);
	scanf("%d",&m);
	A=gsl_matrix_alloc(n,m);   
	B=gsl_matrix_alloc(n,m);   
	gsl_matrix_fscanf (stdin,A);
	gsl_matrix_memcpy(B,A);
        // Use the QR decomp to find the eigenvalues
        e=gsl_vector_alloc(n);
        ans=gsl_vector_alloc(n); 
	gsl_vector_set_zero(e);
	tol=TOL;                      //control the error
	maxite=MAXITE;                //control the iterate
	flag=QREvals(A,e,maxite,tol); //e is the vector contains eigenvalues of A 
	gsl_sort_vector(e);           //sort the eigenvalues into ascending numerical order
	for(i=0;i<n;i++)
	{
		 gsl_vector_set(ans,i,gsl_vector_get(e,n-1-i)); 
	}

	if(flag>0)
	{
		printf("QR Eigenvalues of Matrix A is: \n");
		gsl_vector_fprintf(stdout,ans,"%lf");
	}
        //Compare with gsl build-in function
        gsl_eigen_symmv_workspace * w = gsl_eigen_symmv_alloc (n);
        gsl_vector *evls = gsl_vector_alloc (n);
        gsl_matrix *evec = gsl_matrix_alloc (n, n);
        gsl_eigen_symmv (B, evls, evec, w);
	printf("comparing with eigenvalues calculated by using gsl function:\n");
	gsl_vector_fprintf(stdout,evls,"%lf");
	

	gsl_vector_free(e);
	gsl_vector_free(ans);
	gsl_matrix_free(A);
	gsl_matrix_free(B);
	return 0;
}