试验2.1 多项式插值的振荡现象实验目的:观察多项式插值的振荡现象,了解多项式的次数与逼近效果的关系。
实验内容:问题提出:考虑在一个固定的区间上用插值逼近一个函数。
显然Lagrange 插值中使用的节点越多,插值多项式的次数就越高。
我们自然关心插值多项式的次数增加时,Ln(x)是否也更加靠近被逼近的函数。
Runge 给出的一个例子是极著名并富有启发性的。
设区间[-1,1]上的函数225x11)x (+=f ,考虑区间[-1,1]上的一个等距划分,分点为n2i 1x i +-=,i=0,1,2,…,n则拉格朗日插值多项式为:)x (l 25x11)x (Ln i ni 2i∑=+=,其中的)x (l i ,i=0,1,2,…,n 是n 次拉格朗日插值基函数。
实验要求:1、选择不断增大的分点数目n=2,3,………,画出原函数)x (f 及插值多项式函数)x (Ln 在[-1,1]上的图像,比较并分析试验结果。
2、选择其他的函数,例如定义在区间[-5,5]上的函数4()1x h x x=+,()arctan g x x =,重复上述的实验看其结果如何。
实验步骤及结果分析:1、选择不断增大的分点数目n=2,3,4,5,6,7,8,9,10做)x (f 的拉格朗日插值多项式)x (Ln ,并与原函数值做比较,如下图所示。
观察图像可知:n=2,3时插值函数和原函数差别很大,n=4,5,6时插值函数与原函数的逼近程度相对较好,继续增加插值次数n ,插值函数在插值区域的中间部分收敛,而在这区间外是发散的,此外,n=7,9时在插值中间区域逼近效果不好。
因此,适当提高插值多项式次数,可以提高逼近的精度,但是次数太高反而产生相反的效果。
2、选择其他的函数进行插值。
原函数4()1x h x x=+,区间[-5,5],插值结果如下图:观察图像可知:低次插值时,插值效果不好。
n=7,8,9,10时,在区间[-2,2],插值函数与原函数逼近程度好,但在区间外插值函数发散。
其中,n=8,10插值效果比n=7,9好。
原函数()arctan,区间[-5,5],插值结果如下图:g x x观察图像可知:n=5,6,7,8,9,10时,在区间[-3,3],插值函数与原函数逼近程度好,但在区间外插值函数发散。
其中,n=7,9在插值区间两端发散的更剧烈。
分析在插值区间两端发散的原因:次数越高,计算量就越大,积累误差也大,在整个区间上做高次多项式,当局部插值节点处的值有微小偏差时,可能引起整个区间上函数值的很大变化,使计算不稳定。
Matlab程序如下:function t_charpt2%数值实验二%输入:实验选择,函数式选择,插值节点数%输出:拟合函数及原函数的图形result=inputdlg({'请选择实验,若选2.1,请输入1,否则输入2:'},'charpt_2',1,{'1'}); Nb=str2num(char(result));if(Nb~=1)&(Nb~=2) errordlg('实验选择错误!');return;endpromps={'请选择实验函数,若选f(x),请输入f,若选h(x),请输入h,若选g(x),请输入g:'};result=inputdlg(promps,'charpt 2',1,{'f'});Nb_f=char(result);if(Nb_f~='f'&Nb_f~='h'&Nb_f~='g')errordlg('实验选择错误!');return;endresult=inputdlg({'请输入插值点数N:'},'charpt_2',1,{'10'});Nd=str2num(char(result));if(Nd<1)errordlg('结点输入错误!');return;endswitch Nb_fcase 'f'f=inline('1./(1+25*x.^2)');a=-1;b=1;case 'h'f=inline('x./(1+x.^4)');a=-5;b=5;case 'g'f=inline('atan(x)');a=-5;b=5;endif(Nb==1)x0=linspace(a,b,Nd+1);y0=feval(f,x0);x=a:0.1:b;y=lagrange(x0,y0,x);% clf;把曲线显示在一张图上fplot(f,[a,b],'co');hold on;plot(x,y,'b--');xlabel('x') ;ylabel('y=f(x)o and y=ln(x)--');elseif(Nb==2)x0=linspace(a,b,Nd+1);y0=feval(f,x0);x=a:0.1:b;cs=spline(x0,y0);y=ppval(cs,x);% clf;把曲线显示在一张图上fplot(f,[a,b],'co');hold on;plot(x,y,'k-');xlabel('x');ylabel('y=f(x) o and y=spline(x)-'); end%function y=lagrange(x0,y0,x)%Lagrange插值n=length(x0);m=length(x);for i=1:mz=x(i);s=0.0;for k=1:np=1.0;for j=1:nif (j~=k)p=p*(z-x0(j))/(x0(k)-x0(j));endends=s+p*y0(k);endy(i)=s;end试验3.1 多项式最小二乘拟合实验目的:掌握多项式最小二乘拟合程序,会对数据作多项式最小二乘拟合。
实验内容及要求:编制以函数{}nk k x 0=为基的多项式最小二乘拟合程序,并用于对表中的数据作3次多项式最小拟合。
取权数1≡i w ,求拟合曲线knk kxa∑=**=ϕ中的参数{}k a 、平方误差2δ,并作离散数据{}i i y x ,的拟合函数()x y *=ϕ的图形。
实验步骤及结果分析:拟合函数()x y *=ϕ的图形:拟合曲线knk kxa∑=**=ϕ中的平方误差2δ为2.17619e-005;参数{}k a 为alph 中的四个数值,按降幂排列1.99911,-2.99767,-3.96825e-005,0.549119。
观察图像可知:表中数据均在拟合曲线上,拟合结果很好,平方误差为10-5数量级。
Matlab程序如下:function charpt3%数值实验三:含“实验3.1”和“实验3.2”%子函数调用:dlsa%输入:实验选择%输出:原函数及求得的相应插值多项式的函数的图像以及参数alph和误差rresult=inputdlg({'请选择实验,若选3.1,请输入1,否则输入2:'},'charpt_3',1,{'1'}); Nb=str2num(char(result));if((Nb~=1)&(Nb~=2)) errordlg('实验选择错误!'); return; endx0=-1:0.5:2;y0=[-4.447 -0.452 0.551 0.048 -0.447 0.549 4.552];n=3;%n为拟合阶次if Nb==1alph=polyfit(x0,y0,n);y=polyval(alph,x0);r=(y0-y) * (y0-y)';%平方误差x=-1:0.01:2;y=polyval(alph,x);plot(x,y,'k--');xlabel('x'); ylabel('y0 * and polyfit.y--');hold onplot(x0,y0,'*');title('离散数据的多项式拟合');grid on;elseresult=inputdlg({'请输入权向量w:'},'charpt_3',1,{'[1 1 1 1 1 1 1]'});w=str2num(char(result));[a,b,c,alph,r]=dlsa(x0,y0,w,n);enddisp(['平方误差:',sprintf('%g',r)]);disp(['参数alph:', sprintf('%g\t',alph)])function [a,b,c,alph,r]=dlsa(x,y,w,n)%功能:用正交化方法对离散数据作多项式最小二乘拟合。
%输入:m+1个离散点(x,y,w),x,y,w分别用行向量给出。
% 拟合多项式的次数n,0〈n<m。
%输出:三项递推公式的参数a,b,拟合多项式s(x)的系数c和alph,% 平方误差r=(y-s,y-s),并作离散点列和拟合曲线的图形m=length(x)-1;if (n<1|n>=m) errordlg('错误:n<1或者n>=m!'); return;end%求三项递推公式的参数a,b,拟合多项式s(x)的系数c,其中d(k)=(y,sk);s1=0;s2=ones(1,m+1);v2=sum(w);d(1)=y*w';c(1)=d(1)/v2;for k=1:nxs=x.*s2.^2*w';a(k)=xs/v2;if k==1 b(k)=0;else b(k)=v2/v1;ends3=(x-a(k)).*s2-b(k)*s1;v3=s3.^2*w';d(k+1)=y.*s3*w';c(k+1)=d(k+1)/v3;s1=s2;s2=s3;v1=v2;v2=v3;end%求平方误差rr=y.^2*w'-c*d';%求拟合多项式s(x)的降幂系数alphalph=zeros(1,n+1); T=zeros(n+1,n+2);T(:,2)=ones(n+1,1); T(2,3)=-a(1);if n>=2for k=3:n+1for i=3:k+1T(k,i)=T(k-1,i)-a(k-1)*T(k-1,i-1)-b(k-1)*T(k-2,i-2);endendendfor i=1:n+1for k=i:n+1alph(n+2-i)=alph(n+2-i)+c(k)*T(k,k+2-i);endend%用秦九韶方法计算s(x)的输出序列(t,s)xmin=min(x); xmax=max(x); dx=(xmax-xmin)/(25*m);t=(xmin-dx):dx:(xmax+dx);s=alph(1);for k=2:n+1s=s.*t+alph(k);end%输出点列x-y和拟合曲线t-s的图形plot(t,s,'-',x,y,'*');title('离散数据的正交化多项式拟合');xlabel('x');ylabel('y');grid on;试验5.1 常微分方程性态和R-K 法稳定性试验实验目的:考察下面微分方程右端项中函数y 前面的参数对方程性态的影响(它可使方程为好条件的或坏条件的)和研究计算步长对R-K 法计算稳定性的影响。