数值分析实验二matlab插值法.docx
《数值分析实验二matlab插值法.docx》由会员分享，可在线阅读，更多相关《数值分析实验二matlab插值法.docx（29页珍藏版）》请在冰点文库上搜索。
数值分析实验二matlab插值法
实验二插值法
实验2.1（多项式插值的振荡现象）
实验要求1：
程序：
M文件：
lagrange.m
test1p1.m
Lagrange函数：
%lagrangeinsert
functiony=lagrange（x0,y0,x）
n=length（x0）;m=length（x）;
fori=1:
m
z=x（i）;
s=0.0;
fork=1:
n
p=1.0;
forj=1:
n
ifj~=k
p=p*（z-x0（j））/（x0（k）-x0（j））;
end
end
s=p*y0（k）+s;
end
y（i）=s;
end
主函数：
x1=[-1:
1:
1];y1=1./（1+25*x1.^2）;
x2=[-1:
0.5:
1];y2=1./（1+25*x2.^2）;
x3=[-1:
0.4:
1];y3=1./（1+25*x3.^2）;
x4=[-1:
0.25:
1];y4=1./（1+25*x4.^2）;
x5=[-1:
0.2:
1];y5=1./（1+25*x5.^2）;
x6=[-1:
0.1:
1];y6=1./（1+25*x6.^2）;
x0=[-1:
0.01:
1];
subplot（4,2,1）
y7=1./（1+25*x0.^2）;
plot（x0,y7,'-b'）
y0=lagrange（x1,y1,x0）;
y1=1./（1+25*x0.^2）;
subplot（4,2,2）
plot（x0,y0,'--b'）
y0=lagrange（x2,y2,x0）;
y2=1./（1+25*x0.^2）;
subplot（4,2,3）
plot（x0,y0,'-g'）
y0=lagrange（x3,y3,x0）;
y3=1./（1+25*x0.^2）;
subplot（4,2,4）
plot（x0,y0,'--g'）
y0=lagrange（x4,y4,x0）;
y4=1./（1+25*x0.^2）;
subplot（4,2,5）
plot（x0,y0,'-r'）
y0=lagrange（x5,y5,x0）;
y5=1./（1+25*x0.^2）;
subplot（4,2,6）
plot（x0,y0,'--r'）
y0=lagrange（x6,y6,x0）;
y6=1./（1+25*x0.^2）;
subplot（4,2,7）
plot（x0,y0,'-y'）
数值实验结果及分析：
实验结果：
图一为原函数曲线。
然后图二至图七为插值点分别为2个、4个、5个、8个、10个、20个时的函数图像。
实验分析
在此，我们看得到当采用四个和五个点进行拟合的时候，得到的函数在中间部分拟合得较好。
随着插值点的次数越多，两端出现了明显的震荡。
如图五、图六为插值节点增加到8个、10个的情况，得到的插值函数的图像明显与原函数的相似度很低，并且在两端出现了较大的震荡。
图七增加到20个点时误差震荡极大，曲线基本不可取。
从图中可以明显看到，随着插值节点的增多，虽然在中间的拟合函数拟合得比较好，在两端出现了较为明显的震荡，并且随着插值点的逐渐增多震荡越来越明显。
在这里，我们采用一个很简单的函数对其作出了很直接的验证。
实验要求2：
（1）对：
程序：
M文件：
lagrange.m
test1p2a.m
Lagrange函数：
%lagrangeinsert
functiony=lagrange（x0,y0,x）
n=length（x0）;m=length（x）;
fori=1:
m
z=x（i）;
s=0.0;
fork=1:
n
p=1.0;
forj=1:
n
ifj~=k
p=p*（z-x0（j））/（x0（k）-x0（j））;
end
end
s=p*y0（k）+s;
end
y（i）=s;
end
主函数：
x1=[-5:
2.5:
5];y1=x1./（1+x1.^4）;
x2=[-5:
2:
5];y2=x2./（1+x2.^4）;
x3=[-5:
1:
5];y3=x3./（1+x3.^4）;
x4=[-5:
0.5:
5];y4=x4./（1+x4.^4）;
x5=[-5:
0.4:
5];y5=x5./（1+x5.^4）;
x6=[-5:
0.1:
5];y6=x6./（1+x6.^4）;
x0=[-5:
0.01:
5];
subplot（4,2,1）
y7=x0./（1+x0.^4）;
plot（x0,y7,'-b'）
y0=lagrange（x1,y1,x0）;
y1=x0./（1+x0.^4）;
subplot（4,2,2）
plot（x0,y0,'--b'）
y0=lagrange（x2,y2,x0）;
y2=x0./（1+x0.^4）;
subplot（4,2,3）
plot（x0,y0,'-g'）
y0=lagrange（x3,y3,x0）;
y3=x0./（1+x0.^4）;
subplot（4,2,4）
plot（x0,y0,'--g'）
y0=lagrange（x4,y4,x0）;
y4=x0./（1+x0.^4）;
subplot（4,2,5）
plot（x0,y0,'-r'）
y0=lagrange（x5,y5,x0）;
y5=x0./（1+x0.^4）;
subplot（4,2,6）
plot（x0,y0,'--r'）
y0=lagrange（x6,y6,x0）;
y6=x0./（1+x0.^4）;
subplot（4,2,7）
plot（x0,y0,'-y'）
数值实验结果及分析：
实验结果：
图一为原函数曲线。
然后图二至图七为插值点分别为4个、5个、10个、20个、25个、100个时的函数图像。
实验分析
与第一问得到的拟合函数的结果相类似，随着插值节点愈多，在中间阶段得到的函数图像叫原函数比较相似，但是在两端出现了明显的震荡，并且随着插值点的增多，震荡越来越大且不可控制。
（2）对g（x）=arctanx
程序：
M文件名：
lagrange.m
test1p2b.m
Lagrange函数：
%lagrangeinsert
functiony=lagrange（x0,y0,x）
n=length（x0）;m=length（x）;
fori=1:
m
z=x（i）;
s=0.0;
fork=1:
n
p=1.0;
forj=1:
n
ifj~=k
p=p*（z-x0（j））/（x0（k）-x0（j））;
end
end
s=p*y0（k）+s;
end
y（i）=s;
end
主函数：
x1=[-5:
2.5:
5];y1=atan（x1）;
x2=[-5:
2:
5];y2=atan（x2）;
x3=[-5:
1:
5];y3=atan（x3）;
x4=[-5:
0.5:
5];y4=atan（x4）;
x5=[-5:
0.4:
5];y5=atan（x5）;
x6=[-5:
0.1:
5];y6=atan（x6）;
x0=[-5:
0.001:
5];
subplot（4,2,1）
y7=atan（x0）;
plot（x0,y7,'-b'）
y0=lagrange（x1,y1,x0）;
y1=atan（x1）;
subplot（4,2,2）
plot（x0,y0,'--b'）
y0=lagrange（x2,y2,x0）;
y2=atan（x2）;
subplot（4,2,3）
plot（x0,y0,'-g'）
y0=lagrange（x3,y3,x0）;
y3=atan（x3）;
subplot（4,2,4）
plot（x0,y0,'--'）
y0=lagrange（x4,y4,x0）;
y4=atan（x4）;
subplot（4,2,5）
plot（x0,y0,'-r'）
y0=lagrange（x5,y5,x0）;
y5=atan（x5）;
subplot（4,2,6）
plot（x0,y0,'--r'）
y0=lagrange（x6,y6,x0）;
y6=atan（x6）;
subplot（4,2,7）
plot（x0,y0,'-y'）
数值实验结果及分析：
实验结果：
图一为原函数曲线。
然后图二至图七为插值点分别为2个、4个、5个、8个、10个、20个时的函数图像。
实验分析
当插值节点为五个时，拟合程度不是很好，但是在两端出现的波动较小。
当插值节点为10个时，可以看出，插值点确定出来的函数较原函数比较相识，中间阶段拟合程度比较高，在两端的波动也比较小。
但是当节点增加到20、25、100个时，在插值区间两端出现了较大的波动，也就是龙格现象。
可以看出，等距插值中，增加插值点的个数，并不能使得到的多项式函数与原函数的拟合程度更高。
实验要求3：
程序：
M文件名：
lagrange.m
test1p3.m
Lagrange函数：
%lagrangeinsert
functiony=lagrange（x0,y0,x）
n=length（x0）;m=length（x）;
fori=1:
m
z=x（i）;
s=0.0;
fork=1:
n
p=1.0;
forj=1:
n
ifj~=k
p=p*（z-x0（j））/（x0（k）-x0（j））;
end
end
s=p*y0（k）+s;
end
y（i）=s;
end
主函数：
k=[1:
5];
x1=cos（（2*k-1）*3.1415/12）;
y1=1./（1+25*x1.^2）;
a=[1:
15];
x2=cos（（2*a-1）*3.1415/32）;
y2=1./（1+25*x2.^2）;
b=[1:
50];
x3=cos（（2*b-1）*3.1415/102）;
y3=1./（1+25*x3.^2）;
c=[1:
100];
x4=cos（（2*b-1）*3.1415/102）;
y4=1./（1+25*x4.^2）;
x0=[-1:
0.01:
1];
subplot（3,2,1）
y5=1./（1+25*x0.^2）;
plot（x0,y5,'-b'）
y0=lagrange（x1,y1,x0）;
y1=1./（1+25*x0.^2）;
subplot（3,2,2）
plot（x0,y0,'--b'）
y0=lagrange（x2,y2,x0）;
y2=1./（1+25*x0.^2）;
subplot（3,2,3）
plot（x0,y0,'-g'）
y0=lagrange（x3,y3,x0）;
y3=1./（1+25*x0.^2）;
subplot（3,2,4）
plot（x0,y0,'--g'）;
y0=lagrange（x4,y4,x0）;
y4=1./（1+25*x0.^2）;
subplot（3,2,5）
plot（x0,y0,'--r'）
数值实验结果及分析：
实验结果：
应用上面的lagrange函数以及主函数，分别就
在
区间上做出了图像，如下图的图一。
然后分别采用切比雪夫点进行插值，然后图二至图七为插值点分别为5个、15个、50个、100个时的函数图像。
：
实验分析：
从上面的图中，我们可以看到，采用切比雪夫插值节点来采点，随着采点数目的增多，函数与原函数的拟合程度更相似。
消除了如第一问和第二问随着插值节点的增多，两端出现越来越大且不可控制的波动情况。
故可以认为，采用切比雪夫插值节点来拟合形如
的函数式更为恰当。
实验2.2（样条插值的收敛性）
实验要求
（一）
程序：
M文件名：
test2p1.m
x1=[-1:
1];y1=1./（1+25*x1.^2）;
x2=[-1:
0.5:
1];y2=1./（1+25*x2.^2）;
x3=[-1:
0.4:
1];y3=1./（1+25*x3.^2）;
x4=[-1:
0.25:
1];y4=1./（1+25*x4.^2）;
x5=[-1:
0.2:
1];y5=1./（1+25*x5.^2）;
x6=[-1:
0.1:
1];y6=1./（1+25*x6.^2）;
x0=[-1:
0.01:
1];
subplot（4,2,1）
y7=1./（1+25*x0.^2）;
plot（x0,y7,'-b'）
y0=spline（x1,y1,x0）;
y1=1./（1+25*x0.^2）;
subplot（4,2,2）
plot（x0,y0,'--b'）
y0=spline（x2,y2,x0）;
y2=1./（1+25*x0.^2）;
subplot（4,2,3）
plot（x0,y0,'-g'）
y0=spline（x3,y3,x0）;
y3=1./（1+25*x0.^2）;
subplot（4,2,4）
plot（x0,y0,'--g'）
y0=spline（x4,y4,x0）;
y4=1./（1+25*x0.^2）;
subplot（4,2,5）
plot（x0,y0,'-r'）
y0=spline（x5,y5,x0）;
y5=1./（1+25*x0.^2）;
subplot（4,2,6）
plot（x0,y0,'--r'）
y0=spline（x6,y6,x0）;
y6=1./（1+25*x0.^2）;
subplot（4,2,7）
plot（x0,y0,'-y'）
数值实验结果
应用matlab函数库中提供的spline函数以及原函数表达式函数，分别就
在
区间上做出了原函数图像，如下图的图一。
然后图二至图七为插值点分别为2个、4个、5个、8个、10个、20个时的函数图像。
实验分析
在此看到，当采用五个插值函数的节点时和采用更多节点的三次样条函数对原函数进行插值拟合的时候，得到的函数都是光滑的。
并且随着插值节点的增多，也没有出现如实验2·1所示的随着点的增多而出现龙格现象。
并且可以看出，随着三次样条函数的分段的增多，插值函数的图像与原函数的图像越来越相似。
实验要求
（二）：
程序：
test2p2.m
clear
clc
X=[0,1,2,3,4,5,6,7,8,9,10];
Y=[0.0,0.79,1.53,2.19,2.71,3.03,3.27,2.89,3.06,3.19,3.29];
dy0=0.8;
dyn=0.2;
n=length（X）-1;
d=zeros（n+1,1）;
h=zeros（1,n-1）;
f1=zeros（1,n-1）;
f2=zeros（1,n-2）;
fori=1:
n;h（i）=X（i+1）-X（i）;
f1（i）=（Y（i+1）-Y（i））/h（i）;
end
fori=2:
n;f2（i）=（f1（i）-f1（i-1））/（X（i+1）-X（i-1））;
d（i）=6*f2（i）;
end
d
（1）=6*（f1
（1）-dy0）/h
（1）;
d（n+1）=6*（dyn-f1（n-1））/h（n-1）;A=zeros（n+1,n+1）;
B=zeros（1,n-1）;
C=zeros（1,n-1）;
fori=1:
n-1
B（i）=h（i）/（h（i）+h（i+1））;
C（i）=1-B（i）;
end
A（1,2）=1;
A（n+1,n）=1;
fori=1:
n+1
A（i,i）=2;
end
fori=2:
n
A（i,i-1）=B（i-1）;
A（i,i+1）=C（i-1）;
end
M=A\d;
x=zeros（1,n）;
S=zeros（1,n）;
fori=1:
n
x（i）=X（i）+0.5;
S（i）=Y（i）+（f1（i）-（M（i）/3+M（i+1）/6）*h（i））*（x（i）-X（i））+M（i）/2*（x（i）-X（i））^2+（M（i+1）-M（i））/（6*h（i））*（x（i）-X（i））^3;
end
plot（X,Y,'k'）;holdon;
plot（x,S,'o'）;
title（'三次样条插值效果图'）;
legend（'已知插值节点','三次样条插值'）;
holdoff
返回插值函数：
S（x）=0.8*x-0.001486*x^2-0.008514*x^3（0,1）
S（x）=0.8122*x-0.01365*x^2-0.004458*x^3-0.004055（1,2）
S（x）=0.8218*x-0.01849*x^2-0.003652*x^3-0.01051（2,3）
S（x）=0.317*x^2-0.1847*x-0.04093*x^3+0.9961（3,4）
S（x）=6.934*x-1.463*x^2+0.1074*x^3-8.496（4,5）
S（x）=4.177*x^2-21.26*x-0.2686*x^3+38.5（5,6）
S（x）=53.86*x-8.344*x^2+0.427*x^3-111.7（6,7）
S（x）=6.282*x^2-48.52*x-0.2694*x^3+127.1（7,8）
S（x）=14.88*x-1.643*x^2+0.06076*x^3-41.93（8,9）
S（x）=8.966*x-0.986*x^2+0.03641*x^3-24.18（9,10）
图像：
实验总结：
插值多项式余项公式说明插值节点越多，一般说来误差越小，函数逼近越好，但这也不是绝对的，因为余项的大小既与插值节点的个数有关，也与函数f（x）的高阶导数有关。
换句话说，适当地提高插值多项式的次数，有可能提高计算结果的准确程度，但并非插值多项式的次数越高越好。
当插值节点增多时，不能保证非节点处的插值精度得到改善，有时反而误差更大。
而这个误差就是由于计算机计算里面的舍入误差和截断误差。
所以我们在进行多项式插值的时候，必须考虑这个函数的在插值区间的高阶导数的值是否是存在极大值的。
但是可以考虑的就是，在针对某些函数的时候，可以采用用切比雪夫插值的形式来对原函数进行插值拟合。
这样，随着插值点的增多，插值函数越来越接近于原函数。
在实际的插值过程中，我们不需要一定就用高阶的多项式函数进行插值。
在一些函数中，我们可以采用一次函数分段插值，只要取足够多的点，那么也可以达到很高的拟合程度。
理论上也已经证明，当取点数区域无穷时，得到的函数就等于原函数。
另外我们还可以考虑Hermite插值,三次样条插值等等其他一系列可以考虑的插值函数。
总之，在数值计算过程中，我们必须多方面考虑误差的来源和那些由于计算机内部运算而引起的截断误差和舍入误差，而且在针对不同形式的函数，我们可以采用不同的插值方式来对原函数进行插值。
实验2.3
程序：
M文件名：
test3.m
figure（'position',get（0,'screensize'））
axes（'position',[0011]）
[x,y]=ginput
n=length（x）;
s=（1:
n）';t=（1:
0.05:
n）';
u=interp1（s,x,t,'spline'）;
v=interp1（s,y,t,'spline'）;
clfreset,
plot（x,y,'.',u,v,'-'）,
%xlabel（'X'）,ylabel（'Y'）
运行结果：
x=
0.2017
0.2332
0.2786
0.3122
0.3613
0.3949
0.4059
0.4242
0.4469
0.4740
0.5026
0.5245
0.5553
0.5692
0.5692
0.5553
0.5355
0.5267
0.4857
0.4542
0.4228
0.4067
0.4154
0.4330
0.4513
0.4901
0.5611
0.6182
0.6658
0.6702
0.6636
0.6607
0.6285
0.5875
0.5209
0.4901
0.4579
0.4528
0.4887
0.5135
0.5348
0.5450
0.5714
0.6058
0.6263
0.6431
0.6834
0.7017
0.7053
0.7053
0.6797
0.6321
0.5802
0.5692
0.5253
0.4813
0.4469
0.4616
0.4586
0.4894
0.4967
0.5179
0.5267
0.5362
0.5501
0.5867
0.6021
0.6329
0.6453
0.6651
0.6651
0.6490
0.5684
0.5260
0.4813
0.4425
0.3906
0.3335
0.3730
0.4015
0.4184
0.4323
0.4396
0.4777
0.4872
0.4484
0.4286
0.4037
0.3371
0.3225
0.2903
0.2683
0.2376
0.2171
0.1790
0.1336
0.1600
0.1387
0.1490
0.1468
0.1482
0.1475
0.1490
0.1490
0.1504
y=
0.8374
0.8675
0.8832
0.8947
0.9119
0.9219
0.9219
0.9291
0.9434
0.9606
0.9692
0.9778
0.9778
0.9620
0.9348
0.9176
0.9033
0.8933
0.8746
0.8589
0.8388
0.8288
0.8102
0.8016
0.8087
0.8145
0.8173
0.8489
0.8503
0.8331
0.8116
0.7729
0.7486
0.7457
0.7357
0.7113
0.7085
0.6984
0.6913
0.6884
0.7070
0.7042
0.7056
0.7199
0.7213
0.7314
0.7342
0.7299
0.6970
0.6569
0.6297
0.6211
0.6168
0.6139
0.5924
0.5910
0.5895
0.5752
0.5365
0.5279
0.5251
0.5208
0.5208
0.5193
0.5150
0.4878
0.
升级会员 