龙贝格求积法

#include<iostream>using namespace std;#include"math.h"#define e  1e-5double  f(double x){    double y;    if (x==0)        return y=1.0;//f(0)=1    else y=sin(x)/x;    return y;}void romberg(double a,double b){    int n=1,k=0,i;    double h,T2,S2=0,C2=0,R2=0,T1,C1,S1,R1;    double x;    h=b-a;    T2=h/2*(f(a)+f(b));    do//在要求精度之下    {        x=a+h/2;        R1=R2;        T1=T2;        S1=S2;        C1=C2;        double sum=0;        for(i=0;i<=pow(2,k)-1;i++)           {               sum=sum+f(x);///都是公式               x+=h;           }        T2=(T1+sum*h)/2;        S2=(4*T2-T1)/3;        C2=(16*S2-S1)/15;        R2=(64*C2-C1)/63;        k++;        h=h/2;    }while (fabs((R2-R1))>e);    cout<<R2<<endl;}int main(){    double a,b;    cout<<"输入上下限a,b的值用空格隔开"<<endl;    cin>>a>>b;    cout<<"下限a="<<a<<endl;    cout<<"上限b="<<b<<endl;    cout<<"被积函数为：y=sin(x)/x"<<endl;    romberg(a,b);}