辛普森公式

Simpson算法及其推广形式
摘要：本文研究了辛普森公式的数值积分的计算方法问题，并且更进一步研究了变步长复化的辛普森公式和二重积分的辛普森公式的问题。

首先是对
一维辛普森公式和变步长复化辛普森公式以及二维辛普森公式的推导及
其算法，进行误差分析，并且列举了实例。

然后，对辛普森公式进行改
进，这里的改进最主要是对辛普森公式的代数精度进行提高，从而使辛
普森公式对积分的计算更加精确。

另外，还研究了辛普森公式的推广形
式。

最后，在结论的当中列举了一个例子。

关键词：辛普森公式算法改进推广形式二重积分的辛普森公式
Abstract：This paper first studies the calculation methods of the numerical integration in simpson formula, and then study of the long-simpson
formula and the double integral simpson formula problem. First, study the
algorithm and derived of one-dimensional simpson formula and
step-change in simpson formula, as well as two-dimensional simpson
formula, and then analysis the error. Finally , list the example. In this ,
improve the simpson formula. This improved the most important is to
incre ase the simpson formula’s accuracy of algebra. Besides, we study the
simpson formula’s promotion of forms. At the last, we list a example in
the conclusion.
Key word:The simpson formula, Algorithm, Improve, Promotion of forms, The simpson formula of the two-dimensional integral.
1 引言
辛普森公式主要的研究数值积分（numerical integration）的。

何谓数值积分呢？其是求定积分的近似值的数值方法。

即用被积函数的有限个抽样值的离散或加权平均近似值代替定积分的值。

求某函数的定积分时，在多数情况下，被积函数的原函数很难用初等函数表达出来，因此能够借助微积分学的牛顿-莱布尼兹公式计算定积分的机会是不多的。

另外，许多实际问题中的被积函数往往是列表函数或其他形式的非连续函数，对这类函数的定积分，也不能用不定积分方法求解。

由于以上原因，数值积分的理论与方法一直是计算数学研究的基本课题。

对微积分学做出杰出贡献的数学大师，如I.牛顿、L.欧拉、C.F.高斯等人也在数值积分这个领域做出了各自的贡献，并奠定了它的理论基础。

构造数值积分公式最通常的方法是用积分区间上的n 次插值多项式代替被积函数，由此导出的求积公式称为插值型求积公式。

特别在节点分布等距的情形称为牛顿-柯茨公式，例如梯形公式与抛物线公式就是最基本的近似公式。

但它们的精度较差。

龙贝格算法是在区间逐次分半过程中，对梯形公式的近似值进行加权平均获得准确程度较高的积分近似值的一种方法，它具有公式简练、计算结果准确、使用方便、稳定性好等优点，因此在等距情形宜采用龙贝格求积公式。

当用不等距节点进行计算时，常用高斯型求积公式计算，它在节点数目相同情况下，准确程度较高，稳定性好，而且还可以计算无穷积分。

数值积分还是微分方程数值解法的重要依据。

许多重要公式都可以用数值积分方程导出。

在数值积分的研究辛普森是一位很重要的数学家。

辛普森（Thomas Simpson，公元1710年8月20日─公元1761年5月14日）是英国著名的数学家，他生于英格兰列斯特郡，并卒于同地。

他的父亲是一位纺织工人，所以他主要靠自己力学成材，而他的第一份工作也是纺织。

他对数学的兴趣最初是由一次日蚀所引发起的，他在一位占卜师的指导之下，他学会了算术和基本的代数。

其后他放弃了纺织的工作，而当了一间学校的司阍。

凭借他刻苦和持久的努力，他证明了他在数学方面的能力，以致于1735年他能够解决了数个有关微积分的问题，1737年他便开始撰写有关数学的文章。

在1754年他成为了「淑女日记」一书的编辑，及后他到了伦敦的乌尔威治并出任数学教授一职，直至逝世。

辛普森最为人熟悉的贡献是他在插值法（Interpolation）及数值积分法（Numerical Method of Integration）方面，事实上他在概率方面也有一定的
工作，他在1740年推出了他的「机会的特性和法则」（The Nature and Laws of Chance ），而大部份他在这方面的结果也是建基于棣美弗早期的结果。

另外，当时有一群讲师巡回在伦敦咖啡屋讲学，而辛普森是当中最突出的一位。

他专研有关「误差理论」（Theory of Error ），并且意图证明数算平均数比单一观察较佳。

Simpson 公式就是他的代表的定理。

在辛普森之后，很多后人都在其的基础上，不断完善simpson 公式，使其公式越来越完善。

其中，华罗庚，王元的著作数《值积分及其应用》中，就是用数论的方法研究辛普森公式。

在其研究中，假设()f x 是一个在[],αβ内定义了的函数，以后如果用到几次微商，便假设定()f x 有几次微商。

我们用Euler 求和公式（详细参考该书第一章的内容）来推出普通数值积的Simpson 法。

从而得出下面的内容：
命
12
33
s t r R R R =+ （其中t R ，r R 参照该书第二章的1
和2节） 则得
()1
1
01102246n n s n i i i i R f t dt y y y y n α
β
βα--+==⎛
⎫
-=-
+++ ⎪⎝⎭
∑∑⎰， 而且simpson 公式的余项[由（1）与（2）]为
()()3
22''03
''2201211233633611232n s n
b x b x R f x dx n n b x b x f x dx n n βαβααβαβαα⎛⎫⎛⎫- ⎪ ⎪--⎛⎫⎛⎫⎝⎭ ⎪=-++⨯+ ⎪ ⎪ ⎪⎝⎭⎝
⎭ ⎪
⎝⎭
-⎛⎫-⎛⎫
⎛⎫⎛⎫=+-+ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭⎝
⎭⎰⎰ 其中
()()()1
13
''2012112n l n l n
f t dt y y y n
b x f x dx n n α
β
βαβαβαα-=-⎛⎫
-
++ ⎪⎝⎭
--⎛⎫
⎛
⎫⎛⎫=-+ ⎪ ⎪ ⎪⎝⎭⎝
⎭⎝⎭∑⎰⎰ （1）
()3
1''2
122
11242n r R b x f x dx n n βαβαα--
-⎛-⎫⎛⎫
⎛⎫⎛
⎫=-++- ⎪
⎪ ⎪ ⎪⎝⎭
⎝⎭⎝⎭⎝⎭
⎰
（2） 如果()()IV f x a x b ≤≤存在，由部分积分可得
()()''220'''330122122n
n
b x b x f x dx n b x b x f x dx n
n βααβα
βαα⎛⎫-⎛⎫⎛⎫+-+= ⎪ ⎪ ⎪⎝⎭⎝⎭⎝
⎭-⎛⎫-⎛
⎫⎛⎫-
+-+ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝
⎭⎰⎰
（此处用了()33100,02b b ⎛⎫
== ⎪⎝⎭
）
()()()()()()()'''4
42
440''''''2
44024412212219601222n
o
n IV n IV b x b x f x n
n b x b x f x dx
n n f f n b x b x f x dx
n
n b x b x n βαβααβαβααβααββαβααβα-⎛⎫-⎛
⎫⎛⎫=-
+-+ ⎪ ⎪ ⎪⎝⎭⎝
⎭⎝⎭-⎛⎫-⎛⎫⎛
⎫⎛⎫++-+ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝
⎭⎝⎭-=--⎛⎫-⎛⎫⎛⎫⎛⎫++-+ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝
⎭⎝⎭-⎛⎫
=+ ⎪⎝⎭⎰⎰0112960n
IV f x dx n βαα⎛⎫-⎛⎫⎛⎫--+ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭
⎰
即
()5
440111232960n
IV s R b x b x f x dx n n βαβαα-⎛⎫-⎛⎫
⎛⎫⎛⎫=+--⨯+ ⎪
⎪ ⎪ ⎪⎝⎭
⎝⎭⎝⎭⎝⎭
⎰、 定理1.1 如果()(),IV f x M x αβ≤≤≤，则 ()5
44
1802s
M
R n βα-≤
⋅⋅
证. 由于 ()4
40432
1432
20432
432
1122960111111222(2412247201261236096011111122224122472012612360960n b x b x dx x x x x x x n dx x x x x x x ⎛⎫+-- ⎪⎝⎭
⎛⎫⎛
⎫⎛⎫⎛⎫+++ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭ ⎪=-+-+-+-++ ⎪ ⎪⎝⎭
⎛⎫⎛
⎫⎛⎫⎛⎫+++ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭ +-+-+-+-++ ⎝⎭
⎰⎰121664
)6
63180218021802
dx n n ⎪⎪⎪⎛⎫=+= ⎪⋅⋅⋅⎝⎭⎰
故得定理。

定理1.2 如果
()
''f x 是单调非负递减函数，则
()()3
''3
324s M
R f n βαα-≤
.
证. 由
()1122220
102212x x b x dx dx ⎛⎫
=-+= ⎪⎝⎭
⎰⎰ 及对1
02
η<<
，常有 ()2
20
2
20
2320
12211112212226312224108b x b x dx x x x x dx x x dx η
η
η
ηη⎛⎫⎛⎫+- ⎪ ⎪⎝⎭⎝
⎭⎛⎫
⎛⎫⎛⎫=-+++-++ ⎪ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭
⎛⎫=
-=-≤ ⎪
⎝⎭
⎰
⎰⎰
由
()3
''22011232n
s R b x b x f x dx n n βαβαα-⎛⎫-⎛⎫
⎛⎫⎛⎫=+-+ ⎪
⎪ ⎪ ⎪⎝⎭
⎝⎭⎝⎭⎝
⎭⎰ 以及第二中值公式可得
()()()
()
3
''2203''3
11232324s R f b x b x dx n M
f n
ξβααβαα-⎛⎫⎛⎫⎛
⎫≤+- ⎪ ⎪ ⎪⎝⎭⎝⎭⎝
⎭-≤
⎰
即得定理1.2.
【】参考文献：华罗庚，王元， 数值积分及其应用 科学出版社 1963年 前人的经验都是很多的，这里只是列举了其中一个。

现在就来介绍一下一些预备知识。

首先，来介绍一下代数精度。

代数精度是用来衡量求积公式精确程度的一个概念，代数精度越高则求积公式就越准确。

定义1.1 如果一个一般的求积公式对于所以次数不大于m 的多项式()f x 都准确成立，即
()()0
n
b
k k a
k f x f x λ==∑⎰
而对于次数为m+1的某一个多项式()f x 并不准确成立，则称该求积公式具有m 次代数精度。

根据定义1.1，插值型求积公式()()0
n
b
k k a
k f x f x λ=≈∑⎰至少要具有n 次代数精度。

由于任一m 次多项式
()1110m m m m m p x a x a x a x a --=++⋅⋅⋅++ 是由k
x ()
0,1,2,,k m =⋅⋅⋅叠加而成的，因此可以利用函数
()k f x x =()0,1,2,,k m =⋅⋅⋅来判断求积公式的代数精度。

定理1.3 表明n+1个节点的求积公式是内插求积公式至少具有n 次代数精确度，也可能更高。

判断求积公式代数精度的步骤：
1.在内插求积公式两端分别代入最简单的代数多项式n x x Λ,,1进行计算，即代入f(x)=1, f(x)=x,…,f(x)=x n 计算。

2 .判断两端数值相等时的代数多项式的最高次数，则可以确定出代数精度。

即对于任意不超过n 次的代数多项式都准确成立，而对某一个n ＋1次代数多项式不成立，则其代数精度为n 。

例1. 确定积分公式[])1()(3)(23
1)(211
1
f x f x f dx x f ++≈⎰-中的待定常数，并确
定其代数精度。

解 令f (x )＝1，
左边＝21)(1
11
1==⎰⎰--dx dx x f 右边＝2)132(3
1
=++
令f (x )＝x
左边＝⎰⎰--==11110)(xdx dx x f 右边＝[]13231
21++x x
令f(x)＝x 2， 左边＝321
12=
⎰-dx x 右边＝)132(3
12
2
21++x x 得：
[]⎩⎨⎧=-=⎩⎨⎧-==⎪⎪⎩⎪⎪⎨
⎧=
++=++1266
.06899.05266.02899.032)132(3
1013231
21
21222121x x x x x x x x 或得
令f (x )＝x 3，
左边＝ 右边＝0.6106或 0.3494 ∴左≠右 ∴该求积公式具有2次代数精度。

下面，再来研究一下拉格朗日插值公式。

设连续函数y = f （x ）在[a, b]上对给定n + 1个不同结点：
x 0, x 1, …, x n 分别取函数值
y 0, y 1, …, y n
其中 y i = f (x i ) i = 0, 1, 2,…, n 试构造一个次数不超过n 的插值多项式
n n n x a x a x a a x P ++++=Λ2210)(
使之满足条件
i i n y x P =)( i = 0, 1, 2,…, n
先求n 次多项式l k (x ) k = 0, 1,…, n
使
⎩⎨
⎧≠==i
k i k x l i k ,
0,
1)( （3）
若做出这样的)(x l k ，则P n (x )的次数≤n ，另外，由（3）
i i n
k k k i n y x l y x P ==∑=)()(1 i = 0, 1, 2,…, n
即()n P x 满足插值条件n i y x P i i n ,,2,1,0,)(Λ==。

于是问题归结为具体求
出基本插值多项式l k (x )。

根据（3）式，x k 以外所有的结点都是l k (x )的根，因此令
∏≠=+--=-----=n
k
j j j n k k k x x x x x x x x x x x x x l 01110)
()())(())(()(λλΛΛ
又由l k (x k ) = 1，得：
1
1
3=⎰
-dx x
)
())(())((1
1110n k k k k k k k x x x x x x x x x x -----=
+-ΛΛλ
所以有：
∏
≠=+-+---=----------=
n
k
j j j
k j n k k k k k k k n k k k x x x x x x x x x x x x x x x x x x x x x x x x x l 011101110)
())(())(()())(())(()(ΛΛΛΛ
代入（3）式即得()n P x 的表达式
k n k n
k
j j j
k j k n
k k n y x x x x y x l x P ∑∏
∑=≠==⎪⎪⎪
⎭
⎫
⎝⎛--==000)()( （4）
（4）式称为拉格朗日插多项式。

我们把差f (x )－()n P x 称为用插值多项式()n P x 代替f (x )的余项，误差或
插值余项，记为：
)()()(x P x f x R n n -=
下面讨论余项的形式。

设x 是区间[a, b]中任意固定的数，若x 是插值节点x i ，则0)(=i n x R ，i =
0, 1, 2,…, n ，如果x 不是节点，考虑如下函数：
)()
()
()()()(11x R x t t P t f t n n n n ++-
-=ωωϕ
其中
∏=+-=n
k k n x x x 0
1)()(ω
由插值条件n i y x P i i n ,,2,1,0,)(Λ==知
0)(=i x ϕ i = 0, 1,…, n
而
0)()
()
()()()(11=-
-=++x R x x x P x f x n n n n ωωϕ
所以，ϕ(t )在[a, b]上有n + 2个互异零点。

应用罗尔定理，ϕ’(t )在ϕ(t )的每两个零点间有一个零点，即ϕ’(t ) 在[a, b]上至少有n + 1个零点，进一步对ϕ’(t )应用罗尔定理可知，ϕ" (t )在
[a, b]上至少有n 个零点，继续上述讨论就可推得ϕ(n +1)
(t )在[a, b]内至少有一
个零点，记之为ξ，即
0)()1(=+ξϕn
因为P n (t )为不高于n 次的多项式。

所以
)!1()(0
)()1(1)1(+==+++n t t P n n n n ω
于是
)
()
()!1()()(01)
1()1(x x R n f
n n n n ++++-
==ωξξϕ
移项、整理可得
)()!
1()
()(1)1(x n f x R n n n +++=ωξ
综上所述，得以下余项定理：
定理1.4： 设f (n )(x )在区间[a, b]上连续，f (n + 1)
(x )在[a, b]上存在，x 0, x 1,…,
x n 是[a, b]上互异的数，记插值问题n i y x P i i n ,,2,1,0,
)(Λ==的余项为
)()()(x P x f x R n n -=，那么，当x ∈ [a, b]时，有如下估计
)()!
1()
()(1)1(x n f x R n n n +++=ωξ
其中
∏=+-=n
j j n x x x 1
1)()(ω
2 辛普森公式及其算法
2.1 有限区间上的一维辛普森公式
2.1.1 利用插值多项式导出求积公式
设()f x 为被积函数，[],a b 为积分区间，01,,,n x x x ⋅⋅⋅为[],a b 内的n+1个互异的点，记()n L x 为相应的拉格朗日插值多项式，()n R x 为相应的插值余项，那么我们有
()()()n n f x I f R f =+ （2.1.1）
将上面的（2.1）式两边同事积分得
()()()b b b
n
n a
a
a
f x dx I f dx R f dx =+⎰⎰⎰
如果我们取
()()b
b
n a
a
f x dx I f dx ≈⎰
⎰ （2.1.2）
那么截断误差为
()b
n a
R f dx ⎰
利用拉格朗日插值公式（4），我们有
()()()0n
n k k k I x f x LB x ==∑
其中()k LB x ，0,1,,k n =⋅⋅⋅，是基函数。

从而上面的（2.1.2）式可以改成
()()()0
n
b
b k
k
a
a
k f x dx f x LB x dx =≈∑⎰⎰ （2.1.3）
由于各()k f x 都是常数，()k LB x 是多项式函数，所以从理论上讲，上面的（2.1.3）的右边是可以计算的。

记
()b
k k a w LB x dx =⎰ （2.1.4）
由上面的（2.1.3）即得
()()0
n
b
k
k
a
k f x dx f x w
=≈∑⎰ （2.1.5）
从而得到一种求数值积分的方法。

上面的（2.1.4）和（2.1.5）被称为插值公式导出的求积公式。

回归第十章关于多项式插值的结论可知，由于任意次数不超过n 的多项式与它的任意n+1个基点的插值多项式恒等，再由求积公式的代数精度定义，我们立即得到：由n+1个基点的拉格朗日插值多项式所形成的求积公式的代数精度至少是n 。

为了强调上面的（2.1.5）是利用n+1个基点的插值多项式导出的求积公式，我们把公式中的各k w 改写成(),w n k ，从而有
()(),b
k a w n k LB x dx =⎰ （2.1.6）
()()()0
,n
b
k
a
k f x dx f x w n k =≈∑⎰ （2.1.7）
2.1.2 牛顿-科特斯求积公式的推导
虽然从理论上讲，可以利用上面的（2.1.6）和（2.1.7）完成()f x 在区间[],a b 上的数值积分的近似计算，但是利用（2.1.6）计算各(),w n k 仍然是一件比较麻烦的事。

牛顿-科特斯求积公式的目的就是利用01,,,n x x x ⋅⋅⋅是区间[],a b 的一个等分分划的特点进一步简化（2.1.6）中的各(),w n k 的计算。

假如我们将积分区间[],a b 划分为n 等分，记
()
0,0,1,,k b a h n
x a x a kh k n x a th -⎫
=
⎪⎪=⎬⎪
=+=⋅⋅⋅⎪
=+⎭
（2.1.8）
那么我们有
()(),0,1,,,,0,1,,0
j k j x x t j h j n
x x k j h k j n dx hdt x a t x b t n -=-=⋅⋅⋅⎫
⎪
-=-=⋅⋅⋅⎪
⎪
=⎬⎪=⇒=⎪
=⇒=⎪
⎭ （2.1.9）
再注意到对0,1,,k n =⋅⋅⋅，拉格朗日插值的基函数()k LB x 可表示为
()0n
j k j k j
j k
x x LB x x x =≠-=-∏
所以利用上面（2.1.8）和（2.1.9）给出的记号约定以及相应的积分变换，我们有
()000,n n b n j a j j k j j k
j k
x x t j
w n k dx dx x x k j
==≠≠--==--∏∏⎰⎰
从而可以得到
()()
()000
,n
n
n
j j k
j j k
h
w n k t j dt k j =≠=≠=
--∏⎰∏ （2.1.10）
再注意到
()()()()()()
()011121!!
n
j j k
n k
k j k k n k k n k =≠--=-⋅⋅⋅⨯--⋅⋅⋅-+⎡⎤⎡⎤⎣⎦⎣⎦
=--∏
所以由上面的（2.1.10）可以得到
()()()()00
1,!!n k
n n j j k
h w n k t j dt k n k -=≠-=--∏⎰ 再利用
()
b a h n
-=
立即可以得到
()()()()()00
1,!!n k
n n j j k
w n k b a t j dt nk n k -=≠-=---∏⎰ (2.1.11) 在上面的(2.1.11)中，记
()()()()00
1,!!n k
n n j j k
C n k t j dt nk n k -=≠-=--∏⎰ （2.1.12） 则上面的（2.1.11）可表示为
()()(),,w n k b a C n k =-
从而前面的（2.1.7）可以表示为
()()()()0
,n
b
k
a
k f x dx b a C n k f x ==-∑⎰ （2.1.13）
上面的（2.1.12）和（2.1.13）称为牛顿-科特斯求积公式，利用（2.1.12）得到的各(),C n k 称为牛顿-科特斯系数。

虽然看上去利用（2.1.12）式计算为牛顿-科特斯系数也还是有些麻烦，但是与前面的（2.1.4）对比可知，各(),C n k 与具体问题无关，可以把它们“一次性地”计算出来，并制成表。

这样，人们就可以通过查表得到这些值，从而可以免去这一部分的计算。

这样一来，上面的（2.1.13）就真的称为机械求积公式了。

虽然在现代的计算机条件下牛顿-科特斯求积公式和牛顿-科特斯系数的作用已经很有限，但是牛顿-科特斯求积公式所体现出来的独具创造性的思维方式和大胆的探索精神为科学计算领域中的研究树立了光辉的榜样。

2.1.3 辛普森公式
在牛顿-科特斯求积公式中，如果取n=2时，那么k 可以取0，1，2，此时形成的的求积 公式就是辛普森公式，利用（2.1.10），我们可以得到
()()
()()()()
()()()()
()()2
2
1
1
00
1
0112,01220!2!6
142,10221!1!
6112,20122!0!
6
C t t dt C t t dt C t t dt -=--=
⨯⨯-=--=
⨯⨯-=
--=
⨯⨯⎰⎰⎰ 所以有
()()()()14
16626b
a
a b f x dx b a f a f f b ⎡+⎤⎛⎫=-+
+ ⎪⎢⎥⎝⎭⎣⎦
⎰
（2.1.14） 【】参考文献：甄西丰编著 实用数值计算方法 清华大学出版社 2006年
2.1.4误差分析
应用Newton-Cotes 型求积分公式（2.1.13）计算定积分()()b
a I f f x =⎰dx
时，一方面由于它是由式（2.1.1）去掉余项()n E f 得到的，因而产生了离散误差()n E f ；另一方面，由于计算机的字长是有限的，函数值可能带有误差，并且计算()n I f 还会有舍入误差。

关于离散误差，我们有下面的定理。

定理2.1 设n 为偶数且()f x 在[],a b 上有n+2阶连续导数，则Newton-Cotes 型求积分公式（2.1.13）的离散误差为
()()()
()()()2320
12!
n n n
n h f E f t t t n dt n η++=
-⋅⋅⋅-+⎰ (),a b η∈
若n 为奇数，且()f x 在[],a b 上有n+1阶连续导数，则
()()()
()()()1220
11!
n n n
n h f E f t t t n dt n ξ++=
-⋅⋅⋅-+⎰
(),a b ξ∈
特别，n=2时，Simpson 公式的离散误差为
()()()54290h E f f η=- 2
b a
h -=
由于Newton-Cotes 型求积分公式的误差公式证明比较繁杂，下面给出辛普森公式的误差的证明。

证明：首先考虑构造一个三次差值多项式()3p x 满足条件：
()()()()333,,22a b a b p a f a p b f b p f ++⎛⎫⎛⎫=== ⎪ ⎪⎝⎭⎝⎭和''322a b a b p f ++⎛⎫⎛⎫= ⎪ ⎪⎝⎭⎝⎭。

可以
证明满足上述条件的()3p x 与()f x 的误差为
()()()()()()2
434!2f a b f x p x x a x x b ξ+⎛⎫
-=--- ⎪⎝⎭
（2.1.15）
其中a b ξ≤≤。

对（2.1.15）式两边分别进行积分得
()()()()()()2
434!2b b
b
a
a
a
f a b f x dx p x dx x a x x b dx ξ+⎛⎫
-=--- ⎪⎝⎭
⎰⎰
⎰
又因为()3p x 为三次代数多项式，而辛普森公式的代数精度为3.所以有
()()()()()()()3333141662614
16
626b
a
a b p x dx b a p a p p b a b b a f a f f b ⎡+⎤
⎛⎫=-++ ⎪
⎢⎥⎝⎭⎣⎦
⎡+⎤⎛⎫=-+
+ ⎪⎢⎥⎝⎭⎣⎦
⎰
从而有
()()()()()
()
()()2
414
1
6
6264!
2b
a
b
a
a b f x dx b a f a f f b f
a b x a x x b dx ξ⎡+⎤⎛⎫--+
+ ⎪⎢⎥⎝⎭⎣⎦
+⎛⎫
=--
- ⎪⎝⎭
⎰
⎰
据提设，()
()4f
ξ在[a,b]上连续，而()()2
02a b x a x x b +⎛⎫---≤ ⎪⎝
⎭，由积分中
值定理可知，在[a,b]中存在一个点η，使
(
)
()
()()()()()()()()()
2
42
45
44!
214!22880
b
a
b a f a b x a x x b dx
a b f x a x x b dx b a f ξηη+⎛
⎫--- ⎪⎝
⎭+⎛
⎫=--- ⎪⎝
⎭-=-
⎰
⎰ 因此
()()()()()(
)
()
5
414166262880
b
a
a b f x dx b a f a f f b b a f η⎡+⎤⎛⎫--++ ⎪⎢⎥⎝⎭⎣⎦
-=-
⎰
得证。

【】参考文献：韩国强编著，数值分析，华南理工大学出版社，2005年
2.2 有限区间的变步长复化辛普森公式
2.2.1 复化辛普森公式
假设()012212221,,,,,,,,k k n k x x x x x x x --⋅⋅⋅⋅⋅⋅把积分空间[],a b 划分为2n 等分，我们也可以认为是其中的()02221,,,,,n k x x x x -⋅⋅⋅⋅⋅⋅把区间划分为n 等份，并且21k x -就是第
k 个子空间()221,k k x x -⎡⎤⎣⎦的中点。

记k I ，1,2,,k n =⋅⋅⋅为()f x 在第k 个子空间()221,k k x x -⎡⎤⎣⎦应用辛普森公式所求
得的积分值，则有
()()
()()2122146k k k k b a I f x f x f x n ---⎡
⎤=
++⎣
⎦
记 21
k n
n k k S I ===∑
则有
()()()()1
20221211426n n n n k k k k b a S f x f x f x f x n --==-⎡⎤
=+++⎢⎥⎣⎦
∑∑ （2.2.1）
上面的（2.2.1）式称为复化辛普森公式，虽然我们也可以直接编写2n S 的计算机程序，但是没有必要那么急，而是我们要改进一下变步长复化辛普森公式的性能。

2.2.2 变步长复化辛普森公式
假设()012212221,,,,,,,,k k n k x x x x x x x --⋅⋅⋅⋅⋅⋅把积分空间[],a b 划分为2n 等分，那么我们可以得到下面的3个不同的积分值：
（1） 利用2k x ，0,1,2,,k n =⋅⋅⋅这n+1个点处的函数值和复化梯形公式计算出
n T ；
（2） 利用j x ，0,1,2,,2k n =⋅⋅⋅这2n+1个点处的函数值和复化梯形公式计算出
2n T ；
（3） 利用复化辛普森公式计算出2n S ； 显然，2n S ，n T ，2n T 之间应该存在一定的关系。

在这里，我们应该先知道复化梯形公式
()()()10111
22n n n j j b a T f x f x f x n -=⎡⎤-=++⎢
⎥⎣⎦
∑ （2.2.2）
另外
()()()21202111
222n n n j j b a T f y f y f y n -=⎡⎤-=++⎢
⎥⎣⎦
∑ （2.2.3） 这里的这两个式子的证明就不给出了。

按照（2.2.2）和（2.2.3）式方括号内的表示形式整理（2.2.1）式的方括号内的数学表达式，不难得到
()2126n b a S n
-=∑+∑ （2.2.4）
其中
()()()()1
1022121
1
2244n
n n k k k k f x f x f x f x --==∑=+++∑∑
()()()1
20221
2n n k k f x f x f x -=∑=++∑
上面的1∑，2∑还可以进行进一步调整为
()()()()()()2110211
2022111
42211222n n j j n n k k f x f x f x f x f x f x -=-=⎫⎡⎤∑=++⎪
⎢⎥⎪⎣⎦⎬⎡⎤⎪
∑=++⎢⎥⎪
⎣⎦⎭
∑∑ （2.2.5）
利用前面的（2.2.2）和（2.2.3）式以及上面的（2.2.5）,我们可以得到
12241
,6363
n n b a b a T T --∑=∑= （2.2.6） 再由（2.2.1）和（2.2.6）式即得到
2243n n n T T
S -= （2.2.7）
按照我们习惯记法，取2k n =，则有122k n +=，利用[]1S k +表示由（2.2.7）式所得到的2n S ，利用复化梯形公式得到的梯形序列[][]0,1,,T T ⋅⋅⋅我们有
[][][]()141/3S k T k T k +=+- （2.2.8）
所以由梯形序列[][]0,1,,T T ⋅⋅⋅我们可以得到了一个新的序列式[][]1,2,,S S ⋅⋅⋅称之为辛普森序列。

2.2.3 复化辛普森公式的算法实现 （1）复化辛普森公式的计算步骤
1）确定步长()/h b a N =-（N 为等份数）。

()122,0S f a h S =+=。

2）对1,2,,1k N =⋅⋅⋅-计算
()112,S S f a kh h =+++ ()22S S f a kh =++ 3） ()()12426S h f a S S f b =+++⎡⎤⎣⎦ （2）算法流程图
复化辛普森公式的算法流程图见图1。

（3） 复化辛普森公式的matlab 程序
按照（2.2.1）编写复化辛普森求积函数（函数名：s_quad.m ） function I=S_quad(x,y) %复化辛普森求积公式，其中
%x 为向量，被积函数自变量的等距节点；
%y 为向量，被积函数在节点处的函数值。

n=length(x);m=length(y) if n~=m
error(‘the lengths of X and Y must be equal ’); return; end
if rem(n-1,2)~=0 I=T_quad(x,y); Return; End
N=(n-1)/2;h=(x(n)-x(1))/N;a=zeros(1,n); For k=1:N
a(2*k-1)= a(2*k-1)+1; a(2*k)= a(2*k)+4; a(2*k+1)= a(2*k+1)+1 end
I=h/6*sum(a*y');
下面给一个例子，用复化辛普森公式求积2
1
1x e dx --⎰，在积分区间中点与点之
间的间隔取为0.1。

解：输入
x=-1:0.1:1;y=exp(-x.^2); I=S_quad(x,y) 得到
I=1.4936
【】参考文献：薛毅编著，数值分析与实验，北京化工大学出版社，2005年
2.3 二重积分的辛普森公式
我们上面讨论的辛普森公式的求积分方法，稍经修改就可以用来计算重积分。

(),R
I f x y dxdy =⎰⎰ （2.3.1）
的计算法，此处，假设积分区域为矩形域：
(){},|,R x y a x b c y d =≤≤≤≤。

给定正整数n 和m ，取步长
()(),h b a n k d c m =-=-。

把积分I 写成
(),b
d
a
c
I f x y dxdy =⎰
⎰
对积分 (),d c
f x y dy ⎰
视x 为常数，令,0,1,,2j y c jk j m =+=⋅⋅⋅。

应用复合辛普森公式得
()()()()()()()
()102212114
4
4
,,2,4,,3,,,180
m m d
j j m c
j j k f x y dy f x y f x y f x y f x y d c k f x c d y μμ--==⎡
⎤=+++⎢⎥
⎣⎦
-∂-
∈∂∑∑⎰
于是
()()()()()()
10212121444
2,,334,,33,180
m b b j a a j m b b
j m a a j b
a
k k I f x y dx f x y dx
k k f x y dx f x y dx d c k f x dx y
μ-=-==+++-∂-
∂∑⎰⎰∑⎰⎰⎰
再令,0,1,,2i x a ih i n =+=⋅⋅⋅，则对每一0,1,,2j m =⋅⋅⋅，有
()()()()()()()
()
102212114
4
4
,,2,4,,3,,,180
n n
b
j j i j i j n j a
i i j j j h f x y dy f x y f x y f x y f x y f y b a h
a b x
εε--==⎡⎤
=+++⎢⎥
⎣⎦
∂--
∈∂∑∑⎰
从而
()()()()()()
()()
1002021011
111
2002221
11
1`1
2122211
1
[,2,4,9,2,4,8,2,n n i i i i m m n n j i j j j i m n m i j n j j i j hk
I f x y f x y f x y f x y f x y f x y f x y f x y --==---===---===≈+++++++∑∑∑∑∑∑∑∑
()()
()()
()()()
()1
0212211
11
212122111
1
102222121
1
224,8,16,4,,2,4,,].
m m n j i j j j i m n m
i j n j j i i n n
m i m i m i i n m f x y f x y f x y f x y f x y f x y f x y f x y ---===---===--==++++++++∑∑∑∑∑∑∑∑
其离散误差为
()()
()()
()
()()
44412200441
442121224
4
1
444
,,[2540,,4]
,.180
m j j j m
j j m m j b
a
f y k b a h f y E x x f y f y x x d c k f x dx y εεεεμ-=--=∂--∂=+∂∂∂∂++∂∂-∂-
∂∑∑
⎰
设44f
x
∂∂和44f y ∂∂的R 上连续，则
()()()()()()()()()4444444
4
4444
,,[6]540180,,[].
180
k b a h f d c b a k f E m x y d c b a f f h k x y ημημημη
μ--∂--∂=-
∂∂--∂∂=-
+∂∂%%%%
其中(),ημ，(),ημ%% R ∈。

现在，我们假设积分区域R 为曲边梯形()(),;,R a b c x d x =⎡⎤⎣⎦，它是上下分别由连续曲线()y d x =和()y c x = ()a x b ≤≤所限制，两侧由直线x=a 和x=b 所限制，并且(),x a b ∈时，()()c x d x <。

若(),f x y 在R 上连续，则
()()()
()
()
,,R
b d x a
c x I f x y dxdy f x y dy dx
==⎰⎰⎰
⎰
（2.3.2）
计算积分（2.3.2）的近似值，我们对两个变量x 和y 都应用辛普森公式。

关于变量x ，取步长()h b a =-；关于变量y ，取步长
()()()
2
d x c x k x -=
它是x 的函数。

这样，我们有
()
()()()()()()()()()()()()()()()()()()()()()()()()()()()()()()()[,4,,]3{[,4,,]334[,4,3
,][,4,,]3
b
a k x I f x c x f x c x k x f x d x dx k a h f a c a f a c a k a f a d a k a h f a h c a h f a h c a h k a h f a h d a h k
b f b
c b f b c b k b f b
d b ≈+++≈+++++++++++++++++++⎰
下面，我们给出应用复化辛普森公式计算积分（2.3.2）的一种算法。

例2 应用复化辛普森公式计算积分
()()
()
,b
d x a
c x I dx f x y dy =⎰⎰
的近似值。

输入 断点a ，b;整数m ，n. 输出 I 的近似值S. Step 1 ()2h b a n ←-. Step 2 10;S ← 20;S ← 30.S ←
Step 3 对0,1,,2i n =⋅⋅⋅做
()()()()()()()()123;
2;
,,;0;0;
x a ih g d x c x m K f x c x f x d x K K ←+←-←+←←
对0,1,,21j m =⋅⋅⋅-做
()();,;
y c x jg z f x y ←+←
若j 是偶数，则22,K K z ←+
否则33;K K z ←+ ()12324p K K K g ←++ 若0i =或者2i n = 则11S S p ←+
否则，若i 是偶数，则22S S p ←+ 否则33S S p ←+ Step 4 ()12324 3.S S S S h ←++ Step 5 输出（S ）; 停机。

【】参考文献：林成森编著，数值计算方法 上册，科学出版社，2005年
3 辛普森公式的改进与推广
3.1辛普森公式的改进
3.1.1辛普森公式余项渐进性定理
对于定积分()b a
I f x dx =⎰，在满足有关的条件下，有辛普森公式
()()()46
s b a
I f a f c f b -=
++⎡⎤⎣⎦ （3.1.1） 以及公式的余项为
()()()5
42880
b a r f η-=-
； （3.1.2）
其中，[],,2
b a
c a b η+=
∈。

这里（3.1.2）中的点η为辛普森公式余项的“中间点”。

首先我们就是要证明辛普森公式余项的“中间点”具有当积分区间的长度趋于零时，“中间点”趋于区间的中的位置的性质。

命题1 设[]5,f C a b ∈，则存在[]1234,,,,a b ηηηη∈，满足
（1）()()()()()()5214
2
m
m
f f b f a f a b a b a η=+-+-， （3.1.3）
（2）()()()()4452f b f a f b a η=+-， （3.1.4） （3）()()()
()()()52
34
28
m
m
b a f f
c f
a f a
b a η-=
++
- （3.1.5） （4）()()()
()44542
b a f
c f a f η-=+ （3.1.6）
定理3.1 设[]5,f C a b ∈，若()()50f a ≠，则对辛普森公式余项的“中间点”η，
1
lim
2
b a a
b a η→-=
- （3.1.7）
证 由辛普森公式即可得到：
()()()()()()()5
4462880b
a
b a b a
f x dx f a f c f b f η---++=-⎡⎤⎣
⎦⎰
（3.1.8） 于是，就有
()()()()()
()()5
4288046b a b a f x dx f a f c f b b a f η-⎧⎫--++-=⎡⎤⎨⎬⎣⎦⎩⎭
⎰ （3.1.9）
()()()()()()()()()()()
()
()
5
6
444{288046}/b a b a f x dx f a f c f b f a b a b a f f a b a η-⎧⎫
--++⎡⎤⎨⎬⎣⎦⎩⎭---=--⎰ （3.1.10）
对（3.1.10）式b a →左端时的极限，连续4次运用L.洛比达法则，有： lim b a
→[（3.1.10）式左端]= ()()()lim{884m m b a
f b f c b a →-+--
()()()()()()()()2444
1}/34f b f c f a b a b a ⎡⎤++--⎢⎥⎣⎦
（3.1.11）
借助命题1，可知存在[]1234,,,,a b σσσσ∈，对于（3.1.11）式可有： lim b a
→[（3.1.11）式]= (
)
()()()()
()5551231lim[422
b a
f f f σσσ→-+-
()()()()2
2
544]/3f b a b a σ+-- (3.1.12)
注意到在定理的条件下，b a →必有i a σ→，()1,2,3,4i =则对（3.1.12） 就有
()()()()()
()(
)
()()
()
()
()()
()()
()
()()()
()5551232
2
54555551lim[422
4]/31lim[422
14]/3;2
b a
b a
f f f f b a b a f
a f
a f a f a f a σσσσ→→-+-
+--=-+-+=
（2.1.13） 另外，对（3.1.10）式右端求b a →时的极限。

则根据Lagrange 中值定理，存在
(),a λη∈
使得：
(
)
()()()()()()(
)
()4455lim
lim lim
b a
b a b a
f f a a f b a b a a
f a b a
ηηλη→→→--⎡⎤
=⋅⎢⎥--⎣⎦-=⋅- （3.1.14）
由于（3.1.13）=（3.1.14），就有1
lim
2
b a
a
b a
η→-=
-。

证毕。

3.1.2公式的改进 定义3.1 定义广义阶梯函数
()(),0,x t x t
x t x t
+⎧->-=⎨
≤⎩
（3.1.15） 命题2 设()()()
()56515!x a
R x x t f t dt =
-⎰，则 ()
()()()()46515!x a
R x x t f t dt =-⎰； （3.1.16）
其中[][]6,,,x a b f C a b ∈∈。

证略。

命题3 对[],t a b ∈，
()()()()()()()()6555
1[46!
1]0
4K t b t b a c t b a b t x t b a ++=
-------+--≤ （3.1.17）
定理3.2 设函数[]6,f C a b ∈，则对定积分()b
a
I f x dx =⎰，可有数值积分公式
()()()5
41,2880
s
b a I I f
c -=-
（3.1.18）
以及公式的余项
()()
()5
611120960b a r f η-=-； （3.1.19）
其中[],,2
b a
c a b η+=
∈。

证：依题可有()()()55f x T x R x =+，其中
()()()()()()()55
'
55!
f a T x f a f a x a x a =+-+⋅⋅⋅+- （3.1.20）
()()()
()56515!
x a R x x t f t dt =
-⎰； （3.1.21） 则按peano 核误差的方法，求积余项 ()()()()()()()
()()()()()
()()()()()()(
)
()
15
455555
4555555
4546
2880
46
2880
462880
b
a b a
b
a r I I
b a
f x dx f a f c f b b a f c b a
T x dx T a T c T b b a T c R x dx b a R a R c R b b a R c =--=-
++⎡⎤⎣⎦--
-=-
++⎡⎤⎣⎦-+
+--++⎡⎤⎣⎦-+⎰⎰⎰ 可算得
()()()()
()
()55555
4546
2880
b
a
b a
T x dx T a T c T b b a T c --
++⎡⎤⎣⎦-+
=⎰
由此可知本公式求积精度为5，比辛普森公式代数精度提高了2。

则借助定义3.1及命题2可以求得：
()()()()()()()()()()55555
46
5462880
b
a
b a b a
R x dx R a R c R b b a R c f t K t dt
--
++⎡⎤⎣
⎦-+=⎰
⎰ （3.1.22）
其中，
()()()()()()()()6555
1[46!
1]
4K t b t b a c t b a b t x t b a ++=
-------+-- （3.1.23）
由命题3，()K t 在[a,b]上不变号，则根据广义中值定理就可以得到 ()
()()()()
()5
6611120960
b
a b a r f K t dt f ηη-==-
⎰。

（3.1.24） 证毕。

举例
求()1
500.166666666666667I x dx ==⎰的数值解。

为了便于与辛普森公式，复化辛普森公式比较，有关的结果详见表1，表2.其中(),s s n I I 分别表示用辛普森公式，将积分区间等分为n 个小区间用辛普森公式对例子复化求积计算的结果，1I 表示定理7计算的结果。

表1 s I 与()s n I 的结果比较
表2 当n 取不同值时()s n I 的计算结果
【】参考文献：李毅夫著，辛普森公式“中间点”渐进性定理和辛普森公式的改进， 贵州师范大学学报（自然科学版），2007年11月
3.2 辛普森公式的推广
3.2.1 3次代数精度的条件
在之前我们都可以知道，辛普森公式具有3次代数精度。

辛普森公式实质利
用了,,2
a b
a b +这三点信息，从而我们可以求出一般情况下三点求积公式代数精度具有3次的条件。

已知的3个点01,x x 和2x ，其积分公式如下
()()()()0
1
1
2
2
b
a
f x dx A f x A f x A f x ≈++⎰ （3.2.1）
令()21,,f x x x =使（3.2.1）式准确成立，得
()()012
01222001122
22233
012,1,21.3A A A b a x A x A x A b a x A x A x A b a ⎧
⎪++=-⎪
⎪
++=-⎨⎪
⎪++=-⎪⎩
解方程组得
012012,,,D D D
A A A D D D
=
== 其中
()()01212
220
1201
22
2
2
3
322
11111
1
,2
13
b a D x x x D b a x x x x x b a x x -==--
()()()()02012
22
210
22012
332
22
3
311
1
111,22113
3
b a b a D x b a x D x x b a x b a x x x b a --=-=--- 具有3次代数精度的条件是()3f x x =是公式（3.2.1）准确成立，即满足
()012
333
4
401214
x A x A x A b a ++=- (3.2.2) 把0A ，1A ，2A 代入上式得到
()012
333
4
40121.4
x D x D x D b a D ++=- (3.2.3) 直接展开（2.2.3）式，计算比较复杂，构造下列的行列式，并利用行列式性质展开
()()()()()()()()()
()()()012012
00
120
12000
12333333
2
20122
01
2
2
123
32
3333332
20
122
222
2
212
12
3333
30
1
2
0121020212
201120210201
1
112131111
2
11113
12
x x x b a
x D x D x D x b a x x x x x b a x x x x x x b a x x x b a x x x x x x x x x b a x x x b a x x x x x x x x x b a x x x x x x x x x -++=--=--
-+
-=-----
-++-()()()()()()()0213
301210202113x x x b a x x x x x x x x x --+-++---
代人（3.2.3）式化简得到具有3次代数精度的条件是
()()()()()()2
2012011202335501212
11
34x x x b a x x x x x x b a x x x a b a b -
-+++-++=++（3.2.4）
3.2.2 公式的推广
很显然，辛普森公式就是012,,2
a b
x a x x b +===的情况，是满足（3.2.4）式的特例，
当0x a =时，代入（3.2.4）式得到1x 和2x 满足的条件
()()()22121211
22336
x x a b x x a ab b -+
++=++ （3.2.5） 当2x b =时，代入（3.2.4）式得到1x 和0x 满足的条件
()()()22101011
23236
x x a b x x a ab b -+
++=++ （3.2.6） 当12
a b
x +=
时，代入（3.2.4）式得到2x 和0x 满足的条件 02x x a b +=+ （3.2.7）
此情况比较简单，0x 和2x 是关于12a b x +=
对称，若02
a b
x h +=-，则22
a b
x h +=
+。

下面来推导其求积公式 ()012222b
a
a b a b a b f x dx A f h A f A f h +++⎛⎫⎛⎫⎛⎫
≈-+++ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭
⎰
（3.2.8）
其中012,,A A A 是待定系数，通过待定系数方法可以求得求积公式。

令()21,,f x x x =使（3.2.8）准确成立，得
()()01222
012
22233
012
,1,2
2221.2
223A A A b a a b a b a b h A A h A b a a b a b a b h A A h A b a ⎧
⎪
++=-⎪⎪+++⎪⎛⎫⎛⎫⎛⎫-+++=-⎨ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭⎪⎪+++⎛⎫⎛⎫⎛⎫⎪-+++=- ⎪ ⎪ ⎪⎪⎝⎭⎝⎭⎝⎭⎩ 解得方程组，可以得到012,,A A A 求积公式如下
()
()()()3
3
3
2
22
242122242b
a
b a b a b a a b a b a b f x dx f h b a f f h h h h ⎡⎤---+++⎛⎫
⎛⎫⎛⎫
≈
-+--++⎢⎥ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭
⎢⎥⎣⎦
⎰ （3.2.9）
公式（3.2.9）可以看作辛普森公式的推广，可以验证该推广公式具有3次代
数精度。

当2b a h -=时，推广的公式（3.2.9）就是辛普森公式，也就是说关于
2a b
+对称的两个与2
a b
+这三个点组成的积分公式可以达到 3次代数精度。

3.2.3 推广的简单应用
直接利用公式（3.2.9）可以推导一些求积公式，如求下列求积公式
()1
0120
111424f x dx A f A f
A f ⎛⎫
⎛⎫
⎛⎫≈++ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭
⎰
这里的a=0，b=1，h=
1
4
，代入公式（3.2.9）可以得到 ()3
02
2
2243b a A A h -===，()3
121
123
b a A b a h -=--=-
求积公式如下
()1
0211121343234f x dx f f f ⎛⎫⎛⎫⎛⎫≈-+ ⎪ ⎪ ⎪⎝⎭⎝⎭⎝⎭
⎰ 该公式具有3次代数精度。

【】参考文献：高尚著，辛普森公式的推广， 大学数学， 2007年4月
4 总结
辛普森公式在数值积分的地位是毋庸置疑的，有人曾经说过，数值分析的实际工作有95%的工作量可归结为运用辛普森法则和线性插值。

在这论文里，学习了辛普森公式的初级知识，都是在前人的基础上总结出来的。

在总结的部分，试着利用Excel 这个软件来，运用辛普森公式来计算积分的数值。

本章我们主要是通过计算0sin I xdx π
=⎰来说明辛普森公式的计算过程。

1，启动Excel ，建立新的工作文件。

2，在适当的单元区域先建立11个数据点：A 列的11数据从0等间距地变化到π，可从第二个数据（A5）开始按公式（“=A4+PI()/10”）计算，之后将公式复制到其他单元格，B 列为与A 对应的函数值，即对应x 值的正弦值sinx ，由此可以得到下面图2，图3：
图2 数据点的公式
图3 生成的数据点
3，在D列的对应位置输入公式，计算辛普森方法积分所需要的中间参数，并在适当的单元格内输入积分结果的计算公式，见图4，图5.
图4 辛普森公式的输入公式
图5 输出结果
4，求得的结果为2.00011，与数学方法计算的精确结果是一致的。

当然，利用Excel实现的过程中，还要其他方法，下面就来做一下其他形式的求解方法。

方法一，见图6，图7.
图6 方法一的计算公式
图7 方法一的计算结果
在这种求解方法中，11个数据点（x，y）的生成方法同前面的。

与前面介绍方法有所不同的是，D列为对应y值的一倍，只包括所以数据中两端的数据点；E列为对应y值的两倍，包括所有中间的数据点；F列同样为对应y值的两倍，但仅包括中间数据点中间隔的数据点。

在第16行的小计项，分别计算了D列，E列，F列所有数值的和。

三个累加的总和作为合计项存放在单元格D18内。

最后在单元格D20内输入公式“D18*PI()/10/3”,从而计算出积分结果。

最后的效果图7.与前面的计算结果完全一致。

方法二，见图8，图9.
在方法二中，应用了IF（）函数确定复化辛普森公式中不同函数前面的乘积系数（4或者2）。

分析复化辛普森公式我们可以发现，这个乘积系数是由数据点的次序决定的。

因此，在工作文件的C列为每个数据点进行编号。

在D列格单元格计算公式中，这些编号作为IF（）函数的判别条件。

图8 方法二的计算公式
图9 方法二的计算结果
这种处理形式的计算效果如图9.计算结果与前面的完全一致。

这里要说明一下，函数说明MOD（number，divisor）
功能：返回两个数相除的余数，结果的正负号与除数相同。

参数说明：
Number为被除数；
Divisor为除数。

函数类别：数学与三角函数
由上面的几个方法，根据着辛普森公式的原理，利用Excel方便快捷地求解了积分问题。

与编程程序求解相比，Excel方法更加简便，更加实用。

【】参考文献：谢剑，赵彤编著，Excel在建筑工程中的应用快速计算工具实例，天津大学出版社，2004年
本文通过对各种相关文献和资料的收集和整理，经过认真的阅读和细致的分析，综合分析了辛普森公式的有关知识，通过复化辛普森公式和二重积分的辛普森公式来完善了辛普森公式的内容，但是这也只是辛普森公式的知识的一小部分。

论文虽然经过多次修改，但由于个人能力有限，经验不足，论文中还存在一些不足之处，希望各位老师同学批评指正。