Matlab插值
实验目的
（1）熟悉拉格朗日插值法、分段线性插值、三次样条插值等多项式的应用，了解各方法的使用范围。
（2）熟练掌握三次样条插值不同条件下的使用和多项式插值的高次插值的使用。
实验要求
实验步骤要有模型建立，模型求解、结果分析。
实验要求
（1）在区间[-1.，1]上分别取n=10，n=20用两组等距节点对龙格函数f(x)=1/(1+25x2)做拉格朗日插值及三次样条插值(第一边界条件，端点的一阶导数值），对每个n值，分别画出插值函数及f(x)的图像。
（2）已知数据
|
x |
0 |
1 |
4 |
9 |
16 |
25 |
36 |
49 |
64 |
|
y |
0 |
1 |
2 |
3 |
4 |
5 |
6 |
7 |
8 |
可以得到平方根函数的近似，在区间[0，64]上作图。
(1)用这9个点作8次拉格朗日多项式插值L8（x）.
(2)用三次样条插值（第一边界条件）程序求S（x）。S'(x0)=f0',S'(xn)=fn',从得到结果看在【0,64】上，哪个插值更精确；在区间【0,1】上，n=8两种插值哪个更精确？
实验步骤
（1）解：本报告编写拉格朗日插值法代码得到拉格朗日每个n值的情形，使用MATLAB自带的interp1()函数实现三样条插值。拉格朗日代码如下，
1 function yi=Lagrange(x,y,xi) 2 m=length(x);n=length(y);p=length(xi); 3 if m~=n 4 error('向量x与y的长度必须一致'); 5 end 6 s=0; 7 for k=1:n 8 t=ones(1,p); 9 for j=1:n 10 if j~=k 11 t=t.*(xi-x(j))./(x(k)-x(j)); 12 end 13 end 14 s=s+t.*y(k); 15 end 16 yi=s; 17 end
n=10时代码如下，（n=20只需改第2行的值）
1 %n=10 2 n=10; 3 x=linspace(-1,1,n); 4 y=1./(1+25*x.^2); 5 x1=linspace(-1,1,100); 6 y1=1./(1+25*x1.^2); 7 %调用插值函数 8 p=Lagrange(x,y,x1);%拉格朗日 9 p1=interp1(x,y,x1,'spline');%三样条插值 10 figure,plot(x1,y1,'m','LineWidth',2); 11 hold on 12 plot(x1,p,'g','LineWidth',2); 13 hold on 14 plot(x1,p1,'b','LineWidth',2); 15 legend('龙格函数','拉格朗日插值函数','三样条插值函数'); 16 %grid on; 17 title('n=10的插值函数及原函数图像'); 18 xlabel('x轴'); 19 ylabel('y轴'); 20
运行程序，结果如下图
当n=10时，可见三样条插值较拉格朗日插值效果更好。
当n=20时，运行程序，（为了区分图像，本报告对龙格函数线条做了修改）
由上图可见，对龙格函数插值使用三样条插值较拉格朗日插值效果更好。
（2）解：本报告编写相应的拉格朗日代码和三样条插值代码，并画出[0,64]和[0,1]的这三者的函数图像及误差图。拉格朗日代码
1 function [yy,p]=lagrange1(x1,y1,xx) 2 %本程序为Lagrange1插值，其中x1,y1 3 %为插值节点和节点上的函数值，输出为插值点xx的函数值， 4 %xx可以是向量。 5 n=length(x1); 6 m=length(y1); 7 if m~=n 8 error('输入有误！！'); 9 end 10 if nargin<3 11 xx=linspace(min(x1),max(x1),n*10); 12 syms x 13 for i=1:n 14 t=x1; 15 t(i)=[]; 16 L(i)=prod((x-t)./(x1(i)-t));% L向量用来存放插值基函数 17 end 18 u=sum(L.*y1); 19 p=simplify(u); % p是简化后的Lagrange插值函数（字符串） 20 yy=subs(p,x,xx); %p是以x为自变量的函数，并求xx处的函数值 21 end
求三样条插值S(x)代码
1 %三样条拟合函数 2 %输入：x1,y1插值点 3 %输出：方程与函数值 4 function [yy,s3]=sanci(x1,y1) 5 syms x; 6 n=length(x1); 7 m=length(y1); 8 if n~=m 9 error('输入参数有误'); 10 end 11 p=polyfit(x1,y1,3); 12 s3=p(1)+p(2).*x+p(3).*x.^2+p(4).*x.^3; 13 x2=linspace(min(x1),max(x1),x1(n)); 14 yy=polyval(p,x2); 15 %plot(x2,yy); 16 end
解题1及题2的部分，解题1部分图
L8（x）= -(x*(143*x^7 - 29260*x^6 + 2366546*x^5 - 97191380*x^4 + 2171047879*x^3 -26340674360*x^2 + 166253376432*x - 577880352000))/435891456000。
解题2部分截图
S（x）= (3500*x^3)/6927 + (2762705724454173*x^2)/9007199254740992 - (3397610023959411*x)/576460752303423488 + 6795786867330573/147573952589676412928
作[0,64]的函数图像，代码
1 %9个点 2 x=[0 1 4 9 16 25 36 49 64]; 3 y=[0 1 2 3 4 5 6 7 8]; 4 x1=linspace(0,64,64); 5 %调用插值函数 6 [p,S]=lagrange1(x,y,x1);%拉格朗日 7 [p1,s3]=sanci(x,y);%三样条插值 8 figure,plot(x1,sqrt(x1),'rs','LineWidth',2); 9 hold on 10 plot(x1,p,'g','LineWidth',2); 11 hold on 12 plot(x1,p1,'b','LineWidth',2); 13 %设置图例，命名并设置放在左上角位置 14 legend('平方根函数','拉格朗日插值函数','三样条插值函数','Location','northwest'); 15 grid on; 16 title('[0,64]的n=9的插值函数及原函数图像'); 17 xlabel('x轴'); 18 ylabel('y轴'); 19
运行示例：
在[0,64]的插值中，显然三样条插值更好。那么在[0,1]，n=8条件下的插值呢,对上述程序稍作修改，运行结果如下，
显然，在[0,1]，n=8的条件插值中，拉格朗日的效果更好，代码如下
1 %8个点 2 x=linspace(0,1,8); 3 y=sqrt(x); 4 x1=linspace(0,1,20); 5 %调用插值函数 6 [p,S]=lagrange1(x,y,x1);%拉格朗日 7 [p1,s3]=sanci(x,y);%三样条插值 8 figure,plot(x1,sqrt(x1),'rs','LineWidth',2); 9 hold on 10 plot(x1,p,'g','LineWidth',2); 11 hold on 12 plot(x1,p1,'bo','LineWidth',2); 13 %设置图例，命名并设置放在右下角位置 14 legend('平方根函数','拉格朗日插值函数','三样条插值函数','Location','southeast'); 15 grid on; 16 17 title('[0,1]的n=8的插值函数及原函数图像'); 18 xlabel('x轴'); 19 ylabel('y轴');
小结
在解第2题时，发现似乎在小区间的插值，拉格朗日的插值效果要比三次样条的插值效果更好，到大区间时，三次样条的插值效果要比拉格朗日的插值效果更好。
