三次样条插值

三次样条插值（cubic spline interpolation）在相邻节点之间分别使用次数至多三的多项式，并要求片段在内部节点处的函数值、一阶导数和二阶导数连续。它以多段低次曲线匹配数据，区别于用一个高次多项式穿过全部节点。端点还须给出适当边界条件，才能唯一确定整个样条。

取节点 $(0, 0), (1, 1), (2, 0)$ 。只要求每段通过端点有许多选择；即便加上内部一、二阶导连续，也还需要两个端点条件。本文采用自然边界条件 $S^{″} (0) = S^{″} (2) = 0$ ，由此得到一个可完整手算的自然三次样条。

以节点二阶导数重建片段

设 $x_{0} < \dots < x_{n - 1}$ ，数据值为 $y_{i}$ ，区间长度为 $h_{i} = x_{i + 1} - x_{i} > 0$ ，记节点二阶导数 $M_{i} = S^{″} (x_{i})$ 。由于三次多项式的二阶导是一次函数，在第 $i$ 段上必须有 $S^{″} (x) = M_{i} \frac{x_{i + 1} - x}{h_{i}} + M_{i + 1} \frac{x - x_{i}}{h_{i}} .$ 对它积分两次，再用 $S (x_{i}) = y_{i}, S (x_{i + 1}) = y_{i + 1}$ 确定两个积分常数，得到 $\begin{matrix} S_{i} (x) = & \frac{M_{i} (x_{i + 1} - x)^{3} + M_{i + 1} (x - x_{i})^{3}}{6 h_{i}} \\ + (y_{i} - \frac{M_{i} h_{i}^{2}}{6}) \frac{x_{i + 1} - x}{h_{i}} + (y_{i + 1} - \frac{M_{i + 1} h_{i}^{2}}{6}) \frac{x - x_{i}}{h_{i}} . \end{matrix}$ 代入端点会留下相应 $y_{i}$ ，再求两次导数则回到前一式，所以这份重建同时满足高度和曲率条件。^[1]

一阶连续给出三对角系统

第 $i$ 段右端的一阶导数为 ${S_{i}}^{'} (x_{i + 1}) = \frac{y_{i + 1} - y_{i}}{h_{i}} + \frac{h_{i}}{6} (M_{i} + 2 M_{i + 1}),$ 左端为 $(y_{i + 1} - y_{i}) / h_{i} - h_{i} (2 M_{i} + M_{i + 1}) / 6$ 。在内部节点 $x_{i}$ ，令左片段的右导数等于右片段的左导数，得到 $h_{i - 1} M_{i - 1} + 2 (h_{i - 1} + h_{i}) M_{i} + h_{i} M_{i + 1} = 6 (\frac{y_{i + 1} - y_{i}}{h_{i}} - \frac{y_{i} - y_{i - 1}}{h_{i - 1}}) .$ 自然边界取 $M_{0} = M_{n - 1} = 0$ ，只需求解 $n - 2$ 个内部未知量。每行最多涉及相邻三个未知量，因此是三对角线性系统。

它的系数矩阵实对称，主对角元为 $2 (h_{i - 1} + h_{i})$ ，严格大于该行非对角元绝对值之和。对任意非零向量 $z$ ，利用 $2 z_{i} z_{i + 1} \geq - z_{i}^{2} - z_{i + 1}^{2}$ 估计二次型，得到 $z^{𝖳} A z \geq \sum_{i} (A_{i i} - \sum_{j \neq i} | A_{i j} |) z_{i}^{2} > 0.$ 所以矩阵正定而可逆，节点二阶导数唯一存在，再由重建公式得到自然样条的存在与唯一性。严格递增节点保证这里每个 $h_{i}$ 都正，重复节点不能直接代入同一系统。

三节点的全部计算

主例两段长度均为一，唯一内部方程为 $4 M_{1} = 6 (- 1 - 1) = - 12$ ，故 $M_{1} = - 3$ 。重建为 $S (x) = {\begin{matrix} \frac{3}{2} x - \frac{1}{2} x^{3}, & 0 \leq x \leq 1, \\ \frac{3}{2} (2 - x) - \frac{1}{2} (2 - x)^{3}, & 1 \leq x \leq 2. \end{matrix}$ 两段在一处值均为一，一阶导均为零，二阶导均为负三；两外端二阶导为零。中点查询 $x = 1 / 2$ 给出 $S = 3 / 4 - 1 / 16 = 11 / 16$ 。

自然样条通过零零、一一、二零，节点二阶导数零、负三、零，中点二分之一处值为十六分之十一 — 左右片段在峰值处不仅高度相同，切线和二阶导数也相同。

变式。 数据改成 $(0, 1), (1, 3), (2, 5)$ ，相邻斜率都为二，系统右端为零，故全部 $M_{i} = 0$ ，样条就是直线 $1 + 2 x$ 。自然样条准确重现直线，但不一定重现任意二次多项式：对节点值 $0, 1, 4$ ，自然条件给出 $M = (0, 3, 0)$ ，第一段是 $(x^{3} + x) / 2$ ，在二分之一处为 $5 / 16$ ，并非 $x^{2}$ 的 $1 / 4$ 。真正的二次函数端点二阶导为二，与自然条件不合。

线性时间构造与二分查询

下面实现自然边界条件。输入有限实数、节点严格递增，至少两个节点；返回不可变的节点二阶导元组及求值函数，不改写输入。查询限定在数据区间内，避免把最后一个三次片段的延长误当成已保证的外推。

Python 3

from bisect import bisect_right
from math import isfinite

def natural_cubic_spline(nodes, values):
    nodes = tuple(map(float, nodes))
    values = tuple(map(float, values))
    n = len(nodes)
    if n < 2 or len(values) != n:
        raise ValueError("need equal-length data with at least two nodes")
    if not all(isfinite(z) for z in nodes + values):
        raise ValueError("finite data required")
    h = [nodes[i+1] - nodes[i] for i in range(n-1)]
    if any(step <= 0 for step in h):
        raise ValueError("nodes must increase strictly")
    m = n - 2
    inside = [0.0] * m
    if m:
        diagonal = [2*(h[i-1] + h[i]) for i in range(1, n-1)]
        rhs = [6*((values[i+1]-values[i])/h[i]
                 - (values[i]-values[i-1])/h[i-1])
               for i in range(1, n-1)]
        for k in range(1, m):
            factor = h[k] / diagonal[k-1]
            diagonal[k] -= factor * h[k]
            rhs[k] -= factor * rhs[k-1]
        inside[-1] = rhs[-1] / diagonal[-1]
        for k in range(m-2, -1, -1):
            inside[k] = (rhs[k] - h[k+1]*inside[k+1]) / diagonal[k]
    second = (0.0, *inside, 0.0)
    def evaluate(x):
        x = float(x)
        if not isfinite(x) or x < nodes[0] or x > nodes[-1]:
            raise ValueError("query outside the data interval")
        i = min(bisect_right(nodes, x) - 1, n - 2)
        left, right, step = x - nodes[i], nodes[i+1] - x, h[i]
        return ((second[i]*right**3 + second[i+1]*left**3)/(6*step)
                + (values[i]-second[i]*step**2/6)*right/step
                + (values[i+1]-second[i+1]*step**2/6)*left/step)
    return second, evaluate

second, spline = natural_cubic_spline([0, 1, 2], [0, 1, 0])
assert second == (0.0, -3.0, 0.0)
assert spline(0.5) == 11/16
assert spline(0) == 0 and spline(1) == 1 and spline(2) == 0
_, line = natural_cubic_spline([0, 1, 2], [1, 3, 5])
assert line(0.5) == 2
_, two = natural_cubic_spline([0, 2], [1, 5])
assert two(1) == 3
try:
    natural_cubic_spline([0, 0, 1], [0, 1, 2])
except ValueError:
    pass
else:
    raise AssertionError("repeated node accepted")
print(second, spline(0.5))

消元每行只更新相邻一项，回代也每行一次，构造时间与存储均为 $O (n)$ ；查询通过二分查找所在区间，耗时 $O (\log n)$ ，片段代入为常数成本。精确算术下消元主元为正：首个内部主元大于其右侧系数，递推中减去的项小于下一行的左侧系数，故下一主元仍为正。浮点计算在极端节点尺度下仍需误差控制，正定性不是无限精度的实现保证。

端点条件与形状限制

自然条件是两端二阶导为零，不是两端一阶导为零。若已知端点斜率，可使用夹持样条；周期数据则须施加周期连接条件。它们产生不同的线性系统，也可能产生不同插值曲线。

三次样条的二阶光滑并不保证单调或无过冲。数据本身带噪时，精确穿过每点也不同于平滑拟合；选择插值、形状保持方法还是平滑样条，取决于需要保留哪些条件。

参考文献

↑ Tobin A. Driscoll、Richard J. Braun：Fundamentals of Numerical Computation: Cubic splines，分段三次与连接条件、端点条件。访问日期：2026-10-22。

[1] Tobin A. Driscoll、Richard J. Braun：Fundamentals of Numerical Computation: Cubic splines，分段三次与连接条件、端点条件。访问日期：2026-10-22。

[1]