龙贝格求 积分

1、龙贝格（Romberg）求积法1.算法理论Romberg求积方法是以复化梯形公式为基础，应用Richardson外推法导出的数值求积方法。由复化梯形公式 可以化为 一般地，把区间a,b逐次分半k1次，（k1,2,n）区间长度（步长）为，其中mk2k1。 记由 从而 （1）按Richardson外推思想，可将（1）看成关于，误差为的一个近似公式，因而，复化梯形公式的误差公式为 （2）取有 （3）误差为的误差公式 2.误差及收敛性分析（1）误差，对复化梯形公式误差估计时，是估计出每个子区间上的误差，然后将n个子区间上的误差相加作为整个积分区间上的误差。（2）收敛性，记，由于 上面两个累加式都是积分

2、和，由于在区间上可积可知，只要的分划的最大子区间的长度时，也即时，它们的极限都等于积分值。可见，只要在区间上的可积的条件满足，由复化梯形求积公式计算所得复化梯形值序列都收敛于实际积分值。3.算法：（1）、输入：a，b，epsilon（2）、令ba，计算（3）、令k2，（4）、令，计算（5）、for j2,3,k （6）、if 输出; return;（7）、kk1; ; 跳转（4）（8）、结束4.实例用龙贝格算法计算积分值I=5.龙贝格求积源程序：/Romberg method for Integral/Projectname hu1.cpp/Executable file hu1.exe/Da

3、te:2004.12/By Hu Bentao#include stdafx.h#include stdio.h#include iostream.h#include math.h#include conio.h#include stdlib.h#define N1 20#define N2 20double fun(double x)/被积函数设置return (8/(1+x*x);void Romberg(double a,double b,double(*fun)(double),double epsilon)int i,j,k;double h=b-a,temp;double T202

4、0;FILE *fp;if(fp=fopen(200402157.txt,w)=NULL)/将运算的中间结果和最终结果保存到文件200402157.txtputs(nopen file error!n);return ;fflush(stdin);fprintf(fp,n龙贝格求积计算的中间结果：);printf(n龙贝格求积计算的中间结果：);T11=h*(*fun)(a)+(*fun)(b)/2;fprintf(fp,nT11=%f,T11);printf(nT11=%f,T11);k=2;h/=2;while(1)temp=0;for(i=1;i=pow(2,(k-2);i+)temp+

5、=(*fun)(a+(2*i-1)*h);Tk1=temp*h+Tk-11/2;fprintf(fp,nT%d1=%ft,k,Tk1);printf(nT%d1=%ft,k,Tk1);for(j=2;j=k;j+)Tkj=Tkj-1+(Tkj-1-Tk-1j-1)/(pow(4,j-1)-1);fprintf(fp,T%d%d=%ft,k,j,Tkj);printf(T%d%d=%ft,k,j,Tkj);j-;if(fabs(Tkj-Tkj-1)epsilon)fprintf(fp,n积分结果tI(f)=);/coutendlOK!endlTkj=Tkj;fprintf(fp,%ft,Tkj)

6、;/积分值输出printf(n积分结果tI(f)=);/coutendlOK!endlTkj=Tkj;printf(%ft,Tkj);/积分值输出break;k+=1;h/=2;fclose(fp);/ hu1.cpp : Defines the entry point for the console application.#include stdafx.h#include stdio.h#include iostream.h#include math.h#include conio.h#include stdlib.hextern void Romberg(double a,double b,double(*fun)(double),double epsilon);extern double fun(double x);int main(int argc, char* argv)double a,b,epsil