数学之家

标题: 龙贝格求积法C语言程序 [打印本页]

作者: EMP震荡波 时间: 2011-5-13 16:57
标题: 龙贝格求积法C语言程序
#include<stdio.h>
#include<math.h>
void main()
{
int i,k,m;
double a,b,h,g=0.0,t[20][20],e=0.001;
double f(double x);
scanf("%lf%lf",&a,&b);
h=b-a;
t[0][0]=h*(f(a)+f(b))/2.0;
for(k=1;;k++)
{
  for(i=1;i<=pow(2,k-1);i++)
g+=f(a+(i-1.0/2.0)*h);
  t[k][0]=(t[k-1][0]+h*g)/2.0;
  for(m=1;m<=k;m++)
t[k-m][m]=(pow(4,m)*t[k+1-m][m-1]-t[k-m][m-1])/(pow(4,m)-1);
  if((t[0][k]-t[0][k-1]<=e)&&(t[0][k-1]-t[0][k]<=e))
  {
printf("%f\n",t[0][k]);
break;
  }
  h/=2.0;
}
}
double f(double x)
{
return 1/(x*x);
}

被积函数为\frac{1}{x^2}
令a=1，b=5，输出约1.6，而精确解为0.8，错在哪里？

欢迎光临数学之家 (http://www.2math.cn/)