文章目录

一、实验目的
二、实验设备信息
三、实验内容
- （一）拉格朗日插值多项式
- （二）牛顿插值多项式
四、实验步骤
- （一）拉格朗日插值函数实现
- （二）牛顿插值函数实现
- （三）观察拉格朗日插值和牛顿插值的龙格现象
五、实验结果
- （一）拉格朗日插值
- （二）牛顿插值
- （三）龙格现象
六、总结
七、代码

一、实验目的

（一）以实验的方式，理解高阶插值的病态性，观察拉格朗日插值和牛顿插值的龙格现象。

（二）进一步理解拉格朗日插值和牛顿插值的思想，体验插值问题在程序中的应用过程。

二、实验设备信息

（一）设备信息：

图1 - 设备信息

（二）运行环境 : Python 3.8

（三）IDE : PyCharm Professional 2022.2.2

三、实验内容

（一）拉格朗日插值多项式

1、为了构造通过 n + 1 n+1 n+1 个节点 x 0 < x 1 < ⋯ < x n x_0<x_1<\cdots<x_n x0<x1<⋯<xn 的 n n n 次插值多项式：
L n ( x j ) = y i , j = 0 , 1 , ⋯ , n . ( 1 ) L_n\left(x_j\right)=y_i,\quad j=0,1,\cdots,n.\quad\quad\quad (1) Ln(xj)=yi,j=0,1,⋯,n.(1)
先定义 n n n 次插值基函数。

2、 n n n 次插值基函数定义：若 n n n 次多项式 l j ( x ) ( j = 0 , 1 , ⋯ , n ) l_j\left(x\right)\left(j = 0,1,\cdots,n\right) lj(x)(j=0,1,⋯,n) 在 n + 1 n+1 n+1 个节点 x 0 < x 1 < ⋯ < x n x_0<x_1<\cdots<x_n x0<x1<⋯<xn 上满足条件：
l j ( x k ) = { 1 , k = j , 0 , k ≠ j , j , k = 0 , 1 , ⋯ , n ( 2 ) l_j\left(x_k\right)=\left\{\begin{array}{ll} 1, & k=j, \\ 0, & k \neq j, \end{array} \quad j, k=0,1, \cdots, n\right. \quad\quad\quad (2) lj(xk)={1,0,k=j,k=j,j,k=0,1,⋯,n(2)
就称这 n + 1 n+1 n+1 个 n n n 次多项式 l 0 ( x ) , l 1 ( x ) , ⋯ , l n ( x ) l_0(x),l_1(x),\cdots,l_n(x) l0(x),l1(x),⋯,ln(x) 为节点 x 0 , x 1 , ⋯ , x n x_0,x_1,\cdots,x_n x0,x1,⋯,xn 上的 n n n 次插值基函数。利用类似的推导方法，可得到 n n n 次插值基函数为：
l k ( x ) = ( x − x 0 ) ⋯ ( x − x k − 1 ) ( x − x k + 1 ) ⋯ ( x − x n ) ( x k − x 0 ) ⋯ ( x k − x k − 1 ) ( x k − x k + 1 ) ⋯ ( x k − x n ) , k = 0 , 1 , ⋯ , n ( 3 ) 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)}, \quad k=0,1, \cdots, n \quad\quad\quad (3) lk(x)=(xk−x0)⋯(xk−xk−1)(xk−xk+1)⋯(xk−xn)(x−x0)⋯(x−xk−1)(x−xk+1)⋯(x−xn),k=0,1,⋯,n(3)
显然，（3）满足条件（2）。于是满足条件（1）的插值多项式 L n ( x ) L_n(x) Ln(x) 可以表示为：
L n ( x ) = ∑ k = 0 n y k l k ( x ) . ( 4 ) L_n(x)=\sum_{k=0}^n y_k l_k(x). \quad\quad\quad (4) Ln(x)=k=0∑nyklk(x).(4)

3、拉格朗日插值多项式：由 l k ( x ) l_k\left(x\right) lk(x) 的定义，知：
L n ( x j ) = ∑ k = 0 n y k l k ( x j ) = y j , j = 0 , 1 , ⋯ , n . L_n\left(x_j\right)=\sum_{k=0}^n y_k l_k\left(x_j\right)=y_j, \quad j=0,1, \cdots, n . Ln(xj)=k=0∑nyklk(xj)=yj,j=0,1,⋯,n.
形如（4）式的插值多项式 L n ( x ) L_n\left(x\right) Ln(x) 称为拉格朗日插值多项式。若引入记号：
ω n + 1 ( x ) = ( x − x 0 ) ( x − x 1 ) ⋯ ( x − x n ) , ( 5 ) \omega_{n+1}(x)=\left(x-x_0\right)\left(x-x_1\right) \cdots\left(x-x_n\right), \quad\quad\quad (5) ωn+1(x)=(x−x0)(x−x1)⋯(x−xn),(5)
容易求得：
ω n + 1 ′ ( x k ) = ( x k − x 0 ) ⋯ ( x k − x k − 1 ) ( x k − x k + 1 ) ⋯ ( x k − x n ) . \omega_{n+1}^{\prime}\left(x_k\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). ωn+1′(xk)=(xk−x0)⋯(xk−xk−1)(xk−xk+1)⋯(xk−xn).
于是公式（4）可改写成：
L n ( x ) = ∑ k = 0 n y k ω n + 1 ( x ) ( x − x k ) ω n + 1 ′ ( x k ) . ( 6 ) 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)}. \quad\quad\quad (6) Ln(x)=k=0∑nyk(x−xk)ωn+1′(xk)ωn+1(x).(6)

（二）牛顿插值多项式

1、插值多项式的逐次生成：一般情形已知 f f f 在插值点 x i ( i = 0 , 1 , ⋯ , n ) x_i\left(i=0,1,\cdots,n\right) xi(i=0,1,⋯,n) 上的值为 f ( x i ) ( i = 1 , 2 , ⋯ , n ) f\left(x_i\right)\left(i=1,2,\cdots,n\right) f(xi)(i=1,2,⋯,n) ，要求 n n n 次插值多项式 P n ( x ) P_n\left(x\right) Pn(x) 满足条件：
P n ( x i ) = f ( x i ) , i = 0 , 1 , ⋯ , n , ( 7 ) P_n\left(x_i\right)=f\left(x_i\right), \quad i=0,1, \cdots, n, \quad\quad\quad (7) Pn(xi)=f(xi),i=0,1,⋯,n,(7)
则 P n ( x ) P_n\left(x\right) Pn(x) 可表示为：
P n ( x ) = a 0 + a 1 ( x − x 0 ) + ⋯ + a n ( x − x 0 ) ⋯ ( x − x n − 1 ) , ( 8 ) P_n(x)=a_0+a_1\left(x-x_0\right)+\cdots+a_n\left(x-x_0\right) \cdots\left(x-x_{n-1}\right), \quad\quad\quad (8) Pn(x)=a0+a1(x−x0)+⋯+an(x−x0)⋯(x−xn−1),(8)
其中 a 0 , a 1 , ⋯ , a n a_0,a_1,\cdots,a_n a0,a1,⋯,an 为待定系数，可由条件（7）确定。

2、均差及其性质

（1）均差的定义：称 f [ x 0 , x k ] = f ( x k ) − f ( x 0 ) x k − x 0 f\left[x_0,x_k\right]=\frac{f\left(x_k\right)-f\left(x_0\right)}{x_k - x_0} f[x0,xk]=xk−x0f(xk)−f(x0) 为函数 f ( x ) f\left(x\right) f(x) 关于点 x 0 , x k x_0,x_k x0,xk 的一阶均差。 f [ x 0 , x 1 , x k ] = f [ x 0 , x k ] − f [ x 0 , x 1 ] x k − x 1 f\left[x_0,x_1,x_k\right]=\frac{f\left[x_0,x_k\right]-f\left[x_0,x_1\right]}{x_k - x_1} f[x0,x1,xk]=xk−x1f[x0,xk]−f[x0,x1] 称为 f ( x ) f\left(x\right) f(x) 的二阶均差。一般地，称：
f [ x 0 , x 1 , ⋯ , x k ] = f [ x 0 , ⋯ , x k − 2 , x k ] − f [ x 0 , x 1 , ⋯ , x k − 1 ] x k − x k − 1 ( 9 ) f\left[x_0, x_1, \cdots, x_k\right]=\frac{f\left[x_0, \cdots, x_{k-2}, x_k\right]-f\left[x_0, x_1, \cdots, x_{k-1}\right]}{x_k-x_{k-1}} \quad\quad\quad (9) f[x0,x1,⋯,xk]=xk−xk−1f[x0,⋯,xk−2,xk]−f[x0,x1,⋯,xk−1](9)
为 f ( x ) f\left(x\right) f(x) 的 k k k 阶均差（均差也成为差商）。

（2）均差的基本性质

a. k k k 阶均差可表示为函数值 f ( x 0 ) , f ( x 1 ) , ⋯ , f ( x k ) f\left(x_0\right),f\left(x_1\right),\cdots,f\left(x_k\right) f(x0),f(x1),⋯,f(xk) 的线性组合，即
f [ x 0 , x 1 , ⋯ , x k ] = ∑ j = 0 k f ( x j ) ( x j − x 0 ) ⋯ ( x j − x j − 1 ) ( x j − x j + 1 ) ⋯ ( x j − x k ) . ( 10 ) f\left[x_0, x_1, \cdots, x_k\right]=\sum_{j=0}^k \frac{f\left(x_j\right)}{\left(x_j-x_0\right) \cdots\left(x_j-x_{j-1}\right)\left(x_j-x_{j+1}\right) \cdots\left(x_j-x_k\right)} . \quad\quad\quad (10) f[x0,x1,⋯,xk]=j=0∑k(xj−x0)⋯(xj−xj−1)(xj−xj+1)⋯(xj−xk)f(xj).(10)
这个性质也表明均差与节点的排列次序无关，称为均差的对称性，即：
f [ x 0 , x 1 , ⋯ , x k ] = f [ x 1 , x 0 , x 2 , ⋯ , x k ] = ⋯ = f [ x 1 , ⋯ , x k , x 0 ] . f\left[x_0, x_1, \cdots, x_k\right]=f\left[x_1, x_0, x_2, \cdots, x_k\right]=\cdots=f\left[x_1, \cdots, x_k, x_0\right]. f[x0,x1,⋯,xk]=f[x1,x0,x2,⋯,xk]=⋯=f[x1,⋯,xk,x0].
b. 由性质（a）和式（9）可得：
f [ x 0 , x 1 , ⋯ , x k ] = f [ x 1 , x 2 , ⋯ , x k ] − f [ x 0 , x 1 , ⋯ , x k − 1 ] x k − x 0 . ( 9 ) ′ f\left[x_0, x_1, \cdots, x_k\right]=\frac{f\left[x_1, x_2, \cdots, x_k\right]-f\left[x_0, x_1, \cdots, x_{k-1}\right]}{x_k-x_0} . \quad\quad\quad (9)' f[x0,x1,⋯,xk]=xk−x0f[x1,x2,⋯,xk]−f[x0,x1,⋯,xk−1].(9)′
c. 若 f ( x ) f\left(x\right) f(x) 在 [ a , b ] \left[a,b\right] [a,b] 上存在 n n n 阶导数，且节点 x 0 , x 1 , ⋯ , x n ∈ [ a , b ] x_0,x_1,\cdots,x_n \in \left[a,b\right] x0,x1,⋯,xn∈[a,b] ，则 n n n 阶均差与导数的关系为：
f [ x 0 , x 1 , ⋯ , x n ] = f ( n ) ( ξ ) n ! , ξ ∈ [ a , b ] . ( 10 ) f\left[x_0, x_1, \cdots, x_n\right]=\frac{f^{(n)}(\xi)}{n !}, \quad \xi \in[a, b] . \quad\quad\quad (10) f[x0,x1,⋯,xn]=n!f(n)(ξ),ξ∈[a,b].(10)

（3）均差计算可列均差表如下（表1 - 均差表）

表1 - 均差表

x k x_k xk	f ( x k ) f\left(x_k\right) f(xk)	一阶均差	二阶均差	三阶均差	四阶均差
x 0 x_0 x0	f ( x 0 ) f\left(x_0\right) f(x0)
x 1 x_1 x1	f ( x 1 ) f\left(x_1\right) f(x1)	f [ x 0 , x 1 ] f\left[x_0,x_1\right] f[x0,x1]
x 2 x_2 x2	f ( x 2 ) f\left(x_2\right) f(x2)	f [ x 1 , x 2 ] f\left[x_1,x_2\right] f[x1,x2]	f [ x 0 , x 1 , x 2 ] f\left[x_0,x_1,x_2\right] f[x0,x1,x2]
x 3 x_3 x3	f ( x 3 ) f\left(x_3\right) f(x3)	f [ x 2 , x 3 ] f\left[x_2,x_3\right] f[x2,x3]	f [ x 1 , x 2 , x 3 ] f\left[x_1,x_2,x_3\right] f[x1,x2,x3]	f [ x 0 , x 1 , x 2 , x 3 ] f\left[x_0,x_1,x_2,x_3\right] f[x0,x1,x2,x3]
x 4 x_4 x4	f ( x 4 ) f\left(x_4\right) f(x4)	f [ x 3 , x 4 ] f\left[x_3,x_4\right] f[x3,x4]	f [ x 2 , x 3 , x 4 ] f\left[x_2,x_3,x_4\right] f[x2,x3,x4]	f [ x 1 , x 2 , x 3 , x 4 ] f\left[x_1,x_2,x_3,x_4\right] f[x1,x2,x3,x4]	f [ x 1 , x 2 , x 3 , x 4 , x 5 ] f\left[x_1,x_2,x_3,x_4,x_5\right] f[x1,x2,x3,x4,x5]
⋯ \cdots ⋯	⋯ \cdots ⋯	⋯ \cdots ⋯	⋯ \cdots ⋯	⋯ \cdots ⋯	⋯ \cdots ⋯

3、牛顿插值多项式：根据均差定义，将 x x x 看成 [ a , b ] \left[a,b\right] [a,b] 上一点，可得：
f ( x ) = f ( x 0 ) + f [ x , x 0 ] ( x − x 0 ) , f [ x , x 0 ] = f [ x 0 , x 1 ] + f [ x , x 0 , x 1 ] ( x − x 1 ) , ⋯ f [ x , x 0 , ⋯ , x n ] = f [ x 0 , x 1 , ⋯ , x n ] + f [ x , x 0 , ⋯ , x n ] ( x − x n ) . f\left(x\right)=f\left(x_0\right)+f\left[x,x_0\right]\left(x-x_0\right),\\ f\left[x,x_0\right]=f\left[x_0,x_1\right]+f\left[x,x_0,x_1\right]\left(x-x_1\right),\\ \cdots\\ f\left[x,x_0,\cdots,x_n\right]=f\left[x_0,x_1,\cdots,x_n\right]+f\left[x,x_0,\cdots,x_n\right]\left(x-x_n\right). f(x)=f(x0)+f[x,x0](x−x0),f[x,x0]=f[x0,x1]+f[x,x0,x1](x−x1),⋯f[x,x0,⋯,xn]=f[x0,x1,⋯,xn]+f[x,x0,⋯,xn](x−xn).
只要把最后一步代入前一式，就得到：
f ( x ) = f ( x 0 ) + f [ x 0 , x 1 ] ( x − x 0 ) + f [ x 0 , x 1 , x 2 ] ( x − x 0 ) ( x − x 1 ) + ⋯ + f [ x 0 , x 1 , ⋯ , x n ] ( x − x 0 ) ⋯ ( x − x n − 1 ) + f [ x , x 0 , ⋯ , x n ] ω n + 1 ( x ) = P n ( x ) + R n ( x ) \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) \\ &+f\left[x, x_0, \cdots, x_n\right] \omega_{n+1}(x)=P_n(x)+R_n(x) \end{aligned} f(x)=f(x0)+f[x0,x1](x−x0)+f[x0,x1,x2](x−x0)(x−x1)+⋯+f[x0,x1,⋯,xn](x−x0)⋯(x−xn−1)+f[x,x0,⋯,xn]ωn+1(x)=Pn(x)+Rn(x)
其中，
P n ( x ) = f ( x 0 ) + f [ x 0 , x 1 ] ( x − x 0 ) + f [ x 0 , x 1 , x 2 ] ( x − x 0 ) ( x − x 1 ) + ⋯ + f [ x 0 , x 1 , ⋯ , x n ] ( x − x 0 ) ⋯ ( x − x n − 1 ) ( 11 ) \begin{aligned} P_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) \cdots\left(x-x_{n-1}\right)\quad\quad\quad (11) \\ \end{aligned} Pn(x)=f(x0)+f[x0,x1](x−x0)+f[x0,x1,x2](x−x0)(x−x1)+⋯+f[x0,x1,⋯,xn](x−x0)⋯(x−xn−1)(11)

R n ( x ) = f ( x ) − P n ( x ) = f [ x , x 0 , ⋯ , x n ] ω n + 1 ( x ) ( 12 ) \begin{aligned} R_n(x)=& f(x)-P_n(x)=f\left[x, x_0, \cdots, x_n\right] \omega_{n+1}(x)\quad\quad\quad (12) \end{aligned} Rn(x)=f(x)−Pn(x)=f[x,x0,⋯,xn]ωn+1(x)(12)

其中 ω n + 1 ( x ) \omega_{n+1}\left(x\right) ωn+1(x) 由式（5）定义。由式（11）确定的多项式 P n ( x ) P_n\left(x\right) Pn(x) 显然满足插值条件（7），且次数不超过 n n n ，它就是形如（8）式的多项式，其系数为：
a k = f [ x 0 , x 1 , ⋯ , x k ] , k = 0 , 1 , ⋯ , n . a_k=f\left[x_0,x_1,\cdots,x_k\right],\quad k=0,1,\cdots,n. ak=f[x0,x1,⋯,xk],k=0,1,⋯,n.
我们称 P n ( x ) P_n\left(x\right) Pn(x) 为牛顿均差插值多项式。

四、实验步骤

（一）拉格朗日插值函数实现

1、完成 n n n 次插值基函数的代码实现。

2、完成拉格朗日插值多项式的代码实现。

3、利用数据检验拉格朗日插值多项式的求解正确性。

已知未知函数 y = f ( x ) y=f\left(x\right) y=f(x) 的6个观测点 ( x i , y i ) ( i = 0 , 1 , ⋯ , 5 ) \left(x_i,y_i\right)\left(i=0,1,\cdots,5\right) (xi,yi)(i=0,1,⋯,5) 的值如表2所列，计算拉格朗日插值函数 y = L ( x ) y=L\left(x\right) y=L(x) ，并求出 x = 1.5 , 2.6 x=1.5,2.6 x=1.5,2.6 处函数的估计值。

表2 - 拉格朗日插值观测点数据

x i x_i xi 1 2 3 4 5 6

y i y_i yi 16 18 21 17 15 12

x i x_i xi	1	2	3	4	5	6
y i y_i yi	16	18	21	17	15	12

4、将得出的插值函数可视化并观测结果。

（二）牛顿插值函数实现

1、完成 n n n 阶均差求解的代码实现。

2、完成牛顿插值多项式的代码实现。

3、利用数据检验牛顿插值多项式的求解正确性。

已知未知函数 y = f ( x ) y=f\left(x\right) y=f(x) 的6个观测点 ( x i , y i ) ( i = 0 , 1 , ⋯ , 5 ) \left(x_i,y_i\right)\left(i=0,1,\cdots,5\right) (xi,yi)(i=0,1,⋯,5) 的值如表3所列，计算牛顿插值函数 y = L ( x ) y=L\left(x\right) y=L(x) ，并求出 x = 1.5 , 2.6 x=1.5,2.6 x=1.5,2.6 处函数的估计值。

表3 - 牛顿插值观测点数据

x i x_i xi 1 2 3 4 5 6

y i y_i yi 16 18 21 17 15 12

x i x_i xi	1	2	3	4	5	6
y i y_i yi	16	18	21	17	15	12

4、将得出的插值函数可视化并观测结果。

（三）观察拉格朗日插值和牛顿插值的龙格现象

1、在区间 [ 5 , − 5 ] \left[5,-5\right] [5,−5] 上，用 n + 1 n+1 n+1 个等距节点作拉格朗日插值多项式 L n ( x ) L_n\left(x\right) Ln(x) 以及牛顿插值多项式 P n ( x ) P_n\left(x\right) Pn(x) ，使得它在节点处的值与函数 y = 1 1 + x 2 y=\frac{1}{1+x^2} y=1+x21 在对应节点处的值相等。

2、分析 n = 6 , 8 , 10 n=6,8,10 n=6,8,10 时，多项式的次数与逼近误差的关系。表4，5，6为 n n n 取不同值时， x i x_i xi 的取值以及对应的 y i y_i yi 的取值。

表4 - n=6 时 x 和 y 的取值

x i x_i xi	-3.5	-2.5	-1.5	1.5	2.5	3.5
y i y_i yi	0.08	0.14	0.31	0.31	0.14	0.08

表5 - n=8 时 x 和 y 的取值

x i x_i xi	-4.5	-3.5	-2.5	-1.5	1.5	2.5	3.5	4.5
y i y_i yi	0.05	0.08	0.14	0.31	0.31	0.14	0.08	0.05

表6 - n=10 时 x 和 y 的取值

x i x_i xi	-4.5	-3.5	-2.5	-1.5	-0.5	0.5	1.5	2.5	3.5	4.5
y I y_I yI	0.05	0.08	0.14	0.31	0.8	0.8	0.31	0.14	0.08	0.05

3、求得 n n n 取不同值时的插值函数，可视化插值函数以及原函数 y y y ，根据可视化结果，对比插值函数图像与原函数 y y y 的图像，观察高阶插值函数的龙格现象。

五、实验结果

（一）拉格朗日插值

1、根据表2的6个观测点求得拉格朗日插值多项式后，分别代入 x = 1.5 x=1.5 x=1.5 和 x = 2.6 x=2.6 x=2.6 求得： L 5 ( 1.5 ) = 14.918 L_5\left(1.5\right)=14.918 L5(1.5)=14.918， L 5 ( 2.6 ) = 20.8846 L_5\left(2.6\right)=20.8846 L5(2.6)=20.8846 。

图2 - 拉格朗日插值代码运行结果

2、将求得的拉格朗日插值多项式图像可视化，得到图3所示结果。

图3 - 拉格朗日插值多项式

3、从图2中可以清晰看出，所得到的拉格朗日插值多项式满足条件（1），即 L n ( x j ) = f ( x j ) ( x j = 0 , 1 , ⋯ , n ) L_n\left(x_j\right)=f\left(x_j\right)\left(x_j=0,1,\cdots,n\right) Ln(xj)=f(xj)(xj=0,1,⋯,n) ，从而证明拉格朗日插值多项式求解的正确性。

（二）牛顿插值

1、根据表3的6个观测点求得牛顿插值多项式后，分别代入 x = 1.5 x=1.5 x=1.5 和 x = 2.6 x=2.6 x=2.6 求得： P 5 ( 1.5 ) = 14.918 P_5\left(1.5\right)=14.918 P5(1.5)=14.918， P 5 ( 2.6 ) = 20.8846 P_5\left(2.6\right)=20.8846 P5(2.6)=20.8846 。

图4 - 牛顿插值代码运行结果

2、将求得的牛顿插值多项式图像可视化，得到图5所示结果。

图5 - 牛顿插值多项式

3、从图5中可以清晰看出，所得到的牛顿插值多项式满足条件（7），即 P n ( x i ) = f ( x i ) ( i = 0 , 1 , ⋯ , n ) P_n\left(x_i\right)=f\left(x_i\right)\left(i=0,1,\cdots,n\right) Pn(xi)=f(xi)(i=0,1,⋯,n) ，从而证明牛顿插值多项式求解的正确性。

（三）龙格现象

1、拉格朗日插值的龙格震荡现象

（1）先做出在区间 [ − 5 , 5 ] \left[-5,5\right] [−5,5] 上，函数 y = 1 1 + x 2 y=\frac{1}{1+x^2} y=1+x21 的函数图像。

（2）根据表4、5、6所给出的 x x x 的值以及所对应的 y y y 的值，分别求解出不同的插值点个数 n n n 所对应的拉格朗日插值多项式，并将函数图像可视化，得到图6所示结果。

图6 - 拉格朗日插值龙格现象

（3）根据图6，与原函数 y = 1 1 + x 2 y=\frac{1}{1+x^2} y=1+x21 的图像相比较，当插值节点数 n n n 增加时，由于插值多项式的次数也随之增加，拉格朗日插值多项式逼近函数 y = 1 1 + x 2 y=\frac{1}{1+x^2} y=1+x21 的效果随之下降，插值函数相较于原函数更加失真，振荡现象严重。

2、牛顿插值的龙格震荡现象

（1）先做出在区间 [ − 5 , 5 ] \left[-5,5\right] [−5,5] 上，函数 y = 1 1 + x 2 y=\frac{1}{1+x^2} y=1+x21 的函数图像。

（2）根据表4、5、6所给出的 x x x 的值以及所对应的 y y y 的值，分别求解出不同的插值点个数 n n n 所对应的牛顿插值多项式，并将函数图像可视化，得到图7所示结果。

图7 - 牛顿插值龙格现象

（3）根据图7，与原函数 y = 1 1 + x 2 y=\frac{1}{1+x^2} y=1+x21 的图像相比较，当插值节点数 n n n 增加时，由于插值多项式的次数也随之增加，牛顿插值多项式逼近函数 y = 1 1 + x 2 y=\frac{1}{1+x^2} y=1+x21 的效果随之下降，插值函数相较于原函数更加失真，振荡现象严重。

六、总结

在我们的直观印象中，利用插值多项式对某一函数进行逼近，计算相应的函数值，当多项式的次数越高，插值点的数量越多时，函数值的预测就会越准确。然而，实验的结果表明，高次插值不但计算复杂，而且预测效果反而不理想。从实验当中，我们更加直观且深刻的理解到了高阶插值的病态性，也观测到了随着插值次数的增加，插值函数所出现的龙格现象。

查阅资料后，了解到，出现龙格现象是由于对任意的插值节点，当 n → + ∞ n\to+\infty n→+∞ 时，插值函数不一定收敛于 f ( x ) f\left(x\right) f(x) 。从另一个角度来理解，是因为拉格朗日插值以及牛顿插值反映的是所有插值点的特征，是对于插值点的整体求解得出的多项式，容易受到少量特殊的点的影响，从而造成对总体特征的预测准确度的下降。

高阶插值函数的病态性迫使我们探求其他一些更稳定、更高效的插值格式，如分段插值的方法。在未来学习中，我们也许会对效果更好的插值方法进行探索。

七、代码

# 导入所需要的库
import numpy as np
from pyecharts.charts import Line, Scatter
from pyecharts import options as opts

# n 次插值基函数
def n_basis_function(x, k, x_list):numerator = 1for i in range(0, k):numerator = numerator * (x - x_list[i])for i in range(k + 1, len(x_list)):numerator = numerator * (x - x_list[i])denominator = 1for i in range(0, k):denominator = denominator * (x_list[k] - x_list[i])for i in range(k + 1, len(x_list)):denominator = denominator * (x_list[k] - x_list[i])l_k = numerator / denominatorreturn l_k

# 拉格朗日插值多项式
# x为需要求值的点 x_list为插值点的横坐标 y_list为插值点的纵坐标
def Lagrange(x, x_list, y_list):L_n = 0for i in range(len(x_list)):L_n = L_n + y_list[i] * n_basis_function(x, i, x_list)return L_n

# 拉格朗日插值多项式检验
x_list = [1.0, 2.0, 3.0, 4.0, 5.0, 6.0]
y_list = [16.0, 18.0, 21.0, 17.0, 15.0, 12.0]x = 1.5
print(Lagrange(x, x_list, y_list))
x = 2.6
print(Lagrange(x, x_list, y_list))# 可视化
pre_x_list = []
pre_y_list = []
for i in range(0, 6001, 6):pre_x_list.append(i / 1000)pre_y_list.append(round(Lagrange(i / 1000, x_list, y_list), 4))scatter = Scatter(init_opts = opts.InitOpts(theme = "macarons"))
scatter.add_xaxis(xaxis_data = np.array(x_list))
scatter.add_yaxis(series_name = "实际值",y_axis = np.array(y_list),label_opts = opts.LabelOpts(position = "outside"))
scatter.set_global_opts(title_opts = opts.TitleOpts("拉格朗日插值多项式",pos_left = "40%",pos_top = "5%"),legend_opts = opts.LegendOpts(pos_top = "15%"))line = Line(init_opts = opts.InitOpts(theme = "macarons"))
line.add_xaxis(xaxis_data = pre_x_list)
line.add_yaxis(series_name = "拉格朗日插值函数",y_axis = np.array(pre_y_list),label_opts = opts.LabelOpts(is_show = False),symbol_size = 0,linestyle_opts = opts.LineStyleOpts(width = 2))
line.set_global_opts(xaxis_opts = opts.AxisOpts(is_show = False))
scatter.overlap(line)
scatter.render("拉格朗日插值.html")

# n 阶均差求解 (利用均差表)
def difference_quotient(x_list, y_list):denominator = x_list[-1] - x_list[0]if len(x_list) > 2:return (difference_quotient(x_list[1:], y_list[1:]) - difference_quotient(x_list[:-1], y_list[:-1])) / denominatorelif len(x_list) == 2:return (y_list[-1] - y_list[0]) / denominatorelse:return y_list[0]

# n 阶均差求解 (利用均差的性质)
def difference_quotient(x_list, y_list):if len(x_list) >= 2:d_q = 0for i in range(len(x_list)):denominator = 1.0for j in range(0, i):denominator = denominator * (x_list[i] - x_list[j])for j in range(i + 1, len(x_list)):denominator = denominator * (x_list[i] - x_list[j])d_q = d_q + (y_list[i] / denominator)return d_qelse:return y_list[0]

# 牛顿插值多项式
def func(x, len, x_list):f = 1for i in range(len):f = f * (x - x_list[i])return fdef Newton(x, x_list, y_list):P_n = 0for i in range(len(x_list)):P_n = P_n + difference_quotient(x_list[:i + 1], y_list[:i + 1]) * func(x, i, x_list)return P_n

# 牛顿插值多项式检验
x_list = [1.0, 2.0, 3.0, 4.0, 5.0, 6.0]
y_list = [16.0, 18.0, 21.0, 17.0, 15.0, 12.0]x = 1.5
print(Newton(x, x_list, y_list))
x = 2.6
print(Newton(x, x_list, y_list))# 可视化
pre_x_list = []
pre_y_list = []
for i in range(0, 6001, 6):pre_x_list.append(i / 1000)pre_y_list.append(round(Newton(i / 1000, x_list, y_list), 4))scatter = Scatter(init_opts = opts.InitOpts(theme = "macarons"))
scatter.add_xaxis(xaxis_data = np.array(x_list))
scatter.add_yaxis(series_name = "实际值",y_axis = np.array(y_list),label_opts = opts.LabelOpts(position = "outside"))
scatter.set_global_opts(title_opts = opts.TitleOpts("牛顿插值多项式",pos_left = "40%",pos_top = "5%"),legend_opts = opts.LegendOpts(pos_top = "15%"))line = Line(init_opts = opts.InitOpts(theme = "macarons"))
line.add_xaxis(xaxis_data = pre_x_list)
line.add_yaxis(series_name = "牛顿插值函数",y_axis = np.array(pre_y_list),label_opts = opts.LabelOpts(is_show = False),symbol_size = 0,linestyle_opts = opts.LineStyleOpts(width = 2))
line.set_global_opts(xaxis_opts = opts.AxisOpts(is_show = False))
scatter.overlap(line)
scatter.render("牛顿插值.html")

# 观察龙格现象
# 设置观测点
x_list_6 = [-3.5, -2.5, -1.5, 1.5, 2.5, 3.5]
x_list_8 = [-4.5, -3.5, -2.5, -1.5, 1.5, 2.5, 3.5, 4.5]
x_list_10 = [-4.5, -3.5, -2.5, -1.5, -0.5, 0.5, 1.5, 2.5, 3.5, 4.5]y_list_6 = [0.08, 0.14, 0.31, 0.31, 0.14, 0.08]
y_list_8 = [0.05, 0.08, 0.14, 0.31, 0.31, 0.14, 0.08, 0.05]
y_list_10 = [0.05, 0.08, 0.14, 0.31, 0.8, 0.8, 0.31, 0.14, 0.08, 0.05]x_list_all = [x_list_6, x_list_8, x_list_10]
y_list_all = [y_list_6, y_list_8, y_list_10]n_list = [6, 8, 10]

# 插值点数据可视化
scatter = Scatter(init_opts = opts.InitOpts(theme = "macarons"))
scatter.add_xaxis(xaxis_data = np.array(x_list_all[2]))
scatter.add_yaxis(series_name = "实际值",y_axis = np.array(y_list_all[2]),label_opts = opts.LabelOpts(position = "outside"))
scatter.set_global_opts(legend_opts = opts.LegendOpts(pos_top = "15%"),yaxis_opts = opts.AxisOpts(max_ = 2,min_ = -2))# 函数 y
def y(x):return round(1.0 / (1.0 + x ** 2), 4)pre_x_list = []
pre_y_list = []
for i in range(-6000, 6001, 6):pre_x_list.append(i / 1000)pre_y_list.append(y(i / 1000))line = Line(init_opts = opts.InitOpts(theme = "macarons"))
line.add_xaxis(xaxis_data = pre_x_list)
line.add_yaxis(series_name = "函数 y",y_axis = np.array(pre_y_list),label_opts = opts.LabelOpts(is_show = False),symbol_size = 0,linestyle_opts = opts.LineStyleOpts(width = 3))
line.set_global_opts(xaxis_opts = opts.AxisOpts(is_show = False),yaxis_opts = opts.AxisOpts(max_ = 5,min_ = -5))
scatter.overlap(line)

# 拉格朗日龙格现象
for i in range(len(x_list_all)):pre_x_list = []pre_y_list = []for j in range(-6000, 6001, 6):pre_x_list.append(j / 1000)pre_y_list.append(round(Lagrange(j / 1000, x_list_all[i], y_list_all[i]), 4))line = Line(init_opts = opts.InitOpts(theme = "macarons"))line.add_xaxis(xaxis_data = pre_x_list)line.add_yaxis(series_name = f"节点数 n={n_list[i]} 拉格朗日插值函数",y_axis = np.array(pre_y_list),label_opts = opts.LabelOpts(is_show = False),symbol_size = 0,linestyle_opts = opts.LineStyleOpts(width = 3))line.set_global_opts(xaxis_opts = opts.AxisOpts(is_show = False),yaxis_opts = opts.AxisOpts(max_ = 5,min_ = -5))scatter.overlap(line)scatter.set_global_opts(legend_opts = opts.LegendOpts(orient = "vertical",pos_top = "60%"),title_opts = opts.TitleOpts("拉格朗日插值龙格现象",pos_left = "40%",pos_top = "5%"))
scatter.render("拉格朗日插值龙格现象.html")

# 观测点数据可视化
scatter = Scatter(init_opts = opts.InitOpts(theme = "macarons"))
scatter.add_xaxis(xaxis_data = np.array(x_list_all[2]))
scatter.add_yaxis(series_name = "实际值",y_axis = np.array(y_list_all[2]),label_opts = opts.LabelOpts(position = "outside"))
scatter.set_global_opts(legend_opts = opts.LegendOpts(pos_top = "15%"),yaxis_opts = opts.AxisOpts(max_ = 2,min_ = -2))# 函数 y
def y(x):return round(1.0 / (1.0 + x ** 2), 4)pre_x_list = []
pre_y_list = []
for i in range(-6000, 6001, 6):pre_x_list.append(i / 1000)pre_y_list.append(y(i / 1000))line = Line(init_opts = opts.InitOpts(theme = "macarons"))
line.add_xaxis(xaxis_data = pre_x_list)
line.add_yaxis(series_name = "函数 y",y_axis = np.array(pre_y_list),label_opts = opts.LabelOpts(is_show = False),symbol_size = 0,linestyle_opts = opts.LineStyleOpts(width = 3))
line.set_global_opts(xaxis_opts = opts.AxisOpts(is_show = False),yaxis_opts = opts.AxisOpts(max_ = 5,min_ = -5))
scatter.overlap(line)

# 牛顿插值龙格现象
for i in range(len(x_list_all)):pre_x_list = []pre_y_list = []for j in range(-6000, 6001, 6):pre_x_list.append(j / 1000)pre_y_list.append(round(Newton(j / 1000, x_list_all[i], y_list_all[i]), 4))line = Line(init_opts = opts.InitOpts(theme = "macarons"))line.add_xaxis(xaxis_data = pre_x_list)line.add_yaxis(series_name = f"节点数 n={n_list[i]} 牛顿插值函数",y_axis = np.array(pre_y_list),label_opts = opts.LabelOpts(is_show = False),symbol_size = 0,linestyle_opts = opts.LineStyleOpts(width = 3))line.set_global_opts(xaxis_opts = opts.AxisOpts(is_show = False),yaxis_opts = opts.AxisOpts(max_ = 5,min_ = -5))scatter.overlap(line)scatter.set_global_opts(legend_opts = opts.LegendOpts(orient = "vertical",pos_top = "60%"),title_opts = opts.TitleOpts("牛顿插值龙格现象",pos_left = "40%",pos_top = "5%"))
scatter.render("牛顿插值龙格现象.html")

拉格朗日插值和牛顿插值的龙格现象相关推荐

matlab全域基函数,多项式函数插值：全域多项式插值（一）单项式基插值、拉格朗日插值、牛顿插值 [MATLAB]...
全域多项式插值指的是在整个插值区域内形成一个多项式函数作为插值函数.关于多项式插值的基本知识,见"计算基本理论". 在单项式基插值和牛顿插值形成的表达式中,求该表达式在某一点处的值 ...
Matlab实现线性插值、抛物插值、牛顿插值、拉格朗日插值、分段抛物插值、分段线性插值
目录线性插值原理流程图代码抛物插值原理流程图代码拉格朗日插值代码牛顿插值原理代码分段线性插值代码线性插值原理流程图单个点的线性插值代码 X=[0.2 0.4]; ...
【插值】牛顿插值、拉格朗日插值、三次样条插值的Python代码实现
插值简介插值即根据有限的离散点绘制出穿过所有样本点的曲线,从直观上想象似乎画一条穿过n个特定点的曲线有无数种画法,但从数学意义上来说我们希望画出的曲线能够尽量平滑,震荡幅度尽量小能够在非样本点上符合 ...
1月16日：拉格朗日中值定理，罗尔定理，柯西中值，拉格朗日插值，牛顿插值，重心插值，拉格朗日乘子法的证明
拉格朗日中值定理 https://www.bilibili.com/video/BV117411E7kx?from=search&seid=17921778669593975548 拉格朗日中 ...
数值分析之拉格朗日插值、牛顿插值、分段线性插值实现
1.拉格朗日插值法考虑全局信息的比较经典的插值方法,编程简单,计算量大. #coding=utf-8 from matplotlib import pyplot as pltdef Lg(data, ...
插值问题（拉格朗日插值、牛顿插值）
agui_lagrange.m: function f=agui_lagrange(x0,y0,x) % x0为节点向量,y0为节点上的函数值,x为插值点,f为返回插值 n=length(x0);m= ...
Hermite插值是牛顿插值的极限情形
Hermite插值可以看作牛顿插值的极限状况.为什么可以这么说呢?我们来看一个实例: 构造一个三次多项式 $p_3$ 使得 $p_3(0)=0$,$p_3(1)=1,p_3'(0)=1,p_3'(1) ...
插值法：拉格朗日插值、牛顿插值
拉格朗日插值法 (*以下定义选自维基百科) 算法流程图算法代码 [cpp] view plaincopy #include<iostream> #include<string> ...
【数值分析】插值法：拉格朗日插值、牛顿插值
本科课程参见:<软件学院那些课> 拉格朗日插值法 (*以下定义选自维基百科) 算法流程图算法代码 #include<iostream> #include<string& ...

拉格朗日插值和牛顿插值的龙格现象