龙贝格求积法
龙贝格求 积分
龙贝格(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 。
龙贝格求积公式
2
f( ) 2
T3 , 1 T3 , 2
T3 , 3
T4 , 1 T5 , 1
M
T4 , 2 T5 , 2
M
T4, 3 T5 , 3
M
T4 , 4 T5 , 4
M
用龙贝格方法计算积分的步骤为：
（1）：准备初值，先用梯形公式计算积分近
似 （值 2）：：T1按 变b 2步a[长f (a梯) 形f (公b)]式计算积分近似值：
4 0.918741799 0.916327874 0.916297224 0.916294351
5 0.916905342 0.916293190 0.916290077 0.916290776
T5 , 4 0.916290776
例3：取=0.00001,用龙贝格方法计算积分
1
I
4
dx
01 x2
解：由题意
f(x)=4/(1+x2) a=0 b=1 f(0)=4 f(1)=2 由梯形公式得
T1=1/2[f(0)+f(1)]=3 计算f(1/2)=16/5 用变步长梯形公式得
T2=1/2[T1+f(1/2)]=3.1 由加速公式得
S1=1/3(4T2-T1)=3.133333333
求出f(1/4) f(3/4) 进而求得 T4=1/2{T2+1/2[f(1/4)+f(3/4)]}
Simpson加速公式： Cn
42 S2n Sn 42 1
I
C2n
1 43 1(C2n
Cn )
43C2n Cn 43 1
Cotes加速公式：Rn
43C2n Cn 43 1
Romberg 值序列
计算方法-4.6-4.7龙贝格、高斯求积公式
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 ,
引入 T1 (k 1),令
S n S 2 k 1
即
Sn S2 k 1 T1 ( k 1)
S 2 n T1 ( k )
--------(6)
当然
2016/8/14
10
因此由复合Simpson公式的余项
I S2 n 1 ( S2 n Sn ) 15
Cn
可得 令 即
1次,3次和5次
复合梯形、复合Simpson、复合Cotes公式的收敛阶分别为
2阶、4阶和6阶 无论从代数精度还是收敛速度,复合梯形公式都是较差的 有没有办法改善梯形公式呢?
2016/8/14
2
一、复合梯形公式的递推化
将定积分I f ( x)dx 的积分区间 [a , b]分割为n等份 a ba x a jh , j 0 , 1 , , n 各节点为 h j n 复合梯形(Trapz)公式为
n 1 ba Tn [ f ( a) 2 f ( x j ) f (b)] 2n j 1
b
--------(1)
如果将 [a, b]分割为2n等份，而h (b a) /n不变, 则
n 1 n 1 ba T2 n [ f ( a ) 2 f ( x j ) 2 f ( x 1 ) f (b)] j 4n j 1 j 0 2 --------(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 =+-=-应当是更好的结果。
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 ,
1 h n1 Tn f ( x 1 ) j 2 2 j 0 2 1 b a n 1 1 Tn f (a ( j )h) 2 2n j 0 2 1 b a n 1 ba Tn f ( a ( 2 j 1) ) 2 2 n j 0 2n
Romberg算法的代 数精度为m的两倍 Romberg算法的收敛 阶高达m+1的两倍
T0 (0 ) T0 (1) T0 ( 2 ) T0 ( 3) T1 (0 ) T1 (1) T1 ( 2 ) T2 (0 ) T2 (1)
T3 (0 )
龙贝格积分
例：计算
I T2n 1 I Tn 4
令
--------(7)
即 当然
Cn C2 k 1 T2 (k 1) C 2 n T2 (k )
--------(8)
同样由复合Cotes公式的余项
I C2 n 1 (C2 n Cn ) 63
64 1 64 1 得 I C2 n Cn T2 ( k ) T2 ( k 1) 63 63 63 63
43 C 2n C n Rn 3 4 1
Romberg
T1 = T0( 0 )
T2 = T0( 1 )
<?
算法： T4 = T0( 2 )
T8 = T0
(3)
S1 = T1( 0 ) S2 = T1
数值分析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 Integration Algorithm）是用于数
值积分的一种高效的迭代方法。
它通过连续的二分、四分、八分等等
区间的方式，逐渐逼近最终的积分值，从而提高计算的精度。
该算法的基本思想是利用Richardson外推技术，结合复合梯形
法则，逐渐缩小区间并增加采样点数，以得到更精确的积分值。
下面
我们来介绍龙贝格求积算法的步骤：
1. 将积分区间[a, b]进行二分，得到初始的两个子区间；
2. 对每个子区间应用复合梯形公式进行数值积分，可以得到初始的近似积分值；
3. 利用Richardson外推技术，对不同精度的积分值进行线性组合，得到更高精度的积分值；
4. 重复步骤2和3，将积分区间不断地二分，并逐步增加采样点数，直到达到所需的精度要求。
龙贝格求积算法的主要优点是在保持高精度的能够有效减少计算量。
该算法还可以通过预先计算一些常见函数在一些固定的点上的值，以进一步提高计算速度。
总结起来，龙贝格求积算法通过利用复合梯形法则和Richardson 外推技术，逐渐逼近积分值的精确结果。
它是一种高效且精确的数值
积分方法，广泛应用于科学计算和工程领域。
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
计算方法--龙贝格积分法
结论:从运行结果可见， 当设定误差 eps=1e-5 时， 为了满足误差要求， 运行 n=14 时才能满足， 这样精度比较高，但是代价就是运行的时候速度有点慢。
按照题目要求， 要求尽可能的使结果精确一些， 因此我这里采用的是根据误差要求进行的计 算，计算程序如下： function [R,T,n]=Romberg3(myfunc,a,b,eps) %This function uses Romberg to calculate integration %myfunc is the integrand,a is integral lower bound %b is integral upperbound %eps is tolerance format long; n=5; while 1 T(n,:)=0; T(1,1)=(b-a)*(myfunc(a)+myfunc(b))/2; for k=1:n-1 sum=0; for i=1:2^k sum=sum+myfunc(a+(2*i-1)*(b-a)/(2^(k+1))); end T(k+1,1)=0.5*T(k,1)+(b-a)*sum/(2^(k+1)); end for i=2:n T(i,2)=T(i,1)+1/3*(T(i,1)-T(i-1,1)); end for i=3:n T(i,3)=T(i,2)+1/15*(T(i,2)-T(i-1,2)); end for i=4:n T(i,4)=T(i,3)+1/(4^3-1)*(T(i,3)-T(i-1,3)); end if abs(T(n,4)-T(n-1,4))<eps break; else n=n+1; end end R=T(n,4); 以上便是这次编写的龙贝格积分的程序。 在matlab的command window运行如下： >> myfuc=inline('log(1+x)/(1+x^2)','x'); >> eps=1e-5;
龙贝格求积法
h p
)
Fk
q pk 1
(h)
,k
0,1, 2,L
定义的序列{Fk(h)}有
Fk (h) F (0) an(n)1h pn1
a h (n) pn2 n2
a h (n) pn3 n3
L
,
其中an(n)k (k 1, 2, 3,L )与h 无关，q>1.
Richardson外推法应用非常广泛且有效，下面介绍应用 于数值积分的情形。
16 15
S4
1 15
S2
3.1415946
64 1 R1 63 C2 63 C1 3.141586292
1 sin x
算例结2 果用见R表om4-b5(ekr代g算表法二计分算次I数)。0 计x算值dx.的得误到差的不梯超形过值，计
0.510-6.
表4-5
k
T2k
S k 1 2
§4.4 外推原理与Romberg求积方法
4.4.1 外推原理
在科学与工程计算中，很多算法与步长h有关，特别是数值 积分、数值微分和微分方程数值解的问题。对于这些算法，我 们可以通过外推技巧提高计算精度。
例1 计算的近似值。
数值分析6.3 复化求积公式、龙贝格求积公式讲解
精度不够可将步长逐次分半. 设将区间 [a, b]分为n等
分，共有n+1个分点，如果将求积区间再分一次，则
分点增至2n+1个，我们将二分前后两个积分值联系
起来加以考虑. 注意到每个子区间[xk, xk+1]经过二分
只增加了一个分点
x k 1/ 2
1 ( x k xk 1 ) 2
设hn=(ba)/n, xk=a+khn (k=0,1,,n),在[xk, xk+1]
I f ( x )dx
b a k 0 n 1 xk 1 xk
f ( x )dx
每个子区间[xk, xk+1]上的积分用梯形公式, 得
xk 1 xk
h f ( x )dx [ f ( xk ) f ( xk 1 )] 2
xk 1 xk
I
k 0
6.3 复化求积公式
从求积公式的余项的讨论中我们看到，被积函数
所用的插值多项式次数越高，对函数光滑性的要求也
越高.另一方面，插值节点的增多(n的增大)，在使用
牛顿-柯特斯公式时将导致求积系数出现负数(当n≥8
时, 牛顿-柯特斯求积系数会出现负数)，即牛顿-柯特
斯公式是不稳定的，不可能通过提高阶的方法来提高 求积精度.
b n 1 xk 1 xk a
I f ( x )dx
k 0
f ( x )dx
h n 1 I [ f ( xk ) 4 f ( xk 1/2 ) f ( xk 1 )] 6 k 0
n 1 n 1 h [ f (a ) 4 f ( xk 1/2 ) 2 f ( xk ) f ( b)] 6 k 0 k 1
(龙贝格求积)
针对这类问题的算法技巧是在不同区间上预测被积函数 变化的剧烈程度确定相应的步长.
这种方法称为自适应积分方法.
设给定精度要求 0，计算积分
b
I a f (x)dx
的近似值.先取步长 h b a ，应用辛普森公式有
I b f (x)dx S(a, b) b a ( h )4 f (4) (),
注意到每个子区间等分分为设将区间一次如果将求积区间再二分1个分点共有子区间上的积分值为用复化梯形公式求得该分点经过二分只增加了一个注意到不难导出下列递推公式逐次分半算法据梯形公式计算得9588510求出将区间二等分进一步二分求积区间9445135计算积分值8414709例题1精确至三位有效数字例题25695注意到57例如
例如：对区间[a, b], 令h b a构造梯形值序列{T2k }
h
T1
[ 2
f (a)
f
(b)]
逐次分半算法
把区间二等分，每个小区间长度为 h/2=(b-a）/2，
于是
T2 =T1/2+[h/2]f(a+h/2） 把区间四(22)等分，每个小区间长度为h/2 2 =（b-a）
/4，于是
T4 =T2/2+[h/22][f(a+h/4)+f(a+3h/4)]
················
把[a,b] 2k 等分，分点xi=a+(b-a)/ 2k ·i （i =0， 1，2 ···2k）每个小区间长度为(b-a)/ 2k ，由归
纳法可得
T T 1
2k 2
2 k1
ba 2k
注意到每个子区间
[ xk , xk1 ]经过二分只增加了一个分点
《龙贝格求积》课件
龙贝格求积法的历史背景
01
龙贝格求积法是由瑞典数学 家龙贝格在19世纪末提出的
。
02
它的出现为数值积分的发展 奠定了基础，成为数学领域
的一项重要成果。
03
随着计算机技术的发展，龙 贝格求积法在科学计算、工 程技术和数据分析等领域得
缺点
计算量大
对于大规模问题，龙贝格求积法可 能需要较大的计算资源和时间，因 为需要进行多次迭代和数值逼近。
对初值敏感
该方法对初值的选择比较敏感，如 果初值选择不当，可能会导致算法
不收敛或者收敛到非期望的解。
需要选择合适的参数
龙贝格求积法的精度和稳定性与参 数的选择密切相关，需要仔细选择 合适的参数值。
通过差商的方式计算插值多项式的导数。
导数在数值分析中的应用
用于估计函数的局部变化趋势，提高数值计算的精度。
龙贝格求积公式的推导
03
龙贝格求积公式的定义
龙贝格求积公式的推导过程
龙贝格求积公式的应用
利用插值多项式及其导数构造求积公式， 用于数值积分。
通过构造插值多项式和其导数，利用牛顿 -莱布尼茨公式推导得出。
对未来研究的展望
01
入 ，未来可以对龙贝格求积法进 行进一步的改进和优化，以提 高其计算效率和适用范围。
02
与其他方法的比较研究
可以开展龙贝格求积法与其他 积分方法的比较研究，深入探 讨各种方法的优缺点和应用场 景，为实际应用提供更加全面
的理论支持。
03
扩展应用领域
对不规则区域处理困难
对于不规则积分区域，龙贝格求积 法可能需要进行额外的处理和调整 ，增加了计算的复杂度。
龙贝格求 积分
龙贝格（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++k k h 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 。
数值计算方法 龙贝格求积公式 - 龙贝格求积公式
类似地可验证：
3
S2
4
4
1
T4
4
1
1
T2
龙
贝 格 算 法
Sn
4
4
1
T2
n
4
1
1
Tn
即
Sn
4 3 T2n
1 3
Tn
Sn __ 辛 普 森 积 分 值
注 意 复 化 辛 普 森 求 积 公式 的 余 项
3
Rn
I
Sn
ba 180
h 2
4
f
(4) ( )
O(h4
)
龙 贝 格
R2 n
I
S2n
T4
1 2
T2
1 4
[
f
(
1 4
)
f
(
3 4
)]
0.9445135
典型例题
例2
用复化梯形公式计算积分I 1 dx
0 (1 x ) x
精确至三位有效数字
I
1 dx 0 (1 x )
x t x
1 2dt 0 1 t2
1
g(t )dt
0
T1
1 [g(0) 2
g(1)] 1.5
T2
1 2
[T1
1
g(0.5)] 1.55
T4
1 2
[T2
1 ( g(0.25) 2
g(0.75))] 1.5656
典型例题
T8
1 2 [T4
1 ( g(0.125) 4
g(0.375)
g(0.625)
g(0.875))]
1.5695
注意到:
T8 T4
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
(由上节梯形公式)
h 6
f
n1
a 2
k 1
f
xk
n1
4
k 0
f
x
k
1 2
f
b
Sn
（由上节Simpson公式）
故（4.4.6）式成立.
同理，由上节近似式
I
S2n
1 15
S
2n
Sn
类似推导可得：
C2n
h p
)
Fk
q pk 1
(h)
,k
0,1, 2,L
定义的序列{Fk(h)}有
Fk (h) F (0) an(n)1h pn1
a h (n) pn2 n2
a h (n) pn3 n3
L
,
其中an(n)k (k 1, 2, 3,L )与h 无关，q>1.
Richardson外推法应用非常广泛且有效，下面介绍应用 于数值积分的情形。
)
p
O(hs
),
削去截断误差的主项，得新的算法
q pF( h) F(h)
F1(h)
q qp 1
F (0) O(hs ),
我们称这个计算过程为Richardson外推法。 F1(h)逼近F(0)的截 断误差是 O(hs ).
只要知道F(h)的更加完整的关于h幂的展开式，而无需知 道展开式中各个系数的具体数值，就能重复使用Richardson外 推法，直到截断误差达到容许误差。用归纳法可以证明下面更 一般的定理。
可将上述外推思想推广到一般情况。设F(h)是计算F(0)的一 种近似算式，带截断误差的表达式为
F(h) F(0) aphp O(hs ), s p,
其中，a p 与h无关。如果我们用h和h/q (q>1)两种步长分别计算
F(h)和F(h /q)，则有
F
(
h q
)
F
(0)
ap
(
h q
① T1 ② T2 ④ T4
③ S1 ⑤ S2
⑥ C1
⑨ C2
⑩ R1
⑦ T8
⑾ T16
⑧S4 ⑿ S8
⒀ C4
⒁ R2
Sn
43T2n
1 3 Tn , Cn
16 15
S2n
1 15
Sn , Rn
64 63 C2n
1 63 Cn ,
T1,T2、T4、T6 ……由梯 形公式、复化梯形公式
dx
1
t 2 1
二、Romberg算法
变步长梯形求积法算法简单，但精度较差，收敛速度较慢， 但可以利用梯形法算法简单的优点，形成一个新算法，这就是龙 贝格求积公式。龙贝格公式又称逐次分半加速法。
Romberg算法是在积分区间逐次分半的过程中，对用复化梯 形法产生的近似值进行加权平均，以获得准确度较高的近似值的 一种方法，具有公式简练，使用方便，结果较可靠的优点。
§4.4 外推原理与Romberg求积方法
4.4.1 外推原理
在科学与工程计算中，很多算法与步长h有关，特别是数值 积分、数值微分和微分方程数值解的问题。对于这些算法，我 们可以通过外推技巧提高计算精度。
例1 计算的近似值。
由函数 sin x 的Taylor展开式有
n sin
n
3
3!n2
Tn Sn 即：
Sn
4 3
T2
n
1 3
Tn
4T2n Tn 4 1
(4.4.6)
这表明在收敛缓慢的梯形数列 T2k 的基础上，若将T2n 与 Tn
按(4.4.6)作线性组合就可产生收敛速度较快的Simpson序列：
S2k ：S1 、S2 、S4...
关于(4.4.3.6)的证明：由(4.4.1)可知：
.
可能是更好的结果。
（4.4.4）
1 Tn T2n 3 (T2n Tn )
这就是说，用梯形法二分前后两个积分值 Tn 和T2n作线性组合
即
41 Tn 3 T2n 3 Tn
有可能比
T2n
更好地接近于积分
b
a
f
xdx
的真值
I
(4.4.5)
这种新近似值 Tn实质上又是什么呢？可以验证：
T2
4
1
1 3
T1
3.133333
C1
16 15
S2
1 15
S1
3.142118
S2 3 T4 3 T2 3.141569
T8
1 2 T4
1 8
f
(
1 8
)
f
(
3 8
)
f
(
5 8
)
f
(
7 8
)
3.138989
S4
4 3
T8
1 3
T4
3.141593
C2
C k2 2
R k3 2
O O.9207355 1 O.9397933 2 O.9445135 3 0.9456909
O.9461459 O.9460869 O.9460833
O.9460830 O.9460831
O.9460831
我们看到，这里利用二分3次的数据(它们的精度都很差，只有二三位
是有效数字)，通过三次加速求得 R1 =0.9460831，这个结果的每一位
Romberg序列为止。
2、Tmk 可用二维数组来存放并参加运算，也可用一维数组。
3、 对于积分限为无穷的积分 ，可利用变量代换化成有限区
间的积分然后再进行计算。例如： dx
1 1 x2 x3
令
x1 t
则
dx
1 t2
dt
代入得：
dx
1 1 x2 x3
0
4T2n Tn
3
Tn
2ba n
n k 1
f
a 2k
1
1 2
b
n
a
1 3
1 h
3
2
f
a
n1
2
k 1
f
xk
f
b
n
2h
k 1
f
a
k
1 2
h
值小于事先给定的精度要求。取最后一次的数值积分值作为积分的
近似值。
四、几点说明：
4m
1
1、在上面“加工”过程中的系数4m 1 和4m 1 ，当
m
4
时， 1
4m 1
1 255
而另一个 系数则接近于1，也就是
新公式与原公式差别不大，Tmk Tmk1 但工作量却大增。
因此，在实际计算中常规定 m 3，即计算到出现
16 15
S2n
1 15
Sn
42
S2n Sn 42 1
（4.4.7）
即将Simpson序列 S2k 按(4.4.7) 作线性组合就可产生收敛
速度更快的新序列 C2k : C1,C2,C4,K
------------Cotes序列
由近似式
I
C2n
1 63
(C2
n
Cn )
16 15
S4
1 15
S2
3.1415946
64 1 R1 63 C2 63 C1 3.141586292
1 sin x
算例结2 果用见R表om4-b5(ekr代g算表法二计分算次I数)。0 计x算值dx.的得误到差的不梯超形过值，计
0.510-6.
表4-5
k
T2k
S k 1 2
m=2在Simpson序列中，m=3在Cotes序列中，Tmk 引入上面记号
后，Romberg算法所用到的各个计算公式可统一化为：
T (0) 1
ba
2
f
(a)
f
(b)
T1i
1 2
T1i
1
b 2i
a
2i1 j1
f
(a
2
j 1 2i (b
a) , i
证 T2n 计算结果的误差很小。这样直接用计算结果来估计误差的方法
通常称作误差的事后估计法。
按式(4.4.3),积分近似值与 T2n
的误差大致等于
1 3
T2
n
Tn
因此如果
用这个误差值作为 T2n 的一种补偿，可以期望所得到的
T
T2 n
1 3
(T2n
Tn )
4 3
T2
n
1 3
Tn
f
(
1 2
)
1 3 2
1 16 25
3.1
S1
4 3
T2
1 3
T1
3.133333
T4
1 2
T2
1 4
f