How not to write a loop

[Hume Smith]

ok at quintus.UUCP (Richard A. O'Keefe) wrote in <712 at cresswell.quintus.UUCP>:
>In article <832 at unmvax.unm.edu>, mike at turing.UNM.EDU (Michael I. Bushnell) writes:
>> I see NOTHING which precludes:
>> 	float x;
>> 	for (x = 0; x < 20; x += .5) printf("%f\n", x);
>> The output would, of course, be
>> 	0.0
>> 	0.5
>> 	...
>> 	19.5
>
>You should not be surprised to see as the last value
>	19.99999
>I just tried the very similar loop
>	for (x = 0; x < 21; x += .3)	/* .3 divides 21 exactly */
>By analogy with your "of course", the last output should obviously
>be 20.7.  In fact, when I tried it just now, the last output was
>	20.999994

I would be quite surprised to see 19.99999 come out in the former
loop.  .5 is exactly representable in floating binary, and should cause
no problem as an increment.  .3 is not representable in a finite binary
fraction (.01001 something, it repeats eventually), and cannot be summed
up to get an integral value.

If the machine did its arithmetic in decimal (Texas Intstruments?  Sharp?)
then the .3 increment would again not be a problem, since the representation
is exact again.

This, however, is the point.  The behaviour of a particular increment
is unreliable.

(Not to be insulting at all, but this may be why our Computer Science
department makes its majors take our Numerical Methods Course.  The
very first thing we do is Machine Epsilon and roundoff error.  This
comment is just to make enough text to keep inews happy.)
-- 
  Hume Smith                        UUCP:
Math Department     {uunet|watmath|utai|garfield|mnetor}!
   Acadia U                   dalcs!aucs!850347s