数值分析(15)样条插值课件.ppt

1、 样条是绘图员用于描绘光滑曲线的一种机样条是绘图员用于描绘光滑曲线的一种机械器件械器件,它是一些易弯曲材料制成的窄条或棒条它是一些易弯曲材料制成的窄条或棒条.在绘制需要通过某点的光滑曲线时在绘制需要通过某点的光滑曲线时,对它在这些对它在这些点的位置上点的位置上“压铁压铁”,它就被强制通过或接近图它就被强制通过或接近图表上确定的描绘点表上确定的描绘点.“样条函数样条函数”这个术语意在这个术语意在点出这种函数的图象与机械样条画出的曲线很点出这种函数的图象与机械样条画出的曲线很象象.第五节第五节 样条插值样条插值011(),:()0,1,2,nniyf xa ba baxxxxbf xinx 若若函

2、函数数在在上上连连续续，对对于于区区间间上上的的一一个个分分划划给给定定节节点点上上函函数数值值，。若若函函数数S()S()满满足足21110(1)()0 1(2)(),(3)(),0 11(),.,iiiinnS xyinS xCa bS xx xinS xxxxx ，；，即即在在整整体体上上是是二二阶阶连连续续的的；在在每每一一个个小小区区间间（，）是是三三次次多多项项式式;则则称称为为三三次次样样条条函函数数。称称为为内内节节点点称称为为外外节节点点.一、三次样条插值问题的提法一、三次样条插值问题的提法(1)3(1)42,nnn 共共有有个个条条件件缺缺两两个个条条件件，由由边边界界条条

3、件件给给出出。插插值值条条件件分分析析1,0,1,2,.,iy in 由由（）已已知知节节点点上上函函数数值值。这这是是n n+1 1个个条条件件21112(),()(),()(),()(),1,2,.,1iiiiiiiiiiiiS xCa bSxSxSxS xSxSxin 由由（）隐隐含含着着在在内内节节点点上上应应有有1001112113(),(),;(),;()(),;iinnnS xx xSxxxxSxxxxS xSxxxx 由由（）在在每每个个上上表表达达式式不不同同，故故应应分分段段构构造造：()4S xn故故构构造造需需要要个个条条件件3(1)n 个个条条件件000(0)(0)(

4、30)(0)(0)(0)nnnS xS xSxSxSxSx（）周周期期边边界界条条件件00()(),()1(nnSxfxSxfx常常见见的的边边界界条条边边界界条条件件（）给给定定端端点点一一件件有有以以下下三三种种：称称为为固固定定边边阶阶导导数数值值界界条条件件。000()()()()()0()02nnnSxfxSxfxSxSx ，特特别别，当当（）给给定定端端点点二二阶阶导导数数值值和和，称称为为自自然然边边界界条条件件。()S x构构造造三三次次样样条条插插值值函函数数的的基基本本方方法法(1 1)三三转转角角插插值值法法(2 2)三三弯弯矩矩插插值值法法以以下下面面问问题题为为例例介

5、介绍绍三三转转角角插插值值法法和和三三弯弯矩矩插插值值法法.0,()(0,1,.,)(),()()()(1,2,.,1)(),()()(0,)iiniijjjjyf xinfxfxS xS xyinS xyS xfxjn 对对于于分分划划已已给给相相应应的的函函数数值值以以及及边边界界点点上上的的一一阶阶导导数数值值.求求一一个个三三次次样样条条函函数数使使之之满满足足问问题题二、三转角方程求解法二、三转角方程求解法1()(),(,)iiiiiiiiyfxS xmmymy（）未未知知，但但可可设设三三转转角角法法的的基基本本思思只只是是想想111112,()()()()()1.iiiiiiii

6、iiiiimx xHermiteS xh x yhx yh x mhx mnm （）如如能能求求出出，则则可可在在上上构构造造三三次次插插值值函函数数共共有有个个待待求求4,(0,1,2,.,)imin （）再再由由三三转转角角方方程程边边界界条条件件（补补充充两两个个方方程程）封封闭闭的的方方程程组组，可可求求出出1(3)()(),1,2,.,111iiiiimSxSxinnn 如如何何求求？利利用用在在节节点点上上二二阶阶导导数数连连续续的的条条件件由由导导出出三三转转角角方方程程（个个方方程程要要解解个个未未知知数数）11111(),(,),()()()()()1.iiiiiiiiiii

7、iiiiiiiS xmmymymx xHermiteS xh x yhx yh x mhx m 设设只只是是如如果果已已知知，则则在在每每个个小小区区间间上上，可可构构造造两两点点三三次次插插值值建建立立三三转转角角方方程程11(),(),(),()iiiihxhxhxhxHermite 为为插插值值基基函函数数2211311()2()()12iiiiiiiiiiixxxxxxhxxh xxxxxh2113()(2()()iiiiixxhxxhxh 2211122()()()()()()iiiiiiiixxxxxxxxh xhxhh ，1,iixx 于于是是,在在 上上221113322111

8、22()(2()()(2()()()()()()iiiiiiiiiiiiiiiiiiix xhx xx xhxxS xyyhhx xx xx xx xmmhh 11122113624642()6(2)()iiiiiiiiiiiiiixxxxxxSxmmhhxxxyyh 求求二二阶阶导导数数有有1111122111131624642()6(2)(),iiiiiiiiiiiiiiiixxxxxxSxmmhhxxxyyhxx 同同理理，在在也也可可以以得得到到11,2,-1()()iiiiiinSxSxx （）上上，由由在在内内节节点点11112211111112()3()iiiiiiiiiiiii

9、yyyymmmhhhhhh 111123(,)(1,2,1)iiiiiiiiiiimmmfxxfxxin 化化简简整整理理得得到到三三转转角角方方程程111,1iiiiiiiiihhhhhh 其其中中 0000,()()nnnnmmS xyS xyyy 这这时时，得得到到方方程程组组2 2.加加边边界界条条件件构构造造已已知知两两端端点点一一阶阶导导数数闭闭方方值值和和封封程程组组 1122222211101112102122232322211211122223(,)3(,)3(,)3(,)nnnnnnnnnnnnnnnnnmmmmfxxfxxyfxxfxxfxxfxxfxxfxx 1nny

10、4,(0,1,2,.,)iMin （）再再由由三三弯弯矩矩方方程程边边界界条条件件（补补充充两两个个方方程程）封封闭闭的的方方程程组组，可可求求出出1()(),(,)iiiiiiiiyfxSxMMyMy 三三弯弯矩矩法法的的基基本本思思想想（）未未知知，但但可可设设只只是是2().iiiMMyS x（）如如能能求求出出，则则可可由由和和 构构造造1()(),1,2,.,111iiiiSxS xinnn 可可利利用用在在节节点点上上一一阶阶导导数数连连续续条条件件由由导导出出三三弯弯矩矩方方程程（个个方方程程要要解解个个未未知知数数）(3)iM如如何何求求?三三、三三弯弯矩矩方方程程求求解解法法

11、11(1,),iiiiiiiiiiiS xxa b c dxxMhx在在中中（）建建立立和和之之间间的的关关系系代代入入两两个个端端点点 和和并并记记得得到到132()()()()(0 11),41iiiiiiiiiiiiiixxS xa xxb xxc xxdina b c dn、建建立立三三弯弯矩矩方方程程在在，上上，三三次次样样条条函函数数可可表表示示为为，其其中中是是四四个个待待定定系系数数。（共共有有个个待待求求系系数数）()(0 1)1iiiiiSxMinMynM 设设，共共有有个个待待求求。11323211(),()()()()()1,iiiiiiiiiiiiiiiiiiiiii

12、iiiiiiiiiiiiiS xxxhxxyS xa xxb xxc xxddyS xa hb ha b cc hdMd在在中中代代入入两两个个端端点点 和和并并（）建建立立记记和和之之的的关关系系得得到到间间211()()3()2()()6()2()6()22()62iiiiiiiiiiiiiiiiiiiiiiiiiS xSxa xxb xxcSxa xxbMSxa xxbbMSxa hb 对对求求导导于于是是11111()/6/226,iiiiiiiiiiiiiiiiiiiiiiiiaMMhbMdyyyh Mh MchMMyyabcd 得得到到以以上上得得到到用用和和表表示示系系数数，的的

13、关关系系式式。（2）构造三弯矩方程）构造三弯矩方程2 (),1()3(-)2(-)S xxxiiSxax xbx xciiiiii 利利用用在在内内节节点点上上一一阶阶导导数数连连续续的的条条件件，在在区区间间上上（2）构造三弯矩方程）构造三弯矩方程2 (),1()3(-)2(-)S xxxiiSxax xbx xciiiiii 利利用用在在内内节节点点上上一一阶阶导导数数连连续续的的条条件件，在在区区间间上上2111111111 (-0)(0)()()32,iiiiiiiiiiiiiiiiixSxSxSxSxahb hccabcc 在在节节点点 上上应应有有，即即，得得到到将将前前已已得得到

14、到的的和和 的的表表达达式式代代入入上上式式 有有12111111,()3(-)2(-)iiiiiiiixxSxax xbx xc 在在区区间间上上11-1-1-1-1-1-1-1-1-12-6(-)-2-26iiiiiiiiiiiiiiiiiiiyyh Mh MhMMhyyhMhMMhh 11111126(,)(1,2,1)iiiiiiiiiiiihMhhMh Mf xxf xxgin （+）化化简简后后得得到到三三弯弯矩矩方方程程001101112212221112()2()2()nnnnnnhhhhMghhhhMghhhhMg 11nn 三三弯弯矩矩方方程程组组共共个个方方程程要要解解个

15、个未未知知数数，是是不不封封闭闭的的方方程程组组,要要求求解解必必须须补补充充边边界界条条件件.00()()nnS xyS xy 已已知知两两端端点点一一阶阶导导2 2.加加边边界界条条件件构构造造封封闭闭值值和和方方程程组组 数数 ,012()3(-)2(-)000000-(2)10001()-00006026(,-)0001010 xxSxax xbx xcyyhMMySxchh Mh Mf xxy 在上得补充方程 ,12()3(-)2(-)1111112()32111111nxxnnSxax xbx xcnnnnnnySxahbhcnnnnnnn 在上111126(,)nnnnnnnab

16、chMhMyf xxn-1n-1n-1n-1n-1n-1代代入入系系数数,后后整整理理可可得得补补充充方方程程00010100111111111126(,)2()6(,)(1,2,1)26(,)iiiiiiiiiiiiinnnnnnnnh Mh MfxxyghMhh Mh Mfxxfxxg inhMhMyfxxgM 将将以以上上求求解解的的方方程程合合在在一一起起,得得方方程程组组即即：00001100111122111122()2()2nnnnnnnnnnMghhMghhhhMghhhhhhMg 3200000000,013211111111,12321111111,()()()()(),;

17、()()()(),;()()()()3.()iiiiinnnnnnnnMa b c dS xS xa xxb xxc xxd xx xS xa xxb xxc xxd xx xSxaxxbxxcxdS xx 由由三三次次样样条条插插值值函函数数计计算算，就就可可得得到到分分段段三三次次样样条条插插值值函函数数是是：1,1,;nnxxx 401,()na ba baxxxbxf xC 插插值值区区间间的的 一一个个剖剖分分是是 设设被被插插值值函函数数 S S()是是满满足足第第一一类类或或第第二二类类边边界界条条件件的的三三次次样样条条 插插值值，则则在在插插值值区区间间上上有有估估计计式式

18、()()4(4)()()()(0,1,2)kkkkfxSxc hfxk 10121513(),384248iii nxxccc 0 0其其中中 h h=m ma ax x四、样条插值的收敛性四、样条插值的收敛性-50500.20.40.60.810阶 导 数 精 确 结 果样 条 结 果-505-1-0.500.511阶 导 数 精 确 结 果样 条 结 果-505-2-1.5-1-0.500.512阶 导 数 精 确 结 果样 条 结 果-505-5053阶 导 数 精 确 结 果样 条 结 果Chzh1.m三次样条插值常用命令三次样条插值常用命令(1)y0=interp1(x,y,x0,s

19、pline)y0=interp2(x,y,x0,spline)（2）y0=spline(x,y,x0)(一维一维)（3）pp=spline(x,y);（结构数组）（结构数组）y0=fnval(pp,x0);poly(x0,y0),pp.coefs（系数）（系数）drp=fnder(pp,n),drv=fnval(drp,a)求导数求导数 itp=fnint(pp),itv=fnval(itp,a,b)*-1,1 求积分求积分 fnplt(pp),fnplt(drp,n),fnplt(itp)绘图绘图（4）pp=csape(x,y,边界条件选项边界条件选项);（5）pp=csapi(x,y);演

20、示程序演示程序 sptest.m【例】根据连续时间函数的采样数据，利用【例】根据连续时间函数的采样数据，利用spline重构该连续重构该连续函数，并检查重构误差。函数，并检查重构误差。t=-5:0.5:5;w=exp(-abs(t);N0=length(t);tt=linspace(t(1),t(end),10*N0);ww=spline(t,w,tt);error=max(abs(ww-exp(-abs(tt)plot(tt,ww,b);hold onstem(t,w,filled,r);hold off error=0.0840-50500.10.20.30.40.50.60.70.80.

21、91【例】对于函数【例】对于函数 ，很容易求得，很容易求得 。本例将借此演示样条函数求数值不定积分、导函数。本例将借此演示样条函数求数值不定积分、导函数的能力。的能力。xysin xxxdxxS0cos1sin)(xycos （1）不定积分样条函数、导数样条函数的求取和精度分析。）不定积分样条函数、导数样条函数的求取和精度分析。x=(0:0.1:1)*2*pi;y=sin(x);pp=spline(x,y);int_pp=fnint(pp);der_pp=fnder(pp);%xx=(0:0.01:1)*2*pi;err_yy=max(abs(ppval(pp,xx)-sin(xx)err_i

22、nt=max(abs(ppval(int_pp,xx)-(1-cos(xx)err_der=max(abs(ppval(der_pp,xx)-cos(xx)err_yy=0.0026err_int=0.0010err_der=0.0253（2）不定积分样条函数、导数样条函数的使用不定积分样条函数、导数样条函数的使用%计算计算y(x)在区间在区间1,2上的定积分上的定积分DefiniteIntegral.bySpline=ppval(int_pp,1,2)*-1;1;DefiniteIntegral.byTheory=(1-cos(2)-(1-cos(1);%计算计算dy(3)/dxDeriva

23、tive.bySpline=fnval(der_pp,3);Derivative.byTheory=cos(3);Derivative.byDiference=(sin(3.01)-sin(3)/0.01;DefiniteIntegral,Derivative DefiniteIntegral=bySpline:0.9563 byTheory:0.9564Derivative=bySpline:-0.9895 byTheory:-0.9900 byDiference:-0.9907 （3）绘制三个样条函数的图形）绘制三个样条函数的图形fnplt(pp,b-);hold onfnplt(int_pp,m:),fnplt(der_pp,r-);hold offlegend(y(x),S(x),dy/dx)01234567-1-0.500.511.522.5y(x)S(x)dy/dx习题五习题五 P208-14数值实验题五数值实验题五 P210-4ginput(10)