试验二插值法及数据拟合

1、试验二 插值法一、 实验目的（1） 学会Lagrange 插值和牛顿插值等基本插值方法；（2） 讨论插值的Runge现象，掌握分段线性插值方法（3） 学会Matlab提供的插值函数的使用方法，会用这些函数解决实际问题。二、 实验要求（1） 按照题目要求完成实验内容；（2） 写出相应的Matlab 程序；（3） 给出实验结果(可以用表格展示实验结果)；（4） 分析和讨论实验结果并提出可能的优化实验。（5） 写出实验报告。三、 实验步骤1、用编好的Lagrange 插值法程序计算书本P66 的例1、用牛顿插值法计算P77的例1。2、已知函数在下列各点的值为：0.20.40.60.81.00.980

2、.920.810.640.38试用4次牛顿插值多项式对数据进行插值，根据，画出图形。3、在区间-1,1上分别取用两组等距节点对龙格函数作多项式插值，对不同值，分别画出插值函数及的图形。4、下列数据点的插值01491625364964012345678可以得到平方根函数的近似，在区间0,64上作图。（1） 用这9个点作8次多项式插值。附：试验报告格式样本（正式报告这行可删除）佛山科学技术学院实 验 报 告课程名称 数值分析 实验项目 插值法 专业班级 姓名 学号 指导教师 成 绩 日 期 月 日 一、实验目的1、学会Lagrange 插值、牛顿插值和 分段线性插值等基本插值方法；2、讨论插值的R

3、unge现象，掌握分段线性插值方法3、学会Matlab提供的插值函数的使用方法，会用这些函数解决实际问题。二、实验原理1、拉格朗日插值多项式2、牛顿插值多项式3、分段线性插值三、实验步骤1、用MATLAB编写独立的拉格朗日插值多项式函数2、用MATLAB编写独立的牛顿插值多项式函数3、利用编写好的函数计算本章P66例1、P77例1的结果并比较。4、已知函数在下列各点的值为：0.20.40.60.81.00.980.920.810.640.38试用4次牛顿插值多项式对数据进行插值，根据，画出图形。5、在区间-1,1上分别取用两组等距节点对龙格函数作多项式插值，对不同值，分别画出插值函数及的图形。

4、6、下列数据点的插值01491625364964012345678(1)可以得到平方根函数的近似，在区间0,64上作图。（2） 用这9个点作8次多项式插值。源程序与解答：1、Lagrange 插值多项式源代码function yi=Lagrange(x, y, xi)% Lagrange 插值多项式，其中% x - 向量，全部的插值节点% y - 向量，插值节点处的函数值% xi - 标量，自变量x% yi - xi 处的函数估计值n=length(x); m=length(y);if n=m error('The lengths of X and Y must be equal

5、9;); return;endp=zeros(1,n);% 对向量p赋初值0for k=1:n t=ones(1,n); for j=1:n if j=k if abs(x(k)-x(j)<eps error('the DATA is error!'); return; end t(j)=(xi-x(j)/(x(k)-x(j); end end p(k)=prod(t);endyi=sum(y.*p);（2）function yi=New_Int(x, y, xi)% Newton 基本插值公式，其中% x - 向量，全部的插值节点，按行输入% y - 向量，插值节点处的

6、函数值，按行输入% xi - 标量，自变量x% yi - xi 处的函数估计值n=length(x); m=length(y);if n=m error('The lengths of X and Y must be equal'); return;end% 计算均差表YY=zeros(n); Y(:,1)=y' % Y(:,1)表示矩阵中第一列的元素for k=1:n-1 for i=1:n-k if abs(x(i+k)-x(i)<eps error('the DATA is error!'); return; end Y(i,k+1)=(Y(

7、i+1,k)-Y(i,k)/(x(i+k)-x(i); endend% 计算Newton插值公式 N(xi)yi=0;for i=1:n z=1; for k=1:i-1 z=z*(xi-x(k); end yi=yi+Y(1,i)*z;end3、（1）求书本P66例1的解答x=144,169,225;>> y=12,13,15;>> yi=Lagrange(x,y,175)yi =13.230158730158733（2）书本P77例1的解答x=0.40,0.55,0.65,0.80;>> y=0.41075,0.57815,0.69675,0.88811

8、;>> yi=New_int(x,y,0.596)yi =0.6319144055040004、x=0.2:0.2:1.0; y=0.98,0.92,0.81,0.64,0.38;xx=0.2:0.08:1.0; m=length(xx); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,'o');plot(xx,z,'r*');hold off图一 练习4的图形5、（1）画龙格函数的图形的matlab代码 a=-1;b=1;n=100;h=(b-a)/n;>

9、;> x=a:h:b;y=1./(1+25.*x.2);>> plot(x,y,'k')图二 龙格函数的图形（2）龙格函数的Lagrange（）插值函数画图源程序function Runge(n)% Runge现象% n - 等距离节点a=-1; b=1; h=(b-a)/n;x=a:h:b; y=1./(1+25.*x.2);xx=a:0.01:b; yy=1./(1+25.*xx.2); m=length(xx); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,

10、9;o');plot(xx,z,'r-');hold off图三：Runge(10)的图形图四：Runge(20)的图形6、function L8%平方根函数的8次多项式插值x=0 1 4 9 16 25 36 49 64; y=sqrt(x);m=length(x); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, x(i);endhold onplot(x,y,'o'); plot(x,y,'o');xlabel('x');ylabel('y');plot(x,z,

11、'r-');hold off图五 的图形x=0 1 4 9 16 25 36 49 64;y=0 1 2 3 4 5 6 7 8;xx=0:0.5:64; yy=sqrt(xx);m=length(xx);z=zeros(1,m);for i=1:mz(i)=Lagrange(x,y,xx(i);endplot(x,y,'o',xx,yy,'.',xx,z,'r-'); 图六 的8次Lagrange插值函数图形四、实验结果 MATLAB desktop keyboard shortcuts, such as Ctrl+S, are

12、 now customizable. In addition, many keyboard shortcuts have changed for improved consistency across the desktop. To customize keyboard shortcuts, use Preferences. From there, you can also restore previous default settings by selecting "R2009a Windows Default Set" from the active settings

13、drop-down list. For more information, see Help. Click here if you do not want to see this message again. >> x=144,169,225;y=12,13,15;yi=Lagrange(x,y,175)yi = 13.230158730158733>> x=0.40,0.55,0.65,0.80;y=0.41075,0.57815,0.69675,0.88811;yi=New_int(x,y,0.596)Warning: Could not find an exact

14、 (case-sensitive) match for 'New_int'.C:UsersAdministratorDesktop实验二New_Int.m is a case-insensitive match and will be used instead.You can improve the performance of your code by using exactname matches and we therefore recommend that you update yourusage accordingly. Alternatively, you can

15、disable this warning usingwarning('off','MATLAB:dispatcher:InexactCaseMatch').This warning will become an error in future releases. Y = 0.410750000000000 1.116000000000000 0 0 0.578150000000000 0 0 0 0.696750000000000 0 0 0 0.888110000000000 0 0 0Y = 0.410750000000000 1.1160000000000

16、00 0 0 0.578150000000000 1.186000000000000 0 0 0.696750000000000 0 0 0 0.888110000000000 0 0 0Y = 0.410750000000000 1.116000000000000 0 0 0.578150000000000 1.186000000000000 0 0 0.696750000000000 1.275733333333333 0 0 0.888110000000000 0 0 0Y = 0.410750000000000 1.116000000000000 0.279999999999998 0

17、 0.578150000000000 1.186000000000000 0 0 0.696750000000000 1.275733333333333 0 0 0.888110000000000 0 0 0Y = 0.410750000000000 1.116000000000000 0.279999999999998 0 0.578150000000000 1.186000000000000 0.358933333333334 0 0.696750000000000 1.275733333333333 0 0 0.888110000000000 0 0 0Y = 0.41075000000

18、0000 1.116000000000000 0.279999999999998 0.197333333333340 0.578150000000000 1.186000000000000 0.358933333333334 0 0.696750000000000 1.275733333333333 0 0 0.888110000000000 0 0 0yi = 0.631914405504000>> x=0.2:0.2:1.0; y=0.98,0.92,0.81,0.64,0.38;xx=0.2:0.08:1.0; m=length(xx); z=zeros(1,m);for i

19、=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,'o');plot(xx,z,'r*'); hold off>> a=-1;b=1;n=100;h=(b-a)/n;x=a:h:b;y=1./(1+25.*x.2); plot(x,y,'k')>> a=-1; b=1; h=(b-a)/n;n=10;x=a:h:b; y=1./(1+25.*x.2);xx=a:0.01:b; yy=1./(1+25.*xx.2); m=length(xx); z=zeros(1,m);fo

20、r i=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,'o');plot(xx,z,'r-');hold off>> a=-1; b=1; h=(b-a)/n;n=20;x=a:h:b; y=1./(1+25.*x.2);xx=a:0.01:b; yy=1./(1+25.*xx.2); m=length(xx); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, xx(i);endhold onplot(x,y,'o');plot(xx,z,'r-');hold off>> x=0 1 4 9 16 25 36 49 64; y=sqrt(x);m=length(x); z=zeros(1,m);for i=1:m z(i)=Lagrange(x, y, x(i);endhold onplot(x,y,'o'); plot(x,y,'o');xlabel('x');ylabel('y');plot(x