C++拟合多项式

#include <iostream>
#include <vector>
#include<cmath>
#include <ctime>//eigen核心部分
#include <Eigen/Core>//稠密矩阵的代数运算（逆、特征值等）
#include <Eigen/Dense>
using namespace Eigen;
using namespace std;
/// <summary>
/// 拟合多项式
/// </summary>
/// <param name="degree">多项式的次</param>
/// <param name="x">x方向的数据点</param>
/// <param name="y">y方向的数据点</param>
/// <param name="a">多项式的系数</param>
/// <param name="R2">决定系数</param>
void Polynomial(int degree, std::vector<double> x, std::vector<double> y, std::vector<double> *a,double* R2) {
    // 生成矩阵A
    MatrixXd A(degree + 1, degree + 1);
    A.setZero();
    for (int j = 0; j < degree + 1; j++) {
        for (int i = 0; i < degree + 1; i++) {
            for (int k = 0; k < (int)x.size(); k++) {
                A(i, j) += 2 * pow(x[k], degree - j) * pow(x[k], degree - i);
            }
        }
    }
    // 生成向量B
    VectorXd B(degree + 1);
    B.setZero();
    for (int j = 0; j < degree + 1; j++) {
        for (int k = 0; k < (int)y.size(); k++) {
            B(j) += 2 * pow(x[k], degree - j) * y[k];
        }
    }

    //多项式方程组系数向量 = 矩阵A的逆 * 向量B
    VectorXd X = A.inverse() * B;
    for (int i = 0; i < X.size(); i++) {
        a->push_back(X(i));//把多项式的系数推到向量指针里
    }

    //计算决定系数R2
    double sum_y = 0.0;
    for (int i = 0; i < (int)y.size(); i++) {
        sum_y += y[i];
    }
    double mean_y = sum_y / y.size();//计算平均值
    double ssr = 0, sst = 0;
    std::vector<double> fit_y;//拟合y
    for (int i = 0; i < (int)y.size(); i++) {
        double value = 0;
        for (int n = degree; n >= 0; n--) {
            value += X(degree - n) * pow(x[i], n);//逐个计算拟合后的y
        }
        fit_y.push_back(value);//把拟合后的y推到向量里
    }
    for (int i = 0; i < (int)y.size(); i++) {
        ssr += pow(fit_y[i] - mean_y, 2);
        sst += pow(y[i] - mean_y, 2);
    }
    *R2 = ssr / sst;
}
int main(int argc, char** argv) {
    // 定义多项式的次
    int degree = 3;

    // 定义数据点
    std::vector<double> x = { 1, 2, 3 ,4, 5 };
    std::vector<double> y = { 1.1, 3.2, 6.7, 18.3, 100.3 };
    std::vector<double>* a = new std::vector<double>();;
    double R2 = 0;
    double* ptr = &R2;
    Polynomial(degree,x,y,a, ptr);
    //打印多项式系数 a0 a1... an
    int i = 0;
    for (double elem : *a) {
        cout << "a" << i++ << " = " << elem << endl;
    }
    cout << "R2 = "<< R2 << endl;
}