第二章代数差值_第1页
第二章代数差值_第2页
第二章代数差值_第3页
第二章代数差值_第4页
第二章代数差值_第5页
已阅读5页,还剩127页未读 继续免费阅读

下载本文档

版权说明:本文档由用户提供并上传,收益归属内容提供方,若内容存在侵权,请进行举报或认领

文档简介

第二章代数差值第1页,课件共132页,创作于2023年2月2

某化学反应中,在有限个时刻t(min),测得生成物质量浓度y(10-3g/cm3)的如下数据

那么在时刻t=5分,t=18分时的浓度是多少?ti12346810121416yi4.006.418.018.799.539.8610.3310.4210.5310.61

本章主要讨论利用插值方法寻求函数的近似问题。

用函数来表示变量间的数量关系广泛应用于各学科领域,但在实际问题中,往往是通过实验、观测以及计算等方法,获得函数在一些点上的函数值。如何通过这些离散数据找出函数的一个满足精度要求且便于使用的近似表达式,是经常遇到的问题。例如:第2页,课件共132页,创作于2023年2月§1多项式插值问题

3

设函数y=f(x)在区间[a,b]连续,给定n+1个点已知

f(xk)=yk(k=0,1,…,n),在函数类中寻找一函数φ(x)

作为

f(x)

的近似表达式,使满足这时称

y=f(x)为被插值函数,

φ(x)

称为插值函数,xk称为插值节点,式(2.2)称为插值条件,寻求插值函数φ(x)

的方法称为插值方法.一、插值问题a≤x0<x1<…<xn≤b(2.1)φ(xk)=f(xk)=yk,k=0,1,…,n(2.2)第3页,课件共132页,创作于2023年2月4

在构造插值函数时,函数类的不同选取,P对应各种不同的插值方法,这里我们主要研究函数类是代数多项式,,即所谓的多项式插值问题。

多项式插值,从几何角度看,就是寻求n次代数曲线

y=pn(x)

通过个n+1

点(xk,yk)(k=0,1,…,n)作为

f(x)

的近似(如下图).二、多项式插值问题第4页,课件共132页,创作于2023年2月5(2.4)

Pn

表示所有次数不超过n的多项式函数类,若pn(x)∈Pn

,则

pn(x)=a0+a1x+…+anxn

是由n+1个系数唯一确定的。当满足如下的插值条件时得到关于系数a0、a1、…

、an

的线性方程组pn(xk)=f(xk)=yk,k=0,1,…,n(2.3)第5页,课件共132页,创作于2023年2月6其系数行列式为Vandermonde(范德蒙)行列式只要

i≠j

就有xi≠xj

,因而

V(x0,x1,…,xn)≠0,于是方程组(2.4)有唯一解,也就是说,当节点不重合时,n+1个节点的插值条件能唯一确定一个n次插值多项式,从而有:

定理2.1

给定n+1

个互异节点

x0,x1,…,xn

上的函数值y0,y1,…,yn

,则满足插值条件

pn(xk)=f(xk)=yk,k=0,1,…,n的n次插值多项式pn(x)

存在且唯一。

第6页,课件共132页,创作于2023年2月7这样只要求解方程组拉格朗日插值多项式牛顿插值多项式3.分段线性插值多项式4.三次样条插值多项式便可确定插值多项式pn(x)

,但相对来讲计算较复杂。然而,插值多项式的唯一性保证了无论用什么方法获得满足插值条件的多项式都是同一个多项式pn(x)

,因此可以采用其它更简便的方法来确定多项式。下面就介绍几种常用的方法:(2.4)第7页,课件共132页,创作于2023年2月§2Lagrange插值多项式8的函数值

已知

y=f(x)

在n+1

个点构造n次多项式

pn(x),使得从而得到

f(x)

的近似计算式

第8页,课件共132页,创作于2023年2月一、线性插值(n=1)9求解L1(x)=a1x+a0已知使得f(x)≈L1(x),x∈[x0,x1].根据点斜式得到如果令则称l0(x),l1(x)为一次插值多项式的基函数。这时:并称其为一次Lagrange插值多项式。f(x)≈L1(x)=y0l0(x)+y1l1(x)第9页,课件共132页,创作于2023年2月二、抛物线插值(n=2)10求解L2(x)=a2x2+a1x+a0使得f(x)≈L2(x),x∈[x0,x2].关于二次多项式的构造采用如下方法:令已知并由插值条件得到L2(x)=A(x-x1)(x-x2)+B(x-x0)(x-x2)+C(x-x0)(x-x1)L2(x0)=y0,L2(x1)=y1

,L2(x)=y2第10页,课件共132页,创作于2023年2月11于是得到则有这时:f(x)≈L2(x)=y0l0(x)+y1l1(x)+y2l2(x)…紧凑格式如果令并称其为二次Lagrange插值多项式。第11页,课件共132页,创作于2023年2月12另外,如果再引进记号

则其导数为

而且

于是

第12页,课件共132页,创作于2023年2月13这样,就得到二次拉格朗日插值多项式的三种表示形式…紧凑格式

这样就得到在区间[a,b]上关于

f(x)的近似计算式…基函数表示

………………ω3(x)表示式

下面给出n次拉格朗日插值多项式的构造。第13页,课件共132页,创作于2023年2月三、n次Lagrange插值多项式14已知n+1组离散数据按照二次Lagrange插值多项式的构造方法,令:将插值条件Ln(x0)=y0

代入,得到:同理,由插值条件Ln(x1)=y1

,得到:第14页,课件共132页,创作于2023年2月15由插值条件Ln(xn)=yn,得将代入下式:得到:第15页,课件共132页,创作于2023年2月16也就是第16页,课件共132页,创作于2023年2月17…紧凑格式

…基函数表示

引进基函数则有第17页,课件共132页,创作于2023年2月18还有一种表示式……………ωn+1(x)

表示式

其中ωn+1(x)=(x-x0)(x-x1)…(x-xn).

从而就得到了在区间[a,b]上的关于函数f(x)的n次多项式近似计算式:

关于这样的近似计算,需要在计算机上进行,一般我们利用紧凑格式设计程序,程序的流程图如下:

第18页,课件共132页,创作于2023年2月19开始输入

y=0,j=0

t=1

j=n?N

j=

j+1输出yY结束第19页,课件共132页,创作于2023年2月20

例:某化学反应中,在有限个时刻t(min),测得生成物质量浓度y(10-3g/cm3)的如下数据

那么在时刻t=5分,t=16.4分时的浓度是多少?ti12346810121416yi4.006.418.018.799.539.8610.3310.4210.5310.61解:编写程序进行计算(1).先编写Lagrange插值算法子程序(2).再编写主程序Lagrange2调用子程序。(3).计算在t=5,t=16.4的浓度值,并画出曲线图。第20页,课件共132页,创作于2023年2月21%lagrange.mfunctiony=lagrange(x0,y0,x)n=length(x0);m=length(x);fork=1:mz=x(k);s=0.0;forj=1:np=1.0;

fori=1:n

ifi~=jp=p*(z-x0(i))/(x0(j)-x0(i));

end

ends=p*y0(j)+s;endy(k)=s;end%lagrange2.mx0=[12 346810121416];y0=[4.006.418.018.799.539.8610.3310.4210.5310.61];x=[1:0.1:16];y=lagrange(x0,y0,x);x1=5;x2=16.4;y1=lagrange(x0,y0,x1)y2=lagrange(x0,y0,x2)plot(x0,y0,'.k','markersize',20)holdonplot(x,y,'-r','markersize',30)holdonplot(x1,y1,'*b','markersize',8)holdonplot(x2,y2,'*b','markersize',8)legend('原数值点','Lagrange曲线',2)第21页,课件共132页,创作于2023年2月22计算结果:y(5)=9.2300,y(16.4)=6.5872第22页,课件共132页,创作于2023年2月四、插值余项(误差估计)23

定理2.2:设

f(x)∈Cn[a,b],f(n+1)(x)

在(a,b)内存在,x0,x1,…,xn

[a,b]

中的n+1个互异节点,则当x∈[a,b]

时,n

次Lagrange插值多项式的截断误差为:(2.5)

证明:令x

是[a,b]中任意固定的数,如果

x是插值节点xi,则(2.5)左右两端均为零,等式(2.5)成立。如果x不是节点

xi

,作辅助函数如下:其中,ωn+1(x)=(x-x0)(x-x1)…(x-xn).第23页,课件共132页,创作于2023年2月24其中ωn+1(x)=(x-x0)(x-x1)…(x-xn).可以检验同理可得而且可知

(t)在区间[a,b]内有n+2个零点:x0、x1、…、xn

、x.

由Rolle定理可知,

(t)的导函数

’(t)在(a,b)内有n+1个零点:再利用Rolle中值定理可知

’’(t)

在(a,b)内有n个零点。以此类推,可知

(n+1)(t)在(a,b)内至少有1个零点。第24页,课件共132页,创作于2023年2月25即而根据求得于是从而得到误差估计式:第25页,课件共132页,创作于2023年2月26对于误差估计式如果存在,则可以估计误差限:当n=1时第26页,课件共132页,创作于2023年2月27于是,得到如下Lagrange插值多项式及其误差估计这里f(x)=Ln(x)+Rn(x)当f(x)≈Ln(x)

时,误差为Rn(x)。第27页,课件共132页,创作于2023年2月28例2.1已知分别用线性插值和二次插值计算sin0.3367.解:设(1).取x0,x1

作线性插值于是关于误差,由得到:第28页,课件共132页,创作于2023年2月29(2).取x0,x1,x2

作二次插值得到关于误差,由得到第29页,课件共132页,创作于2023年2月本节(§2)要点30掌握Lagrange插值多项式的构造方法及具体结构掌握Lagrange插值多项式误差分析方法和结果3.编写Lagrange插值多项式计算程序进行实际计算练习:已知函数y=f(x)的如下离散数据x012y2312(1).试用线性插值求函数在x=1.5处的近似值。(2).试用二次插值求函数在x=1.5处的近似值。第30页,课件共132页,创作于2023年2月§3差商及Newton插值多项式一、差商及其性质31

Lagrange

插值多项式的优点是格式整齐规范,但其缺点是:当需要增加节点时,其基函数都要发生变化,需要重新计算,这在实际计算中会影响效率。下面介绍的Newton插值法会弥补这一不足。1.差商的定义设y=f(x)在n+1个互异点

x0,x1,…,xn处的函数值为:则称为f(x)关于点xi,xj

的一阶差商。第31页,课件共132页,创作于2023年2月32并称为f(x)关于点xi,xj,xk

的二阶差商。一般地,称k-1

阶差商的一阶差商为f(x)关于点x0,x1,…xk的k

阶差商。第32页,课件共132页,创作于2023年2月33例如,已知f(x)在的函数值为:可以求得第33页,课件共132页,创作于2023年2月342.差商的性质性质1:k阶差商是由函数值的线性组合而成,即其中证明:以k=2进行证明。由得到第34页,课件共132页,创作于2023年2月35由得到从而第35页,课件共132页,创作于2023年2月36证明:以k=1

为例性质2:差商具有对称性,即k阶差商f[x0,

x1,…,xk-1,xk]

中,任意调换xi,xj

的次序,其值不变。以k=2为例第36页,课件共132页,创作于2023年2月37性质3:若f(x)为一n次多项式,则f[x,x0]为关于x的n-1次多项式。证明:已知由于故类似的可以得到:也就是说,对多项式求一次差商,次数降低一次。第37页,课件共132页,创作于2023年2月383.差商的计算

为构造Newton插值多项式方便起见,计算差商时,采用列表的方式进行。

第38页,课件共132页,创作于2023年2月39

例2.2

已知函数y=f(x)

的如下离散数据(1,0)、(2,2)、(4,12)、(5,20)、(6,70),试求其各阶差商.

解:列差商表计算xy一阶差商二阶差商三阶差商四阶差商1022412520670

2

5

8

50

1

1

21

0

5

1第39页,课件共132页,创作于2023年2月二、Newton插值多项式40

对于区间[a,b]内的离散点及相应的函数值,计算如下差商:可以求得:第40页,课件共132页,创作于2023年2月41依次将下式带入上式得到:第41页,课件共132页,创作于2023年2月42则可以将函数f(x)

表示成:令:容易验证因此Nn(x)满足插值条件,是一个n次插值多项式。第42页,课件共132页,创作于2023年2月43并称为n次Newton插值多项式。如果f(x)≈Nn(x)

则误差为:由Newton插值多项式的结构可以看出,在构造Newton插值多项式时,必须首先计算各阶差商。关于Newton插值多项式,有以下几个特点:第43页,课件共132页,创作于2023年2月441.Newton插值多项式与Lagrange插值多项式的误差相同

这是因为Newton插值多项式与Lagrange插值多项式满足相同的插值条件因此,Newton插值多项式与Lagrange插值多项式的误差相同。这样,由得到这个表达式给出了n+1

阶差商与n+1

阶导数之间的关系式。

Nn(xi)=f(xi),i=0,1,…,n,Ln(xi)=f(xi),i=0,1,…,n为什么?第44页,课件共132页,创作于2023年2月45例2.3已知,试求其如下差商解:由差商与导数的关系式得到第45页,课件共132页,创作于2023年2月462.Newton插值多项式具有递推式由Newton插值多项式可知所以,具有递推公式:由此可知:当求出n-1次插值多项式后,再增加一个节点时,只需要增加一项的计算即可。第46页,课件共132页,创作于2023年2月47

例2.3已知f(x)

的五组数据(1.0)、(2,2)、(3,12)、(4,42)、(5,116),求N4(x)。如果再增加一个节点(6,282),求出N5(x),并计算N4(1.5)、N5(1.5).解:先由前五组数据列差商表10223124425116210307441022240.5得到:如果,再增加一点(6,282),就在上表中增加一行计算差商。628216646810.1第47页,课件共132页,创作于2023年2月48由Newton公式的递推式得到:得到:练习题:已知离散数据(1,0)、(2,2)、(4,12)、(5,20)

求3次牛顿插值多项式,增加一点(6,70)后,再求出4次牛顿插值多项式。第48页,课件共132页,创作于2023年2月本节(§3)要点491.掌握差商及其性质,导数与差商的关系2.掌握Newton

插值多项式的构造方法及具体结构3.掌握Newton插值多项式的误差结果4.编写Newton插值多项式计算程序进行实际计算思考题:如何实现差商表和Newton插值多项式的程序设计。

第49页,课件共132页,创作于2023年2月50for(j=0;j<=n-1;j++)for(i=n;i>=j+1-1;i--)y0[i]=(y0[i]-y0[i-1])/(x0[i]-x0[i-j]);第一步:计算差商Newton插值多项式程序设计xjyj一阶差商二阶差商三阶差商F[j]x0y0F[0]x1y1f[x0,x1]F[1]x2y2f[x1,x2]f[x0,x1,x2]F[2]x3y3f[x2,x3]f[x1,x2,x3]f[x0,x1,x2,,x3]F[3]fori=1:n-1forj=n:-1:i+1y0(j)=(y0(j)-y0(j-1))/(x0(j)-x0(j-i));endend第50页,课件共132页,创作于2023年2月51for(k=0;k<=m;k++){z[k]=y0[0];for(j=1;j<=n;j++){t=1.0;for(i=0;i<=j-1;i++)t=t*(x[k]-x0[i]);z[k]=z[k]+t*y0[j];}第二步:计算m个点的函数值fork=1:mz(k)=y0(1);fori=2:nt=1.0;forj=1:1:i-1t=t*(x(k)-x0(j));endz(k)=z(k)+t*y0(i);endendfor(j=0;j<=n-1;j++)for(i=n;i>=j+1-1;i--)y0[i]=(y0[i]-y0[i-1])/(x0[i]-x0[i-i]);第51页,课件共132页,创作于2023年2月52关于离散数据:构造了lagrange插值多项式:Newton插值多项式:第52页,课件共132页,创作于2023年2月§4分段插值多项式4.1高次插值的龙格现象53

根据插值条件构造的插值多项式作为函数的近似,有很大局限性,插值多项式的次数随着节点个数的增多而升高,但高次插值多项式的近似效果有时并不理想,看下面的例子。取等距节点对于求出可构造Lagrange插值多项式Ln(x).第53页,课件共132页,创作于2023年2月54

n=10时,插值多项式

L10(x)

的计算结果见下图:第54页,课件共132页,创作于2023年2月55

从图中可见,在x=0附近L10(x)

对f(x)

有较好的近似,而点x

距零点越远,近似效果越差,以至于完全失真。实际上,当n→∞

时,在|x|<0.36….

范围内,L10(x)收敛于f(x),而在这个区间外,L10(x)

是不收敛的。这个现象被称为Runge(龙格)现象。Runge现象表明,为减小近似误差,盲目地提高插值多项式次数是不可取的,实际上,很少采用高于7次的插值多项式。因此,只能通过缩小插值区间的办法达到减小误差的目的,这就是下面要讨论的低次分段插值多项式。第55页,课件共132页,创作于2023年2月4.2分段线性插值56为了提高近似程度,可以考虑用分段线性插值来逼近原函数。如下图所示:设

y=f(x)

在节点a=x0<x1<…<xn=b

处的函数值为yi=f(xi),i=0,1,2,…,n.第56页,课件共132页,创作于2023年2月57

这时的插值函数为分段函数:

在区间[xi-1,xi]上的线性函数为误差为:第57页,课件共132页,创作于2023年2月58

易见,S(x)是平面上以点(xi,yi)(i=0,1,2,..,n)

为节点的折线。有如下的特点:(1)S(x)

在[xi-1,xi]上为次数不超过一次的多项式。

若f(x)C2[a,b]

,由线性插值的误差公式:(2)S(x)∈C[a,b];

(3)S(x)∈C1[xi-1,xi];得到第58页,课件共132页,创作于2023年2月59则有:令可以按如下的方式考虑:关于整体误差:第59页,课件共132页,创作于2023年2月60可见,当

h→0

时,分段线性插值S(x)

收敛于f(x)

若记

则对任一

x∈[a,b]都有

于是,我们再考虑在节点处可导的插值多项式的构造。

值的注意的是:分段线性插值函数S(x)虽然有很好的收敛性质,且在整个插值区间[a,b]上是连续的,但插值函数S(x)在节点处却不光滑。也就是说,S(x)在节点处的一阶导数不一定存在(或不连续)!第60页,课件共132页,创作于2023年2月4.3分段Hermite插值

61

分段线性插值多项式S(x),在插值区间[a,b]上只能保证连续性,而不光滑。要想得到在插值区间上光滑的分段插值多项式,可采用分段Hermite插值。

如果已知函数y=f(x)

在节点a=x0<x1<…<xn=b

处的函数值和导数值:则在小区间[xi-1,xi]

上有四个插值条件:故能构造一个三次多项式Hi(x),我们称其为三次埃尔米特(Hermite)插值多项式。

yk=f(xk),yk’=f’(xk),k=0,1,…,nyi-1=f(xi-1),yi=f(xi),y’i-1=f’(xi-1),yi’=f’(xi),第61页,课件共132页,创作于2023年2月62其中Hi(x),x∈[xi-1,xi]满足条件:这时,在整个[a,b]上可以用分段三次Hermite插值多项式来逼近f(x)

。第62页,课件共132页,创作于2023年2月63关于

Hi(x)

的构造,我们可以通过基函数来进行:由插值条件:可以给出基函数满足的条件为:第63页,课件共132页,创作于2023年2月64由于

i-1(x)满足条件:我们称

为三次Hermite插值基函数。所以

i-1(x)

应该具有以下形式:下面来求基函数:再将代入可得而且第64页,课件共132页,创作于2023年2月65将得到类似地有:代入下式:第65页,课件共132页,创作于2023年2月66代入下式得到这样,便求出了分段三次Hermite插值多项式:第66页,课件共132页,创作于2023年2月67

【定理】若f(x)∈C3[a,b],且f(4)(x)在(a,b)内存在,则在[xi-1,xi

]

上有其中Proof:1)当x为xi-1或xi时,上式显然成立;2)当x不为xi-1,xi时,作辅助函数则有关于误差,有如下定理:第67页,课件共132页,创作于2023年2月68

如果记则有:即反复应用Rolle中值定理可知在内至少有一个零点,使得,而将t=

i

代入上式可得结论。定理证毕。第68页,课件共132页,创作于2023年2月69关于整体误差,若

f(x)∈C4[a,b],则可按如下方式考虑:记则有:

于是,当h→0

时,R(x)→0。说明分段三次Hermite插值H(x)

收敛于f(x)

。第69页,课件共132页,创作于2023年2月70

上节中的分段低次插值有很好的收敛性,但其光滑性不够理想,为了提高光滑性,就得增加节点上的导数插值条件,这显然是困难的。

由H(x)的表达式,可以的知它有如下特性:(2)H(x)

C1[a,b];(3)H(x)

C2[xi-1,xi].

(1)H(x)

在[a,b]上是分段函数,且在每个[xi-1,xi]是一个三次多项式;比分段线性插值的光滑性要好一些。

理想的结果是:构造一个多项式,在不增加已知条件的情况下,既有好的逼近效果,又有好的光滑性。下一节讨论这样的问题。第70页,课件共132页,创作于2023年2月

本节(§4)问题712.分段线性插值有何优缺点?如何估计误差?

4.如何分段线性插值算法的程序设计?1.何为高次插值的Runge现象,应如何避免?

3.分段三次Hermite插值有何优缺点,如何估计误差5.如何构造满足以下条件的插值多项式并估计误差?第71页,课件共132页,创作于2023年2月§5三次样条插值72

在解决实际问题时,有时对插值函数曲线的光滑性有较高的要求,如机翼形线、船体放样型值线、汽车外形型值线等,都要求有二阶的光滑度,即具有二阶连续导函数。

2)分段线性插值函数,逼近程度较好,但光滑性差;

然而,1)高阶插值曲线(Lagrange、Newton插值曲线)虽然具有较好的光滑性,但会出现Runge现象;3)分段三次Hermite插值,逼近程度好,光滑性也比较好,但需增加更多的插值条件,有些场合不实用,且节点处的二阶导数也不一定连续。

第72页,课件共132页,创作于2023年2月73

理想的结果是:在不增加更多条件的情况下,构造出的插值多项式,既不出现Runge现象,也有很好的逼近程度,同时还有很好的光滑性。下面讨论的样条函数就具有这样的特性!

所谓三次样条一词来源于一种绘图工具:样条尺,柔软富有弹性,作图时可用它作出不同弯曲程度的光滑曲线。在实际工作中,技术人员利用这种样条尺来绘制光滑曲线,样条尺是用易弯曲的材料做成,如易弯曲的木条等。

在样条尺上加压铁固定,可以强制它通过平面图上的固定点(xi,yi),i=0,1,2,…,n,然后绘图人员再沿着这条变形的样条尺画出曲线,这种曲线称为样条曲线,样条曲线具有很好的光滑性。第73页,课件共132页,创作于2023年2月74

当我们用直尺作曲线时,得到左图的形状;当用样条尺作图时,会得到右图的形状,具有好的光滑性。

下面,我们从理论上对这一曲线进行研究,并着重给出三次样条函数的构造方法。第74页,课件共132页,创作于2023年2月5.1样条(Spline)函数的定义75

已知函数

y=f(x)

在节点

a=x0<x1<…<xn=b

处的函数值为如果函数S(x)

满足条件:

(1)S(x)是一个分段的m次多项式,且S(xk)=yk

;则称

S(x)

是m次样条插值函数。(2)S(x)

在[a,b]具有m-1阶连续导数。1.m次样条函数的定义yk=f(xk)

,k=0,1,2,…,n第75页,课件共132页,创作于2023年2月76

2.

三次样条函数的定义

已知函数y=f(x)

在节点a=x0<x1<…<xn=b

处的函数值为如果函数S(x)

满足条件:

(1)S(x)

是一个分段的三次多项式,且S(xk)=yk;则称S(x)

是三次样条插值函数。(2)S(x)

在[a,b]具有二阶连续导数。yk=f(xk)

,k=0,1,2,…,n第76页,课件共132页,创作于2023年2月77

S(x)的具体形式为:

其中Si

在上[xi-1,xi]

是三次多项式

有四个待定常数。这样S(x)

共有4n个待定常数要确定。

根据已知条件及三次样条函数的性质,一共有已知条件。(n+1)+(n-1)+(n-1)+(n-1)=4n-2也就是说,需要4n个条件来确定这些系数。还缺两个条件,一般通过增加边界条件来补充。第77页,课件共132页,创作于2023年2月78常用的边界条件有以下几种:边界条件1:

边界条件2:

边界条件3:假定函数y=f(x)是以b-a

为周期的周期函数,这时要求S(x)

也是周期函数,即

这样我们便可以求出唯一的三次样条函数来。第78页,课件共132页,创作于2023年2月5.2样条函数的求解1.三转角算法79由n+1组数据

只要求出在区间[xi-1,xi]上的三次多项式

Si(x)

i=1,2,…,n即可。加边界条件构成的三次样条函数为分段函数a=x0<x1<…<xn=b,yk=f(xk)

,k=0,1,2,…,n第79页,课件共132页,创作于2023年2月80已知Si(x)

满足条件

由于

Si(x)

是三次多项式,还缺少两个条件。在这里假设

Si(x)

在两个端点的导数值已知:(2.1)(2.2)(2.1)、(2.2)

满足三次Hermite插值多项式的条件,故可以得到:第80页,课件共132页,创作于2023年2月81上式给出了Si(x)

的表达式,但是表达式中的mi-1和mi还是未知的,还需要求出mi-1、mi

来才具体求出了三次样条函数。

这时,根据三次样条函数的定义:S(x)

在[a,b]上具有二阶连续导函数,可知在节点

xi处第81页,课件共132页,创作于2023年2月82而:于是,由两边求导得:Si(x)Si+1(x)第82页,课件共132页,创作于2023年2月83于是同理可得第83页,课件共132页,创作于2023年2月84根据二阶导数连续性条件

,即及

得到

第84页,课件共132页,创作于2023年2月85用

除等式两侧得到将等式

合并整理得:第85页,课件共132页,创作于2023年2月86令:简化为这是一个方程组,可以进一步写为:第86页,课件共132页,创作于2023年2月87

式(2.3)给出了n+1个未知数

m0,m1,…,mn

的n-1个方程,若要求解出来,还需要补充两个方程,这时可以考虑所给定的边界条件。(2.3)第87页,课件共132页,创作于2023年2月88得到这样,我们只需要求出

m0,m1,…,mn

即可。边界条件1:这时的方程改写为:第88页,课件共132页,创作于2023年2月89或者,表示为矩阵方程(2.4)由于λk+μk=1,方程(2.4)的系数矩阵是严格对角占优矩阵,方程(2.4)有惟一解,并可用追赶法求解。第89页,课件共132页,创作于2023年2月90便得到三次样条函数:解出m0,m1,…,mn

以后,代入下式第90页,课件共132页,创作于2023年2月91也就是:边界条件2:

已知:故可以得到

第91页,课件共132页,创作于2023年2月92将与前面得到的方程组:结合构成一个含n+1个未知量和n+1个方程的方程组,便可求解出m0,m1,…,mn

:第92页,课件共132页,创作于2023年2月93求解此方程组,也可以求出三次样条函数。(2.5)或者表示为矩阵方程第93页,课件共132页,创作于2023年2月94由第一个等式和第二个等式得到

y0=yn,m0=mn边界条件3:周期边界条件由第三个等式说明:再由第94页,课件共132页,创作于2023年2月95即:得到由m0=mn

整理为第95页,课件共132页,创作于2023年2月96在等式两边同除以得到令得到第96页,课件共132页,创作于2023年2月97结合,与将得第97页,课件共132页,创作于2023年2月98并将

表示为矩阵方程:(2.6)

其系数矩阵称作周期三对角矩阵,也是严格对角占优,因而方程组有唯一解。第98页,课件共132页,创作于2023年2月三次样条函数三转角算法的实现流程99Step1:输入节点x0,x1,…,xn

,函数值y0,y1,…,yn、

边界条件及

x.Step3:根据边界条件,求解相应的方程得到

m0,m1,…,mnStep2:计算Step4:判断

x∈[xi-1,xi]

?Step5:计算

y≈si(x)Step6:输出

y第99页,课件共132页,创作于2023年2月100三次样条函数曲线图第100页,课件共132页,创作于2023年2月101设

例2-4已知函数y=f(x)

的如下数据,试求其在区间[0,3]上的三次样条插值函数S(x)。

解这里边界条件是求得第101页,课件共132页,创作于2023年2月102已知由方程组及得到方程组解得第102页,课件共132页,创作于2023年2月103这样便求得代入表达式便得到所求的三次样条函数第103页,课件共132页,创作于2023年2月本节(§5-1、2)要点104什么是三次样条函数?三次样条函数的边界条件是如何给出的?三次样条函数的三转角算法是如何构造的?如何设计程序实现三转角算法?画出流程图。三对角线性方程组和一般线性方程组如何求解?完成P49习题2-12、2-13。第104页,课件共132页,创作于2023年2月2.三弯矩算法105

只要求出在区间[xi-1,xi]上的三次多项式

Si(x)

i=1,2,…,n

即可。在这里我们采用另一种方法进行求解,并称为三弯矩算法。加边界条件构成的三次样条函数为分段函数a=x0<x1<…<xn=b,yk=f(xk)

,k=0,1,2,…,n对于离散点及其上的函数值:第105页,课件共132页,创作于2023年2月106已知Si(x)

满足条件

Si(xi-1)=yi-1,Si(xi)=yi(2.1)在此假设

Si”(xi-1)

=Mi-1,Si”(xi)

=Mi

(2.2)对于三次多项式Si(x)

,其二阶导函数Si(xi-1)为线性函数,故可设积分两次得到(2.3)第106页,课件共132页,创作于2023年2月107可以求得:再由条件

Si(xi-1)=yi-1,Si(xi)=yi(2.1)(2.3)代入(2.3)得到:再根据:Si(x)Si+1(x)及第107页,课件共132页,创作于2023年2月108这时于是根据得到也就是第108页,课件共132页,创作于2023年2月109两边同除以上式得令则有这是含有n-1个方程的线性方程组,具体形式为:第109页,课件共132页,创作于2023年2月110

该线性方程组含有n+1个未知量,要解出唯一的一组解来,还需要补充两个方程,这可以有边界条件来实现。第110页,课件共132页,创作于2023年2月111边界条件1:

也就是由得即第111页,课件共132页,创作于2023年2月112再由得即结合可以得到求解M0、M1、…、Mn的线性方程组。第112页,课件共132页,创作于2023年2月113其中方程组(2.6)的矩阵形式为第113页,课件共132页,创作于2023年2月114边界条件2:

这说明这时方程组直接转化为恰定方程组第114页,课件共132页,创作于2023年2月115矩阵形式为(2.9)第115页,课件共132页,创作于2023年2月116边界条件3:

由得第116页,课件共132页,创作于2023年2月117代入得令第117页,课件共132页,创作于2023年2月118则结合得第118页,课件共132页,创作于2023年2月119矩阵形式为第119页,课件共132页,创作于2023年2月120针对三种边界条件求解相应的方程组(2.9)第120页,课件共132页,创作于2023年2月121得到M0、M1、…、Mn后代入就可以得到三次样条函数并把以上算法称作三弯矩算法。第121页,课件共132页,创作于2023年2月122

例2-5

求三次样条插值函数S(x),满足自然边界条件(S”(x0)=S”(xn)=0),已知离散数据如下:

i01234xi0.250.300.390.450.53yi0.5000.54770.62450.67080.2780

解:已知M0=M4=0,先求出步长h1=0.05,h2=0.09,h3=0.06,h4=0.08,再计算第122页,课件共132页,创作于2023年2月123代入方程组得到解得M0=0,M1=-1.8806,M2=-0.8836,

M3=-1.0261,M4=0.(2.9)第123页,课件共132页,创作于2023年2月124代入得到第124页,课件共132页,创作于2023年2月5.3三次样条插值函数的误差估计1.如果f(x)C[a,b],且划分的网格比125其中一致有界,则当时,S(x)一致收敛于f(x).2.如果f(x)C4[a,b],且S(x)满足边界条件I、II,则即h0时,如果划分比

一致有界,第125页,课件共132页,创作于2023年2月5.4三次样条插值函数的程序设计1.三对角方程求解的追赶法126#include<stdio.h>#include<math.h>intcetrd(doubleb[],intn,intm,doubled[]){intk,j;doubles;

if(m!=(3*n-2)){printf("error\n");return(-2);}

for(k=0;k<=n-2;k++){j=3*k;s=b[j];

if(fabs(s)+1.0==1.0)

{printf("fail\n");return(0);}b[j+1]=b[j+1]/s;d[k]=d[k]/s;b[j+3]=b[j+3]-b[j+2]*b[j+1];d[k+1

温馨提示

  • 1. 本站所有资源如无特殊说明,都需要本地电脑安装OFFICE2007和PDF阅读器。图纸软件为CAD,CAXA,PROE,UG,SolidWorks等.压缩文件请下载最新的WinRAR软件解压。
  • 2. 本站的文档不包含任何第三方提供的附件图纸等,如果需要附件,请联系上传者。文件的所有权益归上传用户所有。
  • 3. 本站RAR压缩包中若带图纸,网页内容里面会有图纸预览,若没有图纸预览就没有图纸。
  • 4. 未经权益所有人同意不得将文件中的内容挪作商业或盈利用途。
  • 5. 人人文库网仅提供信息存储空间,仅对用户上传内容的表现方式做保护处理,对用户上传分享的文档内容本身不做任何修改或编辑,并不能对任何下载内容负责。
  • 6. 下载文件中如有侵权或不适当内容,请与我们联系,我们立即纠正。
  • 7. 本站不保证下载资源的准确性、安全性和完整性, 同时也不承担用户因使用这些下载资源对自己和他人造成任何形式的伤害或损失。

评论

0/150

提交评论