3
$\begingroup$

I encounted the following precision issue. I define c = a + b, but the result for 'c - a - b' does not equal zero. I am not sure how to resolve this.

'''

In[1]:= 
a = N[Exp[3]];
b = N[Exp[4]];
c = a + b;
c - a - b

Out[4]= 7.10543*10^-15

'''

$\endgroup$
6
  • 1
    $\begingroup$ What happens when you do Chop[%] on the final output? does it become zero? see help on Chop or you can use infinite precision. a = N[Exp[3], Infinity]; b = N[Exp[4], Infinity]; and you will get exact zero without doing Chop $\endgroup$
    – Nasser
    Sep 27 at 7:47
  • 1
    $\begingroup$ But note that N[Exp[3], Infinity] just keeps it as Exp[3] and will not make it real number. Once you convert things to approximate numbers, you'll get approximate numbers as result of computation. $\endgroup$
    – Nasser
    Sep 27 at 7:56
  • 1
    $\begingroup$ notice also that N[...] use machine precision (unless you give second argument to N). on same PC, I typed the same thing on Matlab, and got same output: a=exp(3); b=exp(4); c=a+b; c-a-b gives 7.1054e-15 This is the nature of using real numbers. !Mathematica graphics $\endgroup$
    – Nasser
    Sep 27 at 8:05
  • 2
    $\begingroup$ You can avoid some of the problems by using c - (a + b) instead of c - a - b. Now you will get 0. as output. This is because the order makes difference in terms of real numbers computations (large number - small number, vs. difference of two numbers close to each others). see What Every Computer Scientist Should Know About Floating-Point Arithmetic $\endgroup$
    – Nasser
    Sep 27 at 8:11
  • 1
    $\begingroup$ Everything is indeed correct. You can see the same phenomenon even on a scientific calculator. As @Nasser said, this is the nature of floating point computations. $\endgroup$
    – yarchik
    Sep 27 at 9:51
1
$\begingroup$

This phenomenon occurs for simple cases too:

0.3 + 1.1 - 1.4
2.22045*10^-16`

These numbers do not have finite length in binary and so they're truncated when expressed as an inexact number. We can see this with RealDigits:

RealDigits[1.1, 2]
{{1, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0,
  1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 
  0, 0, 1, 1, 0, 0, 1, 1, 0, 1, 0}, 1}

The identity $1.1_{10} = 1.0001100110011_2...$ gets truncated to $53$ digits to the right of the decimal (on a $64$ bit machine). The same happens for $0.3$ and $1.4$ — and in fact only dyadic rationals can be represented exactly in floating point arithmetic.

We can see the error by expressing the truncated reals as rational numbers:

SetPrecision[0.3 + 1.1, ∞] - SetPrecision[1.4, ∞]
1/4503599627370496
N[%]
2.22045*10^-16
$\endgroup$

Not the answer you're looking for? Browse other questions tagged or ask your own question.