-----Original Message----- From: boost-bounces@lists.boost.org [mailto:boost-bounces@lists.boost.org] On Behalf Of Joel Eidsath Sent: 22 December 2007 23:02 To: boost@lists.boost.org Subject: [boost] [random] poisson_distribution problem and fix
The poisson_distribution code will fail for any mean larger than 750 or so (on my x86_64 machine). The problem is a comparison with exp(-_mean) which evaluates as zero right around there.
A far better fix would be to use an algorithm other than Knuth's for this. Knuth's algorithm is simple, but O(n). Here is some code I wrote, converting some old Fortran code that a friend had got from Numerical Recipes once upon a time. It works, but would need actual testing before use in production code:
#ifndef F_POISSON_HEADER #define F_POISSON_HEADER
<snipped>
#endif
This looks too obviously lifted from NR in C (or FORTRAN) to fit the Boost license. But the algorithm should still be usable? (But I am not qualified to judge). ln gamma (and poisson) is available from the Math Toolkit, but this is a random requirement and you may value speed more than accuracy? Paul --- Paul A Bristow Prizet Farmhouse, Kendal, Cumbria UK LA8 8AB +44 1539561830 & SMS, Mobile +44 7714 330204 & SMS pbristow@hetp.u-net.com