*/
// Gaussian Elimination
// Using UTF-8 (Unicode) encoding of Russian characters
//
// Приведение матрицы к ступенчатому виду
// методом Гаусса
//
#include 
#include 

const double EPS = 1e-8;

int gauss(double *a, int m, int n);

void swapRows(double *a, int m, int n, int i1, int i2);
void addRows(double *a, int m, int n, 
    int i, int k, double lambda);
void multRow(double *a, int m, int n, int i, double lambda);

double det(const double *a, int m, int n);
bool readMatrix(const char *path, double **a, int *m, int *n);
bool writeMatrix(
    const char *path, const double *a, int m, int n
);
void printMatrix(const double *a, int m, int n);

// Читаем матрицу из файла "input.txt",
// приводим ее к ступенчатому виду и записываем
// результат в файл "output.txt".
// Формат файла:
//     m n
//     элементы матрицы по строкам
// где m - число строк матрицы, n - число столбцов.
// Например, матрица размера 2x3 записывается как
//     2 3
//     a11 a12 a13
//     a21 a22 a23
// Программа также печатает ранг матрицы и для
// квадратной матрицы ее определитель.
int main() {
    int m, n;
    double *a;
    if (!readMatrix("input.txt", &a, &m, &n)) {
        printf("Cannot read a matrix.\n");
        return (-1);
    }
    printMatrix(a, m, n);
    printf("----\n");  
  
    int rank = gauss(a, m, n);
    printMatrix(a, m, n);
    printf("----\n");  

    printf("rank = %d\n", rank);
    if (n == m) {
        double d = det(a, n, m);
        printf("determinant = %.6f\n", d);
    }
    
    if (!writeMatrix("output.txt", a, m, n)) {
        printf("Cannot write a matrix.\n");
        return (-1);
    }

    return 0;
}

int gauss(double *a, int m, int n) {
    int i = 0;
    int j = 0;
    
    while (i < m && j < n) {
        // 1. Ищем макс. ненулевой эл-т в j-м столбце,
        //    начиная со строки i
        double amax = (-1.);
        int kmax = j;
        for (int k = i; k < m; ++k) {
            if (fabs(a[k*n + j]) > amax) {
                amax = fabs(a[k*n + j]);
                kmax = k;
            }
        }
        if (amax <= EPS) {
            // Нулевой столбец
            for (int k = i; k < m; ++k) {
                a[k*n + j] = 0.;
            }
            ++j;
            continue;
        }
        if (kmax != i) {
            // Меняем местами строки i и kmax
            swapRows(a, m, n, i, kmax);
        }
        
        // Зануляем j-й столбец, начиная со строки i+1
        for (int k = i+1; k < m; ++k) {
            // От строки k вычитаем строку i, умноженную
            // на коэффициент lambda
            double lambda = a[k*n + j]/a[i*n + j];
            addRows(a, m, n, k, i, -lambda);
            a[k*n + j] = 0.;
        }
        ++i; ++j;    
    }
    return i;
}

bool writeMatrix(
    const char *path, const double *a, int m, int n
) {
    FILE *f = fopen(path, "wt");
    if (f == NULL) {
        return false;
    }
    fprintf(f, "%d %d\n", m, n);
    for (int i = 0; i < m; ++i) {
        for (int j = 0; j < n; ++j) {
            if (j > 0)
                fprintf(f, " "); // Разделитель
            fprintf(f, "%16.6f", a[i*n + j]);    
        }
        fprintf(f, "\n");
    }
    fclose(f);
    return true;
}

void printMatrix(const double *a, int m, int n) {
    for (int i = 0; i < m; ++i) {
        for (int j = 0; j < n; ++j) {
            if (j > 0)
                printf(" ");
            printf("%16.6f", a[i*n + j]);
        }
        printf("\n");
    }
}

bool readMatrix(
    const char *path, double **a, int *m, int *n
) {
    FILE *f = fopen(path, "rt");
    if (f == NULL) {
        return false;
    }
    
    if (fscanf(f, "%d%d", m, n) < 2) {
        fclose(f);
        return false;
    }
    
    int rows = *m;
    int cols = *n;
    
    double *matr = new double[rows*cols];
    *a = matr;
    
    for (int i = 0; i < rows; ++i) {
        for (int j = 0; j < cols; ++j) {
            if (fscanf(f, "%lf", &(matr[i*cols + j])) < 1) {
                fclose(f);
                return false;
            }    
        }
    }
    fclose(f);
    return true;
}

double det(const double *a, int m, int n) {
    double d = 1.;
    for (int i = 0; i < n && i < m; ++i) {
        d *= a[i*n + i];
    }
    return d;
}


void swapRows(double *a, int m, int n, int i1, int i2) {
    for (int j = 0; j < n; ++ j) {
        double tmp = a[i1*n + j];
        a[i1*n + j] = a[i2*n + j];
        a[i2*n + j] = (-tmp);
    }
}

void addRows(double *a, int m, int n, 
    int i, int k, double lambda) {
    for (int j = 0; j < n; ++ j) {
        a[i*n + j] += a[k*n + j]*lambda;
    }
}

void multRow(double *a, int m, int n, int i, double lambda) {
    for (int j = 0; j < n; ++ j) {
        a[i*n + j] *= lambda;
    }
}
/*