1、数值分析课程设计报告专业课程设计成绩评定书指导教师评语成 绩： 指导教师 时 间： 答辩小组意见 设计成绩： 答辩组长： 2审定系主任： 题目：用熟悉的计算机语言编程上机完成（1）用 Newton-Cotes 公式计算积分 的近似值，自己设置不同精度要求，102dxe对结果进行比较分析。（2）用 Romberg 积分法计算积分 的近似值，自己设置不同精度要求，102xe对结果进行比较分析；与（1）的结果进行比较分析，谈谈你的体会。（3）记 ，在上面的计算中 只取 4 位有效数字或 7 位有效数字，2)(xef)(f计算结果有什么不同。（4）上面计算精度可达 8-20 位有效数字吗？若可以请说明
2、实现过程，并举例。一摘要：在 matlab 环境下熟悉的运用计算机编程语言并结合 Newton-Cotes和 Romberg 的理论基础对函数求积分，在运行完程序后以及对运行结果做出各方面的分析和比较。二实验设计目地用熟悉的计算机语言编程上机完成 Newton-Cotes 公式计算积分和Romberg 积分法计算积分。三Newton-Cotes 和 Romberg 的理论基础牛顿柯斯特公式：设将积分区间a,b划分为 n 等份，步长 h=（b-a）/n,选取3等距节点 =a+kh 构造出的插值型求积分公式 I ，称为牛kx )(0)kfabnkCn顿柯斯特（Newton-Cotes）公式，式中
3、称为柯特斯系数，引进变换k)(x=a+th,则有 dtjndtjtabhnkjknkjnkC)()!()()/( 00)( 1 龙贝格求积公式：（1）.梯形法的递推化 设将区间a,b分为 n 等份，共有 n+1 个分点，如果将求积区间再二分一次，则分点增至 2n+1 个，我们将二分前后两人积分值联系起来加以考察，注意到每个子区间 经过二分只增加了一个分点xk1,，用复合梯形公式求得该子区间上的积分值为)(211xkk，注意，这里 h =(b-a)/n 代表二分前的步长，将)()412kfffh每个子区间上的积分值相加得 ，从而1021012 )()()4nkknkkk xxTfhff利用上式可
4、导出下列递推分式：022nkknfh（2）.外推技巧 从梯形公式出发，将区间a,b逐次二分可提高求积公式精度，当a,b分为 n 等份 I nabhabfTnn ,),(12若记 ，当区间a,b分为 2n 等份时，则有 ，并且有)(hn )2(TI.)0(),(12)(lim0 ThabITf(3).龙贝格算法 将上述外推技巧得到的公式重新引入记号 ),(0hT)(1hS等从而可将上述公式定成统一形式,2C)(3hR经过 m(m=1.2.3)次加速后，余项便)(121)(4TTmmm 取下列形式： .)( )1(2)1(2hmIh4这方法通常称为理查森外推加速方法。设以 表示二分次后求得的梯形值
5、，且以 表示序列 的 m 次加速Tkm)( Tkm)( k)(值，则依递推公式可得 .21,1)()()( 4kkmk公式也称为龙贝格算法，计算过程如下：.取 k=0,h=b-a,求 )(2)0( bfahT令 （k 记区间a,b二分次数）1.求梯形值 ，即按递推公式计算)(0kbTk)(0.求加速值，按公式 a 逐个求出 T 表的第 k 行其余各元素T)(jk（j=1,2,k）.iv.若|T)0(kT)(1|=8errorendxk=fun(1,:);fk=fun(2,:);a=min(xk);b=max(xk);n=mm-1;endif nargin=4xk=linspace(a,b,n+
6、1);if isa(fun,function_handle)6fx=fun(xk);end endCk=cotescoeff(n);Ak=(b-a)*Ck;y=Ak*fx;function Ck=cotescoeff(n)for i=1:n+1k=i-1;Ck(i)=(-1)(n-k)/factorial(k)/factorial(n-k)/n*quadl(t)intfun(t,n,k),0,n);endfunction f=intfun(t,n,k)f=1;for i=0:k-1,k+1:nf=f.*(t-i);end运行结果 y=mulNewtonCotes(0,1,2,1)xk =0 0
7、.5000 1.0000y =0.73147 vpa(y,7)ans =.7313703 y=mulNewtonCotes(0,1,5,2)xk =0 0.2000 0.4000 0.6000 0.8000 1.0000y =0.7468 vpa(y,7)ans =.746824龙贝格求积分function z=romberg(a,b,e)h=(b-a);f=(x)exp(-x.2);TT(1,1)=h.*(f(b)+f(a)/2;8k=2;TT(1,2)=TT(1,1)./2+h/2.*f(a+h/2);TT(2,1)=TT(1,2).*4/3-TT(1,1)./3;z=TT(2,1);wh
8、ile abs(TT(k,1)-TT(k-1,1)./TT(k,1)=ek=k+1;h=h./2;for j2=1:2.(k-2)ff(1,j2)=f(a+h*(j2-1/2);endfff=sum(ff).*h/2;TT(1,k)=TT(1,k-1)./2+ffffor j1=2:kTT(j1,k-j1+1)=4(j1-1).*TT(j1-1,k-j1+2)./(4(j1-1)-1)-TT(j1-1,k-j1+1)/(4(j1-1)-1);z=TT(j1,k-j1+1);endendfunction r=f()r=(x)exp(-x.2)运算结果 z=romberg(0,1,10e-2)z
9、=0.74729 vpa(z,7)ans =.7471804 z=romberg(0,1,10e-5)TT =0.6839 0.7314 0.74300.7472 0 0TT =0.6839 0.7314 0.7430 0.74590.7472 0.7469 0 00.7468 0 0 0z =0.7468 vpa(z,7)ans =10.7468240 z=romberg(0,1,10e-10)TT =0.6839 0.7314 0.74300.7472 0 0TT =0.6839 0.7314 0.7430 0.74590.7472 0.7469 0 00.7468 0 0 0TT =0.
10、6839 0.7314 0.7430 0.7459 0.74660.7472 0.7469 0.7468 0 00.7468 0.7468 0 0 00.7468 0 0 0 0TT =0.6839 0.7314 0.7430 0.7459 0.7466 0.74680.7472 0.7469 0.7468 0.7468 0 00.7468 0.7468 0.7468 0 0 00.7468 0.7468 0 0 0 00.7468 0 0 0 0 0z =110.7468 vpa(z,7)ans =.7468241五.结果分析根据牛顿柯斯特求积分的原理将所求函数的积分区间分为多个小的积分区间
11、，先求出每个小积分区间上的函数值，然后再将每个小区间上的求积结果加起来就是我们所要求的总函数的积分值，当函数区间所分的小区间的个数越多的时候总函数所计算出来的结果就精确，其原因就是所分的区间数越多，计算时每个小区所带来的误差就越小，当将总的积分值加起来的时候所带来的总的误差也就越小，所以最后的结果的精度越高，而用龙贝格求积和牛顿柯斯特也一样是要将总区间分为很多小的相等的区间，只是他们所用的计算原理不一样，当用此方法求积分时，所设的误差限越小，所求得的积分值就越是精确。六.参考文献1孙祥，徐流美，吴清。Matlabl 7.0 基础教程【M】 。北京：清华大学出版社，2005.2薛毅。数值分析与实验【M】.北京：北京理工大学出版社，2005.3汪卉琴，刘目楼。数值分析【M】 。北京：冶金工业出版社，2004.4丁丽娟，程杞元。数值计算方法【M】 。北京：北京理工大学出版社，2005.5肖伟，刘忠，曾新勇，等。Matlabl 程序设计与应用【M】 。北京：清华大学出版社，2005.6李庆扬，易大义，王能超。现代数值分析。北京：高等教育出版社，1995.
