On 6/7/07, Barton Willis <willisb at unk.edu> wrote:
> (%i65) ff:exp(-x^2/50^2)*sin(x)^2*cos(2*x)$
> (%i66) gg : integrate(ff,x,-4*%pi,4*%pi)$
> (%i67) float(gg);
> (%o67) 0.0
> (%i68) float(expand(gg));
> (%o68) -6.153361529128077
>
> Look at both gg and expand(gg). You should see why (%o67) = 0.0, but
> (%o68) isn't 0.0.
%o67 = 0.0 looks to me like a bug in GCL's floating point math.
There is one term exp(10000) and another one exp(-10000).
The first yields a floating point inf and the second, 0.0.
The product should be NaN (not-a-number).
SBCL could probably be coaxed into yielding NaN for %o67.
(Clisp cannot, to the best of my knowledge.)
About this stuff,
>(%i164) integrate(ff,x,-4*%pi,4*%pi);
>(%o164) -(%e^(-36)*sqrt(%pi)*(3*erf((4*%pi+18*%i)/3)-6*%e^27*erf((4*%pi
>+9*%i)/3)-6*%e^27*erf((4*%pi-9*%i)/3)+3*erf((4*%pi-18*%i)/3)+6*%
>e^36*erf((4*%pi)/3)))/8
>
>(%i165) float(%);
>(%o165) -5.1390589659105944*10^-17*(3.0*erf(0.33333333333333*(18.0*%i
>+12.56637061435917))-3.192289443610792*10^+12*
>erf(0.33333333333333*(9.0*%i+12.56637061435917))-3.192289443610792*10^
>+12*erf(0.33333333333333*(12.56637061435917-9.0*%i))+3.0*
>erf(0.33333333333333*(12.56637061435917-18.0*%i))+2.5867389201337716*10^
>+16)
I wonder if Maxima could output something nicer here.
(1) Could erf applied to a complex float (i.e. float + %i * float)
yield a float or complex float?
(2) I tried abs(%o165) and realpart(%o165) and imagpart(%o165)
to try to get real floats from the %o165 but abs and friends didn't
seem to know how to do something useful with it. I wonder what
we could do to change that. Maybe (2) is moot if (1) is implemented.
FWIW
Robert Dodier