#include <stdio.h>
#include <stdint.h>
#include <stdbool.h>
#include <stdlib.h>

#define IN_FILENAME 	"gauss.in"
//#define IN_FILENAME 	"./teste/grader_test6.in"

#define OUT_FILENAME	"gauss.out"

#define EPS 		0.0000001

bool readInputData(
	double 		**A,
	int32_t 	*n,
	int32_t		*m)
{
	bool bRes = 0;
	FILE *file;
	
	file = fopen(IN_FILENAME, "r");

	if (file==NULL)
		bRes = 1;

	if (bRes == 0)
	{
		fscanf(file, "%d %d", n, m);

		*A = (double*)malloc((*n) * ((*m)+1) * sizeof(double));

		if (*A == NULL)
			bRes = 1;
	}

	if (bRes == 0)
	{			
		int32_t i;
		double *pA;
		pA = *A;
		for (i=0; i< (*n)*((*m)+1); i++)
			fscanf(file, "%lf", pA++);
	}
	
	if (file!=NULL)
		fclose(file);

	return bRes;
}

bool writeOutputData(
	double 		*x,
	int32_t 	m,
	bool 		hasSolution)
{
	bool bRes = 0;
	FILE *file;

	file = fopen(OUT_FILENAME, "w");

	if (file==NULL)
		bRes = 1;

	if (bRes == 0)
	{
		if (hasSolution == true)
		{
			int32_t i;
			double *pX = x;
			for(i=0; i<m; i++)
				fprintf(file, "%15.10lf ", *pX++);
			fprintf(file, "\n");
		}
		else
		{
			fprintf(file, "Imposibil\n");
		}
	}

	if (file!=NULL)
		fclose(file);

	return bRes;
}

void displaySystem(
	double 		*A, 
	int32_t 	n,
	int32_t 	m)
{
	int32_t i, j;
	for(i=0; i< n; i++)
	{
		for(j=0; j<m+1; j++)
		{
			printf("%15.10lf ", A[i * (m+1) +j]);
		}
		printf("\n");
	}
	printf("\n");
}

bool checkSystem(
	double 		*A,
	double 		*x,
	int32_t 	n,
	int32_t 	m)
{
	bool bRes = true;
	int32_t i, j;
	printf("%d\n", n);
	for(i = 0; i < n && bRes == true; i++)
	{
		double val = 0;
		for(j=0; j < m; j++)
		{
			val += A[i * (m+1) + j] * x[j];
		}
		if (abs(val - A[i*(m+1)+m]) > EPS)
			bRes = false;

		printf("i:%4d %15.10lf %15.10lf %d\n", i, val, A[i*(m+1)+m], bRes);
	} 
	return bRes;
}


bool solve(
	double 		*A,
	int32_t 	n,
	int32_t 	m,
	double 		*x,
	bool 		*hasSolution)
{
	bool bRes = 0;
	
	int32_t i, j;
	
	*hasSolution = true;

	i = j = 0; 

	while (i<n && j < m && (*hasSolution) == true)
	{
		int32_t x = i;

		while (x < n && abs(A[x * (m+1) + j]) <= EPS)
			x++;

		if (x==n)
		{
			/* var j is a free variable */
			j++;
		}	
		else
		{
			int32_t l, u;
			/* line i <-> line x */
			if (x != i)
			{
				for(l=0; l<m+1; l++)
				{
					double aux;
					aux = A[i * (m+1) + l];
					A[i * (m+1) + l] = A[x * (m+1) + l];
					A[x * (m+1) + l] = aux;
				}
			}

			for (l = j+1; l < (m+1); l++)	// elements from 0 to j-1 - already 0
			{
				A[i * (m+1) + l] /= A[i * (m+1) + j];
			}
			A[i * (m+1) + j] = 1;
			
			for (u=i+1; u < n; u++)   
			{
				double val = A[u * (m+1) + j];
				for (l=j; l < (m+1); l++ )   // elements from 0 to j-1 - already 0
					A[u * (m+1) + l] -= (A[i * (m+1) + l] * val);
			}
			j++; 
			i++;
		}
		/*displaySystem(A, n, m);*/
	}
	
	/*displaySystem(A, n, m);*/

	for (i=n-1; i>=0; i--) 		/* all eq need to be checked in order to see if there is a solution or not */
	{
		j = 0;
		while (j < m+1 && abs(A[i * (m+1) + j]) <= EPS)
			j++;
		if (j==m)
		{
			*hasSolution = false;
		}
		else if (j<m)
		{
			/* compute variable j */
			int32_t l;
			x[j] = A[i * (m+1) + m];
			for(l = j+1; l < m; l++)
			{
				x[j] -= A[i * (m+1) + l] * x[l];
			} 
			
		}
	}
	

	return bRes;
}


int main()
{
	bool bRes = 0;

	/* n - nEquations, m - nVariables */
	int32_t n, m;
	/* equations matrix */	
	double 	*A = NULL;
	/* variables */
	double 	*x = NULL;
	/* specify if the system has a solution */
	bool hasSolution = true;

	bRes = readInputData(&A, &n, &m);
	
	/*displaySystem(A, n, m);*/

	if (bRes == 0)
	{
		x = (double*)calloc(m, sizeof(double));
		if (x==NULL)
			bRes = 1;
	}
	
	if (bRes == 0)
	{
		bRes = solve(A, n, m, x, &hasSolution);
	}

	/*displaySystem(A, n, m);*/
/*
	if (hasSolution == true)
	{
		bool checks = checkSystem(A, x, n, m);
		if (checks == true)
			printf("check\n");
		else 
			printf("fail\n");
	}
	*/

	if (bRes == 0)
	{
		writeOutputData(x, m, hasSolution);		
	}
	
	free(A);
	free(x);

	return (int)bRes;
}