实验一函数插值方法报告

1、实验一 函数插值方法报告一、问题提出 对于给定的一元函数的n+1个节点值。试用Lagrange公式求其插值多项式或分段二次Lagrange插值多项式。数据如下： （1） 0.4 0.55 0.65 0.80 0.95 1.05 0.41075 0.578150.696750.90 1.00 1.25382 求五次Lagrange多项式，和分段三次插值多项式，计算,的值。（提示：结果为, ） （2） 1 2 3 4 5 6 7 0.368 0.135 0.050 0.018 0.007 0.002 0.001 试构造Lagrange多项式，计算的,值。（提示：结果为, ）二、要求 1、 利用La

2、grange插值公式 编写出插值多项式程序； 2、 给出插值多项式或分段三次插值多项式的表达式； 3、 根据节点选取原则，对问题（2）用三点插值或二点插值，其结果如何； 4、 对此插值问题用Newton插值多项式其结果如何。Newton插值多项式如下： 其中： 三、目的和意义 1、 学会常用的插值方法，求函数的近似表达式，以解决其它实际问题； 2、 明确插值多项式和分段插值多项式各自的优缺点； 3、 熟悉插值方法的程序编制； 4、 如果绘出插值函数的曲线，观察其光滑性。 四、实验学时：2学时五、实验步骤： 1进入C或matlab开发环境；2根据实验内容和要求编写程序；3调试程序；4运行程序；5

3、撰写报告,讨论分析实验结果.解：一、编写插值函数结构程序Lagrange插值多项式M文件：lagrange1.mfunction A1,LN,L1,B1=lagrange1(X,Y)m=length(X); LN=ones(m,m);for k=1: m x1=1; for i=1:m if k=i x1=conv(x1,poly(X(i)/(X(k)-X(i); endendL1(k,:)=x1; B1(k,:)=poly2sym (x1)endA1=Y*L1;LN=Y*B1分段三次艾尔米特插值多项式的M文件：Hermite3.mfunction f,ff = Hermite3(x,y,y1

4、)syms t;f = 0.0;if(length(x) = length(y) if(length(y) = length(y1) n = length(x); else disp('y和y的导数的维数不相等'); return; endelse disp('x和y的维数不相等！ '); return;end for i=1:n h = 1.0; a = 0.0; for j=1:n if( j = i) h = h*(t-x(j)2/(x(i)-x(j)2); a = a + 1/(x(i)-x(j); end end f = f + h*(x(i)-t)*

5、(2*a*y(i)-y1(i)+y(i);end ff = subs(f,'t');（1）、求五次Lagrange多项式，和分段三次插值多项式。在主显示区，输入五次Lagrange多项式程序：>> X=0.4 0.55 0.65 0.80 0.95 1.05;>> Y=0.41075 0.57815 0.69675 0.90 1.00 1.25382;>> A1,LN,L1,B1=lagrange1(X,Y)>> plot(X,A1); >> F=poly2sym(A1) 运行后，输出五次Lagrange多项式的结果：

6、A1 = 121.6264 -422.7503 572.5667 -377.2549 121.9718 -15.0845F = (2139673480305281*x5)/17592186044416 - (1859275536318005*x4)/4398046511104 + (9836621836743*x3)- (414796119737013*x2)/1099511627776 + (2145751274873259*x)/17592186044416 - 1061478972867847/70368744177664拉格朗日插值多项式的图如下：在主显示区，

7、输入分段三次艾尔米特插值多项式的程序：>> x=0.4 0.55 0.65 0.80 0.95 1.05;>> y=0.41075 0.57815 0.69675 0.90 1.00 1.25382;>> y1=2.3440 0.9032 1.4329 0.9903 0.9170 5.1439;>> f,ff = Hermite3(x,y,y1);>> ff运行后，分段三次艾尔米特插值多项式的输出结果：ff = (6400000000*(t - 4/5)2*(t - 11/20)2*(t - 13/20)2*(t - 19/20)2*

8、(t - 21/20)2*(2240245151070481*t)/140737488355328 - 52393133567890089/8796093022208000)/184041 - (16000000*(6348013345609171*t)/140737488355328 - 85523418631741336287/1759218604441600000)*(t - 2/5)2*(t - 4/5)2*(t - 11/20)2*(t - 13/20)2*(t - 19/20)2)/169 + (16000000*(4105617466549689*t)/2814749767106

9、56 - 5238387122042657959/703687441776640000)*(t - 2/5)2*(t - 4/5)2*(t - 13/20)2*(t - 19/20)2*(t - 21/20)2)/9 - (256000000*(35097*t)/10000 - 46347/12500)*(t - 2/5)2*(t - 11/20)2*(t - 13/20)2*(t - 19/20)2*(t - 21/20)2)/81 - (400000000*(13147*t)/20000 - 449611/400000)*(t - 2/5)2*(t - 4/5)2*(t - 11/20)2

10、*(t - 19/20)2*(t - 21/20)2)/81 - (10000000000*(84913*t)/11000 - 1833347/220000)*(t - 2/5)2*(t - 4/5)2*(t - 11/20)2*(t - 13/20)2*(t - 21/20)2)/9801分段三次艾尔米特插值多项式的图如下： (2)、试构造Lagrange多项式。结果为, 在主显示区，输入程序：>> X=1 2 3 4 5 6 7;>> Y=0.368 0.135 0.050 0.018 0.007 0.002 0.001;>> A1,LN,L1,B1=l

11、agrange1(X,Y)>> plot(X,A1);>> F=poly2sym(A1)运行后，输出结果的Lagrange多项式的结果：A1 =0.0001 -0.0016 0.0186 -0.1175 0.4419 -0.9683 0.9950F = (4304240283865561*x6)/73786976294838206464 - (7417128346304051*x5)/4611686018427387904 + (223*x4)/12000 - (2821*x3)/24000 + (994976512675275*x2)/2251799813685248

12、 - (19367*x)/20000 + 199/200Lagrange多项式的图如下： 二、计算函数值计算函数值的主程序：lagrangezhi.mfunction y,R=lagrangezhi(X,Y,x,M)n=length(X); m=length(x);for i=1:m z=x(i);s=0.0; for k=1:n p=1.0; q1=1.0; c1=1.0;for j=1:n if j=kp=p*(z-X(j)/(X(k)-X(j); end q1=abs(q1*(z-X(j);c1=c1*j; end s=p*Y(k)+s; end y(i)=s;endR=M*q1/c1;

13、（1）、计算、 的值。在主显示区，输入程序：>> x=0.596; M=1; X=0.4,0.55,0.65,0.80,0.95,1.05;>> Y=0.41075,0.57815,0.69675,0.90,1.00,1.25382; >> y,R=lagrangezhi(X,Y,x,M)运行结果：y = 0.6257R = 2.2170e-008在主显示区，输入程序：>> x=0.99; M=1; X=0.4,0.55,0.65,0.80,0.95,1.05;>> Y=0.41075,0.57815,0.69675,0.90,1.0

14、0,1.25382; >> y,R=lagrangezhi(X,Y,x,M)运行结果：y = 1.0542R = 5.5901e-008（2）、计算、的值在主显示区，输入程序：>> x=1.8; M=1; X=1,2,3,4,5,6,7;>> Y=0.368,0.135,0.050,0.018,0.007,0.002,0.001;>> y,R=lagrangezhi(X,Y,x,M)运行结果：y = 0.1648R =0.0059在主显示区，输入程序：>> x=6.15; M=1; X=1,2,3,4,5,6,7;>> Y

15、=0.368,0.135,0.050,0.018,0.007,0.002,0.001;>> y,R=lagrangezhi(X,Y,x,M)运行结果：y = 0.0013R =0.0042三、Newton插值多项式Newton插值多项式主程序M文件：Newton.mfunction A,C,L,wcgs,Cw= Newton(X,Y)n=length(X); A=zeros(n,n); A(:,1)=Y's=0.0; p=1.0; q=1.0; c1=1.0;for j=2:nfor i=j:nA(i,j)=(A(i,j-1)- A(i-1,j-1)/(X(i)-X(i-j

16、+1);endb=poly(X(j-1);q1=conv(q,b); c1=c1*j; q=q1;endC=A(n,n); b=poly(X(n); q1=conv(q1,b);for k=(n-1):-1:1C=conv(C,poly(X(k); d=length(C); C(d)=C(d)+A(k,k);endL(k,:)=poly2sym(C); Q=poly2sym(q1);syms Mwcgs=M*Q/c1; Cw=q1/c1;在主显示区，输入的程序:>> x=0.4 0.55 0.65 0.80 0.95 1.05;>> y=0.41075 0.57815

17、0.69675 0.90 1.00 1.25382;>> A,C,L,wcgs,Cw= Newton(x,y)>> syms x;>> ezplot(L,0 1.1);运行结果如下，得到A = 0.4108 0 0 0 0 0 0.5782 1.1160 0 0 0 0 0.6967 1.1860 0.2800 0 0 0 0.9000 1.3550 0.6760 0.9900 0 0 1.0000 0.6667 -2.2944 -7.4261 -15.3020 0 1.2538 2.5382 7.4861 24.4514 63.7551 121.6264C

18、 = 121.6264 -422.7503 572.5667 -377.2549 121.9718 -15.0845 L = (8558693921221117*x5)/70368744177664 - (3718551072636019*x4)/8796093022208 + (5036350380412441*x3)/8796093022208 - (3318368957896111*x2)/8796093022208 + (536437818718315*x)/4398046511104 - 8491831782942691/562949953421312 wcgs = (M*(x6 -

19、 (22*x5)/5 + (1583*x4)/200 - (3721*x3)/500 + (542206127247039*x2)/140737488355328 - (4682696525551953*x)/4503599627370496 + 4111390143022055/36028797018963968)/720 Cw =0.0014 -0.0061 0.0110 -0.0103 0.0054 -0.0014 0.0002牛顿插值多项式的图如下：在主显示区，输入的程序:>> x=1 2 3 4 5 6 7;>> y=0.368 0.135 0.050 0.0

20、18 0.007 0.002 0.001;>> A,C,L,wcgs,Cw= Newton(x,y)>> syms x;>> ezplot(L,0 8);运行结果如下，得到：A = 0.3680 0 0 0 0 0 0 0.1350 -0.2330 0 0 0 0 0 0.0500 -0.0850 0.0740 0 0 0 0 0.0180 -0.0320 0.0265 -0.0158 0 0 0 0.0070 -0.0110 0.0105 -0.0053 0.0026 0 0 0.0020 -0.0050 0.0030 -0.0025 0.0007 -0.0004 0 0.0010 -0.0010 0.0020 -0.0003 0.0005 -0.0000