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

#define N 3
double A[N][N]= {{1,2,3},{4,5,7},{-3,2,8}};
double *P[N];

void matprint(double B[N][N]){
    for(int i=0;i<N;i++){
        for(int j=0;j<N;j++){
            printf("%g ",B[i][j]);
        }
        printf("\n");
    }
}
void matmemprint(double B[N][N]){
    for(int i=0;i<N;i++){
        for(int j=0;j<N;j++){
            printf("%p:%g ",&B[i][j],B[i][j]);
        }
        printf("\n");
    }
}
void pointermatprint(double *P[N]){
    for(int i=0;i<N;i++){
        for(int j=0;j<N;j++){
            printf("%g ",P[i][j]);
        }
        printf("\n");
    }
}
void pointerUprint(double *P[N]){
    for(int i=0;i<N;i++){
        for(int j=0;j<N;j++){
            if(j>=i) printf("%g ",P[i][j]);
            else printf("0 ");
        }
        printf("\n");
    }
}
void pointerPprint(double *P[N]){
    double *pmin=P[0];
    for(int i=1;i<N;i++){
        if(pmin>P[i]) pmin=P[i];
    }
    printf("the real first row was at %p\n",pmin);
    for(int i=0;i<N;i++){
        printf("P[%d]=%p row %ld\n",
            i,P[i],(P[i]-pmin)/N);
    }
}
void pointerLprint(double *P[N]){
    for(int i=0;i<N;i++){
        for(int j=0;j<N;j++){
            if(j<i) printf("%g ",P[i][j]);
            else if(j==i) printf("1 ");
            else printf("0 ");
        }
        printf("\n");
    }
}
int main(){
    matprint(A);
    matmemprint(A);
    for(int i=0;i<N;i++) P[i]=&A[i][0];
    pointermatprint(P);
    for(int k=0;k<N-1;k++){
        double pmax=fabs(P[k][k]); int imax=k;
        for(int i=k+1;i<N;i++){
            if(pmax<fabs(P[i][k])){
                pmax=fabs(P[i][k]); imax=i;
            }
        }
        printf("imax=%d\n",imax);
        double *t=P[k]; P[k]=P[imax]; P[imax]=t;
        printf("after swapping\n");
        pointermatprint(P);
        for(int i=k+1;i<N;i++){
            double alpha=P[i][k]/P[k][k];
            for(int j=k;j<N;j++){
                P[i][j]-=alpha*P[k][j];
            }
            P[i][k]=alpha;
        }
        printf("after some elimination\n");
        pointermatprint(P);
    }
    printf("U=\n"); pointerUprint(P);
    printf("L=\n"); pointerLprint(P);
    pointerPprint(P);
    return 0;
}