0

When I compute the following in Matlab

myeps = abs(3*(4/3-1)-1);
format long e
eps_myeps = [eps ; myeps]

The output is as follows:

eps_myeps =

     2.220446049250313e-16
     2.220446049250313e-16

Why is myeps not 0? Why does this not hold when the base is 3 instead of 2?

Questioner5453
  • 1
  • 1
  • `myeps` is not zero because of floating-point inaccuracies. Same reason why `sin(pi)` does not give zero. More information [here](https://stackoverflow.com/q/686439/2586922). The second part of your question, _Why does this not hold when the base is 3 instead of 2_, is unclear. What do you mean exactly? Please edit the question accordingly. Note that powers of 2 (not too large or too small) do have an exact representation, as described in the link; perhaps that explains your second issue – Luis Mendo Nov 14 '22 at 19:51

2 Answers2

2

Code 4/3 is 4/3 in math. 4/3 is not exactly encodable as a floating point number. Most floating point numbers are dyadic rationals (an integer times some power of 2) and a nearby value is used.*1 Much like we can not write 4/3 exactly in decimal, only 1.3333333 and stop after so many digits.

In this case, The subtraction is expected to be exact as well as the multiplication and final subtraction. Yet the first quotient is not 4/3 and so the final result might not be 0.0.

*1
decimal 1.3333333333333332593184650249895639717578887939453125
hex 0x1.5555555555555

chux - Reinstate Monica
  • 143,097
  • 13
  • 135
  • 256
  • Why does this error only occur in base 2 and 10? Why doesn't it occur in base 3? – Questioner5453 Nov 14 '22 at 21:32
  • 1
    @Questioner5453 _It_ does not occur in base 3 as the division here was 3 - and so exact. _It_ does occur in base 3 if the division here was 5. IOWs, is the divisor divisible by the base? – chux - Reinstate Monica Nov 14 '22 at 21:34
  • 1
    @Questioner5453 Try doing `3*(4/3-1)-1` on paper with using only decimal numbers of finite size, perhaps 5 digits: 4/3 --> 1.3333. 1.3333 - 1.0000 --> 3.3330. 3.0000* 3.3330e-1 --> 9.9990e-1. 9.9990e-1 - 1.0000 --> -1.0000e-4 and not zero. Same thing happens with base 2 floating point numbers - just more digits. Note: 2.220446049250313e-16 is power(2, -52). – chux - Reinstate Monica Nov 14 '22 at 21:42
0

I feel the OP is also searching for an answer for "How could obtain a 0 as an answer", besides "Why is myeps not 0". So here is the answer. You'll need the Symbolic Math Toolbox for MATLAB.

You can create symbolic numbers by using sym. Symbolic numbers are exact representations, unlike floating-point numbers.

So if you type 3*(sym(4)/sym(3)-1)-1 into MATLAB (once you got the toolbox), the answer will be exact

0

Just mind the sym(4)/sym(3) part. If you try sym(4/3), MATLAB will get a floating-point first, then try to convert it to symbolic. That will lose precision and not produce 0 as an answer.

X Zhang
  • 1,081
  • 7
  • 17
  • 1
    Correction: Actually 3*(sym(4/3)-1)-1 will produce exactly 0 in MATLAB because the sym engine has an algorithm to convert floating point values calculated from p/q to exactly p/q for "modest-sized integers p and q." However, one shouldn't rely on this and converting the integers individually first as shown above or using a string for the integer division like '4/3' is the way to go. – James Tursa Nov 15 '22 at 06:14
  • @JamesTursa You are right, then you are right again. By that I mean 3*(sym(4/3)-1)-1 indeed produce 0, but sym(3*(4/3-1)-1) will produce some a non-zero fraction. So, better safe than unpredictable. – X Zhang Nov 15 '22 at 07:01