Kek

#define _CRT_SECURE_NO_WARNINGS

#include <iostream>
#include <vector>
#include <cmath>
#include <algorithm>

using namespace std;

typedef long double ld;
const ld eps = 1e-9;
struct Matrix {
    int n;
    vector<vector<ld> > A;
    Matrix(int n) :n(n) {
        A.assign(n, vector<ld>(n, 0));
    }
    Matrix(const vector<vector<ld> >& a) {
        n = a.size();
        A = a;
    }
    Matrix():n(0) {};
};
pair<Matrix, Matrix> LU_fract(const Matrix& A) {
    int n = A.n;
    vector<vector<ld> > L, U;
    L.assign(n, vector<ld>(n, 0));
    U.assign(n, vector<ld>(n, 0));
    for (int i = 0; i < n; ++i) {
        L[i][i] = 1;
    }
    for (int i = 0; i < n; i++)
    {
        U[0][i] = A.A[0][i];
        L[i][0] = A.A[i][0] / U[0][0];
    }
    for (int i = 1; i < n; i++)
    {
        for (int j = 1; j < n; j++)
        {
            ld sum = 0;
            if (i > j)
            {
                for (int k = 0; k<i; ++k) {
                    sum += (L[i][k] * U[k][j]);
                }
                L[i][j] = (A.A[i][j] - sum) / U[j][j];
            }
            else {
                for (int k = 0; k<j; ++k) {
                    sum += (L[i][k] * U[k][j]);
                }
                U[i][j] = A.A[i][j] - sum;
            }
        }
    }
    return{ Matrix(L),Matrix(U) };
}
vector<ld> Kramer(const Matrix& Start, const vector<ld>& v) {
    ld main_opr = 1;
    pair<Matrix, Matrix> LU = LU_fract(Start);
    int n = Start.n;
    for (int i = 0; i < n; ++i) {
        main_opr *= LU.second.A[i][i];
    }
    vector<ld> ans(n, 1);
    for (int j = 0; j < n; ++j) {
        Matrix A = Start;
        for (int i = 0; i < n; ++i) {
            A.A[i][j] = v[i];
        }
        pair<Matrix, Matrix> next = LU_fract(A);
        for (int i = 0; i < n; ++i) {
            ans[j] *= next.second.A[i][i];
        }
    }
    for (int i = 0; i < n; ++i) {
        ans[i] /= main_opr;
    }
    return ans;
}
vector<ld> Relaxation(const Matrix& A, const vector<ld>& b, const ld& par) {
    vector<ld> cur = b;
    int n = A.n;
    ld checker = eps;
    while (true) {
        vector<ld> old = cur;
        for (int i = 0;i < n; ++i) {
            ld sum = 0;
            for (int j = 0; j < n; ++j) {
                if (j == i) continue;
                sum -= (A.A[i][j] * cur[j]);
            }
            sum += b[i];
            sum *= (par/A.A[i][i]);
            cur[i] = sum+(1-par)*cur[i];
        }
        ld mx = abs(cur[0] - old[0]);
        for (int i = 0; i < n; ++i) {
            mx = max(abs(cur[i] - old[i]), mx);
        }
        if (mx < checker) break;
    }
    return cur;
}
vector<ld> Spinning(const Matrix& A, const vector<ld>& b) {
    Matrix cur = A;
    int n = A.n;
    vector<ld> v = b;
    for (int j = 0; j < n - 1; ++j) {
        for (int i = j+1; i < n; ++i) {
            ld denum = sqrt(cur.A[j][j] * cur.A[j][j] + cur.A[i][j] * cur.A[i][j]);
            ld c = cur.A[j][j]/denum;
            ld s = cur.A[i][j] / denum;
            vector<ld> new_Aj(n,0);
            vector<ld> new_Ai(n,0);
            ld new_bj = c*v[j] + s*v[i];
            ld new_bi = -s*v[j] + c*v[i];
            for (int z = 0; z < n; ++z) {
                new_Aj[z] = c*cur.A[j][z]+s*cur.A[i][z];
                new_Ai[z] = -s*cur.A[j][z] + c*cur.A[i][z];
            }
            cur.A[i] = new_Ai;
            cur.A[j] = new_Aj;
            v[j] = new_bj;
            v[i] = new_bi;
        }
    }
    vector<ld> ans(n, 0);
    for (int i = n - 1; i >= 0; --i) {
        ld sum = v[i];
        for (int k = i + 1; k < n; ++k) {
            sum -= (cur.A[i][k] * ans[k]);
        }
        sum /= cur.A[i][i];
        ans[i] = sum;
    }
    return ans;
}
vector<ld> Precised_Spinning(const Matrix& A, const vector<ld>& b) {
    vector<ld> answ_x0 = Spinning(A, b);
    int n = A.n;
    while (true) {
        vector<ld> e(n);
        for (int i = 0; i < n; ++i) {
            ld sum = 0;
            for (int j = 0; j < n; ++j) {
                sum += (A.A[i][j] * answ_x0[j]);
            }
            e[i] = b[i] - sum;
        }
        ld prec = e[0];
        for (int i = 0; i < n; ++i) {
            prec = max(prec, abs(e[i]));
        }
        if (prec < eps) return answ_x0;
        vector<ld> p = Spinning(A, e);
        for (int i = 0; i < n; ++i) {
            answ_x0[i] += p[i];
        }
    }
}
vector<ld> LU_eigenvalues(const Matrix& Start) {
    Matrix A = Start;
    int n = Start.n;
    for (int t = 0; t < 1000; ++t) {
        pair<Matrix, Matrix> LU = LU_fract(A);
        Matrix next(n);
        for (int i = 0; i < n; ++i) {
            for (int j = 0; j < n; ++j) {
                ld sum = 0;
                for (int z = 0; z < n; ++z) {
                    sum += (LU.second.A[i][z] * LU.first.A[z][j]);
                }
                next.A[i][j] = sum;
            }
        }
        A = next;
    }
    vector<ld> ans(n);
    for (int i = 0; i < n; ++i) {
        ans[i] = A.A[i][i];
    }
    return ans;
}
vector<pair<ld, vector<ld> > > LU_pair_eigen(const Matrix& Start) {
    vector<ld> Eigen = LU_eigenvalues(Start);
    int n = Start.n;
    vector<pair<ld, vector<ld> > > ans(n);
    for (int i = 0; i < n; ++i) {
        Matrix A = Start;
        vector<ld> empty_b(n, 0);
        for (int j = 0; j < n; ++j) {
            A.A[j][j] -= Eigen[i];
        }
        vector<ld> now_e = Kramer(A, empty_b);
        ans[i] = { Eigen[i],now_e };
    }
    return ans;
}
vector<ld> Variation(const Matrix& A, const vector<ld>& b) {
    int n = A.n;
    vector<ld> x = b;
    vector<ld> e(n);
    for (int i = 0; i < n; ++i) {
        ld sum = 0;
        for (int j = 0; j < n; ++j) {
            sum += (A.A[i][j] * x[j]);
        }
        e[i] = b[i] - sum;
    }
    vector<ld> p = e;
    while (true) {
        ld mx = e[0];
        for (int i = 0; i < n; ++i) {
            mx = max(mx, abs(e[i]));
        }
        if (mx < eps) break;
        vector<ld> q(n);
        for (int i = 0; i < n; ++i) {
            ld sum = 0;
            for (int j = 0; j < n; ++j) {
                sum += (A.A[i][j] * p[j]);
            }
            q[i] = sum;
        }
        ld alpha = 0;
        ld num = 0, denum = 0;
        for (int i = 0;i<n;++i){
            num += (e[i] * p[i]);
            denum += (q[i] * p[i]);
        }
        alpha = num / denum;
        vector<ld> new_x(n);
        for (int i = 0; i < n; ++i) {
            new_x[i] = x[i] + alpha*p[i];
        }
        vector<ld> new_e(n);
        for (int i = 0; i < n; ++i) {
            new_e[i] = e[i] - alpha*q[i];
        }
        ld beta = 0;
        num = 0, denum = 0;
        for (int i = 0; i < n; ++i) {
            num += (new_e[i] * q[i]);
            denum += (p[i] * q[i]);
        }
        beta = num / denum;
        vector<ld> new_p(n);
        for (int i = 0; i < n; ++i) {
            new_p[i] = new_e[i] - beta*p[i];
        }
        p = new_p;
        e = new_e;
        x = new_x;
    }
    return x;
}

int main()
{
    freopen("input.txt", "r", stdin);
    freopen("output.txt", "w", stdout);
    int n;
    cin >> n;
    Matrix A(n);
    for (int i = 0; i < n; ++i) {
        for (int j = 0; j < n; ++j) {
            cin >> A.A[i][j];
        }
    }
    vector<ld> v(n);
    for (int i = 0; i < n; ++i) {
        cin >> v[i];
    }
   // vector<ld> Relax_ans = Relaxation(A, v, 1.5).first;
   // vector<ld> Kramer_ans = Kramer(A, v);
    vector<ld> Spinning_ans = Spinning(A, v);
   // vector<ld> Variation_ans = Variation(A, v);
    cout.precision(10);
    //cout << "Relaxation method answer is : \n";
    //for (int i = 0; i < n; ++i) {
    //    cout << fixed << Relax_ans[i] << '\n';
    //}
    //cout << "\nKramer method answer is :\n";
    //for (int i = 0; i < n; ++i) {
    //    cout << fixed << Kramer_ans[i] << '\n';
    //}
    cout << "\nPrecised Spinning method answer is : \n";
    for (int i = 0; i < n; ++i) {
        cout << fixed << Spinning_ans[i] << '\n';
    }
    //cout << "\nVariation method answer is : \n";
    //for (int i = 0; i < n; ++i) {
    //    cout << fixed << Variation_ans[i] << '\n';
    //}
    int n1;
    cin >> n1;
    Matrix Test_Eigenvalues(n1);
    for (int i = 0; i < n1; ++i) {
        for (int j = 0; j < n1; ++j) {
            cin >> Test_Eigenvalues.A[i][j];
        }
    }
    vector<ld> eigen = LU_eigenvalues(Test_Eigenvalues);
    cout << "\nEigenvalues : \n";
    for (int i = 0; i < n1; ++i) {
        cout << fixed << eigen[i];
        cout << '\n';
    }
}