Would someone explain please what's going on here?
cl-user(1): (/ 209 5.5)
38.0
cl-user(2): (/ (* 100 2.09) 5.5)
37.999996
cl-user(3): (* 100 (/ 2.09 5.5))
38.0
cl-user(4):
I've tried ACL, CMUCL the result was identical. Are there any clues
about that? Of course I'd prefer 38.0 instead of 37.999996
--
Vladimir Zolotykh
Vladimir Zolotykh <······@eurocom.od.ua> writes:
> I've tried ACL, CMUCL the result was identical. Are there any clues
> about that? Of course I'd prefer 38.0 instead of 37.999996
http://docs.sun.com/source/806-3568/ncg_goldberg.html
http://en.wikipedia.org/wiki/Floating_point
--
;; Matthew Danish -- user: mrd domain: cmu.edu
;; OpenPGP public key: C24B6010 on keyring.debian.org
Matthew Danish <··········@cmu.edu> writes:
> Vladimir Zolotykh <······@eurocom.od.ua> writes:
>> I've tried ACL, CMUCL the result was identical. Are there any clues
>> about that? Of course I'd prefer 38.0 instead of 37.999996
>
> http://docs.sun.com/source/806-3568/ncg_goldberg.html
> http://en.wikipedia.org/wiki/Floating_point
Beyond these links (which you should read -- the first in particular
should be required knowledge before doing any numerical computation),
it might also help to consider the return values of (rational 2.09)
and (rationalize 2.09).
Christophe
Vladimir Zolotykh <······@eurocom.od.ua> writes:
> Would someone explain please what's going on here?
>
> cl-user(1): (/ 209 5.5)
> 38.0
> cl-user(2): (/ (* 100 2.09) 5.5)
> 37.999996
> cl-user(3): (* 100 (/ 2.09 5.5))
> 38.0
> cl-user(4):
>
> I've tried ACL, CMUCL the result was identical. Are there any clues
> about that? Of course I'd prefer 38.0 instead of 37.999996
Well, try this:
CL-USER> (* 100 2.09)
208.99998
Also of interest might be:
CL-USER> (decode-float (float 100 0.0))
0.78125
7
1.0
CL-USER> (decode-float 2.09)
0.5225
2
1.0
CL-USER> (* (* 0.5225 0.78125) (expt 2 (+ 7 2)))
208.99998
CL-USER> (decode-float (float 209 0.0))
0.81640625
8
1.0
CL-USER> (* 0.81640625 (expt 2 8))
209.0
Beyond that, I don't understand floating point math well enough to
explain it except in the handwavy way to point out that whenever you
do it you can lose little bits of accuracy as when you multiplied 100
by 2.09 getting 208.9998 rather than 209.0.
-Peter
--
Peter Seibel ·····@javamonkey.com
Lisp is the red pill. -- John Fraser, comp.lang.lisp
In article <··············@javamonkey.com>,
Peter Seibel <·····@javamonkey.com> wrote:
> Beyond that, I don't understand floating point math well enough to
> explain it except in the handwavy way to point out that whenever you
> do it you can lose little bits of accuracy as when you multiplied 100
> by 2.09 getting 208.9998 rather than 209.0.
The basic problem is easy to understand: 0.1 (a.k.a. 1/10) is a
repeating decimal when represented in base 2.
It is possible (and actually makes an educational and very doable
exercise) to write floating point math routines that represent numbers
in base 10 which do not suffer from this problem. (Of course, there is
an analogous problem for a base 10 representation when dealing with
fractions that are infinitely repeating decimals in base 10, like 1/3,
1/7, etc.)
rg
Ron Garret <·········@flownet.com> writes:
> The basic problem is easy to understand: 0.1 (a.k.a. 1/10) is a
> repeating decimal when represented in base 2.
>
> It is possible (and actually makes an educational and very doable
> exercise) to write floating point math routines that represent numbers
> in base 10 which do not suffer from this problem. (Of course, there is
> an analogous problem for a base 10 representation when dealing with
> fractions that are infinitely repeating decimals in base 10, like 1/3,
> 1/7, etc.)
I seem to recall that the 80387 math coprocessor had BCD instructions
and represented BCD numbers in 80 bits. I also seem to recall the
Borland C++ 2.0 compiler supporting BCD arithmetic although I don't
recall the declerations you needed to use to make BCD numbers.
Babbage's difference machine used decimal numbers, IIRC.
Would we count in octal if we had no pinkies?
--
An ideal world is left as an excercise to the reader.
--- Paul Graham, On Lisp 8.1
David Steuber <·····@david-steuber.com> writes:
> Ron Garret <·········@flownet.com> writes:
>
> > The basic problem is easy to understand: 0.1 (a.k.a. 1/10) is a
> > repeating decimal when represented in base 2.
> >
> > It is possible (and actually makes an educational and very doable
> > exercise) to write floating point math routines that represent numbers
> > in base 10 which do not suffer from this problem. (Of course, there is
> > an analogous problem for a base 10 representation when dealing with
> > fractions that are infinitely repeating decimals in base 10, like 1/3,
> > 1/7, etc.)
>
> I seem to recall that the 80387 math coprocessor had BCD instructions
> and represented BCD numbers in 80 bits. I also seem to recall the
> Borland C++ 2.0 compiler supporting BCD arithmetic although I don't
> recall the declerations you needed to use to make BCD numbers.
>
> Babbage's difference machine used decimal numbers, IIRC.
And the point is? You still could not divide 10 camels between 3 brothers.
--
__Pascal Bourguignon__ http://www.informatimago.com/
Our enemies are innovative and resourceful, and so are we. They never
stop thinking about new ways to harm our country and our people, and
neither do we.
David Steuber <·····@david-steuber.com> writes:
> I seem to recall that the 80387 math coprocessor had BCD instructions
> and represented BCD numbers in 80 bits. I also seem to recall the
The x87 has instructions to load 80-bit BCD integers. These are
converted to 80-bit/64-bit/32-bit floating point numbers on load
(depends on the setting of the FPU). The 8 FP registers have a width
of 80-bits internally, though.
Regards,
--
____________________________
Julian Stecklina / _________________________/
________________/ /
\_________________/ LISP - truly beautiful
At the risk of applying a band-aid where you need a tourniquet, you
should realize that you are not only using inexact numbers (floating
point), you are almost certainly using single-precision floating point
numbers.
Try 2.09d0, 5.5d0 instead of 2.09 and see if you are any happier.
For further information, browse the Hyperspec on
*read-default-float-format*
http://www.lispworks.com/reference/HyperSpec/Body/v_rd_def.htm#STread-default-float-formatST
double-float
http://www.lispworks.com/reference/HyperSpec/Body/t_short_.htm#double-float
exponent marker
http://www.lispworks.com/reference/HyperSpec/Body/26_glo_e.htm#exponent_marker
and related topics.
Vladimir Zolotykh wrote:
> Would someone explain please what's going on here?
>
> cl-user(1): (/ 209 5.5)
> 38.0
> cl-user(2): (/ (* 100 2.09) 5.5)
> 37.999996
> cl-user(3): (* 100 (/ 2.09 5.5))
> 38.0
> cl-user(4):
>
> I've tried ACL, CMUCL the result was identical. Are there any clues
> about that? Of course I'd prefer 38.0 instead of 37.999996
>
>
The reason that this is surprising is that the problem is hidden. If you
were to use an approximate value in a calculation you wouldn't expect an
'exact' answer:
3 X 1/3 = 1 vs. 3 X 0.3333 = 0.9999 != 1
The computer uses approximations for most floating-point numbers, but it
is not as obvious as the 0.3333 above.
It's true that any rational number can be expressed as an infinitely
repeating decimal:
1/2 = 0.50000... = 0.49999...
1/9 = 0.1111...
However, to manipulate such a representation we have to truncate it at
some point. We can't work with an infinite string of digits. The
computer has to do the same thing in binary. For instance,
2.09 => 10.000101110000101000111101011100001010001111...
But the computer is storing it as a single float this way:
10.00010111000010100011110000000
If you use RATIONAL you see
(rational 2.09) => 8766095/4194304
and
(- 209/100 (rational 2.09)) => 9/104857600
So your 2.09 is an approximation.
Next, when you multiply by 100 you get:
11010000.1111111111111111000000
But 209 in binary is 11010001. The approximation for 2.09 led to an
inexact product.
The number 5.5 can be represented 'nicely' in binary since the
fractional part is a negative integral power of 2:
5.5 = 101.100000000000
So in your first calculation, dividing the exact integer 209 by the
exact 5.5 gives you exactly 38.0:
rational (/ 209 5.5)) => 38
However, the inexact (* 100 2.09) divided by 5.5 produces:
100101.11111111111111111100000
On the other hand, in your third example the initial division is inexact:
(- 209/550 (rational (/ 2.09 5.5))) => 1/209715200
But because of the way the computer's representation is truncated, when
you multiply by 100 you do get exactly 38 again. This is simply an
accident. Don't count on it.
(/ 2.09 5.5) =>
0.01100001010001111010111000000
209/550 =>
0.011000010100011110101110000101000111101...
(rational (* 100 (/ 2.09 5.5))) => 38
This issue is simply a fact of life. Read Goldberg's paper. It's hard
work, but it will answer your questions.
David Sletten
On Mon, 27 Sep 2004 18:55:07 GMT, David Sletten <·····@slytobias.com>
wrote:
> This issue is simply a fact of life. Read Goldberg's paper. It's hard
> work, but it will answer your questions.
I got the point. I should read Goldberg's paper however at first
glance it looked frightening.
May I trouble you a bit more? There is another difficulty:
The following C program invariably prints '=' (GCC 2.95.3).
#include <stdio.h>
#include <stdlib.h>
int main(int argc, char *argv[])
{
float a = 2.09 * 100;
printf("%c\n", (a == 209) ? '=' : '!');
}
After your explanation I'd suppose that it would print '!'? Does GCC
treat floats differently?
--
Vladimir Zolotykh
In C, without specifically using a 'f' or 'l' exponent marker,
floating-point constants are double precision. For *this particular*
choice of numbers, IEEE754 double precision avoids visible round-off
error, while single precision does not.
In Lisp, without an explicit exponent marker, the reader treats
floating-point constants according to the value of
*read-default-float-format*, which will have the initial value
'single-float.
> #include <stdio.h>
> #include <stdlib.h>
>
> int main(int argc, char *argv[])
> {
> float a = 2.09 * 100;
> printf("%c\n", (a == 209) ? '=' : '!');
> }
>
> After your explanation I'd suppose that it would print '!'? Does GCC
> treat floats differently?
When you write float a = 2.09 * 100, what actually happens is that
both 2.09 and 100 are converted to double precision (the internal FPU
registers are double precision), the multiplication is done in double
precision, and then the result is truncated to a float. This should
be true across compilers; it's not gcc specific. This should help
elucidate the difference:
#include <stdio.h>
int main(int argc, char **argv) {
int a = 100;
float b = 2.09;
float c = a * b;
float d = 100 * 2.09;
printf("%f %f\n", c, d);
}
As a side note, you're generally much better off using double
precision for everything unless you're very low on space. Using
double precision for everything is nearly always the right correct
choice.
Cheers,
rif
Many thanks gentlemen!
I've never been so interested and what is more important from this
discussion I've learnt more about floats than after reading some dry
books and rambling misc. manuals.
--
Vladimir Zolotykh
Vladimir Zolotykh <······@eurocom.od.ua> writes:
> May I trouble you a bit more? There is another difficulty:
>
> The following C program invariably prints '=' (GCC 2.95.3).
So does the equivalent (but shorter :) Common Lisp program:
(= (* 2.09 100) 209)
> #include <stdio.h>
> #include <stdlib.h>
>
> int main(int argc, char *argv[])
> {
> float a = 2.09 * 100;
> printf("%c\n", (a == 209) ? '=' : '!');
> }
>
> After your explanation I'd suppose that it would print '!'? Does GCC
> treat floats differently?
Some floats work out correctly in base 10, but that isn't a guarantee
that all floats will exhibit this property.
--
Thomas A. Russ, USC/Information Sciences Institute
Vladimir Zolotykh <······@eurocom.od.ua> writes:
> float a = 2.09 * 100;
> printf("%c\n", (a == 209) ? '=' : '!');
> After your explanation I'd suppose that it would print '!'? Does GCC
> treat floats differently?
Use gcc -S -o f.s f.c to get the assembler source:
main:
pushl %ebp
movl %esp, %ebp
subl $8, %esp
andl $-16, %esp
movl $0, %eax
subl %eax, %esp
movl $0x43510000, -4(%ebp) ;; Oh oh! No multiplication!
subl $8, %esp
flds -4(%ebp)
flds .LC0
fxch %st(1)
fucompp
fnstsw %ax
andb $69, %ah
cmpb $64, %ah
je .L4
jmp .L2
.L4:
movl $61, -8(%ebp)
jmp .L3
.L2:
movl $33, -8(%ebp)
.L3:
pushl -8(%ebp)
pushl $.LC1
call printf
addl $16, %esp
leave
ret
Try this:
#include <stdio.h>
#include <stdlib.h>
int main(int argc, char *argv[])
{
float a = 2.09 * 100;
float b = 2.09;
printf("a: %c\n", (a == 209) ? '=' : '!');
printf("b: %c\n", (b*100 == 209) ? '=' : '!');
}
cd /tmp/
make f ; ./f
cc f.c -o f
a: =
b: !
Or that for some more fun:
#include <stdio.h>
#include <stdlib.h>
int main(int argc, char *argv[])
{
float f=1147433643.0;
int i=1147433600;
printf("%f %c %d \n",f,(f==i)?'=':'!',i);
}
cd /tmp/
make f ; ./f
cc f.c -o f
1147433600.000000 = 1147433600
So, you still want to program in C?
--
__Pascal Bourguignon__ http://www.informatimago.com/
Our enemies are innovative and resourceful, and so are we. They never
stop thinking about new ways to harm our country and our people, and
neither do we.