Hi, all.
I wrote this function and want someone to tell me if it seems ok.
If it doesn't do what it claims, I'd like to know whether what I'm
trying to do is possible in standard common lisp.
(defun single-float-above (f)
"Return the next representable single-float value greater than F."
(check-type f single-float)
(cond
((zerop f) least-positive-single-float)
((= most-positive-single-float f)
(error "~f is most-positive-single-float." f))
(t (multiple-value-bind (significand exponent integer-sign)
(integer-decode-float f)
(let ((result (float (* (1+ significand)
(expt (float-radix f) exponent))
f)))
(if (minusp integer-sign)
(- result)
result))))))
Thanks in advance.
Yuji.
From: Christophe Rhodes
Subject: Re: A function returning the next representable floating point value of an argument
Date:
Message-ID: <sqoeb1ggko.fsf@cam.ac.uk>
"yuji" <········@nifty.ne.jp> writes:
> I wrote this function and want someone to tell me if it seems ok.
> If it doesn't do what it claims, I'd like to know whether what I'm
> trying to do is possible in standard common lisp.
>
> (defun single-float-above (f)
> "Return the next representable single-float value greater than F."
> (check-type f single-float)
> (cond
> ((zerop f) least-positive-single-float)
This clause loses if the implementation supports signed zeros. I'd
write it as
((eql f 0.0) least-positive-single-float)
((eql f -0.0) 0.0)
and if the implementation treats -0.0 the same as 0.0, then the second
branch is never used.
> ((= most-positive-single-float f)
> (error "~f is most-positive-single-float." f))
> (t (multiple-value-bind (significand exponent integer-sign)
> (integer-decode-float f)
> (let ((result (float (* (1+ significand)
> (expt (float-radix f) exponent))
> f)))
Hmm. Without thinking /too/ hard about this, I'd venture that this
isn't going to work in an implementation which represents its floating
point values as IEEE single-floats, for F at the boundary between one
exponent and the next, in round-to-even rounding mode. In this
clause, I'd generate the float as you do, test it for equality against
the input float, and if it is equal, try again with (+ significand 2).
> (if (minusp integer-sign)
> (- result)
This would seem to be generating the next float further away from
zero, not the next greatest float, for negative f. (The logic in my
above paragraph needs to be more complex if you are really trying to
get the next greatest float.)
> result))))))
>
> Thanks in advance.
I hope that helps.
(As for whether what you want to do is possible in portable Common
Lisp: I doubt you can write something that is robust against all
possible implementations of floating point ever; on the other hand, I
think you probably /can/ write something which works on more than just
IEEE754 implementations)
Christophe
Thanks for the help.
Hmm, floating point stuff seems more delicate than I hoped.
I only know basic concepts of floating point arithmetic mostly from
scanning through the famous "What Every Computer Scientist Should Know
About Floating-Point Arithmetic" paper and barely heard of IEEE754.
>From your last comment on the portable implementation of this function,
I conclude that
it's beyond my current capacity.
Thanks again.
Yuji.
From: Bruno Haible
Subject: Re: A function returning the next representable floating point value of an argument
Date:
Message-ID: <d6vvu0$cak$1@laposte.ilog.fr>
Christophe Rhodes wrote:
>> (t (multiple-value-bind (significand exponent integer-sign)
>> (integer-decode-float f)
>> (let ((result (float (* (1+ significand)
>> (expt (float-radix f) exponent))
>> f)))
>
> Hmm. Without thinking /too/ hard about this, I'd venture that this
> isn't going to work in an implementation which represents its floating
> point values as IEEE single-floats, for F at the boundary between one
> exponent and the next, in round-to-even rounding mode.
> ...
> This would seem to be generating the next float further away from
> zero, not the next greatest float, for negative f.
You are right. But this one should work, uniformly producing the next
greater single-float, considering -0.0 and +0.0 as equal and ignoring
the problem of "denormalized floats":
(defun single-float-above (f)
"Return the next representable single-float value greater than F."
(check-type f single-float)
(cond
((zerop f) least-positive-single-float)
((= most-positive-single-float f)
(error "~f is most-positive-single-float." f))
(t (multiple-value-bind (significand exponent integer-sign)
(integer-decode-float f)
(if (minusp integer-sign)
(progn
(when (zerop (logand significand (- significand 1)))
(setq significand (* (float-radix f) significand))
(setq exponent (- exponent 1)))
(float (* (1+ (- significand)) (expt (float-radix f) exponent)) f))
(float (* (1+ significand) (expt (float-radix f) exponent)) f))))))
The 'logand' part probably needs to be changed if (float-radix f) is not 2.
However, I don't know how the significand (mantissa) is normalized in
implementations where (float-radix f) is 4 or 16 or something like this.
P.S.: What do you need this for? You can implement TYPEP and SUBTYPEP
without needing this kind of operation on floating-point numbers.
Bruno
Thank you for the comment.
> P.S.: What do you need this for? You can implement TYPEP and SUBTYPEP
> without needing this kind of operation on floating-point numbers.
Bruno san, are you a psychic or something? Yes, I'm implementing
SUBTYPEP now.
When reducing a number type spec into a simpler form, it helps if you
know the next number of a given one.
the type spec (or (single-float 0.0 0.2503) (single-float 0.2506 1.0))
could be transformed into (single-float 0.0 1.0), if you know there's
no single-float value on a given machine between 0.2503 and 0.2506,
then the two intervals (0.0 0.2503) and (0.2506 1.0) are continuous.
When it comes to integers, this kind of translation is a piece of case.
(or (integer 0 1) (integer 2 3)) => (integer 0 3).
Yuji