Beefy Boxes and Bandwidth Generously Provided by pair Networks
Do you know where your variables are?
 
PerlMonks  

adding numbers and floating point errors

by smeenz (Sexton)
on Nov 22, 2006 at 21:24 UTC ( #585615=perlquestion: print w/ replies, xml ) Need Help??
smeenz has asked for the wisdom of the Perl Monks concerning the following question:

Hi,

Today I discovered that floating point errors seem to be much more of a problem than I had anticipated.

Specifically, it seems that if I set up a simple loop to add 0.05 to a variable and then display the result, it quickly goes wrong.

For example, the following code:

my $c = 0; for ($x = 0; $x <= 1; $x += 0.05){ $c++; print "[$c] x= $x (".sprintf("%20.40f",$x).") \n"; }

produces this:

[1] x= 0 (0.0000000000000000000000000000000000000000) [2] x= 0.05 (0.0500000000000000027755575615628913510591) [3] x= 0.1 (0.1000000000000000055511151231257827021182) [4] x= 0.15 (0.1500000000000000222044604925031308084726) [5] x= 0.2 (0.2000000000000000111022302462515654042363) [6] x= 0.25 (0.2500000000000000000000000000000000000000) [7] x= 0.3 (0.2999999999999999888977697537484345957637) [8] x= 0.35 (0.3499999999999999777955395074968691915274) [9] x= 0.4 (0.3999999999999999666933092612453037872910) [10] x= 0.45 (0.4499999999999999555910790149937383830547) [11] x= 0.5 (0.4999999999999999444888487687421729788184) [12] x= 0.55 (0.5499999999999999333866185224906075745821) [13] x= 0.6 (0.5999999999999999777955395074968691915274) [14] x= 0.65 (0.6500000000000000222044604925031308084726) [15] x= 0.7 (0.7000000000000000666133814775093924254179) [16] x= 0.75 (0.7500000000000001110223024625156540423632) [17] x= 0.8 (0.8000000000000001554312234475219156593084) [18] x= 0.85 (0.8500000000000001998401444325281772762537) [19] x= 0.9 (0.9000000000000002442490654175344388931990) [20] x= 0.95 (0.9500000000000002886579864025407005101442)

which is... okayish.. because the errors only show up when I force sprintf to display to 40dp... the problem is that even though I'm only adding 0.05 (and not 0.00000000000000000000000005), the error still creeps into the results pretty quickly.

If I change the loop to run from 0 to 4, which is 100 interations, perl tells me at iteration 74 that 3.65 + 0.05 = 3.69999999999999.

I'm not sure what to think of this.. I wouldn't have expected that adding 0.05 to a number 74 times would start to produce errors as soon as the 1st decimal place.

Comment on adding numbers and floating point errors
Select or Download Code
Re: adding numbers and floating point errors
by blazar (Canon) on Nov 22, 2006 at 21:39 UTC
      Thanks.. I had tried -q float and -q precision.. hadn't thought of -q decimal.
Re: adding numbers and floating point errors
by imp (Priest) on Nov 22, 2006 at 21:42 UTC
    BrowserUK provided an excellent explanation for why this occurs a few years ago.

    To (mostly) fix the issue for your script you can use bignum.

    use bignum; my $c = 0; for ($x = 0; $x <= 4; $x += 0.05){ $c++; print "[$c] x= $x (".sprintf("%20.40f",$x).") \n"; }
    [78] x= 3.85 (3.8500000000000000888178419700125232338905) [79] x= 3.9 (3.8999999999999999111821580299874767661095) [80] x= 3.95 (3.9500000000000001776356839400250464677811) [81] x= 4 (4.0000000000000000000000000000000000000000)
      Thanks. That explanation does indeed explain the issue.
Re: adding numbers and floating point errors
by philcrow (Priest) on Nov 22, 2006 at 21:50 UTC
    Unless your interval is a fraction of a power of two or can be represented by a small number of them added together, your number can't be represented exactly in binary. Think about what you are asking your poor binary computer to do. .05 is 1/20. The first number that isn't bigger than your number is 1/32. Subtract that from 1/20. Now look for 1/2**n which is smaller than what's left. Continue until it comes out even or you run out of precision.

    Decimal fractions are rarely representable exactly in computing memory.

    To prevent errors from building up, round more often.

    As for

    3.65 + 0.05 = 3.69999999999999
    I wouldn't say the error is in the first place. It's actually in the last place you showed. Any prior rounding would work. If the difference was significant in the first digit, you could round in the second place without getting the correct answer.

    Decimals tie people in knots even when computers aren't involved. I remember a heated discussion in a math class I took, where the teacher was making no headway in convincing some student that .3 repeating is 1/3.

    Phil

      Thanks to everyone.. I see what the problem is now... I had thought that these sorts of numbers would be represented as whole numbers and then the decimal place added back in for the final result, so 0.0005 + 0.0005 would be the same as adding 5+5 (or 00000101 and 00000101), and then shifting the decimal point appropriately at the output stage, but now that you've explained that the computer is trying to approximate the actual fraction using powers of 2, it makes a lot more sense as to what is happening.
Re: adding numbers and floating point errors
by spiritway (Vicar) on Nov 22, 2006 at 22:00 UTC

    Computers don't handle fractional floating point values very well. Anything other than powers of 1/2 (or their multiples) results in a form that repeats endlessly. This results in truncation or rounding errors, and these will tend to accumulate when you iterate.

    You might want to consider keeping your data in rational form (e.g., bigrat) while you are performing your calculations. Doing this can save you lots of problems with rounding or truncation errors. Once you've finished your calculations you can convert your rationals into floating point format, if you wish, for display.

      Anything other than powers of 1/2 (or their multiples) results in a form that repeats endlessly.

      Not quite. Anything other than sums of powers of 2 result in a form that repeats endlessly.

      Sums of powers of 2 do not repeat endlessly. For example, 0.3125 is not a power of 2, yet it does not repeat because it is a sum of powers of 2 (2-2 + 2-4).

      Note that "powers of 1/2" and "powers of 2" are synonymous.

      Update: Nevermind. What you said and what I said are mathematically equivalent. 0.3125 is a multiple of a power of 2 (5 * 2-4).

Re: adding numbers and floating point errors
by ikegami (Pope) on Nov 22, 2006 at 22:18 UTC

    If you really have a loop and you want to avoid the accumulation of error, you can avoid use an integer as the loop counter, and compute the current value of the real from the integer.

    # Loop over [0.00..1.00] by 0.05 for (my $i = 0; $i <= 100; $i += 5) { my $x = $i / 100; printf("[%d] x= %g (%20.40f)\n", ++$c, $x, $x); }
    # Loop over [0.00..1.00] by 0.05 for (my $i = 0; $i <= 20; $i++) { my $x = $i / 20; printf("[%d] x= %g (%20.40f)\n", $i+1, $x, $x); }
    # Loop over [0.00..1.00] by 0.05 for (0..20) { my $x = $_ / 20; printf("[%d] x= %g (%20.40f)\n", $_+1, $x, $x); }

    All of the above output the following:

    [1] x= 0 (0.0000000000000000000000000000000000000000) [2] x= 0.05 (0.0500000000000000030000000000000000000000) [3] x= 0.1 (0.1000000000000000100000000000000000000000) [4] x= 0.15 (0.1499999999999999900000000000000000000000) [5] x= 0.2 (0.2000000000000000100000000000000000000000) [6] x= 0.25 (0.2500000000000000000000000000000000000000) [7] x= 0.3 (0.2999999999999999900000000000000000000000) [8] x= 0.35 (0.3499999999999999800000000000000000000000) [9] x= 0.4 (0.4000000000000000200000000000000000000000) [10] x= 0.45 (0.4500000000000000100000000000000000000000) [11] x= 0.5 (0.5000000000000000000000000000000000000000) [12] x= 0.55 (0.5500000000000000400000000000000000000000) [13] x= 0.6 (0.5999999999999999800000000000000000000000) [14] x= 0.65 (0.6500000000000000200000000000000000000000) [15] x= 0.7 (0.6999999999999999600000000000000000000000) [16] x= 0.75 (0.7500000000000000000000000000000000000000) [17] x= 0.8 (0.8000000000000000400000000000000000000000) [18] x= 0.85 (0.8499999999999999800000000000000000000000) [19] x= 0.9 (0.9000000000000000200000000000000000000000) [20] x= 0.95 (0.9499999999999999600000000000000000000000) [21] x= 1 (1.0000000000000000000000000000000000000000)
Re: adding numbers and floating point errors
by traveler (Parson) on Nov 22, 2006 at 22:24 UTC
Re: adding numbers and floating point errors
by mattr (Curate) on Nov 23, 2006 at 16:57 UTC
    Well, I just saw the exact same question posted somewhere else a couple days ago so it is interesting you came up with the same exact question. Anyway I found some more info about it (see end of my post).

    If you just print it without using sprintf, there is no problem. If you do use sprintf you are intentionally trying to see how the internal math library is working. IIRC 13 significant digits should be fine though your example looks like 15 digits. Maybe that's a 16 bit precision C float. It happens due to the math routines and IIRC that floating point numbers are generally represented in scientific notation, which means it only guarantees a fixed number of digits precision. Also shortcuts are made for speed sometimes. If you needed as many digits as you are requesting you should be using one of the modules that does that for you as other people have noted for arbitrary precision or more precision (I remember BigNum myself), PDL or others might help too.

    If you had a problem showing up without using sprintf I'd be a little more worried. Conceivably a 64bit computer might be using higher precision routines in their base libraries, but I could not say myself.

    Bottom line is, this is a well known artifact common to computers and is not a bug in Perl.

    HOWEVER, You want to be really careful about numbers the computer spits out when you are doing anything important like financial, scientific or statistical calculations. A rounding error or something like this can get ugly. Some interesting info along here.. also I am not sure about whether there is a perl facility to (looks like Math::BigFloat would) provide access to 80 bit extended precision values.

    Wikipedia on floating point addition
    General Decimal Arithmetic website at IBM and the FAQ , especially the part on how many digits precision are needed for decimal arithmetic. (It says business requirements are typically 20-30 digits but interest rate calculations could require over 2000 digits.)
    What Every Computer Scientist Should Know About Floating-Point Arithmetic (html) and pdf (neat.)
    perlnumber (the pod on perl's numerics)

    UPDATED:Boneheaded mistake made unusable links! Sorry. Now fixed.

Re: adding numbers and floating point errors
by Anonymous Monk on Nov 23, 2006 at 19:36 UTC
    I htink this looks horrible, but is not :-)
    You noticed the first time the error gets less than zero! If you calculate the error, you'll see it's around 10E-15 at this point (at last on my pc).
    #!/usr/bin/perl my $x=0.0; my $y=0; my $e_min=1.0; my $e_max=0.0; my $n=0; my $e; for ( 1 .. 999 ) { $x += 0.05; $y += 5; $e = abs($x - $y / 100); $n++ if $e == 0.0; $e_min = $e if $e > 0.0 && $e < $e_min; $e_max = $e if $e > $e_max; printf "%4d x = %.30f e = %E\n", $_, $x, $e; } printf "min %E max %E exact %d\n", $e_min, $e_max, $n; result 998 x = 49.899999999999302247033483581617 e = 6.963319E-13 999 x = 49.949999999999299404862540541217 e = 7.034373E-13 min 2.775558E-17 max 7.034373E-13 exact 16
    You'll find the lines arround 190 interesting:
     189 x = 9.449999999999999289457264239900 e = 0.000000E+00
     190 x = 9.500000000000000000000000000000 e = 0.000000E+00
     191 x = 9.550000000000000710542735760100 e = 0.000000E+00
    

Log In?
Username:
Password:

What's my password?
Create A New User
Node Status?
node history
Node Type: perlquestion [id://585615]
Approved by Corion
Front-paged by tye
help
Chatterbox?
and the web crawler heard nothing...

How do I use this? | Other CB clients
Other Users?
Others examining the Monastery: (13)
As of 2014-09-17 15:01 GMT
Sections?
Information?
Find Nodes?
Leftovers?
    Voting Booth?

    How do you remember the number of days in each month?











    Results (86 votes), past polls