用java写三次样条插值

以下是Java实现三次样条插值的示例代码：

public class CubicSplineInterpolation {
    private double[] x;
    private double[] y;
    private double[] h;
    private double[] a;
    private double[] b;
    private double[] c;
    private double[] d;
    
    public CubicSplineInterpolation(double[] x, double[] y) {
        this.x = x;
        this.y = y;
        int n = x.length;
        
        // 计算h数组
        h = new double[n - 1];
        for (int i = 0; i < n - 1; i++) {
            h[i] = x[i + 1] - x[i];
        }
        
        // 计算a数组
        a = new double[n];
        for (int i = 1; i < n - 1; i++) {
            a[i] = h[i - 1] / (h[i - 1] + h[i]);
        }
        
        // 计算b数组
        b = new double[n];
        for (int i = 1; i < n - 1; i++) {
            b[i] = 1 - a[i];
        }
        
        // 计算c数组
        c = new double[n];
        for (int i = 1; i < n - 1; i++) {
            c[i] = 6 / (h[i - 1] + h[i]) * ((y[i + 1] - y[i]) / h[i] - (y[i] - y[i - 1]) / h[i - 1]);
        }
        
        // 计算d数组
        double[] u = new double[n - 1];
        double[] z = new double[n];
        u[1] = 2 * (a[1] + b[1]);
        z[1] = c[1] / u[1];
        for (int i = 2; i < n - 1; i++) {
            u[i] = 2 * (a[i] + b[i]) - Math.pow(b[i - 1], 2) / u[i - 1];
            z[i] = (c[i] - b[i - 1] * z[i - 1]) / u[i];
        }
        u[n - 2] = a[n - 2] * (2 - b[n - 3]) + u[n - 3];
        z[n - 1] = (c[n - 1] - b[n - 2] * z[n - 2]) / u[n - 2];
        d = new double[n - 1];
        for (int i = n - 2; i >= 1; i--) {
            d[i] = (z[i] - b[i] * d[i + 1]) / a[i];
        }
        
        // 计算c数组
        c[0] = c[n - 1] = 0;
        for (int i = 1; i < n - 1; i++) {
            c[i] = d[i] * h[i - 1] / 6 + d[i + 1] * h[i - 1] / 3 - (y[i + 1] - y[i]) / h[i] + (y[i] - y[i - 1]) / h[i - 1];
        }
    }
    
    public double interpolate(double x0) {
        int n = x.length;
        int i = 0;
        while (i < n - 1 && x[i + 1] < x0) {
            i++;
        }
        double t = (x0 - x[i]) / h[i];
        double y0 = (1 - t) * y[i] + t * y[i + 1] + t * (1 - t) * (a[i] * (1 - t) + b[i] * t) * Math.pow(h[i], 2) / 6 + t * (1 - t) * (a[i] * Math.pow(1 - t, 2) + b[i] * Math.pow(t, 2)) * Math.pow(h[i], 2) / 6 + t * (1 - t) * (1 - a[i]) * y[i] + t * (1 - t) * a[i] * y[i + 1];
        return y0;
    }
}

使用示例：

double[] x = {1, 2, 3, 4, 5};
double[] y = {1, 4, 9, 16, 25};
CubicSplineInterpolation interp = new CubicSplineInterpolation(x, y);
double y0 = interp.interpolate(2.5); // 计算x=2.5处的插值