NIST publication 800-90B, Recommendations for the Entropy Sources Used for Random Bit Generation, contains an interesting continued fraction.
Define a function
where
is the incomplete gamma function.
NIST 800-90B gives [1] a continued fraction implement F, but the continued fraction includes a parameter k for no apparent reason. NIST cites a paper that in turn cites Abramowitz and Stegun formula 6.5.31.
It appears that the k in NIST 800-90B corresponds to a-1 in A&S 6.5.31, where a is the first argument to the incomplete gamma function and the exponent on 1/z, and so a = 3.
The continued fraction in A&S is
However, it’s not entirely clear how to fill in the dots. Apparently the denominators alternate between z and 1, but what about the numerators. The pattern we are give is
1, 1-a, 1, 2-a, 2
That’s as far as A&S goes. My initial guess was that we should ignore the 1 at the beginning of the series. In that case the rest of the series would be
3-a, 3, 4-a, 4, …
and that turns out to be right. If you don’t ignore the first 1, it can be hard to figure out how it fits into the pattern (because it doesn’t!).
How many terms of the continued fraction do we need? It depends on the size of the argument z. For large z, say z > 4, the terms shown above give a good approximation. But for smaller z, the approximation is terrible. For example, F(1) = 5, but the continued fraction evaluated at 1 gives -3/2.
The NIST paper requires evaluating F for values around 1, and it would take some work to figure out how far to carry the continued fraction before it is reliable in that range.
Related posts
[1] Why does NIST include a way to compute F? The incomplete gamma function is difficult to implement well and some mathematical libraries don’t include it. Having a stand-alone implementation of this function might make a testing program more self-contained.
This reminds me of something I noticed the other day about the continued fraction for e. It’s usually given as e = [2; 1, 2, 1, 1, 4, 1, 1, 6, …]. The first 2 seems not to follow the pattern. The number which does seem to follow the pattern is e-1 = [1; 1, 2, 1, 1, 4, 1, 1, 6, …].
But look what happens if we try starting the pattern at an earlier point! Going back one step gives us [0; 1, 1, 2, 1, 1, 4, 1, 1, 6, …], which must be 1/(e-1). So going back one step further gives us [1; 0, 1, 1, 2, 1, 1, 4, 1, 1, 6, …] which is defined as 1 + 1/[0; 1, 1, 2, 1, 1, 4, 1, 1, 6, …], which is 1 + 1/(1/(e-1)), which is 1 + (e-1) which is e!
So instead of using the ordinary continued fraction of e, which has a broken pattern, we could instead use the expression [1; 0, 1, 1, 2, 1, 1, 4, 1, 1, 6, …], in which the pattern is unbroken!
@”it’s not entirely clear how to fill in the dots”
Just looked up several continued fractions in A&S via the INDEX. Their dots fill pattern seems to be: Ignore the first token if necessary and proceed in the most obvious / simple way.