第四章 数值积分
§ 4.1 数值积分的一般概念值型求积公式。的数值积分公式称为插的插值函数！由此得到及函数值上对应于分划是区间的简单函数？？若逼近误差：
易算！数值求积公式：
此时的简单函数数值方法：找一个逼近的原函数不易找到。问题：
)(),(,),(),(
],[)()(
)(
)()()(
)()()(
)()()(
),()(
)(?,)()(
110
110
b
a
nn
nn
b
a
b
a
b
a
xfxfxfxf
bxxxxa
baxpxf
dxxR
pIfIfE r r
dxxppIdxxf
xRxpxf
xpxf
xfdxxffI
次代数精度。称此求积公式具有
）。则不精确成立（次多项式个而对于某都精确成立（的多项式于
）对于任意次数不高公式（所谓代数精度：若求积
。并且代数精度尽量的高我们的任务是希望误差
）（
一般形式为：显然，数值求积公式的
m
pE r rxpm
pE r rxpm
pIfIfE r r
xffI
mm
mm
n
i
ii
0)
~
()(
~
1
),0)()(
1
)()()(
1 )()(
11
0
称为高斯型求积公式。则由此得到的求积公式上的最佳平方逼近函数的简单函数是区间若逼近,],[)( baxf
。时，算法是数值稳定的则当且数值稳定性：若收敛性：
次代数精度。公式具有则求积而若
：数精度可采用如下方法确定一个求积公式的代然联系。的误差小。他们没有必高，也不代表求积公式数精度不代表代数精度高。代）求积公式的误差小，
量指标。是求积公式优良性的度数精度精度的度量标志。而代）求积公式的误差计算注：
n
i
iiiii
n
kk
yxf
nfE r r
xE r rmkxE r r
0
1
,,)()5
.,0)()4
m
,0)(,,,1,0,0)(
)3
2
1
§ 4.2 插值型求积公式
b
a
n
n
n
b
a
ii
n
i
iin
n
n
dxx
n
f
fE r r
dxxlxffIfI
xf
xfxfxxxa
.)(
)!1(
)(
)(:
.)(,)()()(
)(,
),(),(
1
)1(
0
1010
截断误差拉格朗日求积公式：
为逼近函数，得到构成的拉格朗日函数作及对应的型值用分划拉格朗日求积公式一、
但不一定数值稳定。虽然
）代数精度至少为注：
,)2;1
0
ab
n
n
i
i
1
)(
12
)(
)(
))()((
2
)(
1*
5
6)4
,.1
3
3
1
1
0
代数精度为时）梯形公式（
柯特斯求积公式：）几个低次牛顿公式。求积公式。即复合求积段低次插值下的实用的求积公式应是分也难以保证。所以真正果很差，自然求积精度时拉格朗日插值逼近效我们知道数值不稳定。说明对大的但为偶数时代数精度为柯特斯公式，且求积公式称为牛顿）等距节点的拉格朗日
f
ab
fE r r
bfaf
ab
fI
n
n
nnn
n
i
i
3
)(
80
3
)(
3
) ],()2(3)(3)([
8
)(
38/3
3
2
),(
90
)(
)]()
2
(4)([
6
)(
2,*
)4(5
3
3
)4(
5
2
2
代数精度为时）公式（
代数精度为时）辛埔生求积公式（
fhfE r r
ab
h
bfhafhafaf
ab
fI
nS i m p s o n
ab
hf
h
fE r r
bf
ba
faf
ab
fI
nS i m p s o n
二、复化求积公式
1、复化梯形公式
1
1
1
1
1
3
1
1
1
)())()((
2
))()((
2
)(
)(
12
))()((
2
)()(
./)(,:0,],[
1
n
i
n
i
iin
n
i
i
n
i
ii
b
a
n
i
x
x
i
ihafhbfaf
h
xfxf
h
fT
f
h
xfxf
h
dxxfdxxf
nabhniihaxba
i
i
复化梯形求积公式为：
的分划：给定
难的。
的估计是困然而或积分值，需要事先给定的近似满足条件实际计算中，要想得到注：
值稳定的。复合梯形求积公式是数数值稳定性：
收敛性：
。代数精度为截断误差：
)(m a x.
)(
)(
.1:1,,
2
).()(
1
)(
12
)(
)(
12
)(
00
0
2
2
1
3
xfhn
fE r r
abn h
nih
h
hfE r r
f
hab
f
h
fE r r
n
n
i
i
n
i
i
in
n
n
i
in
2、复合辛埔生求积公式
))()(4)((
3
)(
)(
90
))()(4)((
3
)()(
./)(,2,2:0,
],[
1
212)1(22
1
)4(
5
1
212)1(2
1
2
)1(2
m
i
iiim
m
i
i
m
i
iii
b
a
m
i
x
x
i
xfxfxf
h
fS
f
h
xfxfxf
h
dxxfdxxf
nabhmnmiihax
ba
i
i
：复化辛埔生求积公式为的分划：给定
1
11
2
))2((
3
2
))12((
3
4
))()((
3
)(
m
i
m
i
m
hiaf
h
hiaf
h
bfaf
h
fS即
).()(
3
)(
1 80
)(
)(
90
)(
4
2
)4(
4
1
)4(
5
2
hfE r r
f
hab
f
h
fE r r
m
n
i
im
收敛性：
。代数精度为截断误差：
难的。
的估计是困然而或积分值，需要事先给定的近似满足条件实际计算中，要想得到注：
数值稳定的。复合辛埔生求积公式是数值稳定性：
)(m a x.
)(
)(2
.1:1,
3
2;:1,
3
4;
3
)4(
2
2
0
2
0
2120
xfhn
fE r r
abm h
mi
h
mi
hh
m
m
i
i
m
i
i
iin
小结,1）在同一积分问题和同样精度要求下，复合梯形比复合辛埔生计算量大。
2）复合求积公式是实用的方法，但事先确定满足精度要求的
n或 h比较困难。
三、不事先选步长的求积算法
1、区间逐次分半法
( 1 ),,2,1,0
.))12((
2
1
))12(()(
2
))()((
4
)
2
(
2
))()((
4
)())()((
2
).
2
(2],[
2
1
111
2
1
11
12
1
1)2(2
1
12
1
11
1
1
1
m
hiafhTT
hiafhihaf
h
bfaf
h
h
iaf
h
bfaf
h
ihafhbfaf
h
T
ab
hbaT
m
mm
m
m
i
mmmm
i
mm
i
m
mm
i
mmm
i
mm
m
m
m
m
m
m
份的复合梯形求积公式等分为将区间记
( 2 ),4/3)(
))((4/3))(()(
)(
4
1
)(
)
~
()
~
~
(
.)(
2
1
)
~
~
(),
~
~
(
12
)(
)(
2
1
)
~
(
))
~
(
12
(
4
1
)
~
(
12
)(
11
11
1
2
1
2
2
1
1
22
1
1
1
mmm
mmmmm
mm
i
i
m
m
m
i
i
m
mm
m
TTTfI
TfITfITfITT
TfITfIm
ff
m
fff
h
TfI
ff
f
h
f
h
TfI
m
m
）（充分大时，于是：当
。
充分大时，有，当根据平均值稳定性理论又综合上述，得到基于复合梯形求积公式的自动步长求计算法：
否则继续！；分值则停止计算，并输出积若对于
11
2
1
111
1
0
0
,4/3
.))12((
2
1
2
,2,1,0;
));()((
2
mmm
i
mmmm
m
m
TTT
hiafhTT
h
h
m
abh
bfaf
ab
T
m
程序清单：
function [n,t]=tixingfaqiuji(fun_name,a,b,epsl)
%求函数在区间 [a,b]上的积分 t,n=2^m。
h=b-a;
t0=h/2*(feval(fun_name,a)+feval(fun_name,b));
m=0;
while 1
h=h/2;
t=0;
for i=1:2^m
t=t+feval(fun_name,a+(2*i-1)*h);
end
t=t0/2+h*t;
if abs(t-t0)<=3/4*epsl,break,end;
m=m+1;
t0=t;
End
N=2^m;
20 1)( dxexf x演示：求积分执行 [n,t]=tixingfaqiuji(0,2,0.0001)的结果：
n = 64,t = 4.0070.
2、查尔逊（ Richardson)外推技术与龙贝格（ Romberg)积分法查尔逊（ Richardson)外推：;)()()()(
0
,;)(
.0,)(
2211
21
21
*
0
21
21
*
0
*
0
构造新的算法：
有对于且其台劳展开为：设算法
k
k
pk
k
pppp
k
p
k
pp
hqahqahqaFqhF
q
ppp
hahahaFhF
hFhF
kpkpp
p
hahaFq hFqqhFhF 2
1
1
2
*00
1 1
)()()(
以此类推有：
阶算法它是 1
1
*11
.,2,1
1
)()(
)( 1
j
p
jp
j
p
j
j
p
j
haF
q
hFqqhF
hF j
j
j
进一步，再加上步长收缩 0<q<1
.,,2,1;,2,1
1
)()(
)(
1
1
1
mjm
q
hqFqhqF
hqF
j
j
p
jm
j
pjm
jjm
j
。，则只要 *1
1
*
1
*
)()()(
)()()(
);()(),()( 1
FhFhFhF
hohFhF
hoFhFhoFhF
mmm
p
mm
p
m
p
m
m
mm
**将上述技术应用于复合梯形求计算法称之为 龙贝格算法：
14
),1()1,1(4
)
2
1
(1
),1()
2
1
()1,1(
),(
.,,2,1,3,2,1
2/)(),0(
2
2
j
j
j
j
m
m
jmjTjmjT
jmjTjmjT
jmjT
mjm
abhmT
求积公式。则有时的复合梯形是记
)()0,(
./)(
)(
)12
2
2
4
4
2
2
m
m
k
kn
homT
nabh
hahahaIfTT( h )
（
龙贝格算法进程：
T ( m,0) 1)-mT ( 1,m)T ( 0,
T ( 2,0) T ( 1,1) T ( 0,2)
T ( 1,0) T ( 0,1)
T ( 0,0)
龙贝格算法程序清单：
function [n,s]=longbeigesf(fun_name,a,b,epsl)
%计算函数 在区间 [a,b]上的积分值，采用龙贝格算法。
h=b-a;
t(1,1)=h/2*(feval(fun_name,a)+feval(fun_name,b));
m=1;
while 1
h=h/2;
t(1,m+1)=0;
for i=1:2^(m-1)
t(1,m+1)=t(1,m+1)+feval(fun_name,a+(2*i-1)*h);
end
t(1,m+1)=t(1,m)/2+h*t(1,m+1);
for j=1:m
t(j+1,m+1-j)=(4^j*t(j,m-j+2)-t(j,m-j+1))/(4^j-1);
end
if abs(t(m+1,1)-t(m,1))<=epsl,break,end;
m=m+1;
end
s=t(m+1,1);
n=2^m;
程序一：二维数组 T(i,j)的情况
function [n,s]=longbeigesf2(fun_name,a,b,epsl)
%计算函数 在区间 [a,b]上的积分值，采用龙贝格算法。
h=b-a;
t0(1)=h/2*(feval(fun_name,a)+feval(fun_name,b));
m=1;
while 1
h=h/2;
t1(1)=0;
for i=1:2^(m-1)
t1(1)=t1(1)+feval(fun_name,a+(2*i-1)*h);
end
t1(1)=t0(1)/2+h*t1(1);
for j=1:m
t1(j+1)=(4^j*t1(j)-t0(j))/(4^j-1);
end
if abs(t1(m+1)-t0(m))<=epsl,break,end;
m=m+1;
t0=t1;
end
s=t1(m+1);
n=2^m;
程序二：一维数组 t0(n),t1(n)代替 T(I,J)
20 1)( dxexf x演示：求积分执行 [n,s]=longbeigesf(0,2,0.0001)的结果：
n = 8,s = 4.0070
执行 [n,s]=longbeigesf2(0,2,0.0001)的结果：
n = 8,s = 4.0070
插值型算法总结：
*插值型算法适用于有限区间 [a,b]上连续函数的积分问题；
**龙贝格积分算法是插值型算法中最有效的算法；
§ 4.3 高斯型求积公式高斯型求积公式是一种既能求有限区间上函数的积分，
又能求无限区间上函数积分的方法。
.
)()()(
),
),()()()(
1
高代数精度稳定性好，而且具有最使其不仅收敛性和数值求积公式：
的。我们的目的是构造可以有限也可以是无限区间（
上的权函数，为区间，其中问题：求
n
i
ii
b
a
b
a
xfdxxfx
ba
baxdxxfx
那么 一个 n个互异节点的求积公式，最高代数精度是多少？
我们知道：一个 n个互异节点的求积公式，代数精度至少为 n-1.
.2 nx ii 精度达到使上述求积公式的代数和不存在这样的节点结论一：
.2
,0)(,0)()(
)()()()(
1
22
2
2
1
nn
xfxfx
xxxxxxxf
b
a
n
i
ii
n
的代数精度个互异节点的求积公式而
，证明：取
有可能达到 2n-1吗？
.12
)()(
)(
)(
)()(
)(
)(
)()()(
,,,)(
)(),
*
1
21
n
dx
xpxx
xp
x
dx
xxx
x
x
xfdxxfx
xxxnxp
nxba
b
a
ini
n
b
a
ini
n
i
n
i
ii
b
a
nn
的代数精度达到朗日求积公式：
作为节点的拉格个零点的项式次正交多的上关于权选取区间（
结论二：
的多项式。均为次数不超过和
，
的多项式，则，它是次数不超过证明：
1)()(
)()()()(
12)(
nxrxq
xrxxqxp
nxp
n?
n
i
ii
b
a
ii
n
i
ii
b
a
b
a
nnn
n
n
nn
b
a
b
a
n
b
a
xpdxxpx
xrxp
xrdxxrxdxxpx
xpaxp
a
x
xpcxpcxpcxq
dxxrxdxxxqxdxxpx
1
1
111100
)!()()(
)()(
),()()()()(
)(),(
1
)(
)()()()(
)()()()()()()(
，又的首一系数，为
*
定义：对于 n个互异节点的求积公式，如果它的代数精度是 2n-1，
则称它为高斯型求积公式。
是高斯型求积公式。
则求积公式
）
个零点。的次正交多项式的上带权取区间（节点）
型求积公式的方法：结论二明确了构造高斯
n
i
ii
b
a
b
a
ini
n
i
n
n
xfdxxfx
dx
xpxx
xp
x
nxpn
xbaxxx
1
21
)()()(
)()(
)(
)(2
)(
)(),,,,1
高斯型求积公式的截断误差：
的多项式。为次数不超过其中，
项式。则的多条件：
为满足个零点，的次正交多项式的上带权为区间（设
12)(
),(
)!2(
)(
)()(
:1),()(),()(
)()(
)(),,,,
2
)2(
21
nxp
x
n
f
xpxf
nixfxpxfxp
xpnxpn
xbaxxx
n
n
iiii
n
n
b
a
n
n
n
b
a
n
n
b
a
n
nn
i
ii
b
a
n
nn
i
ii
b
a
n
n
b
a
b
a
dxxpx
na
f
dxx
n
f
xfE r r
dxx
n
f
xxf
dxx
n
f
xxp
dxx
n
f
xdxxpxdxxfx
)()(
)!2(
)(
)(
)!2(
)(
)()(
)(
)!2(
)(
)()(
)(
)!2(
)(
)()(
)(
)!2(
)(
)()()()()(
2
2
)2(
2
)2(
2
)2(
1
2
)2(
1
2
)2(
误差：高斯型求积公式的阶段
数值稳定性,
n
i
b
a
i
i
dxx
ni
1
)()2
:1,0)1
结论：
有上述结论看出，高斯型求积公式是数值稳定的。
由于正交多项式随区间合权函数的不同而不同，因此有不同类型的高斯型求积公式。下面介绍几种常用的高斯型求积公式。
12
2
])!2[(
)!(2
)2(
)(
)(*
)
22
(
2
)(*
,:1,
))()(1(
2
*
)*
)()(
1
2
42)2(
1
1
22
1
1
1
nn
n
n
f
fE r r
dtt
abba
f
ab
dxxf
ni
xLx
xxLLe g e n d r e
xfdxxf
l e g e n d r eG a u s s
nn
b
a
ini
i
in
n
i
ii
！
同时也可查表。系数可查表得到；的零点（多项式求积公式、
)(
)!2(
)!(
)(*
,:1,
))((
)!(
*
)*
)()(
2
)2(
2
2
2
1
0
n
ini
i
in
n
i
ii
x
f
n
n
fE r r
ni
xUx
n
xxULag ue r r e
xfdxxfe
l a g u e r r eG a u s s
同时也可查表。系数可查表得到；的零点（多项式求积公式、
)(
)!2(2
)!(
)(*
,:1,
))((
!2
*
)*
)()(
3
)2(
2
1
1
2
n
n
in
n
i
in
n
i
ii
x
f
n
n
fE r r
ni
xH
n
xxHH e r m i t e
xfdxxfe
H e r m i t eG a u s s
同时也可查表。系数可查表得到；的零点（多项式求积公式、
)(
)!2(2
2
)(*
:1,*
2
1)(2
c os
)*
)()(
1
1
4
)2(
2
1
1
1 2
n
n
i
i
n
n
i
ii
f
n
fE r r
ni
n
n
in
x
xTc he by s he v
xfdxxf
x
c he by s he vG aus s
。系数；
的零点（多项式求积公式、
本章评述：
高斯型求积公式使用少节点可得到高精度的结果，这是它明显的优点，另一优点是能计算许多广义积分。
高斯型求积公式也有明显的缺点，这就是节点和求积系数难以得到，需要查表，并且但节点数目增加，原来节点都作废，
先前计算的函数值不能重新使用。
在很多情况下，高斯型求积公式收敛快速性可弥补上述缺点。
从经验来说，对于有限区间的积分问题利用插值型求积公式是比较方便的。
附录：二维区域上的数值积分
),(
),()(
)()(
,),()(
),(),(
m
0j
j
)(
)(
0
)(
)(
)(
)(
ji
xd
xc
ii
b
a
n
i
ii
xd
xc
D
b
a
xd
xc
yxf
dyyxfxg
xgdxxgI
dyyxfxg
dyyxfdxdx dyyxfI
i
i
则令问题：
§ 4.1 数值积分的一般概念值型求积公式。的数值积分公式称为插的插值函数！由此得到及函数值上对应于分划是区间的简单函数？？若逼近误差：
易算！数值求积公式：
此时的简单函数数值方法：找一个逼近的原函数不易找到。问题：
)(),(,),(),(
],[)()(
)(
)()()(
)()()(
)()()(
),()(
)(?,)()(
110
110
b
a
nn
nn
b
a
b
a
b
a
xfxfxfxf
bxxxxa
baxpxf
dxxR
pIfIfE r r
dxxppIdxxf
xRxpxf
xpxf
xfdxxffI
次代数精度。称此求积公式具有
）。则不精确成立（次多项式个而对于某都精确成立（的多项式于
）对于任意次数不高公式（所谓代数精度：若求积
。并且代数精度尽量的高我们的任务是希望误差
）（
一般形式为：显然，数值求积公式的
m
pE r rxpm
pE r rxpm
pIfIfE r r
xffI
mm
mm
n
i
ii
0)
~
()(
~
1
),0)()(
1
)()()(
1 )()(
11
0
称为高斯型求积公式。则由此得到的求积公式上的最佳平方逼近函数的简单函数是区间若逼近,],[)( baxf
。时，算法是数值稳定的则当且数值稳定性：若收敛性：
次代数精度。公式具有则求积而若
：数精度可采用如下方法确定一个求积公式的代然联系。的误差小。他们没有必高，也不代表求积公式数精度不代表代数精度高。代）求积公式的误差小，
量指标。是求积公式优良性的度数精度精度的度量标志。而代）求积公式的误差计算注：
n
i
iiiii
n
kk
yxf
nfE r r
xE r rmkxE r r
0
1
,,)()5
.,0)()4
m
,0)(,,,1,0,0)(
)3
2
1
§ 4.2 插值型求积公式
b
a
n
n
n
b
a
ii
n
i
iin
n
n
dxx
n
f
fE r r
dxxlxffIfI
xf
xfxfxxxa
.)(
)!1(
)(
)(:
.)(,)()()(
)(,
),(),(
1
)1(
0
1010
截断误差拉格朗日求积公式：
为逼近函数，得到构成的拉格朗日函数作及对应的型值用分划拉格朗日求积公式一、
但不一定数值稳定。虽然
）代数精度至少为注：
,)2;1
0
ab
n
n
i
i
1
)(
12
)(
)(
))()((
2
)(
1*
5
6)4
,.1
3
3
1
1
0
代数精度为时）梯形公式（
柯特斯求积公式：）几个低次牛顿公式。求积公式。即复合求积段低次插值下的实用的求积公式应是分也难以保证。所以真正果很差，自然求积精度时拉格朗日插值逼近效我们知道数值不稳定。说明对大的但为偶数时代数精度为柯特斯公式，且求积公式称为牛顿）等距节点的拉格朗日
f
ab
fE r r
bfaf
ab
fI
n
n
nnn
n
i
i
3
)(
80
3
)(
3
) ],()2(3)(3)([
8
)(
38/3
3
2
),(
90
)(
)]()
2
(4)([
6
)(
2,*
)4(5
3
3
)4(
5
2
2
代数精度为时）公式（
代数精度为时）辛埔生求积公式（
fhfE r r
ab
h
bfhafhafaf
ab
fI
nS i m p s o n
ab
hf
h
fE r r
bf
ba
faf
ab
fI
nS i m p s o n
二、复化求积公式
1、复化梯形公式
1
1
1
1
1
3
1
1
1
)())()((
2
))()((
2
)(
)(
12
))()((
2
)()(
./)(,:0,],[
1
n
i
n
i
iin
n
i
i
n
i
ii
b
a
n
i
x
x
i
ihafhbfaf
h
xfxf
h
fT
f
h
xfxf
h
dxxfdxxf
nabhniihaxba
i
i
复化梯形求积公式为：
的分划：给定
难的。
的估计是困然而或积分值，需要事先给定的近似满足条件实际计算中，要想得到注：
值稳定的。复合梯形求积公式是数数值稳定性：
收敛性：
。代数精度为截断误差：
)(m a x.
)(
)(
.1:1,,
2
).()(
1
)(
12
)(
)(
12
)(
00
0
2
2
1
3
xfhn
fE r r
abn h
nih
h
hfE r r
f
hab
f
h
fE r r
n
n
i
i
n
i
i
in
n
n
i
in
2、复合辛埔生求积公式
))()(4)((
3
)(
)(
90
))()(4)((
3
)()(
./)(,2,2:0,
],[
1
212)1(22
1
)4(
5
1
212)1(2
1
2
)1(2
m
i
iiim
m
i
i
m
i
iii
b
a
m
i
x
x
i
xfxfxf
h
fS
f
h
xfxfxf
h
dxxfdxxf
nabhmnmiihax
ba
i
i
：复化辛埔生求积公式为的分划：给定
1
11
2
))2((
3
2
))12((
3
4
))()((
3
)(
m
i
m
i
m
hiaf
h
hiaf
h
bfaf
h
fS即
).()(
3
)(
1 80
)(
)(
90
)(
4
2
)4(
4
1
)4(
5
2
hfE r r
f
hab
f
h
fE r r
m
n
i
im
收敛性：
。代数精度为截断误差：
难的。
的估计是困然而或积分值，需要事先给定的近似满足条件实际计算中，要想得到注：
数值稳定的。复合辛埔生求积公式是数值稳定性：
)(m a x.
)(
)(2
.1:1,
3
2;:1,
3
4;
3
)4(
2
2
0
2
0
2120
xfhn
fE r r
abm h
mi
h
mi
hh
m
m
i
i
m
i
i
iin
小结,1）在同一积分问题和同样精度要求下，复合梯形比复合辛埔生计算量大。
2）复合求积公式是实用的方法，但事先确定满足精度要求的
n或 h比较困难。
三、不事先选步长的求积算法
1、区间逐次分半法
( 1 ),,2,1,0
.))12((
2
1
))12(()(
2
))()((
4
)
2
(
2
))()((
4
)())()((
2
).
2
(2],[
2
1
111
2
1
11
12
1
1)2(2
1
12
1
11
1
1
1
m
hiafhTT
hiafhihaf
h
bfaf
h
h
iaf
h
bfaf
h
ihafhbfaf
h
T
ab
hbaT
m
mm
m
m
i
mmmm
i
mm
i
m
mm
i
mmm
i
mm
m
m
m
m
m
m
份的复合梯形求积公式等分为将区间记
( 2 ),4/3)(
))((4/3))(()(
)(
4
1
)(
)
~
()
~
~
(
.)(
2
1
)
~
~
(),
~
~
(
12
)(
)(
2
1
)
~
(
))
~
(
12
(
4
1
)
~
(
12
)(
11
11
1
2
1
2
2
1
1
22
1
1
1
mmm
mmmmm
mm
i
i
m
m
m
i
i
m
mm
m
TTTfI
TfITfITfITT
TfITfIm
ff
m
fff
h
TfI
ff
f
h
f
h
TfI
m
m
）（充分大时，于是：当
。
充分大时，有，当根据平均值稳定性理论又综合上述，得到基于复合梯形求积公式的自动步长求计算法：
否则继续！；分值则停止计算，并输出积若对于
11
2
1
111
1
0
0
,4/3
.))12((
2
1
2
,2,1,0;
));()((
2
mmm
i
mmmm
m
m
TTT
hiafhTT
h
h
m
abh
bfaf
ab
T
m
程序清单：
function [n,t]=tixingfaqiuji(fun_name,a,b,epsl)
%求函数在区间 [a,b]上的积分 t,n=2^m。
h=b-a;
t0=h/2*(feval(fun_name,a)+feval(fun_name,b));
m=0;
while 1
h=h/2;
t=0;
for i=1:2^m
t=t+feval(fun_name,a+(2*i-1)*h);
end
t=t0/2+h*t;
if abs(t-t0)<=3/4*epsl,break,end;
m=m+1;
t0=t;
End
N=2^m;
20 1)( dxexf x演示：求积分执行 [n,t]=tixingfaqiuji(0,2,0.0001)的结果：
n = 64,t = 4.0070.
2、查尔逊（ Richardson)外推技术与龙贝格（ Romberg)积分法查尔逊（ Richardson)外推：;)()()()(
0
,;)(
.0,)(
2211
21
21
*
0
21
21
*
0
*
0
构造新的算法：
有对于且其台劳展开为：设算法
k
k
pk
k
pppp
k
p
k
pp
hqahqahqaFqhF
q
ppp
hahahaFhF
hFhF
kpkpp
p
hahaFq hFqqhFhF 2
1
1
2
*00
1 1
)()()(
以此类推有：
阶算法它是 1
1
*11
.,2,1
1
)()(
)( 1
j
p
jp
j
p
j
j
p
j
haF
q
hFqqhF
hF j
j
j
进一步，再加上步长收缩 0<q<1
.,,2,1;,2,1
1
)()(
)(
1
1
1
mjm
q
hqFqhqF
hqF
j
j
p
jm
j
pjm
jjm
j
。，则只要 *1
1
*
1
*
)()()(
)()()(
);()(),()( 1
FhFhFhF
hohFhF
hoFhFhoFhF
mmm
p
mm
p
m
p
m
m
mm
**将上述技术应用于复合梯形求计算法称之为 龙贝格算法：
14
),1()1,1(4
)
2
1
(1
),1()
2
1
()1,1(
),(
.,,2,1,3,2,1
2/)(),0(
2
2
j
j
j
j
m
m
jmjTjmjT
jmjTjmjT
jmjT
mjm
abhmT
求积公式。则有时的复合梯形是记
)()0,(
./)(
)(
)12
2
2
4
4
2
2
m
m
k
kn
homT
nabh
hahahaIfTT( h )
（
龙贝格算法进程：
T ( m,0) 1)-mT ( 1,m)T ( 0,
T ( 2,0) T ( 1,1) T ( 0,2)
T ( 1,0) T ( 0,1)
T ( 0,0)
龙贝格算法程序清单：
function [n,s]=longbeigesf(fun_name,a,b,epsl)
%计算函数 在区间 [a,b]上的积分值，采用龙贝格算法。
h=b-a;
t(1,1)=h/2*(feval(fun_name,a)+feval(fun_name,b));
m=1;
while 1
h=h/2;
t(1,m+1)=0;
for i=1:2^(m-1)
t(1,m+1)=t(1,m+1)+feval(fun_name,a+(2*i-1)*h);
end
t(1,m+1)=t(1,m)/2+h*t(1,m+1);
for j=1:m
t(j+1,m+1-j)=(4^j*t(j,m-j+2)-t(j,m-j+1))/(4^j-1);
end
if abs(t(m+1,1)-t(m,1))<=epsl,break,end;
m=m+1;
end
s=t(m+1,1);
n=2^m;
程序一：二维数组 T(i,j)的情况
function [n,s]=longbeigesf2(fun_name,a,b,epsl)
%计算函数 在区间 [a,b]上的积分值，采用龙贝格算法。
h=b-a;
t0(1)=h/2*(feval(fun_name,a)+feval(fun_name,b));
m=1;
while 1
h=h/2;
t1(1)=0;
for i=1:2^(m-1)
t1(1)=t1(1)+feval(fun_name,a+(2*i-1)*h);
end
t1(1)=t0(1)/2+h*t1(1);
for j=1:m
t1(j+1)=(4^j*t1(j)-t0(j))/(4^j-1);
end
if abs(t1(m+1)-t0(m))<=epsl,break,end;
m=m+1;
t0=t1;
end
s=t1(m+1);
n=2^m;
程序二：一维数组 t0(n),t1(n)代替 T(I,J)
20 1)( dxexf x演示：求积分执行 [n,s]=longbeigesf(0,2,0.0001)的结果：
n = 8,s = 4.0070
执行 [n,s]=longbeigesf2(0,2,0.0001)的结果：
n = 8,s = 4.0070
插值型算法总结：
*插值型算法适用于有限区间 [a,b]上连续函数的积分问题；
**龙贝格积分算法是插值型算法中最有效的算法；
§ 4.3 高斯型求积公式高斯型求积公式是一种既能求有限区间上函数的积分，
又能求无限区间上函数积分的方法。
.
)()()(
),
),()()()(
1
高代数精度稳定性好，而且具有最使其不仅收敛性和数值求积公式：
的。我们的目的是构造可以有限也可以是无限区间（
上的权函数，为区间，其中问题：求
n
i
ii
b
a
b
a
xfdxxfx
ba
baxdxxfx
那么 一个 n个互异节点的求积公式，最高代数精度是多少？
我们知道：一个 n个互异节点的求积公式，代数精度至少为 n-1.
.2 nx ii 精度达到使上述求积公式的代数和不存在这样的节点结论一：
.2
,0)(,0)()(
)()()()(
1
22
2
2
1
nn
xfxfx
xxxxxxxf
b
a
n
i
ii
n
的代数精度个互异节点的求积公式而
，证明：取
有可能达到 2n-1吗？
.12
)()(
)(
)(
)()(
)(
)(
)()()(
,,,)(
)(),
*
1
21
n
dx
xpxx
xp
x
dx
xxx
x
x
xfdxxfx
xxxnxp
nxba
b
a
ini
n
b
a
ini
n
i
n
i
ii
b
a
nn
的代数精度达到朗日求积公式：
作为节点的拉格个零点的项式次正交多的上关于权选取区间（
结论二：
的多项式。均为次数不超过和
，
的多项式，则，它是次数不超过证明：
1)()(
)()()()(
12)(
nxrxq
xrxxqxp
nxp
n?
n
i
ii
b
a
ii
n
i
ii
b
a
b
a
nnn
n
n
nn
b
a
b
a
n
b
a
xpdxxpx
xrxp
xrdxxrxdxxpx
xpaxp
a
x
xpcxpcxpcxq
dxxrxdxxxqxdxxpx
1
1
111100
)!()()(
)()(
),()()()()(
)(),(
1
)(
)()()()(
)()()()()()()(
，又的首一系数，为
*
定义：对于 n个互异节点的求积公式，如果它的代数精度是 2n-1，
则称它为高斯型求积公式。
是高斯型求积公式。
则求积公式
）
个零点。的次正交多项式的上带权取区间（节点）
型求积公式的方法：结论二明确了构造高斯
n
i
ii
b
a
b
a
ini
n
i
n
n
xfdxxfx
dx
xpxx
xp
x
nxpn
xbaxxx
1
21
)()()(
)()(
)(
)(2
)(
)(),,,,1
高斯型求积公式的截断误差：
的多项式。为次数不超过其中，
项式。则的多条件：
为满足个零点，的次正交多项式的上带权为区间（设
12)(
),(
)!2(
)(
)()(
:1),()(),()(
)()(
)(),,,,
2
)2(
21
nxp
x
n
f
xpxf
nixfxpxfxp
xpnxpn
xbaxxx
n
n
iiii
n
n
b
a
n
n
n
b
a
n
n
b
a
n
nn
i
ii
b
a
n
nn
i
ii
b
a
n
n
b
a
b
a
dxxpx
na
f
dxx
n
f
xfE r r
dxx
n
f
xxf
dxx
n
f
xxp
dxx
n
f
xdxxpxdxxfx
)()(
)!2(
)(
)(
)!2(
)(
)()(
)(
)!2(
)(
)()(
)(
)!2(
)(
)()(
)(
)!2(
)(
)()()()()(
2
2
)2(
2
)2(
2
)2(
1
2
)2(
1
2
)2(
误差：高斯型求积公式的阶段
数值稳定性,
n
i
b
a
i
i
dxx
ni
1
)()2
:1,0)1
结论：
有上述结论看出，高斯型求积公式是数值稳定的。
由于正交多项式随区间合权函数的不同而不同，因此有不同类型的高斯型求积公式。下面介绍几种常用的高斯型求积公式。
12
2
])!2[(
)!(2
)2(
)(
)(*
)
22
(
2
)(*
,:1,
))()(1(
2
*
)*
)()(
1
2
42)2(
1
1
22
1
1
1
nn
n
n
f
fE r r
dtt
abba
f
ab
dxxf
ni
xLx
xxLLe g e n d r e
xfdxxf
l e g e n d r eG a u s s
nn
b
a
ini
i
in
n
i
ii
！
同时也可查表。系数可查表得到；的零点（多项式求积公式、
)(
)!2(
)!(
)(*
,:1,
))((
)!(
*
)*
)()(
2
)2(
2
2
2
1
0
n
ini
i
in
n
i
ii
x
f
n
n
fE r r
ni
xUx
n
xxULag ue r r e
xfdxxfe
l a g u e r r eG a u s s
同时也可查表。系数可查表得到；的零点（多项式求积公式、
)(
)!2(2
)!(
)(*
,:1,
))((
!2
*
)*
)()(
3
)2(
2
1
1
2
n
n
in
n
i
in
n
i
ii
x
f
n
n
fE r r
ni
xH
n
xxHH e r m i t e
xfdxxfe
H e r m i t eG a u s s
同时也可查表。系数可查表得到；的零点（多项式求积公式、
)(
)!2(2
2
)(*
:1,*
2
1)(2
c os
)*
)()(
1
1
4
)2(
2
1
1
1 2
n
n
i
i
n
n
i
ii
f
n
fE r r
ni
n
n
in
x
xTc he by s he v
xfdxxf
x
c he by s he vG aus s
。系数；
的零点（多项式求积公式、
本章评述：
高斯型求积公式使用少节点可得到高精度的结果，这是它明显的优点，另一优点是能计算许多广义积分。
高斯型求积公式也有明显的缺点，这就是节点和求积系数难以得到，需要查表，并且但节点数目增加，原来节点都作废，
先前计算的函数值不能重新使用。
在很多情况下，高斯型求积公式收敛快速性可弥补上述缺点。
从经验来说，对于有限区间的积分问题利用插值型求积公式是比较方便的。
附录：二维区域上的数值积分
),(
),()(
)()(
,),()(
),(),(
m
0j
j
)(
)(
0
)(
)(
)(
)(
ji
xd
xc
ii
b
a
n
i
ii
xd
xc
D
b
a
xd
xc
yxf
dyyxfxg
xgdxxgI
dyyxfxg
dyyxfdxdx dyyxfI
i
i
则令问题：
