Is there a Bug in Trapzd C Code or not?


Adriana
03-14-2006, 08:38 AM
Hello folks!!
I am very new at this site, and not much expert about c programming (as well as forums), well, I'm just like a child!!
I have a doubt :confused: , when I look at trapzd code :

float trapzd(float (*func)(float), float a, float b, int n) {
float x,tnm,sum,del;
static float s; :eek:
int it,j;
if (n == 1) {
return (s=0.5*(b-a)*(FUNC(a)+FUNC(b)));
} else {
for (it=1,j=1;j<n-1;j++)
it <<= 1;
tnm=it;
del=(b-a)/tnm;
x=a+0.5*del;
for (sum=0.0,j=1;j<=it;j++,x+=del)
sum += FUNC(x);
s=0.5*(s+(b-a)*sum/tnm); :confused: :eek:
return s;
}
}

Why the variable 's' was'nt initialized?? It's updated without having been initialized? Is that right? I suppose I'm quite stupid, but I can't really understand...

hankel
03-22-2006, 07:39 AM
There is no bug in Trapzd, if you use it in an integration routine (such as qsimp(), for instance), you first call trapzd with n=1,
this sets a first value of s, subsequent calls just add corrections
to the previous result. If you use it directly with n=2 you will get a wrong result because you only apply corrections to a result whose first estimate hasn't been computed correctly