插值

插值的是通过已知数据点构建新数据点的过程。其目标是找到一个函数，该函数在给定的数据点处精确匹配这些点，然后用这个函数来估计未知点的值。插值函数会精确通过所有已知数据点，但可能会在数据点之间出现过度摆动（尤其是在使用高阶多项式插值时）。

1 插值方法

1.1 背景

考虑一组点 $(x, y)$ ，

(x_0,y_0),\:(x_1,y_1),\:\ldots,\:(x_n,y_n).

$x$ -值被称为横坐标或节点。假设 $y$ 值可由函数 $f$ 计算得到，即

y_i=f(x_i)

图1：三个横坐标 $(x_0,x_1,x_2)$ 和两个可能的函数 $f(x)$ 的插值多项式表达。对于给定的点， $p(x)$ 可能并不是函数 $f(x)$ (右图) 较好的拟合

多项式插值，简单来说，就是找一条曲线（我们称之为多项式），让它恰好穿过我们预先指定的几个点。假设我们有 $n+1$ 个这样的点，分别是 $(x_0,y_0),\dots,(x_n,y_n)$ 。

那么，我们的目标就是找到一个恰当的多项式 $p$ ，它能满足以下条件：对于每一个我们指定的 $x_i$ ，多项式 $p$ 在这点的值都应该是 $y_i$ 。换句话说，这个多项式在这些特定的点上的高度，应该和我们给定的高度完全一致。

这里有两个重要的保证：

唯一性：如果你找到了这样一条曲线（多项式），那么它是独一无二的。这意味着在相同条件下，你不可能找到另一条完全不同的曲线也能穿过所有这些点。
存在性：无论你给定的点是什么样的，你总能找到这样一条多项式曲线。也就是说，这样的曲线总是存在的。

虽然我们可以用不同的方法来找到这条曲线，但不管用哪一种方法，最终得到的曲线（或者说多项式）都是一样的。这就意味着，所有方法都会指向同一个答案，即那条唯一的、能够穿过所有给定点的多项式曲线。

1.2 一般情况

当我们处理多项式 $p$ 时，我们可以通过一组特殊的函数，我们称之为基 $\phi_0, \cdots, \phi_n$ ，来表达它。这种表达方式让 $p$ 变成了这些基函数的一个组合，每个基函数都乘以一个系数 $c_i$ ，即：

p = c_0\phi_0 + c_1\phi_1 + \cdots + c_n\phi_n.

要想让这个多项式真正符合我们的需要（也就是通过我们指定的点），我们需要找到恰当的系数 $c_i$ ，使得对于每个点 $(x_i, y_i)$ ，多项式在 $x_i$ 处的值正好是 $y_i$ 。这实际上就是解一个有 $n+1$ 个方程和 $n+1$ 个未知数的线性方程组。

当然，为了让这个过程在实际中更加可行，我们需要选取一些特别的基函数 $\phi_i$ 。我们希望这些基函数能让多项式 $p$ 的构建和计算变得简单、高效，同时对计算误差具有一定的抵抗力。更具体来说，我们希望基函数能满足以下条件：

构建起来相对容易。
计算系数 $c_i$ 的过程简单明了。
在数值计算中相对稳定，不会因为小小的计算误差就导致结果大相径庭。
一旦找到了系数 $c_i$ 和基函数 $\phi_i$ ，评估多项式 $p(x)$ 的值就变得很容易。
这些基函数具有一定的灵活性，能适应数据的微小变化。

1.3 单项式插值

最直接、最简单的基函数选择就是单项式集合，即让 $\phi_i$ 就是 $x^i$ 。这样一来，多项式就变成了：

p(x) = c_0 + c_1x + \cdots + c_nx^n.

用这种方法可以建立一个方程组并求解：

如果给定 $n$ 个横纵坐标分别互不相同的点 $(x_i,y_i),i=1,2,\ldots,n$ , 存在这样的一个多项式

p(x)=\sum_{i=0}^{n-1}c_i\:x^i

根据构造条件，有

\begin{cases}p(x_1)=y_1\\p(x_2)=y_2\\\cdots\\p(x_n)=y_n\end{cases}

将上述线性方程组中的 $c_i,i=1,2,\ldots,n$ ,视为未知量，其系数矩阵的行列式 $A$ 恰好为范德蒙行列式，

V(x_1,x_2,\cdots,x_n)=\begin{bmatrix}1&1&\cdots&1\\x_1&x_2&\cdots&x_n\\x_1^2&x_2^2&\cdots&x_n^2\\\vdots&\vdots&&\vdots\\x_1^{n-1}&x_2^{n-1}&\cdots&x_n^{n-1}\end{bmatrix}

故(数学归纳法证明)

det(A)=\prod_{1\leq i<j\leq n}(x_j-x_i)

又 $x_i$ 互不相同， $i=1,2,\ldots,n$ ,因此 $det(A)\neq0$ ,方程组有唯一解

令：

\left\{\begin{array}{ccc}y_0&=&\alpha_0+\alpha_1x_0+\alpha_2x_0^2+\cdots+\alpha_nx_0^n\\y_1&=&\alpha_0+\alpha_1x_1+\alpha_2x_1^2+\cdots+\alpha_nx_1^n\\&&\vdots\\y_n&=&\alpha_0+\alpha_1x_n+\alpha_2x_n^2+\cdots+\alpha_nx_n^n\end{array}\right.

\mathbf{Y}=\left[\begin{array}{c}y_0\\y_1\\\vdots\\y_n\end{array}\right],\quad\mathbf{B}=\left[\begin{array}{ccccc}1&x_0&x_0^2&\ldots&x_0^n\\1&x_1&x_1^2&\ldots&x_1^n\\\vdots&\vdots&\vdots&\ddots&\vdots\\1&x_n&x_n^2&\ldots&x_n^n\end{array}\right],\quad\mathbf{\alpha}=\left[\begin{array}{c}\alpha_0\\\alpha_1\\\vdots\\\alpha_n\end{array}\right]

即：通过求解关于 $\alpha$ 的线性方程组 $\mathbf{Y}=\mathbf{B}\alpha$ 即可得到插值这组点的 $n$ 次多项式曲线

p(x)=\alpha_0+\alpha_1x+\alpha_2x^2+\cdots+\alpha_nx^n

尽管这个方程组的解法并不总是最理想的。对于较小的 $n$ （多项式的次数不是特别高），这个方法还凑活。但随着 $n$ 的增加，这个方法的局限性就开始显现出来了。方程可能会变得复杂难解，而且在理论探索和实际应用中，这种方法的效果也不是特别理想。

1.4 拉格朗日插值

假设我们希望构造一个多项式 $p(x)$ ，而且我们希望这个多项式的系数 $c_i$ 尽可能简单明了。理想的情况是，如果我们能让这个多项式的形式是

p(x) = y_0L_0(x) + \cdots + y_nL_n(x)

这里的 $L_0, \cdots, L_n$ 是一组特殊的基础函数。也就是说，每一个基础函数 $L_i$ 都对应一个 $y_i$ ，这样的表达方式称为多项式的拉格朗日形式。

要构建这样的拉格朗日形式，我们需要找到一组特殊的多项式 $L_i$ ，它们需要满足以下条件：

每个 $L_i$ 的次数都是 $n$ ，这意味着它们都是由 $n$ 个项组成的多项式。
对于每个 $i$ ，多项式 $L_i$ 在 $x_i$ 处的值是 1。
同时，对于不等于 $i$ 的任意 $j$ ， $L_i$ 在 $x_j$ 处的值是 0。

为了满足这些条件，我们可以构造一种特殊形式的多项式，它在除 $x_i$ 之外的所有 $x_j$ 处都有根。这样的多项式可以表示为：

L_i(x) = \prod_{j \neq i} \frac{x - x_j}{x_i - x_j}.

这样构造的 $L_i$ 刚好满足上述所有条件。

根据这些基础函数 $L_i$ ，我们得到了一个重要的结论：任何一组点 $(x_0, y_0), \cdots, (x_n, y_n)$ 的插值多项式都可以用拉格朗日形式表示为：

p(x) = \sum_{i=0}^n y_iL_i(x).

这个结论是由 $L_i$ 的构造方式直接得出的，它保证了多项式 $p$ 不仅是一个次数为 $n$ 的多项式，而且确切地通过了所有给定的点。

举个例子，在线性情况下，我们可以通过两个点 $(x_0, y_0)$ 和 $(x_1, y_1)$ 构造一条直线。对于这种情况，拉格朗日基础函数的表达式很简单。

不过，需要注意的是，尽管拉格朗日形式在理论上很有吸引力，并且在函数值的表示上非常自然，但在实际计算中，计算 $L_i$ 并不总是那么直接和友好。拉格朗日形式主要用于理论分析，因为它的表达式易于操作，特别是在进行微分和积分时。

Tips：尽管拉格朗日形式的直接计算可能不是特别高效，但有一些方法可以用来更高效地构造和评估这种形式的多项式。例如，有一种称为重心公式的方法，可以用来高效地计算拉格朗日形式的多项式。

以一个具体的例子来说明如何构建插值多项式。假设我们有以下三个数据点：

\begin{array}{c|c} x_i & y_i \\ \hline 0 & -1 \\ 1/2 & 2/3 \\ 1 & 8/9 \\ \end{array}

我们的目标是找到一个多项式，它不仅与这些数据点完美吻合，而且其形式尽可能简单。为此，我们采用拉格朗日插值法。下面是根据这些点构建的拉格朗日基础多项式：

对于 $x_0 = 0$ 的点，

L_0(x) = \frac{(x-x_1)(x-x_2)}{(x_0-x_1)(x_0-x_2)} = 2(x-\frac{1}{2})(x-1)

对于 $x_1 = \frac{1}{2}$ 的点，

L_1(x) = -4x(x-1)

对于 $x_2 = 1$ 的点，

L_2(x) = 2x(x-\frac{1}{2})

这些基础多项式帮助我们构建了一个多项式 $p(x)$ ，它精确地穿过所有给定的数据点：

p(x) = -2(x-\frac{1}{2})(x-1) - \frac{2}{3}x(x-1) + \frac{16}{9}x(x-\frac{1}{2})

这就是我们通过给定数据点构建的唯一的二次多项式。

为了直观展示，我们在图中展示了这些基础多项式以及它们构成的插值多项式。

函数 $f(x)$ 和插值多项式 $p(x)$ （虚线）如下所示。

结果表明，尽管我们只有三个数据点，但构建的插值多项式与原始函数 $f(x) = x - 9^{-x}$ 非常接近，这显示了插值多项式在数据恢复方面的强大能力

1.5 牛顿插值多项式

牛顿插值多项式的核心思想是，我们可以一步步地构建这个多项式，每次只增加一个数据点。想象一下，你正在画一条曲线，每次只根据一个新的点来调整曲线的形状。

首先，我们定义 $p_k$ 为一个特殊的多项式，它能准确通过前 $k+1$ 个点，也就是 $(x_0, y_0), \ldots, (x_k, y_k)$ 。

当 $k=0$ ，也就是只有一个点时，事情非常简单， $p_0$ 就是一个常数，值等于那个点的 $y$ 值：

p_0(x) = y_0.

假设我们已经知道 $p_{k-1}$ ，也就是通过前 $k$ 个点的多项式。现在我们要加入第 $k+1$ 个点 $(x_k, y_k)$ ，并稍微调整 $p_{k-1}$ ，让它同时通过这个新点，并且多项式的次数增加一次。我们可以把新的多项式 $p_k$ 写成这样：

p_k = p_{k-1} + c_k(x-x_0)\cdots(x-x_{k-1}).

这里的 $c_k$ 是一个常数。这个新加的项不会影响之前的点，因为对于任意 $i < k$ ， $(x-x_0)\cdots(x-x_{k-1})$ 在 $x_i$ 处的值都是0。所以我们只需要找一个合适的 $c_k$ ，让 $p_k(x_k) = y_k$ 就可以了。

通过这种方法，我们就可以逐步构建出牛顿插值多项式：

p_k = c_0 + c_1(x-x_0) + c_2(x-x_0)(x-x_1) + \cdots + c_k(x-x_0)\cdots(x-x_{k-1})

或者用一个更简洁的表示方式：

p_k = \sum_{j=0}^k c_j \prod_{i=0}^{j-1}(x-x_i)

这里的 $\prod_{i=0}^{-1}\cdots$ 被认为是1。

值得注意的是，每个系数 $c_j$ 只和前 $j+1$ 个点有关。这种表示方式不仅计算效率高，而且我们还有一些巧妙的方法来计算这些系数 $c_j$ ，这些方法既快速又简单！

1.6 差分牛顿插值

我们有时需要用一种特别的方法来理解和处理数据点之间的关系，这种方法叫做差分。想象一下，你有一串珠子，差分就像是观察并记录每两颗珠子之间的距离。

首先，对于单个数据点 $(x_i, y_i)$ ，我们用 $f[x_i]$ 表示它的值，就是 $y_i$ 。

接着，如果我们有两个相邻的数据点 $(x_{i-1}, y_{i-1})$ 和 $(x_i, y_i)$ ，我们可以计算这两个点之间的"差分"：

f[x_{i-1}, x_i] = \frac{f[x_i] - f[x_{i-1}]}{x_i - x_{i-1}}

这其实就是在计算这两个点之间的变化率。

当我们有更多的点时，我们可以用同样的方法，逐步扩展差分的计算。对于一系列的点 $x_i, x_{i+1}, \ldots, x_j$ ，差分的计算可以这样表示：

f[x_i, x_{i+1}, \ldots, x_{j-1}, x_j] = \frac{f[x_{i+1}, \ldots, x_{j-1}, x_j] - f[x_i, x_{i+1}, \ldots, x_{j-1}]}{x_j - x_i}

为了让这些计算看起来不那么复杂，我们用 $\gamma_{j\ell}$ 来表示这些差分，其中 $\ell$ 是距离，或者说是参与计算的点的数量。然后，我们可以用一个简单的公式来计算这些 $\gamma_{j\ell}$ ：

\gamma_{j0} = y_j, \quad \gamma_{j\ell} = \frac{\gamma_{j, \ell-1} - \gamma_{j-1, \ell-1}}{x_j - x_{j-\ell}}

这个方法可以帮助我们构建牛顿插值多项式。这种多项式可以精确地通过我们给定的点 $(x_0, y_0), \ldots, (x_n, y_n)$ 。牛顿插值多项式可以用这样的公式表示：

p_n(x) = f[x_0] + \sum_{j=1}^n f[x_0, x_1, \ldots, x_j] \prod_{i=0}^{j-1}(x - x_i)

虽然这个公式看起来有点复杂，但它实际上提供了一个非常高效的方式来计算和评估插值多项式。

让我们用一个具体的例子来展示这个方法。假设我们有一个函数 $f(x) = x^3 - 2x^2 + 1$ ，我们关注的点是 $x_i = 0, 1, 2, 3$ 。我们可以用一个表格来组织和计算差分，这个表格的每个格子都是基于它左边和上方的格子来计算的。这个过程可以用下面的表格表示：

\begin{array}{c|c|c|c|c} i & x_i & \gamma_{j0} & \gamma_{j1} & \gamma_{j2} & \gamma_{j3} \\ \hline 0 & 0 & \color{red}{f[x_0]} & & & \\ 1 & 1 & f[x_1] & \color{red}{f[x_0,x_1]} & & \\ 2 & 2 & f[x_2] & f[x_1,x_2] & \color{red}{f[x_0,x_1,x_2]} & \\ 3 & 3 & f[x_3] & f[x_2,x_3] & f[x_1,x_2,x_3] & \color{red}{f[x_0,x_1,x_2,x_3]} \\ \end{array}

通过这种方式，我们可以非常直观地看到每个差分是如何计算的。例如，在这个特定的例子中，我们可以计算出：

\gamma_{11} = f[x_0, x_1] = -1, \quad \gamma_{22} = f[x_0, x_1, x_2] = 1, \quad \gamma_{33} = f[x_0, x_1, x_2, x_3] = 1

最后，我们可以根据这些差分构建我们的牛顿插值多项式 $p_3(x)$ ：

p_3(x) = 1 - x + x(x - 1) + x(x - 1)(x - 2)

这个构建的多项式 $p_3(x)$ 和原始函数 $f(x)$ 实际上是完全一样的。这不仅展示了差分方法的准确性，还说明了它在数据分析和数学建模中的强大应用。

1.7 评估插值多项式的效率

在实际应用中，我们通常会先构造一个插值多项式，然后用它来估算大量点的值。这个过程中，如何高效地计算多项式在这些点上的值变得尤为重要。

牛顿形式的插值多项式提供了一种高效的评估方法，这种方法称为霍纳方法。这个方法的思想是从多项式的最高次项开始，逐步向下计算。牛顿形式的多项式可以这样表示：

p_n(x) = c_0 + (x - x_0) \left( c_1 + (x - x_1) \left( c_2 + \cdots + (x - x_{n-1})c_n \right) \cdots \right)

这里的 $c_j = f[x_0, x_1, \cdots, x_j]$ 是之前计算出的差分值。

基于这个形式，我们可以设计一个简单而有效的算法来计算 $p_n(x)$ 。这个算法的输入包括想要评估的点 $x$ ，节点/值对 $(x_0, y_0), \ldots, (x_n, y_n)$ ，以及系数 $(c_0, \ldots, c_n)$ 。算法的步骤如下：

初始化 $y$ 为最高次系数 $c_n$ 。
从 $n-1$ 递减到 0，依次进行如下计算：
$y \gets c_k + (x - x_k) \cdot y$
返回最终计算出的 $y$ 值。

这个算法大约需要 $3n$ 次浮点运算，这几乎是插值评估过程中可能达到的最高效率。尽管拉格朗日形式的插值多项式也能以类似数量的操作进行评估，但牛顿形式在实现上更为简洁直接。

C++示例代码:

#include <iostream>
#include <vector>

//x 一个包含多项式节点的向量。
//c 一个包含差分系数的向量。
//evalPoint 希望评估多项式的点。
double evaluateNewtonPolynomial(const std::vector<double>& x, 
                                const std::vector<double>& c, 
                                double evalPoint) {
    // 获取系数的数量，即n+1
    int n = c.size() - 1;
    
    // 初始化y为最高次系数
    double y = c[n];
    
    // 使用霍纳方法从n-1递减到0计算多项式的值
    for (int k = n - 1; k >= 0; k--) {
        y = c[k] + (evalPoint - x[k]) * y;
    }
    
    // 返回计算结果
    return y;
}

int main() {
    // 示例：节点/值对
    std::vector<double> x = {1, 2, 3};  // x0, x1, x2
    std::vector<double> c = {2, -1, 3}; // c0, c1, c2

    // 计算在x=1.5处的插值多项式的值
    double evalPoint = 1.5;
    double result = evaluateNewtonPolynomial(x, c, evalPoint);
    
    std::cout << "The value of the polynomial at x = " << evalPoint << " is: " << result << std::endl;
    
    return 0;
}

1.8 差分的泛化应用

在数学中，差分的概念可以被广泛地推广和应用。核心思想是，只要每个点 $x_i$ 都是唯一的，我们前面讨论的方法都是适用的。

通常情况下，我们考虑的这些点 $x_i$ 是按照顺序排列的，也就是说，它们是递增或递减的。但实际上，这些点的排列顺序并不重要！我们的方法同样适用于任意顺序排列的点。这个特性特别有用，比如说，如果我们已经有了一些点（比如 0, 0.5, 1），然后突然需要在中间插入一个新的点（比如 0.75），我们完全可以做到这一点，而不需要重新开始。
差分的概念也可以被扩展到包含重复的点。这种情况通常出现在我们想要插值某个函数的导数时。例如，在赫尔米特插值（Hermite interpolation）中，就需要考虑这种情况。
当点 $x_i$ 之间的间隔是相等的，也就是说它们是等间距排列的，我们可以推导出一些特别简洁的公式。在这种情况下，差分的计算和应用会变得更加直接和高效。

差分是一个非常强大的工具，它能够以多种方式适应和解决不同的数学问题。

2 多项式插值的误差分析

当我们使用多项式插值来估算一个函数 $f$ 的值时，一个自然的问题是：这个估算值有多准确？也就是说，我们构造的插值多项式 $p_n$ 在多大程度上接近真正的函数值 $f(x)$ ？这个问题的核心是探讨"插值误差"。

假设我们有一个插值多项式 $p_n$ ，它在若干个数据点 $(x_0, x_1, \ldots, x_n)$ 上与函数 $f$ 的值完全相同，即：

p_n(x_j) = f(x_j), \quad j=0,\ldots,n.

这些点都位于某个区间 $[a, b]$ 内，并且是递增排列的。

如果我们希望 $p_n$ 是 $f$ 的一个好的近似，那么我们关心的是，在某个子区间 $[c, d]$ （这个区间包含在 $[a, b]$ 内）上， $p_n(x)$ 与 $f(x)$ 的接近程度。我们用"插值误差"来描述这种接近程度，它可以用以下公式表示：

\max_{x \in [c, d]} |f(x) - p_n(x)|.

那么，我们面临的关键问题包括：

插值误差的界限是什么？
这个误差界限如何受到 $n$ （插值点的数量）和这些点的分布的影响？
误差界限在不同区间上如何变化？

找到这些问题的答案，可以帮助我们制定策略，在某个特定区间内通过插值来近似函数 $f$ 。值得注意的是，我们通常关注的区间是以插值点 $x_0$ 和 $x_n$ 为端点的区间。

当 $x$ 落在这个区间之外时，用 $p_n(x)$ 来估计 $f(x)$ 被称为外推，而不是插值。我们将会发现，相比插值，外推的准确性通常会更差。因此，除非有迫切的需求，我们应该尽量避免外推。

2.1 插值多项式的误差公式

在探讨插值多项式的误差公式之前，把它和泰勒定理进行对比是非常有帮助的。在泰勒定理的情况下，我们有一个函数 $f$ ，它在某个点 $x_0$ 及其附近是光滑的，我们知道在 $x_0$ 处的函数值和前 $n$ 个导数值：

f(x_0), f'(x_0), \ldots, f^{(n)}(x_0)

基于这些信息，我们可以构造一个多项式 $T_n(x)$ 来近似 $f$ ：

f(x) = T_n(x) + R_n(x)

其中， $T_n(x)$ 是泰勒多项式，可以表示为：

T_n(x) = f(x_0) + f'(x_0)(x-x_0) + \frac{f''(x_0)}{2!}(x-x_0)^2 + \cdots + \frac{f^{(n)}(x_0)}{n!}(x-x_0)^n

而 $R_n(x)$ 是泰勒公式的余项，它表示了多项式近似和真实函数值之间的差异：

R_n(x) = \frac{f^{(n+1)}(\xi_x)}{(n+1)!}(x-x_0)^{n+1}

这里 $\xi_x$ 是 $x_0$ 和 $x$ 之间的某个数。

对于插值多项式，我们有 $n+1$ 个值：

f(x_0), \ldots, f(x_n)

在 $n+1$ 个不同的点上。但与泰勒定理的主要结果相似，我们有拉格朗日误差公式：

定理（拉格朗日误差公式）：假设 $f \in C^{n+1}([a,b])$ ，并且有 $n+1$ 个节点：

x_0 < x_1 < \cdots < x_n

都包含在区间 $[a, b]$ 内。设 $p_n(x)$ 是插值多项式。那么对于每个 $x \in [a, b]$ ，都有：

f(x) = p_n(x) + E(x)

其中误差项 $E(x)$ 为：

E(x) = \frac{f^{(n+1)}(\xi_x)}{(n+1)!} \prod_{j=0}^n (x - x_j)

这里的 $\xi_x$ 是一个依赖于 $x$ 的数，位于区间 $[a, b]$ 内。

注意这个误差公式在形式上和泰勒定理非常相似，除了多项式的形式是 $(x-x_0)\cdots(x-x_n)$ 而不是 $(x-x_0)^{n+1}$ 。它的证明本质上和泰勒定理的证明是一致的。

值得指出的是，这个定理不仅适用于插值（ $x$ 在 $[x_0, x_n]$ 内），也适用于外推（ $x$ 在 $[x_0, x_n]$ 之外）。但要注意，外推时的误差通常会比插值时更大，所以应该谨慎使用。

2.2 插值误差的影响因素

在讨论多项式插值的误差时，我们发现误差的大小通常由三个主要因素决定：

$(n+1)$ 阶导数的最大值：这个因素反映了函数 $f$ 的曲率有多大。通常，我们用 $M$ 来表示函数在区间 $[a, b]$ 上 $(n+1)$ 阶导数的最大绝对值：
$M = \max_{x \in [a, b]} |f^{(n+1)}(x)|.$
由于 $\xi_x$ 的确切值未知，但知道它位于 $[a, b]$ 内，所以这个因素被 $M$ 限制。
节点间距离的影响：这个因素取决于插值节点 $x_j$ 与待估计点 $x$ 之间的距离。如果这些节点相互靠近，且靠近 $x$ ，这个乘积通常会较小；但如果节点间相隔较远，或者 $x$ 离节点群较远（尤其是在外推情况下），这个乘积可能会较大。
因子 $(1/(n+1)!)$ ：这个因子随着 $n$ 的增加而减小，对误差有抑制作用。但与此同时，其他两个因素可能会随着 $n$ 的增加而增大。因此，误差是否随着 $n$ 的增长而减少，往往取决于这些竞争因素之间的相互作用。

基于这些因素，我们可以得到误差的一个上界，即在某个区间 $I$ 内，插值多项式 $p_n(x)$ 和函数 $f(x)$ 之间的最大差异：

\max_{x \in I} |p_n(x) - f(x)| = \max_{x \in I} |E_n(x)| \le \frac{M}{(n+1)!} \max_{x \in I} |\omega_{n+1}(x)|.

需要注意的是，即使我们知道 $x$ 位于某个小区间 $I$ 内，我们也不能假设 $(n+1)$ 阶导数的评估点 $\xi$ 一定位于 $I$ 内，因此我们在评估 $f^{(n+1)}$ 的上界时，仍然被限制在整个区间 $[a, b]$ 内。

另一方面， $\omega$ 的计算只需考虑区间 $I$ 内的情况。所以如果我们只关心在某个小区间（例如 $[x_1, x_2]$ ）内的插值， $\omega$ 的计算就只需要考虑 $[x_1, x_2]$ 内的情况。但是，对于 $f^{(n+1)}$ 的最大值，我们仍然需要考虑整个 $[a, b]$ 区间。

2.3 一个插值误差的示例

以一个具体的例子来理解多项式插值的误差。
考虑函数 $f(x) = e^x$ ，我们知道该函数在点 $x_1 = 0$ 和 $x_2 = 1$ 的值，并想要在区间 [0,1] 内对它进行近似。对这两个点使用线性插值，我们得到插值多项式：

p_1 = 1 + (e - 1)x

对于这个插值多项式，误差的形式为：

E(x) = \frac{f''(\xi_x)}{2!}x(x - 1)

这里的 $\xi_x$ 是一个位于 [0,1] 内的未知数。 $f''(x)$ 的最大值可以通过考虑 $f''(x) = e^x$ 在 [0,1] 内的最大值来确定：

M = \max_{x \in [0,1]}|e^x| \le e

对于其他项，我们可以找到 $x(x - 1)$ 在 [0,1] 内的最大绝对值，这发生在 $x = 1/2$ 时：

\max_{x \in [0,1]}|x(x - 1)| \le \frac{1}{4}

因此，我们得到误差的一个上界：

|p_1(x) - f(x)| \le \frac{e}{8} \approx 0.34，对于 x 在 [0,1] 内

虽然这个近似不是非常理想，但考虑到 $p_1$ 只使用了两个函数值，这样的结果也在预期之中。然而，如果我们用同样的插值 $p_1$ 来近似 $f(x) = e^x$ 在更宽的区间 [0,2] 内，误差会更大。在这种情况下，我们得到一个更大的误差界限：

|p_1(x) - f(x)| \le 2e^2 \approx 15，对于 x 在 [0,1] 内。

现在，回到区间 [0,1]，假设我们使用 10 个等间距的点 $x_j = jh$ ，其中 $j=0, \ldots, 9$ ，并且 $h = 1/9$ 是点之间的间距。此时，误差表达式变为：

E(x) = \frac{f^{(10)}(\xi_x)}{10!} \prod_{j=0}^{9}(x - x_j)

对于 $f^{(10)}(\xi_x)$ ，我们知道：

|f^{(10)}(\xi_x)| \le e

对于乘积项，我们可以做一个粗略的估计，注意到对于所有的 $j$ ，有 $|x - x_j| \le 1$ 。因此，我们得到：

|E(x)| \le \frac{e}{10!} \approx 7.5 \times 10^{-7}

实际上，这个估计可能过于保守，因为一些节点 $x_j$ 可能比其他节点更接近 $x$ 。

2.4 一般误差界限：等间距点

假设我们有以下一系列等间距的点，间距为 $h$ ：

a = x_0 < x_1 < \cdots < x_n = b

x_{i}=a+jh,\quad h=\frac{b-a}{n}.

现在，我们关注一个位于区间 $[a, b]$ 内的特定点 $x$ 。当这个点 $x$ 接近区间的端点时，插值误差通常会达到最大。为什么呢？这是因为插值多项式在接近其定义范围的边缘时通常表现得不太稳定。

具体来说，假设 $x$ 落在 $[x_0, x_1]$ 内。在这种情况下，我们有：

|x - x_j| \le (j + 1)h \text{ 对于 } j = 0, \cdots, n

这意味着 $x$ 与每个点 $x_j$ 之间的距离最多是 $(j+1)h$ 。

因此，当我们确信 $x$ 与区间 $[x_1, b]$ 内所有点的距离不会超过一定值时，我们可以得出以下结论：

|\omega_{n+1}(x)| \le \prod_{j=0}^n |x - x_j| \le (n+1)!h^{n+1} \text{ 对于 } x \in [a, b]

基于这个结果，我们可以得出 $n$ 次插值多项式的误差 $E_n$ 的一个上界。但是，当我们增加更多的插值点时会发生什么呢？因子 $h^{n+1}$ 会随着 $h \to 0$ 而减小，但函数 $f$ 的 $(n+1)$ 阶导数可能会增大。是否增加插值点能改善误差，取决于 $(n+1)$ 阶导数增长的速度。

例如，考虑函数 $f(x) = 1/x$ ，其定义在区间 $[1, 2]$ 上。对于 $x \in [1, 2]$ ：

|f^{(n+1)}(x)| \le (n+1)! \left(\max_{x \in [1, 2]} |1/x|\right) = (n+1)!

在这种情况下， $h^{n+1}$ 的减小速度足以抵消 $(n+1)!$ 的增长：

|E_n(x)| \le (n+1)! \left(\frac{1}{n}\right)^{n+1} \sim \frac{\sqrt{2\pi n}}{e}e^{-n} \text{ 当 } n \to \infty

这里我们使用了斯特林近似：

n! \sim e^{-n}n^n\sqrt{2\pi n} \text{ 当 } n \to \infty

这种抵消效果是典型的，并且通常会导致误差的指数级衰减。然而，也存在一些极端情况，其中 $(n+1)$ 阶导数的增长过快，导致插值误差的急剧增加（龙格现象）。

2.5 龙格现象：插值多项式的一个典型问题

当我们尝试通过增加插值节点来改善插值多项式的精度时，情况并不总是像我们预期的那样。一个经典的例子是，对于看似稳定的函数：

\frac{1}{x^2+25}

我们在区间 [-1,1] 内使用等间距的节点：

-1 = x_0 < \cdots < x_n = 1

来构建插值多项式 $p_n(x)$ 。但是，当我们观察不同 $n$ 值下的插值多项式时，会发现一个有趣而又令人困惑的现象：多项式在端点附近出现剧烈振荡，而且随着 $n$ 的增加，这种振荡变得越来越明显。如图所示

实际上，插值多项式和函数之间的最大误差：

\max_{x \in [-1,1]} |p_n(x) - f(x)|

误差随着 $n \to \infty$ 呈指数级增长。由于多项式需要在每个节点处急转弯以符合函数值，导致其导数急剧增大，一旦出现大的导数，多项式就容易失控，难以收敛。

那么这种现象和误差公式有什么关联呢？关键在于这两个竞争因素：

\frac{f^{(n+1)}(\xi)}{(n+1)!} \quad \text{和} \quad \omega_{n+1}(x) = \prod_{j=0}^n (x - x_j)

可以验证，对于较大的 $n$ ， $f^{(n)} / (n+1)!$ 随 $n$ 增长得相当快，快到足以使得 $\omega_{n+1} / (n+1)!$ 的衰减（即使是指数级）也不足以抵消它。然而，通过选择不同的点使得 $\omega_n$ 更小，我们可以优化 $\omega_n$ 的大小，使得插值多项式效果更好。这正是切比雪夫插值的基础。对于 $n$ 个点，节点的选择是：

x_j = \cos\left(\frac{2j+1}{2n}\pi\right), \quad j = 0, \ldots, n-1

如图所示，对于 $n=10$ 和 $n=50$ ，插值误差现在表现良好！

所以，问题并不在于点数本身的增加，而在于点的分布。合理选择插值点可以显著改善插值多项式的性能。

关键点：这里的教训是，等间距的高阶插值多项式可能会导致灾难性的结果，如果可能，应该避免使用。然而，通过合理选择插值点，即使是高阶的插值多项式也可以发挥良好的效果。

更多背景：一般来说，哪些点集和函数会导致问题，这是一个复杂的问题。事实上，对于任何给定的选择点的方案，都有可能找到一个特定的函数来引发问题，即使是对于切比雪夫节点。然而，对于这些节点导致问题的函数，在实际应用中通常是极端的异常情况，这使得切比雪夫节点成为多数“合理”函数的良好选择。

在龙格现象的例子中，问题在于函数 $f$ 在复平面上有一个非常接近插值区间的奇点，这是等间距点对于函数 $f$ 失效的真正原因。

2.6 分段插值：高阶插值问题的另一种解决方案

对于高阶插值可能带来的问题，分段插值提供了一个完全不同的解决方案。与其在整个区间上使用一个高阶多项式，不如在每个小的子区间上使用低阶（通常是最高三阶）的多项式。这样的做法避免了之前提到的高阶插值中的问题。具体来说，我们有以下节点：

x_0 < x_1 < \cdots < x_n

对于每个子区间 $[x_i, x_{i+1}]$ ，我们可以构造一个线性（或更高阶，但通常不超过三阶）的插值多项式。这种方法的好处在于，增加节点数量并不会引起之前章节中提到的问题。根据拉格朗日误差公式，误差会随着间距 $h = x_{i+1} - x_i$ 趋向于零。

分段线性插值的一个缺点是，它在节点处不可微。换句话说，曲线在这些点上可能有尖锐的转折，这在某些应用中可能不是理想的特性。

对于需要函数值和导数都很平滑的场景，三次样条是一个常见的选择。这种方法的思想是，在每个子区间上使用三次插值多项式，并且利用以下数据：

f(x_i), f'(x_i), f(x_{i+1}), f'(x_{i+1})

这实际上是一种称为赫尔米特插值的方法，它不仅匹配函数值，还匹配导数。然后，我们将这些插值多项式“粘合”在一起，确保在节点处函数值和导数都连续，从而创建一个整体上连续且平滑的插值多项式。

在需要根据数据构建复杂函数表示的场景中，比如在计算机图形学中用于形状轮廓的绘制，三次样条是一种常见的选择，它能够提供既平滑又准确的插值结果。

2.7 插值多项式的证明与分析

在理解插值多项式时，了解其形式及其构建方法至关重要。我们已知，对于 $n+1$ 个特定的点 $(x_0, \ldots, x_n)$ ，存在一个唯一的多项式 $p_n(x)$ ，能够准确地在这些点上插值给定的函数 $f(x)$ ，即：

p_n(x_j) = f(x_j), \quad j=0, \ldots, n

此多项式可以用牛顿形式表示，而牛顿形式中的系数可以通过递归方法定义和计算：

f[x_i] = f(x_i)

f[x_i, x_{i+1}, \ldots, x_{j-1}, x_j] = \frac{f[x_{i+1}, \ldots, x_{j-1}, x_j] - f[x_i, x_{i+1}, \ldots, x_{j-1}]}{x_j - x_i}

在证明牛顿形式的有效性时，我们不是通过递归定义这些系数，而是直接将它们定义为多项式 $(\prod_{i=0}^{j-1}(x-x_i))$ 的系数，并且证明它们满足上述递归关系。

具体地，我们从常数多项式开始，即插值点 $(x_i, f(x_i))$ 的常数多项式，然后逐步构造出插值多项式 $p_{n-1}(x)$ 和 $p_n(x)$ 。我们声称， $p_n(x)$ 可以从 $p_{n-1}(x)$ 通过加上一个调整项得到：

p_n(x) = p_{n-1}(x) + \frac{x - x_0}{x_n - x_0} (q(x) - p_{n-1}(x))

这里的动机是让 $p_{n-1}(x)$ 和 $p_n(x)$ 在点 $x_0, \ldots, x_{n-1}$ 上一致，并确保它们在点 $x_n$ 上也一致。我们通过添加一个项 $g(x)(q(x)-p_{n-1}(x))$ 来实现这一点，其中 $g(x)$ 是一个在 $x_n$ 处为 1，而在 $x_0$ 处为 0 的函数，因此我们选择 $g(x) = \frac{x - x_0}{x_n - x_0}$ 。

通过这种方法，我们可以确保在所有插值点上 $p_n(x)$ 与 $f(x)$ 一致。由于两边都是不超过 $n$ 次的多项式，并且它们在 $n+1$ 个点上相同，所以它们必须是同一个多项式。这就完成了证明。

接下来，我们分析误差。当我们将 $x$ 视为另一个插值点时，我们可以构建一个在点 $x_0, \ldots, x_n, x$ 处插值 $f$ 的多项式，并取其 $n+1$ 阶导数。得到一个在这些点处插值 $f$ 的 $n+1$ 次多项式：

g(t) := p_n(t) + f[x_0, \ldots, x_n, x] \prod_{j=0}^n (t - x_j)

然后，我们考虑 $h(t) := f(t) - g(t)$ 。由于 $g(t)$ 在点 $x_0, \ldots, x_n, x$ 处插值 $f$ ， $h(t)$ 在这 $n+2$ 个点上有零点。根据罗尔定理，如果一个连续函数在区间 $[a, b]$ 内有 $k$ 个零点，那么它的导数在 $(a, b)$ 内至少有 $k-1$ 个零点。因此， $h^{(n+1)}$ 在区间 $[a, b]$ 内至少有一个零点，称为 $\xi_x$ 。因此我们得到：

f[x_0, \ldots, x_n, x] = \frac{f^{(n+1)}(\xi_x)}{(n+1)!}

最终，由于 $g(x)$ 在 $x$ 处插值 $f$ ，我们得到：

f(x) = p_n(x) + \frac{f^{(n+1)}(\xi_x)}{(n+1)!} \prod_{j=0}^n (x - x_j)

这就是误差公式，揭示了插值多项式和实际函数之间的差距。

2.8 赫米特插值

赫米特插值是一种特殊的插值方法，它不仅匹配了节点处的函数值，而且还匹配了节点处函数值的导数。这种方法的目标是构建一个多项式 $p$ ，使得对于每个节点 $x_j$ ， $j = 0, 1, 2, \ldots, n$ ，有：

p(x_j) = f(x_j), \quad p'(x_j) = f'(x_j)

赫米特插值可以扩展到每个点有任意数量的导数情况。关键定理表明，给定一个函数 $f$ 及其在各个 $x_j$ 处的导数，存在一个唯一的多项式 $p$ ，其度数不超过 $2n+1$ ，满足上述条件。

值得注意的是，这个多项式的度数是 $2n+1$ 。通常，这个度数比提供的条件数量少一，因为每个节点提供了两个条件，共 $2n+2$ 个条件。误差公式可以表示为：

f(x) = p(x) + \frac{f^{(2n+2)}(x)}{(2n+2)!} \prod_{j=0}^n (x - x_j)^2

这个公式与包含重复 $n+1$ 个点的 $2n+2$ 个点的拉格朗日公式相似。赫米特基础与拉格朗日基础类似，但构建起来相对繁琐。不过，赫米特插值可以通过差分的方法轻松计算。

我们定义一组 $2n+2$ 个点：

z_{2i} = z_{2i+1} = x_i, \quad 0 \le i \le n

也就是新的 $z$ 点是 $x_0, x_0, x_1, x_1, x_2, x_2, \ldots$ 这样的序列。我们按照之前的方式计算差分 $f[z_i, \ldots, z_{i+j}]$ ，除了在出现除以零的情况时，我们使用导数替换商：

f[z_{2i}, z_{2i+1}] = f'(x_i) \text{ 而不是 } \frac{f[z_{2i+1}] - f[z_{2i}]}{z_{2i+1} - z_{2i}}

这种改变只影响差分表的第一列；其余的差分计算和之前一样进行。结果是赫米特插值多项式。

来看一个实际的例子，设我们有以下数据和它的导数值：

f(-1) = 2, \quad f'(-1) = -1, \quad f(1) = 0, \quad f'(1) = 3

差分表如下，导数值用方框标记，对角线条目像以前一样用于 $p(x)$ ：

\begin{array}{c|c|ccccc}i&z_i&f[z_i]&f[z_i,z_{i-1}]&\cdots&\cdots\\0&-1&2\\1&-1&2&\boxed{\color{red}{-1}}\\2&1&0&-1&\color{red}{0}\\3&1&0&\boxed{3}&2&\color{red}{1}\\\end{array}

因此，我们得到赫米特插值多项式：

p(x) = 2 - (x + 1) + 0 \cdot (x + 1)^2 + 1 \cdot (x + 1)^2(x - 1)

通常，包含导数信息可以大幅提升插值的质量。虽然需要知道导数，但这通常是值得的。如果对于赫米特插值使用适当的点（例如 $n+1$ 个等间距点），那么赫米特插值多项式 $p_{2n+1}$ 确实能够很好地逼近原函数 $f$ 。

然而，即使赫米特插值能提供高精度的结果，高阶插值仍可能不是最佳选择。更常见的实践是分段方法，在点之间的区间内使用三次插值。这种方法在实际中广泛应用，特别是在需要连续导数的场合，如计算机图形学中的形状轮廓绘制。

2.9 三次样条插值

三次样条插值是一种用于通过一组散点构造平滑曲线的方法，不仅在这些点上曲线的函数值与给定数据相吻合，而且在整个区间上曲线的一阶和二阶导数都连续，从而保证曲线的平滑过渡。

三次样条插值的核心是在每个相邻数据点间构造一个三次多项式。假设数据点为 $(x_0, y_0), (x_1, y_1), \ldots, (x_n, y_n)$ ，对于相邻的数据点 $(x_i, y_i)$ 和 $(x_{i+1}, y_{i+1})$ ，我们构造三次多项式 $S_i(x)$ ，其中 $i = 0, 1, \ldots, n-1$ 。每个三次多项式 $S_i(x)$ 的一般形式为：

S_i(x) = a_i(x - x_i)^3 + b_i(x - x_i)^2 + c_i(x - x_i) + d_i

这些三次多项式需要满足以下条件：

函数值匹配： $S_i(x_i) = y_i$ 且 $S_i(x_{i+1}) = y_{i+1}$ 。
一阶导数连续： $S_{i-1}'(x_i) = S_i'(x_i)$ 。
二阶导数连续： $S_{i-1}''(x_i) = S_i''(x_i)$ 。

通常还会对两端的样条函数施加自然边界条件或其他边界条件，以确保整个区间内曲线的平滑性。

示例：

对于给定的三个数据点 $(1, 2)$ , $(2, 3)$ , $(3, 5)$ ，我们有两个三次多项式 $S_0(x)$ 和 $S_1(x)$ 。它们的一般形式为：

S_0(x) = a_0(x - 1)^3 + b_0(x - 1)^2 + c_0(x - 1) + d_0

S_1(x) = a_1(x - 2)^3 + b_1(x - 2)^2 + c_1(x - 2) + d_1

我们需要解决的线性方程组来自以下条件：

函数值条件：
- $S_0(1) = 2$ 即 $d_0 = 2$ 。
- $S_0(2) = S_1(2) = 3$ 即 $a_0 + b_0 + c_0 + d_0 = 3$ 和 $d_1 = 3$ 。
- $S_1(3) = 5$ 即 $a_1 + b_1 + c_1 + d_1 = 5$ 。
一阶导数连续：
- $S_0'(2) = S_1'(2)$ 即 $3a_0 + 2b_0 + c_0 = 3a_1 + c_1$ 。
二阶导数连续：
- $S_0''(2) = S_1''(2)$ 即 $6a_0 + 2b_0 = 6a_1 + 2b_1$ 。
自然边界条件（二阶导数为零）：
- $S_0''(1) = 0$ 即 $2b_0 = 0$ 。
- $S_1''(3) = 0$ 即 $6a_1 + 2b_1 = 0$ 。

将这些条件转换为线性方程组，我们得到：

\begin{aligned} d_0 &= 2 \\ a_0 + b_0 + c_0 + d_0 &= 3 \\ a_1 + b_1 + c_1 + d_1 &= 5 \\ 3a_0 + 2b_0 + c_0 &= 3a_1 + c_1 \\ 6a_0 + 2b_0 &= 6a_1 + 2b_1 \\ b_0 &= 0 \\ 6a_1 + 2b_1 &= 0 \\ \end{aligned}

通过求解这个方程组，我们可以得到 $a_0, b_0, c_0, d_0, a_1, b_1, c_1, d_1$ 的值，并完全确定 $S_0(x)$ 和 $S_1(x)$ 。

最终的样条函数是由这些三次多项式在各自的区间上的表达式组合而成的，它们在整个插值区间上形成一条平滑且连续的曲线。三次样条插值以其出色的平滑性和较高的准确性，在曲线构造、图形设计以及工程计算等领域中有着广泛的应用。它能够生成既符合数据点，又光滑连续的曲线，是一种在实际中非常实用的插值方法。