#include "include/stdlib.h"
#include "include/stdio.h"
#include "include/time.h"      /*Bibliotecas a incluir para el uso de función gettimeofday*/
#include "include/unistd.h"
#include "include/fcntl.h"
#include "include/assert.h"
#include <math.h>
#define MAXSIZE 1000
#define N 40      /* Nº de iteraciones para cada tamaño de arreglo */
#define SWAP(x, y) temp = (x); (x) = (y); (y) = temp



int a,b,e,i,j,k,t=0,d;
long Sumaheap,Sumaquick,Sumaradix,Promedioheap,Promedioquick,Promedioradix,key,antes,despues;
float Desvradix,Desvheap,Desvquick,Eradix,Equick,Eheap,E2radix,E2quick,E2heap,sop;
typedef int keytype;
FILE *fp,*fp1;                 	/*Puntero a archivo para utilización con fprintf()*/

struct timeval timereal;  /*Declaración de estructuras tipo timeval y timezone, estos tipos son predefinidos*/
struct timezone timez;

struct timeval *tv;      /*Declaración de punteros a estructuras para utilizarlos como parámetros en función gettimeofday() */
struct timezone *tz;




main(int argc,char *argv[])

{

fp=fopen( "gnuplot.dat","w");
fp1=fopen( "desv.dat","w");
tv=&timereal;            /*Punteros apuntando a las dos estructuras anteriormente definidas*/
tz=&timez;

if(argc!=2){printf("Revisar parámetros en la llamada al programa\n");exit(1);} 		/* Se verifica número de parámetros entregados al llamar
											la función, en caso de error informa y sale del programa*/	

a=atoi(argv[1]); 					/*Convierte el argumento entregado al llamar la función a valor entero*/

if((a>MAXSIZE)||(a==0)){printf("El tamaño del arreglo debe ser menor o igual a 1000\n");
exit(1);}  /* Se verifica el cumplimiento de las restricciones para el argumento(nº entero entre 1 y 1000, en caso de error informa y sale*/

 d=a/30;
 if(d<=1)b=1;else b=(int)d;  /* La variable b corresponde al tamaño del paso entre puntos consecutivos de la gráfica,aproxima
 				por el entero superior*/
 
 
 
for(i=b;i<=a;i=i+b) /* Se toma 30 (arreglos de distinto tamaño)puntos para cada gráfica, este número no es fijo pues depende de la aproximación de la variable b al transformarla a entero+/ */
		{
		Eheap=Equick=Eradix=E2heap=E2quick=E2radix=0;		
		for(k=0;k<N;k++) /*Se calcula un tiempo promedio en base a 40  repeticiones de cada ordenamiento para arreglos distintos*/
				{
				long A[i];
				for(j=0;j<i;j++)/*Se llena el arreglo en forma aleatoria*/
					{
					A[j]=random();
					}
					e=gettimeofday(tv,tz); /*Se llama a la función gettimeofday para obtener el "tiempo"*/
					if(e==0)antes=(*tv).tv_usec+(*tv).tv_sec*1000000;else exit(1);	/*Toma el valor inicial en microsegundos y sale en caso de error en la llamada a función*/	
                                        heapsort(i,A);
                                        e=gettimeofday(tv,tz);/*Se llama a la función gettimeofday para obtener el "tiempo"*/
					if(e==0)despues=(*tv).tv_usec+(*tv).tv_sec*1000000;else exit(1);	/*Toma el valor final en microsegundos y sale en caso de error en la llamada a función*/	
                                        sop=(float)despues-(float)antes;
                                        Sumaheap+=sop;
                                        Eheap+=(sop/N);
                                        E2heap+=(sop*sop)/N;

                                        e=gettimeofday(tv,tz); /*Se llama a la función gettimeofday para obtener el "tiempo"*/
					if(e==0)antes=(*tv).tv_usec+(*tv).tv_sec*1000000;else exit(1);	/*Toma el valor inicial en microsegundos y sale en caso de error en la llamada a función*/	
                                        quicksort(A,0,i);
                                        e=gettimeofday(tv,tz);/*Se llama a la función gettimeofday para obtener el "tiempo"*/
					if(e==0)despues=(*tv).tv_usec+(*tv).tv_sec*1000000;else exit(1);	/*Toma el valor final en microsegundos y sale en caso de error en la llamada a función*/	
                                        sop=(float)despues-(float)antes;
                                        Sumaquick+=sop;
                                        Equick+=(sop/N);
                                        E2quick+=(sop*sop)/N;


                                        e=gettimeofday(tv,tz); /*Se llama a la función gettimeofday para obtener el "tiempo"*/
					if(e==0)antes=(*tv).tv_usec+(*tv).tv_sec*1000000;else exit(1);	/*Toma el valor inicial en microsegundos y sale en caso de error en la llamada a función*/	
                                        radixsort(A,i);
                                        e=gettimeofday(tv,tz);/*Se llama a la función gettimeofday para obtener el "tiempo"*/
					if(e==0)despues=(*tv).tv_usec+(*tv).tv_sec*1000000;else exit(1);	/*Toma el valor final en microsegundos y sale en caso de error en la llamada a función*/	
                                        sop=(float)despues-(float)antes;
                                        Sumaradix+=sop;
                                        Eradix+=(sop/N);
                                        E2radix+=(sop*sop)/N;

/*********************************************Cálculo de parámetros para calcular tiempo promedio y desviación standard*****************************/
                }
		Promedioquick=Sumaquick/N;
		Promedioradix=Sumaradix/N;
		Promedioheap=Sumaheap/N;
		
		fprintf( fp,"%d\t",i);                  /*Escritura a archivo que contiene los tiempos promedio de ejecución*/
		fprintf( fp,"%d\t",Promedioheap);
		fprintf( fp,"%d\t",Promedioquick);
		fprintf( fp,"%d\n",Promedioradix);
		
		Desvheap=sqrt(E2heap-(Eheap*Eheap));
		Desvquick=sqrt(E2quick-(Equick*Equick));
		Desvradix=sqrt(E2radix-(Eradix*Eradix));
		fprintf( fp1,"%d\t",i);                 /*Escritura a archivo que contiene los valores obtenidos de desv. standard para los distintos tamaños de arreglo y los diferentes algoritmos de ordenamiento*/
		fprintf( fp1,"%f\t",Desvheap);
		fprintf( fp1,"%f\t",Desvquick);
		fprintf( fp1,"%f\n",Desvradix);
		
				}
		
	fclose(fp);
	fclose(fp1);
	system("gnuplot plot.gnu"); /* Se grafica en Gnuplot los distintos valores obtenidos*/
	system("gnuplot plot1.gnu");
	exit(1);

    }

                       /*Funciones de ordenamiento*/
/*************************************************************************************/
int quicksort(r,lo,up)
  int lo,up,r[up];
  {
  int i, j,tempr;
     while (up>lo) {
          i = lo;
          j = up;
          tempr = r[lo];
          while (i<j){
               for (;r[j]>tempr;j-- );
               for ( r[i]=r[j];i<j&&r[i]<=tempr;i++);
               r[j]=r[i];
               }
          r[i]=tempr;
                 if(i-lo<up-i){quicksort(r,lo,i-1);lo=i+1;}
               else {quicksort(r,i+1,up);up=i-1;}
          }
     }
/*************************************************************************************/
int heapsort(int n, keytype a[])
{
	int i, j, k;
	keytype x;

	for (k =n/2;k>=1;k--) {
		i = k; x=a[i];
		while ((j=2*i)<=n){
			if (j<n&&a[j]<a[j+1])j++;
			if (x>=a[j])break;     /* En caso de cumplirse la condición se "rompe" el ciclo*/
			a[i]=a[j];i=j;
		}
		a[i]=x;
	}
	while (n>1){
		x=a[n];a[n]=a[1];n--;
		i=1;
		while((j=2*i)<=n){
			if (j<n&&a[j]<a[j+1])j++;
			if (x>=a[j])break;     /* En caso de cumplirse la condición se "rompe" el ciclo*/
			a[i]=a[j];i=j;
		}
		a[i]=x;
	}
}

/***************************************************************************************/
#define         DIGITOS         10
#define         SHIFT           1
#define         BASE            (1<<SHIFT)
int radixsort(int ar[],int n)  /* Radix consta de 2 partes  radix propiamente tal y la subrutina Lineal*/
        {       int i;
        for(i=0;i<DIGITOS;++i)
                {Lineal(ar,n,i);  }
        }

int Lineal(int ar[],int n,int digit)
        {       int cuenta[BASE];  int i;  int *w;  int sum;
        for(i=0;i<BASE;++i)
                {  cuenta[i]=0; }

        for(i=0;i<n;++i)
                {cuenta[(ar[i]>>(SHIFT*digit))%(BASE)]++;}

        sum=cuenta[0]; cuenta[0]=0;
        for(i=1;i<BASE;++i)
                {int tmp;
                tmp=cuenta[i]; cuenta[i]=sum; sum+=tmp;
                }


        w=(int *)malloc(sizeof(int)*n); assert(w!=0);   /*Asigna memoria y deja a w apuntando a la memoria asignada*/

        for(i=0;i<n;++i)
                {w[cuenta[(ar[i]>>(SHIFT*digit))%BASE]++]=ar[i]; }
        for(i=0;i<n;++i)
                {ar[i]=w[i];}

        free(w); /*Libera la memoria asignada y apuntada por w*/
        }