实验次样条插值实验

1、数值分析实验报告姓名忘川学号1205025106系另U数学系班级12级主讲教师指导教师实验日期2014/6/25专业信息与计算科学专业课程名称数值分析同组实验者无一、实验名称：实验四、三次样条插值实验二、实验目的：1.掌握三次样条插值的运用；2.了解拉格朗日插值在高次上的误差。、实验内容及要求:一、一1一一，_、给正函数f(x)=,x?5,5,下点a=-5+k(k=0,1,L,10),x+1(1) 编写三次样条程序求Si0(x)，取Sio(Xk)=f(xQ(k=0,1,L,10),斗好-5)=f(-5),Sio(5)=f(5).(2) 作原函数f(x)、Langrage插值函数L10(x)和三

2、次样条插值函数§0(x)的图像，并比较它们的区别。附：算法描述：三次样条插值(CubicSplineInterpolation)PURPOSE:TofindapiecewisecubicsplinefunctionMi3Mi13yiMiy1Mi1S(x)=(x+-x)+(x-x)+(h)(xix)+(3hi)(xx).6hi6hih6h6(Xxx1,i=0,1,，n-1).whereSxi=Mii=0,1,n,andS''(x)=Mjx*x+Mi+JxxWxi,xz,i=0,1,，n1.hhiINPUT:interpolatedpoints(xi,yi);clampe

3、dboundarycondition*(冷)=)：,SgNy：.OUTPUT:Mii=0,1,nStep1Fori=1,2",n-1hihihi1hi-0二=1)hihi1_6N1fyi-yii=(一)hihi1hi1hiStep2Byclampedboundarycondition,Setn。=6(yTy。_y。),n=!(yny二加)hih加加Step3Solvetri-diagonalsystem一2ft1-M01飞1菁2RM1%2P2M2役-11-OtnA2PnXMn:.«n2JMn_1InStep4OutputMii=0,1,n-Stop四、实验步骤(或记录)(1

4、)编写三次样条程序求Sio(x)，取Sio(Xk)=f(Xk)(k=0,1,L,10),So¥-5)=f(-5),Si0(5)=f(5).程序如下：functions,m=selfspline(x0,y0,df,x,conds)n=length(x0);h=diff(x0);b=ones(1,n)*2;forj=2:n-1a(j)=h(j-1)/(h(j-1)+h(j);c(j)=h(j)/(h(j-1)+h(j);d(j)=6*(h(j-1)*(y0(j+1)-y0(j)-h(j)*(y0(j)-y0(j-1)/(h(j)*h(j-1)*(h(j)+h(j-1);enda(1:n-

5、2)=a(2:n-1);switchcondscase1a(n-1)=1;c(1)=1;d(1)=6*(y0(2)-y0(1)/(x0(2)-x0(1)-df(1)/(x0(2)-x0(1);d(n)=6*(df(2)-(y0(n)-y0(n-1)/(x0(n)-x0(n-1)/(x0(n)-x0(n-1);case2a(n-1)=0;c=0;d(1)=2*df(1);d(n)=2*df(2);otherwiseerror('conds值错误，conds只能为1或2')endfori=2:nr=a(i-1)/b(i-1);b(i)=b(i)-r*c(i-1);d(i)=d(i)

6、-r*d(i-1);endm(n)=d(n)/b(n);fori=n-1:-1:1m(i)=(d(i)-c(i)*m(i+1)/b(i);endforj=1:length(x)fori=1:n-1ifx(j)>=x0(i)&x(j)<=x0(i+1)s(j)=m(i)*(x0(i+1)-x(j)A3/(6*h(i)+m(i+1)*(x(j)-x0(i)A3/(6*h(i)+(y0(i)-m(i)*h(i)A2/6)*(x0(i+1)-x(j)/h(i)+(y0(i+1)-m(i+1)*h(i)A2/6)*(x(j)-x0(i)/h(i);endendend在matlab的命

7、令窗口输入:>>k=0:10;>>x0=-5+k;>>y0=1./(x.A2+1);>>x=-5:0.5:5;>>y=1./(x.A2+1);>>df=diff(y)得到如下数据：df=0.00860.01180.01660.02450.03790.06210.10770.19230.30000.2000-0.2000-0.3000-0.1923-0.1077-0.0621-0.0379-0.0245-0.0166-0.0118-0.0086(2)作原函数f(x)、Langrage插值函数L10(x)和三次样条插值函数S,

8、0(x)的图像，并比较它们的区别。先写拉格朗日插值的程序，如下：functiony=lagr1(x0,y0,x)n=length(x0);m=length(x);fori=1:mz=x(i);s=0.0;fork=1:np=1.0;forj=1:nifj=kp=p*(z-x0(j)/(x0(k)-x0(j);endends=p*y0(k)+s;endy(i)=s;endend然后直接在matlab的命令窗口输入：x0=-5-4-3-2-1012345;y0=1./(x0.A2+1);x=-5:0.01:5;df=diff(y0);conds=1;y=selfspline(x0,y0,df,x,conds);y1=1./(x.A2+1);y2=lagr1(x0,y0,x);plot(x,y1,'b',x,y,'y-',x,y2,'r:')legend('原函数，三次样条插值，'Lagrange插值)即可得到以下图形：1.50.50原函数三次样条插值Lagrange插值-0