计算方法-数值积分b

1、5.3 复合求积复合求积 /* Composite Quadrature */高次插值有高次插值有Runge 现象现象，故采用分段低次插值，故采用分段低次插值 分段低次合成的分段低次合成的 Newton-Cotes 复合复合求积公式。求积公式。 复合梯形公式复合梯形公式在每个在每个 上用梯形公式：上用梯形公式：,1kkxx mkxfxfxxdxxfkkkxkxkk,., 1,)()(2)(11111)()(2)(2mkkbfxfafhbamkkkxfxfhdxxf11)()(2)(= Tm将将 a, b区间区间m等分，等分，步长步长h=(b-a)/m，分点分点xk=a+kh, k=0,1, m

2、。用。用低阶低阶牛顿牛顿柯特斯公式求柯特斯公式求子区间子区间xk, xk+1上的积分值，再累加得到积分的近似值。上的积分值，再累加得到积分的近似值。思思路路复合公式的余项复合公式的余项？5.3 Composite Quadrature 复化复化 Simpson 公式公式)()(4)(6)(1211 kkkxxxfxfxfhdxxfkkkx21 kx1 kx44444 )()(2)(4)(6)(1020121mkmkkkbabfxfxfafhdxxf= Sm327)(7)()(12)(32)(90)(14/32/14/11kkkkknxnxxfxfxfxfxfhxf 复化复化 Cotes 公式公

3、式327)()(32)(12)()(14904/32/111104/1bfxfxfxfxff(a)hkkmkmkkkCm 若若f(x)在积分区间在积分区间a, b上分别具有二阶、四阶和六上分别具有二阶、四阶和六阶阶连续导数连续导数，则复化积分公式的，则复化积分公式的余项余项分别是分别是 定理定理)(fhabCdx)x(f)(fhabSdx)x(f)(fhabTdx)x(f)(mba)(mbamba 66442)2(9452)2(18012其中，其中，a,b，且当，且当h充分小时，又有充分小时，又有)2(9452)2(180112155642)a(f)b(fhCdx)x(f)a(f)b(fhSd

4、x)x(f)a(f)b(fhTdx)x(f)()(mbambamba 5.3 Composite Quadrature证证 这里考察这里考察复化梯形公式复化梯形公式的余项的余项公式，其余由同学们自己完成。公式，其余由同学们自己完成。因此因此 )()()(112)()()(121)(212213mmmbafffmhabfffhTdxxf )(12)(2fhabTdxxfmba )(122fhab 中值定理中值定理)()(121121 1211020afbfdx)x(fh)(flimhTdx)x(flimbamkkhmbah 另有另有)x(f 由于由于 在在a, b上连续，故每个上连续，故每个小区

5、间小区间上的积分使用梯形公上的积分使用梯形公式时，有误差为式时，有误差为 kkkkxxfh,),(12113 )()(121)(2afbfhTdxxfmba5.3 Composite Quadrature 收敛速度与误差估计：收敛速度与误差估计： 若一个积分公式的误差满足若一个积分公式的误差满足 且且C 0，则，则称该公式是称该公式是 p 阶收敛阶收敛的。的。 ChfRphlim0)(,)(,)(642hOChOShOTnnn例：例：计算计算dxx 10142 解：解： )1()(2)0(161718fxffTkk8kxk 其中其中= 3.138988494)1(24)0(2413102/14

6、fxfxffSkkkk)()(34kxk其中其中= 3.141592502运算量基本运算量基本相同相同定义定义5.3 Composite QuadratureQ: 给定精度给定精度 ，如何取，如何取 m ?例如：要求例如：要求 ，如何判断，如何判断 m = ？|mTI)()(122afbfhfR上例上例中若要求中若要求 ，则，则610|mTI622106| )0() 1 (|12| |hffhfRm00244949.0 h即：取即：取 m = 409通常采取将区间通常采取将区间不断二分不断二分的方法，即取的方法，即取 m = 2k2/)()(,0hbfafTabh,a,binput 01TT

7、2/)()()(, 2/111hbfxfafThhmkk10TT 1 Toutput需计算导数需计算导数不实用不实用重复计算函重复计算函数值不经济数值不经济5.3 Composite Quadrature5.4 龙贝格龙贝格积分积分 /* Romberg Integration */ 梯形法的递推化梯形法的递推化11101)()(22)()(2mkkmkkkmbfxff(a)hxfxfhT积分区间积分区间a, b m等分的复化梯形公式是等分的复化梯形公式是 如果把区间如果把区间2m等分，即在原来的小等分，即在原来的小区间区间xk, xk+1上增加分点上增加分点xk+1/2=(xk +xk+1)

8、/2，变为两个小区间。于是有，变为两个小区间。于是有 )()(2)(4)(12/11kkkkxkxxfxfxfhdxxf102/11011012/12)(2)()(4)()(2)(4mkkmkkkmkkkkmxfhxfxfhxfxfxfhTmabhxfhTTmkkmm ,102/12)(221把把(T(T2m2m- T- Tm m)/3)/3作为误差的作为误差的修正值修正值加到加到T T2m2m上去，得到上去，得到 龙贝格算法龙贝格算法5.4 Romberg Integration 322mmmTTTI )4(3132222mmmmmmTTTTTT 上式的精度上式的精度完全有可能完全有可能比比

9、T T2m2m好。好。412 mmTITI考察考察例：例：计算，检验上述论断。计算，检验上述论断。 dxx 10142 mhTmTmSm113.13.133333320.53.13117653.13333333.141568640.253.13898853.14156863.141592580.1253.14094163.1415925306253.14142993.14159265S Sm m T T2m2m)(mmTT2431= =?T2m2m1012/1)()(2)(4mkkkkxfxfxfh1011012/12)()(6)()(2)(3)4(31mkkkmkk

10、kkmmxfxfhxfxfxfhTT证明证明)4(312mmmTTS 1012/1)()(4)(6mkkkkxfxfxfh101)()(2mkkkxfxfhTm m= = Sm m同理，考察同理，考察1612 mmSISI1522mmmSSSI )6(mmmmmmSSSSSS2222115115= = Cm m因此还有因此还有)46(6312mmmCCR Romberg 公式公式)(6hOCm)(4hOSm6412mmCICI5.4 Romberg Integration mmmSTT1442mmmCSS144222mmmRCC144323 Romberg 算法：算法： R4 T1 T8 T4

11、 T2 S1 S2 C1 R1 C2 S4 mabh)x(fhTTmk/kmm ,2211021211 T16 14 R2 13 C4 12 S8 梯梯形形递递推推化化公公式式Romberg公式公式m2mRR 5.4 Romberg Integration 5.5 高斯型高斯型积分积分 /* Gaussian Quadrature */bankkkxfAdxxf0)()(用用n+1个节点构造具有个节点构造具有2n+1次代数精度的求积公式次代数精度的求积公式将节点将节点 x0 xn 以及系数以及系数 A0 An 都作为待定系数。都作为待定系数。令令 f (x) = 1, x, x2, , x2n

12、+1 代入可求解，得到的公式代入可求解，得到的公式具有具有2n+1 次代数精度。这样的节点称为次代数精度。这样的节点称为Gauss 点点，公式称为公式称为Gauss 型求积公式型求积公式。例：例：求求 的的 2 点点 求积求积 公式公式11)(dxxf有有 次代数精度。次代数精度。)x(fA)x(fAdx)x(f110011 解：解：限定求积节点限定求积节点x0=-1, x1=1，得到，得到插值型求积公式插值型求积公式)(f)(fdx)x(f1111 1用解非线性方程用解非线性方程组求高斯点的方组求高斯点的方法很困难！法很困难！如果设如果设 ，我们对式中的系数，我们对式中的系数A0, A1和节

13、点和节点x0, x1不事先加以限制，而是适当地选其值，可使所不事先加以限制，而是适当地选其值，可使所得的公式有得的公式有 3 次代数精度。次代数精度。 1-11100)()()(xfAxfAdxxf代入代入 f (x) = 1, x, x2, x3,要求准确成立，得到要求准确成立，得到 5.5 Gaussian Quadrature 03202311300211200110010 xAxA/xAxAxAxAAA31322020102021/)(xxAAxx从而有从而有)(f)(fdx)x(f333311 称之为称之为2点点Gauss公式公式，具有，具有3次代数精度。次代数精度。 bankkkx

14、fAdxxf0)()(构造具有构造具有2n+1次代数精度的次代数精度的n+1点点Gauss公式公式?133331010AAxx/如果式中如果式中x0 xn 为为 Gauss 点点, 则公式则公式 至少有至少有 2n+1 次代数精度。次代数精度。bankkkxfAdxxf0)()(对任意次数对任意次数不大于不大于n 的多项式的多项式 P(x)， P(x) w(x)的次数的次数不不大于大于2n+1，则代入公式应则代入公式应精确成立精确成立：nkkkkbaxwxPAdxxwxP0)()()()(0= 0 x0 xn 为为 Gauss 点的充要条件是点的充要条件是 与与任意次数不大于任意次数不大于n

15、的多项式的多项式 P(x) 正交正交，即成立，即成立 nkkxxxw0)()(定理定理利用区间利用区间a, b上的上的n+1次正交多项式次正交多项式确定确定Gauss点点；然后；然后利用利用代数精度代数精度确定确定求积系数求积系数。作一个作一个n+1次的多项式次的多项式 )()()()(100nnkkxxxxxxxxxw0)()(badxxwxP求求 Gauss 点点 求求w(x)5.5 Gaussian Quadrature 高斯高斯勒让德求积公式勒让德求积公式/*Gauss-Legendre 公式公式Legendre 多项式族：多项式族： 定义在定义在 1, 1上上nnnnnxdxdnxP

16、)1(!21)(2)33035(81)()35(21)()13(21)()(1)(244232210 xxxPxxxPxxPxxPxP递推公式递推公式 )x(nP)x(xP)n()x(P)n(nnn11121 Legendre多项式多项式Pn(x)对于任意对于任意n-1次的多项式在次的多项式在-1, 1上上正交正交。 定理定理5.5 Gaussian QuadraturennnnnxdxdnxP)1(!21)(2 Legendre多项式多项式Pn(x)对于任意对于任意n-1次的多项式在次的多项式在-1, 1上上正交正交。 定理定理证明证明： 令令n)x()x(g12 有有110 01 n,.,

17、k,)x(gx)k(11)(11)()(!21)()(dxxgxQndxxQxPnnn11)1(11)1()()()()(!21dxxgxQxgxQnnnn11)1()()(!21dxxgxQnnn 11)2(11)2()()()()(!21dxxgxQxgxQnnnn11)()()(!2)1(dxxgxQnnnn00= 0k次勒让德多项次勒让德多项式式的的k个零点个零点就就是是k k个高斯点个高斯点，用来构造用来构造k k点高点高斯斯勒让德求积勒让德求积公式公式。 详见教详见教材材P.140P.140表表5.15.1若若Q(x)是次数小于是次数小于n的多的多项式，则恒有项式，则恒有Q(n)(

18、x)=05.5 Gaussian Quadrature高斯高斯勒让德求积公勒让德求积公式式nnkkkGxfAdxxf110)()(n=0: 一点公式一点公式011)0(2)(Gfdxxf中矩形公式中矩形公式badxxf)()(knkktababfAab2220?9460831. 02121sin21sin3010kkkkttAdxxx9460831. 0sin21sin21sin301110kkkkxxAdxxxdxxx解解：作变量替换：作变量替换 x=1/2+t/2本题可以不作变量替换本题可以不作变量替换: : 10dxxxsin例例 用用四点四点高斯高斯勒让德公式计算勒让德公式计算(I=0.9460831)n=1: 两点公式两点公式n=3: 四点公式四点公式)33998104. 0()33998104. 0(65214515. 0)86113631. 0()86113631. 0(34785485. 03ffffG )