AGLA 2

#include <iostream>
#include <iomanip>
#include <math.h>

using namespace std;

class Matrix
{
private:
    int size;
    double** matrix, ** inverse;

public:
    Matrix(int size) {
        (*this).size = size;
        matrix = new double* [size];
        inverse = new double* [size];
        for (int i = 0; i < size; i++) {
            matrix[i] = new double[size];
            inverse[i] = new double[size];
        }

        for (int i = 0; i < size; i++)
            for (int j = 0; j < size; j++)
                inverse[i][j] = (i == j) ? inverse[i][j] = 1 : inverse[i][j] = 0;
    }

    double getElement(int i, int j) { return matrix[i][j]; }
    int getSize() { return size; }

    void permutation(int r1, int r2) {
        double* temp = *(matrix + r1);
        *(matrix + r1) = *(matrix + r2);
        *(matrix + r2) = temp;

        temp = *(inverse + r1);
        *(inverse + r1) = *(inverse + r2);
        *(inverse + r2) = temp;
    }

    void elimination(int upperLine, int bottomLine) {
        double constant = matrix[bottomLine][upperLine] / matrix[upperLine][upperLine];

        if (matrix[bottomLine][upperLine] != 0) {

            for (int j = 0; j < size; j++) {
                matrix[bottomLine][j] -= matrix[upperLine][j] * constant;
                inverse[bottomLine][j] -= inverse[upperLine][j] * constant;
            }
        }
    }

    void inverseMatrix() {
        for (int i = 0; i < size; i++) {
            for (int j = 0; j < size; j++) {
                cout << fixed << setprecision(2) << inverse[i][j];
                if (j != size - 1)
                    cout << " ";
            }
            cout << endl;
        }
    }

    void normalization() {
        for (int i = 0; i < size; i++) {
            double constant = matrix[i][i];
            for (int j = 0; j < size; j++) {
                matrix[i][j] /= constant;
                inverse[i][j] /= constant;
            }
        }
    }

    friend istream& operator>>(istream& cin, Matrix base) {
        for (int i = 0; i < base.size; i++)
            for (int j = 0; j < base.size; j++)
                cin >> base.matrix[i][j];
        return cin;
    }
    friend ostream& operator<<(ostream& cout, Matrix base) {
        for (int i = 0; i < base.size; i++) {
            for (int j = 0; j < 2 * base.size; j++) {
                if (j < base.size) {
                    if (base.matrix[i][j] < 0.005 && base.matrix[i][j] > -0.005)
                        cout << "0.00";
                    else
                        cout << fixed << setprecision(2) << base.matrix[i][j];
                }
                else {
                    if (base.inverse[i][j - base.size] < 0.005 && base.inverse[i][j - base.size] > -0.005)
                        cout << "0.00";
                    else
                        cout << fixed << setprecision(2) << base.inverse[i][j - base.size];
                }

                if (j != 2 * base.size - 1)
                    cout << " ";
            }
            cout << endl;
        }
        return cout;
    }
};

void output(int& step, const Matrix& base, bool permutation) {
    cout << "step #" << step << ": ";
    if (permutation)
        cout << "permutation";
    else {
        cout << "elimination";
        ++step;
    }
    cout << endl << base;
}

void calculusFirst(Matrix& base, int line, int& step) {
    int maxLine = line;
    for (int i = line + 1; i < base.getSize(); i++) {
        if (abs(base.getElement(i, line)) > abs(base.getElement(maxLine, line)))
            maxLine = i;
    }

    if (maxLine != line) {
        base.permutation(line, maxLine);
        output(step, base, true);
    }

    for (int i = line + 1; i < base.getSize(); i++) {
        base.elimination(line, i);
        output(step, base, false);
    }
}

void calculusSecond(Matrix& base, int line, int& step) {

    for (int i = line - 1; i >= 0; i--) {
        base.elimination(line, i);
        output(step, base, false);
    }
}

int main()
{
    int n, step = 1;
    cin >> n;
    Matrix base(n);
    cin >> base;
    cout << "step #0: Augmented Matrix\n" << base;

    cout << "Direct way:\n";
    for (int i = 0; i < n; i++)
        calculusFirst(base, i, step);

    cout << "Way back:\n";
    for (int i = n - 1; i > 0; i--)
        calculusSecond(base, i, step);

    base.normalization();
    cout << "Diagonal normalization:" << endl << base;
    cout << "result:" << endl;
    base.inverseMatrix();

    return 0;
}