三次样条:

三次样条是一种使用满足给定 m 个控制点的三次多项式的样条。为了推导三次样条的解，我们假设端点处的二阶导数为 0，这进而提供了一个边界条件，在 m-2 个方程的基础上增加了两个方程，使它们变得可解。一维三次样条的方程组可以表示为：

Si (x) = ai + bi x + ci x^{2} + d_i x^{3}

我们取函数 y = f(x) 的一组点 [xi, yi]，i = 0, 1, …, n。三次样条插值是一条分段连续的曲线，它会经过表中的每一个值。

• 以下是 K=3 阶样条的条件：
• s 的定义域在区间 [a, b] 内。
• S, S‘, S‘‘ 都是 [a,b] 上的连续函数。

S(x) =\begin{bmatrix} S0 (x), x \epsilon [x0, x1] \\ S1 (x), x \epsilon [x1, x2] \\ … \\ … \\ … \\ S{n-1} (x), x \epsilon [x1, x_2] \end{bmatrix}

这里 Si(x) 是将在子区间 [xi, xi+1] 上使用的三次多项式。

因为有 n 个区间，且每个方程有 4 个系数，为了解这个样条，我们需要 4n 个参数。我们可以根据每个三次样条方程应满足两端值这一事实得到 2n 个方程：

Si(xi) = yi \\ Si(x{i+1}) = y{i+1}

上述三次样条方程不仅应该是连续且可微的，而且在控制点处的一阶和二阶导数也必须是连续的。

S{i-1}^{‘}(xi) = S{i}^{‘}(xi)

以及

S{i-1}^{‘‘}(xi) = S{i}^{‘‘}(xi)

对于 1, 2, 3…n-1，这提供了 2n -2 个方程约束。所以，我们需要再增加 2 个方程来求解上述三次样条。为此，我们将使用一些自然边界条件。

自然三次样条:

在自然三次样条中，我们假设样条在边界点处的二阶导数为 0：

S^{‘‘}(x0) = 0 \\ S^{‘‘}(xn) = 0

现在，由于 S(x) 是一个三阶多项式，我们知道 S‘‘(x) 是一个进行插值的线性样条。因此，首先我们构造 S‘‘(x)，然后对其积分两次得到 S(x)。

现在，让我们假设 ti = xi (i= 0, 1,…n)，以及 zi = S^{‘‘}(xi) (i=0,1,2..n)，并且根据自然边界条件有 z0= zn =0。对三次样条进行两次微分会得到一个线性样条，可以写成：

Si^{‘‘}(x) = zi \frac{x- t{i+1}}{ti – t{i+1}} + z{i+1} \frac{x – ti}{t{i+1} – t_i}

其中，

hi = t{i+1} – t_i ; t=0,1, 2,…,n

现在，方程变为：

S‘‘(x) = z{i+1} \frac{x- ti}{hi} + z{i} \frac{t{i+1} – x}{hi}

对该方程积分两次以获得三次样条：

S(x) =\frac{z{i+1}}{6hi}(x – ti)^{3} + \frac{z{i}}{6hi}( t{i+1} – x)^{3} + Ci(xi -t) + Di(t{i+1} -x)

其中，

Ci = \frac{y{i+1}}{hi} – \frac{z{i+1} hi}{6} \\ Di = \frac{y{i}}{hi} – \frac{z{i} hi}{6}

现在，为了检查导数在 ti 处的连续性，即 S^{‘}{i}(ti) = S^{‘}{i-1}(t_i)。我们首先需要求导并应用该条件：

S‘(x) = \frac{z{i+1}}{2hi}(x – ti)^{2} – \frac{z{i}}{2hi}( t{i+1} – x)^{2} + \frac{1}{hi}(y{i+1}-yi)- \frac{hi}{6}(z{i+1} -zi)

其中，

bi = \frac{1}{hi}(y{i+1}-yi)

将上述方程代入连续性条件并求解，可以得到以下方程：

6(bi -b{i-1}) = h{i-1}z{i-1} + 2(h{i-1} + hi)zi + h{i}z_{i+1}

令 vi = 6(bi – b_{i-1})，上述方程可以写成矩阵形式：

\begin{bmatrix} 1& 0 & 0& .& .& .& .& .& .& .& .& 0& \\ h0& 2(h0 + h1)& h1& & & & & & & & & .& \\ .& h1 & 2(h1 + h2)& h2& & & & & & & & .& \\ .& .& .& .& & & & & & & & .& \\ .& & & .& .& .& & & & & & .& \\ .& & & & .& .& & & & & & .& \\ .& & & & & .& .& & & & & .& \\ .& & & & & & .& .& & & & .& \\ .& & & & & & & & .& h{n-2}& 2(h{n-1} + h{n-2})& h{n-1}& \\ 0& .& .& .& .& .& .& .& .& 0& 0& 1& \\ \end{bmatrix} \begin{bmatrix} z0 \\ z1\\ z2\\ .\\ .\\ .\\ .\\ .\\ z{n-1}\\ z{n} \end{bmatrix} = \begin{bmatrix} 0\\ v1 \\ .\\ .\\ .\\ .\\ .\\ .\\ v_{n-1}\\ 0 \end{bmatrix}

实现示例

• 在这个实现示例中，我们将为函数 f(x) = 1/x 在 2 到 10 之间的点进行样条插值，使用满足自然边界条件的三次样条。

Python3


“

#imports

import matplotlib.pyplot as plt

import numpy as np

from scipy.interpolate import CubicSpline, interp1d

plt.rcParams[‘figure.figsize‘] =(12,8)

x = np.arange(2,10)

y = 1/(x)

apply cubic spline interpolation

cs = CubicSpline(x, y, extrapolate=True)

apply natural cubic spline interpolation

ns = CubicSpline(x, y,bc_type=‘natural‘, extrapolate=True)

Apply Linear interpolation

linear_int = interp1d(x,y)

xs = np.arange(2, 9, 0.1)

ys = linear_int(xs)

plot linear interpolation

plt.plot(x, y,‘o‘, label=‘data‘)

plt.plot(xs,ys, label="interpolation", color=‘green‘)

plt.legend(loc=‘