龙贝格积分实验报告.docx

二、Romberg积分法
1.变步长Romberg积分法的原理
复化求积方法对于提高精度是行之有效的方法，但复化公式的一个主要缺点在于要事先估计出部长。若步长过大，则精度难于保证；若步长过小，则计算量又不会太大。而用复化公式的截断误差来估计步长，其结果是步长往往过小，而且和在区间上的上界的估计是较为困难的。在实际计算中通常采用变步长的方法，即把步长逐次分半（也就是把步长二等分），直到达到某种精度为止，这种方法就是Romberg积分法的思想。
在步长的逐步分半过程中，要解决两个问题：
1. 在计算出后，如何计算,即导出和之间的递推公式；
2. 在计算出后，如何估计其误差，即算法的终止的准则是什么。
首先推导梯形值的递推公式，在计算时，需要计算个点处的函数值在计算出后，在计算时，需将每个子区间再做二等分，共新增个节点。为了避免重复计算，计算时，将已计算的个点的数值保留下来，只计算新增个节点处的值。为此，把表示成两部分之和，即

由此得到梯形值递推公式

因此

由复化梯形公式的截断误差有

若变化不大时，即，则有

式（2）表明，用作为定积分的近似值，其误差大致为，因此其终止条件为

其中是预先给定的精度。
2.Romberg积分公式
将上述方法不断推广下去，可以得到一个求积分的序列，而且这个序列很快收敛到所求的定积分。记
，将区间等分的梯形值。，将区间等分的Simpson
，将区间等分的Cotes。，将区间等分的Romberg。
由其可构造一个序列，次序列称为Romberg序列，并满足如下递推关系：

以上递推公式就是Romberg积分递推公式。
3.Romberg积分程序
1. 置，精度要求，；
2. 计算；
3. 置，并计算；
4. 置
5. 计算；
6. 若 ，则转（7）；否则置，转（5）；
7. 若，则停止计算（输出），否则转（3）。
4.Romberg积分法的应用
function [T,n] = romb(f,a,b,eps)
double R;
if nargin<4,eps=1e-8;
end
h=b-a;R(1,1)=(h/2)*(feval(f,a)+feval(f,b));
n=1;J=0;err=1;
while (err>eps)
J=J+1;h=h/2;S=0;
for i=1:n
x=a+h*(2*i-1);
S=S+feval(f,x);
end
R(J+1,1)=R(J,1)/2+h*S;
for k=1:J
R(J+1,k+1)=(4^k*R(J+1,k)-R(J,k))/(4^k-1);
end
err=abs(R(J+1,J+1)-R(J+1,J));
n=2*n;
end
R;
T=R(J+1,J+1)
End
其中输入项：f为被积函数，ab为积分区间的端点值，ep为积分精度；输出项：T是逐次积分表值，n是迭代次数，R是最后积分值。
4.1程序调用
可以将被积分函数编成函数文件，也可以直接使用内联函数来表示被积分函数，示例如下:

>>f=inline('1/(1+x.^2)','x');
>> [T,n,R]=romb(f,2,9,1e-9)
运行后得出其迭代次数，最终积分结果以及龙贝格积分矩阵如表2-1所示，
迭代次数N=64，最终的积分值R=0.3530.
0.7427
0.0000
0.0000
0.0000
0.0000
0.0000
0.0000
0.4833
0.3969
0.0000
0.0000
0.0000
0.0000
0.0000
0.3905
0.3596
0.3571
0.0000
0.0000
0.0000
0.0000
0.3628
0.3536
0.3532
0.3532
0.0000
0.0000
0.0000
0.3555
0.3530
0.3530
0.3530
0.3530
0.0000
0.0000
0.3536
0.3530
0.3530
0.3530
0.3530
0.3530
0.0000
0.3531
0.3530
0.3530
0.3530
0.3530
0.3530
0.3530
表2-1 龙贝格积分矩阵
3.课本例题求解

1 当迭代精度ep=1e-9的条件下，迭代次数N=32，迭代结果R=0.6931
表2-2 式1对应的龙贝格积分矩阵
0.7500
0.0000
0.0000
0.0000
0.0000
0.0000
0.7083
0.6944
0.0000
0.0000
0.0000
0.0000
0.6970
0.6933
0.6932
0.0000
0.0000
0.0000
0.6941
0.6932
0.6931
0.6931
0.0000
0.0000
0.6934
0.6931
0.6931
0.6931
0.6931
0.0000
0.6932
0.6931
0.6931
0.6931
0.6931
0.6931
2 当迭代精度ep=1e-9的条件下，迭代次数N=32，迭代结果R=0.2722.
表2-3 式2对应的龙贝格积分矩阵
0.1733
0.0000
0.0000
0.0000
0.0000
0.0000
0.2488
0.2740
0.0000
0.0000
0.0000
0.0000
0.2665
0.2723
0.2722
0.0000
0.0000
0.0000
0.2708
0.2722
0.2722
0.2722
0.0000
0.0000
0.2718
0.2722
0.2722
0.2722
0.2722
0.0000
0.2721
0.2722
0.2722
0.2722
0.2722
0.2722
3对于积分，由于积分下限0为其奇点，理论上无法进行数值积分，本题中近似取下限为1*10-9来进行计算。当迭代精度ep=1e-9的条件下，迭代次数N=16，迭代结果R=0.2722.

表2-4 式3对应的龙贝格积分矩阵
0.8466
0.0000
0.0000
0.0000
0.0000
0.8288
0.8228
0.0000
0.0000
0.0000
0.8241
0.8225
0.8225
0.0000
0.0000
0.8229
0.8225
0.8225
0.8225
0.0000
0.8226
0.8225
0.8225
0.8225
0.8225
4.对于积分，同样积分下限0为积分函数的奇点，理论上无法进行数值积分运算，本题中仍取积分下限近似为1*10-9进行计算。当迭代精度ep=1e-9的条件下，迭代次数N=16，迭代结果R=1.3708.
表2-5 式4对应的龙贝格积分矩阵
1.2854
0.0000
0.0000
0.0000
0.0000
1.3498
1.3713
0.0000
0.0000
0.0000
1.3655
1.3708
1.3708
0.0000
0.0000
1.3695
1.3708
1.3708
1.3708
0.0000
1.3704
1.3708
1.3708
1.3708
1.3708