数值分析——插值(一)、Runge现象
在数值分析中,插值方法是基础且重要的。本文将介绍Lagrange插值公式与Newton插值公式。此外,针对Runge现象,本文给出了稍稍详细的讨论。
一、Lagrange插值公式
假设函数 \(y=f(x)\) 在取定的\(n+1\)个互异的基点 \(x_0,x_1,\cdots,x_n\) 处的值已知分别为 \(y_0=f(x_0),y_1=f(x_1),\cdots,y_n=f(x_n)\),现在要寻找多项式 \(p(x)\) 使得 $$p(x_k)=f(x_k), \quad k=0,1,\cdots,n$$
记
\]
令
\]
这就是满足条件的多项式。
例 利用等距基点计算函数 \(f(x)=\frac{\sin{e^{-x^2}}}{\sqrt{x+2}}\) 在区间 \([-1,1]\) 上的值。
import math
def fun1(t):
return math.sin(math.exp(-t**2))/math.sqrt(t+2)
def LagIn(X,P,x):
Pol=0
length=len(P)
count1=0
while count1 <= length-1:
lpoly=P[count1]
count2=0
while count2 <= length-1:
if count1 != count2:
lpoly=lpoly*(x-X[count2])/(X[count1]-X[count2])
else:
lpoly=lpoly
count2=count2+1
Pol=Pol+lpoly
count1=count1+1
return Pol
n=int(input("输入基点个数:"))
A=[]
K=[]
Result=[]
count=0
while count <= n:
A.append(fun1(-1+count*2/n))
K.append(-1+count*2/n)
count=count+1
count=0
print(LagIn(K,A,0.1))
考虑误差问题。记
\]
若\(f(x)\)在包含\(n+1\)个插值基点\(x_0,x_1,\cdots,x_n\)的区间\([a,b]\)上有\(n\)阶连续导数,且在\((a,b)\)内存在\(n+1\)阶有界导数,那么对\([a,b]\)上的每一点\(x\)必存在一点\(\eta \in (a,b)\)使得
\]
二、均差与Newton插值公式
仍然使用上一节的记号。考虑形如
\]
的多项式,满足条件
\]
由线性代数的知识容易知道,\(a_0,a_1,\cdots,a_n\)是唯一确定的。考虑这种形式的多项式的好处是,当需要增加基点时,可以减少计算量。
下面引入均差的概念来确定该插值多项式的系数。
\(f(x)\)关于基点\(x_0,x_1,\cdots,x_n\)的n阶均差定义为
\]
规定
\]
是\(f(x)\)关于基点\(x_i\)的零阶均差。
利用数学归纳法可以证明:n阶均差可以写成
\]
与利用Lagrange插值法得到的多项式比较可知,\(N_n(x)\)的系数
\]
三、Runge现象
考虑函数
\]
对该函数在等距基点上插值,分别取插值基点个数为5,10,15,20,30,作出误差图像,如下所示:
import math
import matplotlib.pyplot as mpl
def fun1(t):
return (1+t**2)**(-1)
def LagIn(X,P,x):
Pol=0
length=len(P)
count1=0
while count1 <= length-1:
lpoly=P[count1]
count2=0
while count2 <= length-1:
if count1 != count2:
lpoly=lpoly*(x-X[count2])/(X[count1]-X[count2])
else:
lpoly=lpoly
count2=count2+1
Pol=Pol+lpoly
count1=count1+1
return Pol
n=int(input("输入基点个数:"))
A=[]
K=[]
Result=[]
count=0
while count <= n:
A.append(fun1(-5+count*10/n))
K.append(-5+count*10/n)
count=count+1
count=0
B=[]
F=[]
while count <= 10000:
B.append(-5.01+count*10.02/10000)
F.append((fun1(B[count])-LagIn(K,A,B[count])))
count=count+1
mpl.plot(B,F)
mpl.show()
可以看到,插值多项式在端点附近有较大的误差。事实上可以证明
\]
这就是所谓的Runge现象。下面将进行仔细讨论。
在区间\([a,b]\)上,记插值基点为\(\{x_j^{(n)}\},\quad j=0,1,\cdots,n\),插值多项式为\(p_n(x)\)。首先估计误差项\(f(x)-p_n(x)\)。
由Newton插值法,在区间\([a,b]\)上取与插值基点相异的一个基点 \(x=t\),得到
\]
由于\(p_{n+1}(t)=f(t)\),因此
\]
下面记
\]
下面在复平面上考虑问题,插值基点视作实轴上的点,仍保留原来的记号。
由留数定理
\]
其中
\]
C是曲线,使得插值基点在它内部,\(f\)在其内部解析。类似可以求得\(f[z,x_0,x_1,\cdots,x_n]\)的表达式。从而得到
\]
为了进一步估计,需要先证明以下引理:
引理 假设\(\{x_j^{(n)}\}\)是区间\([a,b]\)上的等距基点,定义
\]
则
\]
证明:
\]
由黎曼积分的定义知
\]
引理成立。
考虑曲线
\]
命题 \(\{x_j^{(n)}\}\)是包含在曲线\(C(\rho)\)中的等距插值基点,若\(f\)在\(C(\rho)\)中解析,则\(p_n\)在\(C(\rho'),\rho' < \rho\)上一致收敛于\(f\)。
证明:假设\(z\)在曲线\(C(\rho')\)上,显然\(C(\rho)\)和\(C(\rho')\)是不交的,因此有以下估计:
\]
由于\(|f(\eta)|\)和\(|\eta-z|\)在\(C(\rho)\)上分别由上界和下界,因此
\]
考虑引理中的结果,由\(\lim_{n \to \infty}\sigma_n(z)=\rho'\)及\(\lim_{n \to \infty}\sigma_n(\eta)=\rho\),对\(\epsilon < \frac{|p-p'|}{3}\)能够取到\(n\)使得
\]
从而可以证明一致收敛。
考虑\(f\)在\(C(\rho)\)中并不解析的情况。仅仅讨论一种特殊情况:若只有一个一阶极点\(z^*\),此时
\]
其中\(C^*=C(\rho^*)\)是围绕极点的一条路径。由上面的讨论知,第一项趋于零。利用留数定理,第二项可以写成
\]
其中
\]
是有界的。
\]
此时\(f(z)-p_n(z)\)并不收敛。
上面的讨论并没有给出Runge现象的充要条件。当插值区间过大时,可能会出现Runge现象;缩小插值区间,可以避免Runge现象,因此可以考虑分段插值的方法。此外也可以采用Chebyshev基点进行插值,这就不在本文的讨论范围内了。
参考
[1]Epperson J. On the Runge Example[J]. The American Mathematical Monthly, 1987, 94(4):329-341.
[2]林成森. 数值计算方法[M]. 科学出版社, 2005.
[3]作者本人的实验报告.
最新文章
- ABP文档 - EntityFramework 集成
- 使用navigator.geolocation来获取用户的地理位置信息
- H264与RTP
- bounds的剖析
- 【QT】ui转代码
- 简单几何(判断矩形的位置) UVALive 7070 The E-pang Palace(14广州B)
- Python 正则表达式指南 zz
- C Primer Plus之C预处理器和C库
- Hadoop概念学习系列之Hadoop HA进一步深入(二十八)
- HDFS原理讲解
- [wikioi]乌龟棋
- C# 开发系列(二)
- Java Web - HTML 常用标签和符号
- JQuery官方学习资料(译):CSS
- 系统IO
- ORACLE数据库管理员的职责
- 【IOS学习】【Swift语言】
- dialog提交表单
- linux时间
- aused by: org.springframework.beans.factory.BeanCreationException: Error creating bean with name &#39;sessionFactory&#39; defined in class path resource [applicationContext.xml]: Invocation of init method fai