CompiledFunction returns machine numbers smaller than $MinMachineNumber

Question

When thinking on the workaround for this LogLogPlot bug suggested by halirutan I noticed that CompiledFunction actually can return machine numbers smaller than $MinMachineNumber. Consider:

f = Compile[{{t, _Real}}, Exp[-9 t^2]];
num = f[8.872]
num < $MinMachineNumber
MachineNumberQ[num]
MachineNumberQ[num2 = ToExpression[ToString[num, StandardForm]]]

2.191864698767832`*^-308

True

True

False

As you see, num is smaller than $MinMachineNumber BUT MachineNumberQ returns True. Conversion of this number into String and backward returns arbitrary-precision number.

Moreover, CompiledFunction accepts such numbers and works with them correctly:

f2 = Compile[{{x, _Real}}, Sqrt@x];
f2[num]

1.480494747970364`*^-154

How this behavior could be explained? Does it mean that CompiledFunction actually CAN work with arbitrary-precision numbers?

P.S. I use Mathematica 8.0.4 under Windows 7 x64.

Answer 1

Congratulations! You find one of subnormal positive double :)

Another example

f = Compile[{{t, _Real}}, 2.0^t];

f[-1074]
f[-1075]

5.*10^-324
0

MachineNumberQ@f[-1074]

True

This doesn't mean that CompiledFunction can work with arbitrary-precision numbers.

Update

Normally Mathematica prevents such numbers

2.^-1074.
% // MachineNumberQ

4.940656458413*10^-324
False

But after a "hack" we can create subnormal positive double manually

SetSystemOptions["CatchMachineUnderflow" -> False];
2.^-1074.
% // MachineNumberQ

5.*10^-324
True

See also Improving speed of Exp[].

Note that Mathematica displays only the first digit of this number (a little bug). All digits:

RealDigits[%%]

{{4, 9, 4, 0, 6, 5, 6, 4, 5, 8, 4, 1, 2, 4, 6, 5}, -323}

Answer 2

In addition to @ybeltukov's excellent answer, I thought it would be worth noting the behaviour of RuntimeOptions, when compiling to either the Wolfram Virtual Machine (WVM) or to C, for these subnormal positive doubles.

f = Compile[{{t, _Real}}, Exp[-9 t^2]];

Needs["GeneralUtilities`"]
f[1.]    (* = 0.00012341 *)
f[8.872] (* = 2.191*10^-308 *)

Do[f[1.], {i, 10000}] // AccurateTiming
(* 0.0066 seconds *)

Do[f[8.872], {i, 10000}] // AccurateTiming
(* 0.0652 seconds *)

The difference in timings is an interesting one because of the subnormal result of f[8.872] - take a look at this plot to see.

range = {6., 7., 8., 8.871, 8.872, 9., 10., 11., 12.};
times = AccurateTiming[Do[f[#], {i, 10000}]] & /@ range;
ListPlot[Transpose[{range,times}]]

The default for "RuntimeOptions" is to not catch machine underflow. This can be turned on e.g.

f2 = Compile[{{t, _Real}}, Exp[-9 t^2],
                RuntimeOptions -> {"CatchMachineUnderflow" -> True}];
f2[1.]    (* all fine as before *)
f2[8.872] (* uh-oh! CompiledFunction::cfne *)

CompiledFunction::cfne: Numerical error encountered; proceeding with uncompiled evaluation.

This error has an associated performance hit, as it jumps out to uncompiled evaluation.

Do[f2[1.], {i, 10000}] // AccurateTiming
(* 0.0063 seconds *)

Do[f2[8.872], {i, 10000}] // AccurateTiming
(* 0.3324 seconds *)

Now look at the timing performance when the underflow is caught. Again, the sub-normal f2[8.872] is a tricky one to handle compared to subsequent results.