龙贝格求+积分.doc

龙贝格(Romberg),应用Richardson外推法导出的数值求积方法。由复化梯形公式   可以化为=一般地,把区间[a,b]逐次分半k-1次,(k=1,2,……,n)区间长度(步长)为,其中mk=2k-1。 记=由=  从而=-  (1)按Richardson外推思想,可将(1)看成关于,误差为的一个近似公式,因而,复化梯形公式的误差公式为-==   (2)取=有  -=  (3)误差为的误差公式  =+(1)误差,对复化梯形公式误差估计时,是估计出每个子区间上的误差,然后将n个子区间上的误差相加作为整个积分区间上的误差。(2)收敛性,记,由于=上面两个累加式都是积分和,由于在区间上可积可知,只要的分划的最大子区间的长度时,也即时,它们的极限都等于积分值。可见,只要在区间上的可积的条件满足,由复化梯形求积公式计算所得复化梯形值序列都收敛于实际积分值。:(1)、输入:a,b,epsilon(2)、令=b-a,计算=(3)、令k=2,=(4)、令=,计算(5)、for j=2,3,……,k(6)、if   {输出; return;}(7)、k=k+1;  =;  跳转(4)(8)、=://RombergmethodforIntegral//////Date://ByHuBentao#include""#include""#include""#include""#include""#include""#defineN120#defineN220doublefun(doublex){//被积函数设置return(8/(1+x*x));}voidRomberg(doublea,doubleb,double(*fun)(double),doubleepsilon){inti,j,k;doubleh=b-a,temp;doubleT[20][20];FILE*fp;if((fp=fopen("","w"))==NULL){//("\nopenfileerror!\n");return;}fflush(stdin);fprintf(fp,"\n龙贝格求积计算的中间结果:");printf("\n龙贝格求积计算的中间结果:");T[1][1]=h*((*fun)(a)+(*fun)(b))/2;fprintf(fp,"\nT[1][1]=%f",T[1][1]);printf("\nT[1][1]=%f",T[1][1]);k=2;h/=2;while(1){temp=0;for(i=1;i<=pow(2,(k-2));i++){temp+=(*fun)(a+(2*i-1)*h);}T[k][1]=temp*h+T[k-1][1]/2;fprintf(fp,"\nT[%d][1]=%f\t",k,T[k][1]);printf("\nT[%d][1]=%f\t",k,T[k][1]);for(j=2;j<=k;j++){T[k][j]=T[k][j-1]+(T[k][j-1]-T[k-1][j-1])/(pow(4,j-1)-1);fprintf(fp,"T[%d][%d]=%f\t",k,j,T[k][j]);printf("T[%d][%d]=%f\t",k,j,T[k][j]);}j--;if((fabs(T[k][j]-T[k][j-1]))<epsilon){fprintf(fp,"\n积分结果\tI(f)=");//cout<<endl<<"OK!!"<<endl<<"T"<<k<<j<<"="<<T[k][j];fprintf(fp,"%f\t",T[k][j]);//积分值输出printf("\n积分结果\tI(f)=");//cout<<endl<<"OK!!"<<endl<<"T"<<k<<j<<"="<<T[k][j];printf("%f\t",T[k][j]);//积分值输出break;}k+=1;h/=2;}fclose(fp);}//:Definestheentrypointfortheconsoleapplication.#include""#include""#include""#include""#includ