拉格朗日插值实验报告

实验名称： 实验一 拉格朗日插值

1 引言

我们在生产生活中常常会遇到这样的问题：某个实际问题中，函数f(x)在区间[a,b]上存在且连续，但却很难找到其表达式，只能通过实验和观测得到有限点上的函数表。显然，根据这些点的函数值来求其它点的函数值是非常困难的。有些情况虽然可以写出表达式，但结构复杂，使用不方便。所以我们总是希望根据已有的数据点（或函数表）来构造某个简单函数P(x)作为f(x)的近似值。插值法是解决此类问题的一种比较古老的、但却很常用的方法。它不仅直接广泛地应用于生产实际和科学研究中，而且也是进一步学习数值计算方法的基础。

2 实验目的和要求

运用Matlab编写三个.m文件，定义三种插值函数，要求一次性输入整函数表，并利用计算机选择在插值计算中所需的节点。分别通过分段线性插值、分段二次插值和全区间上拉格朗日插值计算f(0.15)，f(0.31)，f(0.47)的近似值。已知函数表如下：

x f(x) 0.0 0.39894 0.1 0.39695 0.195 0.39142 0.3 0.38138 0.401 0.36812 0.5 0.35206 3 算法原理与流程图

（1）原理

设函数y=在插值区间[a,b]上连续，且在n+1个不同的插值节点a≤x0,x1,…,xn≤b上分别取值y0,y1,…,yn。目的是要在一个性质优良、便于计算的插值函数类Φ中，求一简单函数P(x)，满足插值条件P(xi)=yi(i=0,1,…,n)，而在其他点x≠xi上，作为f(x)近似值。求插值函数P(x)的方法称为插值法。在本实验中，采用拉格朗日插值法。 ①分段低次插值

当给定了n+1个点x0f(x)P1(x)yi1xxixxi1yi

xi1xixixi1这种分段低次插值叫分段线性插值，又称折线插值。

类似地，我们可以选取距离x最近的三个节点xi-1，xi与xi+1，然后进行二次插值，即得

xxji1f(x)P2(x)yk1xkxjki1jjiki1 这种分段低次插值叫分段二次插值，又称分段抛物线插值。

②全区间上拉格朗日插值

对节点xi(i=0,1,…,n)中任一点xk(0≤k≤n)，作一n次多项式lk(x)，使它在该点上的取值为

. ... .c

... ...

1，在其余点xi(i=0,1,…,k-1,k+1,…,n)上取值为零。对应于每一节点xk(k=0,1,…,n)，都能写出一个满足此条件的多项式，这样写出了n+1个多项式l0(x),l1(x),…,ln(x)，其中

lk(x)Ak(xx0)(xx1)(xxk1)(xxk1)由条件lk(xk)1可得

•(xxn)；

Ak(xkx0)1(xkxk1)(xkxk1)(xkxn)

于是我们可以得出如下的拉格朗日n次插值多项式（对于全区间上的插值，n取函数表的长度）

Pn(x)y0l0(x)y1l1(x)k0ynln(x)ykn(xx0)(xxk1)(xxk1)(xxn)(xkx0)(xkxk1)(xkxk1)(xkxn)

（2）流程图

分段线性插值

分段二次插值 全区间拉格朗日插值

4 程序代码及注释

1、分段线性插值 %分段线性插值 function y=piece_linear(x0,y0,x) % x0，y0为已知点，x为待求点 n=length(x0);p=length(y0);m=length(x); . ... .c

... ...

% n，p，m分别为x0，y0，x长度 if n~=p fprintf('Error! Please input again!\\n'); % x0和y0长度不等时，报错 else for i=1:m z=x(i); sum=0.0; l=0; %给l赋初值,根据x的值确定l if zx0(n) fprintf('Error!x(%d) is out of range!\\n',i); break; end %当插值点超出围时，报错 for j=2:n if z ... ...

function y=piece_square(x0,y0,x) % x0，y0为已知点，x为待求点 n=length(x0);p=length(y0);m=length(x); % n，p，m分别为x0，y0，x长度 if n~=p fprintf('Error! Please input again!\\n'); % x0和y0长度不等时，报错 else for i=1:m z=x(i); sum=0.0; l=0; %给l赋初值,根据x的值确定l if zx0(n) fprintf('Error!x(%d) is out of range!\\n',i); break; end %当插值点超出围时，报错 for j=1:n-2 p=0.5*(x0(j)+x0(j+1)); if z

... ...

%输出插值结果与所需节点 end end end 3、拉格朗日全区间插值 %拉格朗日全区间插值 function y=lagrange(x0,y0,x) % x0，y0为已知点，x为待求点 n=length(x0);p=length(y0);m=length(x); %n，p，m分别为x0，y0，x长度 if n~=p fprintf('Error! Please input again!\\n'); %x0和y0长度不等时，报错 else for i=1:m z=x(i); s=0.0; if zx0(n) fprintf('Error!x(%d) is out of range!\\n',i); break; end %当插值点超出围时，报错 for k=1:n p=1.0; for j=1:n if j~=k p=p*(z-x0(j))/(x0(k)-x0(j)); end end s=p*y0(k)+s; end y(i)=s; fprintf('y(%d)=%.5f\\n',i,y(i)); %输出插值结果 end end end 5 算例分析

1、 测试示例

. ... .c

... ...

>> x=[1 2 3 4]; >> y=[2 3 4]; >> y2=lagrange(x,y,x0) Error! Please input again! >> x=[1 2 3 4]; >> y=[2 3 4 5]; >> x0=[0.5 5.5]; >> y2=lagrange(x,y,x0) Error!x(1) is out of range! >> x=[1 2 3 4]; >> y=[2 3 4 5]; >> x0=[1.5 5.5]; >> y2=lagrange(x,y,x0) y(1)=2.50000 Error!x(2) is out of range! y2 = 2.500000000000000 2、首先输入函数变及待求点 >> x=[0.0 0.1 0.195 0.3 0.401 0.5]; >> y=[0.39894 0.39695 0.39142 0.38138 0.36812 0.35206]; >> x0=[0.15 0.31 0.47]; 注：保证在matlab工作目录中有三个.m文件 3、分段线性插值

y0=piece_linear(x,y,x0) y(1)=0.394039 x1=0.100 y1=0.39695,x2=0.195 y2=0.39142 y(2)=0.380067 x1=0.300 y1=0.38138,x2=0.401 y2=0.36812 y(3)=0.356927 x1=0.401 y1=0.36812,x2=0.500 y2=0.35206 y0 = 0.394039473684211 0.380067128712871 0.356926666666667 . ... .c

... ...

4、分段二次插值

>> y1=piece_square(x,y,x0) y(1)=0.394460 x1=0.100 y1=0.39695 x2=0.195 y2=0.39142 x3=0.300 y3=0.38138 y(2)=0.380225 x1=0.195 y1=0.39142 x2=0.300 y2=0.38138 x3=0.401 y3=0.36812 y(3)=0.357247 x1=0.300 y1=0.38138 x2=0.401 y2=0.36812 x3=0.500 y3=0.35206 y1 = 0.394460319548872 0.380224691595373 0.357246844884488 5、全区间拉格朗日插值

>> y2=lagrange(x,y,x0) y(1)=0.39447 y(2)=0.38022 y(3)=0.35722 y2 = 0.394472803878061 0.380219062454732 0.357222112339485 6 讨论与结论

1、使用tic,toc函数计算下列四种方法计算上述问题所运行的时间 Function 运行时间(s) lagrange(x0,y0,x) 0.000272 piece_linear(x0,y0,x) 0.000375 piece_square(x0,y0,x) 0.000272 从三次实验结果可知，三个程序的运行时间都很短。 2、程序优化

由分段线性插值和分段二次插值的原理，x取值在函数表围时，插值结果有意义，而当x

. ... .c

... ...

取值在函数表围以外，利用分段线性插值公式仍可以进行运算并得到一个值，但其结果不准确；分段二次插值则无法找到三个合适的点以求插值，不予以输出结果；若输入的函数表x与y的长度不相等，则无法插值。所以加入以下判断以提高插值的准确性 n=length(x0);p=length(y0);m=length(x); if n~=p fprintf('Error! Please input again!\\n'); if zx0(n) fprintf('Error!x(%d) is out of range!\\n',i); break; end 3、作图比较

上图为三种方法的插值曲线，其中x取0到0.5，步长为0.001，由图可得，三种曲线非常接近，这说明我们用拉格朗日插值计算所给点函数值的近似值时，引起的误差还是比较小的。

参考文献

[1] 易大义，云宝，有法. 计算方法(第2版)，大学. p.29-53.

[2] 琨 高思超 毕靖 编著 MATLAB2010从入门到精通 电子工业

. ... .c