第11章实验十数值积分.docVIP

第11章 实验十数值积分 实验目的：理解Newton-Cotes求积公式，会用复化梯形公式、复化辛普生公式计算 11.1 Newton-Cotes求积公式 设，为等距节点，且，前四个求积公式为： 当充分可微，则Newton-Cotes求积公式的余项是一更高阶导数项。梯形公式的精度为n=1，若，则： 辛普生公式的精度为n=3，若，则： 辛普生3/8公式的精度为n=3，若，则： 布尔公式的精度为n=5，若，则： 例11.1 考虑函数，等距节点为 ， 对应的函数值为 利用Newton-Cotes求积前四个公式，步长为h=0.5，计算得： 例11.2. 考虑函数，在固定区间[a,b]=[0,1]内的积分，使用Newton-Cotes求积前四个公式。 解：对梯形公式，h=1，且有 对辛普生公式，h=1/2，且有 对辛普生3/8公式，h=1/3，且有 对布尔公式，h=1/4，且有 该积分的真解为：=1.3082506046426…. 11.2 复化梯形公式 通过的n+1个等步长节点， 逼近积分。这里，。 程序11.1 function I=trapez_v(f,h) I=h*(sum(f)-(f(1)+f(length(f)))/2); 功能：用复化求积公式进行函数积分 调用格式：I=trapez_v(f,h) f：等距节点上的函数值序列 h：步长 例11.3 已知积分精确值=4.006994,分析用复化梯形法计算下述积分时，剖分区间数对误差的影响。 解：程序如下： clear Iexact=4.006994; a=0;b=2; fprintf(\n Extended Trapezoidal Rule\n); fprintf( n I Error\n); n=1; for k=1:6, n=2*n; h=(b-a)/n; i=1:n+1; x=a+(i-1)*h; f=sqrt(1+exp(x)); I=trapez_v(f,h); I=h*(sum(f)-(f(1)+f(length(f)))/2); fprintf(%3.0f %10.5f %10.5f\n, n, I, Iexact-I); end 执行程序的结果为： Extended Trapezoidal Rule n I Error 2 4.08358 -0.07659 4 4.02619 -0.01919 8 4.01180 -0.00480 16 4.00819 -0.00120 32 4.00729 -0.00030 64 4.00707 -0.00008 程序11.2 function I=trapez_n(f_name,a,b,n) h=(b-a)/n; x=a+(0:n)*h; f=feval(f_name,x); I=trapez_v(f,h) function I=trapez_g(f_name,a,b,n) n=n;hold off h=(b-a)/n; x=a+(0:n)*h; f=feval(f_name,x); I=h/2*(f(1)+f(n+1)); if n1 I=I+h*sum(f(2:n));end h2=(b-a)/100; xc=a+(0:100)*h2; fc=feval(f_name,xc); plot(xc,fc,r); hold on title(Trapezoidal Rule); xlabel(x);ylabel(y); plot(x,f) plot(x,zeros(size(x))) for i=1:n+1; plot([x(i),x(i),],[0,f(i)]);end 功能：用复化求积公式进行函数积分，trapez_n不画出图形，而trapez_g画出图形 调用格式：trapez_v(‘f_name’,a,b,n)或trapez_g(‘f_name’,a,b,n) f_name：被积函数的文件名 a：x的上限 b：x的上限 n：剖分区间数 例11.4 已知积分精确值=4.006994,用trapez_g画出图形并计算出值