Open Zemyla opened 4 years ago
its a bug/problem in current random
a better algorithm with decent perf can be found here https://github.com/cartazio/old-random/blob/v1.3.0/src/Data/Distribution/FloatingInterval.hs
you raised a very good point @Zemyla :)
i've been in a rabbit hole on some ghc patching experiments for vector, i'll dust this stuff off post haste :)
@Zemyla could you please compare your algorithm to random-1.2.0
? I believe we generate more than 24/53 bits now.
This approach does not cover denormalized floats (@curiousleo correct me, if I'm wrong, please), which is discussed in #53, but AFAIU yours does not as well.
I don't claim to fully understand @Zemyla's code, but my impression is that it is based on ideas similar to http://allendowney.com/research/rand/ -- the countLeadingZeros
here can be used to read a uniform bitstring as a geometric distribution, which is how the exponent of a uniformly random IEEE float ought to be generated.
I would like to see something along those lines in a future version of random
, either as the default way to generate floating point numbers or as an optional "more accurate but slower" method, depending on how it performs.
However, code like this needs to be tested thoroughly. In addition, while generating numbers in the unit interval is nice, translating the unit interval into an arbitrary target interval via addition and multiplication leads to a loss in precision. So ideally, the uniform floating point generation method would be able to directly generate numbers in an arbitrary interval.
This is what https://gitlab.com/christoph-conrads/rademacher-fpl/ does. Note how thoroughly tested it is - which is necessary, because the code gets pretty complex.
@Zemyla if this is something you're interested in, you may also be interested in the discussions we had in the run-up to the 1.2.0 release: https://github.com/idontgetoutmuch/random/issues/113#issuecomment-624041080 (and the whole surrounding issue) as well as https://github.com/idontgetoutmuch/random/issues/105. I am experimenting with a Haskell implementation of a "uniform floating point number in an arbitrary range generator" here: https://github.com/curiousleo/random-float.
It seems like we could exploit the fact that those types have more precision closer to 0. The algorithm would look like this:
And then randomRFloating could be defined to take advantage of the increased precision near 0: