8
$\begingroup$

I want to understand how small numerical errors propagate through a computation and may drastically change the final result. To this end, I consider a toy world in which every number is limited to exactly 5 significant digits.

The output of

$MaxPrecision = 5;
$MinPrecision = 5;

x = SetPrecision[10^(-10), 5];
x // FullForm
y = 1 + x;
y // FullForm
y - 1 // FullForm

is

1.`5.*^-10

1.`5.

1.`5.*^-10

The last line puzzles me. The FullForm of y shows that the Precision of 5 significant digits is not sufficient to resolve the small deviation x from 1. However, y - 1 still yields the exact result, instead of y - 1 = 0.

If one writes explicitly

1.`5. - 1 // FullForm

this yields, of course,

0

as it should.

Does this mean that FullForm[y] hides something?

$\endgroup$
  • 2
    $\begingroup$ I think you need NumericalMath`$NumberBits[y] in order to see the guard bits that I guess are otherwise hidden in e.g. FullForm and InputForm. Maybe there is a cleaner way, I'm not sure. $\endgroup$ – Daniel Lichtblau Feb 2 '17 at 23:40
  • $\begingroup$ @Daniel Why not post an answer? $\endgroup$ – Mr.Wizard Feb 2 '17 at 23:52
  • $\begingroup$ @Mr.Wizard It's at best an answer/2 ... $\endgroup$ – Daniel Lichtblau Feb 3 '17 at 0:08
  • $\begingroup$ You might want to look at the Computer Arithmetic package, if you're interested in this sort of thing. $\endgroup$ – J. M. is away Feb 3 '17 at 4:22
8
$\begingroup$

This is at best a partial answer.

I don't know in full detail what is going on with FullForm or InputForm but I will venture to guess they are influenced in some way by the $XXXPrecision settings.

$MaxPrecision = 5;
$MinPrecision = 5;
InputForm[x = SetPrecision[10^(-10), 5]]
(* Out[6]//InputForm=
1.`5.*^-10 *)
InputForm[y = 1 + x]
(* Out[7]//InputForm=
1.`5. *)
InputForm[y - 1]
(* Out[8]//InputForm=
1.`5.*^-10 *)

There do seem to be guard bits retained, and they can be found using NumericalMath`$NumberBits.

NumericalMath`$NumberBits[y]

(* Out[9]= {1, {1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 
  0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 0, 1,
   1, 1, 1, 1, 0, 0, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 0, 0, 1, 
  1, 1}, 1} *)

NumericalMath`$NumberBits[y - 1]

(* Out[10]= {1, {1, 1, 0, 1, 1, 0, 1, 1, 1, 1, 1, 0, 0, 1, 1, 0, 1}, {1, 
  1, 1, 1, 1, 1, 0, 1, 1, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
   0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
  0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
   0, 0, 0, 0, 0, 0, 0, 0, 0}, -33} *)

Why are bits considered "bad" in y now regarded as "good" in y-1? I do not know.

If you make the perturbation value sufficiently small, the precision tracking will discard it.

x = SetPrecision[10^(-20), 5];
y = 1 + x;
InputForm[{x, y, y - 1}]

(* Out[45]//InputForm=
{1.`5.*^-20, 1.`5., 0} *)

Let me add two quotes from Jerry Keiper, "Mathematica Numerics: Controlling the Effects of Numerical Errors in Computation", 1992 Mathematica Developer Conference:

When SetPrecision is applied to an exact number (i.e., an integer or a rational number), a bignum is created with a few more bits than are requested and the unknown bits are so indicated. These explicit, hidden bits are in fact used in the subsequent calculations, but the boundary between the known and unknown digits is calculated for each result.

and

By changing the values of MinPrecision and MaxPrecision to be the same, say 100 you can get fixed, high-precision arithmetic. You should be aware that there is a difference between this arithmetic and machine arithmetic (other than the precision). All bignum arithmetic in Mathematica uses hidden digits, i.e., digits that are included in the calculations but are neither displayed nor considered significant. When loss of digits occurs due to cancelation of leading digits a positive value of MinPrecision may cause some of the hidden digits to become exposed as the least significant digits of the result. If there are already at least MinPrecision digits in the number no hidden digits will be exposed.

$\endgroup$
  • 3
    $\begingroup$ To give credit where due: The "Let me add..." with quotes from Jerry Keiper was added by another user. I only wish I had recalled that reference. $\endgroup$ – Daniel Lichtblau Feb 3 '17 at 16:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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