拉格朗日插值

拉格朗日插值绘制龙格现象

一、问题叙述

龙格反例1/(1+x^2)说明高次代数插值会导致误差很大。在区间[-5,5]上取等距结点构造10次拉格朗日插值多项式用计算机绘制图形显示龙格现象。 二、理论分析

1. 拉格朗日插值：假设有（n+1）个拉格朗日插值结点x0?x1??xn ，已知函数值

y0?f(x0),y1?f(x1),,yn?f(xn)

求n次多项式Ln(x)使其满足插值条件f(xj)?yj(j?0,1,,n)

类似于二次插值方法，根据插值结点构造(n+1)个拉格朗日插值基函数

lk(x)?(x?x0)?(x?xk?1)(x?xk?1)?(x?xn)

(xk?x0)?(xk?xk?1)(xk?xk?1)?(xk?xn)?1j?k每一个基函数都是零点多项式lk(xj)??,(j?0,1n)

0j?k?Ln(x)满足插值条件 Ln(xj)?f(xj)拉格朗日插值基函数：lk(x)??j?0j?kn(j?0,1,,n)

(x?xj)(xk?xj)拉格朗日插值多项式：Ln??lj(x)yj

j?0n2. 切比雪夫插值：n阶切比雪夫多项式定义为

Tn(x)?cos(narccosx),(n?0,1,),x?[?1,1]

若令??arccosx ，则有Tn?cosn? 。由余弦函数性质，有

cos(n?1)??cos(n?1)??2cos?cosn?

所以有递推关系 Tn?1(x)?2xTn(x)?Tn?1(x),(n?1,2比雪夫多项式零点为xk?cos 3. Hermite插值

如果f(x)在区间[a,b]上连续可导，x0,x1,L xn?[a,b] 是互异的，那么存在唯一的多项式H2n?1(x)满足多项式在这些点上的值与函数f(x)的值相等、多项式在这些点的一阶导数值与函数的一阶导数值相等。 这个多项式可以表示为

)，又有T0?1,T1?x，所以n阶切

2k?1?,(k?0,1,2n,n?1) 。

H2n?1(x)??f(xi)[1?2(x?xi)li'(x)li2(x)

i?0n??f'(xi)(x?x)li2(x)i?0n

其中 li(x)??j?0j?knn(x?xj)(xk?xj),i?0,1,L,n

li'(x)??j?0j?k1,i?0,1,L,n xi?xj三、算法MATLAB实现

（1）拉格朗日插值描绘龙格现象，代码如下：

function f = Language(x,y,x0) syms t l;

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); end 主程序 x1=-5:1:5; x2=-5:2:5; x3=-5:2/3:5;

y11=1./(1+x1.^2);次拉格朗日插值 y12=1./(1+x2.^2);%5次拉格朗日插值 y13=1./(1+x3.^2);次拉格朗日插值 x0=-5:0.001:5; %调用拉格朗日函数

y1=Language(x1,y11,x0); y2=Language(x2,y12,x0); y3=Language(x3,y13,x0); y0=1./(1+x0.^2); plot(x0,y0); hold on

plot(x0,y2,'r'); hold on

plot(x0,y1,'k'); hold on

plot(x0,y3,'g'); hold on xlabel('x'); ylabel('y');

title('原函数f(x)=1/(1+x^2) 等距拉格朗日插值');

legend('原函数','5次拉格朗日插值','10次拉格朗日插值','15次拉格朗日插值'); grid on

插值效果如下：

图1 等距拉格朗日插值

图1分别进行5次、10次、15次拉格朗日插值，可以看出5次插值拟合效果不是很好，10次插值在插值区间的边界处出现很大波动，明显偏离原函数，15次插值看到波动情况加强，故得出结论：拉格朗日插值次数不宜过高。 高次插值边界出现这种波动现象叫做龙格现象。

为避免上述现象，分别采用切比雪夫插值，埃米特插值和样条插值来避免龙格现象 （2）切比雪夫插值：用切比雪夫多项式零点代替等距结点，其他基本不变

x1=-5:1:5; y11=1./(1+x1.^2); x0=-5:0.001:5;

y1=Language(x1,y11,x0); y0=1./(1+x0.^2); k=0:1:10;

xx=5*cos((2*k+1)*pi/22); % 用切比雪夫多项式零点代替等距结点 yy=1./(1+xx.^2); y4=Language(xx,yy,x0); plot(x0,y0,'k'); hold on

plot(x0,y1,'r'); hold on

plot(x0,y4,'b',xx,yy,'ob');

插值效果如下：

图2 切比雪夫插值

由图2可以看出10次切比雪夫插值很好的和原函数拟合，并且波动较小，说明伪振荡现象得到控制。

（3）样条插值：样条插值直接采用MATLAB插值函数spline()，代码如下： x=-5:1:5;

y=1./(1+x.^2); xx=-5:0.01:5;

yy=spline(x,y,xx); %使用样条插值方法 x0=-5:0.001:5; y0=1./(1+x0.^2);

plot(x,y,'ok',xx,yy,'r',x0,y0,'b');

插值效果如下：

图3 样条插值

由图3可以看出样条插值拟合效果最好，与原函数基本一致，但样条插值处理1/（1+25*x^2） 函数时，还是有波动的，故不能完全说明样条是最好的插值方式 （4）Hermite插值：根据埃米特插值表达式写出代码如下：

function yy=Hermite(x,y,dy,xx) n=length(y); m=length(x); l=length(dy); k=length(xx);

if m~=n,error('向量长度不一致'); end;

if n~=l,error('向量长度不一致'); end;

z=zeros(1,k); for j=1:k s=0; for t=1:m; a=0;b=1; for i=1:n; if x(t)~=x(i)

a=a+1/(x(t)-x(i));

b=b*((xx(j)-x(i))/(x(t)-x(i))); end end

s=s+(y(t)*(1-2*(xx(j)-x(t))*a)*b^2+dy(t)*(xx(j)-x(t))*b^2); end

z(j)=s; end yy=z;

插值效果如下：

图4 Hermite插值

由图4可以看出，10个Hermite插值结点在区间内与原函数拟合程度非常好，基本能解决伪振荡现象。

综上所述可以得出结论：等距拉格朗日插值会产生伪振荡现象，插值次数不宜过高，切比雪夫插值，埃米特插值，样条插值可以很好的解决龙格函数的伪振荡现象。

本文来源：https://www.bwwdw.com/article/t4rf.html