数值分析课程设计baogao

1、数值分析课程设计报告专业课程设计成绩评定书 指导教师评语 成 绩： 指导教师 时 间： 答辩小组意见设计成绩： 答辩组长： 审定 系主任： 题目：用熟悉的计算机语言编程上机完成（1）用Newton-Cotes公式计算积分的近似值，自己设置不同精度要求，对结果进行比较分析。（2）用Romberg积分法计算积分的近似值，自己设置不同精度要求，对结果进行比较分析；与（1）的结果进行比较分析，谈谈你的体会。（3）记，在上面的计算中只取4位有效数字或7位有效数字，计算结果有什么不同。（4）上面计算精度可达8-20位有效数字吗？若可以请说明实现过程，并举例。一摘要：在matlab环境下熟悉的运用计算机编程

2、语言并结合Newton-Cotes和Romberg的理论基础对函数求积分，在运行完程序后以及对运行结果做出各方面的分析和比较。二实验设计目地用熟悉的计算机语言编程上机完成Newton-Cotes公式计算积分和Romberg积分法计算积分。三Newton-Cotes和Romberg的理论基础牛顿柯斯特公式：设将积分区间a,b划分为n等份，步长h=（b-a）/n,选取等距节点=a+kh构造出的插值型求积分公式I，称为牛顿柯斯特（Newton-Cotes）公式，式中称为柯特斯系数，引进变换x=a+th,则有龙贝格求积公式：（1）.梯形法的递推化 设将区间a,b分为n等份，共有n+1个分点，如果将求积

3、区间再二分一次，则分点增至2n+1个，我们将二分前后两人积分值联系起来加以考察，注意到每个子区间经过二分只增加了一个分点，用复合梯形公式求得该子区间上的积分值为，注意，这里h =(b-a)/n代表二分前的步长，将每个子区间上的积分值相加得，从而利用上式可导出下列递推分式：（2）.外推技巧 从梯形公式出发，将区间a,b逐次二分可提高求积公式精度，当a,b分为n等份I 若记，当区间a,b分为2n等份时，则有，并且有I.(3).龙贝格算法 将上述外推技巧得到的公式重新引入记号 等从而可将上述公式定成统一形式经过m(m=1.2.3)次加速后，余项便取下列形式：这方法通常称为理查森外推加速方法。设以表示

4、二分次后求得的梯形值，且以表示序列的m次加速值，则依递推公式可得公式也称为龙贝格算法，计算过程如下：.取k=0,h=b-a,求令（k记区间a,b二分次数）.求梯形值，即按递推公式计算.求加速值，按公式a逐个求出T表的第k行其余各元素T（j=1,2,k）.iv.若|TT|<(预先设定的精度)，则终止计算，并取TI；否则令k+1k转（2）继续计算.kh0b-a1234 T 表表指出了计算过程，第2列h=给出了子区间长度，i表示第i步外推.可以证明，如果f(x)充分光滑，那么T表每一列的元素及对角线元素均收敛到所要求的积分值I.四程序代码及运行结果牛顿柯斯特求积分：function y=mul

5、NewtonCotes(a,b,m,n)fun=(x)exp(-x.2);xk=linspace(a,b,m+1)for i=1:m s(i)=NewtonCotes(fun,xk(i),xk(i+1),n);endy=sum(s);function y,Ck,Ak=NewtonCotes(fun,a,b,n)%fun=(x)exp(-x.2)if nargin=1 mm,nn=size(fun); if mm>=8 error end xk=fun(1,:); fk=fun(2,:); a=min(xk); b=max(xk); n=mm-1;endif nargin=4 xk=lin

6、space(a,b,n+1); if isa(fun,'function_handle') fx=fun(xk); end end Ck=cotescoeff(n); Ak=(b-a)*Ck; y=Ak*fx'function Ck=cotescoeff(n) for i=1:n+1 k=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:n f=f.*(t-i);en

7、d运行结果>> y=mulNewtonCotes(0,1,2,1)xk = 0 0.5000 1.0000y = 0.7314>> 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

8、)/2;k=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);while abs(TT(k,1)-TT(k-1,1)./TT(k,1)>=e k=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:k TT(j1,k-j1+1)=4(j1-1).*TT(j1-1,k-j1+2)./(4(j1-1)-1)-TT(j1-1,k-j1

9、+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 = 0.7472>> vpa(z,7) ans = .7471804 >> z=romberg(0,1,10e-5)TT = 0.6839 0.7314 0.7430 0.7472 0 0TT = 0.6839 0.7314 0.7430 0.7459 0.7472 0.7469 0 0 0.7468 0 0 0z = 0.7468>> vpa(z,7) a

10、ns = .7468240 >> z=romberg(0,1,10e-10)TT = 0.6839 0.7314 0.7430 0.7472 0 0TT = 0.6839 0.7314 0.7430 0.7459 0.7472 0.7469 0 0 0.7468 0 0 0TT = 0.6839 0.7314 0.7430 0.7459 0.7466 0.7472 0.7469 0.7468 0 0 0.7468 0.7468 0 0 0 0.7468 0 0 0 0TT = 0.6839 0.7314 0.7430 0.7459 0.7466 0.7468 0.7472 0.74

11、69 0.7468 0.7468 0 0 0.7468 0.7468 0.7468 0 0 0 0.7468 0.7468 0 0 0 0 0.7468 0 0 0 0 0z = 0.7468>> vpa(z,7) ans = .7468241五.结果分析根据牛顿柯斯特求积分的原理将所求函数的积分区间分为多个小的积分区间，先求出每个小积分区间上的函数值，然后再将每个小区间上的求积结果加起来就是我们所要求的总函数的积分值，当函数区间所分的小区间的个数越多的时候总函数所计算出来的结果就精确，其原因就是所分的区间数越多，计算时每个小区所带来的误差就越小，当将总的积分值加起来的时候所带来的总的误差也就越小，所以最后的结果的精度越高，而用龙贝格求积和牛顿柯斯特也一样是要将总区间分为很多小的相等的区间，只是他们所用的计算原理不一样，当用此方法求积分时，所设的误差限越小，所求得的积分值就越是精确。六.参考文献1孙祥，徐流美，吴清。Matlabl 7.0基础教程【M】。北京：清华大学出版社，2005.2薛毅。