C语言实现龙贝格求积

#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;

}