插值

1. 拉格朗日插值
2. 均差与牛顿插值多项式
3. 埃尔米特插值
4. 分段低次插值
5. 三次样条插值

---

有些老师如同大山一样压抑着学生。

华中师范大学真是太多好老师了，可惜这辈子，在本科期间是无缘再见了。

• 埃米尔特插值
• 张晓平老师的主页
• 潘建瑜老师的主页
• 朱老师主页
• 图论的课件

潘老师的学习建议

重基础，多练习，勤思考
I hear and I forget, I see and I remember, I do and I understand.

• 注意掌握各种数值方法的思想和原理，而不是死记硬背
• 注意数值方法的处理技巧，以及与计算机编程的结合（编程实践）
• 注重必要的数值计算训练（必要的手算推导能力）
• 要重视算法的误差分析、收敛性和稳定性

纸上得来终觉浅，绝知此事要躬行。

问题的提出

白话：根据有限数据来推测未知点的值

插值函数：

$$P(x_i)=y_i\qquad (i=0,1,\cdots,n)$$

插值多项式：

$$P(x)=a_0+a_1x+\cdots+a_nx^n$$

确定曲线（插值函数）使其通过插值点，用它近似表示数据的理想函数。

多项式插值

给定$n+1$个点，求次数不超过$n$的多项式。 $a_i$为未知数

$$\begin{cases} a_0+a_1x_0+\cdots+a_nx_0^n=y_0\\ a_0+a_1x_1+\cdots+a_nx^n_1=y_1\\ ...\\ a_0+a_1x_n+\cdots+a_nx^n_n=y_n \end{cases}$$

化成范德蒙矩阵，其行列式值不为零，解唯一。

行列式值与解的关系 待补充

拉格朗日插值

线性插值

只有两个节（数据）点$(x_k,y_k),(x_{k+1},y_{k+1})$

$$P(x)=a_0+a_1x=L_1(x)\tag{线性插值}$$

要求插值函数上的值与节点值严格一致，在图像上显示出来就是一条直线,其斜率为：

$$k=\frac{y_{k+1}-y_{k}}{x_{k+1}-x_{k}}$$

转化成：

$$L_1(x)=y_kl_k(x)+y_{k+1}l_{k+1}(x)$$

其中：

$$l_k(x)=\frac{x-x_{k+1}}{x_k-x_{k+1}}\qquad l_{k+1}(x)=\frac{x-x_{k}}{x_{k+1}-x_{k}}$$

$l_k(x),l_{k+1}(x)$为线性插值基函数，在节点上满足：

$$l_k(x_k)=1\qquad l_k(x_{k+1})=0 \qquad l_{k+1}(x_k)=0 \qquad l_{k+1}(x_{k+1})=1\qquad$$

习题1： (1,0) (2,4)

$$l_1=\frac{x-1}{2-1}=x-1\qquad l_2=\frac{x-2}{1-2}=2-x$$

$$l_1(x)=0+4x-4=4x-4$$

抛物插值

有三个节点时：

$$P(x)=a_0+a_1x+a_2x^2=L_2(x)\tag{抛物插值}$$

进而化简得到基函数：下面节点减去所有其他的点，上面用x替换节点

$$L_2(x)=y_{k-1}l_{k-1}(x)+y_{k}l_{k}(x)+y_{k+1}l_{k+1}(x)$$

其中：

$$l_{k-1}(x)=\frac{(x-x_k)(x-x_{k+1})}{(x_{k-1}-x_{k})(x_{k-1}-x_{k+1})}$$

基函数序数与x序数的关系：相等为1，不相等为0。

$$l_{k-1}(x_{k-1})=1\qquad l_{k-1}(x_{k})=0\qquad l_{k-1}(x_{k+1})=0\qquad ...$$

习题2： (0,1) (1,2) (2,3)

$$l_1(x)=\frac{(x-1)( x-2)}{(0-1)(0-2)}...$$

$$L_2(x)=l_1(x)y_1+...$$

拉格朗日插值多项式

当节点数n比3大时

$$L_n(x_j)=\sum_{k=0}^n y_kl_k(x_j)=y_j\qquad (j=0,1,\cdots,n)$$

n次插值基函数性质：

$$l_j(x_k)= \begin{cases} 1 & k=j\\ 0 & k≠j\\ \end{cases}$$

$$l_{k}(x)=\frac{\left(x-x_{0}\right) \cdots\left(x-x_{k-1}\right)\left(x-x_{k+1}\right) \cdots\left(x-x_{n}\right)}{\left(x_{k}-x_{0}\right) \cdots\left(x_{k}-x_{k-1}\right)\left(x_{k}-x_{k+1}\right) \cdots\left(x_{k}-x_{n}\right)}$$

引入记号：

$$\omega_{n+1}(x)=(x-x_{0})(x-x_{1})\cdots (x-x_{n})$$

所以：

$$\omega^{'}_{n+1}(x_k)=(x_k-x_{0})\cdots (x_k-x_{k-1})(x_k-x_{k+1})\cdots (x_k-x_{n})$$

于是公式(3.5)可以化简为：

$$L_{n}(x)=\sum_{k=0}^{n} y_{k} \frac{\omega_{n+1}(x)}{\left(x-x_{k}\right) \omega_{n+1}^{\prime}\left(x_{k}\right)}$$

特别篇- 拉格朗日插值的C++ 实现

题目：求以下数据的以多项式为基底、以拉格朗日插值多项式为基底、以牛顿插值多项式为基底 的二次插值多项式：

$f(x)$	0	-3	4
x	1	-1	2

解答：

1. 多项式插值：

   直接设未知数，代入解方程即可，可以求得：

$$F(x)= \frac{-20+9x+11x^2}{6}$$

2. 拉格朗日插值多项式

$$\begin{aligned} L_0(x)&=\frac{1}{2}-\frac{x^2}{2}\\ L_1(x)&=\frac{2-3x+x^2}{6}\\ L_2(x)&=\frac{-1+x^2}{3}\\ \therefore L2(x)&=\frac{-14+9x+5x^2}{6} \end{aligned}$$

3. 牛顿插值多项式

   设此二次多项式为：

   $$P_{2}(x)= f\left(x_{0}\right)+f\left[x_{0}, x_{1}\right]\left(x-x_{0}\right)+f\left[x_{0}, x_{1}, x_{2}\right]\left(x-x_{0}\right)\left(x-x_{1}\right)$$

   差商表如下：

   $x_i$	$f(x_i)$	一阶差商	二阶差商
   1	0
   -1	-3	$\frac{3}{2}$
   2	4	$\frac{7}{3}$	$\frac{5}{6}$

   代入可以求得$P_2(x)$为：

   $$P_2(x)=\frac{-14+9x+5x^2}{6}$$

   可见在最高阶次为2的情形下，拉格朗日插值与牛顿插值 得到的结果形式上一致

4. 编程实现：

01_Lagrange.c

// 拉格朗日插值 @抛物插值 
// Monday, September 26, 2022 6:43:15 PM CST
// 田欣洋
// modified from  https://blog.csdn.net/landcruiser007/article/details/79314887
// 采用g++ 编译，由于gcc不支持 new, delete 命令。
#include<stdio.h>

//定义函数 lagrange(坐标x, 坐标y, 插值点x, 最高阶次)
float lagrange(float x[],float y[],float xx,int n)
{
    int i,j;
    float *a,yy=0;
    a=new float[n];
    for(i=0;i<=n-1;i++)
    {
        a[i]=y[i];
        for(j=0;j<=n-1;j++)
            if(j!=i)a[i]*=(xx-x[j])/(x[i]-x[j]);
        yy+=a[i];
    }
    delete a;
    return yy;
}

// 主函数
int main()
{
    //输入所给的坐标点
    float y[3]={0, -3, 4};
    float x[3]={1, -1, 2};
    //输入要求的或要验证的 坐标x的值
    float xx1=1,xx2=0,xx3=2,yy1,yy2,yy3;
    yy1=lagrange(x,y,xx1,3);
    yy2=lagrange(x,y,xx2,3);
    yy3=lagrange(x,y,xx3,3);
    printf("x1=%-20f,y1=%f\n",xx1,yy1);
    printf("x2=%-20f,y2=%f\n",xx2,yy2);
    printf("x3=%-20f,y3=%f\n",xx3,yy3);
    return 0;
}

输出结果

x1=1.000000            ,y1=0.000000
x2=0.000000            ,y2=-2.333333
x3=2.000000            ,y3=4.000000

可见与理论结果是很相符合的。

插值余项与误差估计

罗尔定理

余项定义：

$$R_n(x)=f(x)-L_n(x)$$

$$R_{n}(x)=f(x)-L_{n}(x)=\frac{f^{(n+1)}(\xi)}{(n+1) !} \omega_{n+1}(x)$$

牛顿插值

为何要用牛顿插值法？

拉格朗日插值与插值节点的依赖性太大，用牛顿插值法可以逐次进行插值。

将拉格朗日线性插值看作零次插值的修正；线性插值看作抛物插值的修正；…

差商（微商的离散形式，也叫均差）：

• 点 $x_0, x_1$ 互异，定义函数 $f(x)$ 的一阶均差

$$f\left[x_0, x_1\right]=\frac{f\left(x_1\right)-f\left(x_0\right)}{x_1-x_0}$$

• 点 $x_0, x_1, x_2$ 互异，定义 $f(x)$ 的二阶均差

$$f\left[x_0, x_1, x_2\right]=\frac{f\left[x_1, x_2\right]-f\left[x_0, x_1\right]}{x_2-x_0}$$

• 点 $x_0, x_1, \ldots, x_k$ 互异，定义 $f(x)$ 的 $n$ 阶均差定义为

$$f\left[x_0, x_1, \ldots, x_k\right]=\frac{f\left[x_1, \ldots, x_k\right]-f\left[x_0, \ldots, x_{k-1}\right]}{x_k-x_0}$$

$$f\left[x_{i}, x_{i+1}, \cdots, x_{i+k}\right]=\frac{f\left[x_{i+1}, x_{i+2}, \cdots, x_{i+k}\right]-f\left[x_{i}, x_{i+1}, \cdots, x_{i+k-1}\right]}{x_{i+k}-x_{i}}$$

差商的性质

1. k阶差商

   $$f[x_0, x_1, · · · , x_k] =\sum_{j=0}^k\frac{f(x_j)}{\omega^{'}_{k+1}(x_j)}$$

2. 差商与其所含节点的排列次序无关，即

   $$f[x_i , x_{i+1}] = f[x_{i+1}, x_i]$$

3. 设$f(x)$在包含互异节点$x_0, x_1, · · · , x_n$的闭区间$[a, b]$上有n阶导数，则

   $$f[x_0, x_1, · · · , x_n] =\frac{f^{n}(\xi)}{n!}\quad \xi \in(a,b)$$

牛顿插值多项式

2021-10-16 12:31:19

不学了！！！我要去收拾东西，晚上再坐上32个小时的火车，去向不可知的远方。楼下母亲快把饭做好了，饿！饿！饿！

回到学校再来过！！！

依次迭代可以得到，牛顿插值多项式的形式。把 $x$ 看成一个点，由均差定义得

$$\begin{aligned} &f(x)=f\left(x_0\right)+\left(x-x_0\right) f\left[x, x_0\right] \\ &f\left[x, x_0\right]=f\left[x_0, x_1\right]+\left(x-x_1\right) f\left[x, x_0, x_1\right] . \\ &f\left[x, x_0, x_1\right]=f\left[x_0, x_1, x_2\right]+\left(x-x_2\right) f\left[x, x_0, x_1, x_2\right] \\ &\ldots \ldots \\ &f\left[x, x_0, \ldots, x_{n-1}\right]=f\left[x_0, \ldots, x_n\right]+\left(x-x_n\right) f\left[x, x_0, \ldots, x_n\right] \end{aligned}$$

把后一式依次代入前一式得

$$\begin{aligned} &f(x)=f\left(x_0\right)+f\left[x_0, x_1\right]\left(x-x_0\right)+f\left[x_0, x_1, x_2\right]\left(x-x_0\right)\left(x-x_1\right)+\cdots\\ &+f\left[x_0, x_1, \cdots, x_n\right]\left(x-x_0\right) \cdots\left(x-x_{n-1}\right)\\ &+\underline{f\left[x, x_0, \cdots, x_n\right] \omega_{n+1}(x)}=\underline{N_n(x)}+\underline{R_n(x)}\\ \text{其中 }&\omega_{n+1}(x)=\left(x-x_0\right)\left(x-x_1\right) \cdots\left(x-x_n\right) . \end{aligned}$$

简化一下，于是得到牛顿插值多项式：

$$\begin{aligned} N_{n}(x)=& f\left(x_{0}\right)+f\left[x_{0}, x_{1}\right]\left(x-x_{0}\right)+f\left[x_{0}, x_{1}, x_{2}\right]\left(x-x_{0}\right)\left(x-x_{1}\right)+\cdots \\ &+f\left[x_{0}, x_{1}, \cdots, x_{n}\right]\left(x-x_{0}\right)\left(x-x_{1}\right) \cdots\left(x-x_{n-1}\right) \\ \tilde{R}_{n}(x)=& f\left[x, x_{0}, x_{1}, \cdots, x_{n}\right]\left(x-x_{0}\right)\left(x-x_{1}\right) \cdots\left(x-x_{n}\right) \end{aligned}$$

称$N_n(x)$为$n$次牛顿插值多项式,$R̃_n (x)$为牛顿型插值余项。

• 插值多项式存在且唯一：$P_n(x)≡ L_n(x)$

差分形式

步长：等距节点$x_k=x_0+kh\ (k=0,1,\dots,n)$

称$\Delta f_k=f_{k+1}-f_k$为一阶差分

$$\Delta f_{k}=f_{k+1}-f_{k}$$

• 差分与差商的关系：

$$\begin{aligned} &f\left[x_{0}, x_{1}, \cdots, x_{k}\right]=\frac{\Delta^{k} y_{0}}{k ! h^{k}}, \quad k=1,2, \cdots, n \\ &f\left[x_{n}, x_{n-1}, \cdots, x_{n-k}\right]=\frac{\nabla^{k} y_{n}}{k ! h^{k}}, \quad k=1,2, \cdots, n \end{aligned}$$

• 差分与导数的关系

  $$\Delta^n y_o=h^nf^{(n)}(\xi)\ \xi \in (x_0,x_n)$$

• 牛顿向前插值多项式

  $$P_{n}\left(x_{0}+t h\right)=f\left(x_{0}\right)+t \Delta y_{0}+\frac{t(t-1)}{2 !} \Delta^{2} y_{0}+\cdots+\frac{t(t-1) \cdots(t-n+1)}{n !} \Delta^{n} y_{0}$$

  其余项为：

  $$R_n(x_0+th)=\frac{t(t-1) \cdots(t-n+1)}{n !}h^{n+1}f^{(n+1)}(\xi)\quad\xi \in (x_0,x_n)$$

课后作业

第一题

已知数据：

$$\begin{array}{c|c|c|c|c|c} \hline x_i & 1 & 2 & 2.5 & 4 & 6 \\ \hline u_i & 11.5 & 22.2 & 22.0 & 38.3 & 51.1 \\ \hline \end{array}$$

使用拉格朗日和一般多项式方法进行插值。

解：

1. 一般多项式插值由题可得多项式为：

   $$$$

a_0+a_1+a_2+a_3+a_4=11.5 \
a_0+a_1 2+a_2 4+a_{33} 8+a_4 \times 16=22.4 \
a_{10}+a_1 2.5+a_2 .6 .25+a_3 \times 15.625+a_4 \times 38.0625=22.0 \
a_0+a_1 \times 4+a_2 \times 16+a_3 \times 64+a_4 \times 2156=38.3 \
a_0+a_1 \times 6+a_2 \times 36+a_3 \times 216+a_4 \times 1286=51.1
\end{array}\right.

$$其矩阵形式为：$$

\left[\begin{array}{ccccc}
1 & 1 & 1 & 1 & 1 \
1 & 2 & 4 & 8 & 16 \
1 & 25 & 6.25 & 15.625 & 31.0625 \
1 & 4 & 16 & 64 & 256 \
1 & 6 & 36 & 216 & 1286
\end{array}\right]\left[\begin{array}{l}
a_0 \
a_1 \
a_2 \
a_3 \
a_4
\end{array}\right]=\left[\begin{array}{c}
11.5 \
22.2 \
22.0 \
38.3 \
51.1
\end{array}\right]

$$求解可得：$$

\left{\begin{array}{l}
a_0=-60.26 \
a_1=10 9.205 \
a_2=-43.3825 \
a_3=7.1275 \
a_4=-1.215
\end{array}\right.

$$所以一般多项式的插值公式为：$$

u_1(x)=\sum_n^R a_i x=-1.215 x^4+7.1275 x^3-43.8325 x^2+10 9. 205 x-60.26

2. 拉格朗日插值 设该插值多项式为$u(x)$则

\begin{aligned}
\mu(x)&=\sum_{i=1}^n u P_i(x)\=
&11.5 \times \frac{(x-2)(x-2.5)(x-4)(x-6)}{(-1) \times(-1.5) \times(-3) \times(-5)}+22.2 \times \frac{(x-1)(x-2.5)(x-4)(x-6)}{1 \times(-0.5) \times(-2) \times(-4)}\&+22 \times \frac{(x-1)(x-2)(x-4)(x-6)}{1.5 \times 0.5 \times(-1.5) \times(-3.5)} +38.3 \times \frac{(x-1)(x-2)(x-2.5)(x-6)}{3 \times 2 \times 1.5 \times(-2)}\ &+51.1 \times \frac{(x-1)(x-2)(x-2.5)(x-4)}{5 \times 4 \times 3.5 \times 2}\=
&0.51 x(x-2)(x-2.5)(x-4)(x-6)+(-5.55)(x-1)(x-2.5)(x-4)(x-6)\ &+5.58(x-1)(x-2)(x-4)(x+6)
+(-2.12)(x-1)(x-2)(x-2.5)(x-6)\ & +0.365(x-1)(x-2)(x-2.5)(x-4)
\end

### 第二题 用雅可比迭代求解方程组。（初始向量：$\boldsymbol{x}^{(0)}=\left[\begin{array}{lll}0 & 0 & 0\end{array}\right]^{\mathrm{T}}$）

\left[\begin{array}{ccc}
1 & 2 & -2 \
1 & 1 & 1 \
2 & 2 & 1
\end{array}\right]\left[\begin{array}{l}
x_1 \
x_2 \
x_3
\end{array}\right]=\left[\begin{array}{l}
1 \
3 \
5
\end{array}\right]

$$解：由题可得：$$

\left{ \begin{aligned}
&x_1^{(k+1)}=1-2 x_2^{(k)}+2 x_3^{(k)} \
&x_2^{{(k+1)}=3-x_1}-x_3^{(k)} \
&x_3^{(k+1)}=5-2 x^{(k)}-2 x_2^{(k)}
\end{aligned}\right.

取初始向量：$\boldsymbol{x}^{(0)}=\left[\begin{array}{lll}0 & 0 & 0\end{array}\right]^{\mathrm{T}}$，则有：

\left{ \begin{aligned}
&x^{(1)}=\left[\begin{array}{lll}
1 & 3 & 5
\end{array}\right]^{\top} \
&x^{(2)}=\left[\begin{array}{lll}
5 & -3 & -3
\end{array}\right]^{\top} \
&x^{(3)}=\left[\begin{array}{lll}
1 & 1 & 1
\end{array}\right]^{\top} \
&x^{(4)}=\left[\begin{array}{lll}
1 & 1 & 1
\end{array}\right]^{\top} \
\end{aligned}\right.

因为：$\left\|x^{(3)}-x^{(4)}\right\|<0.01$ 在误差允许范围内所以$[1,1,1]$为最终解。 ### 第三题 试实现高斯-赛德尔迭代算法的程序化，并利用程序求解方程组。

\left[\begin{array}{ccc}
8 & -3 & 2 \
4 & 11 & -1 \
6 & 3 & 12
\end{array}\right]\left[\begin{array}{l}
x_1 \
x_2 \
x_3
\end{array}\right]=\left[\begin{array}{l}
20 \
33 \
36
\end{array}\right]

$$解：根据算法得到如下程序$$

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52

#include<stdio.h> #include<math.h> //设置最大迭代次数 #define N 20 //设置目标精度 #define exc 1e-5 //定义高斯塞德尔迭代函数，*a为行列式矩阵 void gauss(float *a,int n,float x[]) { int i,j,k=1; float d,dx,eps; for(i=0;i<=n-1;i++)x[i]=0.0;//定义初值为：[0,0,0] while(1) { eps=0; for(i=0;i<=n-1;i++)//精度循环 { d=0; for(j=0;j<=n-1;j++)//核心循环 { if(j==i)continue; d+=*(a+i*(n+1)+j)*x[j]; } dx=(*(a+i*(n+1)+n)-d)/(*(a+i*(n+1)+i)); eps+=fabs(dx-x[i]);//求代数精度的绝对值 x[i]=dx; } if(eps<exc)//比较代数精度 {printf("迭代次数为：%d\n",k);return;} if(k>N) { printf("false\n"); return; } k++; } } //主程序 int main() { int i; float a[3][4]= { 8,-3, 2,20, 4,11,-1,33, 6, 3,12,36, }; float x[3]; gauss(a[0],3,x); for(i=0;i<3;i++)printf("x[%d]=%f\n",i+1,x[i]); }

输出结果为

迭代次数为：3
x[1]=3.009814
x[2]=1.996807
x[3]=0.995891