积分程序
辛普森法
#include <stdio.h>
#include <math.h>
#define F(X) (4.0/(1+X*X))
static float EPS = 1.0E-14;
static int COUNT=1;
static double a=0.0, b=1.0;
static double M[32],T[32],S[32],C[32],R[32],E[32];
double GETM(int K)
{ unsigned long j,n=1;
double x=0,y=0,step=0;
for(j=0;j<K;j++) n*=2;
step = (b-a)/n;
x = a+step/2;
for(j=0;j<n;j++){y+=F(x);x+=step;}
return (y*step);
}
Simpson(float EPS)
{ int k=0;
T[0]= (b-a)*(F(a)+F(b))/2.0;
S[0]=T[0];
E[0]=1.0E10;
for(k=0;k<20;k++)
{ M[k] = GETM(k);
T[k+1]=(M[k]+T[k])/2.0;
S[k+1]=(4.0*T[k+1]-T[k])/3.0;
COUNT++;
E[k+1]=fabs(S[k+1]-S[k]);
if(k<3)continue;
if(E[k+1]<EPS)break;
}
return;
}
ShawResult()
{ int k;
system("cls");
printf("\n K M[K] T[K] S[K]");
printf(" E[K]");
printf("\n-------------------------------------------");
printf("-----------------------------------");
for(k=0;k<COUNT;k++)
printf("\n%2d %20.15lf%20.15lf%20.15lf%12.4e",k,M[k],T[k],S[k],E[k]);
printf("\n--------------------------------------------");
printf("-----------------------------------");
getch();
system("cls");
return;
}
SaveResult()
{ int k;
FILE * fp;
fp=fopen("Simpson.htm","w");
if(!fp) return;
fprintf(fp,"<html><head></head>");
fprintf(fp,"<body bgcolor = #006699 text = #FFff00>");
fprintf(fp,"<pre><font size=\"6\">");
fprintf(fp,"\n K M[K] T[K] S[K]");
fprintf(fp," E[K]");
fprintf(fp,"\n----------------------------------------------");
fprintf(fp,"----------------------------------");
for(k=0;k<COUNT;k++)
fprintf(fp,"\n%2d %20.15f%20.15f%20.15f%12.4le",k,M[k],T[k],S[k],E[k]);
fprintf(fp,"\n----------------------------------------------");
fprintf(fp,"----------------------------------");
fprintf(fp,"</pre></font></body>");
fclose(fp);
return;
}
main()
{ Simpson(EPS);
ShawResult();
SaveResult();
}