c++最小二乘法拟合多项式

C++实现最小二乘法拟合多项式的步骤如下：

1. 定义自变量x、因变量y和多项式的次数m

double x[n], y[n];
int m;

2. 输入自变量x和因变量y的数据

for (int i = 0; i < n; i++) {
    cin >> x[i] >> y[i];
}

3. 输入多项式的次数m

cin >> m;

4. 定义多项式的系数c

vector<double> c(m + 1, 0.0);

5. 计算多项式的系数c

vector<double> v_x(n), v_y(n);
for (int i = 0; i < n; i++) {
    v_x[i] = 1.0;
    v_y[i] = y[i];
    for (int j = 1; j <= m; j++) {
        v_x[i] *= x[i];
    }
    for (int j = 0; j <= m; j++) {
        c[j] += v_x[i] * v_y[i];
        if (j < m) {
            v_x[i] *= x[i];
        }
    }
}

vector<vector<double>> a(m + 1, vector<double>(m + 1, 0.0));
for (int i = 0; i <= m; i++) {
    for (int j = 0; j <= m; j++) {
        for (int k = 0; k < n; k++) {
            if (i == 0 && j == 0) {
                a[i][j] += 1.0;
            } else {
                a[i][j] += pow(x[k], i + j);
            }
        }
    }
}

vector<double> b(m + 1, 0.0);
for (int i = 0; i <= m; i++) {
    for (int j = 0; j < n; j++) {
        b[i] += pow(x[j], i) * y[j];
    }
}

for (int i = 0; i <= m; i++) {
    for (int j = i + 1; j <= m; j++) {
        double ratio = a[j][i] / a[i][i];
        for (int k = i; k <= m; k++) {
            a[j][k] -= ratio * a[i][k];
        }
        b[j] -= ratio * b[i];
    }
}

for (int i = m; i >= 0; i--) {
    for (int j = i - 1; j >= 0; j--) {
        double ratio = a[j][i] / a[i][i];
        for (int k = i; k >= j; k--) {
            a[j][k] -= ratio * a[i][k];
        }
        b[j] -= ratio * b[i];
    }
}

for (int i = 0; i <= m; i++) {
    c[i] = b[i] / a[i][i];
}

6. 输出多项式的系数c

cout << "The coefficients of the polynomial: ";
for (int i = 0; i <= m; i++) {
    cout << c[i] << " ";
}
cout << endl;

完整代码如下：

#include <iostream>
#include <vector>
#include <cmath>

using namespace std;

int main() {
    int n, m;
    cin >> n;

    double x[n], y[n];
    for (int i = 0; i < n; i++) {
        cin >> x[i] >> y[i];
    }

    cin >> m;

    vector<double> c(m + 1, 0.0);

    vector<double> v_x(n), v_y(n);
    for (int i = 0; i < n; i++) {
        v_x[i] = 1.0;
        v_y[i] = y[i];
        for (int j = 1; j <= m; j++) {
            v_x[i] *= x[i];
        }
        for (int j = 0; j <= m; j++) {
            c[j] += v_x[i] * v_y[i];
            if (j < m) {
                v_x[i] *= x[i];
            }
        }
    }

    vector<vector<double>> a(m + 1, vector<double>(m + 1, 0.0));
    for (int i = 0; i <= m; i++) {
        for (int j = 0; j <= m; j++) {
            for (int k = 0; k < n; k++) {
                if (i == 0 && j == 0) {
                    a[i][j] += 1.0;
                } else {
                    a[i][j] += pow(x[k], i + j);
                }
            }
        }
    }

    vector<double> b(m + 1, 0.0);
    for (int i = 0; i <= m; i++) {
        for (int j = 0; j < n; j++) {
            b[i] += pow(x[j], i) * y[j];
        }
    }

    for (int i = 0; i <= m; i++) {
        for (int j = i + 1; j <= m; j++) {
            double ratio = a[j][i] / a[i][i];
            for (int k = i; k <= m; k++) {
                a[j][k] -= ratio * a[i][k];
            }
            b[j] -= ratio * b[i];
        }
    }

    for (int i = m; i >= 0; i--) {
        for (int j = i - 1; j >= 0; j--) {
            double ratio = a[j][i] / a[i][i];
            for (int k = i; k >= j; k--) {
                a[j][k] -= ratio * a[i][k];
            }
            b[j] -= ratio * b[i];
        }
    }

    for (int i = 0; i <= m; i++) {
        c[i] = b[i] / a[i][i];
    }

    cout << "The coefficients of the polynomial: ";
    for (int i = 0; i <= m; i++) {
        cout << c[i] << " ";
    }
    cout << endl;

    return 0;
}