Ronberg算法计算积分

语法：

result=integral(double a,double b);

参数：

a：积分上限

b：积分下限

function：

f（）积分函数

返回值：

f在a、b之间的返回值

注意：

function f(x)需要自行修改。程序中用的是 x*x；

需要 math.h;

默认精度要求是

1e-5

源程序：

#include<iostream>
#include<math.h>
#define epsilon 0.00001
#define COUNT 100
using namespace std;

double fun(double x)
{
    return x*x;
}

double Romberg(double a,double b)
{
    int m ,n;
    double h,x,s,q,ep;
    double p,*R =new double[COUNT];

            h=b-a;
            R[0]= h*(fun(a)+ fun(b))/2.0;
            m=1;
            n=1;
            ep=epsilon+1.0;
            while ((ep >= epsilon)&& (m <COUNT))
            {
                p = 0.0;
            {
                for(int i=0;i<n;i++)
                {
                    x = a+ (i+0.5)*h ;
                    p= p + fun(x);
                }
                    p= (R[0]+ h*p)/2.0;
                    s = 1.0;
                    for(int k=1;k<=m;k++)
                    {
                        s = 4.0*s;
                        q= (s*p-R[k-1])/(s-1.0);
                            R[k-1]= p;
                            p =q;
                    }
                    p=fabs(q -R[m-1]);
                    m =m + 1;
                    R[m-1]= q;
                    n = n + n;
                    h = h/2.0;
                }
                return (q);
            }
}

int main()
{
    double a,b;
    cout<<"Input a，b：a为下限，b为上限"<<endl;
    cin>>a>>b;
    cout<<"Romberg="<<Romberg(a,b)<<endl;
    system("pause");
    return 0;
}