#include <stdio.h>
#include <math.h> 

#define N 4

int main(void)
{
	double A[N][N]={{2,1,4,1},{3,4,-1,-1},{1,-4,1,5},{2,-2,1,3}}, b[N]={-4,3,9,7}, l[N][N]={0}, y[N]={0}, x[N]={0};
	int i,j,k;

	printf("Antes do escalonamento: \n");
	for (i=0;i<N;i++)
	{
		for (j=0;j<N;j++)
		{
			printf("%8.3lf",A[i][j]);

			if (j==(N-1))
				printf("\n");
			else 
				printf(" ");
		}
	}

	for (i=0;i<N;i++)
	{

		l[i][i] = 1.0;

		for (j=i+1;j<N;j++)
		{
			l[j][i]=A[j][i]/A[i][i];

			for (k=i;k<N;k++)
			{
				A[j][k] = A[j][k] - l[j][i]*A[i][k];
			}
		}
	}

	printf("U= \n");
	for (i=0;i<N;i++)
	{
		for (j=0;j<N;j++)
		{
			printf("%8.3lf",A[i][j]);

			if (j==(N-1))
				printf("\n");
			else 
				printf(" ");
		}
	}

	printf("L= \n");
	for (i=0;i<N;i++)
	{
		for (j=0;j<N;j++)
		{
			printf("%8.3lf",l[i][j]);

			if (j==(N-1))
				printf("\n");
			else 
				printf(" ");
		}
	}

	for (i=0;i<N;i++) // SUBSTITUIÇÃO PROGRESSIVA
	{
		double soma = b[i];

		for (j=0;j<i;j++)
		{
			soma-=l[i][j]*y[j];
		}

		y[i]=soma;
	}

	for (i=N-1;i>=0;i--) // SUBSTITUIÇÃO REGRESSIVA
	{
		double soma = y[i];

		for (j=i+1;j<N;j++)
		{
			soma-=A[i][j]*x[j];
		}

		x[i]=soma/A[i][i];
	}


	printf("Solução: \n");
	
	for (i=0;i<N;i++)
	{
		printf("%8.3lf",x[i]);

		if (i==(N-1))
			printf("\n");
		else 
			printf(" ");	
	}

	return 0;
}