Why does this Gamma function approximation fail?

158 Views Asked by At

I've spent a lot of time poring over this equation from Wikipedia: Gamma approximation

I was curious about its accuracy for low $N$, so I plotted it on Desmos, against my reference curve which I described using the infinite product form of the Gamma function. Very strangely, for high $x$, or $T$ as I've named it in my Desmos plot, the plot grows more accurate... and then drops to the x axis, at height zero, and then promptly disappears from Desmos altogether, even though it's marked as if Desmos has finished calculating its curve. Even more strangely, I notice that for higher $N$, the threshold at which a "high" $T$ causes the plot to disappear actually lowers!

For example, at $N=225$, $T$ can equal $20$ and an accurate, visible curve is shown. For $T$ any much higher than $20$, it disappears! For $N=1000$, the threshold for $T$ is ... $2$!

I know that Desmos has finished computing all the sums and products since its curves are marked as ready, and I don't believe that there's much scope for numerical error here, so what's going wrong?

The second image shows what it looks like when $T$ is greater than this mysterious threshold of convergence.

Good plot: Desmos plot

Weird plot: Desmos plot with missing curve

1

There are 1 best solutions below

4
On BEST ANSWER

This is a machine precision problem.

In short, an IEEE 754 64-bit double has a maximum value of about $2 \cdot 10^{308}$. In your expression, for example, the single expression $T^N$ when $T = 23$ and $N = 225$ is about $2.44 \cdot 10^{306}$, but when $T = 24$ it exceeds the maximum value of a double.

I don't know the internals of Desmos, but this seems like pretty convincing evidence to me that they are using standard 64-bit double precision arithmetic, presumably in javascript, and so when $T = 24$ and $N = 225$, the term $T^N$ is probably returning infinity. Then I would suppose that Desmos voids the rest of the calculation since NaN and infinity are starting to appear.

This sort of problem happens all over the place with computer precision. The order in which operations occur can have big effects on the machine accuracy of the result.