>>>>> "Raymond" == Raymond Toy <toy@rtp.ericsson.se> writes:
>>>>> "Raymond" == Raymond Toy <toy@rtp.ericsson.se> writes:
Raymond> A bug in FPPLUS perhaps? However, the code for float.lisp hasn't
Raymond> changed between 5.5 and 5.9rc1, and this isn't a problem with 5.5.
Raymond> In any case, here is a test case that demonstrates the issue.
Raymond> ans : -1.05b-1;
Raymond> inc : -2.5186b-221;
Raymond> sum : ans+ inc;
Raymond> returns -1.05b-1.
Raymond> But if we look at ans, inc, and the sum, we see:
Raymond> :lisp $ans
Raymond> ((BIGFLOAT SIMP 56) -60528378991859469 -3)
Raymond> :lisp $|inc|
Raymond> ((BIGFLOAT SIMP 56) -41001125856082033 -732)
Raymond> :lisp $|sum|
Raymond> ((BIGFLOAT SIMP 56) -60528378991859470 -3)
Raymond> So, surprisingly, the sum is now less than ans.
I think I found the problem. In rev 1.2, biglsh was replaced by ash
in fpshift. However, ash isn't the same as biglsh for negative
numbers, as the comments for biglsh indicate. It seems as if
(biglsh x n)
is
(- (ash (- x) n))
for negative x.
When I reinstate the original version of fpshift, the test case
changes, and log(.1b0) works again. Cool!
HOWEVER, the comment for the 1.2 change says:
BIGLSH is replaced by ASH since former seems to be buggy and
gives wrong result for 10.0b0-1.
And indeed 10.0b0 - 1 gives 9.5b-1. (!!!!)
After some experimentation, here is a possible fix to make both of
these work again:
(defun fpshift (l n)
(cond ((null *decfp)
(cond ((and (minusp n) (minus l))
;; Left shift of negative number requires some
;; care? If the shift of the positive number would
;; have returned 0, return 0. Otherwise, return the
;; negative of the shifted number. (That is,
;; (truncate l (expt 2 n))).
(let ((s (ash (- l) n)))
(if (zerop s)
0
(- s))))
(t
(ash l n))))
((greaterp n 0.) (times l (expt 10. n)))
((lessp n 0.) (quotient l (expt 10. (minus n))))
(t l)))
Or perhaps the fix is to fpplus so that if we right shift a negative
number so that all that's left is the sign, we replace it with zero
(which is what the above does, I think.)
Someone more familiar with how maxima's bigfloat works may want to
look at this or look at how fpplus works.
Ray