Comment by art4711 | Hacker News Reader

Comment by art4711 | Hacker News Reader

Parent Full thread

I'm sorry, but if you mean Makoto Matsumoto who seems to be the name mentioned in the SFMT library, that library is bogus. It generates numbers in [1,2) and then maps that to [0,1) by subtracting 1. This can be trivially shown to lose information, in fact I added a little hack to my github repo two days ago showing why it's bogus when that specific library was mentioned to me (again). Unless you want to lose one bit of entropy in the first conversion every number gets. And it still doesn't get us anywhere close to arbitrary ranges (which is the hard problem).

I also just looked at stuff that Agner Fog has published, thank you for mentioning him. Same thing. Maybe my understanding of floating point is bad, but his is worse. One generator uses the classic rng()/((double)range_of_rng) approach which is the whole reason why the C++ standard generators aren't uniform (because it will only work you use a non-standard rounding mode for floating points and still doesn't get us anywhere close to arbitrary ranges which are a completely different problem). Another uses [1,2) - 1.0 which loses information. And the third one is "combines two different random number generators for improved randomness", which I guess sounds good, how does he do it? He generates two [0,1) numbers, adds them together then subtracts 1.0 when the result is bigger than 1.0 which both loses precision and compounds the rounding problems.

Literally the first thing I did when trying to figure this out months ago was to randomize the bits in the mantissa with a 0 exponent to get a [1,2) number then subtracted 1.0 from it, then it took me 10 minutes to test that this doesn't actually work.

So, maybe my understanding of floating point numbers is bad, but the authorities you cite are worse.

The fast generators as described by Matsumoto and Fog generate the numbers this way:

"If the 52-bit fraction part is uniformly randomly chosen, then the represented number is uniformly randomly distributed over [1, 2) with 52-bit precision."

Then for [0,1) 1 is subtracted.

I guess you can be confused because you don't understand the binary floating point? Specifically you can have wrong assumption about the numbers that are to be generated by the given algorithms. Let me explain you on the smaller example: let's generate the floats instead of doubles. If all possible input integer values were generated by Matsumoto's algorithm i e [0..8388607] and then from them generated the floats Xi in the range [0..1) then multiplying Xi with the constant 8388608f should get the results: 0f, 1f, 2f, ... 8388605f, 8388606f, 8388607f. And that's what is produced. Doubles have more bits but function the same: the constant which makes the tests easy is 4503599627370496.0 and the last number from the algorithm after the multiplication with that constant should be 4503599627370495, and the previous 4503599627370494 etc. The first three should of course be 0, 1 and 2.

What does your "better" algorithm produce under the given assumptions? Or do you actually have some other assumptions? If you do, it doesn't mean their algorithms are wrong, these do what they are designed to do and obviously correctly.

I suggest you you test your algorithm on floats, using absolutely all inputs and calculating all the outputs. Then you can write here what you achieved, that is, if you started from other goals, did you really reach them for all the inputs? Note that you don't need any random number generator to prove you reached your goal, just show that for given input of all ascending integers you get the desired output -- every result should also appear only once (I see your "tests" on github aren't made so).

But please don't write about C++ or use some notation from the latest C++ standards. Just write the inputs, outputs and describe your algorithms mathematically.

Got pointed back here that there was another response. Sure. If you subtract 1 from [1,2) you will get 52 bit precision (or 23 for floats as your magical constant indicates). This means that the number just before 1.0 will not be picked, then the number before that will be, the next one won't be, etc. You don't pick every other number in the range where you should have problems covering every number without introducing any irregularity to the result. So what? Map your random numbers to the set [0,1,2,3] divide by 4 and you also get a uniformly distributed number in [0,1).

With my method we'd pick the magical constant 16777216.0f and we'd also get 0, 1, 2, ... 16777216.0f is bigger than 8388608f, therefore it's better.

This [1,2)-1.0 method has 0 probability of picking the first number less than 1.0, 2^-52 probability to pick the next one, 0 probability to pick the next one, 2^-52 probability of picking the next one, etc. My algorithm has 2^-53 probability of picking the first number less than 1.0, 2^-53 probability of picking the next one, etc. Down to .5 I can pick every representable number. After that there is a choice. I chose to have the same probability of picking every smaller number (including 0.0) as the probability of picking numbers in the upper half of the range, so I have holes that I can't pick, but at least in the upper half of the range I can pick every number I can represent.

That's why I think that [1,2)-1.0 is flawed. Because it adds an arbitrary limitation to its output, without any good reasons except that it was convenient to implement and they didn't have to think hard.

> With my method we'd pick the magical constant 16777216.0f and we'd also get 0, 1, 2, ... 16777216.0f is bigger than 8388608f, therefore it's better.

The constants I've mentioned (23 bits of input and 8388608.0f) aren't "magic" but are the essence of floats when generating the uniform distribution for floats in the range 0..1.

That's what I want you to actually try in your code: to write it with the inputs of 24 bit integers 0..16777215 and to demonstrate that when you produce the float values and multiply these with the 16777216.0f you get every number from the 0..16777215.0f exactly once. Unless you show me that result (repeatable by anybody) I'll claim that your code doesn't work what it claims to have had achieved.

So once again, to be really clear: don't input any "random" bits. Input the integer numbers 0..16777215 (always only the lower 24 bits of them!) and demonstrate that your code works for floats (as explained above, every resulting float multiplied by the 16777216.0f must then be in the set of 0..16777215.0f and happen only once!). Only then I can believe you the equivalent solution works for doubles too. Before that, I don't believe you really understand what you want to do and what you're actually producing.

Reply on news.ycombinator.com