LU

#include < stdio.h > 
#include < math.h > 
#define N 4
#define FO "%-10.5g"
void scanfm(float * a, int m, int n, char aa) {
    int i,
    j;
    for (i = 0; i < m; i++) for (j = 0; j < n; j++) {
        printf("%c[%d][%d]=", aa, i + 1, j + 1);
        scanf("%f", a + i * n + j);
    }
}
void printm(float * a, int m, int n, char aa) {
    int i,
    j,
    k;
    printf("%c=\n", aa);
    for (i = 0; i < m; i++) {
        for (j = 0; j < n; j++) {
            k = i * n + j;
            if (fabs(a[k]) <= 1e-6) a[k] = 0.0;
            printf(FO, a[k]);
        }
        printf("\n");
    }
}
int lu(float * a, const int n, float * l, float * u) {
    int i,
    j,
    k;
    int p;
    float r;
    for (i = 0; i < n; i++) {
        for (j = i; j < n; j++) {
            p = i * n + j;
            r = 0.0;
            for (k = 0; k <= i - 1; k++) r += l * u[k * n + j];
            u[p] = a[p] - r;
        }
        for (j = i + 1; j < n; j++) {
            p = j * n + i;
            r = 0.0;
            for (k = 0; k <= i - 1; k++) r += l[j * n + k] * u[k * n + i];
            if (fabs(u) <= 1e-5) return 0;
            l[p] = (a[p] - r) / u;
        }
    }
    return 1;
}
void main(void) {
    static float A[N][N],
    L[N][N],
    U[N][N];
    int i,
    j;
    for (i = 0; i < N; i++) L = 1.0;
    scanfm(A[0], N, N, 'A');
    printm(A[0], N, N, 'A');
    i = lu(A[0], N, L[0], U[0]);
    if (i == 0) printf("A can't be resoled!");
    else {
        printm(L[0], N, N, 'L');
        printm(U[0], N, N, 'U');
    }
}