1、毕业论文文献综述 信息与计算科学 定积分的数值计算方法 一、 前言部分 在科学与工程计算中,经常要计算定积分 ( ) ( ) ( ) .baI f f x d x a b ( 1 1) 这个积分的计算似乎很简单,只要求出 f 的原函数 F 就可以得出积分( 1 1)的值,即 ( ) ( ) ( ).I f F b F a ( 1 2) 如果原函数 F 非常简单又便于使用,那么式( 1 2)就提供 了计算起来最快的积分法但是,积分过程往往将导出新的超越函数,例如,简单积分 1dxx 可引出对数函数,它已不是代数函数了;而积分 2xe dx ,将引出一个无法用有限个代数运算、对数运算或指数运算组合
2、表示的函数有些积分虽然容易求解,并且原函数仍然是一个初等函数,但可能过于复杂,以致于人们采用( 1 2)来计算之前还得三思而 行 1 例如 242 21 1 2 1 2 1a r c ta n ( ) l n112 2 4 2 2 1x x xd x Cxx xx , ( 1 3) 采用式( 1 3)这种“精确”表达式时,所需运算次数是个根本问题由式( 1 3)看出,需计算对数和反正切,因此只能计算到一定的近似程度 因此可以 看出,这类表面上是“精确”的方法,实际上也是近似的因此,我们常常需要探讨一些近似计算定积分的数值方法2 通过人们的研究和发现,得出了很多数值计算的方法,比如利用牛顿 -科
3、茨求积公式,复合求积公式,龙 贝格积分法,高斯求积公式,切比雪夫求积法等来解决定积分的数值计算问题 构造数值积分公式最通常的方法是用积分区间上的 n 次插值多项式代替被积函数,由此 导出的求积公式称为插值型求积公式特别在节点分布等距的情形称为牛顿-柯茨公式,例如梯形公式与抛物线公式就是最基本的近似公式但它们的精度较差龙贝格算法是在区间逐次分半过程中,对梯形公式的近似值进行加权平均获得准确程度较高的积分近似值的一种方法,它具有公式简练、计算结果准确、使用方便、稳定性好等优点,因此在等距情形宜采用龙贝格求积公式当用不等距节点进行计算时,常用高斯型求积公式计算,它在节点数目相同情况下,准确程度较高,
4、稳定性好,而且还可以计算无穷积分 3 二、 主题部分 2 1 牛顿 -科茨求积 公 式 4 2 1 1 公式的一般形式 4 将积分( 1 1)中的积分区间 ,ab 分成 n 等分,其节点 kx 为 1, ( )kx a kh h b an ( 0,1, , )kn 对于给定的函数 f ,在节点 kx ( 0,1, , )kn 上的值 ()kfx 为已知那么 f 在 n+1 个节点 01, , , nx x x 上的 n 次代数插值多项式为 0 0( ) ( ) .nn jnkk j kjjkxxp x f xxx 如果记 x a th ,则上式可以写为 0 0( ) ( ) .nnnkk jj
5、ktjp x f xkj ( 2 1) 在积分( 1 1)中的被积函 数 f 用其 n+1 个节点的代数插值多项式 ()npx来代替,可 得 ( ) ( ) ( ) ( )bbnnaaI f f x d x I f p x d x 多项式的积分是容易求出的,因此把上式写为 0( ) ( ) ( )nn n kkI f I f A f x , ( 2 2) 其中 ()0 0 ( ) ,nn nkkjjkb a t jA d t b a cn k j (2 3) ()0 0( 1 ) ( ) .! ( ) !nk nnnk jjkc t j d tk n k n (2 4) 公式( 2 2)称为
6、牛顿 -科茨求积公式 或称为 等距节点求积公式 , kA 称为求积公式系数,()nkc 称为科茨求积系数 牛顿 -科茨求积公式的误差估计 ()nEf ( ) ( )nI f I f,由下面定理给出 定理 2 1 (1) 如果 n 为偶数, ( 2)nf 在 ,ab 上连续,则有 3 ( 2 )( ) ( ) , ,nnnnE f c h f a b, ( 2 5) 其中 201 ( 1 ) ( 2 ) ( )( 2 ) ! nnc t t t t n d tn (2) 如果 n 为奇数, ( 1)nf 在 ,ab 上连续,则有 2 ( 1 )( ) ( ) , ,nnnnE f c h f a
7、 b, ( 2 6) 其中 01 ( 1 ) ( 2 ) ( )( 1 ) ! nnc t t t t n d tn 定义 2 1 如果求积公式0( ) ( )nbkka kf x dx A f x 对所有次数不高于 n 的代数多项式等式精确成立,但存在 n+1 次的代数多项式使等式不成立,则称上式求积公式具有 n 次代数精度 由定理 2 1 可知,牛顿 -科茨求积公式( 2 2)的代数 精度至少是 n 次,而当 n 是偶数时,( 2 2)的代数精度可达 n+1 次 2 1 2 梯形公式 5 在牛顿 -科茨公式( 2 2)中,取 n=1时 (1) (1)011 ,2cc所以有 1( ) ( )
8、 ( ) ( ) .2baI f I f f a f b ( 2 7) 公式( 2 7)称为 梯形公式 ,如果用 连接 , ( )a f a 和 , ( )b f b 的直线来逼近 f ,并对这线性函数进行积分可得到 1()If再用 1()If来逼近 ()If 定理 2 2 若 2 ,f C a b ,则梯形公式( 2 7)的误差为 311 1( ) ( ) ( ) ( ) ( ) , , .12E f I f I f b a f a b 2 1 3 辛普森公式 6 在牛顿 -科茨公式( 2 2)中 ,取 n=2,则有 220 011( 1 ) ( 2 ) ,46c t t d t 221 0
9、14( 2 ) ,26c t t dt 222 011( 1) ,46c t t dt 有此得到 2( ) ( ) ( ) 4 ( ) ( ) .32h a bI f I f f a f f b ( 2 8) 其中 1()2h b a式( 2 8)称为辛普森公式 定理 2 3 若 4 ,f C a b ,则辛普森 公式( 2 8)的误差为 5 ( 4 )22 1( ) ( ) ( ) ( ) , , .90E f I f I f h f a b 2 2 复化求积公式 7 上面已经给出了计算积分 ( ) ( )baI f f x dx的 3 个基本的求积公式:梯形公式,辛普森公式,牛顿 -科茨公
10、式,并给出了它们误差的表达式由这些表达式可知其截断误差依赖于求积区间的长度若积分区间的长度是小量的话,则这些求积公式的截断误差是该长度的高阶小量但若积分区间的长度比较大,直接使用这些公式,则精度难以保证为了提高计算积分的精度,可把积分区间分 为若干个小区间, ()If 等于这些小区间上的积分和,然后对每个小区间上的积分应用上述求积公式,并把每个小区间上的结果累加,所得到的求积公式称为 复化求积公式 将积分区间 ,ab 作 n 等分,并记 , , 0 , 1 , ,kbah x a k h k nn ,于是 110( ) ( )kkn xxkI f f x dx 2 2 1 复化梯形求积公式 8
11、 如果需要求出一个已知函数 ()fx在一个很大区间 ,ab 上的积分,那么我们可以把区间分成 n 个长度为 xh 的小区间,对每一个小区间用梯形法则,然后再把这些小区间上的积分值相加于是就得到了计算定积分的复化梯形公式: 11 0 1 2 10( ) ( ) ( 2 2 2 )22nbi i n na i hhf x d x f f f f f f f (2 9) 整体积分误差等于 n 个小区间上的积分误差之和: 整体误差 = 312( ) ( ) ( )12 nh f f f , 其中 i 是第 i 个小区间上的某一点如果 ()fx在区间 ,ab 上连续,那么由连续函数的性质可知,在区间 ,
12、ab 上存在点 使得 ( )if 的平均值等于 ()f 于是由于 nh b a ,有 整体误差 = 3 22( ) ( ) ( )1 2 1 2n h b af h f O h , 局部误差是 3()Oh ,整体误差是 2()Oh 2.2.2 复化辛普森求积公式 9 对于积分 ()ba f xdx,将 ,ab 等分,每个小区间长度 bah n ,节点记为 ( 0 ,1 , 2 , , )kx a kh k n ,第 k 个小区间记为 1 , ( 1, 2 , , )kkx x k n 记 1kkxx 的中点为112 1 ()2 kkkx x x ,则复化辛普森公式为 111 2( ) ( )
13、( ) 4 ( ) ( )6nbkka kkhf x dx S h f x f x f x 2 3 龙贝格积分 10 现在要介绍用龙贝格( Romberg) 命名 的一个算法,龙贝格首先给出了这种算法的递推形式,假设需要积分 ()baI f x dx ( 2 10) 的近似值在讨论过程中函数 ()fx和区间 ,ab 将保持不变 2 3 1 递推梯形法则 10 设 ()Tn 表示在长度是 ( )/h b a n 的 n 个子区间上积分 I 的梯形法则根据0( ) ( )nba if x d x h f a ih, 我们有 00( ) ( ) ( ) ( )nnniib a b aT h f a
14、ih f a inn , ( 2 11) 这里求和符号中的两撇表示和式中第一项和最后一项减半 2 3 2 龙贝格算法 10 在龙贝格算法中使用上述公式设 ( ,0)Rn 表示具有 2n 个子区间的梯形估计,我们有 1211( 0 , 0 ) ( ) ( ) ( )21( , 0 ) ( 1 , 0 ) ( ( 2 1 ) )2nnniR b a f a f bR n R n h f a i h , ( 2 12) 对于一个适度的 M 值,计算 ( 0 , 0 ) , (1 , 0 ) , ( 2 , 0 ) , , ( , 0 )R R R R M,并且其中没有重复的函数值的计算在龙贝格算法的
15、其余部分中,还要计算附加值 ( , )Rnm 所有这些都可以被理解为积分 I 的估计计算出 ( ,0)RM 后,不再需要被积函数 f 值的计算根据公式 1( , ) ( , 1 ) ( , 1 ) ( 1 , 1 )41mR n m R n m R n m R n m , ( 2 13) 对于 1n 和 1m 构造 R 阵列的各列 定理 2 4(龙贝格算法收敛性定理) 10若 ,f C a b ,则龙贝格阵列中每一列都收敛于 f 的积分因此,对每个 m , lim ( , ) ( )baR n m f x dx 2 4 高斯求积 11 前面研究的求积公式都是事先确定了 n 个节点,然后按使求积
16、公式阶数达到最大的原 则选取最佳权由于自由参数为 n 个,所以阶数一般为 n-1,但如果节点的位置也自由选择,则自由参数的个数将变为 2n,因此求积公式的阶数可达到 2n-1高斯求积公式就是通过选择最佳的节点和权,使求积公式的阶数最大化一般地,对每个 n, n 点高斯公式都是唯一的 ,而且阶数为 2n-1因而,对一定的节点个数,高斯求积公式的精度是最高的但它的求得比牛顿 柯特斯公式要困难得多虽然它的节点和权也可由待定系数法确定,但得到的方程是非线性的 2 4 1 高斯求积公式 11 为说明高斯求积公式,推导区间 1,1 上的两点公式 1 1 1 2 2 21( ) ( ) ( ) ( ) (
17、)I f f x d x w f x w f x G f , 其中的节点 1x 、 2x 及权 1w 、 2w 按使求积公式阶数最大化的原则选取令公式对前四个单项式精确成立,得力矩方程组 112 111 1 2 2 112 2 21 1 2 2 113 3 31 1 2 2 11 2 ,0,2,30.w w d xw x w x xd xw x w x x d xw x w x x d x 这个非线性方程组的一个解为 1 2 1 21 3 , 1 3 , 1 , 1 ,x x w w 另 一 个 解 可 通 过 改 变 1x , 2x 的 符 号 而 得 到 这 样 , 两 点 高 斯 求 积
18、 公 式 为2 ( ) ( 1 3 ) (1 3 )G f f f ,阶数为 3 另外,高斯求积公式的节点也可以由正交多项式得到若 p 是 n 次多项式,且满足( ) 0 , 0 , , 1 ,b ka p x x d x k n 则 p 与 ,ab 区间上所有次数小于 n 的多项式正交,容易证明: 1 p 的 n 个零点都是实的、单的,且位于开区间 (, )ab 2 区间 ,ab 上以 p 的零点为节点的 n 点插值型求积公式的阶数为 2n-1,是唯一的 n点高斯公式 定义 2 212 如果 1n 个节点的求积公式 0( ) ( ) ( )nbkka kx f x d x A f x ( 2
19、 14) 的代数精度达到 21n ,则称式( 2 14)为 高斯型求积公式 ,此时称节点 kx 为高斯点,系数 kA 称为高斯系数 定理 2 512 以 01, , , nx x x 为高斯点的插值型求积公式具有 21n 次代数精确度的充要条件是以这些节点为零点的多项式 1 0 1( ) ( ) ( ) ( )nnx x x x x x x 与任意次数不超过 n 的多项式 ()px 带权 ()x 均在区间 ,ab 上正交,即 1( ) ( ) ( ) 0b na x p x x dx ( 2 14) 定理 2 6 高斯公式 0( ) ( )nbiia if x dx A f x ( 2 15)
20、 的求积系数 kA 全为正,且 2( ) ( ) , 0 , 1 , ,bbk k kaaA l x d x l x d x k n ( 2 16) 定理 2 7 对于高斯公式( 2 14),其余项为 ( 2 2 ) 211( ) ( ) ( ) ( )( 2 2 ) ! bn naR f f x x d xn , ( 2 17) 其中 1 0 1, , ( ) ( ) ( ) ( ) .nna b x x x x x x x 2 4 2 高斯 勒让德( Gauss-Legendre)公式 13 对于任意求积区间 ,ab ,通过变换 22a b b axt,可化为区间 1,1 ,这时 11(
21、) ( )2 2 2ba b a a b b af x d x f t d t 因此,不失一般性,可取 1, 1ab ,考查区间 1,1 上的高斯公式 11 0( ) ( )niiif x dx A f x ( 4 5) 我们知道,勒让德多项式 1 211 111( ) ( 1 )2 ( 1 ) !n nn nn dL x xn d x , ( 4 6) 是区间 1,1 上的正交多项式,因此, 1()nLx 的 n+1 个零点就是高斯公式( 4 5)的 n+1个节点特别地,称 1()nLx 的零点为高斯点,形如( 4 5)的高斯公式称为高斯 勒让德公式 以上这些公式中的节点和求积系数可查表得到
22、 2 4 3 高斯 哈 米特求积公式( Gauss-Hermite) 14 Gauss-Hermite 求积公式 2 ()0( ) ( )nxnkkke f x d x f x , ( 4 7) 其余项为 ( 2 2 )1( 1 ) ! ( ) .2 ( 2 2) ! nn nnR f fn ( 4 8) 2 4 4 高斯 切比雪夫( Gauss-Chebyshev)求积公式 15 区间为 1,1 ,权函数21() 1x x 的 Gauss 型求积公式,其节点 kx 是 Chebyshev 多项式 1()nTx 的零点,即 21c os ( 0 , 1 , , )2( 1 )k kx k nn
23、 ,而 ( 0 ,1, , )1kA k nn 于是得到 121 0( ) 2 1c o s1 2 ( 1 )1 nkf x kd x fnnx ( 4 9) 称为 Gauss-Chebyshev 求积公式,公式的余项为 ( 2 2 )2 ( 1 )2( ) ( ) , ( 1 , 1 )2 ( 2 2 ) ! nn nR f fn , ( 4 10) 这种求积公式可用于计算奇异积分 2 5 递推型高斯求积 10 高斯求积公式不具有递推性:当节点个数一定时,如果自由选择所有的节点和权以达到最高的阶数,则节点个数不同的公式一般没有公共节点,这意味着与一组节点对应的积分值,在用另外一组节点计算积分
24、值时不能被利用 Kronrod 求积公式避免了这种工作量的增加,这类公式是对称的, n 点高斯公式 nG 与2n+1 个点 Kronrod 公式 21nK 对应 21nK 节点的约束条件为:以 nG 的节点作为 21nK 的节点,按求积公式达到最高阶数的要求确定 21nK 中剩下的 n+1 个节点及 2n+1 个权(其中包括 nG 的节点的权)这样,求积公式的阶数可达到 3n+1,而真正 2n+1 个点高斯公式应该是4n+1 阶的,所以精度和效率是一对矛盾 使用两个节点个数不同的求积公式的主要原因是可以用它们的差估计积分近似值的误差使用 Gauss-Kronrod 公式对时,若以 21nK 的值作为积分的近似值,则一半基于理论,一半基于经验,可以得到关于误差的保守估计: 1.521(20 0 )nnGK Gauss-Kronrod 公式不仅有效地提供了较高的精度,还给出了可靠误差估计,所以它被认为是最有效的求积公 式之一,并且构成了主要软件库中求积程序的基础,特别地,公式