插值法概述
核心思想
插值 （Interpolation）就是在给定有限个点的函数值的情况下，寻找一个解析形式的插值函数 φ ( x ) \varphi(x) φ ( x ) f ( x ) f(x) f ( x )
代数插值问题的形式化定义
给定 n + 1 n+1 n + 1
x x 0 x 1 ⋯ x n y y 0 y 1 ⋯ y n \begin{array}{c|ccccc}
x & x_0 & x_1 & \cdots & x_n \\ \hline
y & y_0 & y_1 & \cdots & y_n
\end{array} x y x 0 y 0 x 1 y 1 ⋯ ⋯ x n y n 要构造插值多项式 φ n ( x ) \varphi_n(x) φ n ( x )
φ n ( x ) \varphi_n(x) φ n ( x ) n n n φ n ( x j ) = f ( x j ) = y j \varphi_n(x_j) = f(x_j) = y_j φ n ( x j ) = f ( x j ) = y j j = 0 , 1 , ⋯ , n j = 0, 1, \cdots, n j = 0 , 1 , ⋯ , n
插值多项式的存在唯一性
将 φ n ( x ) = a 0 + a 1 x + ⋯ + a n x n \varphi_n(x) = a_0 + a_1 x + \cdots + a_n x^n φ n ( x ) = a 0 + a 1 x + ⋯ + a n x n
∣ 1 x 0 x 0 2 ⋯ x 0 n 1 x 1 x 1 2 ⋯ x 1 n ⋮ ⋮ ⋮ ⋱ ⋮ 1 x n x n 2 ⋯ x n n ∣ = ∏ 0 ≤ j < i ≤ n ( x i − x j ) \begin{vmatrix}
1 & x_0 & x_0^2 & \cdots & x_0^n \\
1 & x_1 & x_1^2 & \cdots & x_1^n \\
\vdots & \vdots & \vdots & \ddots & \vdots \\
1 & x_n & x_n^2 & \cdots & x_n^n
\end{vmatrix}
= \prod_{0 \le j < i \le n}(x_i - x_j) 1 1 ⋮ 1 x 0 x 1 ⋮ x n x 0 2 x 1 2 ⋮ x n 2 ⋯ ⋯ ⋱ ⋯ x 0 n x 1 n ⋮ x n n = 0 ≤ j < i ≤ n ∏ ( x i − x j ) 故当 x 0 , x 1 , ⋯ , x n x_0, x_1, \cdots, x_n x 0 , x 1 , ⋯ , x n 两两不相等 时，Vandermonde 行列式不为零，方程组有唯一解，即插值多项式存在且唯一 。
Taylor 展开 vs 插值展开
插值余项与 Taylor 余项在结构上高度相似，区别仅在于单点展开 还是多点展开 ：
对比维度 Taylor 展开 插值展开 展开基点 单点 x 0 x_0 x 0 多点 x 0 , x 1 , … , x n x_0, x_1, \dots, x_n x 0 , x 1 , … , x n 余项/误差 R n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ( x − x 0 ) n + 1 R_n(x) = \dfrac{f^{(n+1)}(\xi)}{(n+1)!}(x-x_0)^{n+1} R n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ( x − x 0 ) n + 1 R n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ∏ i = 0 n ( x − x i ) R_n(x) = \dfrac{f^{(n+1)}(\xi)}{(n+1)!}\displaystyle\prod_{i=0}^{n}(x-x_i) R n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) i = 0 ∏ n ( x − x i ) 核心区别 ( x − x 0 ) n + 1 (x-x_0)^{n+1} ( x − x 0 ) n + 1 ω n + 1 ( x ) = ∏ i = 0 n ( x − x i ) \omega_{n+1}(x) = \displaystyle\prod_{i=0}^{n}(x-x_i) ω n + 1 ( x ) = i = 0 ∏ n ( x − x i ) 所需信息 需要 f f f x 0 x_0 x 0 只需要 f f f 近似性质 局部 近似：在 x 0 x_0 x 0 全局 近似：在节点处精确相等，在区间内整体拟合
插值余项
定义
插值多项式 φ n ( x ) \varphi_n(x) φ n ( x ) f ( x ) f(x) f ( x )
R n ( x ) = f ( x ) − φ n ( x ) R_n(x) = f(x) - \varphi_n(x) R n ( x ) = f ( x ) − φ n ( x ) 称为截断误差 ，又称插值余项 。
两种表示形式
对同一插值问题，余项有两种等价写法，分别对应 Lagrange 与 Newton 视角：
形式 表达式 所需条件 特点 Lagrange 型（导数形式） R n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ω n + 1 ( x ) R_n(x)=\dfrac{f^{(n+1)}(\xi)}{(n+1)!}\,\omega_{n+1}(x) R n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ω n + 1 ( x ) f ∈ C n + 1 [ a , b ] f\in C^{n+1}[a,b] f ∈ C n + 1 [ a , b ] 依赖光滑性，便于误差估计 Newton 型（差商形式） R n ( x ) = f [ x 0 , … , x n , x ] ω n + 1 ( x ) R_n(x)=f[x_0,\dots,x_n,x]\,\omega_{n+1}(x) R n ( x ) = f [ x 0 , … , x n , x ] ω n + 1 ( x ) 仅需函数值存在 无需可导假设
其中 ω n + 1 ( x ) = ∏ i = 0 n ( x − x i ) \omega_{n+1}(x)=\displaystyle\prod_{i=0}^{n}(x-x_i) ω n + 1 ( x ) = i = 0 ∏ n ( x − x i )
两种余项的等价性
将 差商与导数的关系 定理应用于 Newton 余项中的 n + 1 n+1 n + 1 f [ x 0 , … , x n , x ] f[x_0,\dots,x_n,x] f [ x 0 , … , x n , x ]
f [ x 0 , … , x n , x ] = f ( n + 1 ) ( ξ ) ( n + 1 ) ! f[x_0,\dots,x_n,x]=\frac{f^{(n+1)}(\xi)}{(n+1)!} f [ x 0 , … , x n , x ] = ( n + 1 )! f ( n + 1 ) ( ξ ) 代回 Newton 余项，即恢复 Lagrange 型余项：
R n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ⋅ ω n + 1 ( x ) R_n(x)=\frac{f^{(n+1)}(\xi)}{(n+1)!}\cdot\omega_{n+1}(x) R n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ⋅ ω n + 1 ( x ) Lagrange 插值
一般形式
对 n + 1 n+1 n + 1 第 i i i 为：
l i ( x ) = ∏ j = 0 j ≠ i n x − x j x i − x j l_i(x) = \prod_{\substack{j=0 \\ j \ne i}}^{n} \frac{x - x_j}{x_i - x_j} l i ( x ) = j = 0 j = i ∏ n x i − x j x − x j 插值基函数次数恰好为 n n n
l i ( x j ) = δ i j = { 1 , i = j 0 , i ≠ j l_i(x_j) = \delta_{ij} = \begin{cases} 1, & i = j \\ 0, & i \ne j \end{cases} l i ( x j ) = δ ij = { 1 , 0 , i = j i = j n n n
φ n ( x ) = ∑ i = 0 n y i l i ( x ) \varphi_n(x) = \sum_{i=0}^{n} y_i\, l_i(x) φ n ( x ) = i = 0 ∑ n y i l i ( x ) 显然 φ n ( x j ) = ∑ i = 0 n y i δ i j = y j \varphi_n(x_j) = \sum_{i=0}^{n} y_i\, \delta_{ij} = y_j φ n ( x j ) = ∑ i = 0 n y i δ ij = y j
Lagrange 插值余项
设 x 0 , x 1 , ⋯ , x n x_0, x_1, \cdots, x_n x 0 , x 1 , ⋯ , x n [ a , b ] [a,b] [ a , b ] φ n ( x ) \varphi_n(x) φ n ( x ) n n n f ( x ) ∈ C n + 1 [ a , b ] f(x) \in C^{n+1}[a,b] f ( x ) ∈ C n + 1 [ a , b ] [ a , b ] [a,b] [ a , b ] x x x
R n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ω n + 1 ( x ) , ξ ∈ ( min { x i , x } , max { x i , x } ) R_n(x) = \frac{f^{(n+1)}(\xi)}{(n+1)!}\,\omega_{n+1}(x), \quad \xi \in \bigl(\min\{x_i,x\},\max\{x_i,x\}\bigr) R n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ω n + 1 ( x ) , ξ ∈ ( min { x i , x } , max { x i , x } ) 其中 ω n + 1 ( x ) = ∏ i = 0 n ( x − x i ) \omega_{n+1}(x) = \displaystyle\prod_{i=0}^{n}(x - x_i) ω n + 1 ( x ) = i = 0 ∏ n ( x − x i ) x x x ξ \xi ξ
证明 ：
选取一个非节点的求值点 x x x
ψ ( t ) = f ( t ) − φ n ( t ) − K ω n + 1 ( t ) \psi(t) = f(t) - \varphi_n(t) - K\,\omega_{n+1}(t) ψ ( t ) = f ( t ) − φ n ( t ) − K ω n + 1 ( t ) 其中常数 K K K
K = R n ( x ) ω n + 1 ( x ) = f ( x ) − φ n ( x ) ω n + 1 ( x ) K=\frac{R_n(x)}{\omega_{n+1}(x)}=\frac{f(x)-\varphi_n(x)}{\omega_{n+1}(x)} K = ω n + 1 ( x ) R n ( x ) = ω n + 1 ( x ) f ( x ) − φ n ( x ) 使得 ψ ( x ) = 0 \psi(x)=0 ψ ( x ) = 0
步骤 1 ：由插值条件 φ n ( x i ) = f ( x i ) \varphi_n(x_i)=f(x_i) φ n ( x i ) = f ( x i ) ω n + 1 ( x i ) = 0 \omega_{n+1}(x_i)=0 ω n + 1 ( x i ) = 0
ψ ( x i ) = 0 , i = 0 , 1 , … , n \psi(x_i)=0,\qquad i=0,1,\dots,n ψ ( x i ) = 0 , i = 0 , 1 , … , n 又 ψ ( x ) = 0 \psi(x)=0 ψ ( x ) = 0 ψ ( t ) \psi(t) ψ ( t ) [ min { x 0 , … , x n , x } , max { x 0 , … , x n , x } ] \bigl[\min\{x_0,\dots,x_n,x\},\max\{x_0,\dots,x_n,x\}\bigr] [ min { x 0 , … , x n , x } , max { x 0 , … , x n , x } ] n + 2 n+2 n + 2
步骤 2 ：反复应用 Rolle 中值定理 ：
ψ ( t ) \psi(t) ψ ( t ) n + 2 n+2 n + 2 ⇒ \Rightarrow ⇒ ψ ′ ( t ) \psi'(t) ψ ′ ( t ) n + 1 n+1 n + 1 ⋯ \cdots ⋯ ψ ( n + 1 ) ( t ) \psi^{(n+1)}(t) ψ ( n + 1 ) ( t ) 1 1 1 ξ \xi ξ ψ ( n + 1 ) ( ξ ) = 0 \psi^{(n+1)}(\xi)=0 ψ ( n + 1 ) ( ξ ) = 0
步骤 3 ：计算 ( n + 1 ) (n+1) ( n + 1 )
φ n ( n + 1 ) ( t ) ≡ 0 \varphi_n^{(n+1)}(t)\equiv 0 φ n ( n + 1 ) ( t ) ≡ 0 ω n + 1 ( n + 1 ) ( t ) ≡ ( n + 1 ) ! \omega_{n+1}^{(n+1)}(t)\equiv (n+1)! ω n + 1 ( n + 1 ) ( t ) ≡ ( n + 1 )!
代入 ψ ( n + 1 ) ( ξ ) = 0 \psi^{(n+1)}(\xi)=0 ψ ( n + 1 ) ( ξ ) = 0
0 = f ( n + 1 ) ( ξ ) − K ⋅ ( n + 1 ) ! ⟹ K = f ( n + 1 ) ( ξ ) ( n + 1 ) ! 0 = f^{(n+1)}(\xi) - K\cdot(n+1)!
\quad\Longrightarrow\quad
K=\frac{f^{(n+1)}(\xi)}{(n+1)!} 0 = f ( n + 1 ) ( ξ ) − K ⋅ ( n + 1 )! ⟹ K = ( n + 1 )! f ( n + 1 ) ( ξ ) 因此
R n ( x ) = K ⋅ ω n + 1 ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ω n + 1 ( x ) R_n(x)=K\cdot\omega_{n+1}(x)=\frac{f^{(n+1)}(\xi)}{(n+1)!}\,\omega_{n+1}(x) R n ( x ) = K ⋅ ω n + 1 ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ω n + 1 ( x ) 注 ：对比 Taylor 余项 R n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ( x − x 0 ) n + 1 R_n(x)=\dfrac{f^{(n+1)}(\xi)}{(n+1)!}(x-x_0)^{n+1} R n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ( x − x 0 ) n + 1 ( x − x 0 ) n + 1 (x-x_0)^{n+1} ( x − x 0 ) n + 1 ω n + 1 ( x ) \omega_{n+1}(x) ω n + 1 ( x )
Lagrange 插值算法实现
Lagrange 插值单插值点
def lagrange_single (x, y, xt): n = len (x) - 1 phi = 0.0 for i in range (n + 1 ): # 计算第 i 个基函数值 li = 1.0 for j in range (n + 1 ): if j != i: li *= (xt - x[j]) / (x[i] - x[j]) # 累加插值结果 phi += y[i] * li return phi 复杂度为 O ( n 2 ) \mathcal{O}(n^2) O ( n 2 )
Lagrange 插值多插值点
利用 ω n + 1 ′ ( x i ) = ∏ j ≠ i ( x i − x j ) \omega_{n+1}'(x_i) = \displaystyle\prod_{j \ne i}(x_i - x_j) ω n + 1 ′ ( x i ) = j = i ∏ ( x i − x j )
l i ( x ) = ω n + 1 ( x ) / ( x − x i ) ω n + 1 ′ ( x i ) φ n ( x ) = ω n + 1 ( x ) ∑ i = 0 n y i ( x − x i ) ω n + 1 ′ ( x i ) \begin{aligned}
l_i(x) &= \frac{\omega_{n+1}(x)/(x - x_i)}{\omega_{n+1}'(x_i)} \\
\varphi_n(x) &= \omega_{n+1}(x) \sum_{i=0}^{n} \frac{y_i}{(x-x_i)\,\omega_{n+1}'(x_i)}
\end{aligned} l i ( x ) φ n ( x ) = ω n + 1 ′ ( x i ) ω n + 1 ( x ) / ( x − x i ) = ω n + 1 ( x ) i = 0 ∑ n ( x − x i ) ω n + 1 ′ ( x i ) y i def lagrange_multi (x, y, xt_list): n = len (x) - 1 # 预计算权重 1/omega'(x_i)，复杂度 O(n^2) wt = [ 1.0 ] * (n + 1 ) for i in range (n + 1 ): for j in range (n + 1 ): if j != i: wt[i] *= (x[i] - x[j]) # 先求积 for i in range (n + 1 ): wt[i] = 1.0 / wt[i] results = [] for xt in xt_list: # 若 xt 恰好等于某个 x_i，直接返回 y_i 避免出现 0/0 flag = False for i in range (n + 1 ): if abs (xt - x[i]) < EPS : results.append(y[i]) flag = True break if flag: continue # 计算 omega_{n+1}(xt)，复杂度 O(n) omega = 1.0 for i in range (n + 1 ): omega *= (xt - x[i]) # 累加，复杂度 O(n) phi = 0.0 for i in range (n + 1 ): phi += y[i] * omega / (xt - x[i]) * wt[i] results.append(phi) return results 复杂度为 O ( n 2 + m n ) \mathcal{O}(n^2 + mn) O ( n 2 + mn ) m m m
Lagrange 插值向量化版
import numpy as np def lagrange_vectorized (x, y, xt_list): # 计算权重 wt[i] = 1 / prod_{j!=i}(x_i - x_j) diff = x[:, None ] - x[ None , :] # (n+1, n+1) 差矩阵 np.fill_diagonal(diff, 1.0 ) # 对角线置 1 避免除零 wt = 1.0 / np.prod(diff, axis = 1 ) # 每行连乘 phi = np.zeros( len (xt)) # 解决 0/0 问题 hits = np.abs(xt[:, None ] - x[ None , :]) < EPS # (m, n+1) hit_mask = np.any(hits, axis = 1 ) # (m,) if np.any(hit_mask): hit_idx = np.argmax(hits, axis = 1 ) # 命中的节点下标 phi[hit_mask] = y[hit_idx[hit_mask]] # 直接返回对应 y 值 # 计算其余位置 xt_valid = xt[ ~ hit_mask] if len (xt_valid) > 0 : denom = xt_valid[:, None ] - x[ None , :] # (m', n+1) omega = np.prod(denom, axis = 1 ) # (m',) phi[ ~ hit_mask] = omega * np.sum(y * wt / denom, axis = 1 ) return phi 复杂度分析
指标 单插值点 m m m 时间复杂度 O ( n 2 ) O(n^2) O ( n 2 ) O ( n 2 + m n ) O(n^2 + mn) O ( n 2 + mn ) 空间复杂度 O ( 1 ) O(1) O ( 1 ) O ( n ) O(n) O ( n ) 新增一点时 全部基函数重算 O ( n 2 ) O(n^2) O ( n 2 ) 全部权重重算 O ( n 2 ) O(n^2) O ( n 2 )
注 ：Lagrange 插值的主要缺点是：每新增一个插值节点，所有基函数均需重新计算 ，无法复用之前的结果。Newton 插值解决了这一问题。
Newton 插值
均差（差商）的定义
零阶均差 ：
f [ x i ] = f ( x i ) f[x_i] = f(x_i) f [ x i ] = f ( x i ) 一阶均差 ：
f [ x i , x j ] = f ( x j ) − f ( x i ) x j − x i f[x_i, x_j] = \frac{f(x_j) - f(x_i)}{x_j - x_i} f [ x i , x j ] = x j − x i f ( x j ) − f ( x i ) n n n
f [ x 0 , x 1 , ⋯ , x n ] = f [ x 1 , x 2 , ⋯ , x n ] − f [ x 0 , x 1 , ⋯ , x n − 1 ] x n − x 0 = ∑ k = 0 n f ( x k ) ∏ j = 0 j ≠ k n ( x k − x j ) = ∑ k = 0 n f ( x k ) ω n + 1 ′ ( x k ) \begin{aligned}
f[x_0, x_1, \cdots, x_n]
&= \frac{f[x_1, x_2, \cdots, x_n] - f[x_0, x_1, \cdots, x_{n-1}]}{x_n - x_0} \\
&= \sum_{k=0}^n \frac{f(x_k)}{\prod\limits_{\substack{j=0 \\ j \ne k}}^n (x_k - x_j)} \\
&= \sum_{k=0}^n \frac{f(x_k)}{\omega_{n+1}'(x_k)} \\
\end{aligned} f [ x 0 , x 1 , ⋯ , x n ] = x n − x 0 f [ x 1 , x 2 , ⋯ , x n ] − f [ x 0 , x 1 , ⋯ , x n − 1 ] = k = 0 ∑ n j = 0 j = k ∏ n ( x k − x j ) f ( x k ) = k = 0 ∑ n ω n + 1 ′ ( x k ) f ( x k ) 均差的性质
线性性 ：
若 f ( x ) = α φ ( x ) + β ψ ( x ) f(x) = \alpha\varphi(x) + \beta\psi(x) f ( x ) = α φ ( x ) + β ψ ( x )
f [ x 0 , ⋯ , x n ] = α φ [ x 0 , ⋯ , x n ] + β ψ [ x 0 , ⋯ , x n ] f[x_0,\cdots,x_n] = \alpha\varphi[x_0,\cdots,x_n] + \beta\psi[x_0,\cdots,x_n] f [ x 0 , ⋯ , x n ] = α φ [ x 0 , ⋯ , x n ] + β ψ [ x 0 , ⋯ , x n ]
线性组合，系数和为 0 ：
f [ x 0 , x 1 , ⋯ , x n ] = ∑ i = 0 n f ( x i ) ω n + 1 ′ ( x i ) f[x_0, x_1, \cdots, x_n] = \sum_{i=0}^{n} \frac{f(x_i)}{\omega_{n+1}'(x_i)} f [ x 0 , x 1 , ⋯ , x n ] = i = 0 ∑ n ω n + 1 ′ ( x i ) f ( x i ) 显然 f [ x 0 , x 1 , ⋯ , x n ] f[x_0, x_1, \cdots, x_n] f [ x 0 , x 1 , ⋯ , x n ] f ( x 0 ) , ⋯ f ( x n ) f(x_0), \cdots f(x_n) f ( x 0 ) , ⋯ f ( x n ) ∑ i = 0 n 1 / ω n + 1 ′ ( x i ) = 0 \sum_{i=0}^{n} 1/\omega_{n+1}'(x_i) = 0 ∑ i = 0 n 1/ ω n + 1 ′ ( x i ) = 0
证明 ：
Lagrange 基函数之和为：
∑ i = 0 n l i ( x ) = ∑ i = 0 n ω n + 1 ( x ) ω n + 1 ′ ( x i ) ( x − x i ) = ∑ i = 0 n ∏ j ≠ i ( x − x j ) ω n + 1 ′ ( x i ) \sum_{i=0}^{n}l_i(x)
=\sum_{i=0}^{n}\frac{\omega_{n+1}(x)}{\omega_{n+1}'(x_i)(x-x_i)}=\sum_{i=0}^{n}\frac{\prod_{j\ne i}(x-x_j)}{\omega_{n+1}'(x_i)} i = 0 ∑ n l i ( x ) = i = 0 ∑ n ω n + 1 ′ ( x i ) ( x − x i ) ω n + 1 ( x ) = i = 0 ∑ n ω n + 1 ′ ( x i ) ∏ j = i ( x − x j ) 这是次数 ≤ n \le n ≤ n n + 1 n+1 n + 1 x k x_k x k 1 1 1 1 1 1 x n x^n x n
对称性 ：
f [ x i , x j ] = f [ x j , x i ] , f [ x i , x j , x k ] = f [ x j , x i , x k ] f[x_i, x_j] = f[x_j, x_i], \quad f[x_i, x_j, x_k] = f[x_j, x_i, x_k] f [ x i , x j ] = f [ x j , x i ] , f [ x i , x j , x k ] = f [ x j , x i , x k ]
多项式求导关系 ：
若 f ( x ) f(x) f ( x ) n n n f [ x 0 , ⋯ , x k , x ] f[x_0, \cdots, x_k, x] f [ x 0 , ⋯ , x k , x ] n − k − 1 n-k-1 n − k − 1
差商与导数的关系 ：
若 f ∈ C n [ a , b ] f\in C^{n}[a,b] f ∈ C n [ a , b ] x 0 , x 1 , … , x n x_0,x_1,\dots,x_n x 0 , x 1 , … , x n ξ \xi ξ
f [ x 0 , x 1 , … , x n ] = f ( n ) ( ξ ) n ! f[x_0,x_1,\dots,x_n]=\frac{f^{(n)}(\xi)}{n!} f [ x 0 , x 1 , … , x n ] = n ! f ( n ) ( ξ ) 证明 ：
设 P n − 1 ( t ) P_{n-1}(t) P n − 1 ( t ) f ( t ) f(t) f ( t ) x 0 , … , x n − 1 x_0,\dots,x_{n-1} x 0 , … , x n − 1 n − 1 n-1 n − 1 x n x_n x n
φ ( t ) = f ( t ) − P n − 1 ( t ) − K ⋅ ω n ( t ) , ω n ( t ) = ( t − x 0 ) ⋯ ( t − x n − 1 ) \varphi(t)=f(t)-P_{n-1}(t)-K\cdot\omega_n(t),\qquad \omega_n(t)=(t-x_0)\cdots(t-x_{n-1}) φ ( t ) = f ( t ) − P n − 1 ( t ) − K ⋅ ω n ( t ) , ω n ( t ) = ( t − x 0 ) ⋯ ( t − x n − 1 ) 取 K = f ( x n ) − P n − 1 ( x n ) ω n ( x n ) K=\dfrac{f(x_n)-P_{n-1}(x_n)}{\omega_n(x_n)} K = ω n ( x n ) f ( x n ) − P n − 1 ( x n ) φ ( x n ) = 0 \varphi(x_n)=0 φ ( x n ) = 0
由插值条件，φ ( x i ) = 0 ( i = 0 , … , n − 1 ) \varphi(x_i)=0\ (i=0,\dots,n-1) φ ( x i ) = 0 ( i = 0 , … , n − 1 ) φ \varphi φ n + 1 n+1 n + 1 φ ( n ) \varphi^{(n)} φ ( n ) ξ \xi ξ
计算得 φ ( n ) ( t ) = f ( n ) ( t ) − K ⋅ n ! \varphi^{(n)}(t)=f^{(n)}(t)-K\cdot n! φ ( n ) ( t ) = f ( n ) ( t ) − K ⋅ n ! ξ \xi ξ
0 = f ( n ) ( ξ ) − K ⋅ n ! ⟹ K = f ( n ) ( ξ ) n ! 0=f^{(n)}(\xi)-K\cdot n!
\quad\Longrightarrow\quad
K=\frac{f^{(n)}(\xi)}{n!} 0 = f ( n ) ( ξ ) − K ⋅ n ! ⟹ K = n ! f ( n ) ( ξ ) 又因为 K = f ( x n ) − P n − 1 ( x n ) ω n ( x n ) = f [ x 0 , … , x n − 1 , x n ] K=\dfrac{f(x_n)-P_{n-1}(x_n)}{\omega_n(x_n)}=f[x_0,\dots,x_{n-1},x_n] K = ω n ( x n ) f ( x n ) − P n − 1 ( x n ) = f [ x 0 , … , x n − 1 , x n ]
f [ x 0 , … , x n ] = f ( n ) ( ξ ) n ! f[x_0,\dots,x_n]=\frac{f^{(n)}(\xi)}{n!} f [ x 0 , … , x n ] = n ! f ( n ) ( ξ ) 推论 ：
当节点全部重合时，差商退化为导数：
f [ x 0 , … , x 0 ⏟ n + 1 个 ] = f ( n ) ( x 0 ) n ! f[\underbrace{x_0,\dots,x_0}_{n+1\text{个}}]=\frac{f^{(n)}(x_0)}{n!} f [ n + 1 个 x 0 , … , x 0 ] = n ! f ( n ) ( x 0 )
差商对活跃节点的求导 ：
若 f ∈ C n + 2 [ a , b ] f\in C^{n+2}[a,b] f ∈ C n + 2 [ a , b ] x 0 , … , x n x_0,\dots,x_n x 0 , … , x n x x x
d d x f [ x 0 , … , x n , x ] = f [ x 0 , … , x n , x , x ] \frac{\mathrm{d}}{\mathrm{d}x}f[x_0,\dots,x_n,x]=f[x_0,\dots,x_n,x,x] d x d f [ x 0 , … , x n , x ] = f [ x 0 , … , x n , x , x ] 证明 ：
记 g ( x ) = f [ x 0 , … , x n , x ] g(x)=f[x_0,\dots,x_n,x] g ( x ) = f [ x 0 , … , x n , x ] y ≠ x y\neq x y = x
f [ x 0 , … , x n , x , y ] = f [ x 0 , … , x n , y ] − f [ x 0 , … , x n , x ] y − x = g ( y ) − g ( x ) y − x f[x_0,\dots,x_n,x,y]=\frac{f[x_0,\dots,x_n,y]-f[x_0,\dots,x_n,x]}{y-x}=\frac{g(y)-g(x)}{y-x} f [ x 0 , … , x n , x , y ] = y − x f [ x 0 , … , x n , y ] − f [ x 0 , … , x n , x ] = y − x g ( y ) − g ( x ) 令 y → x y\to x y → x f ∈ C n + 2 f\in C^{n+2} f ∈ C n + 2 n + 2 n+2 n + 2 f [ x 0 , … , x n , x , x ] f[x_0,\dots,x_n,x,x] f [ x 0 , … , x n , x , x ] g ( x ) g(x) g ( x )
d d x f [ x 0 , … , x n , x ] = f [ x 0 , … , x n , x , x ] \frac{\mathrm{d}}{\mathrm{d}x}f[x_0,\dots,x_n,x]=f[x_0,\dots,x_n,x,x] d x d f [ x 0 , … , x n , x ] = f [ x 0 , … , x n , x , x ] 注 ：
结合性质 5，进一步可得：
f [ x 0 , … , x n , x , x ] = f ( n + 2 ) ( ξ ) ( n + 2 ) ! f[x_0,\dots,x_n,x,x]=\frac{f^{(n+2)}(\xi)}{(n+2)!} f [ x 0 , … , x n , x , x ] = ( n + 2 )! f ( n + 2 ) ( ξ ) 其中 ξ \xi ξ x 0 , … , x n , x x_0,\dots,x_n,x x 0 , … , x n , x
Newton 插值公式
定义Newton 基函数 ：
N 0 ( x ) ≡ 1 , N k ( x ) = ∏ j = 0 k − 1 ( x − x j ) ( k ≥ 1 ) N_0(x)\equiv 1,\qquad N_k(x)=\prod_{j=0}^{k-1}(x-x_j)\;(k\ge 1) N 0 ( x ) ≡ 1 , N k ( x ) = j = 0 ∏ k − 1 ( x − x j ) ( k ≥ 1 ) 显然任意小于等于 n n n φ n ( x ) \varphi_n(x) φ n ( x ) N 0 , ⋯ N n N_0, \cdots N_n N 0 , ⋯ N n
φ n ( x ) = ∑ k = 0 n a k N k ( x ) \varphi_n(x)=\displaystyle\sum_{k=0}^{n}a_k N_k(x) φ n ( x ) = k = 0 ∑ n a k N k ( x ) 则在插值条件 φ n ( x i ) = f ( x i ) ( i = 0 , 1 , … , n ) \varphi_n(x_i)=f(x_i)\;(i=0,1,\dots,n) φ n ( x i ) = f ( x i ) ( i = 0 , 1 , … , n )
a k = f [ x 0 , ⋯ , x k ] , k = 0 , 1 , … , n a_k = f[x_0, \cdots, x_k], \quad k=0,1,\dots,n a k = f [ x 0 , ⋯ , x k ] , k = 0 , 1 , … , n 递推关系：
φ n ( x ) = φ n − 1 ( x ) + N n ( x ) ⋅ f [ x 0 , x 1 , ⋯ , x n ] \varphi_n(x) = \varphi_{n-1}(x) + N_n(x) \cdot f[x_0,x_1,\cdots,x_n] φ n ( x ) = φ n − 1 ( x ) + N n ( x ) ⋅ f [ x 0 , x 1 , ⋯ , x n ] 证明思路 ：
设命题对 k ≤ n − 1 k \le n-1 k ≤ n − 1 n + 1 n+1 n + 1 φ n ( x n ) = f ( x n ) \varphi_n(x_n) = f(x_n) φ n ( x n ) = f ( x n )
φ n ( x n ) = φ n − 1 ( x n ) + ∏ j = 0 n − 1 ( x n − x j ) ⋅ f [ x 0 , ⋯ , x n − 2 , x n ] − f [ x 0 , ⋯ , x n − 2 , x n − 1 ] x n − x n − 1 = φ n − 1 ( x n ) + ∏ j = 0 n − 2 ( x n − x j ) ⋅ ( f [ x 0 , ⋯ , x n − 2 , x n ] − f [ x 0 , ⋯ , x n − 2 , x n − 1 ] ) = φ n − 1 ( x n ) + [ φ ~ n − 1 ( x n ) − φ n − 2 ( x n ) ] − [ φ n − 1 ( x n ) − φ n − 2 ( x n ) ] = φ ~ n − 1 ( x n ) = f ( x n ) \begin{aligned}
\varphi_n(x_n)
&= \varphi_{n-1}(x_n) + \prod_{j=0}^{n-1}(x_n-x_j)
\cdot\frac{f[x_0,\cdots,x_{n-2},x_n]-f[x_0,\cdots,x_{n-2},x_{n-1}]}{x_n-x_{n-1}} \\
&= \varphi_{n-1}(x_n) + \prod_{j=0}^{n-2}(x_n-x_j)
\cdot(f[x_0,\cdots,x_{n-2},x_n]-f[x_0,\cdots,x_{n-2},x_{n-1}]) \\
&= \varphi_{n-1}(x_n) + [\tilde\varphi_{n-1}(x_n) - \varphi_{n-2}(x_n)] - [\varphi_{n-1}(x_n) - \varphi_{n-2}(x_n)] \\
&= \tilde\varphi_{n-1}(x_n) \\
&= f(x_n)
\end{aligned} φ n ( x n ) = φ n − 1 ( x n ) + j = 0 ∏ n − 1 ( x n − x j ) ⋅ x n − x n − 1 f [ x 0 , ⋯ , x n − 2 , x n ] − f [ x 0 , ⋯ , x n − 2 , x n − 1 ] = φ n − 1 ( x n ) + j = 0 ∏ n − 2 ( x n − x j ) ⋅ ( f [ x 0 , ⋯ , x n − 2 , x n ] − f [ x 0 , ⋯ , x n − 2 , x n − 1 ]) = φ n − 1 ( x n ) + [ φ ~ n − 1 ( x n ) − φ n − 2 ( x n )] − [ φ n − 1 ( x n ) − φ n − 2 ( x n )] = φ ~ n − 1 ( x n ) = f ( x n ) 其中 φ ~ n − 1 ( x ) \tilde\varphi_{n-1}(x) φ ~ n − 1 ( x ) x 0 , ⋯ , x n − 2 , x n x_0, \cdots, x_{n-2}, x_n x 0 , ⋯ , x n − 2 , x n
再由插值多项式的唯一性即证。当然直接由待定的系数猛猛求解也是可以的。
注 ：对比 Taylor 展开
对比项 Taylor 展开 Newton 插值 基函数 ( x − x 0 ) k (x-x_0)^k ( x − x 0 ) k ∏ j = 0 k − 1 ( x − x j ) \prod_{j=0}^{k-1}(x-x_j) ∏ j = 0 k − 1 ( x − x j ) 系数 f ( k ) ( x 0 ) k ! \frac{f^{(k)}(x_0)}{k!} k ! f ( k ) ( x 0 ) f [ x 0 , … , x k ] f[x_0,\dots,x_k] f [ x 0 , … , x k ] 信息来源 单点 x 0 x_0 x 0 多点 x 0 , … , x k x_0,\dots,x_k x 0 , … , x k
当 x 1 , … , x k → x 0 x_1,\dots,x_k\to x_0 x 1 , … , x k → x 0 → \to → → \to →
Newton 插值余项
Newton 插值多项式 P n ( t ) P_n(t) P n ( t ) x 0 , … , x n x_0,\dots,x_n x 0 , … , x n f ( t ) f(t) f ( t ) x x x n + 2 n+2 n + 2 x 0 , … , x n , x x_0,\dots,x_n,x x 0 , … , x n , x n + 1 n+1 n + 1 P n + 1 ( t ) P_{n+1}(t) P n + 1 ( t )
根据 Newton 插值的逐次生成性质，增加一个节点只需在末尾追加一项：
P n + 1 ( t ) = P n ( t ) + f [ x 0 , … , x n , x ] ⋅ ( t − x 0 ) ( t − x 1 ) ⋯ ( t − x n ) P_{n+1}(t)=P_n(t)+f[x_0,\dots,x_n,x]\cdot(t-x_0)(t-x_1)\cdots(t-x_n) P n + 1 ( t ) = P n ( t ) + f [ x 0 , … , x n , x ] ⋅ ( t − x 0 ) ( t − x 1 ) ⋯ ( t − x n ) 记 ω n + 1 ( t ) = ( t − x 0 ) ⋯ ( t − x n ) \omega_{n+1}(t)=(t-x_0)\cdots(t-x_n) ω n + 1 ( t ) = ( t − x 0 ) ⋯ ( t − x n )
P n + 1 ( t ) = P n ( t ) + f [ x 0 , … , x n , x ] ⋅ ω n + 1 ( t ) P_{n+1}(t)=P_n(t)+f[x_0,\dots,x_n,x]\cdot\omega_{n+1}(t) P n + 1 ( t ) = P n ( t ) + f [ x 0 , … , x n , x ] ⋅ ω n + 1 ( t ) 因为 x x x P n + 1 ( t ) P_{n+1}(t) P n + 1 ( t ) P n + 1 ( x ) = f ( x ) P_{n+1}(x)=f(x) P n + 1 ( x ) = f ( x ) t = x t=x t = x
f ( x ) = P n ( x ) + f [ x 0 , … , x n , x ] ⋅ ω n + 1 ( x ) f(x)=P_n(x)+f[x_0,\dots,x_n,x]\cdot\omega_{n+1}(x) f ( x ) = P n ( x ) + f [ x 0 , … , x n , x ] ⋅ ω n + 1 ( x ) 移项即得 Newton 型余项 ：
R n ( x ) = f ( x ) − P n ( x ) = f [ x 0 , x 1 , … , x n , x ] ⋅ ω n + 1 ( x ) R_n(x)=f(x)-P_n(x)=f[x_0,x_1,\dots,x_n,x]\cdot\omega_{n+1}(x) R n ( x ) = f ( x ) − P n ( x ) = f [ x 0 , x 1 , … , x n , x ] ⋅ ω n + 1 ( x ) 特点 ：此推导不需要 f f f ，仅依赖函数值与差商定义，本质上是恒等式。
均差表的建立
根据函数表逐列填写：
x i x_i x i f [ x i ] f[x_i] f [ x i ] 一阶差商 二阶差商 三阶差商 四阶差商 x 0 x_0 x 0 f [ x 0 ] f[x_0] f [ x 0 ] f [ x 0 , x 1 ] f[x_0,x_1] f [ x 0 , x 1 ] x 1 x_1 x 1 f [ x 1 ] f[x_1] f [ x 1 ] f [ x 0 , x 1 , x 2 ] f[x_0,x_1,x_2] f [ x 0 , x 1 , x 2 ] f [ x 1 , x 2 ] f[x_1,x_2] f [ x 1 , x 2 ] f [ x 0 , x 1 , x 2 , x 3 ] f[x_0,x_1,x_2,x_3] f [ x 0 , x 1 , x 2 , x 3 ] x 2 x_2 x 2 f [ x 2 ] f[x_2] f [ x 2 ] f [ x 1 , x 2 , x 3 ] f[x_1,x_2,x_3] f [ x 1 , x 2 , x 3 ] f [ x 0 , ⋯ , x 4 ] f[x_0,\cdots,x_4] f [ x 0 , ⋯ , x 4 ] f [ x 2 , x 3 ] f[x_2,x_3] f [ x 2 , x 3 ] f [ x 1 , x 2 , x 3 , x 4 ] f[x_1,x_2,x_3,x_4] f [ x 1 , x 2 , x 3 , x 4 ] x 3 x_3 x 3 f [ x 3 ] f[x_3] f [ x 3 ] f [ x 2 , x 3 , x 4 ] f[x_2,x_3,x_4] f [ x 2 , x 3 , x 4 ] f [ x 3 , x 4 ] f[x_3,x_4] f [ x 3 , x 4 ] x 4 x_4 x 4 f [ x 4 ] f[x_4] f [ x 4 ]
Newton 插值算法实现
Newton 插值标准算法
采用原位存储 。
def newton_interp (x, y, xt): n = len (x) - 1 f = list (y) # 原位存储差商 phi = f[ 0 ] # 初始化 xt 处值 omega = 1.0 # 初始化 N_0(xt) for i in range ( 1 , n + 1 ): omega *= (xt - x[i - 1 ]) # 求 N_i(xt) # 从后往前更新差商，避免覆盖还未用到的值 for k in range (n, i - 1 , - 1 ): f[k] = (f[k - 1 ] - f[k]) / (x[k - i] - x[k]) # 求 f[x_{k-i}, \cdots, x_k] phi += f[i] * omega # 累加当前阶贡献 return phi 复杂度为 O ( n 2 ) \mathcal{O}(n^2) O ( n 2 )
Newton 插值向量化版
import numpy as np def newton_vectorized (x, y, xt_list): x = np.array(x, dtype = float ) f = np.array(y, dtype = float ) n = len (x) - 1 phi = np.full( len (xt_list), f[ 0 ]) omega = np.ones( len (xt_list)) for i in range ( 1 , n + 1 ): omega *= (xt_list - x[i - 1 ]) for k in range (n, i - 1 , - 1 ): f[k] = (f[k - 1 ] - f[k]) / (x[k - i] - x[k]) phi += f[i] * omega return phi Lagrange 与 Newton 插值对比
对比维度 Lagrange 插值 Newton 插值 公式结构 基函数线性组合 均差递推展开 新增节点 全部基函数重算 只增加一项，保留 已有结果 时间复杂度（单点） O ( n 2 ) O(n^2) O ( n 2 ) O ( n 2 ) O(n^2) O ( n 2 ) 适用场景 节点固定，一次性计算 节点逐步增加，在线插值 数值稳定性 可能因大 n n n 同等条件下相同
Runge 现象
现象描述
对于 n n n 节点数量越多，插值结果不一定越接近原函数 。
以 Runge 函数为例：
f ( x ) = 1 1 + 25 x 2 , − 1 ≤ x ≤ 1 f(x) = \frac{1}{1 + 25x^2}, \quad -1 \le x \le 1 f ( x ) = 1 + 25 x 2 1 , − 1 ≤ x ≤ 1 在 [ − 1 , 1 ] [-1,1] [ − 1 , 1 ] 剧烈振荡 。
图 1：Runge 现象
原因分析
根据插值余项公式：
R n ( x ) = f ( n + 1 ) ( ξ ) ( n + 1 ) ! ω n + 1 ( x ) R_n(x) = \dfrac{f^{(n+1)}(\xi)}{(n+1)!}\,\omega_{n+1}(x) R n ( x ) = ( n + 1 )! f ( n + 1 ) ( ξ ) ω n + 1 ( x ) 等距节点时，∣ ω n + 1 ( x ) ∣ |\omega_{n+1}(x)| ∣ ω n + 1 ( x ) ∣ ( n + 1 ) ! (n+1)! ( n + 1 )! f ( x ) f(x) f ( x )
改进方法
方法 思路 Chebyshev 节点 使 ∥ ω n + 1 ∥ ∞ \|\omega_{n+1}\|_\infty ∥ ω n + 1 ∥ ∞ 分段插值 用低次多项式分段拟合，下一讲介绍 样条插值 分段三次多项式，额外要求导数连续
综合复杂度对比
指标 Lagrange（单点） Lagrange（m m m Newton（单点） Newton（m m m 时间复杂度 O ( n 2 ) O(n^2) O ( n 2 ) O ( n 2 + m n ) O(n^2 + mn) O ( n 2 + mn ) O ( n 2 ) O(n^2) O ( n 2 ) O ( n 2 + m n ) O(n^2 + mn) O ( n 2 + mn ) 空间复杂度 O ( 1 ) O(1) O ( 1 ) O ( n ) O(n) O ( n ) O ( n ) O(n) O ( n ) O ( n ) O(n) O ( n ) 新增节点代价 O ( n 2 ) O(n^2) O ( n 2 ) O ( n 2 ) O(n^2) O ( n 2 ) O ( n ) O(n) O ( n ) O ( n + m ) O(n + m) O ( n + m )
