n 1 b a1 ' ' ' ' 利 用 h 和 f ( ) f ( ) a , b k n n k 0


得 到 复 化 梯 形 公 式 的 截 断 误 差 是 ： b a 2 '' 2 R ( T ) hf ( ) O ( h) n 1 2
复化Simpson公式数值实验结果
n=8, n=16, n=32, n=64, n=128, n=256, n=512, n=1024, n=2048, n=4096, n=8192, sp=1.85195370083255 sp=1.85193808500560 sp=1.85193711642985 sp=1.85193705600861 sp=1.85193705223407 sp=1.85193705199819 sp=1.85193705198345 sp=1.85193705198253 sp=1.85193705198247 sp=1.85193705198246 sp=1.85193705198247
数值试验结果

用复化Simpson公式求解，当区间分成64份 时，计算结果： sp=1.85193705600861 比用复化梯形公式，当区间分成8192份 时的计算结果： tp=1.85193704808136 还要精确。
第 二 节 复化求积公式
一、复化求积公式 复化求积公式的基本思想： 将区间[a , b] 分为若干个小子区间，在每个 小子区间上使用低阶的Newton-Cotes公式。然后
把它们加起来，作为整个区间上的求积公式。
1、复化梯形公式
将 区 间 等 分 , a,bn ba h ,x h , (k 0 ,1 , ,n ), k ak n 在 每 个 小 区 间 k 0 ,1 , ,n1 ） xk, xk1,（ 上 用 梯 形 公 式 :
1 0.9 0.8 f(x)=sin(x)/x 0.7 0.6 0.5 0.4 0.3 0.2
x=eps:0.01:pi; y=sin(x)./x; plot(x,y); legend('f(x)=sin(x)/x');
0.1 0
0
0.5
1
1.5
2
2.5
3
3.5
数值试验

复化梯形公式Matlab程序
复化Simpson公式的截断误差为
( b a ) 4( 4 ) 4 R () S h f( ) O () h a , b n 2 8 8 0
Example 1
Approximate the integral sin x dx x 0 Using the Composite Trapezoidal rule and Composite Simpson’s rule
h T (( f x ) f ( x ) ) k 0 , 1 ,, n 1 k k k 1 2
复化梯形公式为
n 1 h T T ( fa () fb () ) h fx (k ) n k 2 k 0 k 1 n 1
截断误差分析：
3 h ' ' 在 区 间 x , x 上 ， R f () , k x , x k 1 kk 1 k k k 1 2 3 n 1 n 1 h ' ' 整 体 误 差 为 R R ( ) f( ) n k k 2 k 0 k 0 1
2、复化Simpson公式
在 每 个 小 区 间 ,x 上 用 S i m p s o n 公 式 x k k 1 h S f( x ) 4 f( x 1) f( x ) ) k ( k k 1 k 6 2
复化Simpson公式为
1 1 h 2n 1n S S f( a ) f( b ) ) h f( x 1) h f( x ) n k ( k k 6 3k 3k k 0 0 1 2 n 1 1 2 T , 其 中 H h f( x 1) n H n n k 3 3 k 0 2 n 1
数值试验
在Matlab命令窗口中，进行如下操作：

>> format long >>f =inline( 'sin(x)/x' ) ; >> a = eps; b =pi; n = 8; % eps 是Matlab最小正数 >> sp = simpson(f, a, b, n)
数值试验
function rs= trapezoid(f,a,b,n) h = (b-a)/n; r= (feval(f,a)+feval(f,b))/2 ; for j = 1:n-1 x=a+j*h ; r= r+ feval(f,x); end rs = r* h; 将此程序存于work 目录中

数值试验
数值试验

复化Simpson公式Matlab程序
function rs= simpson(s,a,b,n) h = (b-a)/n; r= feval(s,a)+feval(s,b); for j = 1:2:n-1 x=a+j*h ; r= r+ 4*feval(s,x); end for j = 2:2:n-2 x=a+j*h ; r= r+ 2*feval(s,x); end 将此程序存于work目录中 rs = r*h/3;
复化梯形公式数值实验结果



ቤተ መጻሕፍቲ ባይዱ



n=8, n=16, n=32, n=64, n=128, n=256, n=512, n=1024, n=2048, n=4096, n=8192,
tp=1.84784230644461 tp=1.85091414036536 tp=1.85168137241373 tp=1.85187313510989 tp=1.85192107295303 tp=1.85193305723690 tp=1.85193605329681 tp=1.85193680231110 tp=1.85193698956462 tp=1.85193703637800 tp=1.85193704808136