C++拟合多项式过零，拟合优度与相关系数求解

void fitPolynomialThroughOrigin(int n, double *x, double *y, int degree, double *coe, double *r_square, double *r)
{
    std::vector<double> coefficients(degree + 1, 0.0);
    std::vector<std::vector<double> > A(degree + 1, std::vector<double>(degree + 1, 0.0));
    std::vector<double> B(degree + 1, 0.0);

    for (int i = 0; i < n; i++)
    {
        for (int j = 0; j <= degree; j++)
        {
            for (int k = 0; k <= degree; k++)
            {
                A[j][k] += pow(x[i], j + k);
            }
            B[j] += y[i] * pow(x[i], j);
        }
    }

    for (int i = 0; i <= degree; i++)
    {
        A[0][i] = 0.0;
    }
    B[0] = 0.0;

    for (int i = 1; i <= degree; i++)
    {
        for (int j = i + 1; j <= degree; j++)
        {
            double factor = A[j][i] / A[i][i];
            for (int k = i; k <= degree; k++)
            {
                A[j][k] -= factor * A[i][k];
            }
            B[j] -= factor * B[i];
        }
    }

    for (int i = degree; i > 0; i--)
    {
        coefficients[i] = B[i] / A[i][i];
    }

    double x_mean = 0.0;
    double y_mean = 0.0;
    for (int i = 0; i < n; i++)
    {
        x_mean += x[i];
        y_mean += y[i];
    }
    x_mean /= n;
    y_mean /= n;

    double res = 0.0;
    double sum_ = 0.0;
    double total = 0.0;
    for (int j = 0; j < n; j++)
    {
        for (int i = 0; i < coefficients.size(); i++)
        {
            sum_ += coefficients[i] * pow(x[j], i);
        }
        res += pow(sum_ - y[j], 2);
        total += pow(y[j] - y_mean, 2);
        sum_ = 0.0;
    }

    double r_squared = 1.0 - res / total;

    double up = 0.0;
    double down1 = 0.0;
    double down2 = 0.0;
    for (int i = 0; i < n; i++)
    {
        up += (x[i] - x_mean) * (y[i] - y_mean);
        down1 += pow(x[i] - x_mean, 2);
        down2 += pow(y[i] - y_mean, 2);
    }

    for (int i = 0; i < degree + 1; i++)
    {
        coe[i] = coefficients[i];
    }
    *r_square = r_squared;
    *r = up / (sqrt(down1) * sqrt(down2));
}

// 多次函数拟合函数
std::vector<double> poly::fitPolynomial(std::vector<double> x, std::vector<double> y, int degree) {
    int n = x.size();

    // 创建一个系数向量，初始化为0
    std::vector<double> coefficients(degree + 1, 0.0);

    // 创建矩阵 A 和向量 B，用于最小二乘法
    std::vector<std::vector<double>> A(degree + 1, std::vector<double>(degree + 1, 0.0));
    std::vector<double> B(degree + 1, 0.0);

    // 填充 A 和 B
    for (int i = 0; i < n; i++) {
        for (int j = 0; j <= degree; j++) {
            for (int k = 0; k <= degree; k++) {
                A[j][k] += std::pow(x[i], j + k);
            }
            B[j] += y[i] * std::pow(x[i], j);
        }
    }

    // 使用高斯消元法求解线性方程组
    for (int i = 0; i <= degree; i++) {
        for (int j = i + 1; j <= degree; j++) {
            double factor = A[j][i] / A[i][i];
            for (int k = i; k <= degree; k++) {
                A[j][k] -= factor * A[i][k];
            }
            B[j] -= factor * B[i];
        }
    }

    // 回代求解系数
    for (int i = degree; i >= 0; i--) {
        coefficients[i] = B[i];
        for (int j = i + 1; j <= degree; j++) {
            coefficients[i] -= A[i][j] * coefficients[j];
        }
        coefficients[i] /= A[i][i];
    }

    return coefficients;
}

// 多次函数拟合函数--过零
std::vector<double> poly::fitPolynomial_through_origin(std::vector<double> x, std::vector<double>y, int degree) {
    int n = x.size();

    // 创建一个系数向量，初始化为0
    std::vector<double> coefficients(degree + 1, 0.0);

    // 创建矩阵 A 和向量 B，用于最小二乘法
    std::vector<std::vector<double>> A(degree + 1, std::vector<double>(degree + 1, 0.0));
    std::vector<double> B(degree + 1, 0.0);

    // 填充 A 和 B
    for (int i = 0; i < n; i++) {
        for (int j = 0; j <= degree; j++) {
            for (int k = 0; k <= degree; k++) {
                A[j][k] += std::pow(x[i], j + k);
            }
            B[j] += y[i] * std::pow(x[i], j);
            cout<<"111B[j]"<<B[j]<<endl;
        }
    }
    // 修改 A 和 B确保零截距
    for (int i = 0; i < degree+1; i++) {
        A[0][i] = 0;
    }
    B[0] = 0;

    //    // 使用高斯消元法求解线性方程组
    //    for (int i = 0; i <= degree; i++) {
    //        for (int j = i + 1; j <= degree; j++) {
    //            double factor = A[j][i] / A[i][i];
    //            for (int k = i; k <= degree; k++) {
    //                A[j][k] -= factor * A[i][k];
    //            }
    //            B[j] -= factor * B[i];
    //            cout<<"222B[j]"<<B[j]<<endl;
    //        }
    //    }

    //    // 回代求解系数
    //    for (int i = degree; i >= 0; i--) {
    //        coefficients[i] = B[i];
    //        for (int j = i + 1; j <= degree; j++) {
    //            coefficients[i] -= A[i][j] * coefficients[j];
    //        }
    //        coefficients[i] /= A[i][i];
    //        cout<<" coefficients[i]"<< coefficients[i]<<endl;
    //    }
    // 使用高斯消元法求解线性方程组
    for (int i = 1; i <= degree; i++) {
        for (int j = i + 1; j <= degree; j++) {
            double factor = A[j][i] / A[i][i];
            for (int k = i; k <= degree; k++) {
                A[j][k] -= factor * A[i][k];
            }
            B[j] -= factor * B[i];
            cout << "222B[j]" << B[j] << endl;
        }
    }

    // 回代求解系数
    for (int i = degree; i >= 1; i--) {
        coefficients[i] = B[i];
        for (int j = i + 1; j <= degree; j++) {
            coefficients[i] -= A[i][j] * coefficients[j];
        }
        coefficients[i] /= A[i][i];
        cout << " coefficients[i]" << coefficients[i] << endl;
    }

    return coefficients;
}