3. 三次样条插值：三次样条的定义、边界条件、三对角矩阵求解
各位同行，咱们接着聊样条插补。前面讲了线性插值和二次插值，说白了就是拿直线和抛物线去拟合曲线。但实际加工中，尤其是模具曲面，曲率变化大，直线段会留下折痕，抛物线又不够灵活。这时候，三次样条就该登场了。
我个人觉得，三次样条是数控插补里最实用的工具之一。它既能保证曲线光滑，又不会像高次多项式那样产生龙格现象——就是那种两头翘起来的振荡。嗯，这里要注意，龙格现象在模具加工里很要命，会让刀具路径出现不必要的波动。
3.1 三次样条的定义
三次样条，说白了就是分段的三次多项式。我们把数据点分成若干段，每段用一个三次函数去拟合，同时保证段与段之间连接处的一阶导数和二阶导数都连续。
数学上，给定 n+1 个数据点 (x₀, y₀), (x₁, y₁), ..., (xₙ, yₙ)，在每个子区间 [xᵢ, xᵢ₊₁] 上，我们定义一个三次多项式：
Sᵢ(x) = aᵢ + bᵢ(x - xᵢ) + cᵢ(x - xᵢ)² + dᵢ(x - xᵢ)³其中 i = 0, 1, ..., n-1。每个区间有 4 个未知数，总共 4n 个未知数。我们需要 4n 个方程才能解出来。
这些方程从哪里来？
- 插值条件：每个数据点必须落在曲线上，共 2n 个方程
- 一阶导数连续：内部节点处左右导数相等，共 n-1 个方程
- 二阶导数连续：内部节点处左右二阶导数相等，共 n-1 个方程
你看，加起来是 4n-2 个方程，还差两个。这两个方程就由边界条件来提供。
核心要点：三次样条的本质是「分段三次、整体二阶连续」。它比线性插值光滑，比高次多项式稳定。我在做汽车覆盖件模具时，经常用三次样条来拟合A级曲面，效果比B样条差不了太多，但计算量小得多。
3.2 边界条件
边界条件就是给首尾两端加上约束。常见的边界条件有三种：
| 类型 | 数学表达 | 适用场景 |
|---|---|---|
| 自然边界 | S''(x₀) = 0, S''(xₙ) = 0 | 曲线两端自由，无约束 |
| 固定边界 | S'(x₀) = f'(x₀), S'(xₙ) = f'(xₙ) | 已知端点切线方向 |
| 周期边界 | S'(x₀) = S'(xₙ), S''(x₀) = S''(xₙ) | 闭合曲线 |
我个人习惯用自然边界条件，因为大多数模具加工中，刀具路径的起点和终点不需要额外约束。但如果你做的是闭合轮廓的加工，比如圆形型腔，那就必须用周期边界条件。
避坑指南：我曾经在加工一个复杂曲面时，用了固定边界条件，但给定的端点导数不准确，结果曲线在两端出现了明显的扭曲。后来我改用自然边界条件，问题就解决了。所以，除非你非常确定端点的导数，否则优先用自然边界。
3.3 三对角矩阵求解
好了，现在我们来解决核心问题——怎么求解这 4n 个未知数？
其实不用直接解 4n 个方程。我们可以利用二阶导数连续的条件，把问题简化为求解二阶导数 mᵢ = S''(xᵢ)。
推导过程我就不细写了，直接给结果。对于内部节点，我们有：
λᵢ · mᵢ₋₁ + 2 · mᵢ + μᵢ · mᵢ₊₁ = 6 · f[xᵢ₋₁, xᵢ, xᵢ₊₁]其中：
- λᵢ = hᵢ / (hᵢ₋₁ + hᵢ)
- μᵢ = hᵢ₋₁ / (hᵢ₋₁ + hᵢ)
- hᵢ = xᵢ₊₁ - xᵢ
- f[xᵢ₋₁, xᵢ, xᵢ₊₁] 是二阶差商
加上边界条件，就构成了一个三对角线性方程组：
| 2 μ₀ 0 ... 0 | | m₀ | | d₀ |
| λ₁ 2 μ₁ ... 0 | | m₁ | | d₁ |
| 0 λ₂ 2 ... 0 | | m₂ | = | d₂ |
| ... ... ... ... ... | | ... | | ... |
| 0 0 0 λₙ 2 | | mₙ | | dₙ |这个矩阵只有主对角线和两条副对角线有非零元素，所以叫三对角矩阵。求解它，用追赶法（Thomas算法）最合适，时间复杂度只有 O(n)。
追赶法核心步骤：
- 追：从上到下，消去下对角线元素，同时更新主对角线和右端项
- 赶：从下到上，回代求解未知数
说白了就是高斯消元法的简化版，专门针对三对角矩阵设计的。
下面是我常用的 C 语言实现：
void thomas_algorithm(double *a, double *b, double *c,
double *d, double *x, int n) {
// a: 下对角线, b: 主对角线, c: 上对角线
// d: 右端项, x: 解向量, n: 方程个数
// 追的过程
for (int i = 1; i < n; i++) {
double w = a[i] / b[i-1];
b[i] -= w * c[i-1];
d[i] -= w * d[i-1];
}
// 赶的过程
x[n-1] = d[n-1] / b[n-1];
for (int i = n-2; i >= 0; i--) {
x[i] = (d[i] - c[i] * x[i+1]) / b[i];
}
}注意：追赶法要求主对角线元素占优，即 |bᵢ| ≥ |aᵢ| + |cᵢ|。好在三次样条生成的矩阵天然满足这个条件，所以不用担心数值稳定性问题。
3.4 知识体系结构图
为了让大家更直观地理解三次样条插值的整体逻辑，我画了一张结构图：
这张图把三次样条插值的三个核心环节串起来了：从输入数据点出发，经过分段三次多项式的定义，加上边界条件的约束，最后通过三对角矩阵求解得到光滑曲线。
3.5 实际应用中的注意事项
最后，我结合自己的经验，给大家几点实用建议：
- 数据点间距要均匀：如果 x 坐标分布不均匀，会导致某些区间曲率变化剧烈。我一般会先做等距重采样，再做样条插值。
- 避免过冲：三次样条在数据点之间可能会出现小幅过冲。如果加工精度要求极高，可以考虑用单调三次样条（PCHIP）来抑制过冲。
- 实时性考虑：追赶法虽然快，但如果在 CNC 系统里做实时插补，建议预计算好系数，运行时只做多项式求值。
个人经验：我在做五轴联动加工时，经常用三次样条来插补刀轴矢量。这时候边界条件要特别小心——如果刀轴在起点和终点有突变，加工出来的表面会有明显的刀痕。我一般会在两端各加一个虚拟点，让刀轴变化更平缓。
好了，三次样条插值就讲到这里。内容不少，但核心就三件事：分段三次多项式、选对边界条件、用追赶法解三对角矩阵。你想想看，是不是比想象中简单？
公众号:蓝海资料掘金营，微信deep3321