资源描述
第五章第五章 数值积分数值积分 5.0 引言引言 5.1 机械求积公式机械求积公式 5.2 Newton-Cotes公式公式5.3 变步长求积公式及其加速收敛技巧变步长求积公式及其加速收敛技巧5.4 Gauss公式公式5.5 小结小结 15.0 引 言 引 言 1.定积分的计算可用著名的牛顿-莱布尼兹公式来计算:()()()baf x dxF bF a=其中 F(x)是 f(x)的原函数之一，可用不定积分求得。然而在实际问题中，往往碰到以下问题：(a)被积函数 f(x)是用函数表格提供的；(b)被积函数表达式极为复杂，求不出原函数，或求出原函数后，由于形式复杂不利于计算；(c)大量函数的原函数不容易或根本无法求出，例如 210 xedx，概率积分 10sinxdxx,正弦型积分 2222204()1sinIrxH xdrxr =回路磁场强度公式 等根本无法用初等函数来表示其原函数，因而也就无法精确计算其定积分，只能运用数值积分。2 所谓数值积分就是求积分近似值的方法。而数值积分只需计算()f x 在节点(1,2,ixin=?)上的值，计算方便且适合于在计算机上机械地实现。25.1 机械求积公式 1 数值积分的基本思想数值积分的基本思想 区间a,b上的定积分，就是在区间a,b内取 n+1 个点()baf x dx01,nx xx?，利用被积函数 f(x)在这 n+1 个点的函数值的某一种线性组合来近似作为待求定积分的值，即 0()()nbkkakf x dxA f x=右端公式称为左边定积分的某个数值积分公式。其中,xk称为积分节点，Ak称为求积系数。因此，一个数值积分公式关键在于积分节点xk的选取和积分系数Ak的决定，其中Ak与被积函数f(x)无关。称为机械求积公式。1.1 简单算例说明简单算例说明 例例 1 求积分 10()xxf x dx 此积分的几何意义相当于如下图所示的曲边梯形的面积。解解：（1）用 f(x)的零次多项式 00()()yLxf x=来近似代替,于是，()f x11000()()xxxxf x dxf xdx 001()(f xxx=（为左矩公式）3 推广推广:（为右矩公式）11001110()()()()xxxxf x dxf xdxf xxx=1100010110()()2()(2xxxxxxf x dxfdxxx)fxx+=+=（为中矩公式）（2）用 f(x)的一次多项式 01100110()()()xxxxLx1f xf xxxxx=+=+来近似代替,于是，()f x11001()()xxxxf x dxLx dx 10010101101001()()1()()()2xxxxxxf xf xdxxxxxxxf xf x=+=+=+=+（为梯形公式）4（3）用 f(x)的二次插值多项式,其中01xxx 011200010101101()()()()()()()()()()()()()()()()xxxxxxxxLxf xxxxxxxxxxxxxf xxxxxf x=+=+来近似代替,于是，()f x11002()()xxxxf x dxL x dx 特别地：当011(2xxx=+=+)时，有 10100101()()()4()()62xxxxxxf x dxf xff x+5 （为 Simpson 公式）2 代数精确度代数精确度 定义：若积分的数值积分公式对于任意()baf x dx 0()()nbkkakf xdxA f x=一个次数不高于 m 次的多项式都精确成立，且存在一个 m+1 次多项式使之不精确成立，则称该数值积分公式的代数精确度为 m。对于代数精确度为 m 的求积公式，若 f(x)是不超过 m 次的代数多项式，求积公式是精确成立的。2.1 算例算例 例 1:有积分公式:11112012()()()()f x dxfff+求该积分公式的代数精确度。这个求积公式的几何意义是曲边梯形的面积近似地用两个梯形面积来代替。-1 0 1-1 0 1 X Y X Y f(-1)f(1)f(-1)f(1)f(0)f(0)6解：(1)取 f(x)=1，定积分1112 dx=,而数值积分11 2 122+=+=，两端相等；(2)取 f(x)=x，定积分110 xdx=,而数值积分112 0 102()+=+=，两端相等；(3)取，定积分2()f xx=12123x dx=,而数值积分22112 0 12()+=+=1，两端不相等；只要取 f(x)=1，f(x)=x 验证了上述求积公式精确成立，就意味着对于任意一个一次多项式，求积公式都是精确成立的；而取时求积公式不精确成立，也就是存在一个二次多项式使求积公式不精确成立；故该求积公式的代数精确度为 1。2()f xx=例 2：在如下求积公式中，求积分节点12,x x和相应的求积系数使其12,A A代数精确度尽可能高。111221()()()f x dxA f xA f x=+=+解：(1)f(x)=1,而数值积分为11 1 2dx=12AA+；得到方程；122AA+=7(2)f(x)=x，而数值积分为11 0 x dx=1 12 2AxA x+；得到方程；1 12 20AxA x+=(3)，2()f xx=1212 3x dx=，而数值积分为221 12 2AxAx+；得到方程221 12223AxA x+=；(4)3()f xx=，131 0 x dx=，而数值积分为331 122AxA x+；得到方程；331 1220AxA x+=综合上述方程：121 122221 122331 1222 (1)0 (2)2 (3)30 (4)AAAxA xAxA xAxA x+=+=+=+=?解得:1211,33xx=121AA=。于是我们得到积分公式 1111()()()33f x dxff=+。再取，有4()f xx=1412 5x dx=，8而数值积分为4411933+=2，两式不相等，求积公式不精确成立了。所以，该积分公式的代数精确度为 3。95.2 Newton-Cotes公式 1 公式的推导公式的推导 Newton-Cotes 公式是由拉格朗日插值公式推导出来的数值积分公式。将区间a,b等分 n 等份,记bahn=,分点为kxakh=+,k=0,1,.,n,这 n+1 个节点上的函数值为(),0,1,kf xkn=?,从而区间a,b上的拉格朗日插值多项式为 0()()()nnkkL xf x l x=k 其中,()kl x为插值基本多项式,与函数 f(x)无关,k=0,1,.,n。()()bbnaaf x dxL x dx 0()()nbk kakf x l xdx=00()()()nbkkaknkkkl x dx f xA f x=由于插值结点是等距节点，故插值多项式可以进一步化简:因为 ，kxak=+hxath=+故 ，()kjxxkj h=()jxxtj h=10011011()()()()()()()()(nkkknkkkkkkx xx xx xx xl x)xxxxxxxx+=?(1)(1)(1)()(1)(1)(1)()t ttktktnk kkkkkkn+=+?0,(1)()!()!nn kjj ktjk n k=因()bkkaAl x dx=，作积分变量代换，xath=+badxhdtdtn=，当时，t=0；当 x=b 时，t=n；xa=故 00,()1(1)()!()!nnjj kn kkbatj dtAn k n k=()(),0,1,nkkA记ba Ckn=?，我们称为柯特斯（Cotes）系数，它不仅与函数 f(x)无关，而且与积 分区间a,b无关。()nkC例如：当 n=1 时 (梯形积分公式中的系数梯形积分公式中的系数)11(1)00(1)1(1)0!1!2Ctd=t=,01(1)10(1)1(0)1!0!2Ctd=t=;11当 n=2 时 (抛物线积分公式中的系数抛物线积分公式中的系数)22(2)00(1)1(1)(2)0!2!26Ctt=dt=,12(2)10(1)2(0)(2)1!1!23Ctt=dt=,02(2)20(1)1(0)(1)2!0!26Ctt=dt=;当 n=3 时 (3/8 积分公式中的系数积分公式中的系数)33(3)00(1)1(1)(2)(3)0!3!38Cttt=dt=，23(3)10(1)3(0)(2)(3)1!2!38Cttt=dt=，13(3)20(1)3(0)(1)(3)2!1!38Cttt=dt=，03(3)30(1)1(0)(1)(2)3!0!38Cttt=dt=；于是，由柯特斯（Cotes）系数公式出发，我们得到 n 阶 NewtonCotes公式：()nkC()0()()()nbnkkakf x dxbaCf x=。2 低阶公式及其余项低阶公式及其余项 12常用的常用的 NewtonCotes 公式公式 a)梯形公式 n=1 时，积分节点为，0 xa=1xb=，则数值积分公式为：()()()2babaf x dxf af b+Y X Y X f(a)f(b)f(a)f(b)f(x)f(x)a b a b O O 其几何意义是曲边梯形的面积 ()baf x dx近似地用梯形面积()()2baf af b+来代替。其余项31()()()(,)12baR ffa b=b)抛物线公式（辛浦生 Simpson 公式）n=2 时，积分节点为x0=a，12abx+=，x2=b；柯特斯系数为(2)(2)(2)02112,63CCC=；则数值积分公式为:13 ()()4()()62babaabf x dxf aff b+其几何意义是曲边梯形的面积 ()baf x dx近似地用由抛物线形成的曲边梯形面积来代替。其余项 5(4)2()()()(,)2880baR ffa b=c)柯特斯公式 n=4 时，积分节点为，0 xa=4xb=，,1,4kbaxakh hk=+=2,3；柯特斯系数为4404790()()CC=,44413232129090()()(),CCC=；则数值积分公式为：()baf x dx 012347()32()12()32()7()90b af xf xf xf xf x+其余项 7(6)48()()()(,)9454baR ffa b=综上所述，Newton-Cotes 数值积分公式具有如下特点：(1)建立在等距积分节点上，(2)是封闭型的，即两个端点 a,b 也是积分节点，14(3)是由拉格朗日插值公式推导而得到的。2.1 Cotes 系数的性质系数的性质()nkC引理引理:n 阶 NewtonCotes 公式()0()()()nbnkkakf x dxbaCf x=的代数精确度至少是 n。证明:如果()f x是一个次数不超过 n 次的多项式，则(1)()0nfx+其拉格朗日插值公式的插值余项为：(1)1()()()()()0(1)!nnnnR xf xL xfxn+=+故()()nf xL x=，这是对一切 x 均相等，精确成立。所以，()0()()()()nbbnnkaakkf x dxL x dxbaCf x=即，数值积分公式的值精确地等于定积分的值，故 n 阶 NewtonCotes 公式的代数精确度至少是 n。性质性质 1:归一公式：()01,1,2,3,nnkkCn=?15 证明：由于数值积分公式的代数精确度至少为 n，故对于，数值积分公式是精确成立的：()1f x=()1()bbaaf x dxdxba=，而()()00()()()nnnnkkkkkbaCf xbaC=由上述两式相等，得到：()01,1,2,3,nnkkCn=?性质 2:对称性：。()()nnknCC=k3 复合求积公式复合求积公式 随着 n 的增加可以减少积分误差，但高阶 N-C 公式又会造成数值不稳定，因而采用复合求积公式。16 3.1 复合梯形公式复合梯形公式 将区间等分 n 等份,bahn=,分点是xk=x0+kh,(k=0,1,.,n),其中0,naxbx=。在每个子区间上用梯形公式 1,kkx x+1 1()()()2kkxkkxhf x dxf xf x+则 1 1 0()()()2nakkbkhf x dxf xf x+=+17 011()2()2()()2nnhf xf xf xf x=+?1()2()()2nkkhf af xf=+b 此公式就是复合梯形公式复合梯形公式。其几何意义是曲边梯形面积近似地用许多小的细条梯形来代替（如图）从图中可以看出，n 越大，则 h 越小，实际面积与近似面积的差，即求积误差也就越小。这与分段插值相类似，所不同的是分段插值函数是不光滑的，而数值积分公式是对一个数的近似，不存在光滑和不光滑的问题。3.2 复合辛浦生公式复合辛浦生公式 将区间a,b等分 n 等份，n 为偶数（n=2n），bahn=，分点是,0122,max x xxb=?1kkxxh+=，简记()kkf xf=，则在每两个子区间(222,kkxx+0,1,2,1km=?)上用辛浦生公式:18222 2222212()46kkxkkkkkxxxf x dxfff+2 221243kkkhfff+=+2 则 1 2212 0()43makkkbkhf x dxfff+=+2()012322212424243mmhmfffffff=+?此公式就是复合辛浦生公式复合辛浦生公式。3.3 算例算例 分别利用梯形公式和 Simpson 公式计算积分：10sinxIdxx=步长 h=1/8。解解：设,0,1,8ixihi=?由复合梯形公式有：1017sin(0)2()2()(1)2 0.94569086xIdxxhff xf xf=+=?由复合 Simpson 公式有：19 10322120sin43 0.94608331kkkkxIdxxhfff+=+=2 积分的相对精确值为 10sin0.94608309xIdxx=。205.3 变步长求积公式及其加速收敛技巧 1.基本思想：基本思想：在计算机上自动选择步长的变步长的方法和加速收敛的技巧在计算机上自动选择步长的变步长的方法和加速收敛的技巧（Richardson 外推）外推）a)定义：定义：若一个积分公式的误差满足 0 limphR fCh=且 C 0，则称该公式是 p 阶收敛阶收敛的。b)复合梯形公式：复合梯形公式：在每个上用梯形公式梯形公式：1,kkxx 1 1 0()()()2nakkbkhf x dxf xf x+=+011()2()2()()2nnhf xf xf xf x=+?1()2()()2nkkhf af xf=+b 误差：21321122()()()1212()(),(,)12()()12nknkkkfhhR ffbanhba fa bhf bfa=c)复合复合 Simpson 公式：公式：1121()()4()()6kkxkkkxhf x dxf xf xf x+1211100()()4()2()()6nnbkkakkhf x dxf af xf xf b+=+误差：误差：4(4)()1802ba hR ff=d)综上所述综上所述:TO 24(),()nnhSO h1.1 算例：算例：21410 xdx+=例：计算 22(1)7811(0)2()(1)16kkTff xf=+=3.138988494(2)4oddeven1(0)4()2()(1)24kkSff xf xf=+=3.141592502 其中其中8kkx=。2 变步长梯形求积公式：变步长梯形求积公式：复合梯形求积时，通常采取将区间不断对分不断对分(一分为二一分为二)的方法，即取 n=2k,注意到区间再次对分时，有:()(2)baf af b+0 b)(aTIf x dx=101()222babaTafb2()()2()()22abaabIf x dxf afTf b+=+=+3 b3221212211(21)()()222iibaibaTf a=()()()2()2)2abai baIf x dxfTaf af b=+=+如此类推,可得 变步长的梯形公式为变步长的梯形公式为(1):递推公式 232111(21)()2221,2,.KKKkbaibaf aK=1b )(KKaTTIf x dx+=(2):步长的自动选取 注意到区间再次对分时注意到区间再次对分时:2111()()1224KKhRff bfaRf+114KKITIT 11()3KKKITTT 因此，若给定精度,则以递推公式计算积分近似值，直至1KKTT终止计算，并以当前值为近似值。KT 3 加速收敛技巧加速收敛技巧 3.1 基本思想基本思想 龙贝格积分龙贝格积分 Romberg Integration，主要采用主要采用 Richardson 外推。外推。变步长梯形公式中：注意到 114KKITIT 即：1441KKTTI 由此，可得梯形积分值外推一次的计算公式为：24141,2,3,.4 1KKKTTSK=一般地，外推二次的计算公式为：21242,3,.41KKKSSCK=外推三次的计算公式为：31343,4,5,.41KKKCCRK=外推 m 次的计算公式为：14,1,2,.41mKKKmddeKm mm=+.重复上述步骤即可获得得系列逼近值：3.2 算例算例 用龙贝格积分法计算定积分（用龙贝格积分法计算定积分（1e6=）：）：12111 25Idxx=+255.4 Gauss公式 1 引言引言 牛顿柯特斯型求积公式是封闭型的（区间a,b的两端点 a,b 均是求积节点）而且要求求积节点是等距的,受此限制，牛顿柯特斯型求积公式的代数精确度只能是(为奇数)或n n1n+(为偶数)。而如果对求积节点也适当的选取,即在求积公式中不仅nkA而且也加以选取,这就可以增加自由度，从而可提高求积公式的代数精确度。kx2高斯求积公式和高斯点 例:10011 1()()()f x dxA f xA f x+其中,固定在01,x x1,1,0,1A A可以适当选取,只有两个自由度,得到的是梯形公式,其代数精确度只有 1。如果对求积节点也进行适当选取,将有四个自由度,得到如下公式:01,x x 1 111()33f x dxff+这个积分公式的代数精确度为 3,这就是高斯型求积公式,上面的求积节点13称为高斯点。定义 1：高斯型求积公式和高斯点,kkA x)n(0,1,k=?的求积公式:对于含有个参数22n+0()()nbkkakf x dxA f x=,适当选取这22n+个参数,可以使得数值积分公式的代数精确度达到,21n+我们称这一类求积公式为高斯型求积公式高斯型求积公式,26称这类求积公式的积分节点为高斯点高斯点。定义 2:如果个求积节点的求积公式的代数精确度为,则这1n+2n+11n+个求积节点称为高斯点高斯点。3.高斯点的特征高斯点的特征 定理定理:设是个相异点,以这01,nx xx?1n+1n+个点为零点的次多项式为，则是高斯点的充要条件是对于任意不超过次的多项式,成立。1n+101()()()()nnxxxxxxx+=?01,nx xx?n()q x 1()()0bnaq xx dx+=证明:(必要性)设为高斯点,则对任意不超过次的多项式01,nx xx?2n+1()f x,均有 0()()nbkkakf x dxA f x=,则对于任意次数不超过次多项式,是次数不超过n()q x1()()nq xx+21n+次的多项式,且注意到,故1()0,(0,1,)nkxk+=?nk 11 0()()()()0nbnkknakq xx dxA q xx+=(充分性)设对于一切次数不超过次的多项式,n()q x成立 1()()0bnaq xx dx+=又设()f x是次数不超过次的多项式,用去除2n+11()nx+()f x,商,余,即()q x()r x1()()()()nf xq xxr x+=+,可知,()q x和均是不超过次的多项式,从而()r xn 1 ()()()()bbbnaaaf x dxq xx dxr x dx+=+()bar x dx=，27又因求积公式是插值多项式的构造导出的,由01,nA AA?的选取,其代数精确度可以达到,而是次数不超过次的多项式,因此成立 n()r xn 0()()nbkkakr x dxA r x=因,所以1()0nkx+=()()kkf xr x=,0,1,kn=?故而也即 0()()nbkkakr x dxA f x=0()()nbkkakf x dxA f x=由于()f x是次数不超过次的多项式,因此该积分公式的代数精确度至少为,因而节点是高斯点。2n+112n+01,nx xx?对于关系 1()()0bnaq xx dx+=,我们称之为正交性,即与任意次多项式正交,而这样的多项式类称为正交多项式。1()nx+n3 高斯勒让德求积公式 Legendre 多项式 21()(1),0,1,2,22!nnnnndP xxndx=?称为勒让德(Legendre)多项式。其具有前面提到的正交性质,即对于任意次数不超过的多项式,成立 n()q x 11 1()()0nq x Px dx+=。因此,多项式1()nPx+的零点就是相应的高斯勒让德求积公式的高斯点。勒让德多项式的前几项如下:0()1P x=,()211()12dP xxdx=x=,222222211()(1)(31),22!2dP xxxdx=28()323333311()(1)53,2 3!2dP xxxdx=x()4241()353038P xxx=+()5351()657015,8P xxxx=+?例 例 用四点(n=3)的Gauss求积公式计算 解 解 先将区间变换为-1，1，令 其中 (精确解=0.467 401)这结果与用n=8 的Romberg求积相当.295.5 总 结 二、NewtonCotes 求积公式 将区间a,b等分n等份，求积节点为xk=a+kh,k=0,1,.,n，其中节点间距bahn=。()0()()()nbnkkakf x dxbaCf x=其中，由n次拉格朗日插值公式的插值基函数l()nkCk(x)决定:。()()()bnkkaba Cl x dx=当 n=1 时,得到梯形公式:()()()2babaf x dxf af b+误差为3()()12baf；当 n=2 时，得到抛物线（辛浦生）公式：()()4()()62babaabf x dxf aff b+误差为5(4)()()2880baf；当 n=4 时，得到柯特斯公式：()baf x dx01347()32()12()32()7()90b a5f xf xf xf xf x+误差为7(6)()()30240baf；30三、复合公式 1.复合梯形公式 1 1 0()()()2nakkbkhf x dxf xf x+=+011()2()2()()2nnhf xf xf xf x+?误差 2()()12baRh f=；2.复合辛浦生公式 b ahn=,n=2m 为偶数,1 0()3mabkhf x dx=()22 12 20123222124424243kkkmmhmffffffffff+?误差 4(4)()()180baRh f=。例：用复合公式计算，精确值 1 04 arctgx dx 12 04ln(1)2ln2 1.7552983xxarctgxax=+=由()4f xarctg=x，将区间0,1分成 4 等份，n=4，h=0.25，得数表：3100.250.50.751()0 0.9799146 1.8545904 2.5740044 3.1415927xf x 用复合梯形公式:401234)(2222hTfffff=+1.7448264=，误差R=0.0104719；用复合抛物线公式：40123(4243hSfffff=+4)=1.7553750，误差 R=0.0000767。四、高斯型求积公式 1.Legendre 多项式 21()(1),0,1,2,22!nnnnndP xxndx=?对于任意次数不超过 n 的多项式 q(x)，成立 11 1()()0nq x Px dx+=。2.高斯勒让德求积公式 勒让德多项式1()nPx+的零点是高斯勒让德求积公式的求积节点。0()()nbkkakf x dxA f x=当 n=1 时 1 111()33f x dxff+，其代数精确度为 3。当 n=2 时,()1 113()580595f x dxfff35+其代数精确度为 5。一般地，n+1 个求积节点的高斯型求积公式的代数精确度为 2n+1。32对于一般的有界区间a,b上的定积分，可以通过变量代换转化为区间-1,1上的定积分。即在积分 ()baf x dx中，令22babaxt+=+。例：用 n=2 的高斯勒让德公式计算。1 04 arctgx dx精确值。12 04ln(1)1.7552983xxarctgxax=+=此处，a=0，b=1，作变换1122xt=+，则 1 1 0 1114 2()22arctgx dxarctgtdt=+于是 11()2()22f tarctgt=+，故 f(0)=0.9272652，3()0.22445625f=，3()1.45150625f=，所以原式()13580595fff+35=1.7553526，误差 R=0.0000543。三个点的高斯型公式比五个点的复合抛物线公式还要精确，误差还要小。33
展开阅读全文