Re: faster log factorials

[James Theiler <jt at lanl dot gov>]

On 12 May 2003, Gerard Jungman wrote:

] On Mon, 2003-05-12 at 08:02, James Theiler wrote:
]
]
] > A still-open question:  If we provide
] > pre-computed values, how many should we provide?  For the straight (no
] > logarithm) values, there is a natural cutoff at 170 since 170! (or is
] > it 171!?) is the largest value that is a valid IEEE double precison
] > number.  But for the logs, we can assume the higher the cutoff the
] > more often we'll be able to provide a fast precomputed value.  We
] > could easily provide thousands, and I think most computers nowadays
] > would not begrudge the memory.  But there may be other issues that I
] > am not considering.
]
] The natural cutoff is where Stirling's formula becomes good enough.

good suggestion.

] Computing that way will also be better than a memory access,
] I imagine.
]
] If we start using Stirling at 170 then it should be good
] enough to take up to sixth order in the series correction.
] So we can probably use the same arrays, etc.
]
] I don't remember how the logic works in the factorial functions.
] It might need a little tweaking to get optimal performance,
] so it takes the right branch without too much fuss.

For large arguments, seems to be gsl_sf_lnfact_e -> gsl_sf_lngamma_e
-> lngamma_lanczos which implements a Stirling-like formula.
However, Lanczos requires two log evals (cf Stirling requires one),
and Lanczos and Stirling both require several terms of a polynomial
(Lanczos has eight terms). Going straight to Stirling from the ln_fact
code could save a few branches, and one logarithm, and possibly a few
polynomial terms.  My intuition is that that won't beat a table
lookup, but it could be a substantial improvement. (And we have to
switch over from table lookup to large n formula at *some* point.) If
I get a chance, I'll give it a try.

jt

---------------------------------------------
James Theiler                     jt@lanl.gov
MS-B244, NIS-2, LANL        tel: 505/665-5682
Los Alamos, NM 87545        fax: 505/665-4414
----- Space and Remote Sensing Sciences -----