Robert Dodier wrote:
> On Wed, Mar 18, 2009 at 2:47 PM, Edwin Woollett <woollett at charter.net> wrote:
>
>
>> I suppose a workaround is to define my own exponential function with
>> graceful underflow behavior, but isn't this something which should be
>> part of the bfloat code??
>>
>
> It appears the problem is due to failing to detect the limits
> of the bigfloats that can be handled by the machine
> (so it hangs or falls over instead of printing a Maxima
> error message). I guess that's a bug.
>
> 1b<foo> attempts to call (EXPT 10 <foo>) so when <foo>
> is a large integer (positive or negative), it is trying to compute
> a very large integer or 1/(very large integer). I'm guessing that
> when Maxima appears to hang, it is actually running EXPT and
> just taking a very long time.
But it doesn't have to do it this way. It could compute 10^foo using
bigfloat operations instead. Potentially slower, and certainly more
roundoff.
In any case, I found the problem. After tracing fpplus, I can see the
issue is adding two numbers with one having a very large negative
exponent. The problem is in haipart (clmacs.lisp). It's called with
that very negative exponent and wants to return the least significant n
bits, and does this by masking. So gcl (and other lisps!) will compute
a huge number to mask out a relatively small number.
Solution is to check the length of the mantissa and if it's smaller than
the desired number of bits, just return the number.
This allows gcl to compute g(40,40) without problem. g(1000,1000) works too.
I'll check in the change soon.
Ray