using namespace std;

#include <iostream>
#include <cstdio>
#include <vector>
#include <cmath>
#include <iomanip>

void print(double matrix[4][5]) {
    for (int i = 0; i < 4; i++) {
        for (int j = 0; j < 5; j++) {
           cout << setw(5) << matrix[i][j];
        }
        cout << endl;
    }
}

void CompletePivoting(double a[4][5]) {
    cout << "Complete pivoting method" << endl;
    print(a);

    int cl[4];
    vector<vector<double>> t(4), matrix(4);
    for (int i = 0; i < 4; i++) {
        matrix[i].resize(5);
        for (int j = 0; j < 5; j++) {
            matrix[i][j] = a[i][j];
        }
    }

    for (int k = 0; k < 3; k++) {
        int p = 0, q = 0;

        double a = matrix[0][0];
        for (int i = 0; i < 4 - k; i++) {
            for (int j = 0; j < 4 - k; j++) {
                if (fabs(matrix[i][j]) > fabs(a)) {
                    a = matrix[i][j];
                    p = i;
                    q = j;
                }
            }
        }

        for (int i = 0; i < 4 - k; i++) {
            if (i != p) {
                double m = -matrix[i][q] / a;
                for (int j = 0; j < 5 - k; j++) {
                    matrix[i][j] += matrix[p][j] * m;
                }
            }
        }

        cl[3 - k] = q;
        t[3 - k] = matrix[p];

        matrix.erase(matrix.begin() + p);

        for (int i = 0; i < 3 - k; i++) {
            matrix[i].erase(matrix[i].begin() + q);
        }
    }

    cl[0] = 0;
    t[0] = matrix[0];

    double x[4];
    for (int i = 0; i < 4; i++) {
        double s = 0;
        for (int j = 0, index = 0; j < t[i].size() - 1; j++) {
            if (j == cl[i]) continue;
            s += x[index++] * t[i][j];
        }
        x[i] = (t[i][t[i].size() - 1] - s) / t[i][cl[i]];
    }

    int counter[4] {0};
    double X[4] {0};
    for (int i = 0; i < 4; i++) {
        X[(counter[cl[i]] += cl[i] + 1) - 1] = x[i];
    }

    for (int i = 0; i < 4; i++) {
        cout << "x" << (i + 1) << " = " << X[i] << endl;
    }


}

void SimpleIteration(double matrix[4][5]) {
    cout << endl << "Simple iteration method" << endl;
    print(matrix);

    double q = 0, result[4], x0[4];

    for (int i = 0; i < 4; i++) {
        double a = matrix[i][i];
        matrix[i][4] /= a;
        for (int j = 0; j < 4; j++) {
            if (j != i) {
                matrix[i][j] /= a * -1;
            } else {
                matrix[i][j] = 0;
            }
        }
    }
    for (int i = 0; i < 4; i++) {
        double a = 0;
        for (int j = 0; j < 4; j++) {
            a += abs(matrix[i][j]);
        }
        if (q < a) {
            q = a;
        }
    }
    for (int i = 0; i < 4; i++) {
        x0[i] = matrix[i][4];
    }

    printf("|   x1   |   x2   |   x3   |   x4   |\n");

    while (true) {
        printf("|");
        for (int i = 0; i < 4; i++) {
            result[i] = matrix[i][4];
            for (int j = 0; j < 4; j++) {
                result[i] += matrix[i][j] * x0[j];
            }
            printf("%7.3f |", result[i]);
        }
        printf("\n");
        double a = abs(result[0] - x0[0]);
        for (int i = 1; i < 4; i++) {
            double t = abs(result[i] - x0[i]);
            if (t > a) a = t;
        }

        if (a <= abs((1 - q) * 1e-2 / q)) break;

        for (int i = 0; i < 4; i++) {
            x0[i] = result[i];
        }
    }
}

int main() {
    double AB[4][5] = {
        {  3, 19, 11,  8, 149 },
        {  9, 31,  3, 18, 257 },
        { 11,  7, 32, 13, 143 },
        { 12, 19, 12,  5, 144 }
    };

    CompletePivoting(AB);

    double matrix[4][5] = {
        {  9,   0,  1, -3,  -5 }, //[3] - [0]
        {  3,  19, 11,  8, 149 }, //[0]
        { -1, -12, 20,  8,  -1 }, //[3] * -1 + [2]
        {  5,   0, 12, 18, 107 }  //[1] + [2] - [3] - [0]
    };

    SimpleIteration(matrix);
}