ben@catnip.UUCP (12/29/86)
I have noticed what I believe is incorrect behavior on the part of the
floating point routines on my machine. The problem is best demonstrated
with a small program:
#include <stdio.h>
main(){
double x;
for (x=0.;x <= 3.;x += .2)
printf("%f\n",x);
printf("Final value: %f less 3 is %e\n",x,x-3.);
}
When this program is run, it produces:
0.000000
0.200000
0.400000
0.600000
0.800000
1.000000
1.200000
1.400000
1.600000
1.800000
2.000000
2.200000
2.400000
2.600000
2.800000
Final value: 3.000000 less 3 is 4.440892e-16
It seems to be going through the loop one time less than I expected.
I don't understand this behavior. 3.0 can certainly be represented
precisely in binary, so only the .2 should subject to roundoff error.
The crazy thing is, it appears that the floating point routines are
rounding this number up instead of down! It is the only explanation
I can think of to explain this output. Is roundup considered acceptable
behavior for floating point routines? Or am I somehow misinterpreting
the problem?
--
Ben Broder
{ihnp4,decvax} !hjuxa!catnip!ben
{houxm,clyde}/henry@utzoo.UUCP (Henry Spencer) (12/31/86)
> ... Is roundup considered acceptable > behavior for floating point routines? Or am I somehow misinterpreting > the problem? The precise details of rounding are very machine-dependent and you are not entitled to make assumptions about them. Don't forget that there is rounding in the additions as well as in the representation of "0.2". It is a grave mistake to use floating-point arithmetic to control a loop which has to iterate an exact number of times. -- Henry Spencer @ U of Toronto Zoology {allegra,ihnp4,decvax,pyramid}!utzoo!henry
herman@ti-csl.CSNET (Herman Schuurman) (12/31/86)
> .... Is roundup considered acceptable > behavior for floating point routines? Or am I somehow misinterpreting > the problem? > > Ben Broder Your floating point routines seem to be using the standard IEEE double precision floating point format, with the rounding mode set to "round to nearest". This means that they are rounding to the nearest representable result, which is quite different from rounding up all the time. The accumulated rounding error shows up in the final bit of the result, which makes it 3.000000000000004e0, rather than 3.0. The big problem is that .2 is not exactly representable in a binary floating point format. It therefore always suffers from rounding errors. If you printed the results of your example more precisely, you would come up with the following results: 0.0e0 0.2e0 0.4e0 0.6000000000000001e0 0.8e0 1.0e0 1.2e0 1.4e0 1.5999999999999999e0 1.7999999999999999e0 1.9999999999999999e0 2.1999999999999999e0 2.4e0 2.6e0 2.8000000000000003e0 3.0000000000000004e0 This clearly shows the rounding errors during addition. Notice that some of the results are too low, whereas other results are to high, indicating a "round to nearest", rather than a "round to positive infinity" (as IEEE calls rounding up). -- Herman Schuurman ARPA: herman%TI-CSL@CSNET-RELAY.ARPA Texas Instruments Inc. CSNET: herman@TI-CSL PO Box 226015 M/S 238 USENET: {ut-sally,convex!smu,texsun,rice}!ti-csl!herman Dallas, Texas 75266 VOICE: (214) 995-0845
ark@alice.UUCP (12/31/86)
In article <442@catnip.UUCP>, ben@catnip.UUCP writes: > > I have noticed what I believe is incorrect behavior on the part of the > floating point routines on my machine. [he includes a sample program that increments a variable, initially 0, by 0.2 until it passes 3; then prints the variable minus 3] > When this program is run, it produces: > 0.000000 > 0.200000 [and so on] > 2.600000 > 2.800000 > Final value: 3.000000 less 3 is 4.440892e-16 > > It seems to be going through the loop one time less than I expected. > I don't understand this behavior. 3.0 can certainly be represented > precisely in binary, so only the .2 should subject to roundoff error. > The crazy thing is, it appears that the floating point routines are > rounding this number up instead of down! It is the only explanation > I can think of to explain this output. Is roundup considered acceptable > behavior for floating point routines? Or am I somehow misinterpreting > the problem? You're not misinterpreting it very much. You are quite right: 0.2 cannot be represented exactly in floating-point. So what number you get is determined by the exact characteristics of your system. You will probably get one of the two nearest representable numbers; from your results it looks like you are getting the one slightly larger. Thus as you add multiple instances of that value to 0, you tend to get numbers that are slightly larger than multiples of 0.2 When you think you're at 3, you're really at 3 plus a small number. This might be easier to understand if you print the value of (0.2 * 5) - 1 This "should" be zero, but it won't.
jsdy@hadron.UUCP (Joseph S. D. Yao) (01/02/87)
One of the few values of teaching Fortran to programmers
was that invariably, eventually, students would be taught
about round-off error. In C, we rarely think that someone
might use a float as a loop control variable. If one needs
a float value iterated, one should use the old Fortran trick:
int i;
float f;
for (i = 0; i < 30; i += 2) {
f = i / 10.0;
...
}
[Do n o t flame me for not using 0:15:1;;f=i/5.0;. This
valuation makes it more obvious to the casual reader what
conversion was done; and the arithmetic is obvious!]
--
Joe Yao hadron!jsdy@seismo.{CSS.GOV,ARPA,UUCP}
jsdy@hadron.COM (not yet domainised)jc@piaget.UUCP (John Cornelius) (01/02/87)
The value .2 is an irrational fraction in binary. Adding the representation of (irrational) .2 to (rational) +3. results in a number slightly higher than (rational) +3.0, thereby terminating the loop. Henry Spencer is quite right about using floating point numbers as loop indices, it is a _bad_ practice and has been for nearly 30 years. -- John Cornelius (...!sdcsvax!piaget!jc)