实验次样条插值实验

作者：香雪oοО | 发布时间：2020-12-07 18:57:39 收藏本文下载本文

数值分析实验报告姓名忘川学号 1205025106 系别数学系班级 12 级主讲教师指导教师实验日期 2014/6/25 专业信息与计算科学专业课程名称数值分析同组实验者无一、实验名称：

实验四、三次样条插值实验二、实验目的：

1．掌握三次样条插值的运用； 2．了解拉格朗日插值在高次上的误差。

三、实验内容及要求：

给定函数21(), [ 5,5],1f x xx= ?+节点 5(0,1, ,10),kx k k =-+ = L(1)编写三次样条程序求10()S x，取10 10 10()()(0,1, ,10),(5)(5),(5)(5).k kS x f x k S f S fⅱⅱ = =-=-= L(2)作原函数()f x、Langrage 插值函数10()L x 和三次样条插值函数10()S x 的图像，并比较它们的区别。

附：

算法描述：

三次样条插值（Cubic Spline Interpolation）PURPOSE: To find a piecewise cubic spline function 3 31 1 11()(-)(-)()()()()6 6 6 6i i i i i ii i i i i ii i i iM M y M y MS x x x x x h x x h x xh h h h         .(1 i ix x x    ,0,1 i ,… 1 ,  n).where     n i M x Si i, , 1 , 0      , and ""11()i ii ii ix x x xS x M Mh h   1[ , ],i ix x x  0,1 i ,… 1 ,  n.INPUT: interpolated points  i iy x ,;clamped boundary condition0 0(),()i n nS x y S x y      .OUTPUT:   n i M i , , 1 , 0   Step 1 For 1 , , 2 , 1   n i  1111 11 1(1)6()iii iii i ii ii i i iii i i ihh hhh hy y y yh h h h        

Step 2 By clamped boundary condition, Set 01n    ,1 00 01 16()y yyh h  ,16()n nn nn ny yyh h   Step 3 Solve tri-diagonal system 0 0 01 1 1 12 2 2 21 1 1 122222n n n nn n nMMMMM                                                    Step 4 Output   n i M i , , 1 , 0  .Stop 四、实验步骤（或记录）(1)编写三次样条程序求10()S x，取10 10 10()()(0,1, ,10),(5)(5),(5)(5).k kS x f x k S f S fⅱⅱ = =-=-= L 程序如下：

function [s,m]=selfspline(x0,y0,df,x,conds)n=length(x0);h=diff(x0);b=ones(1,n)*2;for j=2:n-1 a(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)));end a(1:n-2)=a(2:n-1);switch conds case 1 a(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));case 2 a(n-1)=0;c(1)=0;d(1)=2*df(1);d(n)=2*df(2);otherwise error("conds值错误，conds只能为1或2")end for i=2:n r=a(i-1)/b(i-1);

b(i)=b(i)-r*c(i-1);d(i)=d(i)-r*d(i-1);end m(n)=d(n)/b(n);for i=n-1:-1:1 m(i)=(d(i)-c(i)*m(i+1))/b(i);end for j=1:length(x)for i=1:n-1 if x(j)>=x0(i)&x(j)<=x0(i+1)s(j)=m(i)*(x0(i+1)-x(j))^3/(6*h(i))+m(i+1)*(x(j)-x0(i))^3/(6*h(i))+(y0(i)-m(i)*h(i)^2/6)*(x0(i+1)-x(j))/h(i)+(y0(i+1)-m(i+1)*h(i)^2/6)*(x(j)-x0(i))/h(i);end end end 在 matlab 的命令窗口输入：

>> k=0:10;>> x0=-5+k;>> y0=1./(x.^2+1);>> x=-5：0.5:5;>> y=1./(x.^2+1);>> df=diff(y)得到如下数据：

df = 0.0086 0.0118 0.0166 0.0245 0.0379 0.0621 0.1077 0.1923 0.3000 0.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 插值函数10()L x 和三次样条插值函数10()S x 的图像，并比较它们的区别。

先写拉格朗日插值的程序，如下：

function y=lagr1(x0,y0,x)n=length(x0);m=length(x);for i=1:m z=x(i);s=0.0;for k=1:n p=1.0;for j=1:n if j~=k p=p*(z-x0(j))/(x0(k)-x0(j));end end s=p*y0(k)+s;end y(i)=s;end end 然后直接在 matlab 的命令窗口输入：

x0=[-5-4-3-2-1 0 1 2 3 4 5 ];y0=1./(x0.^2+1);x=-5:0.01:5;df=diff(y0);conds=1;y=selfspline(x0,y0,df,x,conds);y1=1./(x.^2+1);y2=lagr1(x0,y0,x);plot(x,y1,"b",x,y,"y-",x,y2,"r:")legend(" 原函数"," 三次样条插值","Lagrange 插值")即可得到以下图形：

-5-4-3-2-1 0 1 2 3 4 5-0.500.511.52 原函数三次样条插值Lagrange插值(3)分析结果：

根据图像可知，用拉格朗日插值函数进行插值时，在次数较高时会出现龙格现象，导致其值与原函数相差甚大，在一定范围内插值拟合效果不好，而三次样条插值函数可以很好地接近原函数且光滑性较好，是比较理想的插值函数。

五、教师评语（或成绩）教师签字：

年月日