matlab定积分的近似计算 实验报告二

《matlab与数学实验》实验报告实验序号：实验二日期：2015年05月09日班级132132002姓名高馨学号1321320041实验名称定积分的近似计算问题背景描述定积分的很多应用问题中，被积函数甚至没有解析表达式，可能只是一条实验记录曲线，或者是一组离散的采样值，这时只能应用近似方法去计算相应的定积分．实验目的本实验将主要研究定积分的三种近似计算算法：矩形法、梯形法、抛物线法．对于定积分的近似数值计算实验原理与数学模型MATLAB7.11.0主要内容（要点）实现实验内容中的例子，即分别采用矩形法、梯形法、抛物线法计算，取，并比较三种方法的精确程度．分别用梯形法与抛物线法，计算，取．并尝试直接使用函数trapz()、quad()进行计算求解，比较结果的差异．将的近似计算结果与Matlab中各命令的计算结果相比较，试猜测Matlab中的数值积分命令最可能采用了哪一种近似计算方法？并找出其他例子支持你的观点．学习fulu2sum.m的程序设计方法，尝试用函数sum改写附录1和附录3的程序，避免for循环。实验过程记录（含基本步骤、主要程序清单及异常情况记录等）实现实验内容中的例子，即分别采用矩形法、梯形法、抛物线法计算，取，并比较三种方法的精确程度．矩形法（中点法）程序：clearformatlong;a=0;b=1;n=258;h=(b-a)/n;C=0;fori=1:n;xi=(i-1)*h;xj=i*h;xz=(xi+xj)/2;C=C+h*(1/(1+xz*xz));enddisp(C)E=abs((C-pi/4)/(pi/4));fprintf('TherelativeerrorbetweenCandreal-valueisabout:%d\n',E)答案：0.785398476379441TherelativeerrorbetweenCandreal-valueisabout:3.985010e-007（二）梯形法程序：clearformatlong;a=0;b=1;n=258;h=(b-a)/n;B=0;x1=0;y1=1/(1+x1*x1);B=h*y1/2;fori=1:n-1xi=i*h;fxi=1/(1+xi*xi);B=B+h*fxi;endxn=1;yn=1/(1+xn*xn);B=B+h*yn/2;disp(B)E=abs((B-pi/4)/(pi/4));fprintf('TherelativeerrorbetweenBandreal-valueisabout:%d\n',E)答案：0.785397537433464TherelativeerrorbetweenBandreal-valueisabout:7.970021e-007（三）抛物线法程序：clearformatlonga=0;b=1;n=258;h=(b-a)/n;A=0;fori=1:nxj=(i-1)*h;xi=i*h;xk=(xi+xj)/2;fxi=1/(1+xi*xi);fxj=1/(1+xj*xj);fxk=1/(1+xk*xk);A=A+(h/6)*(fxj+4*fxk+fxi);enddisp(A)E=abs((B-pi/4)/(pi/4));fprintf('TherelativeerrorbetweenAandreal-valueisabout:%d\n',E)答案：0.785398163397449TherelativeerrorbetweenAandreal-valueisabout:2.827160e-016从他们的相对误差值，我们可以看出，抛物线法精确程度最高，其次是矩形法，最后是梯形法。二、分别用梯形法与抛物线法，计算，取．并尝试直接使用函数trapz()、quad()进行计算求解，比较结果的差异．梯形法：程序clearformatlong;a=1;b=2;n=120;h=(b-a)/n;B=0;x1=1;y1=1/x1;B=h*y1/2;fori=1:n-1xi=1+i*h;fxi=1/xi;B=B+h*fxi;endxn=2;yn=1/xn;B=B+h*yn/2;disp(B)答案0.684887057990131（2）trapzx=1:1/120:2;y=1./x;trapz(x,y)答案：ans=0.693151520800048（二）抛物线法：（1）程序clearformatlonga=1;b=2;n=120;h=(b-a)/n;A=0;fori=1:nxj=1+(i-1)*h;xi=1+i*h;xk=(xi+xj)/2;fxi=1/xi;fxj=1/xj;fxk=1/xk;A=A+(h/6)*(fxj+4*fxk+fxi);enddisp(A)答案：0.684848377754341quadf=inline('1./x');I=quad(f,1,2)答案：I=0.693147199862970从他们的最终结果可以看出，用梯形法与抛物线法求解的结果比直接使用函数trapz()、quad()进行计算求解的结果都要来得小。使用函数trapz()与梯形法命令的结果接近，说明函数trapz()后台计算时用的是梯形法的算法。。。。。。同理）三、将的近似计算结果与Matlab中各命令的计算结果相比较，试猜测Matlab中的数值积分命令最可能采用了哪一种近似计算方法？并找出其他例子支持你的观点．Matlab中的数值积分命令最可能采用了中点法近似计算。如下的梯形法中的左点法、中点法与右点法中中点法更接近Matlab中各命令的计算结果。（1）clearformatlong;a=0;b=1;n=258;h=(b-a)/n;A=0;B=0;C=0;fori=1:n;xi=(i-1)*h;xj=i*h;xz=(xi+xj)/2;A=A+h*(1/(1+xi*xi));B=B+h*(1/(1+xj*xj));C=C+h*(1/(1+xz*xz));enddisp(A)disp(B)disp(C)答案：0.7863665296815260.7844285451854020.785398476379441（2）x=0:1/100:1;y=1./(1+x.^2);trapz(x,y)答案：ans=0.785393996730782可以看出中点法与答案更接近。（算出误差）4、学习fulu2sum.m的程序设计方法，尝试用函数sum改写附录1和附录3的程序，避免for循环．附录1formatlongn=100;a=0;b=1;inum1=0;inum2=0;inum3=0;symsxfxfx=1/(1+x^2);i=1:n；（不然会把结果直接输出来）xj=a+(i-1)*(b-a)/n;%左点xi=a+i*(b-a)/n;%右点fxj=subs(fx,'x',xj);%左点值fxi=subs(fx,'x',xi);%右点值fxij=subs(fx,'x',(xi+xj)/2);%中点值e=fxj*(b-a)/n;f=fxi*(b-a)/n;g=fxij*(b-a)/n;inum1=sum(e);inum2=sum(f);inum3=sum(g);integrate=int(fx,0,1)integrate=double(integrate)fprintf('Therelativeerrorbetweeninum1andreal-valueisabout:%d\n\n',...abs((inum1-integrate)/integrate))fprintf('Therelativeerrorbetweeninum2andreal-valueisabout:%d\n\n',...abs((inum2-integrate)/integrate))fprintf('Therelativeerrorbetweeninum3andreal-valueisabout:%d\n\n',...abs((inum3-integrate)/integrate))答案：integrate=pi/4integrate=0.785398163397448Therelativeerrorbetweeninum1andreal-valueisabout:3.177794e-003Therelativeerrorbetweeninum2andreal-valueisabout:3.188404e-003Therelativeerrorbetweeninum3andreal-valueisabout:2.652582e-006附录3clearformatlongn=100;a=0;b=1;inum=0;ie=0symsxfxfx=1/(1+x^2);fori=1:nxj=a+(i-1)*(b-a)/n;%左点xi=a+i*(b-a)/n;%右点xk=(xi+xj)/2;%中点fxj=subs(fx,'x',xj);fxi=subs(fx,'x',xi);fxk=subs(fx,'x',xk);ie=ie+(fxj+4*fxk+fxi)*(b-a)/(6*n);endinum=sum(ie);integrate=int(fx,0,1)integrate=double(integrate)fprintf('Therelativeerrorbetweeninumandreal-valueisabout:%d\n\n',...abs((inum-integrate)/integrate))答案：ie=0integrate=pi/4integrate=0.785398163397448Therelativeerrorbetweeninumandreal-valueisabout:2.827160e-016实验结果报告与实验总结其实最重要的就是弄清楚第一题