// Note to self - Counted out! is followed by a core dump if encountered.
#include "main.h"

float *X,*Y,*Z;
sparse_matrix_arr* L;
sparse_vector **svb;
float *_vertexarr, *orgnl_va, precision, *eigs;
int  *_facearr, *orgnl_fa, _nvertex, orgnl_nv, _nface, orgnl_nf;int memorized_thingies;


int main(int argc, char* argv[])
{
	//int i,j,c;
//	int count = 0;
	
	srand( (unsigned)time( NULL ) );
	memorized_thingies = 0; //todo: default inits?
	
	if ( argc != 3 )
	{
		printf("KAKA!\n");
		exit(1);
	}
		
	parse_from_file(argv[1]);
	precision = (atof(argv[2]));
	L = Laplacian(orgnl_fa, _nface, _nvertex);
	matrix_shifting(L);

	X = (float*) calloc( L->n , sizeof(float) );
	Y = (float*) calloc( L->n , sizeof(float) );
	Z = (float*) calloc( L->n , sizeof(float) );
	
	svb = (sparse_vector**) calloc( L->n , sizeof(sparse_vector*) );
	eigs = (float*) calloc( L->n , sizeof(float) );

	printf("\n");	
	
	startLib();
	
	return 0;
}


void parse_from_file(char *address)
{
	char ch;
	int i;
	FILE *fp;
	
	fp = fopen(address,"r");
	
	free(orgnl_va); free(orgnl_va); // Could this be removed?
	free(_vertexarr); free(_facearr);
	
	while ( (ch = getc(fp)) != '\n' )
		;
		
	fscanf(fp,"%d ",&orgnl_nv);
	fscanf(fp,"%d ",&orgnl_nf);
	fscanf(fp,"%d\n",&i);
	
	orgnl_va = (float*) calloc( orgnl_nv * 3 , sizeof(float) );
	orgnl_fa = (int*)   calloc( orgnl_nf * 3 , sizeof(int) );
	
	for ( i = 0 ; i < (orgnl_nv * 3) ; i += 3 )
	{
		fscanf(fp,"%f ",  &(orgnl_va[i + 0]));
		fscanf(fp,"%f ",  &(orgnl_va[i + 1]));
		fscanf(fp,"%f\n", &(orgnl_va[i + 2]));
	}
	
	
	for ( i = 0 ; i < (orgnl_nf * 3) ; i += 3 )
	{
		
		fscanf(fp,"%c ",&ch); // Expected to be three.
		fscanf(fp,"%d ",  &(orgnl_fa[i + 0]));
		fscanf(fp,"%d ",  &(orgnl_fa[i + 1]));
		fscanf(fp,"%d\n", &(orgnl_fa[i + 2]));
	}
	
	_nvertex = orgnl_nv;
	_nface = orgnl_nf;
	_vertexarr = orgnl_va;
	_facearr = orgnl_fa;
	
	
	fclose(fp);
}




/*
void parse_from_stdin()
{
	char ch;
	int i;
	
	free(_vertexarr); free(_facearr);

	while ( (ch = getchar()) != '\n' )
		;
		
	scanf("%d ",&_nvertex);
	scanf("%d ",&_nface);
	scanf("%d\n",&i);
	
	_vertexarr = (float*) calloc( _nvertex * 3 , sizeof(float) );
	_facearr = (int*)   calloc( _nface * 3 , sizeof(int) );
	
	for ( i = 0 ; i < (_nvertex * 3) ; i += 3 )
	{
		scanf("%f ",  &(_vertexarr[i]));
		scanf("%f ",  &(_vertexarr[i + 1]));
		scanf("%f\n", &(_vertexarr[i + 2]));
	}
	
	
	for ( i = 0 ; i < (_nface * 3) ; i += 3 )
	{
		
		scanf("%c ",&ch); // Expected to be three.
		scanf("%d ",  &(_facearr[i]));
		scanf("%d ",  &(_facearr[i + 1]));
		scanf("%d\n", &(_facearr[i + 2]));
	}	
}

*/



void charge_arrays(float **vertexarr, int *nvertex, 
                   int **facearr, int *nface)
{
	*nvertex = _nvertex;
	*nface = _nface;
	
	*vertexarr = _vertexarr;
	*facearr = _facearr;
}


void myexit()
{
  free(_vertexarr);
  free(_facearr);
}


void  compute(int nbasis)
{
	clock_t timer = clock();
	int i,j,t;
	
	if ( nbasis < 1 || nbasis > orgnl_nv )
	{
		printf("Throwing a tantrum - compute\n");
		return;
	}
	
	else if ( memorized_thingies < nbasis )
	{
		our_gram_schmidt(L, precision , svb, eigs, nbasis, memorized_thingies);
		
		for ( i = memorized_thingies ; i < nbasis ; i++ )
			for ( j = 0 ; j < svb[i]->nnz ; j++ )
			{
				t = (svb[i]->indices)[j];
				X[i] += orgnl_va[3*t+0] * (svb[i]->values)[j];
				Y[i] += orgnl_va[3*t+1] * (svb[i]->values)[j];
				Z[i] += orgnl_va[3*t+2] * (svb[i]->values)[j];
			}
			
		memorized_thingies = nbasis;
	}
	
	printf("nbasis was %d, the time compute needed was %d. Note to selves - should be in miliseonds.\n" 
		, nbasis , (int) (clock() - timer));
}


void  compress(int ncoeff)
{
	int i,j,k,t;
	
	if ( ncoeff < 1 || ncoeff > memorized_thingies )
	{
		printf("Throwing a tantrum - compress\n");
		return;
	}
	
	if ( _vertexarr != orgnl_va )
	{
		free(_vertexarr);
		_vertexarr = (float*) calloc( 3 * _nvertex , sizeof(float) );
	}
	else
		for ( i = 0 ; i < 3 * _nvertex ; i++ )
			_vertexarr[i] = 0.0f;
	
	for ( i = 0 ; i < ncoeff ; i++ )
	{
		for ( j = 0 ; j < svb[i]->nnz ; j++ )
		{
			//TODO: Rewrite!
			t = (svb[i]->indices)[j];
			_vertexarr[3*t+0] += X[i] * (svb[i]->values)[j];
			_vertexarr[3*t+1] += Y[i] * (svb[i]->values)[j];
			_vertexarr[3*t+2] += Z[i] * (svb[i]->values)[j];
		}
	}
}