龙贝格求积

#include <stdio.h>
#include <stdlib.h>
#include<conio.h>
#include<math.h>
#include<string.h>

double f(double x)
{
return (4*x*x/(1+x*x));
}

double Romberg(double top,double bottom,double precision)
{
int k=1;
double S,x,T1,T2,S1,S2,C1,C2,R1,R2,h=bottom-top;

/*
转载请注明出处：去转盘网www.quzhuanpan.com
*/

//S梯形公式，T梯型变步长，S梯形加速，C幸普森加速，R龙贝格求积
T1=h*(f(top)+f(bottom))/2;//梯形公式
while(1)
{
S=0;
x=top+h/2;
do
{
S+=f(x);
x+=h;

}while(x<bottom);
T2=(T1+h*S)/2.0;
if(fabs(T2-T1)<precision)
{
return T2;
}
S2=T2+(T2-T1)/3;//梯形加速
if(k==1)
{
T1=T2;
S1=S2;//几下原来的S2
h/=2;
k+=1;
continue;
}
C2=S2+(S2-S1)/15;//新的S2减原来的S2，即使S1,幸普森加速
if(k==2)
{
C1=C2;
T1=T2;
S1=S2;
h/=2;
k+=1;
continue;
}
R2=C2+(C2-C1)/63.0;
if(k==3)
{
R1=R2;
C1=C2;
T1=T2;
S1=S2;
h/=2;
k+=1;
continue;
}
if(fabs(S2-S1)<precision)
{
return S2;
}
C1=C2;
T1=T2;
S1=S2;
h/=2;
k+=1;
if(fabs(R2-R1)<precision);
{
return R2;
}
}
}

int main()
{
double top,bottom,precision,S;
printf("Pletopse input the begin: ");
scanf("%lf",&top);
printf("Pletopse input the end: ");
scanf("%lf",&bottom);
printf("Pletopse input the precision:");
scanf("%lf",&precision);
S=Romberg(top,bottom,precision);
printf("The result is:%lf",S);
getch();
return 0;
}