# How does InternalCompareNumeric work?

In this answer, Oleksandr R. mentioned an undocumented function InternalCompareNumeric and briefly explained its usage as follows:

InternalCompareNumeric[prec, a, b] returns -1, 0, or 1 according to whether a is less, equal, or greater than b when compared at the precision of a or b (whichever is less) minus prec decimal digits of "tolerance". It is the fundamental operation underlying Less, Equal, Greater, LessEqual etc. for finite-precision numeric types.

…To be honest, I failed to understand this description. In my eyes, it seems to suggest that, InternalCompareNumeric will compare

With[{minPre = Precision /@ {a, b} // Min}, N[{a, b}, minPre - prec]]


but it's not true, because

InternalCompareNumeric[1, 1.12, 1.22]
(* -1 *)


while the output should be 0 in my understanding.

Can someone explain the usage of InternalCompareNumeric in a more detailed way?

• my guess the comparison is actually made directly on the binary representation and the "decimal digits" count only approximately translates to a specific number of base 2 digits. Mar 20, 2017 at 22:03
• @george2079 Hmm… but InternalCompareNumeric[1, 112, 122] still gives -1. Mar 21, 2017 at 3:13
• FYI, the split point of -1 and 0 is around7.3456427360042 Mar 26, 2017 at 9:10
• I was discussing this with yode the other day. InternalCompareNumeric[n,a,b] basically answers the question, "how do a and b related ignoring n points of precision". So if you have two numbers and drop $MachinePrecision points of precision it's basically comparing them at the integer level. e.g. InternalCompareNumeric[IntegerPart@$MachinePrecision, 2.1, 2]  (because my $MachinePrecision isn't an int) gives 0, but InternalCompareNumeric[IntegerPart@$MachinePrecision - 1, 2.1, 2] gives 1. Mar 29, 2017 at 6:13
• Additional question: Why does InternalCompareNumeric give an error message when prec is set smaller than $MinPrecision? Mar 29, 2017 at 14:21 ## 4 Answers This is just to highlight another aspect of the other answer given. I would have expected a linear progression from the definition. I.e., if Precision[n]=p then we'd expect the switch point for CompareNumeric[x, n, n+epsilon] to be p. And this holds once you get past $MachinePrecision. Probably happens when Mathematica has to go from 64 bit numbers to whatever its internal infinite-precision numeric is. We can see this here:

pr[x_] :=
SelectFirst[
Range[1, 100, .1],
InternalCompareNumeric[#, SetPrecision[1.1, x], 1.2100] == 0 &
];

Show[
Plot[x, {x, 0, 25}, PlotStyle -> Gray],
{#, pr@#} & /@ Range[1, 25, .01] // ListPlot
]


I looked at the part for n < $MachinePrecision and thought: "oh, a log-type function (probably log2)". But no simple scaling on Log2 worked for me. Show[ListLinePlot@Table[Table[x*N@Log2[n], {n, 15}], {x, 5}], ListPlot[{#, pr@#} & /@ Range[1, 15, .01]]]  So maybe someone else can tackle what sort of pre-form that is. Once we have that we can probably piece back how Mathematica is interpreting things for its sub-machine precision numbers. ### Update As found by happy fish, $MachinePrecision - $MachinePrecision/x fits this really well: Plot[{ pr[x],$MachinePrecision - $MachinePrecision/x }, {x, 1,$MachinePrecision}
]


Noting that $MachinePrecision is 53Log10[2] by simple log rearrangements we have Log10[2^(53*(1 - 1/x))]. Since $MachinePrecision comes from how a 64 bit number is stored (there are 53 bits of precision), this could be saying it's using 53*(1 - 1/x) bits of precision for these numbers with Precision specified below $MachinePrecision. Why that would be, I have no idea. But it's a possibility. ### Update 2 Shadowray points out this is even odder when you have fewer points of precision available before they would be equivalent For example, make pr multivariate: pr[x_, y_: 1] := SelectFirst[ Range[1, 100, .1], With[{ a = SetPrecision[1. + (10^(-y)), x], b = SetPrecision[1. + (2*10^(-y)), 100] }, InternalCompareNumeric[#, a, b] == 0 & ] ];  Now plotting this for different amounts of reduction: Show[ Plot[x, {x, 0, 25}, PlotStyle -> Gray], Table[ {#, pr[#, n]} & /@ Range[1, 25, .05], {n, 17} ] // ListPlot ]  Note that the first 6 curves actually continue. They just jump up to y=x after $MachinePrecision (as I would naively have predicted).

If we use a near-continuous change:

Show[
Plot[x, {x, 0, 25}, PlotStyle -> Gray],
Table[
{#, pr[#, n]} & /@ Range[1, 25, .05],
{n, Range[1, 17, .1]}
] // ListPlot
]


The pre-discontinuity part can still be easily formulated:

Show[
Plot[x, {x, 0, 25}, PlotStyle -> Gray],
Plot[
Evaluate@Table[$MachinePrecision (1 - n/x), {n, 16}], {x, 1,$MachinePrecision}],
Table[
{#, pr[#, n]} & /@ Range[1, 25, .1],
{n, 17}
] // ListPlot
]


What's happening at the discontinuities is somewhat less clear. Note that the oddity at 17 is from the fact that 10^-17 is smaller than 10^-$MachinePrecision and so the system just treats the addition as doing nothing. The jumps need more explanation though. In any case, if we return to the storage argument, now we have 53 n/x bits being dropped, assuming the interpreting the pre-jump part of the curves as showing Log10[2^(53 (1 - n/x))] is valid. • @xzczd made the pr much more sensitive. Shows the real discontinuity at $MachinePrecision now. Mar 29, 2017 at 8:06
• I got similar plots. Although I used FindFormula to get very close fits, but those forms and coefficients do not make sense... Mar 29, 2017 at 8:13
• @happyfish It honestly looks like a -1/r thing (but shifted vertically). Simple (16+-1/r) is way too steep though. Mar 29, 2017 at 8:21
• $MachinePrecision -17.217894/x no idea Mar 29, 2017 at 8:27 • @happyfish Wow. That fits perfectly. Who knows if it's right. Maybe someone brighter than myself will see that and be able to explain why it works. Maybe it's supposed to be $MachinePrecision - $MachinePrecision / x? That's just as good a fit, I think. Mar 29, 2017 at 8:29 After some experiments I found some of the "rules" with which InternalCompareNumeric[c, a, b] operates. If you find any exceptions please notify me. Lets denote minpr = Min[Precision[a],Precision[b]]. If c is greater than pr[minpr] then InternalCompareNumeric[c, a, b] kind of throw away anything after the point. The function pr depends on a and b, the only symmetry I was able to detect is that if minpr>15, pr[minpr]=minpr (maybe connected to MachinePrecision?). The following figure is a possible pr[]. • What is pr[]? Mar 29, 2017 at 11:28 • @MichaelE2 An assumed function relationship between minpr and the threshold? Mar 30, 2017 at 9:25 • @xzczd There's no definition of pr[] in your OP, nor code for it. I notice MB1965 has a code for it (but no definition other than the code). Is it the same as used here? Since the ordering of answers changes over time (this was the top answer when I read it), it would be helpful to point to the definition of pr[] or give the code, if it's different. In the answer, "The following figure is a possible pr[]" suggests to me that there might be various ways to define pr[] and the plot is a possible one. Mar 30, 2017 at 11:03 • @MichaelE2 In my understanding, pr[] is more of a guess in this answer,. (BTW this is indeed the earlest answer. ) Anyway let's wait for Most Perfect Username's clarification :) Mar 30, 2017 at 11:10 • @xzczd A guess at what? There's no def. But yes, I think the OP can easily clarify. Mar 30, 2017 at 11:13 In this answer I will present an overall simulation of InternalCompareNumeric's behavior. Throughout this answer, there are four variables: InternalCompareNumeric[tol, a, b], first=Min[Abs[a],Abs[b]] diff=Abs[a-b] and prec is the minimum precision of the two specified. Here's my findings so far • For exact arithmetic, InternalCompareNumeric gives correct answer regardless of tol. • For machine precision integers (1., 2.), the split point is always $MachinePrecision.
• The order and sign of a, b does not matter, the minimum precision is used for calculation(this can be implied from my definition above).
• The smallest prec that makes the two numbers comparable is $-\log_{10}(\frac{diff}{first+diff})$ (this number comes from the definition of Precision), i.e. when $prec>-\log_{10}(\frac{diff}{first+diff})$, there exist a tol, such that the function returns -1 or 1; when prec is less than this number, the function returns 0 for all tol settings.
• When $prec\in [-\log_{10}(\frac{diff}{first+diff}),\$MachinePrecision]$, there exist a split point$y$, such that when$tol<y$, the function returns -1 or 1, when$tol>y$, the function returns 0. The analytical formula with first, diff, prec as parameters is $$y(first,diff,prec)=\MachinePrecision (1+\frac{1}{prec}\log_{10}(\frac{diff}{first+diff}))$$ • For prec>$MachinePrecision, the split point of tol is exactly the same as prec.
• For prec>$MachinePrecision, the formula above is incorrect. Some experiments show that when the ratio of first and diff is big enough, the comparison is converted to exact arithmetic(independent of tol). There should be two such points when diff/first is small and big, I haven't worked on it yet. A quick fixed of the formula is y[prec,first,diff]=prec-If[#>8,#,0]&[-Log10[diff/(first+diff)]], I will make careful research later when I have time. The following is the code version of the findings above: fishCompareNumeric[tol_, a_, b_] := With[{first = Min[Abs@a, Abs@b], diff = Abs[a - b], prec = Precision@{a, b}}, Module[{split}, Catch[ split = Piecewise[{ {Throw@Sign[a - b], prec === Infinity}, {$MachinePrecision, prec === MachinePrecision},
{Throw@0, prec < -Quiet@Log10[diff/(first + diff)]},
{$MachinePrecision (1 + Log10[diff/(first + diff)]/prec), -Log10[diff/( first + diff)] <= prec <=$MachinePrecision},
{prec, prec > $MachinePrecision} }]; Piecewise[{ {Sign[a - b], tol < split} }]]]]  Here is a plot of my binary search for approximate values, with parameters(OP's example) first=1.1, diff=0.1, horizontal axis being prec and vertical axis being split point (of tol). You may notice the two vertical lines at around 1.07 and 16, they are introduced by the sudden value changes, at the two points in my formula. Here is a plot of my exact formula with same parameters Code for binary search approximation: splitPointBinary[{{val1_, val2_}, {low_, up_}}] := With[{mid = InternalCompareNumeric[(low + up)/2, val1, val2]}, If[InternalCompareNumeric[low, val1, val2] == mid, {{val1, val2}, {(low + up)/2, up}}, {{val1, val2}, {low, (low + up)/2}}]] splitPointApprox[first_, diff_, prec_] := FixedPoint[ splitPointBinary, {{SetPrecision[first, prec], SetPrecision[first + diff, prec]}, {0.1, 1000}}][[2, 1]]  Code for my formula: splitPoint[first_, diff_, prec_] := Piecewise[{{$MachinePrecision (1 + Log10[diff/(first + diff)]/prec),
-Log10[diff/(first + diff)] < prec < $MachinePrecision}, {prec, prec >$MachinePrecision}}]


For OP's example,

splitPointApprox[1.1, 0.1, 2]
(*7.34564*)
splitPoint[1.1, 0.1, 2]
(*7.34564*)


Another one,

splitPointApprox[2.2, 0.5, 6]
(*14.0071*)
splitPoint[2.2, 0.5, 6]
(*14.0071*)

• This looks strange for splitPointApprox[1.1, 0.0000000000001, 2] Mar 29, 2017 at 11:14
• @xzczd answer updated. Mar 31, 2017 at 6:26
• I took the liberty to add a function interpreted from your findings to your answer, feel free to edit or roll back if you don't like it :) . Apr 2, 2017 at 14:03
• @xzczd I like it :) Apr 2, 2017 at 14:41
• Since I haven't think out a systematic enough way to judge the answer, I'll leave it unaccepted for a while, but this is undoubtfully the most plausible answer I received so far, +51 :) . Apr 3, 2017 at 3:21

Let's consider a slightly different approach than above. For InternalCompareNumeric[prec, a, a+ε], lets fix prec and a and find the threshold value for ε for which this switches from 0 to -1. For this purpose I define the following function

tf[p_?NumericQ, y_?NumericQ, z_?NumericQ] :=
InternalCompareNumeric[
p,
SetPrecision[1, y],
SetPrecision[1, y] + 10^-SetPrecision[z, 100]
]


That is we take a = SetPrecision[1, y] and ε = 10^-SetPrecision[p, 100]. (The SetPrecision in the definition of ε is there to avoid that the precision of z in the input can affect the outcome in anyway.)

We can then find the threshold value using simple bisection. For example abusing FindRoot:

tf2[p_?NumericQ, y_?NumericQ] := (z /.
FindRoot[tf[p, y, z] == -1/2, {z, -10, p + y + 1},
WorkingPrecision -> 100, MaxIterations -> 20])


Plotting the threshold value as a function of y for values of p we get

The top curve (p=0) is simply given by z = y. More generally the behavior of tf2 is modelled precisely by

Max[
Piecewise[
{
{y - p, y > MP},
{((MP - p)/MP) y, y < MP}
}
],
MP/2 UnitStep[y - Max[MP, p]]
]


where

MP=\$MachinePrecision


What does this tell us about the behavior of InternalCompareNumeric[prec,a,b]? Well, if the precision of a and b is set higher than machine precision and prec is not too large, the behavior is exactly as described by Oleksandr R.

It becomes a bit weirder when prec is large (compared to the precision of a and b), in which case there is an apparent upperbound of on the threshold difference of ε = 10^(-MP/2).

I don't quite understand the behavior of the threshold when the precision of a and b is smaller than machine precision. Maybe somebody else can shed some light.

Update:

Based on the above we can venture to propose a mock version of InternalCompareNumeric[prec, a, b]:

MockCompareNumeric[p_,a_,b_]:=With[
{ip = Min[Precision[a], Precision[b]]},
With[
{ep = Which[
ip === MachinePrecision, MP-p
ip >= Max[MP, p], Max[ip - p, MP/2],
True, Max[(1 - p/MP) ip, 0]
]
},
If[N[ Abs[a - b]/Max[Abs[a], Abs[b]]] < N[10^-ep], 0, Sign[a-b]]
]
];


I have tested the behavior of this mock function with the actual InternalCompareNumeric over a large range of prec, a b, and they match. (edit 3.4.2017: adjusted the mock function to get the correct behavior for MachinePrecision input.)

Ofcourse, the actual InternalCompareNumeric would not work like this. Instead it probably uses some low level functions working directly on the internal representation of the arbitrary precision numbers. (It certainly is a lot faster) However, the mock function gives us a fairly good idea of the behavior of InternalCompareNumeric[prec,a,b].
In short it determines some effective precision (ep in the mock function) from prec and the minimal precision of the inputs a and b and then determines whether a > b within a precision tolerance set by ep.
The difference between the behavior of ep dependent on whether the input precision is larger or smaller than machine precision probably arises due to differences in implementation leveraging the fact that for precision smaller than machine precision, Mathematica represents its arbitrary precision numbers simply as a pair of machine numbers representing the value and the precision of the number, in which case it can leverage various machine operations to do the comparison. I suspect that similarly the lower bounds on ep arise some checks that can be down efficiently using low level functions. The scaling of ep when ip < MP presumably is some design choice relating to the desired functionality, possibly related with wanting the effective precision to be larger than zero for reasonable values of prec.
• @ J.M. Ah, I see. Well, in this case we are not interested in finding the threshold where the two numbers are no longer considered the "same" by InternalCompareNumeric which is a slightly different question, than find the threshold where the internal representation becomes the same. (Arbitrary precision numbers in Mathematica can have many more digits in their internal representation than their actual precision value.) Mar 29, 2017 at 14:06
• What if one evaluates MockCompareNumeric[5, -0.92528723334935225612., 0.92896486862054583613.] ? Mar 30, 2017 at 9:36