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

#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;
}

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

    // прямой ход
    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;

            swap(b[i], b[m]);

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

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

        // нормируем строку
        b[i] /= 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];

            b[k] -= b[i] * A[k][i];
        }
    }

    vec x(n);

    // обратный ход
    for (int i = n - 1; i >= 0; --i)
    {
        x[i] = b[i];

        for (int j = i + 1; j < n; ++j)
            x[i] -= A[i][j] * x[j];
    }

    return x;
}

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

    cout << "Матрица:\n";

    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)
        {
            res[i][j] = v + i;

            if (i != j) res[i][j] /= 100.0;
            else
            {
                // int p = rand() % 2;
                // if (p)
                res[i][j] = 0;
            }
        }

    return res;
}

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

    for (int i = 0; i < n; ++i)
        for (int j = 0; j < n; ++j)
        {
            res[i][j] = interval_double(gen);

            if (i == j)
            {
                int p = rand() % 2;
                if (p) res[i][j] = 0;
            }
        }

    return res;
}

double solve(int n, mat& A, vec& b, vec& x)
{
    vec checkx = Gauss(A, b);

    double res = 0.0;

    for (int i = 0; i < n; ++i)
        res = max(res, fabs(x[i] - checkx[i]));

    return res;
}

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

    vector <int> sizes = { 10, 50, 100, 250, 500, 1000, 2000, 3000, 4000, 5000 };

    cout << "+-" << "------" << "+-" << "------------" << "+-" << "------------+\n";
    cout << "| " << "n     " << "| " << "Time        " << "| " << "Accuracy    |\n";
    cout << "+-" << "------" << "+-" << "------------" << "+-" << "------------+\n";

    for (auto n : sizes)
    {
        mat A = task_matrix(n);

        vec x(n);

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

        vec b = A * x;

        system_clock::time_point start = system_clock::now();

        double acc = solve(n, A, b, x);

        system_clock::time_point end = system_clock::now();

        duration <double> delta = end - start;

        cout << left << "| " << setw(6) << n << "| " << setw(12) << delta.count() << "| " << setw(12) << acc << "|\n";
    }

    cout << "+" << "------" << "+" << "------------" << "+" << "------------+\n";

    return 0;
}