数值分析实验报告

分数：评卷人：研究生（数值分析）实验报告学号M201772118姓名刘良专业光学工程指导教师刘小英院（系、所）光学与电子信息学院2018年4月30日数值实验2.1多项式差值的振荡现象问题提出：考虑在一个固定的区间上用插值逼近一个函数，显然Lagrange插值中使用的节点越多，插值多项式的次数就越高，我们自然关心插值多项式的次数增加时，(x)nL是否也更加靠近被逼近的函数。Runge给出的一个例子是极其著名并富有启发性的，设区间[-1,1]上函数21(x)125fx实验内容：考虑区间[-1,1]的一个等距划分，分点为2x1,0,1,2,...,iiinn则拉格朗日插值多项式为201(x)(x).125nniiiLlx其中(x),i0,1,2,...,nil是n次Lagrange插值基函数。实验要求：（1）选择不断增大的分点数目n=2,3,···，画出原函数f(x)插值多项式函数nL(x)在[-1,1]上的图像，比较并分析实验结果。（2）选择其他的函数，例如定义在区间[-5,5]上的函数4(x),(x)arctanx,1xhgx重复上述的实验看其结果如何。1、MATLAB实验程序prompt={'如选择f(x),请输入f,如选择h(x),请输入h,如选择g(x),请输入g:'};word='请选择试验函数';answer=inputdlg(prompt,word,1,{'f'});choose=char(answer);prompt={'请输入插值多项式的次数：'};word='插值次数';answer1=inputdlg(prompt,word,1,{'10'});numb=str2num(char(answer1));if(numb1)errordlg('请输入大于0的自然数！')ruturn;endswitchchoosecase'f'f=@(x)1./(1+25*x.^2);a=-1;b=1;case'h'f=@(x)x./(1+x.^4);a=-5;b=5;case'g'f=@(x)atan(x);a=-5;b=5;endx0=linspace(a,b,numb+1);y0=f(x0);x=a:0.1:b;y=lagrange(x0,y0,x);clf;fplot(f,[ab],'r-');holdon;plot(x,y,'b--');xlabel('x');ylabel('精确函数（红色），插值函数（蓝色）');str1=strcat('插值次数为',answer1);str2=strcat('试验函数为',answer)title([str1,str2]);%lagrange插值函数functiony=lagrange(x0,y0,x)n=length(x0);m=length(x);fork=1:mz=x(k);s=0.0;fori=1:np=1.0;forj=1:nif(j~=i)p=p*(z-x0(j))/(x0(i)-x0(j));endends=s+p*y0(i);endy(k)=s;endend2、实验结果对f(x)，g(x)，h(x)，依次输入分点数目为n=2,3,5,8,10,13,15,20,25,其原函数与其插值多项式函数Ln的图像如下图所示：当实验函数选择为f(x)时，取不同分点数目的原函数和插值多项式函数的图像：当实验函数选择为h(x)时，取不同分点数目的原函数和插值多项式函数的图像：当实验函数选择为g(x)时，取不同分点数目的原函数和插值多项式函数的图像：3、结果分析（1）对于f(x)的实验结果分析可知：首先观察到，无论插值次数选择为奇数还是偶数，其插值函数都为偶函数，这是因为原函数是偶函数，在计算插值函数L(x)的过程中，其奇数次的项全为0。随着分点数目的增加，插值函数的次数和项数也在增加，其复杂性也在增加，同时计算次数也在增加，由于计算机位数限制而引入的计算误差也在增加。其次，观察插值次数由高到低的插值函数拟合情况，当插值次数较低（n=2,3,5）时，拉格朗日插值函数不能很好的与原函数曲线吻合，随着次数的增加插值函数越来越逼近与原函数，次数越高插值函数越吻合原函数，但是这只局限于插值区间的中间部分。后面可以观察到插值次数增加的时候，在靠近插值区间端点会出现极大的波动。当n=13时，其偏离程度为0.5左右，n=15时达到0.7，n=20时拟合情况最好，n=25时在端点的偏离程度达到4.5。因此，在拉格朗日基函数的高次插值中，在插值区间的边界部分插值函数会出现很大波动，明显偏离原函数，所以插值次数不宜过高，容易出现龙格现象，但是也不能选择的过低，不然拉格朗日插值函数在区间中部的拟合情况不好。（2）对于h(x)、g(x)的实验结果分析可知:h(x)和g(x)均为奇函数，插值函数也为奇函数。与f(x)的实验分析结果类似，在低次插值时，插值函数不能很好的拟合，在选择过高的插值次数时，同样出现龙格振荡现象。在等距节点的拉格朗日插值中，插值次数不易过高，在高次插值中龙格现象无法避免，而且超过一定次数后，龙格现象会越来越严重。为解决龙格现象可以采用分段低次插值或采用切比雪夫零点为插值点。数值实验3函数逼近与曲线拟合（1）编制以函数0{}knkx为基的多项式最小二乘拟合程序，并用与对表3.11中的数据作3次多项式最小二乘拟合。ix-1.0-0.50.00.51.01.52.0yi-4.447-0.4520.5510.048-0.4470.5494.552取权数i1，求拟合曲线n**0=kkkx中的参数{k}、平方误差2，并做离散数据{,iixy}的拟合函数*y(x)的图形。（3）编制正交化多项式最小二乘拟合程序，并用于求解上题中的3次多项式最小二乘拟合问题，作拟合曲线的图形，计算平方误差，并与上题结果进行比较。1、MATLAB实验程序answer=inputdlg({'请输入1选择实验3.1，否则输入2'},'实验选择窗口',1,{'1'});nb=str2num(char(answer));x0=-1.0:0.5:2.0;y0=[-4.447-0.4520.5510.048-0.4470.5494.552];n=3;%拟合次数w=[1111111];if(nb==1)alph=wlsf(x0,y0,w,n);%按照降次排列的系数矩阵y=polyval(alph,x0);%拟合多项式在x0处的值r=(y0-y)*(y0-y)';%平方误差x=-1:0.01:2;y=polyval(alph,x);plot(x,y,'r-');xlabel('x');ylabel('拟合曲线（红色）离散数据（*）');holdonplot(x0,y0,'*');title('离散数据的多项式拟合');gridon;else[a,b,c,alph,r]=oplsf(x0,y0,w,n);enddisp(['平方误差：',sprintf('%g',r)]);disp(['系数为(按降次)：',sprintf('%g\t',alph)])functionp=wlsf(x0,y0,w,n)%带权重的最小二乘拟合子函数formatshortA=zeros(n+1,n+1);fori=0:nforj=0:nA(i+1,j+1)=sum((x0.^(i+j))*w');endb(i+1)=sum((x0.^i.*y0)*w');enda=A\b';p=fliplr(a');end%g1:初始值为0矩阵；g2:初始值设置为一个长度1*m的单位一矩阵；g3:根据递推关系%v1:初始值：g(k-1)的内积；v2:g(k)的内积；v3:作为一个中间变量。%d(k):y与ψ的内积；c(k)=存储的是正交多项式的系数。function[a,b,c,alph,r]=oplsf(x,y,w,n)%正交化多项式最小二乘拟合子函数m=length(x);g1=0;g2=ones(1,m);%g2表示默认（1,1,1,1,1···，m+1）v2=g2.^2*w';d(1)=y*w';c(1)=d(1)/v2;fork=1:nxs=x.*g2.^2*w';a(k)=xs/v2;%表示gk（s2）的内积与x的乘积，v2表示g2的内积if(k==1)b(k)=0;elseb(k)=v2/v1;endg3=(x-a(k)).*g2-b(k).*g1;v3=g3.^2*w';d(k+1)=y.*g3*w';c(k+1)=d(k+1)/v3;g1=g2;g2=g3;v1=v2;v2=v3;%改变递推值endr=y.*y*w'-c*d';%平方误差alph=zeros(1,n+1);t=zeros(n+1,n+2);t(:,2)=ones(n+1,1);t(2,3)=-a(1);if(n=2)fork=3:n+1fori=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);endendendfori=1:n+1fork=i:n+1alph(n+2-i)=alph(n+2-i)+c(k)*t(k,k+2-i);endendxmin=min(x);xmax=max(x);dx=(xmax-xmin)/(25*m);t=(xmin-dx):dx:(xmax+dx);s=alph(1);fork=2:n+1s=s.*t+alph(k);endplot(x,y,'*',t,s,'-');title('正交多项式拟合');xlabel('x');ylabel('y');gridon;end2、实验结果对表中的数据作3次最小二乘拟合，其离散数据的拟合函数的图形如图3.1所示，输出的平方误差为：2.17619e-05，系数为(按降次)：1.99911-2.99767-3.96825e-050.549119。对表中的数据通过正交多项式3次最小二乘拟合，其离散数据的拟合函数的图形如图3.2所示，输出平方误差：2.17619e-05。系数为(按降次)：1.99911-2.99767-3.96825e-050.549119当拟合次数n取不同的值时，其平方误差得到如下表格：n23456平方误差13.4882.17619e-051.6671e-051.63734e-051.2374e-253、结果分析最小二乘拟合是在离散情况下的最佳平方逼近。从图中我们可以看到，在拟合次数为3次的情况下，各个离散点均能落在拟合曲线上，拟合效果很好，拟合的平方误差为10-5量级。利用多项式进行最小二乘拟合和利用正交多项式进行最小二乘拟合的结果对比发现，两种算法算出的拟合函数系数一致，平方误差一致，所得的结果完全一致，实际上这也是必然的。在算法原理上，用一般多项式进行最小二乘拟合时，其系数矩阵是一个正规方程组，而利用正交多项式的解法时，无非是为了避免多项式拟合时求解病态法方程组，通过寻找一组正交的基函数，使得法方程的系数矩阵变成一个对角阵，使得法方程组得到简化，利于求解。本质上还是同一种算法，所得的结果必然是一致的。由于题中给出了7个离散点，通过依次改变拟合次数n，观察计算得到平方误差的表，在拟合次数越高时（小于离散点数目），其平方误差越小，拟合效果也越好。数值实验4.1数值积分实验目的：复化积分公式计算定积分。实验题目：数值计算下列各式右端定积分的近似值。3221(1)ln2ln32;1dxx1201(2)41dxx102(3)3ln3xdx221(4)exxedx实验要求：（1）若用复化梯形公式、复化Simpson公式和复化Gauss-legendreI型公式做计算，要求绝对误差限为-71=102，分别利用他们的余项对每种算法做出步长的事前估计。（2）分别用复化梯形公式、复化Simpson公式和复化Gauss-legendreI型公式做计算。（3）将计算结果与精确解作比较，并比较各种算法的计算量。1、MATLAB实验程序answer=inputdlg({'