Problem concerning floatin point representation.

I'm trying to fully understand how Matlab represents numbers. Therefore, I compared some numbers using "==". I was able to predict its outcome for several cases. The problem I have is that when I compare 0.1*5 and 0.5 I get that they are the same. I can't understand this because 0,1 can't be represented exactly, but 5 and 1/2 can. How can these two be equal?

PS: I tried it also with 0.1*1, 0.1*2, 0.1*3, 0.1*4 and 0.1, 0.2, 0.3, 0.4. My results where equal, equal, not equal and equal.

Answer 1

[John D'Errico]

There are tiny somewhat unpredictable errors in the least significant bits of a number. Not random, but unpredictable.

As you found out, sometimes such computations yield the same results, sometimes not. Hmm, just what does unpredictable mean? That sounds about right.

So never trust the least significant bits of a floating point number. At least not until you know fully how and EXACTLY when you can trust those results. You are not there yet. There are occasions when you can trust the least significant bits. For example, 1/2 is exactly 0.5. Think about it. Why?

Essentially, what you need to understand is why 0.1 is not representable as a double. It is the same fundamental reason why 1/3 is not representable as a decimal. Again, think about it. Remember that double precision numbers are represented internally using a variation of BINARY. Numbers are not stored in decomal form in the standard IEEE form used by most (but not all) computer languages. Now, think about why 1/10 is NOT representable exactly, yet 1/2 is. Once you fully understand this concept, you will be at least started on your journey.

So consider this infinite sum:

1/10 = sum(2.^[-4 -5 -8 -9 -12 -13 -16 -17 -20 -21 -24 -25 -28 -29 -32 -33 -36 -37 -40 -41 -44 -45 -48 -49 -52 -53 -55 ...])

Note the ... at the end there. There are infinitely many terms in that series. But as a double precision number, MATLAB essentially stops at the 2^-55 term. It must stop somewhere, after all.

Now, how many terms are needed in a similar representation of 0.5, written in terms of a sum of powers of 2?

[Iron Degryse]

1/2 is equal to .1 in binary.

0,1: .0001100..

0,2: .001100.

0,4: .01100..

So given that Matlab can exactly calculate multiplications of two-powers, it's logical that:

error(0.1)=1*error(0,1)

error(0.2)=2*error(0,1) one power difference

error(0.4)=4*error(0,1) two power difference

For multiplication with 3 and 5 I would have to take an extra multplication error into account and it just doesn't make sense that these errors 'dissapear' for 5*0,1.

[John D'Errico]

Re-read what I wrote. Think again.

0.1 is not representable in binary. You get a garbage result down at the end.

sprintf('%0.55f',0.1)
ans =
    '0.1000000000000000055511151231257827021181583404541015625'

This is effectively the number stored in MATLAB when you type 0.1.

You still think you can trust those least significant bits.

sprintf('%0.55f',0.3)
ans =
    '0.2999999999999999888977697537484345957636833190917968750'
>> sprintf('%0.55f',3*0.1)
ans =
    '0.3000000000000000444089209850062616169452667236328125000'

So, both numbers are 0.3. Or are they?

Use a tolerance whenever you compare floating point numbers.

[Iron Degryse]

I understand what you are saying and I already knew that you preferably use a tolerance, but the calculations Matlab makes when calculating 0,1*5 are:

0.1000000000000000055511151231257827021181583404541015625*5.0000000000000000000000000000000000000000000000000000000

This clearly doesn't equal 0.5 but Matlab says it does. Why?

[Stephen23]

"This clearly doesn't equal 0.5 but Matlab says it does. Why?"

You are are focussing so much on the input value that you are totally forgetting about the value of the result. You assume that the result will be exactly five times that input value: i.e. that there exists a binary floating point value that represents exactly five times that input value and that this must be the output value. Clearly your test indicates that your assumption is incorrect.

Below is an example that might help (and before anyone complains, this is not presented as an exact representation of the floating point error handling of the exact operations being performed on those exact floating point numbers... it is just intended to illustrate the behavior).

Take a look at the input value, two binary floating point representable values around 0.5, and five times the pesky trailing non-zero digits that you are confused about:

>> fprintf('%.40f\n',0.1) % input value
0.1000000000000000055511151231257827021182
>> Y = 0.0000000000000000055511151231257827021182; % just the trailing part
>> fprintf('%.40f\n',[0.5,0.5+eps(0.5),5*Y])
0.5000000000000000000000000000000000000000
0.5000000000000001110223024625156540423632
0.0000000000000000277555756156289135105908

Ignoring the 0.5, which of those two binary floating point representable values would five times that trailing part be closest to? Which one do you think the result of 0.1*5 would be? (hint: the closest). Does exactly five times the input value even exist as a floating point number, even though you assumed that it must? (hint: no)

[Iron Degryse]

Okay, so I conclude that in general the calculations include new errors and I won't be able to predict in which direction these errors will accumulate. Thanks for your instructing remarks