Методы вычислений 24.11.22

#define _USE_MATH_DEFINES
#include <iostream>
#include <iomanip>
#include <vector>
#include <random>
#include <chrono>
#include <cmath>
#include <ctime>
using namespace std;
using namespace chrono;

typedef vector <double> vec;
typedef vector <vec> mat;

mt19937 gen(time(NULL));
uniform_real_distribution<double> interval_double(-10, 10);

const double v = 10;

vec operator * (const mat& A, const vec& x)
{
    int n = x.size();

    vec b(n);

    for (int i = 0; i < n; ++i)
        for (int j = 0; j < n; ++j)
            b[i] += A[i][j] * x[j];

    return b;
}

double det(mat A)
{
    int n = A.size();

    double res = 1;

    // прямой ход
    for (int i = 0; i < n; ++i)
    {
        if (fabs(A[i][i]) < 1e-5)
        {
            int m = i;

            for (int k = i + 1; k < n; ++k)
                if (fabs(A[k][i]) > fabs(A[m][i])) m = k;

            for (int j = i; j < n; ++j)
                swap(A[i][j], A[m][j]);

            res *= (-1);
        }

        res *= A[i][i];

        // нормируем строку
        for (int j = i + 1; j < n; ++j)
            A[i][j] /= A[i][i];

        // нормируем строку
        A[i][i] = 1;

        // вычитаем из нижних строк текущую строку
        // получаем в i-ом столбце нули
        for (int k = i + 1; k < n; ++k)
        {
            for (int j = i + 1; j < n; ++j)
                A[k][j] -= A[i][j] * A[k][i];
        }
    }

    return res;
}

void print(mat& A)
{
    int n = A.size();

    for (int i = 0; i < n; ++i, cout << "\n")
        for (int j = 0; j < n; ++j)
            cout << left << setw(10) << A[i][j];
}

void print(vec& a)
{
    int n = a.size();

    cout << "Вектор: (";

    for (int i = 0; i < n - 1; ++i)
        cout << a[i] << ", ";

    cout << a[n - 1] << ")\n\n";
}

void print(mat& A, vec& b)
{
    int n = b.size();

    cout << "Матрица СЛУ имеет вид:\n";

    for (int i = 0; i < n; ++i, cout << "\n")
    {
        for (int j = 0; j < n; ++j)
            cout << left << setw(8) << A[i][j];

        cout << "| " << b[i];
    }

    cout << "\n";
}

mat task_matrix(int n)
{
    mat res(n, vec(n));

    for (int i = 0; i < n; ++i)
        for (int j = 0; j < n; ++j)
        {
            if (i == j)
            {
                if ((rand() % 5) == 0) res[i][j] = 0;
                else res[i][j] = v + i;
            }
            else
            {
                if ((rand() % 10) == 0) res[i][j] = 0;
                else res[i][j] = (v + i) / 100;
            }
        }

    return res;
}

bool check(vec& x, vec& check_x)
{
    const double EPS = 1e-5;

    bool res = true;

    for (int i = 0; i < x.size() && res; ++i)
        res = fabs(x[i] - check_x[i]) < EPS;

    return res;
}

int main()
{
    srand(time(NULL));
    setlocale(LC_ALL, "Russian");
    cout << fixed << setprecision(5);

    int n = 100;
    mat A = task_matrix(n);

    vec x(n);

    for (int i = 0; i < n; ++i)
        x[i] = rand();

    vec b = A * x;

    // print(A, b);

    vec check_x(n);

    double det_system = det(A);

    for (int j = 0; j < n; ++j)
    {
        mat cur = A;

        for (int i = 0; i < n; ++i)
            cur[i][j] = b[i];

        check_x[j] = det(cur) / det_system;
    }

    // print(x);
    // print(check_x);

    if (check(x, check_x)) cout << "n = " << n << " : CORRECT\n";
    else cout << "n = " << n << " : INCORRECT\n";

    return 0;
}