1

我的大胆主张是三次样条的 Octave 实现，与Wolfram 的 Mathworldinterp1(..., "spline")中概述的“自然三次样条”算法不同。我已经编写了自己的后者实现，并将其与函数的输出进行了比较，结果如下：interp1(..., "spline")

样条比较

我发现，当我尝试与 4 个点进行相同的比较时，解也不同，而且，Octave 解与将单个三次多项式拟合到所有四个点相同（实际上并没有为三个间隔生成分段样条） .

我还尝试深入了解 Octave 的样条曲线实现，发现它太迟钝，无法在 5 分钟内阅读和理解。

我知道在实现三次样条时可以选择一些边界条件选项（“自然”与“钳制”）。我的实现使用“自然”边界条件（其中两个外部点的二阶导数设置为零）。

如果Octave 的三次样条确实与标准三次样条不同，那么它到底是什么？

编辑：

上面的比较图中显示的两个解决方案的二阶差绘制在此处：

二阶有限差分

首先，在 Octave 的情况下似乎只有两个三次多项式：一个适合前两个区间，一个适合后两个区间。其次，它们显然没有使用“自然”样条曲线，因为极端的二阶导数不会趋于零。

此外，我认为我在中间（即第三个）点的实现的二阶差为零只是一个巧合，而不是算法所要求的。对一组不同的点重复此测试将确认/反驳这一点。

4

1 回答 1

2

不同的结束条件解释了您的实现与 Octave 之间的差异。Octave 使用非结条件（取决于输入）

看help spline

为了解释您的观察：由于非结条件，三阶导数在第二次和 (n-1) 次中断处是连续的，这就是为什么 Octave 的二阶导数看起来“中断”较少的原因，因为它是前两段和后两段上的连续直线。如果你看三阶导数，你可以更清楚地看到效果——三阶导数只在第三个断点（中间）是不连续的

x = 1:5;
y = rand(1,5);
xx = linspace(1,5);
pp = interp1(x, y, 'spline', 'pp');
yy = ppval(pp, xx);
dyy = ppval(ppder(pp, 3), xx);
plot(xx, yy, xx, dyy);

数据结构也pp看起来像这样

pp =

  scalar structure containing the fields:

    form = pp
    breaks =

       1   2   3   4   5

    coefs =

       0.427823  -1.767499   1.994444   0.240388
       0.427823  -0.484030  -0.257085   0.895156
      -0.442232   0.799439   0.058324   0.581864
      -0.442232  -0.527258   0.330506   0.997395

    pieces =  4
    order =  4
    dim =  1
    orient = first
于 2016-10-29T05:14:40.310 回答