龙贝格算法就积分的c语言程序

C++实现如下：

#include<iostream>困态

#include<cmath>

using namespace std

const int MAXRepeat = 100 //最大允许重毕尺顷复

double function(double x)//被积函数,根据自己的需要手工输入

{

double s

s = 1.0 / (1 + x)

return s

}

void Romberg(double a, double b, double epsion, double f(double x))

{

int m = 1

int n = 1

int k

double h

double ep

double p

double xk

double s

double q

double T[MAXRepeat]

h = b - a

T[0] = 0.5 * h * (f(a) + f(b))

ep = epsion + 1.0

while ((ep >= epsion) &&(m <MAXRepeat))

{

p = 0.0

for (k = 0k <nk++)

{

xk = a + (k + 0.5) * h// n-1

p = p + f(xk) //计算∑f(xk+h/2),T

} // k=0

p = (T[0] + h * p) /手陆 2.0 //T`m`(h/2),变步长梯形求积公式

s = 1.0

for (k = 1k <= mk++)

{

s = 4.0 * s //[pwww.hbbz08.com ow(4,m)T`m`(h/2)-T`m`(h)]/[pow(4,m)-1],2m阶牛顿柯斯特公式，即龙贝格公式

q = (s * p - T[k - 1]) / (s - 1.0)

T[k-1] = p

p = q

}

ep = fabs(q - T[m - 1])

m++

T[m - 1] = q

n++// 2 4 8 16

h /= 2.0

}

for (int i = 0i <mi++)

{

int j

if (!(i % j))

{

cout<<T[i]<<endl

}

else

{

cout<<T[i]<<" "

}

j++

这是我在乱纤vc++6.0上测试你的程序输出的结果，我不知道这个结果对不对，但至少不是弊陪岁0，我建议你在定义double T[limit][limit]；把T数组初始租睁化一下。

#include<iostream.h>

#include<math.h>

# define Precision 0.000001//粗察积分精度要求

# define e 2.71828183

#define MAXRepeat 10 //最大允许重复

double function(double x)//被积函数

{

double s

s=pow(e,x)*cos(x)

return s

}

double Romberg(double a,double b,double f(double x))

{

int m,n,k

double y[MAXRepeat],h,ep,p,xk,s,q

h=b-a

y[0]=h*(f(a)+f(b))/2.0//计算T`1`(h)=1/2(b-a)(f(a)+f(b))

m=1

n=1

ep=Precision+1

while((ep>=Precision)&&(m<MAXRepeat))

{

p=0.0

for(k=0k<nk++)

{

xk=a+(k+0.5)*h// n-1

p=p+f(xk) //计算∑f(xk+h/2),T

} // k=0

p=(y[0]+h*p)/2.0 //T`m`(h/2),变步长梯形求积信纳公式

s=1.0

for(k=1k<=mk++)

{

s=4.0*s// pow(4,m)

q=(s*p-y[k-1])/(s-1.0)//[pow(4,m)T`m`(h/2)-T`m`(h)]/[pow(4,m)-1],2m阶牛顿柯斯特公式，即龙贝格公岩坦茄式

y[k-1]=p

p=q

}

ep=fabs(q-y[m-1])//前后两步计算结果比较求精度

m=m+1

y[m-1]=q

n=n+n // 2 4 8 16

h=h/2.0//二倍分割区间

}

return q

}

main()

{

double a,b,Result

cout<<"请输入积分下限:"<<endl

cin>>a

cout<<"请输入积分上限:"<<endl

cin>>b

Result=Romberg( a, b, function)

cout<<"龙贝格积分结果："<<Result<<endl

return 0

}

本文来自CSDN博客，转载请标明出处：http://blog.csdn.net/liuhongyi666/archive/2008/12/23/3589002.aspx

---