Normal/Gaussian random number generation for D

[Joseph Rushton Wakeling]

On 10/23/2012 04:36 PM, jerro wrote:
> I have an implementation of the Ziggurat algorithm at
> https://github.com/jerro/phobos/blob/master/std/random.d#L2035.

OK, so I had a reasonable (not totally in-depth!) look through.  Looks good! 
And of course reminds of the other benefit of Ziggurat, that it can be used for 
multiple different random number distributions.

It looks like it should be readily possible to integrate the core Ziggurat 
functionality and to convert the normal() function in your code to be a 
NormalRandomNumberEngine struct -- so assuming you like my own architecture 
proposal, let's do this (I'm happy to hear suggestions for alternative designs, 
as I don't feel particularly confident in my API-designing abilities).

For the other distributions, my feeling is that in some cases there's a value in 
also having this "engine" approach, e.g. for exponentially-distributed numbers 
one could use Ziggurat or one could use the approach

     T u = uniform!("[)", T, T, UniformRandomNumberGenerator)(0.0, 1.0, urng);
     return -log(1 - u)/lambda;

... which is not as fast but has a much lower memory footprint.

> It modified the Ziggurat algorithm a bit, so that it doesn't need as many layers
> to work well, which reduces the memory consumption and makes initialization faster.

Can you expand on this, and maybe provide a reference?  I don't doubt your 
code's effectiveness but I think where RNGs are concerned we really need to be 
able to justify our algorithmic choices.  There's too much literature out there 
showing how commonly-used algorithms actually carry statistical flaws.

Bigger picture on my approach to non-uniform random number distributions.  The 
goal is to have the following:

     * Where useful, it should be possible to define and use multiple different
       internal "engines" for generating random numbers from the given
       distribution

     * For each distribution, there should be a function interface and a struct
       interface.

     * The struct implementation should store the distribution parameters and an
       instance of the internal engine (if any).

     * The function implementation should have 2 versions: one which allows the
       user to pass an engine of choice as input, one which contains a static
       instance of the specified engine (hence, thread-safe, distinguished
       according to both engine type and underlying uniform RNG type).

     * ... unless there's no call for distinct underlying engines, in which case
       the function version just takes parameters and uniform RNG :-)

     * The struct version should be useful to couple with an RNG instance to
       create an arbitrary random-number range à la Boost's variate_generator
       class.

So, a nice way to expand on our respective approaches might be to incorporate 
your Ziggurat, adapt it to my Normal engine API, and for me to write a similar 
setup for exponentially-distributed random numbers that uses the simple approach 
above and can also use your Ziggurat implementation.

What do you think?

Best wishes,

     -- Joe