//MÉTODO DA FATORAÇÃO LU EM C++

//ELIEZER BATISTA DA COSTA NETO - 359606 - ENGENHARIA MECÂNICA

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

int main(){

// DECLARAÇÃO DE VARIÁVEIS

    int i,k,w,j;
    double X[6], Y[6], R[6], soma;
    double A[6][6], L[6][6], U[6][6], V[6][6];

    R[1]=120; R[2]=0; R[3]=350; R[4]=0; R[5]=0;

    A[1][1]=9; A[1][2]=0; A[1][3]=-3; A[1][4]=0; A[1][5]=0;
    A[2][1]=1; A[2][2]=-1; A[2][3]=0; A[2][4]=0; A[2][5]=0;
    A[3][1]=0; A[3][2]=-2; A[3][3]=9; A[3][4]=0; A[3][5]=0;
    A[4][1]=0; A[4][2]=1; A[4][3]=6; A[4][4]=-9; A[4][5]=2;
    A[5][1]=5; A[5][2]=1; A[5][3]=0; A[5][4]=0; A[5][5]=-6;

// CÁLCULO DA MATRIZ L e U

    for(i=1;i<=5;i++){
        for(j=1;j<=5;j++){
            if(i!=j){
                L[i][j] = 0;
            }
            if(i==j){
                L[i][j] = 1;
            }
        }
    }

    k=1;
    w=1;

    while(k<=4){
        for(i=1; i<=4; i++){
            V[k+i][k] =(A[k+i][w]/A[k][w]);

            L[k+i][k] = -(V[k+i][k]);

            if(fabs(L[k+i][k])<0.001){
                    L[k+i][k] = 0;
            }

            for(j=1; j<=5; j++){
                    A[k+i][j]= A[k+i][j] - (V[k+i][k])*(A[k][j]);

                    if(fabs(A[k+i][j])<0.001){
                        A[k+i][j] = 0;
                    }
            }
        }
        k = k + 1;
        w = w + 1;
    }

    for(i=1;i<=5;i++){
        for(j=1;j<=5;j++){
            U[i][j] = A[i][j];
        }
    }

    for(i=1;i<=5;i++){
        for(j=1;j<=5;j++){
            printf("%f  ", U[i][j]);
        }
        printf("\n\n");
    }

    printf("\n\n\n\n");

    for(i=1;i<=5;i++){
        for(j=1;j<=5;j++){
            printf("%f  ", L[i][j]);
        }
        printf("\n\n");
    }

    printf("\n\n\n\n");


//CÁLCULO DAS INCÓGNITAS: SABE-SE QUE LUX = R E QUE UX = Y, ASSIM LY=R

    for(i=1;i<=5;i++){
        soma = 0;
        for(j=1; j<=(i-1); j++){
            soma = soma + L[i][j]*Y[j];
        }
        Y[i] = R[i] - soma;
    }

    printf("%f  %f  %f  %f  %f \n\n", Y[1],Y[2],Y[3],Y[4],Y[5]);

    for (i=5; i>=1; i--){
        soma = 0;
        for(j=(i+1); j<=5; j++){
            soma = soma + U[i][j]*X[j];
        }

        X[i] = (Y[i] - soma)/U[i][i];
    }

    printf("%f  %f  %f  %f  %f \n\n", X[1],X[2],X[3],X[4],X[5]);

    return(0);
}