strange behavior with float value multiplication

[Yonghee Kim]

I wrote simple code like this

----------------------

a = 0.2

for i in 1:10

  println(a * i)

end

---------------------------

and what I got is not 0.2, 0.4, 0.6, 0.8, 1.0 .... 

but this.

0.2
0.4
0.6000000000000001
0.8
1.0
1.2000000000000002
1.4000000000000001
1.6
1.8
2.0

println(0.2 * 3) does the same thing. not 0.6 but 0.6000000000000001

does anyone know why this happens? 

[Pranit Bauva]

Dear Yonghee Kim,

Floating point calculations are not exact because the decimal point is
represented by only a fixed number of binary digits. The more is the
bits allocated more will be the precision. This explanation would be
clearer if you see how to convert floating point integers to binary.
http://stackoverflow.com/questions/3954498/how-to-convert-float-number-to-binary

The series will be more and more accurate if more number of terms are
added and will each added term it error will decrease by 2^-(no of
terms). Furthermore this is the problem with almost every programming
language.

Regards,
Pranit Bauva

[Yonghee Kim]

Dear Pranit Bauva 

Thanks for your kind explanation!

your answer helped me a lot. 

[Scott Jones]

That is simply because you are using binary floating point (Float64 in Julia), instead of decimal floating point.  There is a very nice package, DecFP.jl, (kudos to Steven Johnson, @stevengj on GitHub), that wraps the Intel decimal arithmetic library.

You can do `Pkg.add("DecFP") ; using DecFP` (you may get a lot of deprecation warnings, esp. on v0.5, you can also do a `Pkg.checkout("DecFP") ; Pkg.build("DecFP")` to eliminate all but one warning on v0.5), and then you can do:

julia> a = d"0.2"
+2E-1

julia> for i in 1:10 ; println(a * i) ; end
+2E-1
+4E-1
+6E-1
+8E-1
+10E-1
+12E-1
+14E-1
+16E-1
+18E-1
+20E-1

The output is a bit funny, but the answers are exact, unlike most of the answers in the binary floating point case you showed.

Note: the most heavily used language in the world (at least up to the early 2000s, I'm not sure what the stats would be like now), which was the 2nd ANSI standardized language, Cobol,

originally only had decimal arithmetic (Modern Cobols also support IEEE binary floating point).

A number of other heavily used (but not much talked about) languages, such as M/MUMPS, which happened to be the 3rd language standardized by ANSI, even before ANSI C.

The most used version of M/MUMPS these days is Caché ObjectScript - which adds many things, including objects and IEEE binary floats.

The Pick/MultiValue family of database languages also use decimal floating point.

There is now the IEEE 754-2008 standard for decimal floating point, which the DecFP.jl package conforms to (implementing 3 of the formats - there are some more formats which I'm trying to write a package for, which are also implemented in hardware on IBM POWER architecture processors).

There are a number of cases where it really is better to use decimal arithmetic, in order to not get in trouble with the tax man! (Lots of currency operations have this same issue)

Here's a real simple example: you have a store, where you sell something for 70 cents, and there is a 5% sales tax.  How much do tax to you charge the customer?

julia>tax = .7 * 0.05
0.034999999999999996
julia> round(tax,2)
0.03
julia> tax = d".70" * d"0.05"
+350E-4
julia> round(tax,2)
+4E-2

As you can see, with binary floating point, you get the wrong answer (0.03), but with decimal floating point, you get the correct answer (the one the tax authorities will want you to give them).

Depending on your particular use case, you can either use the default IEEE binary floating point numbers in Julia, or the DecFP.jl package.

Scott

[Stefan Karpinski]

I would recommend against using decimal floating-point unless you really must. On almost any hardware you're likely to use, binary floating-point arithmetic will be an order of magnitude faster than decimal floating-point arithmetic. Moreover, if you use decimal floating-point, you will chronically find yourself unable to use standard numerical libraries without converting your data first.

Instead, I would strongly recommend learning a bit about how binary floating-point works. The basic fact you need to know is that the Float64 type can only represent real values of the form:

n*2^k where |n| < 2^53 and -1074 ≤ k ≤ 971 (both integers).

Some basic but potentially unintuitive facts that follow from this:

1. Only rational values whose denominators are powers of two can be represented as binary floating-point values exactly. Most rational numbers are not of this form, of course. When you write `0.1` in Julia (and most programming languages), it does not mean 1/10. Instead, it means the number of the above form that is closest to 1/10, namely:

3602879701896397*2.0^-55 =

0.1000000000000000055511151231257827021181583404541015625.

This value is printed as "0.1" by convention since that's the shortest decimal input that gives you this value, and therefore in some sense the most likely way to input it. However, there are many other decimal inputs that produce the same Float64 value, such as 0.0999999999999999987 and 0.10000000000000001.

2. Unlike rational numbers, which are dense in all intervals of the real number line, floating-point values are not at all dense: there is a gap between each floating-point value and the next representable value. The eps(x) function gives the distance to the next representable value:

eps(1.0) = 2.220446049250313e-16

eps(12.34) = 1.7763568394002505e-15

eps(1e25) = 2.147483648e9.

The gap between moderately sized values like 1.0 and 12.34 is very small, but the gaps between large numbers like 1e25 is very large. In general, the gap between values of size 10^n is around 10^(n-16). That is because Float64s have 53 bits of precision and log10(2^53) ≈ 16.

[Scott Jones]

On Monday, December 28, 2015 at 1:10:25 PM UTC-5, Stefan Karpinski wrote:

I would recommend against using decimal floating-point unless you really must. On almost any hardware you're likely to use, binary floating-point arithmetic will be an order of magnitude faster than decimal floating-point arithmetic. Moreover, if you use decimal floating-point, you will chronically find yourself unable to use standard numerical libraries without converting your data first.

True - I was not advocating using decimal floating point *unless* you must, however, there are a number of use cases where that you do need to use decimal floating point,

and you might find that using decimal floating point can actually be faster than using larger precision binary floating point if you are trying to avoid some of the conversion issues that exist with binary floating point.  For example, in some benchmarks I ran comparing `DecFP.jl` to `BigFloat`, DecFP was about 2-3 times faster than 256-bit BigFloat.

(DecFP is very nice, in that it uses bitstypes, but BigFloat causes a number of allocations for every operation - which may indicate that some work could be done on Julia's BigFloat wrapping, to use bitstypes, instead of having to have both an allocation for the type, and the array of words).

julia> f1() = (a = big"0.0" ; for i = 1:10000 ; a += big"0.1" ; end ; a)

f1 (generic function with 1 method)

julia> f1()

9.999999999999999999999999999999999999999999999999999999999999999999999999859301e+02

julia> f2() = (a = d"0.0" ; for i = 1:10000 ; a += d"0.1" ; end ; a)

f2 (generic function with 1 method)

julia> f2()

+10000E-1

julia> @benchmark f1()

================ Benchmark Results ========================

     Time per evaluation: 2.55 ms [139.77 μs, 4.95 ms]

Proportion of time in GC: 4.52% [0.00%, 16.41%]

        Memory allocated: 1015.63 kb

   Number of allocations: 30000 allocations

       Number of samples: 100

   Number of evaluations: 100

 Time spent benchmarking: 0.31 s

julia> @benchmark f2()

================ Benchmark Results ========================

     Time per evaluation: 792.35 μs [703.65 μs, 881.05 μs]

Proportion of time in GC: 0.00% [0.00%, 0.00%]

        Memory allocated: 0.00 bytes

   Number of allocations: 0 allocations

       Number of samples: 100

   Number of evaluations: 100

 Time spent benchmarking: 0.10 s

julia> f1m() = (a = big"0.0" ; for i = 1:10000 ; a *= big"0.1" ; end ; a)

f1m (generic function with 1 method)

julia> f2m() = (a = d"0.0" ; for i = 1:10000 ; a *= d"0.1" ; end ; a)

f2m (generic function with 1 method)

julia> @benchmark f1m()

================ Benchmark Results ========================

     Time per evaluation: 2.03 ms [0.00 ns, 4.06 ms]

Proportion of time in GC: 3.73% [0.00%, 14.73%]

        Memory allocated: 1015.63 kb

   Number of allocations: 30000 allocations

       Number of samples: 100

   Number of evaluations: 100

 Time spent benchmarking: 0.23 s

julia> @benchmark f2m()

================ Benchmark Results ========================

     Time per evaluation: 756.59 μs [678.94 μs, 834.24 μs]

Proportion of time in GC: 0.00% [0.00%, 0.00%]

        Memory allocated: 0.00 bytes

   Number of allocations: 0 allocations

       Number of samples: 100

   Number of evaluations: 100

 Time spent benchmarking: 0.10 s

[Ismael Venegas Castelló]

You can see the binary representation of those `Float64`s with the `bits` function:

julia> for i in 1:10

         @show bits(a * i)

       end

bits(a * i) = "0011111111001001100110011001100110011001100110011001100110011010"

bits(a * i) = "0011111111011001100110011001100110011001100110011001100110011010"

bits(a * i) = "0011111111100011001100110011001100110011001100110011001100110100"

bits(a * i) = "0011111111101001100110011001100110011001100110011001100110011010"

bits(a * i) = "0011111111110000000000000000000000000000000000000000000000000000"

bits(a * i) = "0011111111110011001100110011001100110011001100110011001100110100"

bits(a * i) = "0011111111110110011001100110011001100110011001100110011001100111"

bits(a * i) = "0011111111111001100110011001100110011001100110011001100110011010"

bits(a * i) = "0011111111111100110011001100110011001100110011001100110011001101"

bits(a * i) = "0100000000000000000000000000000000000000000000000000000000000000"

[Erik Schnetter]

You can use rational numbers instead of floating-point numbers. For example writing 7//10 * 5//100 yields the exact value 7//200.

-erik

On Tue, Dec 29, 2015 at 12:20 AM, Yonghee Kim <maste...@gmail.com> wrote:

what I've been doing is kind of similar to tax problem (game economy simulation)

for now,turning floating number to integer is enough to solve all my headaches. 

But I will remember your recommendation when I ran into more complicated case.

Thank you 

2015년 12월 28일 월요일 오후 10시 38분 21초 UTC+9, Scott Jones 님의 말:

--

Erik Schnetter <schn...@gmail.com> http://www.perimeterinstitute.ca/personal/eschnetter/