写一段三次样条插值函数的c++代码

以下是三次样条插值函数的C++代码：

#include <iostream>
#include <vector>

using namespace std;

// 三次样条插值函数
double splineInterpolation(vector<double>& x, vector<double>& y, double xi) {
    // 计算步长
    int n = x.size();
    vector<double> h(n - 1);
    for (int i = 0; i < n - 1; i++) {
        h[i] = x[i + 1] - x[i];
    }

    // 计算三对角矩阵的系数
    vector<double> a(n - 2), b(n - 1), c(n - 2), d(n - 1);
    for (int i = 0; i < n - 2; i++) {
        a[i] = h[i];
        b[i] = 2 * (h[i] + h[i + 1]);
        c[i] = h[i + 1];
        d[i] = 3 * ((y[i + 2] - y[i + 1]) / h[i + 1] - (y[i + 1] - y[i]) / h[i]);
    }
    b[n - 2] = 2 * h[n - 2];
    d[n - 2] = 3 * ((y[n - 1] - y[n - 2]) / h[n - 2] - (y[n - 2] - y[n - 3]) / h[n - 3]);

    // 追赶法求解三对角矩阵
    vector<double> g(n - 1), p(n - 1);
    g[0] = b[0];
    p[0] = d[0] / g[0];
    for (int i = 1; i < n - 1; i++) {
        g[i] = b[i] - a[i - 1] * c[i - 1] / g[i - 1];
        p[i] = (d[i] - a[i - 1] * p[i - 1]) / g[i];
    }

    // 求解插值函数的系数
    vector<double> c2(n - 1), c3(n - 1);
    c2[n - 2] = p[n - 2];
    for (int i = n - 3; i >= 0; i--) {
        c2[i] = p[i] - c[i] * c2[i + 1] / g[i];
    }
    for (int i = 0; i < n - 1; i++) {
        c3[i] = (c2[i + 1] - c2[i]) / (3 * h[i]);
    }
    vector<double> c1(n - 1);
    for (int i = 0; i < n - 1; i++) {
        c1[i] = (y[i + 1] - y[i]) / h[i] - h[i] * (c2[i + 1] + 2 * c2[i]) / 3;
    }

    // 计算插值函数的值
    int k = 0;
    for (int i = 0; i < n - 1; i++) {
        if (xi >= x[i] && xi <= x[i + 1]) {
            k = i;
            break;
        }
    }
    double dx = xi - x[k];
    return y[k] + c1[k] * dx + c2[k] * dx * dx + c3[k] * dx * dx * dx;
}

int main() {
    vector<double> x = { 0, 1, 2, 3, 4 };
    vector<double> y = { 0, 1, 8, 27, 64 };
    double xi = 2.5;
    double yi = splineInterpolation(x, y, xi);
    cout << "spline interpolation of (" << xi << ", " << yi << ")" << endl;
    return 0;
}

此代码使用向量存储输入数据，并使用三对角矩阵的追赶法求解三次样条插值函数的系数。最终，它计算出插值函数在指定点的值，并输出结果