科学计算与数学建模第三章课件.ppt

1、科学计算与数学建模中南大学数学科学与计算技术学院中南大学数学科学与计算技术学院湘江水流量估计的实际意义湘江水流量估计的实际意义 3.1数值求积的数值求积的Newton-CotesNewton-Cotes（牛顿（牛顿-柯特斯）方法柯特斯）方法 3.2 RombergRomberg（龙贝格）算法（龙贝格）算法3.3GaussGauss（高斯）型求积公式与测量位置的优化选取（高斯）型求积公式与测量位置的优化选取 3.4第三章第三章 湘江流量估计模型湘江流量估计模型3.1 3.1 湘江水流量估计的实际意义湘江水流量估计的实际意义 水流量是水文特征值的一个重要指标，而水文特征值对于水资水流量是水文特征值

2、的一个重要指标，而水文特征值对于水资源的合理利用，防洪以及抗旱具有指导性的作用，因此湘江水流量源的合理利用，防洪以及抗旱具有指导性的作用，因此湘江水流量估计对于湘江流域的社会经济和人民生活具有重大的影响。现根据估计对于湘江流域的社会经济和人民生活具有重大的影响。现根据实际测量得到湘江某处河宽实际测量得到湘江某处河宽700m700m，其横截面不同位置某一时刻的水，其横截面不同位置某一时刻的水深如表深如表3.1.13.1.1所示。若此刻湘江的流速为所示。若此刻湘江的流速为0.5m/s0.5m/s，试估计湘江此刻，试估计湘江此刻的流量。要计算湘江水流量就需要知道其横截面面积，如果知道此的流量。要计算

3、湘江水流量就需要知道其横截面面积，如果知道此处江的水深曲线函数，则其横截面面积为处江的水深曲线函数，则其横截面面积为 。但是在实际中。但是在实际中是不可能精确得到的，那么怎样求出足够高精度的横截面面积的近是不可能精确得到的，那么怎样求出足够高精度的横截面面积的近似值。似值。()bah x dx表表3.1.1 湘江某处横截面不同位置的水深数据湘江某处横截面不同位置的水深数据 单位：单位：m x050100150200250300350400450500550600650700h(x)4.25.95.85.24.55.755.54.85.94.15.14.65.7,4.7 3.1.1 3.1.1

4、数值求积数值求积的必要性的必要性 在高等代数中，曾用牛顿在高等代数中，曾用牛顿-莱布尼兹（莱布尼兹（Newton-Leibniz）公式：）公式：（其中（其中F(X)是是f(x)的一个原函数）来计算定积分。但是，在工程技术的一个原函数）来计算定积分。但是，在工程技术和科学研究中，常常遇到如下情况：和科学研究中，常常遇到如下情况：v f(x)的结构复杂，求原函数困难；的结构复杂，求原函数困难；v f(x)的原函数不能用初等函数表示；的原函数不能用初等函数表示；v f(x)的精确表达式不知道，只给出了一张由实验提供的函数表。的精确表达式不知道，只给出了一张由实验提供的函数表。对于这些情况，要计算积分

5、的精确值都是十分困难的，这就要求建对于这些情况，要计算积分的精确值都是十分困难的，这就要求建立积分的近似计算方法。此外，积分的近似计算又为其它一些数值立积分的近似计算方法。此外，积分的近似计算又为其它一些数值计算，例如微分方程数值解、积分方程数值解等提供了必须的基础。计算，例如微分方程数值解、积分方程数值解等提供了必须的基础。()()()()bbaaf x dxf xF bF a 其中其中 为为n次插值基函数。用次插值基函数。用 近似代替被近似代替被积函数积函数f(x)，则得：，则得：3.1.2 3.1.2 构造数值求积公式的基本方法构造数值求积公式的基本方法v 可以从不同的角度出发，通过各种

6、途径来构造数值求积公式。但常用的一可以从不同的角度出发，通过各种途径来构造数值求积公式。但常用的一个方法是，利用插值多项式来构造数值求积公式。具体做法如下：个方法是，利用插值多项式来构造数值求积公式。具体做法如下：在积分区间在积分区间a,ba,b上取一组点：上取一组点：作作f(x)f(x)的的n n次插值多项式：次插值多项式：0()()()nnkkkLxfxlx01naxxxb()(0,1,)klxkn0()()()()nbbbnkkaaakf x dxL x dxf xl x dx（3.1.1）()nL x011011()()()()()()()()()bbkknkkaakkkkkknx x

7、x xx xx xAl xdxdxxxxxxxxx若记若记 得数值求积公式：得数值求积公式：0()()nbkkakf x dxA f x这样的求积公式称为这样的求积公式称为机械求积公式机械求积公式。（3.1.2）（3.1.3）其中其中 称为求积节点，称为求积节点，称为求积系数。若求积公式（称为求积系数。若求积公式（3.1.33.1.3）中的求积系数中的求积系数 是由（是由（3.1.23.1.2）确定的，则称该求积公式为插值型）确定的，则称该求积公式为插值型求积公式。求积公式。本章主要讨论插值型求积公式。本章主要讨论插值型求积公式。kxkAkA3.1.3 3.1.3 求积公式的余项求积公式的余项

8、 积分积分 的真值与由某求积公式给出的近似之差，称为该的真值与由某求积公式给出的近似之差，称为该求积公式的余项，记作求积公式的余项，记作 。例例3.1.1 3.1.1 求积公式（求积公式（3.1.33.1.3）的余项为：）的余项为：()baf x dx R f 0()()nbkkakR ff x dxA f x()()()()bbbnnaaaR ff xdxL xdxf x L x dx(1)1()()(1)!nbnafR fwx dxn如果求积公式（如果求积公式（3.1.33.1.3）是插值型的，则由上知：）是插值型的，则由上知：101()()()(),(,)nnxx x x xx xab

9、（3.1.43.1.4）为了使一个求积公式能对更多的积分具有良好的实际计算意义，就为了使一个求积公式能对更多的积分具有良好的实际计算意义，就应该要求它对尽可能多的被积函数都准确地成立。在计算方法中，常用应该要求它对尽可能多的被积函数都准确地成立。在计算方法中，常用代数精度这个概念来描述。代数精度这个概念来描述。定义定义3.1.1 3.1.1 若求积公式：若求积公式：对任意不高于对任意不高于m m次的代数多项式都准确成立，而对于次的代数多项式都准确成立，而对于 却不能准却不能准 确成立，则称该公式的代数精度为确成立，则称该公式的代数精度为m m。0()()nbkkakfx dxA fx1mx 例

10、例3.1.2 3.1.2 梯形梯形公式公式 的代数精度的代数精度m=1。一个求积公式的代数精度越高，它就越能对更多的被积函数一个求积公式的代数精度越高，它就越能对更多的被积函数f（x）准确（或较准确）地成立，从而具有更好的实际计算意义。由插值型求准确（或较准确）地成立，从而具有更好的实际计算意义。由插值型求 积公式的余项（积公式的余项（3.1.4）易得：）易得：定理定理3.1.1 3.1.1 含有含有n+1n+1个节点个节点 的插值型求积公式的插值型求积公式（3.1.33.1.3）的代数精度至少为）的代数精度至少为n n。()()()2babafx dxfaf b（3.1.53.1.5）3.1

11、.4 3.1.4 求积公式的代数精度求积公式的代数精度 (0,1,)kxkn3.2 3.2 数值求积的数值求积的Newton-CotesNewton-Cotes（牛顿（牛顿-柯特斯）方法柯特斯）方法 在在3.1 3.1 中，介绍了插值型求积公式及其构造方法。在实际应用中，介绍了插值型求积公式及其构造方法。在实际应用时，考虑到计算的方便，常将积分区间等分之，并取分点为求积分节时，考虑到计算的方便，常将积分区间等分之，并取分点为求积分节点。这样构造出来的插值型求积公式就称为牛顿点。这样构造出来的插值型求积公式就称为牛顿-柯特斯（柯特斯（Newton-Newton-CotesCotes）公式。）公式

12、。本节在介绍一般牛顿本节在介绍一般牛顿-柯特斯公式的基础上，介绍几个常用的牛顿柯特斯公式的基础上，介绍几个常用的牛顿-柯特斯公式以及这些公式在实际计算时的用法。柯特斯公式以及这些公式在实际计算时的用法。3.2.1 Newton-Cotes3.2.1 Newton-Cotes（牛顿（牛顿-柯特斯）公式柯特斯）公式 若将积分区间若将积分区间a,bna,bn等分，取分点作为求积节点，并作变量替换等分，取分点作为求积节点，并作变量替换 ,那么插值型求积公式那么插值型求积公式(3.1.3)(3.1.3)的系数由（的系数由（3.1.23.1.2）可得：）可得：(;0,1,)kb axakh hknnxat

13、h00(1)(1)(1)()!(1)()!(1)(1)(1)(1)()!()!nkn kn knt ttktktnAhdtknkbat ttktktn dtnk nk记记()0(1)(1)(1)(1)()!()!n knnkCt ttktktn dtn knk（3.2.1）则则于是，由（于是，由（3.1.3）就可写出相应的插值型求积公式：）就可写出相应的插值型求积公式：()()nkkAba C()0()()()nbnkkakfx dxbaCfx（3.2.2）这就是一般的牛顿这就是一般的牛顿柯特斯公式，其中柯特斯公式，其中 称为柯特斯系数。称为柯特斯系数。从柯特斯系数的算式（从柯特斯系数的算式（

14、3.2.13.2.1）可以看出，其值与积分区间）可以看出，其值与积分区间a,ba,b及及被积函数被积函数f(x)f(x)都无关，只要给出了积分区间的等分数都无关，只要给出了积分区间的等分数n n，就能毫无困难，就能毫无困难地算出地算出 为了便于应用，部分柯特斯系数列见表为了便于应用，部分柯特斯系数列见表3-13-1()nkC()()()01,nnnnCCC例例3.2.2 3.2.2 当当n=1n=1时有两点公式：时有两点公式：当当n=2n=2时有三点公式时有三点公式当当n=4n=4时有五点公式：时有五点公式：()()()2babaf x dxf af b()()4()()62babaabf x

15、 dxf aff b（3.2.4）01234()7()32()12()32()7()90bab af xdxf xf xf xf xf x（3.2.5）(0,1,)4kbaxakkn其中其中v 求积公式（求积公式（3.2.3）就是梯形公式。）就是梯形公式。v 求积公式（求积公式（3.2.4）称为辛普生（）称为辛普生（Simpson）公式其几何意义就是通）公式其几何意义就是通过三点的抛物线围成的曲边梯形面积近似地代替原曲边梯形面积（见过三点的抛物线围成的曲边梯形面积近似地代替原曲边梯形面积（见图图3.2.1）。）。因此，求积公式（因此，求积公式（3.2.4）又名抛物线公式。求积公式）又名抛物线公

16、式。求积公式（3.2.5）称为柯特斯公式。）称为柯特斯公式。v 梯形公式、辛普生公式和柯特斯公式，是三个最基本、最常用的等梯形公式、辛普生公式和柯特斯公式，是三个最基本、最常用的等距节点下的求积公式。距节点下的求积公式。v 下述定理给出了这些求积公式的余项。下述定理给出了这些求积公式的余项。定理定理3.2.1 3.2.1 若若 在在a,ba,b上连续，则梯形公式（上连续，则梯形公式（3.2.33.2.3）的余项为：）的余项为：若若 在上连续，则辛普生公式（在上连续，则辛普生公式（3.2.43.2.4）的余项为：）的余项为：若若 在上连续，则柯特斯公式在上连续，则柯特斯公式(3.2.5)(3.2

17、.5)的余项为：的余项为：其中其中()fx（3.2.6）31()12baRff(4)()fx（3.2.7）4521()902baRff(6)()fx 6748()9454baRff（3.2.8）,a b3.2.2 3.2.2 复合复合Newton-CotesNewton-Cotes（牛顿（牛顿-柯特斯）公式柯特斯）公式 由定理由定理3.2.1知，当积分区间较大时，直接使用牛顿知，当积分区间较大时，直接使用牛顿-柯特斯公式柯特斯公式所得积分近似值的精度是很难得到保证的。因此在实际应用中，为了所得积分近似值的精度是很难得到保证的。因此在实际应用中，为了既能提高结果的精度，又使算法简便且易在电子计算

18、机上实现，往往既能提高结果的精度，又使算法简便且易在电子计算机上实现，往往采用复合求积的方法。所谓复合求积，就是先将积分区间分成几个小采用复合求积的方法。所谓复合求积，就是先将积分区间分成几个小区间，并在每个小区间上用低阶牛顿区间，并在每个小区间上用低阶牛顿柯特斯公式计算积分的近似值，柯特斯公式计算积分的近似值，然后对这些近似值求和，从而得到所求积分的近似值。由此得到的一然后对这些近似值求和，从而得到所求积分的近似值。由此得到的一些具有更大实用价值的数值求积公式，统称为复合求积公式。些具有更大实用价值的数值求积公式，统称为复合求积公式。v 例例3.2.3 先将区间先将区间a,bn等分，记分点为

19、等分，记分点为 ，其中其中 ，称为步长，然后在每个小区间称为步长，然后在每个小区间 上应用梯形公式（上应用梯形公式（3.2.3），即：），即：就可以导出复合梯形公式：就可以导出复合梯形公式：,(0,1,)kxa kh kn b ahn1,iixx11()(1,2,)2kkxkkxhf x dxfxf xkn1111()()()()2kknnbxkkaxkkhf x dxf x dxf xf x若将所得积分近似值记成若将所得积分近似值记成 ，并注意到，并注意到 ，则上式即为：，则上式即为：仿上，可得复合辛普生公式：仿上，可得复合辛普生公式：和复合柯特斯公式：和复合柯特斯公式：nT0,nxa xb

20、11()()2()()2nbnkakhf x dxTf af xf b（3.2.9）111012()4()2()()6nnbnkakkkhf x dxSf af xf xf b（3.2.10）其中其中定理定理3.2.2 若若f(x)在积分区间在积分区间a,b上分别具有二阶、四阶和六阶连续导上分别具有二阶、四阶和六阶连续导 数数,则符合求积公式则符合求积公式(3.2.9)、(3.2.10)和和(3.2.11)的余项分别为：的余项分别为：11041111300024()732()9012()32()14()7()nbnakknnnkkkkkkhf x dxCf af xf xf xf xf b（3

21、.2.11）113424113,424kkkkkkxxhxxhxxh其中其中 ，且当，且当h充分小时，又有：充分小时，又有：2()()12bnabaf x dxTh f 4(4)()()()1802bnaba hf x dxSf（3.2.14）662()()()()9452bnabahf x dxCf(3.2.12)（3.2.13）,a b21()()()12bnaf x dxThfbfa 41()()()()180 2bnahf x dxSfbfa 6(5)(5)2()()()()945 4bnahf x dx Cfbfa(3.2.15)（3.2.16）（3.2.17）定义定义3.2.1 对

22、于复合求积公式对于复合求积公式 ，若当若当h0时有时有 则称则称 是是p阶收敛的。阶收敛的。定理定理3.2.3 复合求积公式（复合求积公式（3.2.9）、（）、（3.2.10）和（）和（3.2.11）分别具有二）分别具有二阶、四阶和六阶收敛性。阶、四阶和六阶收敛性。证明证明 由收敛性的定义，从（由收敛性的定义，从（3.2.19）可以看出，复合梯形公式（）可以看出，复合梯形公式（3.2.9）具有二阶收敛性。同样，可证明复合辛普生公式（具有二阶收敛性。同样，可证明复合辛普生公式（3.2.10）和复合柯）和复合柯斯特公式（斯特公式（3.2.11）分别具有四阶和六阶收敛性。）分别具有四阶和六阶收敛性。

23、()bnaf xI()(0)bnapf x dxIc chnIv 3.2.3 3.2.3 误差的事后估计与步长的自动选择误差的事后估计与步长的自动选择 虽然可用余项公式（虽然可用余项公式（2.122.12）（）（2.172.17）来估计近似值的误差，也可以根据精）来估计近似值的误差，也可以根据精度要求用这些公式来确定积分区间的等分数，即确定步长度要求用这些公式来确定积分区间的等分数，即确定步长h h。但由于余项公式中包。但由于余项公式中包含被积函数含被积函数f(x)f(x)的高阶导数，在具体计算时往往会遇到困难。因此，在实际应用的高阶导数，在具体计算时往往会遇到困难。因此，在实际应用时，常常利

24、用误差的事后估计法来估计近似值的误差或步长时，常常利用误差的事后估计法来估计近似值的误差或步长h h。该方法的大致做法。该方法的大致做法是：是：将积分区间逐次分半，每分一次就用同一复合求积公式算出相应的积分近似将积分区间逐次分半，每分一次就用同一复合求积公式算出相应的积分近似值，并用前后两次计算结果来判断误差的大小。值，并用前后两次计算结果来判断误差的大小。其原理和具体做法是：其原理和具体做法是：对于复合梯形公式（对于复合梯形公式（3.2.93.2.9），由余项公式（），由余项公式（3.2.123.2.12）或（）或（3.2.153.2.15）可以看出，）可以看出，当当f(x)f(x)在积分区

25、间上变化不大或积分区间在积分区间上变化不大或积分区间a,ba,b的等分数的等分数n n较大（即步长较大（即步长h h较小）较小）时，若将时，若将a,ba,b的等分数改为的等分数改为2n2n（即将步长缩小到原步长（即将步长缩小到原步长h h的一半），则新近的一半），则新近 似值似值T2nT2n的余项约为原近似值余项的的余项约为原近似值余项的1/41/4，即：，即：214nnITIT 其中其中I表示积分表示积分 的真值。的真值。对对I求解得：求解得：此式表明，若用此式表明，若用T2n作为积分真值作为积分真值I的近似值，则其误差约为的近似值，则其误差约为()baf x221()3nnnITTT(3.

26、2.20)21()3nnTT故在将区间逐次分半进行计算的过程中，可以用前后计算结果和来估故在将区间逐次分半进行计算的过程中，可以用前后计算结果和来估计误差与确定步长。计误差与确定步长。先算出先算出TnTn和和T2nT2n，若，若 （为计算结果的允许误（为计算结果的允许误差），则停止计算，并取差），则停止计算，并取T2nT2n作为积分的近似值；否则，将区间再次作为积分的近似值；否则，将区间再次分半后算出新近似值分半后算出新近似值T4nT4n，并检查不等式，并检查不等式 是否成立，是否成立，直到得到满足精度要求的结果为止。直到得到满足精度要求的结果为止。对于复合辛普生公式（对于复合辛普生公式（3.

27、2.103.2.10）和复合柯特斯公式（）和复合柯特斯公式（3.2.113.2.11），），当所涉及的高阶导数在积分区间上当所涉及的高阶导数在积分区间上23 nnTT42nnTT变化不大或积分区间的等分数较大时，由相应的余项公式可以看出：变化不大或积分区间的等分数较大时，由相应的余项公式可以看出：分别对分别对I求解得求解得 2116nnISIS2164nnICIC221()15nnnISSS221()63nnnICCC(3.2.21)(3.2.22)因此，估计新近似值和的误差，并判断计算过程是否需要继续进因此，估计新近似值和的误差，并判断计算过程是否需要继续进行，也可以像使用复合梯形法求积分近

28、似值那样，在将积分区间逐次行，也可以像使用复合梯形法求积分近似值那样，在将积分区间逐次分半进行计算的过程中下去分半进行计算的过程中下去 .3.2.4 3.2.4 复合梯形法的递推算式复合梯形法的递推算式 上段介绍的变步长的计算方案，虽然提供了估计误差与选取步长上段介绍的变步长的计算方案，虽然提供了估计误差与选取步长的简便方法，但还没有考虑到避免在同一节点上重复计算函数值的问的简便方法，但还没有考虑到避免在同一节点上重复计算函数值的问题，故有进一步改进的余地。题，故有进一步改进的余地。先看复合梯形公式。先看复合梯形公式。在利用（在利用（3.2.93.2.9）计算）计算TnTn时，需要计算时，需要

29、计算n+1n+1个点（它们是积分区间个点（它们是积分区间a,bna,bn等分点的分点，不妨简称为等分点的分点，不妨简称为“n n分点分点”）上的函数值。当）上的函数值。当TnTn不满足精度要求时，根据上面提供的计算方案，就应将各个小不满足精度要求时，根据上面提供的计算方案，就应将各个小区间分半，计算出新近似值。若利用（区间分半，计算出新近似值。若利用（3.2.93.2.9）进行计算）进行计算T2nT2n，就，就需要求出需要求出2n+12n+1个点（它们是个点（它们是“2n2n分点分点”）上的函数值。而实际上，）上的函数值。而实际上，在这在这2n+12n+1个个2n2n分点中，包含有分点中，包含

30、有n+1n+1个个n n分点，对应的函数值在计算分点，对应的函数值在计算TnTn时早已算出。为了避免这种重复计算，下面分析近似值时早已算出。为了避免这种重复计算，下面分析近似值T2nT2n与与原有近似值原有近似值TnTn之间的联系。之间的联系。由复合梯形公（由复合梯形公（3.2.9）知式：）知式：若注意到在若注意到在2n分点：分点：中，当中，当k取偶数时是取偶数时是n分点，当分点，当k取奇数时，才是新增加的分点。将新取奇数时，才是新增加的分点。将新增加的分点处的函数值从求和记号中分离出来，就有：增加的分点处的函数值从求和记号中分离出来，就有：2121()2()()42nnkbabaTf af

31、akf bnn(1,2,21)2kbaxakknn即即1211111()2(2)2(21)()422()2()()(21)422nnnkknnkkbababaTf af akf akf bnnnbabababaf af af bf aknnnn（3.2.23）211(21)222nnnkbabaTTf aknn由递推公式（由递推公式（3.2.233.2.23）可以看出，在已经算出）可以看出，在已经算出TnTn的基础上再计算的基础上再计算T2nT2n时，只要计算时，只要计算n n个新分点上的函数值就行了。与直接利用复合梯个新分点上的函数值就行了。与直接利用复合梯形公式（形公式（3.2.93.2.

32、9）求）求T2nT2n相比较，计算工作量几乎节省了一半相比较，计算工作量几乎节省了一半。3.3 Romberg3.3 Romberg（龙贝格）算法（龙贝格）算法龙贝格（龙贝格（Romberg)Romberg)算法是在积分区间逐次分半的过程中，对用复合算法是在积分区间逐次分半的过程中，对用复合梯形产生的近似值进行加权平均，以获得准确程度较高的一种方法，梯形产生的近似值进行加权平均，以获得准确程度较高的一种方法，具有公式简练、使用方便、结果较可靠等优点。本节介绍它的基本原具有公式简练、使用方便、结果较可靠等优点。本节介绍它的基本原理和应用方法。理和应用方法。v 3.3.1 Romberg3.3.1

33、 Romberg（龙贝格）算法的基本原理（龙贝格）算法的基本原理 上节中介绍的递推公式（上节中介绍的递推公式（3.2.233.2.23）或（）或（3.2.243.2.24），虽然具有结构简），虽然具有结构简单，易在电子计算机上实现等优点，但是由它产生的梯形序列单，易在电子计算机上实现等优点，但是由它产生的梯形序列 ，其收敛速度却是非常缓慢的。其收敛速度却是非常缓慢的。例例3.3.1 3.3.1 用此法计算用此法计算 的近似值时，要一直算到的近似值时，要一直算到 才获得误差不超过才获得误差不超过 的近似值（见例的近似值（见例3.2.53.2.5）。因此，用这种方法计）。因此，用这种方法计算更复杂

34、的高精度要求的积分近似值显然是费时、费力甚至是不可能的。算更复杂的高精度要求的积分近似值显然是费时、费力甚至是不可能的。2kT12041dxx512T610如何提高收敛速度，以节约计算工作量，自然是人们极为关心的课题。如何提高收敛速度，以节约计算工作量，自然是人们极为关心的课题。由近似等式（由近似等式（3.2.20），用），用T2n作为积分真值作为积分真值I的近似值，其误差约的近似值，其误差约为为 。因此，如果用因此，如果用 作为作为T2n的一种补偿，可以期望所得到的新的一种补偿，可以期望所得到的新近似值：近似值：有可能有可能 比更好地接近于比更好地接近于 积分的真值。积分的真值。21()3n

35、nTT21()3nnTT221()3nnnnTTTT2nT()baf x dx如在例如在例3.2.53.2.5中，中，和和 是两个精度很差的是两个精度很差的近似值，但如果将它们按（近似值，但如果将它们按（3.3.13.3.1）作线性组合，所得到的近似值：）作线性组合，所得到的近似值：却具有七位有效数字，其准确程度比却具有七位有效数字，其准确程度比 还要高，而计算还要高，而计算 只涉及求九只涉及求九个点上的函数值，其计算工作量仅为计算个点上的函数值，其计算工作量仅为计算 的的 。那么，按（那么，按（3.3.13.3.1）式作线性组合得到的新近似值）式作线性组合得到的新近似值 ，其实质又是什，其实

36、质又是什么呢？通过直接验证，易知么呢？通过直接验证，易知 ，亦即亦即 ：43.13117647T 83.13898849T 484413.1415925033TTT512T4T512T157nTnnTSnS22441334 1nnnnnTTSTT（3.3.2）这表明在收敛速度缓慢的梯形序列这表明在收敛速度缓慢的梯形序列 基础上，若将基础上，若将 与与 按按（3.3.2）作线性组合，就可产生收敛速度较快的辛普生序列）作线性组合，就可产生收敛速度较快的辛普生序列 ：同理，从近似等式（同理，从近似等式（3.2.21）出发，通过类似的分析，可以得到：）出发，通过类似的分析，可以得到：故在辛普生序列故在

37、辛普生序列 的基础上，将的基础上，将 与与 按（按（3.3.3）作线性组合，）作线性组合，就可产生收敛速度更快的柯特斯序列就可产生收敛速度更快的柯特斯序列 ：这种加速过程还可以继续下去。这种加速过程还可以继续下去。2kTnT2nT2kS124,S SS nS2nS2kC124,.C C C 2kS22224161151541nnnnnSSCSS（3.3.3）例例3.3.2 3.3.2 通过通过 与与 的线性组合，可以在柯特斯序列的线性组合，可以在柯特斯序列 的基础上，的基础上，产生一个称为龙贝格序列的新序列产生一个称为龙贝格序列的新序列 ，即：，即：经过进一步的分析，可以证明，当经过进一步的分

38、析，可以证明，当 满足一定条件时，龙贝格序列满足一定条件时，龙贝格序列 比柯特斯序列比柯特斯序列 更快地收敛到积分更快地收敛到积分 的真值的真值 。综上可知，可以在积分区间逐次分半的过程中利用公式（综上可知，可以在积分区间逐次分半的过程中利用公式（3.3.23.3.2）、）、（3.3.33.3.3）和（）和（3.3.43.3.4），将粗糙的近似值，将粗糙的近似值 逐步逐步“加工加工”成越来越精成越来越精确的近似值确的近似值 。也就是说，将收敛速度缓慢的梯形序列。也就是说，将收敛速度缓慢的梯形序列 逐步地逐步地“加工加工”成收敛速度越来越快的新序列：成收敛速度越来越快的新序列：。这种。这种加速的

39、方法就称为龙贝格算法。其加工过程如图加速的方法就称为龙贝格算法。其加工过程如图3.3.13.3.1，其中圆圈中的号码，其中圆圈中的号码表示计算顺序。表示计算顺序。nC2nC2kC2kR32234641636341nnnnnCCRCC（3.3.4）()f x2kR2kC()baf x dxInT,nnnS CR 2kT 222,kkkSCR图3.3.1v 3.3.2 Romberg3.3.2 Romberg（龙贝格）算法计算公式的简化（龙贝格）算法计算公式的简化 为了便于上机计算，引用记号为了便于上机计算，引用记号 来表示各近似值，其中来表示各近似值，其中k k仍代表积分区间的二分次数，而下标仍

40、代表积分区间的二分次数，而下标m m则指出了近似值所在序列的性则指出了近似值所在序列的性质：质：当当m=0m=0时，在梯形序列中；当时，在梯形序列中；当m=1m=1时，在辛普生序列中；当时，在辛普生序列中；当m=2m=2时时在柯特斯序列中在柯特斯序列中 .()kmT例例3.3.4 表表3.3.1中的各近似值，若用中的各近似值，若用 记号表示，则如表记号表示，则如表3.3.2所所示。示。()kmT引入上面的记号后，龙贝格算法所用到的各个计算公式可以统一为：引入上面的记号后，龙贝格算法所用到的各个计算公式可以统一为：（3.3.5）10021001111()()21211,2,3,22240,1,2

41、,;1,2,3,41lllllikmkkmmmmbaTf af bbabaTTfailTTTkm相应的程序框图见图相应的程序框图见图3.3.2。其中。其中 为最大二分次数。最后指出下列几点：为最大二分次数。最后指出下列几点：（1）当）当 较大时，由（较大时，由（3.3.5）第三式知）第三式知 。因此，在实际计。因此，在实际计算中，常规定算中，常规定 ，即在计算到出现龙贝格序列为止。在这种情况下，即在计算到出现龙贝格序列为止。在这种情况下，程序框图程序框图3.3.2应做相应的修改，需将应做相应的修改，需将“按式（按式（3.3.5）计算）计算”改为改为“按式（按式（3.3.5）计算）计算”，并将精

42、度控法，并将精度控法 改为改为 （2）为防止假收敛，可设置最小二分次数）为防止假收敛，可设置最小二分次数 。当。当 时，跳时，跳过精度判别而继续运算；过精度判别而继续运算；（3）可以用二维数组来存放与参加运算，也可用一维数组。可以用二维数组来存放与参加运算，也可用一维数组。0km()(1)1kkmmTT3m()(1)(0)01,kkkTTT()(1)(2)(3)0123,kkkkTTTT(0)(0)1kkTT3433kkTT（）（）minKminKK()kmT图3.3.23.4 Gauss(3.4 Gauss(高斯高斯)型求积公式与测量位置的优化选取型求积公式与测量位置的优化选取 下面介绍一种

43、高精度的求积公式下面介绍一种高精度的求积公式高斯（高斯（GaussGauss）型求积公式。）型求积公式。3.4.1 Gauss3.4.1 Gauss（高斯）型求积公式的定义（高斯）型求积公式的定义 在在3.23.2中，限定把积分区间的等分点作为求积节点，从而构造出中，限定把积分区间的等分点作为求积节点，从而构造出一类特殊的插值型求积公式，即牛顿一类特殊的插值型求积公式，即牛顿-柯特斯公式。这种做法虽然简柯特斯公式。这种做法虽然简化了计算，但却降低了所得公式的代数精度。化了计算，但却降低了所得公式的代数精度。例例3.4.1 在构造形如在构造形如 的两点公式时，如果限定求积节点的两点公式时，如果限

44、定求积节点 ，那么所得插值，那么所得插值型求值公式：型求值公式：的代数精度仅为的代数精度仅为1。100111()()()f x dxA f xA f x(3.4.1)01x 11x 11()(1)(1)f x dxff(3.4.2)但是，如果我们对（但是，如果我们对（3.4.1）中的系数）中的系数A0,A1和节点和节点X0,X1都不加限制，都不加限制，那么就可以适当选取那么就可以适当选取A0,A1和和X0,X1，使所得公式的代数精度，使所得公式的代数精度m1。事实上，若只要求求积公式（事实上，若只要求求积公式（3.4.1）对函数）对函数 都准确成立，只要都准确成立，只要 和和 满足方程组满足方

45、程组 23()1,f xx xx01,x x01,A A解之得解之得01001 122001 133001 120230AAA xAxA xAxA xAx（3.4.3）011AA033x 133x 代入（代入（3.4.1）即得：）即得：容易验证，所得公式（容易验证，所得公式（3.4.4）是代数精度的插值型求积公式。）是代数精度的插值型求积公式。同理，对于一般求积公式：同理，对于一般求积公式：1133()()()33f x dxff（3.4.4）0()()nbkkakf x dxA f x（3.4.5）只要适当选择只要适当选择2n+2个待定参数个待定参数 和和 使它的代数使它的代数精度达到精度达

46、到2n+1也是完全可能的。也是完全可能的。定义定义3.4.1 若形如（若形如（3.4.5）的求积公式代数精度达到了，则称它为高）的求积公式代数精度达到了，则称它为高斯型求积公式，并称相应的求积节点为高斯点。斯型求积公式，并称相应的求积节点为高斯点。kx(0,1,)kAkn 3.4.2 Gauss（高斯）求积公式的构造与应用（高斯）求积公式的构造与应用 可以像构造两点高斯型求积公式（可以像构造两点高斯型求积公式（3.4.4）那样，通过解形如）那样，通过解形如（3.4.3）的方程组来确定高斯点）的方程组来确定高斯点 和和 ，从而，从而构造构造n+1点高斯型求积公式。点高斯型求积公式。但是，这种做法

47、要解一个包含有但是，这种做法要解一个包含有2n+2个未知数的非线性方程组，个未知数的非线性方程组，其计算工量是想当大的。一个比较简单的方法是：其计算工量是想当大的。一个比较简单的方法是：kx(0,1,)kAkn（1）先用区间）先用区间a,b上的上的n+1次正交多项式确定高斯点次正交多项式确定高斯点（2）然后利用高斯点确定求积系数）然后利用高斯点确定求积系数 当积分区间是当积分区间是-1,1时，两点至五点高斯型求积公式的节点、系时，两点至五点高斯型求积公式的节点、系数数T和余项见表和余项见表3.4.1，其中，其中 利用表利用表3.4.1，可以方便地写出相应的高斯型求积公式。例如，可以方便地写出相

48、应的高斯型求积公式。例如，当当n=2时，由表时，由表3.4.1知：知：,(0,1,)kxa b kn(0,1,)kA kn 1,1 0,10.57735027x 0,11A表表3.4.1 故得两点高斯型求积公式：故得两点高斯型求积公式：又如，当又如，当n=3时，由表时，由表3.4.1可以查出三个求积节点和对应的三个可以查出三个求积节点和对应的三个系数（注意，系数系数（注意，系数0.55555556应连用两次），从而得到三点高斯型应连用两次），从而得到三点高斯型求积公式：求积公式：对于一般区间对于一般区间a,b上的积分，也可以利用表上的积分，也可以利用表3.4.1写出高斯型求写出高斯型求积公式。

49、积公式。11()(0.57735027)(0.57735027)f x dxff11()0.55555556(0.77459667)0.88888889(0)0.55555556(0.77459667)f x dxfff先作变量替换，令先作变量替换，令将区间上的积分转化为区间将区间上的积分转化为区间-1,1上的积分：上的积分：记记则等式（则等式（3.4.6）右端的积分为）右端的积分为 22abbaxt11()()222babaabbaf x dxft dt()()22abbag tft（3.4.6）11()g t dt利用表利用表3.4.1，对于给定的，对于给定的n=1,2,3,4，可以写出高

50、斯型求积公式：，可以写出高斯型求积公式：即即 代入（代入（3.4.6）式得：）式得：110()()nkkkg t dtA g t110()()2222nkkkabbaabbafdtA ft0()()222nakkbkbaabbaf x dxA ft（3.4.8）其中系数与节点可在表其中系数与节点可在表3.4.1中查得。中查得。由变量替换式由变量替换式 容易看出，由于求积公式（容易看出，由于求积公式（3.4.7）对变量）对变量t不高于不高于2n+1次的多项次的多项式准确成立，从而求积公式（式准确成立，从而求积公式（3.4.8）对变量）对变量x不高于不高于2n+1次的多项式次的多项式也准确成立，即