# Lagrange插值
標(biāo)簽（空格分隔）： 數(shù)值方法 插值
---
## 插值
在實(shí)際問(wèn)題中，往往只能觀察到某物理量在有限個(gè)點(diǎn)處的值。如何從這有限個(gè)點(diǎn)值將物理量的函數(shù)恢復(fù)出來(lái)，進(jìn)而研究該物理量的特性呢？
另外，在計(jì)算機(jī)中，只能處理有限、離散的信息。它可以通過(guò)保存一些系數(shù)來(lái)表達(dá)特定的函數(shù)，如：通過(guò)存貯$a_0,a_1,\cdots,a_n$來(lái)表示多項(xiàng)式$a_0+a_1x+\cdots+a_nx^n$，
或者存貯$a,b,c$來(lái)表示函數(shù)$a\sin(bx+c)$。但對(duì)于一個(gè)任意的函數(shù)$f(x)$，計(jì)算機(jī)就沒(méi)辦法來(lái)表示。一種可選的方法是存貯函數(shù)$f(x)$在有限個(gè)點(diǎn)處的函數(shù)值，然后在需要的時(shí)候，盡可能的將函數(shù)$f(x)$恢復(fù)出來(lái)。
它們的共同的**數(shù)學(xué)模型**是：已知函數(shù)$f(x)$在若干個(gè)點(diǎn)處的函數(shù)值，如何恢復(fù)這個(gè)原來(lái)的函數(shù)$f(x)$？
1. 通過(guò)有限個(gè)點(diǎn)的函數(shù)有無(wú)窮多個(gè)。因此，想要得到原來(lái)的函數(shù)是幾乎不可能的；
2. 只能找一個(gè)原來(lái)函數(shù)的近似函數(shù)。因此，可以找一個(gè)性質(zhì)“好”的，形式上熟悉的近似函數(shù)；
$\{x_i\}_{i=0}^n$是$n+1$個(gè)互不相同的點(diǎn)，函數(shù)$f(x)$在這些節(jié)點(diǎn)上的函數(shù)值為$f(x_i)$，$\Phi$是已知的函數(shù)類(lèi)(函數(shù)空間)。若$\Phi$中存在一個(gè)函數(shù)$g(x)$，使得
$$
g(x_i)=f(x_i), i=0,1,\cdots,n
$$
則稱(chēng)$g(x)$是$f(x)$的一個(gè)*插值函數(shù)*。$f(x)$稱(chēng)為*原函數(shù)*。$\{x_i\}$稱(chēng)為*插值節(jié)點(diǎn)*。
找一個(gè)性質(zhì)“好”的近似函數(shù)，可以通過(guò)使用熟悉的函數(shù)類(lèi)來(lái)得到。如：
1. 多項(xiàng)式 $P_n=span\{1,x,x^2,\cdots,x^n\}$；
2. 三角函數(shù) $S_n=span\{1,\cos x,\sin x,\cos 2x, \sin 2x,\cdots,\cos nx,\sin nx\}$
## 插值函數(shù)的存在唯一性
怎么得到$g(x)$？
設(shè)$\Phi=span\{\phi_0(x),\cdots,\phi_m(x)\}$是一個(gè)$m+1$維的線性空間，則
$$
g(x)=a_0\phi_0(x)+\cdots+a_m\phi_m(x)
$$
想辦法求出$a_0,\cdots,a_m$即找到了$g(x)$。由插值的要求，$g(x)$滿足
$$
g(x_i)=a_0\phi_0(x_i)+\cdots+a_m\phi_m(x_i)=f(x_i), i=0,1,\cdots,n
$$
寫(xiě)成矩陣形式，有
$$
\begin{pmatrix}
\phi_0(x_0) & \phi_1(x_0) & \cdots & \phi_m(x_0) \\
\phi_0(x_1) & \phi_1(x_1) & \cdots & \phi_m(x_1) \\
\vdots & \vdots & & \vdots \\
\phi_0(x_n) & \phi_1(x_n) & \cdots & \phi_m(x_n) \\
\end{pmatrix}
\begin{pmatrix}
a_0 \\ a_1 \\ \vdots \\ a_m
\end{pmatrix}
=\begin{pmatrix}
f(x_0) \\ f(x_1) \\ \vdots \\ f(x_n)
\end{pmatrix}
$$
這是一個(gè)關(guān)于系數(shù)$a_0, \cdots, a_m$的線性方程組。 利用線性代數(shù)相關(guān)知識(shí)，可以分析出這個(gè)方程的解結(jié)構(gòu)。特別地，有
**定理**
若$m=n$，且
$$
\left|\begin{matrix}
\phi_0(x_0) & \phi_1(x_0) & \cdots & \phi_n(x_0) \\
\phi_0(x_1) & \phi_1(x_1) & \cdots & \phi_n(x_1) \\
\vdots & \vdots & & \vdots \\
\phi_0(x_n) & \phi_1(x_n) & \cdots & \phi_n(x_n) \\
\end{matrix}\right|\neq 0
$$
則插值函數(shù)$g(x)$存在唯一。
## 插值多項(xiàng)式
在已知的函數(shù)形態(tài)中，多項(xiàng)式是最早接觸到、形式簡(jiǎn)單的一類(lèi)函數(shù)。它特別適合計(jì)算機(jī)使用
1. 計(jì)算一個(gè)多項(xiàng)式只需要四則運(yùn)算，方便編程
2. 一個(gè)多項(xiàng)式的積分與導(dǎo)數(shù)仍然是一個(gè)多項(xiàng)式，且可以方便的得到任意階導(dǎo)數(shù)和任意次的積分
多項(xiàng)式還有一些特性，如：
1. $f(x)$是一個(gè)$n$次多項(xiàng)式，若$f(a)=0$，則有$f(x)=g(x)(x-a)$，且$g(x)$是一個(gè)$n-1$次多項(xiàng)式
2. 若$f(a)=0$，且$f'(a)=0$，則有$f(x)=g(x)(x-a)^2$，且$g(x)$是一個(gè)$n-2$次多項(xiàng)式
由于多項(xiàng)式的諸多優(yōu)秀的特性，在以后的學(xué)習(xí)中，一直使用多項(xiàng)式來(lái)近似。
回到插值過(guò)程，若節(jié)點(diǎn)$\{x_i\}_{i=0}^n$互不相同，在$n$次多項(xiàng)式空間$P_n=span\{1,x, x^2,\cdots, x_n\}$中找插值函數(shù)。由線性代數(shù)知識(shí)，有
$$
\left|\begin{matrix}
1 & x_0 & \cdots & x_0^n \\
1 & x_1 & \cdots & x_1^n \\
\vdots & \vdots & & \vdots \\
1 & x_n & \cdots & x_n^n
\end{matrix}\right|=\prod_{i,j=0, \\ j>i}^n(x_j-x_i)\neq 0
$$
這是一個(gè)Vandoormod行列式。因此，有結(jié)論
**定理**
$n+1$個(gè)互不相同節(jié)點(diǎn)上的$n$次插值多項(xiàng)式存在且唯一。
## 待定系數(shù)法得到插值多項(xiàng)式
前面已經(jīng)給出了求解插值多項(xiàng)式的一種方法--待定系數(shù)。
設(shè)插值多項(xiàng)式是$g(x)=a_0+a_1x+\cdots+a_n x^n$，由插值條件，可以得到線性方程組
$$
\begin{cases}
a_0+a_1 x_0+\cdots+a_n x_0^n & = f(x_0) \\
a_0+a_1 x_1+\cdots+a_n x_1^n & = f(x_1) \\
\cdots \\
a_0+a_1 x_n+\cdots+a_n x_n^n & = f(x_n) \\
\end{cases}
$$
解出系數(shù)$a_0, \cdots, a_n$即可。
**例**
求過(guò)點(diǎn)$(1, 3)$, $(2, 2)$, $(3, 4)$的插值多項(xiàng)式。
**解**
總共3個(gè)點(diǎn)，因此使用2次多項(xiàng)式。設(shè)$g(x)=a_0+a_1 x+a_2 x^2$，則有
$$
\begin{cases}
a_0+ a_1 + a_2? &=3 \\
a_0+2a_1 + 4 a_1&=2 \\
a_0+3a_1 + 9 a_2&=4 \\
\end{cases}
$$
**方法的不足**：需要解一個(gè)線性方程組，而且Vandoormod矩陣是一個(gè)病態(tài)矩陣。當(dāng)$n$比較大($>10$)時(shí)，數(shù)值上的小擾動(dòng)，會(huì)帶來(lái)明顯的計(jì)算誤差。
## Lagrange插值
設(shè)插值多項(xiàng)式是
$$
g(x)=f(x_0)l_0(x)+f(x_1)l_1(x)+\cdots+f(x_n)l_n(x)
$$
其中$l_i(x), i=0,\cdots,n$是次數(shù)不超過(guò)$n$次的多項(xiàng)式。下面，想辦法救出$l_i(x)$。
由
$$
f(x_0)=g(x_0)=f(x_0)l_0(x_0)+f(x_1)l_1(x_0)+\cdots+f(x_n)l_n(x_0)
$$
因此，可以假定成立
$$
l_0(x_0)=1, l_i(x_0)=0, i=1,2,\cdots,n
$$
類(lèi)似地，可以得到
$$
l_0(x_1)=0, l_1(x_1)=1, l_i(x_1)=0, i=2,3,\cdots,n
$$
對(duì)$i=0,1,\cdots,n$，都有
$$
l_i(x_i)=1, l_j(x_i)=0, j=0,1,\cdots,j-1,j+1,\cdots,n
$$
可以得到如下的表
-? | $l_0(x)$ | $l_1(x)$ | $\cdots$ | $l_n(x)$
- | :-: | :-: | :-: | -: |
| $x_0$? | 1 | 0 | $\cdots$ | 0
| $x_1$ | 0 | 1 | $\cdots$ | 0
| $\vdots$ | $\vdots$
|$x_n$ | 0 | 0 | $\cdots$ |1
可以看出，每個(gè)函數(shù)$l_i(x)$滿足
$$
l_i(x_j)=\delta_{ij}=\begin{cases} 1, & i=j \\ 0, & i\neq j \end{cases}
$$
$l_i(x)$有零點(diǎn)$x_0$, $x_1$, $\cdots$, $x_{i-1}$, $x_{i+1}$, $\cdots$, $x_n$， 因此有
$$
l_i(x)=a_i(x-x_0)\cdots(x-x_{i-1})(x-x_{i+1})\cdots(x-x_n)
$$
其中$a_i$是某個(gè)實(shí)數(shù)。又$l_i(x_i)=1$，即
$$
a_i(x_i-x_0)\cdots(x_i-x_{i-1})(x_i-x_{i+1})\cdots(x_i-x_n)=1
$$
解得
$$
a_i=\frac1{(x_i-x_0)\cdots(x_i-x_{i-1})(x_i-x_{i+1})\cdots(x_i-x_n)}
$$
得到
$$
l_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)}
$$
稱(chēng)$l_i(x)$為*Lagrange基函數(shù)*。插值函數(shù)
$$
g(x)=f(x_0)l_0(x)+f(x_1)l_1(x)+\cdots+f(x_n)l_n(x)
$$
為*Lagrange插值*，記為$L_n(x)$。
**例**
求過(guò)點(diǎn)$(1, 3)$, $(2, 2)$, $(3, 4)$的插值多項(xiàng)式。
**解**
插值節(jié)點(diǎn)是$1,2,3$，因此Lagrange插值基函數(shù)是
$$
\begin{aligned}
l_0(x)=&\frac{(x-2)(x-3)}{(1-2)(1-3)} \\
l_1(x)=&\frac{(x-1)(x-3)}{(2-1)(2-3)} \\
l_2(x)=&\frac{(x-1)(x-2)}{(3-1)(3-2)} \\
\end{aligned}
$$
因此，插值函數(shù)為
$$
g(x)=3\frac{(x-2)(x-3)}{(1-2)(1-3)}+2\frac{(x-1)(x-3)}{(2-1)(2-3)}
+4\frac{(x-1)(x-2)}{(3-1)(3-2)}
$$
可以看到$L_n(x)$是一個(gè)次數(shù)不超過(guò)$n$次的多項(xiàng)式，由插值多項(xiàng)式的存在唯一性知，$L_n(x)$與待定系數(shù)法得到的多項(xiàng)式是一樣的，而它避開(kāi)了求解一個(gè)病態(tài)的線性方程組。
## Lagrange插值的偽代碼
``` c
fx=0.0;
for i=0,n:
? li=1.0;
? for j=0,i-1:
? ? li=li*(t-x[j])/(x[i]-x[j]);
? end for
? for j=i+1,n:
? ? li=li*(t-x[j])(x[i]-x[j])
? end for
? fx = fx+y[i]*li;
end for
// fx 即為 插值函數(shù) 在 t 處的值
```
## 誤差分析
插值多項(xiàng)式當(dāng)然跟原函數(shù)有誤差，下面來(lái)分析一下這個(gè)誤差。記
$$
R_n(x)=f(x)-L_n(x)
$$
則，由插值的定義知
$$
R_n(x_i)=f(x_i)-L_n(x_i)=0, i=0,1,\cdots,n
$$
因此，可以設(shè)
$$
R_n(x)=K_n(x)(x-x_0)(x-x_1)\cdots(x-x_n)
$$
其中$K_n(x)$是一個(gè)未知的函數(shù)。
下面，想辦法把$K_n(x)$表達(dá)出來(lái)。對(duì)$f(x)$定義域內(nèi)任意的數(shù)$a$，構(gòu)造輔助函數(shù)
$$
\phi(t)=f(t)-L_n(t)-K_n(a)(t-x_0)\cdots(t-x_n)
$$
為$t$的函數(shù)，則有
$$
\phi(x_i)=f(x_i)-L_n(x_i)=0, i=0,1,\cdots,n
$$
及
$$
\phi(a)=f(a)-L_n(a)-R_n(a)=0
$$
也就是說(shuō)，點(diǎn)$x_i$, $i=0,\cdots,n$及$a$均為函數(shù)$\phi(t)$的零點(diǎn)。設(shè)所有這些點(diǎn)在區(qū)間$[b,c]$內(nèi)，則
* $\phi(t)$具有$n+2$個(gè)互不相同的零點(diǎn)，若$\phi(t)$可導(dǎo)，則由Rolle定理，至少存在$n+1$個(gè)互不相同的點(diǎn)是$\phi'(t)$的零點(diǎn)，這些零點(diǎn)在區(qū)間$(b,c)$內(nèi)。
* 同樣的，若$\phi'(t)$可導(dǎo)，則至少存在$(b,c)$內(nèi)的$n$個(gè)互不相同的$\phi''(t)$的零點(diǎn)。
* 因此，若$\phi(t)$具有$n+1$階導(dǎo)數(shù)，則至少存在一個(gè)點(diǎn)$\xi\in(b,c)$滿足$\phi^{(n+1)}(\xi)=0$。
* 從$\phi(t)$的表達(dá)式中可以看到，$\phi(t)$的光滑性與$f(t)$的光滑性相關(guān)。若$f(x)$足夠光滑（比如：有到$m$階的連續(xù)導(dǎo)數(shù)），則$\phi(t)$也具有$m$階連續(xù)導(dǎo)數(shù)。
因此，若$f(x)$具有$n+1$階導(dǎo)數(shù)，則存在$\xi$滿足
$$
0=\phi^{(n+1)}(\xi)=f^{(n+1)}(\xi)-K_n(a)(n+1)!
$$
即有
$$
K_n(a)=\frac{f^{(n+1)}(\xi)}{(n+1)!}
$$
由$a$的任意性，可以得到誤差表達(dá)式
$$
R_n(x)=f(x)-L_n(x)=\frac{f^{(n+1)}(\xi)}{(n+1)!}\prod_{i=0}^n(x-x_i)
$$
**例**
已知$\sin\frac{\pi}6=\frac12$, $\sin\frac{\pi}4=\frac{\sqrt 2}2$, $\sin\frac{\pi}3=\frac{\sqrt 3}2$，用1次和2次Lagrange插值近似$\sin50^\circ$，并估計(jì)誤差
**解**
求解位置$x=50^\circ=\frac{5\pi}{18}$。記$x_0=\frac{\pi}6$, $x_1=\frac{\pi}4$, $x_2=\frac{\pi}3$。
1. 先用$x_0$, $x_1$做線性插值，近似計(jì)算$x$處的值
$$
L_1(x)=\frac12\frac{x-\frac{\pi}4}{\frac{\pi}6-\frac{\pi}4}
+\frac{\sqrt2}2\frac{x-\frac{\pi}6}{\frac{\pi}4-\frac{\pi}6}
$$
近似得到$\sin50^\circ\approx L_1(\frac{5\pi}{18})\dot=0.77614$。根據(jù)誤差表達(dá)式
$$
R_1(x)=\frac{\sin^{(2)}(\xi)}{2!}(x-\frac{\pi}6)(x-\frac{\pi}4),\quad? \xi\in(\frac{\pi}6,\frac{\pi}3)
$$
又
$$
\sin^{(2)}(\xi)=-\sin(\xi)\in(-\frac{\sqrt3}2,-\frac12), \quad \xi\in(\frac{\pi}6,\frac{\pi}3)
$$
得到
$$
R_1(\frac{5\pi}{18})=-\sin(\xi)\frac12\frac{\pi}{18}\frac{4\pi}{18}
\in(-0.01319,-0.00762)
$$
比較真解$\sin50^\circ=0.7660444\cdots$，得到誤差約為$-0.01001$，落在估計(jì)的區(qū)間內(nèi)。
2. 類(lèi)似地，用$x_1$, $x_2$做插值，得到估計(jì)值和誤差$\tilde R_1$
$$
\sin50^\circ\approx0.76008, \quad \tilde R_1\in(0.00538, 0.00660)
$$
實(shí)際誤差約為$0.00596$
3. 用$x_0$, $x_1$, $x_2$三個(gè)點(diǎn)構(gòu)造二次插值多項(xiàng)式
$$
L_2(x)=\frac12\frac{(x-\frac{\pi}4)(x-\frac{\pi}3)}{(\frac{\pi}6-\frac{\pi}4)(\frac{\pi}6-\frac{\pi}3)}
+\frac{\sqrt2}2\frac{(x-\frac{\pi}6)(x-\frac{\pi}3)}{(\frac{\pi}4-\frac{\pi}6)(\frac{\pi}4-\frac{\pi}3)}
+\frac{\sqrt3}2\frac{(x-\frac{\pi}6)(x-\frac{\pi}4)}{(\frac{\pi}3-\frac{\pi}6)(\frac{\pi}3-\frac{\pi}4)}
$$
及誤差
$$
R_2(x)=\frac{\sin^{(3)}(\xi)}{3!}(x-\frac{\pi}6)(x-\frac{\pi}4)(x-\frac{\pi}3), \xi\in(\frac{\pi}6,\frac{\pi}3)
$$
得到
$$
L_2(\frac{5\pi}{18})\approx0.76543, \quad R_2(\frac{5\pi}{18})\in(0.00044,0.00077)
$$
實(shí)際誤差約為$0.00061$
* 在插值中，所求點(diǎn)$x$落在節(jié)點(diǎn)$\{x_i\}_{i=0}^n$所組成的區(qū)間外部的插值，又稱(chēng)為*外插*。從上面的例子可以看到，通常內(nèi)插的誤差要小于外插的誤差。
* 上面的例子還可以看到，高階插值的誤差要小于低階插值的誤差。
## 事后誤差估計(jì)
在誤差表達(dá)式中，$f^{(n+1)}(\xi)$在實(shí)際應(yīng)用中，是不可計(jì)算的。
1. 若$f(x)$已知，則可以通過(guò)$\xi$的范圍來(lái)得到$f^{(n+1)}(\xi)$的取值范圍，進(jìn)而給出誤差的取值范圍來(lái)。
2. 若$f(x)$也未知，則沒(méi)辦法得到誤差的范圍。
在實(shí)際應(yīng)用中，**事后誤差估計(jì)**方法可以對(duì)誤差進(jìn)行估計(jì)。
設(shè)節(jié)點(diǎn)$\{x_0,x_1,\cdots,x_n\}$上的插值多項(xiàng)式是$L_n(x)$，則有誤差
$$
f(x)-L_n(x)=\frac{f^{(n+1)}(\xi)}{(n+1)!}(x-x_0)\cdots(x-x_n)
$$
替換一個(gè)節(jié)點(diǎn)，節(jié)點(diǎn)組$\{x_1,x_2,\cdots,x_{n+1}\}$上的插值多項(xiàng)式是$\tilde L_n(x)$，誤差為
$$
f(x)-\tilde L_n(x)=\frac{f^{(n+1)}(\eta)}{(n+1)!}(x-x_1)\cdots(x-x_{n+1})
$$
這里，做一個(gè)假設(shè)$f^{(n+1)}(\eta)\approx f^{(n+1)}(\xi)$，則有
$$
\frac{f(x)-L_n(x)}{f(x)-\tilde L_n(x)}\approx\frac{x-x_0}{x-x_{n+1}}
$$
得到
$$
f(x)\approx\frac{x-x_{n+1}}{x_0-x_{n+1}}L_n(x)+\frac{x-x_0}{x_{n+1}-x_0}\tilde L_n(x)
$$
進(jìn)而，$L_n(x)$的誤差可以估計(jì)為
$$
f(x)-L_n(x)\approx\frac{x-x_0}{x_{n+1}-x_0}\left(\tilde L_n(x)-L_n(x)\right)
$$
這樣，得到了**事后誤差估計(jì)**。
**注意**： 推導(dǎo)過(guò)程中，假設(shè)$f^{(n+1)}(\eta)\approx f^{(n+1)}(\xi)$ 是沒(méi)有數(shù)學(xué)依據(jù)的， 因此得到的誤差估計(jì)即可能比真實(shí)誤差大，也可能比真實(shí)誤差小。從數(shù)學(xué)上說(shuō)，這個(gè)數(shù)值沒(méi)有意義。
但實(shí)際使用中，效果挺好。
