#1264 Floating point imprecission

closed
nobody
5
2007-09-22
2007-09-03
Abel Pascual
No

I don't understand why Maxima rounds down the number...
I'm currently running Ubuntu 7.04.

Maxima 5.13.0 http://maxima.sourceforge.net
Using Lisp GNU Common Lisp (GCL) GCL 2.6.7 (aka GCL)
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
This is a development version of Maxima. The function bug_report()
provides bug reporting information.
(%i1) 4.6*10^8;
(%o1) 4.5999999999999994E+8

Discussion

  • Abel Pascual
    Abel Pascual
    2007-09-03

    Logged In: YES
    user_id=1881850
    Originator: YES

    Greetings

     
  • Harald Geyer
    Harald Geyer
    2007-09-03

    Logged In: YES
    user_id=929336
    Originator: NO

    Hi Abel!

    Thanks for you interest in maxima. The behaviour, that you see is a basic property of fixed precision floating point numbers (maxima acutally is using maschine doubles IIRC): The number 4.6e8 is not
    representable exactly in base 2 arithmetic, thus maxima rounds the
    number to the next number that can be represented as double float.

    If you need higher precision, you can use bfloats: 4.6b8

    If you want exact numbers, you can use rational arithmetic. In this case enter the number as: 46*10^7 (that is: without decimal dot)

    If you still have questions, you might ask on the Mailinglist. I'm closing your bug report now, because this is nothing we can fix in maxima.

     
  • Harald Geyer
    Harald Geyer
    2007-09-03

    • status: open --> closed
     
  • Abel Pascual
    Abel Pascual
    2007-09-03

    Logged In: YES
    user_id=1881850
    Originator: YES

    Ok, sorry for my mistake, and thanks for your explanation...

     
  • Raymond Toy
    Raymond Toy
    2007-09-04

    • status: closed --> open
     
  • Raymond Toy
    Raymond Toy
    2007-09-04

    Logged In: YES
    user_id=28849
    Originator: NO

    Actually 4.6e8 is an integer: 460000000. This should be exactly representable as a double-float. I think it is an error that maxima prints something else out.

    Reopening this bug.

     
  • Abel Pascual
    Abel Pascual
    2007-09-04

    Logged In: YES
    user_id=1881850
    Originator: YES

    I don't know if there is any relation, but in GCL I'm having the same issue:

    >(* 4.6 (expt 10 8))

    4.5999999999999994E8

    But otherwise in CLISP this is not happening:

    [1]> (* 4.6 (expt 10 8))
    4.6E8

     
  • Raymond Toy
    Raymond Toy
    2007-09-04

    Logged In: YES
    user_id=28849
    Originator: NO

    The issues you see with GCL and CLISP are somewhat related, but outside the scope of this bug.

    If you want to discuss it further, I'd be happy to. Just not here. Please mail me directly or, even better, use the maxima mailing list.

     
  • Logged In: YES
    user_id=1179910
    Originator: NO

    I don't think 4.6*10^8 is parsed as one number. First 4.6 is parsed, then 10^8, and then they are multiplied. To enter 4.6e8 you should enter it as 4.6e8.

    OTOH, 4.6e8 is not printed correctly on clisp and gcl. The reason is that floats are printed with

    (format nil "~ve" 21 4.6e8)

    I think 21 is arbitrarily chosen - if it was 19, all lisps print it correctly, but then maybe some other example would be printed incorrectly. Why don't we just use (format nil "~ve" 16 4.6e8).

    Andrej

     
  • Raymond Toy
    Raymond Toy
    2007-09-05

    Logged In: YES
    user_id=28849
    Originator: NO

    Good question. I wasn't sure, so I checked the code. Look at make-number and readlist in nparse.lisp. If the number is not a bigfloat, readlist is called, which basically calls Lisp read-from-string to read "4.6e8". So the result depends on the underlying Lisp. All Lisps should say (floor 4.6d8) is 46000000 and 0.0d0. Clisp and cmucl do. gcl does not. We could implement our own reader.

    However, for a bigfloat, maxima sees 4.6, converts it to 46/10 and then multiplies by 10^8 to get a rational number. This is then converted to a bigfloat.

    I think 21 is because we want upto 16 digits printed, and the 5 is to make sure there's space to print the exponent marker, a sign, and the exponent itself (3 digits). If 16 is used, the number of significant digits would be too small.

    I think it is a bug in gcl and clisp that it prints something other than 4.6d8.

     
  • Robert Dodier
    Robert Dodier
    2007-09-07

    Logged In: YES
    user_id=501686
    Originator: NO

    (1) I've reported GCL's failure to parse 4.6e8 correctly to the GCL bug tracker as # 20992.
    https://savannah.gnu.org/bugs/index.php?20992
    So as for this 4.6e8 in Maxima + GCL, we can close it because the problem is not in Maxima.

    (2) About Clisp printing 46e7 as 4.5999999999999996E+8, it appears that Clisp is interpreting "~ve" to indicate a coercion to single-float before formatting it. Changing the "~ve" to "~vd" in EXPLODEN in src/commac.lisp appears to tell Clisp to format 46e7 as 4.6E8. SBCL and GCL are happy with "~vd".

    Leaving this report open to gather comments on item (2). If nobody is opposed, I'll modify EXPLODEN as described above.

     
  • Raymond Toy
    Raymond Toy
    2007-09-07

    Logged In: YES
    user_id=28849
    Originator: NO

    We should probably file a bug report for clisp for item 2.

    Changing ~e to ~d is ok, but ~d is really ~a, so we should just use ~a.

    There is a difference, however. ~e pads numbers on the left and ~a pads on the right.t. DDoDont' know if that matters or not.

    In addition, cmucl has a known bug that ~e prints numbers differently from ~a (princ). ~a is more accurate.

     
  • Raymond Toy
    Raymond Toy
    2007-09-07

    Logged In: YES
    user_id=28849
    Originator: NO

    I'm mistaken. Using ~a won't work (which implies ~d won't work either). Consider if fpprintprec is 2. Then
    (format t "~7a" 1.2345678d9) -> 1.2345678d9
    But
    (format t "~7e" 1.2345678d9) -> 1.23d8

    We need to do something else.

     
  • Robert Dodier
    Robert Dodier
    2007-09-07

    Logged In: YES
    user_id=501686
    Originator: NO

    Ray, I think you are correct on both counts: changing e to a or d shouldn't help, the format indicator should remain e; and the Clisp output in this case (format nil "~ve" 21 46d7) indicates a bug in Clisp.

    You or I can report the Clisp bug, whoever gets to it first, and then we can close this report, I think.

     
  • Robert Dodier
    Robert Dodier
    2007-09-08

    Logged In: YES
    user_id=501686
    Originator: NO

    Reported the bug (format nil "~ve" 21 46d7) => 4.5999999999999996d+8 as SF bug report # 1790496 for Clisp.
    http://sourceforge.net/tracker/index.php?func=detail&aid=1790496&group_id=1355&atid=101355

    I think we can close this Maxima bug report now since both of problems observed are problems in Lisp implementations. Marking this report as "won't fix" and "pending" (to be closed automatically in 2 weeks) accordingly.

     
  • Robert Dodier
    Robert Dodier
    2007-09-08

    • status: open --> pending
    • labels: 887077 --> Problem not in Maxima
     
  • Raymond Toy
    Raymond Toy
    2007-09-08

    Logged In: YES
    user_id=28849
    Originator: NO

    There is one thing we could do: If fpprintprec is 0, we want to print everything, so in that case ,we might use ~a. If fpprintprec is something else, we use ~e.

    This will probablly fix the immediate issue. I assume most people leave fpprintprec defaulted.

     
  • Robert Dodier
    Robert Dodier
    2007-09-08

    Logged In: YES
    user_id=501686
    Originator: NO

    Ray, it's OK by me if you want to modify EXPLODEN to treat fpprintprec = 0 as a special case.

     
  • Raymond Toy
    Raymond Toy
    2007-09-10

    Logged In: YES
    user_id=28849
    Originator: NO

    Aargh! I completely misread the problem. This whole discussion has been about 4.6e8. But the original problem is 4.6*10^8, which is completely different. 4.6 isn't exactly representable as a float The printed result is, as best as I can tell, the correct answer.

    However, the issue with clisp printing 4.6e8 incorrectly is still valid. Same with the bug in gcl.

    Perhaps we should just leave exploden as is.

     
  • Logged In: YES
    user_id=1312539
    Originator: NO

    This Tracker item was closed automatically by the system. It was
    previously set to a Pending status, and the original submitter
    did not respond within 14 days (the time period specified by
    the administrator of this Tracker).

     
    • status: pending --> closed