利用Matlab实现Romberg数值积分算法
一、内容摘要
针对于某些多项式积分，利用Newton—Leibniz积分公式求解时有困难，可以采用数值积分的方法，求解指定精度的近似解，本文利用Matlab中的.m文件编写了复化梯形公式与Romberg的数值积分算法的程序，求解多项式的数值积分，比较两者的收敛速度。
二、数值积分公式
1.复化梯形公式求解数值积分的基础是将区间一等分时的Newton—Cotes
求积公式:
I=?f(x)dx?abb?a[f(a)?f(b)] 2其几何意义是，利用区间端点的函数值、与端点构成的梯形面积来近似f(x)在区间[a,b]上的积分值，截断误差为：
?(b?a)3\f(?) ??(a,b) 12具有一次的代数精度，很明显，这样的近似求解精度很难满足计算的要求，因而，可以采用将积分区间不停地对分，当区间足够小的时候，利用梯形公式求解每一个小区间的积分近似值，然后将所有的区间加起来，作为被求函数的积分，可以根据计算精度的要求，划分对分的区间个数，得到复化梯形公式：
n?1(b?a)k(b?a)f(x)dx?[f(a)?f(b)?2?f(a?)]
2nnk?1I=?ba其截断误差为：
R??(b?a)2\hf(?) ??(a,b) 12数值积分算法
使用复化的梯形公式计算的数值积分，其收敛速度比减慢，为此，采用Romberg数值积分。其思想主要是，根据I的近似值T2n加上I与T2n的近似误差，作为新的I的近视，反复迭代，求出满足计算精度的近似解。
用T2n近似I所产生的误差可用下式进行估算：
1??I?T2n?(T2n?T2n?1)
3新的I的近似值：
jj=（0 1 2 ….） T2n?1???Tn 2Romberg数值积分算法计算顺序
i=0 i=1 i=2 i=3 i=4
(1) (2) (4) (7) (11) …
T200 T201
(3) (5) (8) (12) …
T210 T211
(6) (9) (13) …
T220 T221
(10) (14) …
T230 T231
T202 T203 T204
T212 T213
T222
其中，第一列是二阶收敛的，第二列是四阶收敛的，第三列是六阶收敛的，第四列是八阶收敛的，即Romberg序列。
三、复化梯形法以及Romberg算法程序流程图
输入被积函数、积分区间、收敛条件计算Ti0计算Ti+10差值是否满足收敛条件是否输出结果 图1 复化梯形法程序流程图
输入被积函数、积分区间、收敛条件计算并存储Ti0计算并存储Ti+10是差值是否满足收敛条件否构建并存储T1i+1否i是否大于等于2是差值是否满足收敛条件否构建并存储T2i否i是否大于等于3是差值是否满足收敛条件否构建并存储T3i否i是否大于等于4是否差值是否满足收敛条件是是是输出结果
图2 Romberg算法程序流程图
四、计算实例
依据上文所述的流程图，编写复化梯形程序以及Romberg算法程序，并且利用实例验证程序的正确性，示例如下(计算精度)：
???4dx 01?x21表2 计算结果
计算精度 复化梯形 算法 Romberg 算法
×10^-5 时间 近似值 时间 近似值
3. 3.
×10^-7
3. 3.
×10^-9
3. 3.
从上表中可以看出，当要求的计算精度不高时，复化梯形算法与Romberg算法计算时间相差不太大，但是Romberg算法是要快于复化梯形算法的；当要求的计算精度更高的时候，Romberg算法是明显快于复化梯形算法。
本文所编写的程序适用于多项式的数值积分，且对于积分区间内，被积函数在每一点必须有定义，在以后的学习中进一步改进。
利用Matlab实现Romberg数值积分方法
