When toying around with these functions:
(defun almost-equal (a b &optional (delta .001))
(<= (abs (- a b)) delta))
(defun compose-or (&rest functions)
#'(lambda(&rest args)
(if (null functions) nil
(or (apply (first functions) args)
(apply (apply #'compose-or (rest functions)) args)))))
(funcall (compose-or #'< #'almost-equal) 1 .999)
I would have expected the last expression to evaluate to T.
Alas:
? (funcall (compose-or #'< #'almost-equal) 1 .999)
NIL
Then I tried to find out what the problem could be. I evaluated:
? (- 1 .999)
0.0010000000000000009
?
So 1 minus .999 is 0.0010000000000000009 according to my Lisp
implementation (MCL 5.0 on Mac OSX 10.3.8, G4 400 Mhz PPC , 1.19 GB
SDRAM).
Can someone explain this?
Btw: also comments on the code of compose-or are welcome, if there is a
"better" way of doing this kind of thing in CL.
Greetings,
Michiel
"Borkdude" <··············@gmail.com> writes:
> ? (- 1 .999)
> 0.0010000000000000009
That's how floating-point arithmetic works. Maybe what you meant was:
(- 1 999/1000)
? You'll find that many Lispers lack patience with questions from
people who don't understand floating-point math because we *do* have
proper rationals. If you want rational arithmetic, use rationals
(sort of like how if you want arrays, use arrays, not lists).
"Borkdude" <··············@gmail.com> writes:
> So 1 minus .999 is 0.0010000000000000009 according to my Lisp
> implementation (MCL 5.0 on Mac OSX 10.3.8, G4 400 Mhz PPC , 1.19 GB
> SDRAM).
>
> Can someone explain this?
999/1000 cannot be represented exactly in floating point, so your Lisp
uses the closest floating point approximation:
8998192055486251/9007188254740992
which is just a tad smaller than 999/1000
When you subtract this from 1,
9007188254740992 8998192055486251 9007188254741
---------------- - ---------------- = ----------------
9007188254740992 9007188254740992 9007188254740992
you therefore get a number that is just a tad larger than 1/1000
When your lisp prints this number, it prints the shortest decimal
approximation that can unambiguously be read back in as the same
floating point number. That happens to be
0.0010000000000000009
I now have the following version of almost-equal:
1.
(defun almost-equal (a b &optional (delta .001))
(<= (abs (- a b)) delta))
2.
(defun almost-equal (a b &optional (delta .001))
(<= (abs (- a b)) (rationalize delta)))
3.
(defun almost-equal (a b &optional (delta .001))
(<= (abs (- (rationalize a) (rationalize b))) (rationalize delta)))
4.
(defun almost-equal (a b &optional (delta .001))
(<= (abs (- a b)) (+ delta (* 2 SINGLE-FLOAT-EPSILON))))
5.
(defun almost-equal (a b &optional (delta .001))
(<= (abs (- (rationalize a) (rationalize b))) (+ (rationalize delta)
(* 2 SINGLE-FLOAT-EPSILON))))
Version 1 and 2 do not work as I want when applied to the numbers 1 and
.999. Versions 3 and 4 and of course 5 do. In version 3 every number is
rationalized. In version 4 the epsilon is taken into account. In
version 5 both methods are combined.
A few questions:
1) Which version of {3,4,5} would you prefer to use now? And more
important: why? All seem to work, but version 5 might be a little
exaggerated.
2) Are there situations in which SINGLE-FLOAT-EPSILON is not enough and
should be replaced by e.g., DOUBLE-FLOAT-EPSILON? Is there a way to do
this easily without having to check the type of the numbers you deal
with and writing it out for every type?
I see what you're trying to do in 4-5 with two epsilons, but it's
hopeless if you want to use it with substantial computations. If x and
y are epsilon apart, then x+x+x and y+y+y are roughly 3*epsilon apart.
Multiplication is worse: because (x + epsilon)^3 = x^3 + 3*x^2*epsilon
+ lower-order-terms, the difference between
(100+epsilon)^3 and 100^3 is around 3*100^2 = 30,000 epsilons! Imagine
what happens after 1000 lines of code with lots of arithmetic and
unpredictable if-branching. The subject of numerical analysis is all
about these difficulties.
Using rationalize runs into the same problems, no more and no less. It
tries to replace x with a rational that's within epsilon of x. But
then the rational version of x+x+x is still off by 3 epsilons, etc.
The basic answer is, use 1, and choose delta after experimentation with
*your* specific code and *your* specific problem.
A good introduction to these ideas is Numerical Recipes in C [or
Fortran, etc.]
These books are wonderful: smart, snide, funny, and supremely
practical.
"Borkdude" <··············@gmail.com> writes:
> I now have the following version of almost-equal:
>
> 1.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- a b)) delta))
>
> 2.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- a b)) (rationalize delta)))
>
> 3.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- (rationalize a) (rationalize b))) (rationalize delta)))
>
> 4.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- a b)) (+ delta (* 2 SINGLE-FLOAT-EPSILON))))
>
> 5.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- (rationalize a) (rationalize b))) (+ (rationalize delta)
> (* 2 SINGLE-FLOAT-EPSILON))))
>
> Version 1 and 2 do not work as I want when applied to the numbers 1 and
> .999.
You *can't* apply them to the numbers 1 and .999 because you can't
represent .999
> Versions 3 and 4 and of course 5 do. In version 3 every number is
> rationalized. In version 4 the epsilon is taken into account. In
> version 5 both methods are combined.
>
> A few questions:
>
> 1) Which version of {3,4,5} would you prefer to use now? And more
> important: why? All seem to work, but version 5 might be a little
> exaggerated.
What are you trying to accomplish?
There is no foolproof way to deal with rounding because there is no
way of knowing how much rounding has occurred (if there were, it'd be
a simple matter to correct for it, right?)
"Borkdude" <··············@gmail.com> writes:
> 1) Which version of {3,4,5} would you prefer to use now? And more
> important: why? All seem to work, but version 5 might be a little
> exaggerated.
Actually, #6:
6. (defun almost-equal (a b &optional (delta .001))
(<= (abs (- a b)) (+ delta (* 2 SINGLE-FLOAT-EPSILON))))
> 2) Are there situations in which SINGLE-FLOAT-EPSILON is not enough and
> should be replaced by e.g., DOUBLE-FLOAT-EPSILON? Is there a way to do
> this easily without having to check the type of the numbers you deal
> with and writing it out for every type?
Well, there aren't THAT many floating point types. And single float
epsilon should always be at least as big as double float epsilon. But
if it were truly necessary, you could implement this using really ugly
cascading typecase statements:
7. (defun almost-equal (a b &optional (delta .001))
(typecase a
(SINGLE-FLOAT
(typecase b
(SINGLE-FLOAT
(<= (abs (- a b)) (+ delta (* 2 SINGLE-FLOAT-EPSILON))))
(DOUBLE-FLOAT
(<= (abs (- a b)) (+ delta (* 2 DOUBLE-FLOAT-EPSILON))))
...))
(DOUBLE-FLOAT
(typecase b
...))))
but for a truly general solution you might need to consider short and
long floats as well as single and double. I suppose an alterate
implementaiton might be able to use some sort of lookup table keyed by
taking TYPE-OF each of the operands.
The the fundamental problem of floating point imprecision will remain,
because there is the tendency of errors to grow with the number of
operations that are applied to the data. If a and b, for example, are
the products of some computations, rather than being just read by the
Lisp reader, then there isn't any guarantee that the error won't
overwhelm twice the epsilon for the particular floating point format
anyway.
If you really care about this, you will operate using rational numbers
exclusively. Convert all input using RATIONALIZE and then convert the
output using FLOAT.
--
Thomas A. Russ, USC/Information Sciences Institute
"Borkdude" <··············@gmail.com> writes:
> I now have the following version of almost-equal:
>
> 1.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- a b)) delta))
>
> 2.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- a b)) (rationalize delta)))
>
> 3.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- (rationalize a) (rationalize b))) (rationalize delta)))
>
> 4.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- a b)) (+ delta (* 2 SINGLE-FLOAT-EPSILON))))
>
> 5.
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- (rationalize a) (rationalize b))) (+ (rationalize delta)
> (* 2 SINGLE-FLOAT-EPSILON))))
>
> Version 1 and 2 do not work as I want when applied to the numbers 1 and
> .999. Versions 3 and 4 and of course 5 do. In version 3 every number is
> rationalized. In version 4 the epsilon is taken into account. In
> version 5 both methods are combined.
>
> A few questions:
>
> 1) Which version of {3,4,5} would you prefer to use now? And more
> important: why? All seem to work, but version 5 might be a little
> exaggerated.
> 2) Are there situations in which SINGLE-FLOAT-EPSILON is not enough and
> should be replaced by e.g., DOUBLE-FLOAT-EPSILON? Is there a way to do
> this easily without having to check the type of the numbers you deal
> with and writing it out for every type?
Usually, (< DOUBLE-FLOAT-EPSILON SINGLE-FLOAT-EPSILON)
Note that these epsilon values are valid near 1. Since you're
interested with absolute values close to 0, I deemed them close enough
to 1. More sophisticated mathematical analysis would be in order to
get a more precise test.
--
__Pascal Bourguignon__ http://www.informatimago.com/
Until real software engineering is developed, the next best practice
is to develop with a dynamic system that has extreme late binding in
all aspects. The first system to really do this in an important way
is Lisp. -- Alan Kay
Yes, that must explain it. How could I forget, that what I am
programming has still something to do with the physical architecture on
which it is run! :)
Borkdude schrieb:
> Yes, that must explain it. How could I forget, that what I am
> programming has still something to do with the physical architecture on
> which it is run! :)
It depends on your implementation.
One Lisp I tried (IIRC it was LispWorks on Windows) gave me:
CL-USER> (+ 0.000001 0.000001)
0.000002
While CLISP says:
CL-USER> (+ 0.000001 0.000001)
2.0E-6
Anyway, is there a good reason why you said 0.001 instead of using the
rational 1/1000?
If one does not really need ultimate speed, why should one /not/ use
rationals?
Andr�
--
Andr� Thieme schrieb:
> If one does not really need ultimate speed, why should one /not/ use
> rationals?
I remember I was so pissed about possibly losing information by using
the floating-point unit I wrote a little function float2rational. About
30 seconds when I was done with writing and testing it I found
(RATIONALIZE) ;-)
(funcall (compose-or #'< #'almost-equal) 1 (rationalize 0.999))
and
(defun almost-equal (a b &optional (delta 1/1000))
(<= (abs (- a b)) delta))
Andr�
--
"Borkdude" <··············@gmail.com> writes:
> So 1 minus .999 is 0.0010000000000000009 according to my Lisp
> implementation (MCL 5.0 on Mac OSX 10.3.8, G4 400 Mhz PPC , 1.19 GB
> SDRAM).
>
> Can someone explain this?
I think (*) that's within the single-float-epsilon of MCL, isn|t it?
(*) sorry, forgot to install MCL on the mac next to me...
--
(espen)
Borkdude schrieb:
> I would have expected the last expression to evaluate to T.
>
> ? (funcall (compose-or #'< #'almost-equal) 1 .999)
> NIL
Tried it with "Lisp in a Box" under windows:
CL-USER> (defun almost-equal (a b &optional (delta .001))
(<= (abs (- a b)) delta))
ALMOST-EQUAL
CL-USER> (defun compose-or (&rest functions)
#'(lambda(&rest args)
(if (null functions) nil
(or (apply (first functions) args)
(apply (apply #'compose-or (rest functions))
args)))))
COMPOSE-OR
CL-USER> (funcall (compose-or #'< #'almost-equal) 1 .999)
T
Does it help if you do (delta 1/1000) ?
Andr�
--
"Borkdude" <··············@gmail.com> writes:
> When toying around with these functions:
>
> (defun almost-equal (a b &optional (delta .001))
> (<= (abs (- a b)) delta))
You should take into account epsilon
http://www.lispworks.com/reference/HyperSpec/Body/v_short_.htm
In a first approximation:
(defun almost-equal (a b &optional (delta .001))
(<= (abs (- a b)) (+ delta (* 2 SINGLE-FLOAT-EPSILON))))
should give better results.
--
__Pascal Bourguignon__ http://www.informatimago.com/
In a World without Walls and Fences,
who needs Windows and Gates?
"Borkdude" <··············@gmail.com> writes:
> (defun compose-or (&rest functions)
> #'(lambda(&rest args)
> (if (null functions) nil
> (or (apply (first functions) args)
> (apply (apply #'compose-or (rest functions)) args)))))
Here are some variations that cons less when the returned
function is eventually called.
Tested on SBCL 0.8.17.20 with:
(declaim (optimize (speed 1) (space 2) (debug 0) (safety 1)))
(let ((x (compose-or #'consp #'arrayp #'symbolp)))
(time (loop repeat 100000 do (funcall x 1))))
Your version:
;; Evaluation took:
;; 6.623 seconds of real time
;; 5.78 seconds of user run time
;; 0.7 seconds of system run time
;; 0 page faults and
;; 328,374,280 bytes consed.
(defun compose-or (&rest functions)
#'(lambda (&rest args)
(some #'(lambda (function) (apply function args))
functions)))
;; Evaluation took:
;; 5.13 seconds of real time
;; 4.37 seconds of user run time
;; 0.44 seconds of system run time
;; 0 page faults and
;; 256,307,752 bytes consed.
(defun compose-or (&rest functions)
(reduce #'(lambda (left right)
#'(lambda (&rest args)
(or (apply left args)
(apply right args))))
functions :initial-value (constantly nil)))
;; Evaluation took:
;; 3.316 seconds of real time
;; 2.86 seconds of user run time
;; 0.28 seconds of system run time
;; 0 page faults and
;; 152,191,216 bytes consed.
I think the following version is what I would use by default.
(defun compose-or (&rest functions)
#'(lambda (&rest args)
(loop for function in functions
thereis (apply function args))))
;; Evaluation took:
;; 1.995 seconds of real time
;; 1.88 seconds of user run time
;; 0.09 seconds of system run time
;; 0 page faults and
;; 40,065,008 bytes consed.
(defun compose-or (&rest functions)
(compile nil `(lambda (&rest args)
(or ,@(loop for function in functions
collect `(apply ',function args))))))
;; Evaluation took:
;; 0.705 seconds of real time
;; 0.5 seconds of user run time
;; 0.01 seconds of system run time
;; 0 page faults and
;; 8,003,584 bytes consed.
> ? (- 1 .999)
> 0.0010000000000000009
> ?
That's a roundoff error due to the binary floating-point representation.
If you want exact results, use ratios:
* (- 1 999/1000)
1/1000