#include <iostream>
#include <fstream>
#include <vector>
#include <cmath>
#include <limits>
#include <iomanip>
#include <string>
#include <sstream>
#include <cstdlib>

std::vector<double> polynomialFit(const std::vector<double>& x,
                                  const std::vector<double>& y,
                                  int degree) 
{
    int N = (int)x.size();
    int n = degree;

    std::vector<double> coeffs(n + 1, 0.0);

    // Формируем матрицу (n+1) x (n+2) (последний столбец - правая часть)
    std::vector<std::vector<double>> A(n + 1, std::vector<double>(n + 2, 0.0));

    // Заполнение матрицы нормальных уравнений
    for (int i = 0; i <= n; ++i) {
        for (int j = 0; j <= n; ++j) {
            double sum = 0.0;
            for (int k = 0; k < N; ++k) {
                sum += std::pow(x[k], i + j);
            }
            A[i][j] = sum;
        }
        double sum = 0.0;
        for (int k = 0; k < N; ++k) {
            sum += y[k] * std::pow(x[k], i);
        }
        A[i][n + 1] = sum;
    }

    // Прямой ход метода Гаусса
    for (int i = 0; i <= n; ++i) {
        double maxElement = std::fabs(A[i][i]);
        int maxRow = i;
        for (int k = i + 1; k <= n; ++k) {
            if (std::fabs(A[k][i]) > maxElement) {
                maxElement = std::fabs(A[k][i]);
                maxRow = k;
            }
        }
        for (int k = i; k <= n + 1; ++k) {
            std::swap(A[maxRow][k], A[i][k]);
        }
        for (int k = i + 1; k <= n; ++k) {
            double factor = -A[k][i] / A[i][i];
            for (int j = i; j <= n + 1; ++j) {
                if (i == j) {
                    A[k][j] = 0;
                } else {
                    A[k][j] += factor * A[i][j];
                }
            }
        }
    }

    // Обратный ход
    for (int i = n; i >= 0; --i) {
        coeffs[i] = A[i][n + 1] / A[i][i];
        for (int k = i - 1; k >= 0; --k) {
            A[k][n + 1] -= A[k][i] * coeffs[i];
        }
    }

    return coeffs;
}

int main(int argc, char* argv[]) {
    if (argc < 2) {
        std::cerr << "Ошибка: не указана степень полинома.\n"
                  << "Использование: " << argv[0] << " degree\n";
        return 1;
    }

    int degree = std::atoi(argv[1]);
    if (degree < 0) {
        std::cerr << "Ошибка: степень полинома должна быть >= 0.\n";
        return 1;
    }

    std::ifstream infile("data.txt");
    if (!infile.is_open()) {
        std::cerr << "Не удалось открыть файл data.txt.\n";
        return 1;
    }

    std::vector<double> x;
    std::vector<double> y;
    std::string line;
    while (std::getline(infile, line)) {
        std::stringstream ss(line);
        double xi, yi;
        if (ss >> xi >> yi) {
            x.push_back(xi);
            y.push_back(yi);
        }
    }
    infile.close();

    if (x.size() != y.size() || x.empty()) {
        std::cerr << "Ошибка: некорректные входные данные." << std::endl;
        return 1;
    }

    std::vector<double> coeffs = polynomialFit(x, y, degree);

    std::ofstream outfile("out.txt");
    if (!outfile.is_open()) {
        std::cerr << "Не удалось открыть файл out.txt.\n";
        return 1;
    }

    outfile << std::fixed << std::setprecision(25);

    for (size_t i = 0; i < coeffs.size(); ++i) {
        outfile << "a[" << i << "] = " << coeffs[i] << "\n";
    }

    outfile.close();

    std::cout << "Результаты сохранены в out.txt\n";
    std::cout << "DEGREE " << degree << ", количество коэффициентов: " 
              << coeffs.size() << std::endl;

    return 0;
}
