MATLAB实现拉格朗日插值.docx
《MATLAB实现拉格朗日插值.docx》由会员分享，可在线阅读，更多相关《MATLAB实现拉格朗日插值.docx（9页珍藏版）》请在冰点文库上搜索。
MATLAB实现拉格朗日插值
数值分析上机报告
题目：
插值法
学号：
201014924
姓名：
靳会有
一、调用MATLAB内带函数插值
1、MATLAB内带插值函数列举如下：
interp1
interpft
interp2
interp3
interpn
spline
meshgrid
ndgrid
griddata
一维数据内插（查表法）
使用FFT方法的一维数据内插
二维数据内插（查表法）
三维数据内插（查表法）
多维数据内插（查表法）
三次样条内插
为三维绘图产生X和Y阵
为多维函数和内插产生阵列
数据网格
2、取其中的一维数据内插函数（interp1）为例，程序如下：
其调用格式为：
yi=interp1（x,y,xi）
yi=interp1（x,y,xi,method）
举例如下：
x=0:
10:
100
y=[40444652657680828892110];
xi=0:
1:
100
yi=interp1（x,y,xi,'spline'）
3、其他内带函数调用格式为：
Interpft函数：
y=interpft（x,n）
y=interpft（x,n,dim）
interp2函数：
ZI=interp2（X,Y,Z,XI,YI），ZI=imerp2（Z,ntimes）
ZI=interp2（Z,XI,YI），ZI=interp2（X,Y,Z,XI,YI,method）interp3函数：
VI=interp3（X,Y,Z,V,XI,YI,ZI）VI=interp3（V,ntimes）
VI=interp3（V,XI,YI,ZI）VI=interp3（…,method）
Interpn函数：
VI=interpn（X1,X2,X3,…,V,Y1,Y2,Y3,…）
VI=interpn（V,ntimes）
VI=interpn（V,Yl,Y2,Y3,…） VI=interpn（…,method）
Spline函数：
yi=spline（x,y,xi）
pp=spline（x,y）
meshgrid函数：
[X,Y]=meshgrid（x,y）
[X,Y]=meshgrid（x）
[X,Y,Z]=meshgrid（x,y,z）
Ndgrid函数：
[X1,X2,X3,…]=ndgrid（x1,x2,x3,…）
[X1,X2,X3,…]=ndgrid（x）
Griddata函数：
ZI=griddata（x,y,z,XI,YI）
[XI,YI,ZI]=griddata（x,y,z,xi,yi）
[…]=griddata（…method）
二、自编函数插值
1、拉格朗日插值法：
建立M文件：
functionf=Language（x,y,x0）
symstl;
if（length（x）==length（y））
n=length（x）;
else
disp（'x和y的维数不相等！
'）;
return;%检错
end
h=sym（0）;
for（i=1:
n）
l=sym（y（i））;
for（j=1:
i-1）
l=l*（t-x（j））/（x（i）-x（j））;
end;
for（j=i+1:
n）
l=l*（t-x（j））/（x（i）-x（j））;
end;
h=h+l;
end
simplify（h）;
if（nargin==3）
f=subs（h,'t',x0）;%计算插值点的函数值
else
f=collect（h）;
f=vpa（f,6）;%将插值多项式的系数化成6位精度的小数
end
在MATLAB中输入：
x=[18316668707270;]
y=[23335251434046];
f=Language（x,y）
plot（x,y）
结果为：
f=Inf+（-t）*Inf-54329.8*t^2+1503.75*t^3-22.2065*t^4+0.16789*t^5-0.000512106*t^6
图形如下：
MATLAB实现拉格朗日插值
建立如下拉格朗日插值函数：
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
画图程序如下：
x=[-5:
1:
5];
y=1./（1+x.^2）;
x0=[-5:
0.001:
5];
y0=lagrange（x,y,x0）;
y1=1./（1+x0.^2）;
plot（x0,y0,'r'）
holdon
plot（x0,y1,'g'）
注：
画出的图形为n=10的图形
得到图形如下：
n=10的图像
牛顿K次插值多项式
一、实验目的：
1、掌握牛顿插值法的基本思路和步骤。
2、培养编程与上机调试能力。
二、牛顿插值法基本思路与计算步骤：
给定插值点序列（
。
构造牛顿插值多项式
。
输入要计算的函数点
并计算
的值，利用牛顿插值公式，当增加一个节点时，只需在后面多计算一项，而前面的计算仍有用；另一方面
的各项系数恰好又是各阶均差，而各阶均差可用均差公式来计算。
为
的一阶均差。
为
的k阶均差。
均差表：
零阶均差
一阶均差
二阶均差
三阶均差
X0
f（X0）
X1
f（X1）
f[X0,X1]
X2
f（X2）
f[X1,X2]
f[X0,X1,X2]
X3
f（X3）
f[X2,X3]
f[X1,X2,X3]
f[X0,X1,X2X3]
M
M
M
M
M
牛顿插值法计算步骤：
1．输入
值及（
；要计算的函数点
。
2．对给定的
由
计算
的值。
3．输出
。
程序清单：
function[c,d]=newpoly（x,y）
%牛顿插值的MATLAB实现
%这里x为n个节点的横坐标所组成的向量，y为纵坐标所组成的向量。
%c为所求的牛顿插值多项式的系数构成的向量。
n=length（x）;%取x的个数。
d=zeros（n,n）;%构造nXn的空数组。
d（:
1）=y';
forj=2:
n
fork=j:
n
d（k,j）=（d（k,j-1）-d（k-1,j-1））/（x（k）-x（k-j+1））;
end
end
c=d（n,n）;
fork=（n-1）:
-1:
1
c=conv（c,poly（x（k）））;%conv求积，poly（x）将该多项式的系数赋给向量。
m=length（c）;
c（m）=c（m）+d（k,k）;
end
五、测试数据与结果：
测试数据：
（第三章习题第三题第2题）
f（x）=lnx的数值如表所示，构造牛顿插值多项式并求ln0.53的值。
X
0.4
0.5
0.6
0.7
0.8
lnx
-0.916291
-0.693147
-0.510826
-0.357765
-0.223144
解：
由表可知x0=0.4,x1=0.5,x2=0.6,x3=0.7,x4=0.7，函数值：
Y0=-0.916291,y1=-0.693147,y2=-0.510826,y3=-0.357765,y4=-0.223144
建立一个主程序np.m
clc
clear
newpoly（[0.4,0.5,0.6,0.7,0.8],[-0.916291,-0.693147,-0.510826,-0.357765,-0.223144]）
计算结果如下：
ans=
-0.30962.6083-5.48615.6921-2.4744
由此看出所求的牛顿多项式为：
P（x）=-0.3096x4+2.6083x3-5.4861x2+5.6921x-2.4744
P（0.53）=-0.6347。
升级会员 