Hermite插值多项式
前几节中介绍的拉格朗日和牛顿插值多项式， 虽然能够用一个高次的多项式穿过所有的采样点。但在工程实际中，很少有人这样构建一个全局的高次多项式，而是将它拆分成一段段的局部曲线，每段曲线都用低次的多项式来表示，最后再拼接成完整的曲线。 之所以这样做，主要是因为龙格现象。在等距采样的情况下，高次的多项式会在两端剧烈的震荡。虽然改变采样点的分布可以缓解这一现象，但并不能根治，而且很多实际情况我们也只能等距采样。 所以人们希望通过分段来规避这一问题。
分段的低次多项式又会带来另外一个问题，在连接处不够平滑。传统的拉格朗日和牛顿插值多项式，只能保证曲线一定经过采样点，但并没有对一阶导数做约束， 所以两端多项式在连接处的一阶导数可能不连续。如果用这样的曲线控制机械臂，就会在连接点处出现突然的加减速，不利于控制也很危险。 Hermite 插值多项式不仅要求插值多项式在采样点上与原函数相同，还要求它们的一阶导数也相同，即所谓的 \(C^1\) 连续性。这样在连接点处，曲线就会丝滑过渡。
1. Hermite插值多项式
定理 1 :如果函数 \(f(x) \in C^1[a,b]\)，即在闭区间 \([a,b]\) 上有连续的一阶导数。\(\{x_0, \cdots, x_n\}\) 是 \([a,b]\) 上 \(n+1\) 个不同的采样点。 那么存在一个次数不超过 \(2n+1\) 的多项式 \(H_{2n+1}(x)\)，在各采样点上的插值和一阶导数与原函数相同，即 \(H_{2n+1}(x_i) = f(x_i), H_{2n+1}'(x_i) = f'(x_i)\)。 并且该多项式是惟一的。
在拉格朗日多项式中，我们曾构造了一个多项式 \(L_{n,i}(x)\)
$$ L_{n,i}(x) = \frac{(x - x_0)\cdots(x - x_{i-1})(x - x_{i+1}) \cdots (x - x_n)}{(x_i - x_0)\cdots(x_i - x_{i-1})(x_i - x_{i+1}) \cdots (x_i - x_n)} \Longrightarrow L_{n,i}(x_j) = \begin{cases} 0, \text{if} \quad i \neq j \\ 1, \text{if} \quad i = j \end{cases} $$在此基础上，我们构造如下两个函数:
$$ \begin{equation}\label{f1} H_{n,i}(x) = \left[1 - 2(x - x_i) L_{n,i}'(x_i) \right] L_{n,i}^2(x), \qquad \hat{H}_{n,i}(x) = (x - x_i)L_{n,i}^2(x) \end{equation} $$显然 \(i \neq j\) 时，\(H_{n,i}(x_j) = 0, \hat{H}_{n,i}(x_j) = 0\)。由于 \(L_{n,i}(x_i) = 1\)，所以有
$$ H_{n,i}(x_i) = [1 - 2(x_i - x_i)L_{n,i}'(x_i)] \cdot 1 = 1, \qquad \hat{H}_{n,i}(x_i) = (x_i - x_i) \cdot 1^2 = 0 $$那么由此构造的多项式:
$$ \begin{equation}\label{f2} H_{2n+1}(x) = \sum_{i=0}^n \left[ f(x_i)H_{n,i}(x) \right] + \sum_{i=0}^n \left[f'(x_i)\hat{H}_{n,i}(x)\right] \end{equation} $$一定满足 \(H_{2n+1}(x_i) = f(x_i), i \in \{0,1,\cdots,n\}\)。
我们再来证明 \(H_{2n+1}'(x_i) = f'(x_i)\)。对于 \(H_{n,i}(x)\) 根据求导法则，有:
$$ H_{n,i}'(x) = -2L_{n,i}'(x_i) \cdot L_{n,i}^2(x) + \left[1 - 2(x - x_i) L_{n,i}'(x_i) \right] \cdot 2 L_{n,i}(x) L_{n,i}'(x) $$因为 \(L_{n,i}(x_j) = 0, i \neq j\)，所以 \(H_{n,i}'(x_j) = 0\)。而 \(H_{n,i}'(x_i) = -2L_{n,i}'(x_i) + 2L_{n,i}'(x_i) = 0\)。
对于 \(\hat{H}_{n,i}(x)\) 根据求导法则，有:
$$ \hat{H}_{n,i}'(x) = L_{n,j}^2(x) + (x - x_i) \cdot 2 L_{n,i}(x) L_{n,i}'(x) = L_{n,i}(x) \cdot \left[L_{n,j} + 2(x - x_i)L_{n,i}'(x)\right] $$显然有 \(\hat{H}_{n,i}'(x_j) = 0, i \neq j\)，并且 \(\hat{H}_{n,i}'(x_i) = 1\)。综合 \(H_{n,i}'(x), \hat{H}_{n,i}'(x)\)，有
$$ H_{2n+1}'(x) = \sum_{i=0}^n \left[f(x_i) H_{n,i}'(x)\right] + \sum_{i=0}^n \left[f'(x_i) \hat{H}_{n,i}(x) \right] = \sum_{i=0}^n \left[f'(x_i) \hat{H}_{n,i}(x) \right] $$所以有 \(H_{2n+1}'(x_i) = f'(x_i), i \in \{0,1,\cdots,n\}\)。式\((\ref{f1})(\ref{f2})\)被称为Hermite 插值多项式。
我们再来通过反证法证明 Hermite 插值多项式的唯一性。假设存在两个多项式 \(P(x), Q(x)\)，它们的次数都 \(≤ 2n + 1\)。并且都完全满足 \(2n + 2\) 个插值条件， 即，对于 \(\{x_0, x_1, \cdots, x_n\}\) 共 \(n+1\) 个不相同的采样点，都有:
$$ \begin{aligned} P(x_i) = f(x_i) & \qquad P'(x_i) = f'(x_i) \\ Q(x_i) = f(x_i) & \qquad Q'(x_i) = f'(x_i) \\ \end{aligned} $$两个多项式相减，得到 \(R(x) = P(x) - Q(x)\)。由于 \(P(x),Q(x)\) 的次数都 \(≤ 2n + 1\)，所以它们的差的次数也满足:
$$ \text{deq}(R) ≤ 2n + 1 $$对于 \(n+1\) 个采样点，\(R(x)\) 有如下关系:
$$ \begin{array}{l} R(x_i) & = P(x_i) - Q(x_i) & = f(x_i) - f(x_i) & = 0 \\ R'(x_i) & = P'(x_i) - Q'(x_i) & = f'(x_i) - f'(x_i) & = 0 \\ \end{array} $$这意味着，对于每一个采样点 \(x_{i}\)，它不仅是 \(R(x) = 0\) 的根，同时也是 \(R'(x) = 0\) 的根。根据代数相关定理， 若一函数在某点的函数值与导数值皆为 0，则该点必为该函数的重根。所以 \(R(x)\) 至少有 \(n + 1\) 个二重根 \(x_0, \cdots, x_n\)， 其总根数为 \((n+1) \times 2\)。但是一个次数不超过 \(2n+1\) 的多项式，最多只有 \(2n+1\) 个根。这是矛盾的，只可能 \(R(x) \equiv 0\)， 那么有 \(P(x) = Q(x)\)。所以 Hermite 插值多项式是唯一的。
2. 误差余项
定理 2 :如果函数 \(f(x) \in C^{2n+2}[a,b]\)，即在闭区间 \([a,b]\) 上有连续的\(2n+2\)阶导数。\(\{x_0, \cdots, x_n\}\) 是 \([a,b]\) 上 \(n+1\) 个不同的采样点。 那么对于闭区间 \([a,b]\) 上的任意一个 \(x\)，都存在一个数 \(\varepsilon(x)\) 使得在开区间 \((a,b)\) 上有: $$ \begin{equation}\label{f3} f(x) = H_{2n+1}(x) + \frac{(x-x_0)^2\cdots(x-x_n)^2}{(2n+2)!}f^{2n+2}(\varepsilon(x)) \end{equation} $$ 其中 \(H_{2n+1}(x)\) 为上式 \((\ref{f1})(\ref{f2})\) 中定义的 Hermite 插值多项式。
该定理的证明与拉格朗日多项式的误差余项如出一辙，也是通过反复使用中值定理得出。 首先构造一个辅助函数 \(g(t)\)
$$ \begin{equation}\label{f4} g(t)=f(t)-H_{2n+1}(t)-K\cdot \prod _{j=0}^{n}(t-x_{j})^{2} \end{equation} $$对于任意固定点 \(x\)，都可以凑出一个常数 \(K = \frac{f(x) - H_{2n+1}(x)}{\prod_{j=0}^{n} (x-x_j)^2}\) 使得 \(g(x) = 0\)。 显然 \(g(t)\) 在区间 \([a,b]\) 上至少有 \(n + 2\) 个不同的零点 \(x, x_0, \cdots, x_n\)。根据中值定理可以推定，其一阶导数 \(g'(t)\) 在这些零点中间至少有 \(n+1\) 个零点。
根据 Hermite 插值多项式的约束有，\(f'(x_i) = H'_{2n+1}(x_i)\)，代入式 \((\ref{f4})\) 有:
$$ g'(x_i) = f'(x_i) - H'_{2n+1}(x_i) - K \cdot 2(x_i - x_i) \cdot \prod_{i\neq j}(x_i - x_j)^2 - (x_i - x_i)^2\Omega(x_i) = 0 $$上式中 \(\Omega(x_i)\) 是 \(K\cdot \prod _{j=0}^{n}(t-x_{j})^{2}\) 除去 \((x - x_i)^2\) 后的导数项，由于其系数为 \((x_i - x_i)^2 = 0\)，所以我们不用关心 \(\Omega(x_i)\) 的具体形式。 因此 \(x_i, i \in \{0,\cdots,n\}\) 这 \(n+1\) 个采样点也是 \(g'(t)\) 的零点。所以 \(g'(t)\) 共有 \(2n+2\) 个不同的零点。
对 \(g'(t)\) 应用中值定理，可以推定 \(g''(t)\) 有 \(2n+1\) 个零点。如此重复求导，知道 \(2n+2\) 阶导数，\(g^{(2n+2)}(t)\) 在区间 \((a,b)\)内至少还剩 1 个零点， 我们把这个点记为 \(\varepsilon \)，即 \(g^{(2n+2)}(\varepsilon) = 0\)。
由于 \(H_{2n+1}(t)\) 的次数不超过 \(2n+1\)，所以求导 \(2n+2\) 此后直接为 0。\(\prod _{j=0}^{n}(t-x_{j})^{2}\) 则是一个 \(2n+2\) 次的多项式，并且最高次项系数为 1， \(2n+2\) 此求导之后得到常数 \((2n + 2)!\)。整理之后有:
$$ f^{(2n+2)}(\varepsilon) - 0 - K \cdot (2n+2)! = 0 \Longrightarrow K = \frac{f^{(2n+2)}(\varepsilon)}{(2n+2)!} $$最后将 \(K\) 代入
$$ \frac{f(x) - H_{2n+1}(x)}{\prod_{j=0}^{n} (x-x_j)^2} = \frac{f^{(2n+2)}(\varepsilon)}{(2n+2)!} \Longrightarrow f(x) - H_{2n+1}(x) = \frac{(x-x_0)^2\cdots(x-x_n)^2}{(2n+2)!}f^{2n+2}(\varepsilon(x)) $$得到的就是 Hermite 插值多项式的误差余项。同样的 Hermite 的逼近误差也主要有两个方面的因素。其一是 \(f^{(2n+2)}(\varepsilon)\) 表示原函数越平滑，高阶导数越小，误差就越小。 其二是描述采样点位置分布的 \(\prod_{i=0}^{n}(x-x_{i})^{2}\)。
3. Hermite 插值多项式的实现
我们严格按照式\((\ref{f1})(\ref{f2})\)实现了一个朴素的 Hermite 插值。 由于 \(L_{n,i}(x_i), L_{n,i}'(x_i)\) 需要 \(O(n)\) 的乘法运算，\(H_{2n+1}(x)\) 需要将各项类加起来，需要 \(O(n)\) 的加法运算。 所以我们要在两层循环中完成，整体的复杂度是 \(O(n^2)\) 的。暴力展开的过程过于复杂，这里不再提供了。
引入重合采样点，Hermite 插值也可以写成差商的形式。 设有 \(n+1\) 个采样点 \(x_0, x_1, \dots, x_n\)，每个采样点已知函数值\(f(x_i)\)和一阶导数\(f'(x_i)\)。 我们把每个采样点拷贝一份，就可以写出一个长度为 \(2n+2\) 的序列:
$$ z_{0}=x_{0},\;z_{1}=x_{0},\;z_{2}=x_{1},\;z_{3}=x_{1},\;\dots ,\;z_{2n}=x_{n},\;z_{2n+1}=x_{n} $$此时，Hermite 插值多项式可以完全写成差商的形式:
$$ \begin{aligned} H_{2n+1}(x) = &f[z_{0}]+f[z_{0},z_{1}](x-z_{0})+f[z_{0},z_{1},z_{2}](x-z_{0})(x-z_{1})+\dots \\ &+f[z_{0},z_{1},\dots ,z_{2n+1}](x-z_{0})(x-z_{1})\cdots (x-z_{2n}) \end{aligned} $$其中 \(f[z_0, z_1, \dots]\) 是广义差商。说它广义，主要体现在两个重复采样点的差商计算上。因为 \(z_0 = z_1 = x_0\)，所以直接按照上一节介绍的差商公式计算， 就会出现 \(\frac{f(x_0) - f(x_0)}{x_0 - x_0} = \frac{0}{0}\) 的情况。但是根据极限和导数的定义，有:
$$ \lim_{z_1 \to z_0} \frac{f(z_1) - f(z_0)}{z_1 - z_0} = f'(x_0) $$所以，我们定义这种重复采样点的差商为 \(f[x_i, x_i] = f'(x_i)\)。如此就可以写出一个差商表来:
$$ \begin{array}{l} z_0 = x_0: \quad & f[z_0] = f(x_0) & f[z_0, z_1] = f'(x_0) \quad & f[z_0, z_1, z_2] = \frac{f[z_1,z_2]-f[z_0,z_1]}{z_2 - z_0} & f[z_0,z_1,z_2,z_3] \quad & f[z_0,z_1,z_2,z_3,z_4] & f[z_0,z_1,z_2,z_3,z_4,z_5] \\ z_1 = x_0: \quad & f[z_1] = f(x_0) & f[z_1, z_2] = \frac{f[z_2]-f[z_1]}{z_2 - z_1} \quad & f[z_1, z_2, z_3] = \frac{f[z_2,z_3]-f[z_1,z_2]}{z_3 - z_1} & f[z_1,z_2,z_3,z_4] \quad & f[z_1,z_2,z_3,z_3,z_5] \\ z_2 = x_1: \quad & f[z_2] = f(x_1) & f[z_2, z_3] = f'(x_1) \quad & f[z_2, z_3, z_4] = \frac{f[z_3,z_4]-f[z_2,z_3]}{z_4 - z_2} & f[z_2,z_3,z_4,z_5] \quad & \\ z_3 = x_1: \quad & f[z_3] = f(x_1) & f[z_3, z_4] = \frac{f[z_4]-f[z_3]}{z_4 - z_3} \quad & f[z_3, z_4, z_5] = \frac{f[z_4,z_5]-f[z_3,z_4]}{z_5 - z_3} \\ z_4 = x_2: \quad & f[z_4] = f(x_2) & f[z_4, z_5] = f'(x_2) \quad \\ z_5 = x_2: \quad & f[z_5] = f(x_2) \end{array} $$如此，我们就可以实现一个 差商形式的 Hermite 插值， 经过 \(O(n^2)\) 的初始化构造差商表之后，就可以在 \(O(n)\) 的时间复杂度内完成插值计算。