数值分析63 复化求积公式龙贝格求积公式讲解
数值分析第五版李庆扬王能超课件第3章(2)
复化求积公式
h2 h2 6 上例中若要求 | I Tn | 10 ，则 | Rn [ f ] | | f (1) f (0) | 10 12 6
6
h 0.00244949 即：取 n = 409
通常采取将区间不断对分的方法，即取 n = 2k 可用来判断迭代 上例中2k 409 k = 9 时，T512 = 3.14159202 是否停止。 2 1 h 1 注意到区间再次对分时 R2 n [ f ] [ f (b) f (a )] Rn [ f ]
Romberg
T1 = T0( 0 )
<?
算法： T4 = T0( 2 )
T8 = T0
(3)
T2 = T0( 1 )
S1 = T1( 0 )
S2 = T1 S4 = T1
(1) (2)
<?
C1 = T2 C2 = T2
(0) (1)
<?
………………
R1 = T3
第二讲
§1. 复化求积公式
§2. 龙贝格求积公式
高次插值有Runge 现象，故采用分段低 次插值
分段低次合成的 Newton-Cotes 复 合求积公式。
§ 1. 复化求积公式 § 1.拉格朗日插值
2.1 复化梯形公式 1.1 拉格朗日插值
1.2 复化辛普森公式
1.1 复化梯形公式
ba 复合梯形公式： h , xk a k h n
4T2 n Tn 4 1 T2 n Tn 来计算 I 效果是否好些？ 4 1 3 3 Romberg 序列
4 1 T8 T4 = 3.141592502 = S 4 3 3 4T2 n Tn 42 S2n Sn Sn 一般有： Cn 2 41 4 1
龙贝格求 积分
龙贝格(Romberg ）求积法1.算法理论Romberg 求积方法是以复化梯形公式为基础，应用Richardson 外推法导出的数值求积方法。
由复化梯形公式 )]()(2)([2222b f h a f a f h T +++=可以化为)]()]()([2[212112h a f h b f a f hT +++=＝)]([21211h a f h T ++一般地，把区间［a,b ］逐次分半k －1次，（k ＝1，2,……,n)区间长度（步长)为kk m a b h -=，其中mk ＝2k －1。
记k T ＝)1(k T由)1(k T ＝]))12(([21211)1(1∑=---++km j k k k h j a f h T 从而⎰badxx f )(＝)1(kT－)(''122k f h a b ξ- (1)按Richardson 外推思想，可将（1）看成关于k h ，误差为)(2k h O 的一个近似公式，因而，复化梯形公式的误差公式为⎰badxx f )(－)1(k T ＝......4221++kkh K h K ＝∑∞=12i i k i h K （2）取1+k h ＝k h 21有 ⎰ba dx x f )(－)1(1+k T ＝∑∞=+121221i ik ii hK （3）误差为)(2jh O 的误差公式 )(j kT＝)1(-j kT＋141)1(1)1(------j j k j k T T2。
误差及收敛性分析（1）误差，对复化梯形公式误差估计时，是估计出每个子区间上的误差,然后将n 个子区间上的误差相加作为整个积分区间上的误差。
(2）收敛性,记h x i =∆，由于∑=++=ni i i n x f x f h f T 01))]()([2)(＝))()((21101∑∑-==∆+∆n i ni i i i i x x f x x f上面两个累加式都是积分和，由于)(x f 在区间],[b a 上可积可知，只要],[b a 的分划的最大子区间的长度0→λ时,也即∞→n 时,它们的极限都等于积分值)(f I 。
4.4龙贝格求积公式
4 1 T1 ( k − 1) = T0 ( k ) − T0 ( k − 1) 3 3 16 1 T2 ( k − 1) = T1 ( k ) − T1 ( k − 1) 15 15 64 1 T3 ( k − 1) = T2 ( k ) − T2 ( k − 1) 63 63
k = 1 ,2 ,L
因此有如下递推公式 b−a [ f ( a ) + f (b )] T0 (0) = 2
1 T0 (k ) = T0 (k − 1) + hk 2
2 k −1 −1 j =0
∑ f (a + (2 j + 1)h )
k
k = 1, 2 ,L
上式称为递推的梯形公式
思考
递推梯形公式加上一个控制精度,即 可成为自动选取步长的复合梯形公式
带权）正交。 不大于n 不大于 的多项式 P(x) （带权）正交。
k =0
n
x 证明： 证明： “⇒”0 … xn 为 Gauss 点, 则公式 ∫a ρ( x) f ( x)dx ≈ ∑Ak f ( xk ) k=0 次代数精度。 至少有 2n+1 次代数精度。 对任意次数不大于 的多项式 不大于n 对任意次数不大于 ⇔ 求w(x) Pm(x)， Pm(x) w(x)的次数 ， 的次数 求 Gauss 点 不大于2n+1，则代入公式应精确成立： 精确成立： 不大于 ，则代入公式应精确成立 n 0 b ∫ ρ ( x ) Pm ( x ) w ( x )dx = ∑ Ak Pm ( x k ) w ( x k ) = 0
外推加速公式
由复合梯形公式的余项公式
I − T2 n 1 ≈ I − Tn 4 1 I − T2 n ≈ (T2 n − Tn ) 3 4 1 I ≈ T2 n − Tn 3 3 1 f ( x 1 )) − Tn j+ 3 2
龙贝格算法
龙贝格积分1. 算法原理采用复化求积公式计算时，为使截断误差不超过ε，需要估计被积函数高阶导数的最大值，从而确定把积分区间[]b a ,分成等长子区间的个数n 。
首先在整个区间[]b a ,上应用梯形公式，算出积分近似值T1；然后将[]b a ,分半，对 应用复化梯形公式算出T2；再将每个小区间分半，一般地，每次总是在前一次的基础上再将小区间分半，然后利用递推公式进行计算，直至相邻两个值之差小于允许误差为止。
实际计算中，常用ε≤-n n T T 2作为判别计算终止的条件。
若满足，则取n T f I 2][≈；否则将区间再分半进行计算，知道满足精度要求为止。
又经过推导可知，∑=-++=ni i i n n x x f h T T 112)2(221，在实际计算中，取kn 2=，则k a b h 2-=，112)1*2(2++--+=+k i i ab i a x x 。
所以，上式可以写为∑=++--+-+=+kk i k k ab i a f a b T T 211122)2)12((2211k开始计算时，取())()(21b f a f ab T +-=龙贝格算法是由递推算法得来的。
由梯形公式得出辛普森公式得出柯特斯公式最后得到龙贝格公式。
根据梯形法的误差公式，积分值n T 的截断误差大致与2h 成正比，因此步长减半后误差将减至四分之一，即有21114n n T T -≈-将上式移项整理，知2211()3n n n T T T -≈-由此可见，只要二分前后两个积分值n T 和2n T 相当接近，就可以保证计算保证结果计算结果2n T 的误差很小，这种直接用计算结果来估计误差的方法称作误差的事后估计法。
按上式，积分值2n T 的误差大致等于21()3n n T T -，如果用这个误差值作为2n T 的一种补偿，可以期望，所得的()222141333n n n n n T T T T T T =+-=-应当是更好的结果。
数值分析63 复化求积公式龙贝格求积公式讲解
只增加了一个分点
1 xk?1/ 2 ? 2 ( xk ? xk?1)
设hn=(b? a)/n, xk=a+kh n (k=0,1,? ,n),在[xk, xk+1] 上用梯形公式得
T1 ?
hn 2
?f
(
xk
)
?
f ? ( xk ? 1 )
复化求积的基本想法 :
将积分区间 [a, b]n等分, 步长
h?
b
? n
a
,
分点为
xk=a+kh (k=0,1,…,n) , 则由定积分性质知
? ? ? I ?
b
n?1
f ( x )dx ?
xk?1 f ( x )d x
a
k ? 0 xk
每个子区间 上的积分
?xk?1 f ( x )dx xk
用低阶求积公式 , 然后把所有区间的 计算结果求和 ,
注2: 同样也可用 | S4m-S2m |<ε 来控制计算的精度 . 这就是下面要介绍的 龙贝格求 积公式 .
6.4 龙贝格求积公式
6.4.1 梯形公式的递推化
复化求积方法可提高求积精度，实际计算时若
精度不够可将步长逐次分半 . 设将区间 [a, b]分为n等
分，共有 n+1个分点，如果将求积区间再分一次，则 分点增至 2n+1个，我们将二分 前后两个积分值 联系
果T8=0.9456909 只有2位有效数字，而应用复化辛普 森公式计算的结果 S4= 0.9460832 却有6位有效数字 .
注：为了利用余项公式估计误差，要求 f(x)=sin x/x 的高阶导数，由于
龙贝格(Romberg)求积法
龙贝格（Romberg）求积法
复化求积方法对于提高计算精度是行之有效的 方法，但复化公式的一个主要缺点在于要先估计出 步长。若步长太大，则难以保证计算精度，若步长 太小，则计算量太大，并且积累误差也会增大。在 实际计算中通常采用变步长的方法，即把步长逐次 分半，直至达到某种精度为止。
1.1变步长的梯形公式 变步长复化求积法的基本思想是在求积过程中，
)
0.9397933
进一步二分求积区间,并计算新分点上的函数值
f
(
1 4
)
0.9896158
,
f
(
3 4
)
0.9088516
有
T4
1 2 T2
1 4
f
(
1 4
)
f
(
3 4
)
0.9445135
这样不断二分下去，计算结果如P139列表所 示。积分的准确值为0.9460831，从表中可
看出用变步长二分10次可得此结果。
的考积察分T值与nS等n 。份辛卜生公式 S n之间的关系。将
复化梯形公式
Tn
h 2
f
(a)
2
n1 k 1
f
(xk ) f (b)
梯形变步长公式
T2 n
Tn 2
h n1
2 k 0
f (xk1 ) 2
代入(6.9) T 表达式得
h
n1
n1
T
6 f (a) 4k0
f
(
x k
1
)
2
（
2
输入 a,b,ε
）
变
b-ah,
h 2
f(a)+f(b)T1
6b复合求积公式龙贝格算法
步长折半：[xi , xi+1/2] ， [xi +1/2 , xi+1]
n1
xi xi +1/2 xi +1
h T2 n f ( xi ) f ( xi 1 2 ) f ( xi 1 2 ) f ( xi 1 ) i 0 4 n1 h f ( xi ) 2 f ( xi 1 2 ) f ( xi 1 ) i 0 4 h n1 h n1 1 h n1 f ( xi ) f ( xi 1 ) f ( xi 1 2 ) Tn f ( xi 1 2 ) 4 i 0 2 i 0 2 2 i 0 13
1 I T2 n (T2 n Tn ) 3
I Tn 4( I T2n )
3I 4T2n Tn
1 3
验后误差估计式 I T2 n (T2 n Tn )
当
T2n Tn 时，T2n即为所求的近似值。
1 (T2 n Tn ) 3
是T2n 的修正项，它与T2n 之和比T2n 、 Tn更接近与真值，即它是一种补偿。
|| T T2-T|< 2-T1|<
输出T2
16
举例
计算精度满足 | T2n Tn | 107
I [ f ]=0.946083070367
例：用梯形法的递推公式计算定积分 解：
1
0
sin( x ) dx , 要求 x
k
0 1 2 3 4 5 6 7 8 9 10
T (k)
梯形法递推公式
1 h n1 1 h n1 T2 n Tn f ( xi 1 2 ) Tn f ( a ih 0.5h) 2 2 i 0 2 2 i 0
第五讲 复化求积公式
四、自动选取积分步长
事前确定步长的问题 (1) 高阶导数的估计往往是很困难的； (2) 这种估计往往是很保守的，得到的n往往偏大。 为了改正上述缺点，实际常采用“事后估计法” “事后估计法”的基本思想是 (1) 求数值积分时，将区间逐次分半； (2) 利用前后两次的计算结果来判断误差是否满足精度要求， 从而确定n. 下面以复化梯形公式为例来介绍这种步长逐次减半求积法
1 h n1 T f (x ), n k1 2 2 2k0
如何根据Tn和T2n来确定误差是否满足要求？
(ba ) 2 I Tn ( h f ) 1 2 ba h 2 I T2n ( ( ) f ) 1 2 2
则有
如果二阶导数在区 间[a,b]上变化不大
n 1
R (Tn )
复化simpson公式的截断误差
（ 4 ） 若 函 数 f ( x )[ 在 a ,] b 上 连 续 ， 则
ba 4 (4) h5 (4) hf ( ) R ( S n ) f ( k ) I Sn 2 8 8 0 8 8 0 k0 2
0 . 9 4 6 0 8 3 2
1 1 C2 [ 7 f( 0 ) [ 3 2 f( x 1) 1 2 f( x 2) 3 2 f( x 3) ] k k k 1 8 0 k 0 4 4 4
1 4 f( x 7 f( 1 ) ] k)
k 1
1
0 . 9 4 6 0 8 3 0
n 1 h [ 7 f ( x ) 3 2 f ( x ) 1 2 f ( x ) 3 2 f ( x ) 7 f ( x ) ] k 1 2 3 k 1 k k k 9 0 k 0 4 4 4
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
将积分区间 [a, b] 划分为2n等分, 即将每一个区间
[xk, xk+1]经过二等分增加了一个分点
x k ?1/2
?
1 2
(
x
k
?
x k?1)
在每个子区间 [xk, xk+1]上的积分用 辛普森公式 , 得
?x k ? 1 xk
f (x)dx ?
h 6
[
f
(
xk
)
?
4
f
(
xk
? 1/2
)
?
f (xk?1)]
?
(b ? a)5 2880n4
f
(4) (? )
?6.3.2 复化辛普森公式
将积分区间 [a, b] 划分为2n等分, 则
? ? ? I ?
b
n?1
f ( x )dx ?
x2k? 2 f ( x )d x
a
k ? 0 x2k
每个子区间 [x2k, x2k+2]上的积分用 辛普森公式 , 得
?x2 k ? 2 x2k
称为复化辛普森公式 . 记
? ? h
n?1
n?1
Sn ?
[f 6
(a) ?
4
k?0
f
( x k ? 1/2 ) ?
2
k ?1
f
(xk ) ?
f (b)]
若 f(x)? C 4[a,b], 其求积余项 为 h ? b ? a
n
Rn ( f ?
b ? a ( h )4 180 2
f (4) (? ) ?
f (b)]
称为复化梯形公式 . 记
? h
n?1
Tn
?
[ f (a) ? 2
2
k ?1
f (xk ) ?
f (b)]
若 f(x)? C2[a,b], 其求积余项 Rn(f )为(p239)
? ? Rn ( f
)?
I
? Tn
?
n?1
[?
k?0
h3 12
f
?(?? k )] ?
?
h2 ?b ? a n?1 12 n k ? 0
为了提高精度，通常在实际应用中往往采用 将积 分区间划分成若干个小区间，在各小区间上采用低次 的求积公式 (梯形公式或抛物形公式 )，然后再利用积 分的可加性，把各区间上的积分加起来，便得到新的 求积公式，这就是复化求积公式的基本思想 . 我们仅 讨论各小区间均采用同一低次的求积公式的复化求积 公式.
f
( x2k
)
?
4
f
(
x2k?1 )
?
f ( x 2k? 2 )]
? ? h
n?1
n?1
?
[ f (a) ? 3
2
k?1
f (x2k ) ? 4
k?0
f ( x 2k?1 ) ? f (b)]
称为复化辛普森公式 . 记
? ? h
n?1
n?1
Sn
?
[ f (a) ? 3
2
k ?1
f (x2k ) ? 4
f ?(??
k)
Rn (
f
)
?
?
b? a 12
h2
f
??(?
),
? ? [a, b], h ? b ? a .
n
可以看出误差是 h2 阶，且由误差公式得到，当
f(x)? C2[a, b] 时，则有
b
? lim
n? ?
Tn
?
f ( x )dx .
a
即复化梯形公式是 收敛的. 事实上只要 f(x)? C[a, b], 则可得到收敛，因为只要把 Tn改写为
k?0
f (x2k?1) ? f
(b)]
若 f(x)? C 4[a,b], 其求积余项 为 h ? b ? a 2n
Rn( f
)?
I
? Sn
?
?
b ? a ( 2h )4 180 2
f (4) (? ) ?
?
(b ? a)5 2880n4
f (4) (? )
例1 对于函数 f(x)=sin x/x，给出n=8的函数表，
f (x)dx
?
2h 6
[
f
(
x2k
)
?
4
f
(
x2k
?
1
)
?
f ( x2k?2 )]
? I
?
2h 6
n?1
[
k?0
f
( x2k
)
?
4
f
( x2k ?1 )
?
f ( x2k ? 2 )]
? ?
h 3
n?1
[
k?0
f
(
x2k
)
?
4
f
(
x2k
?1)
?
f ( x2k ? 2 )]
? I
?
h 3
n?1
[
k?0
?
h 2
[
f
(
x
k
)
?
f ( x k ?1 )]
? ? ? n?1
I?
k?0
xk ?1 xk
f ( x )dx
?
n?1 k?0
h[
2
f
(
xk
)
?
f ( x k ?1 )]
? ? n?1 h
h
n?1
I
?
k?0 2[f (xk ) ?
f (xk?1)] ?
[ f (a) ? 2
2
k?1
f (xk) ?
就得到整个区间上积分 I 的近似值。
6.3.1 复化梯形公式
将积分区间 [a, b]划分为 n 等分, 则
? ? ? I ?
b
n?1
f ( x )dx ?
xk?1 f ( x )dx
a
k ? 0 xk
每个子区间 [xk, xk+1]上的积分用 梯形公式 , 得
?x k ? 1 xk
f ( x )dx
? ? ? I ?
b
n?1
f ( x )dx ?
xk ?1 f ( x )dx
a
k?0 xk
?h n?1
I
?
6
[
k?0
f
(
xk
)
?
4
f
( x k ? 1/2
)
?
f ( x k?1 )]
? ? h
n?1
n?1
?
[f 6
(a) ?
4
k?0
f
( x k ? 1/2 ) ?
2
k?1
f
(xk ) ?
f (b)]
6.3 复化求积公式
从求积公式的余项的讨论中我们看到，被积函数 所用的插值多项式次数越高，对函数光滑性的要求也 越高.另一方面，插值节点的增多 (n的增大)，在使用 牛顿-柯特斯公式时将导致求积系数出现负数 (当n≥8 时, 牛顿-柯特斯求积系数会出现负数 )，即牛顿 -柯特 斯公式是不稳定的， 不可能通过提高阶的方法来提高 求积精度 .
试用复化梯形公式和复化辛普森公式 计算积分
1 sin x
I ? ?0 x dx .
解 将积分区间 [0,1]划分为8 等分，用复化梯形公式求得
x
f(x)
0
1
1/8 0.9973978
1/4 0.9896158
T8 ? 0.9456909. 而将积分区间 [0, 1]划分为2×4等 分，用复化辛普森公式求得
复化求积的基本想法 :
将积分区间 [a, b]n等分, 步长
h?
b
? n
a
,
分点为
xk=a+kh (k=0,1,…,n) , 则由定积分性质知
? ? ? I ?
b
n?1
f ( x )dx ?
xk?1 f ( x )d x
a
k ? 0 xk
每个子区间 上的积分
?xk?1 f ( x )dx xk
用低阶求积公式 , 然后把所有区间的 计算结果求和 ,
? ? Tn
?
1 ?b ? a n?1
2
? ?
n
k?0
f (xk )?
b? a n
n k ?1
? f ( x k )?.
?
当n →∞时，上式右端括号内的两个和式均收敛到函数
的积分，所以复化梯形公式收敛 . 此外，Tn 的求积系 数均为正，由定理 2知复化梯形公式是稳定的 .
6.3.2 复化辛普森公式