Search code examples
cprng

Getting a random real number in a certain range using WELL512


I'm using the WELL512 pseudorandom number generator function described in this paper. The function returns a random unsigned long value.

How do I use this return value to produce a random real number within a certain range - like a float between 340.92491 and 859812.53198 inclusive.

The documentation for the C rand() function seems to warn against using mod.


Solution

  • Well, mathematically it's just:

    min_value + (max_value - min_value) * (my_random() / (long double)ULONG_MAX)
    

    (Assuming my_random() returns a uniformly distributed number between 0 and ULONG_MAX)

    However, depending on the exact values of min_value, max_value, and ULONG_MAX, some floating point numbers will almost certainly be more likely than others.

    Each possible random unsigned long maps to a float by this formula. But since the number of distinct floating point numbers between min_value and max_value is almost certainly not exactly ULONG_MAX, some unsigned longs will map to the same floating point number or some floating point numbers will have no unsigned long map to them or both.

    Fixing this to make the result truly uniform is... non-trivial, I think. Maybe someone better read than I can cite a paper.

    [edit]

    Or see the answer to this question:

    Generating random floating-point values based on random bit stream

    That answer depends on the internals of the IEEE double representation. I am also not sure I fully understand how it works.

    [edit 2]

    OK now I understand how it works. The idea is to pick a random floating point representation between the min and the max, and then to throw it out with probability inversely proportional to its scale as represented by its exponent. Because for a uniform distribution, numbers between (say) 1/2 and 1 need to be half as likely as those between 1 and 2, but the number of floating point representations in those ranges is the same.

    I think you could make that code more efficient by first picking the exponent on a logarithmic scale -- say, by using ffs on a randomly-chosen integer -- and then picking the mantissa at random. Hm...