CENTRAL SOUTH UNIVERSITY 数值分析实验报告
变步长的梯形积分方法的应用
一、问题背景
实际问题当中常常需要计算积分，有些数值方法，如微分方程和积分方程的求解，也都和积分计算相关联。
依据人们所熟知的微积分基本定理，对于积分
()dx x f I a ⎰=b
， 只要找到被积分函数()x f 的原函数()x F ，便有下列牛顿-莱布尼茨（Newton-Leibniz ）公式：
()()()a F b F dx x f b
-=⎰a . 但实际使用这种求积分方法往往有困难，因为大量的被积函数，诸如
()0sin ≠x x
x ，2x e -等，其原函数不能用初等函数表达，故不能用上述公式计算。
即使能求得原函数的积分，有时计算也十分困难。
例如对于被积函数
()6
11x x f +=，其原函数 ()C x x x x x x x x F ++-+++⎪⎭⎫ ⎝⎛-+=1
313ln 3411arctan 61arctan 3122， 计算()a F ，()b F 仍然很困难，另外，当()x f 是由测量或数值计算给出的一张数据表时，牛顿-莱布尼茨公式也不能直接运用。
因此有必要研究积分的数值计算问题。
二、数学模型
由于牛顿-科特斯积分公式在8≥n 时不具有稳定性，故不能通过提高阶数的方法来提高求积精度。
为了提高精度通常可以把积分区间划分成若干的子区间（通常是等分），再在每个子区间上用低阶求积公式。
这种方法称为复合求积法。
复合梯形法虽然方法简单，但是却不能估计积分精度，这有时候是很不方便的。
要想控制积分精度，可以采用如下的方法，设积分区间已经划分为n 个子区间，这时再把区间划分更细，给出新的积分结果，如果前后两次积分的差比给定的误差容限小的话，则停止细华否则继续增加积分区间。
这种方法原理很简单也 容易实现，但是实际计算中一般采用的比较少，因为这种方法比较机械效率不是太高，实际上采用比较多的通常是Romberg 方法。
三、算法及流程
给定义误差容限小量TOL ，对于()dx x f b
a ⎰，有复合梯形公式
()()()()⎥⎦
⎤⎢⎣⎡++≈∑⎰-=1122n k k b
a x f
b f a f h dx x f 如果前后两次的划分的积分计算结果大于给定的误差TOL ，则增加划分区间，如果满足精度，则停止细化，并输出结果。
MATLAB 实现过程：
%可控精度复合梯形法计算积分问题
function [jifen,num]=kong_tixing(a,b,tol)
%a,b 为积分区间
%tol 为积分精度，默认为10的-3次方
if (nargin==3)
eps=1.0e-3;
end
n=1;
h=(b-a)/2;
int_1=0;
%调用方程函数
int_2=(kong_t_f(a)+kong_t_f(b))/h;
% 如果前后两次误差小于给定的精度，则停止细化积分区间
while abs(int_2-int_1)>tol
n=n+1;
h=(b-a)/n;
int_1=int_2;
int_2=0;
for i=0:n-1
x=a+h*i;
x1=x+h;
int_2=int_2+(h/2)*(kong_t_f(x)+kong_t_f(x1));
end
end
% 积分结果
jifen=int_2;
%区间划分细度
num=n;
将文件以文件名kong_tixing.m 保存。
四、计算结果及分析
计算定积分
dx e T x ⎰-=102
要求输出精度为10-4。
打开Editor 编写如下程序，并将文件以文件名kong_t_f.m 保存。
function f=kong_t_f(x)
f=exp(-x^2);
打开Editor编写如下程序，并将文件以文件名kong_t_main.m保存。
% 可控精度的梯形积分方法
% 精度为0.1
[s_1,num_1]=kong_tixing(0,1,1e-1)
%精度为0.01
[s_2,num_2]=kong_tixing(0,1,1e-4)
% 精度为0.001
[s_1,num_3]=kong_tixing(0,1,1e-7)
% 画出积分图形
x=0:0.02:1;
y=exp(-x.^2);
area(x,y)
在MATLAB命令窗口输入输入
>>kong_t_main
点击Enter键后得到：
s_1 = 0.739986475276682
num_1 =3
s_2 =0.746398247893441
num_2 = 12
s_1 =0.746818876175303
num_3 =108
从输出的结果可以看出，要达到10-4精度，需要把区间划分为12个子区间，而要达到10-7精度，则要把区间划分为108个子区间。
事实上，函数的积分总区间跨度不是很大，所以在划分为108个子区间后已经是对函数取点很密集了。
下图给出了积分的几何意义，积分的结果即图中蓝色区域面积，从输出结果可以看出蓝色区域面积约为0.746818876175303。