数值积分ppt课件 (《计算方法》).ppt

《数值积分ppt课件 (《计算方法》).ppt》由会员分享，可在线阅读，更多相关《数值积分ppt课件 (《计算方法》).ppt（81页珍藏版）》请在三一办公上搜索。

1、2022/12/25,1,第7章 数值积分,1 插值型求积公式 2 复化求积公式3 龙贝格(Romberg)求积方法,2022/12/25,2,1 插值型求积公式,在一元函数的积分学中,我们已经熟知,若函数f(x)在区间a, b上连续且其原函数为F(x) ,则可用牛顿莱布尼兹公式,(71),来求定积分。,2022/12/25,3,公式(71)虽然在理论上或在解决实际问题中都起了很大的作用,但它并不能完全解决定积分的计算问题。因为定积分的计算常常会碰到以下三种情况: (1)被积函数f(x)的原函数F(x)不易找到。许多很简单的函数,例如,其原函数都不能用初等函数表示成有限形式。,2022/12/

2、25,4,(2)被积函数f(x)没有具体的解析表达式。其函数关系由表格或图形表示,无法求出原函数。,其被积函数 的原函数就比较复杂,从数值计算角度来看,计算量太大。,(3)尽管f(x)的原函数能表示成有限形式但其表达式相当复杂。例如定积分,2022/12/25,5,图 7.1,如图7.1,若用左矩形近似地代替曲边梯形,则得到左矩形公式,(72),2022/12/25,6,同样可得到右矩形公式：,(73),2022/12/25,7,图 7.2,如图7.2,若用梯形的面积近似地代替曲边梯形的面积,则得到计算定积分的梯形公式,(74),2022/12/25,8,如图7.3,若用抛物线代替曲线f(x)

3、,则可得到抛物线公式(或辛普生公式),(75),图7.3,2022/12/25,9,此外,众所周知的梯形公式: I(f)(b-a)f(a)+f(b)/2和 Simpson公式: I(f)(b-a)f(a)+4f(a+b)/2)+f(b)/6则分别可以看作用 a, b, c=(a+b)/2, 三点高度的加权平均值 f(a)+f(b)/2 和 f(a)+4f(c)+f(b)/6作为平均高度f()的近似值.,2022/12/25,10,更一般地,取区间a,b内n+1个点 xi,(i=0,1, 2,n) 处的高度f(xi) (i=0,1,n)通过加权平均的方法近似地得出平均高度f(),这类求积方法称为

4、机械求积:,2022/12/25,11,或写成:,数值积分公式,求积系数,求积节点,(1),2022/12/25,12,记,称(2)为数值求积公式,(3)为求积公式余项(误差).,构造或确定一个求积公式，要讨论解决的问题有,(i) 确定求积系数Ak和求积节点xk ；(ii)求积公式的误差估计和收敛性,为了构造形如式(2)的求积公式,需要提供一种判定求积方法精度高低准则,2022/12/25,13,求积公式的代数精度定义1 称求积公式(2)具有m次代数精度,如果它满足如下两个条件: (i)对所有次数 m次的多项式 ,有 (ii)存在m+1次多项式 ,使得,定义1中的条件(i),(ii)等价于:,

5、2022/12/25,14,插值型求积公式,在积分区间a,b 上取n+1个节点xi，i=0,1,2,n,作f(x)的n次代数插值多项式（拉格朗日插值公式）:则有 为插值余项于是有,2022/12/25,15,取称(4)式为插值型求积公式,其中求积系数Ak由(5) 式确定.,(4),(5),2022/12/25,16,2022/12/25,17,推论1 求积系数满足:,误 差,定理1 形如 的求积公式至少有 n 次代数精度 该公式为插值型（即： ）,2022/12/25,18,现用第六章介绍的插值多项式Pn(x)来代替被积函数f(x),即有,1.1 牛顿柯特斯公式(NewtonCotes),取节

6、点为等距,即 a=x0 x1xn=b,建立数值积分公式最基本的思想是选取一个既简单又有足够精度的函数(x), 用(x)代替被积函数f(x),于是有,2022/12/25,19,利用拉格朗日插值多项式,(76),其中,(77),2022/12/25,20,这里yi=f(xi),对式(76)两边积分得,2022/12/25,21,为牛顿柯特斯(Newton-Cotes)求积公式,Rn(f)为牛顿柯特斯求积公式的余项。,我们称,2022/12/25,22,令 x=x0+sh , 0sn dx=hds=(b-a)/nds,(711),2022/12/25,23,Newton-Cotes公式的误差为:,

7、与x有关,注意:由(7-11)式确定的Cotes系数只与i和n有关,与f(x)和积分区间a,b无关,且满足:,(7-9),2022/12/25,24,称Ci(n)为柯特斯求积系数。很显然,当n=1时,可算得,此时式(710)为,(712),这是梯形公式。,2022/12/25,25,当n=2时,可得,于是,(713),这是抛物线(Simpson)公式。,2022/12/25,26,当n=3时,代入(710)式得到求积公式,2022/12/25,27,类似地可分别求出n=4,5,时的柯特斯系数,从而建立相应的求积公式。具体结果见表71。,从表中可以看出,当n7时,柯特斯系数为正;从n8开始,柯特

8、斯系数有正有负。因此,当n8时,误差有可能传播扩大,牛顿 柯特斯求积公式不宜采用。,柯特斯系数Ci (n) 仅与n和i有关,与被积函数f(x)无关,且满足,(715),事实上,式(710)对f(x)=1是准确成立的。,2022/12/25,28,表 71,2022/12/25,29,定理 当阶数n为偶数时, Newton-Cotes公式(8)至少具有n+1次代数精度.,证明 只需验证当n为偶数时,Newton-Cotes公式对f(x)=xn+1的余项为零.,由于f(x)=xn+1,所以f(n+1)(x)=(n+1)! .由式(7-9)得,引进变换t=u+n/2,因为n为偶数,故n/2为整数,于

9、是有,据此可断定R(f)=0,因为上述被积函数是个奇函数.,2022/12/25,30,Newton-Cotes公式的数值稳定性,现在讨论舍入误差对计算结果产生的影响.设用公式 近似计算积分 时,其中计算函数值f(xj)有误差j (j=0,1,2,n).设计算Cj(n)没有误差,中间计算过程中的舍入误差也不考虑,则在式(10 )的计算中,由j引起的误差为,(10),2022/12/25,31,如果Cj(n)都是正数,并设,故en是有界的,即由j引起的误差受到控制,不超过的(b-a)倍,保证了数值计算的稳定性.而当n7时,Cj(n)将出现负数,保证数值稳定性. 因此高阶公式不宜采用,有实用价值的

10、仅仅是几种低阶的求积公式.,将随n增大,因而不能,则有,2022/12/25,32,解 利用梯形公式,利用抛物线公式,例1 试分别用梯形公式和抛物线公式计算积分:,原积分的准确值,2022/12/25,33,现对牛顿柯特斯求积公式所产生的误差作一个分析。由式(79),牛顿柯特斯求积公式的余项为,1.2 误差估计,易知,牛顿柯特斯求积公式(710)对任何不高于n次的多项式是准确成立的。这是因为 f(n+1)()0 故 Rn(f)0,2022/12/25,34,牛顿柯特斯求积公式的代数精确度至少为n。通常在基点个数相等的情况下,代数精确度愈高,求积公式就愈精确。,一般说来,若某个求积公式对于次数不

11、高于m的多项式都准确成立(即Rn(f)0),而对于某一次数为m+1的多项式并不准确成立,则称这一求积公式的代数精确度为m。,定理1 (梯形公式的误差)设f(x)在区间a, b上具有连续的二阶导数,则梯形求积公式的误差为,2022/12/25,35,由于 1(x)=(x-a)(x-b),证 由式(79)知,梯形公式的余项为,1(x)在区间(a, b)内不变号,f()是x的函数且在a,b上连续,故根据积分第二中值定理知,存在某一(a, b)使,2022/12/25,36,定理2 (抛物线公式的误差)设f(x)在a, b上有连续的四阶导数,则抛物线公式的误差为,(717),2022/12/25,37

12、,n = 1:,Trapezoidal Rule,/* 令 x = a+th, h = ba, 用中值定理 */,代数精度 = 1,n = 2:,Simpsons Rule,代数精度 = 3,n = 4: Cotes Rule, 代数精度 = 5,2022/12/25,38,复合求积公式,高次插值有Runge 现象,高阶Newton-Cotes公式会出现数值不稳定,低阶Newton-Cotes公式有时又不能满足精度要求.解决这个矛盾的办法是将积分区间a,b分成若干小区间,在每个小区间上用低阶求积公式计算,然后将它们加起来,这就是复合求积方法.,2022/12/25,39,2 复合求积公式,2.

13、1 复合梯形公式 对于定积分(71),将积分区间a, b分成n个相等的子区间xi,x i+1,这里步长,在每一个子区间 xi,x i+1 上使用梯形公式,则,2022/12/25,40,复化梯形公式积分法,2022/12/25,41,相加后得,(718),(719),若f(x)在a,b上连续,由连续函数的介值定理,存在某一(a,b)使得,2022/12/25,42,因而,于是得到复合梯形公式,(721),其余项为,2022/12/25,43,例2 若用复合梯形公式计算积分,则 当0 x1时,有,因为,又,解 由余项(721)式,问积分区间要等分多少才能保证有五位有效数字?,2022/12/25

14、,44,由于原积分的准确值具有一位整数,因此要使近似积分值有五位有效数字,只需取n满足,两边取对数得,整理后得到,取 n=68.,2022/12/25,45,类似复合梯形公式的做法,把区间a,b分成n个相等的子区间x2i,x2i+2(i=0,1,n-1),设每个子区间上的中点为x2i+1(i=0,1,n-1),且,(722),2.2 复合抛物线公式,在每一个子区间 x2i ,x2i+2 上利用抛物线公式得,2022/12/25,46,复化Simpson公式积分法,2022/12/25,47,相加后得,(723),2022/12/25,48,若f(4)(x)在a, b上连续,则,从而得到复合抛物

15、线公式,(724),其余项为,(725),2022/12/25,49,图 7.4 复合抛物线公式框图,2022/12/25,50,的数据表,分别用复合梯形公式和复合抛物线公式计算,例3 已知函数,x f (x)0 11/8 0.99739782/8 0.98961583/8 0.97672674/8 0.95885105/8 0.93615566/8 0.90885167/8 0.87719251 0.8414709,2022/12/25,51,解 用复合梯形公式,这里,2022/12/25,52,用复合抛物线公式可得,比较上面两个结果T8和S4,它们都需要提供9个点上的函数值工作量基本相同,

16、然而精度却差别很大.,同积分的准确值I(f)=0.9460831比较,复化梯形法的结果T8=0.9456909只有两位有效数字, 而复化,Simpson法的结果S4=0.9460832却有六位有效数字.,2022/12/25,53, 复化梯形公式：,在每个 上用梯形公式：,= Tn,/*中值定理*/,2022/12/25,54, 复化 Simpson 公式：,2022/12/25,55,2.3 变步长公式,前面介绍的复合梯形公式和复合抛物线公式的步长都是预先确定的。它的主要缺点是事先很难估计出n的大小(或步长h的大小),使结果达到预先给定的精度。在实际计算中,我们常常借助于计算机来完成积分步长

17、h的自动选择,即采用变步长求积公式。具体地讲,就是将步长逐次折半,反复利用复合求积公式,直到满足精度要求为止。,2022/12/25,56,下面介绍变步长复合抛物线公式(变步长复合梯形公式留给读者作为练习)。,(726),其中,再把每个子区间分成两半,用,逐次将区间a,b分成2,4,2m等分,并按复合抛物线公式逐次计算积分得到S1,S2,Sm,而,2022/12/25,57,作步长,按复合抛物线公式计算出积分的近似值S2m。对于相邻两次的积分近似值Sm、S2m,考察,设预先给定的精度为,若 d 则以S2m作为所要求的积分近似值,否则继续将区间分半,利用复合抛物线公式求积分,直到满足预给的精度为

18、止。,2022/12/25,58,图 7.5 变步长复合抛物线公式,2022/12/25,59,图7.5 变步长复合抛物线公式,2022/12/25,60,3 龙贝格(Romberg)积分方法,我们已经知道,当被积函数f(x)在区间a,b上连续时,要使得复合梯形公式或复合抛物线公式比较精确地代替定积分,可将分点(即基点)加密,也就是将区间a,b细分,然后利用复合梯形公式或复合抛物线公式求积。,2022/12/25,61,若用Tm表示把a,b作m等分并按复合梯形公式求积的结果,将每一小段再对分,令新的小段的长h=h/2,则T2m与Tm之间有如下关系:,(528),其中,2022/12/25,62

19、,另外,若用Sm表示把a,b分成m(偶数)个小段按复合抛物线公式计算的结果,那么只要把Sm中的m改为2m,h改为h就有,从Tm的定义可得到关系式,(529),2022/12/25,63,我们再举一个计算上半单位圆面积的例子(它的准确面积为/2)。现用内接正多边形的逼近方法来计算。,如图5.6,图(a) 、 图(b)是用同样的内接正多边形计算上半单位圆的面积。图(a)是用梯形方法计算其面积,图(b)是用三角形方法计算其面积。,2022/12/25,64,图 5.6,2022/12/25,65,设正多边形边数为n=2k,则由图(b)利用三角形公式算得面积为,同理,2022/12/25,66,如果组

20、合一下,就会得到更精确的结果,即,同理,2022/12/25,67,再以类似方法组合得,这样继续下去,其值越来越接近上半单位圆面积/2。这种方法可以用到计算定积分,2022/12/25,68,为了推广公式(529)和上述计算上半单位圆面积的组合方法,我们引进龙贝格求积算法。 龙贝格求积算法本来是利用所谓外推法构造出的一种计算积分的方法。为了避免从外推引入而带来理论上的麻烦,我们将直接从构造一个T数表开始。 首先将a,b依次作20,21,22,等分,记,2022/12/25,69,按复合梯形公式(520)算得的值相应地记为T(k)0(k=0,1,2,);把按式(529)算得的S2m依次记为T(k

21、)1(k=0,1,2,崐),而这每一个S2m又理解为由T2m与Tm的线性组合得到的改进值,即,我们可按照类似的方法继续进行改进,也即由S2m与Sm的线性组合得到改进值,依次记为T(k)2(k=0,1,2,),即,2022/12/25,70,这样就可构造出一个数表,(5-30),2022/12/25,71,其中除第0列(即最左一列)的T(k)0是按复合梯形公式计算外,其余各列都按下述规则(对m),(531),递推地计算出来。箭头表示计算流程。其计算步骤为: (1)将区间a,b等分为20,用梯形公式计算T(0)0,即,2022/12/25,72,(2)将区间a,b等分为21,用梯形公式算出T(1)

22、0,即,再由T(0)0,T(1)0根据公式(531)算出T(0)1,即,若 T(0)1-T(0)0, (为预给的精度) 则停止计算;否则继续往下计算;,2022/12/25,73,(3)依次分别算出T(2)0,T(1)1,T(0)2,这一行地往下推算,每一行算完,就得验证T(0)m(m=1,2,)是否满足预给的精度,即若,则停止计算;否则继续进行下一行。为了便于在计算机上实现,可运用下列公式编制程序:,2022/12/25,74,2022/12/25,75,图 5.7,2022/12/25,76,图 5.7,2022/12/25,77,例 4 计算积分,精确到10-4。 解,2022/12/25,78,2022/12/25,79,2022/12/25,80,2022/12/25,81,于是,由于,实际上,