Langrage和Newton插值法的matlab实现
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
仅供参考
1．已知数据如下：
（1）用MATLAB语言编写按Langrage插值法和Newton插值法计算插值的程序，对以上数据进行插值；（2）利用MATLAB在第一个图中画出离散数据及插值函数曲线。
(1.1)langrage插值法编程实现
syms x
x0=[0.2,0.4,0.6,0.8,1.0];
y0=[0.98,0.92,0.81,0.64,0.38];
for i=1:5
a=1;
for j=1:5
if j~=i
a=expand(a*(x-x0(j)));
end
end
b=1;
for k=1:5
if k~=i
b=b*(x0(i)-x0(k));
end
end
A(i)=expand(a/b);
end
L=0;
for p=1:5
L=L+y0(p)*A(p);
end
L
L =
-25/48*x^4+5/6*x^3-53/48*x^2+23/120*x+49/50
（1.2）Newton插值程序实现
clear all
clc
syms x
x0=[0.2,0.4,0.6,0.8,1.0];
y0=[0.98,0.92,0.81,0.64,0.38];
for k=1:5
for i=1:k
a=1;
b=0;
for j=1:k
if j~=i
a=a*(x0(i)-x0(j));
end
end
b=b+y0(i)/a;
end
A(k)=b;
end
B=[1,(x-x0(1)),(x-x0(1))*(x-x0(2)),(x-x0(1))*(x-x0(2))*(x-x0(3)),(x-x 0(1))*(x-x0(2))*(x-x0(3))*(x-x0(4))];
L1=A.*B;
l=0;
for m=1:5
l=l+L1(m);
end
L=expand(l)
L =
61/100+13/30*x+383/48*x^2-155/24*x^3+475/48*x^4
（2）画图
x0=[0.2,0.4,0.6,0.8,1.0];
y0=[0.98,0.92,0.81,0.64,0.38];
subplot(1,2,1);
plot(x0(1),y0(1),'+r',x0(2),y0(2),'+r',x0(3),y0(3),'+r',x0(4),y0(4),' +r',x0(5),y0(5),'+r')
x=0:0.05:1;
y=-25/48.*x.^4+5/6.*x.^3-53/48.*x.^2+23/120.*x+49/50;
subplot(1,2,2);
plot(x,y)
2．给定函数21(),[1,1]125f x x x ，利用上题编好的Langrage 插值程序（或Newton 插值程序），分别取3个，5个、9个、11个等距节点作多项式插值，分别画出插值函数及原函数()f x 的图形，以验证Runge 现象、分析插值多项式的收敛性。
取三个节点如下
clear
clc
x0=0:0.5:1;
y0=1./(1+25.*x0.^2);
syms x
for i=1:3
a=1;
for j=1:3
if j~=i
a=expand(a*(x-x0(j)));
end
end
b=1;
for k=1:3
if k~=i
b=b*(x0(i)-x0(k));
end
end
A(i)=expand(a/b);
end
L=0;
for p=1:3
L=L+y0(p)*A(p);
end
L
L =
575/377*x^2-1875/754*x+1
x1=0:0.0001:1;
y1=1./(1+25.*x1.^2);
y2=575/377.*x1.^2-1875/754.*x1+1; plot(x1,y1,'+r')
hold on
plot(x1,y2,'*k')
取五个节点如下
clear
clc
x0=0:0.25:1;
y0=1./(1+25.*x0.^2);
syms x
for i=1:5
a=1;
for j=1:5
if j~=i
a=expand(a*(x-x0(j)));
end
end
b=1;
for k=1:5
if k~=i
b=b*(x0(i)-x0(k));
end
end
A(i)=expand(a/b);
end
L=0;
for p=1:5
L=L+y0(p)*A(p);
end
L
L =
1570000/3725137*x^4-9375000/3725137*x^3+16996575/3725137*x^2-25546875 /7450274*x+1
x1=0:0.0001:1;
y1=1./(1+25.*x1.^2);
y2=1570000/3725137.*x1.^4-9375000/3725137.*x1.^3+16996575/3725137.*x1 .^2-25546875/7450274.*x1+1 ;
plot(x1,y1,'+r')
hold on
plot(x1,y2,'*k')