<html><head></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; ">Brian,<div><br></div><div>I must reset my mental model of how conversions to and from float work. I wasn't too surprised by your examples until I came to (= 0.010000000000000009 0.01) which according to my previous understanding should have returned T. I.e., I was under the impression that 0.01 would have an internal representation identical to the other and therefore the comparison should return true.</div><div><br></div><div>Many thanks for the links to reading material.</div><div><br></div><div>Cheers</div><div>Sudhir</div><div><br><div><div>On Jun 4, 2012, at 7:44 AM, Brian Hayes wrote:</div><br class="Apple-interchange-newline"><blockquote type="cite"><div style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; "><div><br></div><div>This isn't really a Lisp issue; it's an IEEE-754 floating-point issue. And you needn't invoke anything as complicated as the floor function to see this behavior. Mere addition and subtraction show it clearly enough:</div><div><br>(setf *read-default-float-format* 'double-float)<br><br>? (- 0.59 0.58)<br>0.010000000000000009<br>? (+ 0.58 0.010000000000000009)<br>0.59<br>? (+ 0.58 0.01)<br>0.59<br>? (= (+ 0.58 0.010000000000000009) (+ 0.58 0.01))<br>t<br>? (= 0.010000000000000009 0.01)<br>nil<br><br>Two excellent papers --- one by Will Clinger, the other by Guy Steele and JonL White --- discuss the decimal-binary conversion problem at the root of this weirdness. The papers are in the 1990 ACM conference on Programming Language Design and Implementation (available through the ACM Digital Library, and probably elsewhere on the net if you go looking).<br><br>Another eye-opening paper (not Lisp-oriented): <br><br> Monniaux, David. 2008. The pitfalls of verifying floating-point computations.<br> <a href="http://hal.archives-ouvertes.fr/docs/00/28/14/29/PDF/floating-point-article.pdf">http://hal.archives-ouvertes.fr/docs/00/28/14/29/PDF/floating-point-article.pdf</a><br><br>The 2008 version of the IEEE-754 standard for floating point includes a decimal floating-point format. This would fix the problem in the specific case (- 0.59 0.58), but of course there are infinitely many real numbers that have no exact decimal representation, so it's not magic.<br><br>On the other hand, exact rationals are heavenly:<br><br>? (floor 59/100 1/100)<br>59<br>0<br><br>All those quaint notions like the associative law and the distributive law --- they actually work! (If only we could ignore the diagonals of squares and the circumferences of circles.)<br><br>Brian Hayes<br><br></div><div><br></div><br><div><div>On 3 Jun 2012, at 5:58 PM, Sudhir Shenoy wrote:</div><br class="Apple-interchange-newline"><blockquote type="cite"><div style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; ">Hi all,<div><br></div><div>There was a recent discussion on the Lispworks mailing list that I thought may be of interest for Clozure CL users. Essentially, someone was complaining about a possible bug in Lispworks where (floor 0.59d0 0.01d0) returned 0.010000000000000009D0. The reasonable expectation (even given the vagaries of inexact floating point arithmetic) that the remainder is always less than the divisor is violated.</div><div><br></div><div>Joe Marshall pointed out that this is probably due to the way floor is being computed (see below) and CLISP seems to do it the "right" way. CCL, which defines floor in terms of truncate gives</div><div><div><br></div><div>? (truncate 0.59d0 0.01d0)</div><div>58</div><div>0.010000000000000009D0</div><div><br></div><div>which does seem strange ...</div><div><br></div><div>Cheers</div><div>Sudhir</div><div><br><div>Begin forwarded message:</div><br class="Apple-interchange-newline"><blockquote type="cite"><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px;"><span style="font-family:'Helvetica'; font-size:medium; color:rgba(0, 0, 0, 1.0);"><b>From: </b></span><span style="font-family:'Helvetica'; font-size:medium;">Joe Marshall <<a href="mailto:jmarshall@alum.mit.edu">jmarshall@alum.mit.edu</a>><br></span></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px;"><span style="font-family:'Helvetica'; font-size:medium; color:rgba(0, 0, 0, 1.0);"><b>Subject: </b></span><span style="font-family:'Helvetica'; font-size:medium;"><b>Re: Double precision Floor bug?</b><br></span></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px;"><span style="font-family:'Helvetica'; font-size:medium; color:rgba(0, 0, 0, 1.0);"><b>Date: </b></span><span style="font-family:'Helvetica'; font-size:medium;">June 4, 2012 4:27:15 AM GMT+09:00<br></span></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px;"><span style="font-family:'Helvetica'; font-size:medium; color:rgba(0, 0, 0, 1.0);"><b>To: </b></span><span style="font-family:'Helvetica'; font-size:medium;">Guy Footring <<a href="mailto:Guy@footring.net">Guy@footring.net</a>><br></span></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px;"><span style="font-family:'Helvetica'; font-size:medium; color:rgba(0, 0, 0, 1.0);"><b>Cc: </b></span><span style="font-family:'Helvetica'; font-size:medium;"><a href="mailto:lisp-hug@lispworks.com">lisp-hug@lispworks.com</a><br></span></div><div style="margin-top: 0px; margin-right: 0px; margin-bottom: 0px; margin-left: 0px;"><span style="font-family:'Helvetica'; font-size:medium; color:rgba(0, 0, 0, 1.0);"><b>Reply-To: </b></span><span style="font-family:'Helvetica'; font-size:medium;">Joe Marshall <<a href="mailto:jmarshall@alum.mit.edu">jmarshall@alum.mit.edu</a>><br></span></div><br>This was bothering me this weekend so I did some more playing around.<div>Computing the quotient part of the answer (the first result) is straightforward:</div><div>you just divide and drop the fractional part. It is the remainder that is</div>
<div>the problem. There are several ways to compute the remainder, and</div><div>it seems that Lispworks does it this way:</div><div><br></div><div><div>(defun floor (n divisor)</div><div> (let* ((quotient (/ n divisor))</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> (integer-part (truncate quotient))</div><div> (remainder (- n (* integer-part divisor))))</div><div> (values integer-part remainder)))</div>
</div><div><br></div><div>But you could use the fractional part of the quotient instead like this:</div><div><br></div><div><div>(defun floor (n divisor)</div><div> (let* ((quotient (/ n divisor))</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> (integer-part (truncate quotient))</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> (remainder (* (- quotient integer-part) divisor)))</div><div> (values integer-part remainder)))</div></div><div><br></div><div>If the arithmetic were exact, these would produce the same answers,</div>
<div>but in floating point, the first produces</div><div><span style="">58 0.010000000000000009d0</span></div><div><span style="">but the second produces</span></div><div><span style="">58 </span><font color="#222222" face="arial, sans-serif">0.009999999999999929d0</font></div>
<div><font color="#222222" face="arial, sans-serif"><br></font></div><div><font color="#222222" face="arial, sans-serif">This second way of computing the remainder satisfies your expectations</font></div><div><font color="#222222" face="arial, sans-serif">about the constraints.</font></div>
<div><font color="#222222" face="arial, sans-serif"> </font></div><div><font color="#222222" face="arial, sans-serif">Incidentally, clisp appears to do it the second way.</font></div><div><font color="#222222" face="arial, sans-serif"><div>
[1]> (floor .59d0 .01d0)</div><div>58 ;</div><div>0.009999999999999929d0</div><div><br></div><div><br></div></font></div><div>-- <br>~jrm<br>
</div>
</blockquote></div><br></div></div>_______________________________________________<br>Openmcl-devel mailing list<br><a href="mailto:Openmcl-devel@clozure.com">Openmcl-devel@clozure.com</a><br><a href="http://clozure.com/mailman/listinfo/openmcl-devel">http://clozure.com/mailman/listinfo/openmcl-devel</a><br></blockquote></div><br></div></blockquote></div><br></div></body></html>