math/rand/v2: a new API for math/rand and a first v2 for std

[rsc]

Based on earlier discussions in #26263 and #21835 as well as discussions with @robpike, I suggest adding a new version of math/rand, imported as math/rand/v2, to the standard library. This GitHub Discussion is meant to gather feedback before moving to a proposal.

The math/rand/v2 API would use math/rand as a starting point and then make the following backwards incompatible changes:

1. Remove Rand.Read and top-level Read. It is almost always a mistake to pretend that a pseudo-random generator is a good source of arbitrarily long byte sequences. The kinds of simulations and non-deterministic algorithms for which math/rand is a good choice almost never need byte sequences. Read is the only common piece of API between math/rand and crypto/rand, and code should essentially always use crypto/rand.Read instead. (math/rand.Read and crypto/rand.Read are problematic because they have the the same signature; math/rand.Int and crypto/rand.Int also both exist, but with different signatures, meaning code never accidentally mistakes one for the other.)

2. Remove Source.Seed, Rand.Seed, and top-level Seed. Top-level Seed is deprecated as of Go 1.20. Source.Seed and Rand.Seed assume that the underlying source can be seeded by a single int64, which is only true of a limited number of sources. Specific source implementations can provide Seed methods with appropriate signatures, or none at all for generators that cannot be reseeded; the details of seeding do not belong in the general interface.

   Note that removing top-level Seed means that the top-level functions like Int will always be randomly seeded rather than deterministically seeded. math/rand/v2 will not pay attention to the randautoseed GODEBUG setting that math/rand does; auto-seeding for top-level functions is the only mode. This means in turn that the specific PRNG algorithm used by the top-level functions is unspecified and can change from release to release without breaking any existing code.

3. Change the Source interface to have a single Uint64() uint64 method, replacing Int63() int64. The latter was overfitting to the original Mitchell & Reeds LFSR generator. Modern generators can provide uint64s, and uint64s have fewer special cases at call sites.

4. Remove Source64, which is unnecessary now that Source provides the Uint64 method.

5. Use a more straightforward implementation in Float32 and Float64. Taking Float64 as an example, it originally used float64(r.Int63()) / (1<<63), but this has the problem of occasionally rounding up to 1.0, which Float64 must not. We tried changing it to float64(r.Int63n(1<<53) / (1<<53), which avoids the rounding problem, but we decided it was a backwards compatibile change to break the value streams that rand.Rand generators, and instead added a retry loop for the rare 1.0 case. Now we can make the breaking change, which is simpler and faster.

   Note that some people have observed that the simple division does not make use of all the possible float64 values in the range [0, 1). For example values like 1/(1<<54), 1/(1<<55), 3/(1<<55), and so on are not generated at all, both in today's math/rand and in this simpler algorithm. Only the values 0 and 1/(1<<53) are, not the ones in between. It is possible to introduce even more complex algorithms to spread out the low values more while preserving something that can be deemed a uniform distribution, but these algorithms seem like overkill. The simple division should continue to suffice.

6. Implement Rand.Perm in terms of Rand.Shuffle. Shuffle is a bit more efficient, and then we have only one implementation. It has been suggested elsewhere to remove Rand.Perm instead, but keeping it costs little (especially if implemented in terms of Shuffle) and avoids unnecessary churn for users.

7. Rename Int31, Int31n, Int63, Int64n to Int32, Int32n, Int64, Int64n. The names are unnecessarily pedantic and confusing.

8. Add Uint32, Uint32n, Uint64, Uint64n, Uint, Uintn, both as top-level functions and methods on Rand. At the least, we need Uint64 to provide access to the widest values that the Source can provide. But we may as well complete the set while we are here.

9. Use Lemire's algorithm in Intn, Uintn, Int32n, Uint32n, Int64n, Uint64n. Preliminary benchmarks show a 40% savings compared to v1 Int31n and a 75% savings compared to v1 Int63n. (Like with Float64, since this changes the values returned, it is a breaking change that can only be applied in math/rand/v2.)

10. Add a new Source implementation, PCG-DXSM, with this API:

    func NewPCG(seed1, seed2 uint64) *PCG
    type PCG struct { ... }
    func (p *PCG) Uint64() uint64
    func (p *PCG) Seed(seed1, seed2 uint64)
    

    PCG is a simple, efficient algorithm with good statistical randomness properties. The DXSM variant was introduced by the author specifically to correct a rare, obscure shortcoming in the original (PCG-XSLRR) and is now the default generator in Numpy.

    PCG is provided for people who are writing simulations and need a seeded, deterministic source of randomness suitable for most purposes. Note that PCG is not assumed in any code. In particular the top-level functions need not use PCG (and in the current prototype do not). If in the future we add another algorithm, it will sit alongside PCG as a true peer.

11. Remove the Mitchell & Reeds LFSR generator and NewSource. An alternative to removing it would be to give it a name like NewLFSR or NewMitchellReeds, but that would serve no purpose. The obvious purpose would be to provide the original math/rand source so that code that needs to reproduce the exact math/rand outputs can do so. But the other changes to Rand, in routines like Intn and Float64, change the derived values from Rand's methods for a given Source input stream. Since preserving the math/rand Source stream would not suffice to preserve the math/rand Rand stream, preserving the math/rand Source is not worthwhile. Programs that need exactly the math/rand value streams can continue to use that package; it will not be removed.

The most direct motivation for this proposal is to clean up math/rand and fix these many lingering problems, especially the use of an outdated generator, slow algorithms, and the unfortunate collision with crypto/rand.Read.

A more indirect but equally important motivation is to set an example for other v2 packages in the standard library. Creating math/rand/v2 will let us work out tooling issues (support for v2 packages in gopls, goimports, and so on) in a relatively rarely used package with fairly low stakes attached before moving on to more commonly used, higher-stakes packages like sync/v2 or encoding/json/v2.

Once math/rand/v2 has shipped, we would tag and delete x/exp/rand. This would keep programs that already use x/exp/rand working (by reference to the tagged version or older versions of x/exp) but allow us to delete the code from the main branch of x/exp, making clear that development on it has ended.

A prototype of these math/rand/v2 changes can be found at https://go.dev/cl/502506 and the CLs below it in the stack.

[iand]

This could be an opportunity to rename the package to pseudorand which conveys very directly that these are not cryptographically secure random number generators.

[smyrman]

I am just going to jump in here, since in the end, names matters.

Also there is already an established name for this. If you see math/rand and math/rand/v2, it is clear what their relationship is and which one you should prefer for new code. If you see math/rand and math/pseudorand, not so much. (We're not deleting math/rand.)

Agree, the v2 makes the relationship obvious. However, based on precedent, packages can be deprecated, and users pointed to the new implementation by both docs and linters. So if there is a good enough naming suggestion, I would argue that a rename is possible, even if we don't get to try out the v2 suffix.

The argument for a rename, is aiming at helping new users not to pick the wrong package for their random use-case; for me personally, or any other experienced Go users, it doesn't matter.

I am not saying math/rand/v2 can't work, but if it was going to be something else, then math/arbitrary might work? It's not as short as rand, but shorter than pseudorand, and IMO a fairly accurate description of the package intent.

roll := 1 + arbitrary.IntN(5) // this game is clearly rigged

[earthboundkid]

“arbitrary” would be a good name for a third party package, where it’s good to have a name that is clever and a subtle reference to what it’s about, but for the standard library, the names of packages should be very boring and predictable.

[Merovius]

"arbitrary" is definitely incorrect for what the package does. return 42 is "arbitrary". But it would be wrong for a package implementing a PRNG.

FWIW in a parallel universe, we might chosen the names prng and csprng - they are inscrutable on first sight, but at least they are terse, accurate and distinct. But really, rand seems to be the name people would generally expect, distinguishing them by being crypto or not crypto seems good enough and as pointed out, their API will no longer overlap so you won't accidentally import the wrong one. I think the name is good as proposed.

[smyrman]

"arbitrary" is definitely incorrect for what the package does.

Alright, I suppose I shouldn't actually have suggested a name myself; I am not a native English speaker. My main point is really just that with proper deprecation of math/rand, then I don't really see a package rename as that problematic. Of course, if there is a better name available. But I suppose before that, the real question is wether the current name is problematic, and how well deprecation warnings in the docks works in terms of discouraging usage in new code.

Does anyone have any tools to gather numbers on this, or know if there have been some numbers collected in the past?

I see previous issues such as #11871 as an indication that there at least have been some concerns in the past regarding (security related) misuse of the math/rand package. Back then, the issue was solved by documentation.

[Merovius]

To be clear, if I would build a new PL from scratch today, I'd likely name its equivalent of math/rand also "rand". I'd be too concerned that people wouldn't find prng. So it's not really that a rename would be problematic, IMO. It's that I think it's an appropriate name for the concept.

And the structural fix for #11871 were #54880 and #56319, even if it took a while. And removing Seed and Read is part of the proposal, so the APIs of math/rand and crypto/rand become disjoint, so you do really have to confront that difference - you can't accidentally mix them up.

And beyond that, I really don't know what we could do. I don't think anything short of naming it insecure_rand would stop people from using it for crypto (and FWIW I'd got into a lot of arguments with people who used math/rand for crypto even after being told explicitly that it's a bad idea). And IMO if we'd name it insecure_rand, we'd also have to use insecure_json and insecure_os. Like, there is nothing that makes math/rand less secure than any other package outside the crypto/ tree.

I'm not saying there is not a great, descriptive and scrutable name. If someone wants to suggest something, great. But I wouldn't try to bend over backwards to avoid rand.

[bcmills]

Add Uint32, Uint32n, Uint64, Uint64n, Uint, Uintn, both as top-level functions and methods on Rand. At the least, we need Uint64 to provide access to the widest values that the Source can provide. But we may as well complete the set while we are here.

Now that we have generics, could we provide two or four top level-functions (maybe Int, Uint, IntN, UintN?) parameterized on the return type, instead of individually-sized variations?

Combined with better type inference, that could also be a significant noise reduction for common cases (such as choosing a random time.Duration).

[bcmills]

rand.Intn(s, 6) doesn't seem that bad to me.

(Note that the word rand need not be repeated between the package and function name — since this is math/rand/v2, nothing requires us to spend the good Intn name on a function with the bad “implicitly acquire a globally-shared source” semantics. 🙃)

I certainly find

	d := rand.Intn(s, 30*time.Second)

much clearer than

	d := time.Duration(r.Int64n(int64(30*time.Second)))

despite the transformation from a receiver to a free function. 😅

[zephyrtronium]

(Note that the word rand need not be repeated between the package and function name — since this is math/rand/v2, nothing requires us to spend the good Intn name on a function with the bad “implicitly acquire a globally-shared source” semantics. 🙃)

On the other hand, I think it's important for most users to have easy access to a source that doesn't need inputs to start using. This has traditionally been the global mutexed source. In Go 1.21, it will become runtime.fastrand64 by default. If we change the package-level functions to take a Source as input, then I'd like to also see a package-level func Default() Source, for example, which returns a Source that is ready to use.

[jimmyfrasche]

An adapter type could help for repeated calls though oneoffs would be a bit messy. Sketch below: Names are most likely wrong, just first thing I thought of. Function bodies elided.

type LessThan[T constraint] struct { r *Rand }
func NewLessThan[T constraint](*Rand) LessThan[T]
func (lt LessThan[T]) Draw() T
type Slice[T any] struct { r *Rand }
func NewSlice[T any](*Rand) Slice[T]
func (s Slice[T]) Shuffle(slice []T)
func (s Slice[T]) Choice(slice []T) T
func (s Slice[T]) Sample(n int, slice []T) []T

[neild]

All true, but r.Intn(6) is still a nicer API than rand.RandIntn(s, 6). If we were forced away from the latter because of concerns about generics, I would argue we are using generics wrong or generics are wrong.

I agree that r.Intn(6) is a nicer API.

I'm not sure what the best way to fit generic functions into the API is, but I'd also argue that if we don't have a good way to put generic functions in the new package then we are designing the package wrong or generics are wrong. Maybe generics are wrong and we need to figure out generic methods.

[rsc]

Generic methods would be nice but they are fundamentally incompatible with reflection and interfaces. Everything has limitations, including generics. We should not redesign generics just because it has limitations.

[raggi]

We had a code base that was affected by the performance of the prior API. We ended up with a sync.Pool wrapping instances that certainly could be replaced by a v2 implementation.

I did a casual study of the performance of the current API, and the exp/ API. https://gist.github.com/raggi/25035bf3fc5cd7c1a3b02dd3987a87fa

In our use case what is really desired is regularly doing this:

func veryRegularlyCalled(userGroup T) {
   seed := deriveSeed(userGroup)
   r := rand.New()
   r.Seed(seed)
   r.Intn(...)
}

That is, the cardinality of userGroup is very high, and the frequency with which this method is called is very high. The goal of the code is to project groups of objects into repeatable but otherwise random groups. There are many ways to do this, but a seeded random function is one potentially reasonable way.

I'm really curious if it is necessary, and/or what the motivation is to maintain the concept of a "Source" that is separate from the random number generator. This apparently flexibility in the API was unable to solve the problem at hand, and so it is not clear if it would solve for future problems either. The API mandates that the values escape to the heap, which adds a GC pressure cost to these kinds of use cases. Unless std itself provides alternative sources, this may be better handled by users defining an interface, and offering alternative implementations behind that interface. It should be possible to offer re-usable range-draw helpers such as Uint32 by way of something that depends on an interface instead, for example randv2.Uint32(r interface { Int63() int64 }) rather than baking them on an object with a strict internal representation.

I would love to see an API with PCG that gives me convenient access (offers a variety of drawn types as the Rand API does) while also being able to be constructed, seeded, and discarded entirely on the stack.

[raggi]

In the use case more than a single int is drawn. Perhaps my summary didn't do a good enough job, but the code in question is not public. I appreciate your explanation, but the benchmark you sent demonstrates the same thing already demonstrated in the gist I sent in my original post, which also demonstrated drawing more than a single element.

Unfortunately I think this is very distracting, I'd love to see more time put into the question of whether forcing the allocation in the API design is a useful and necessary element of the design. Perhaps I'm being unimaginative, but I don't really see it.

[Merovius]

If we would (as suggestion elsewhere in the thread) remove rand.Rand as a concept and instead have top-level functions that take a Source, it would be possible to add a Seed(seed1, seed2 uint64) method to PCG and then this wouldn't allocate, I believe:

s := new(rand.PCG)
s.Seed(42, 23) // or whatever
for i := 0; i < 100; i++ {
    rand.Intn(s, 100)
}

This would also mean you could put a PCG as a field in a struct without an extra indirection.

[zephyrtronium]

If I'm not mistaken, the allocation in the math/rand/v2 benchmark comes from the fact that the PCG is placed in an interface-typed field of Rand, which prevents devirtualization. Even without a Seed method, I would expect NewPCG not to allocate if its return value is only passed to functions.

[DeedleFake]

I'd like to see a way to decrease the allocations and other overhead, too. A while back I had a need to write an image.Image implementation that produced a random, but consistent color pattern. In other words, given an instance pat that was initialized with some seed, pat.At(x, y) should always return the same color.Color for the same values of x and y. I started out by simply reseeding the global math/rand functions but that was way too slow because of the mutex and rendered the concurrent implementation that I had mostly useless. The next thing I did was create new rand.Rand instances and seed them manually for every call to At(), but that also had a ton of overhead. I wound up directly using a custom copy of the PCG source in exp/rand that was designed not to allocate without using rand.Rand at all simply because there was no way to get rid of the overhead.

[rsc]

The API does not force an allocation at all. The current compiler happens not to split up the stack-allocated rand.Rand into variables, because if it did then it would devirtualize the source too. Compilers can improve, but APIs are forever. We should pick the right API. Especially since there is a clear path to this one allocation going away.

[icholy]

I'd argue in favor of using v2/math/rand instead of math/rand/v2. This approach aligns better with semantic import versioning if one views the stdlib as one large module with no prefix. To be clear, the v2 directory would not contain copies of all existing packages, only those which required a breaking change.

[mikeschinkel]

@rsc — Strictly from an IDE user perspective where the IDE auto-selects the prefix by right-most segment of import path, I would personally really appreciate it if the import path were /math/v2/rand so that I don’t end up either having to search and replace them after doing a routine refactoring, or worse forgetting to and the accidentally submitting a pull request with v2. aliases.

We already have to deal with so much of this when working with Kubernetes packages, it would suck to have to deal with it even more.

However, while I will understand if you consider that an IDE problem so not something for Go to consider, it doesn’t mean I won’t still be sad about it. 😢

P.S. I came here after reading about your proposal in GoLang Weekly to address this, so I am not bike-shedding merely after reading the OP’s comment.

[Merovius]

@mikeschinkel I strongly believe IDEs (and goimports) should handle this correctly. In particular, note that while the import path would be math/rand/v2, the package name would still be rand, so selectors will still be rand.Intn etc. And in my opinion, goimports (and any IDE doing their own thing) should default to choosing the last known major version, if they are looking for package rand. If they aren't currently, then we should fix that. Because that's a concern for third party modules and packages as well.

[mikeschinkel]

@Merovius — I agree with the "should" aspect of your argument.

However, as a pragmatist, I realize that "should" and "does" are often not the same and that real-world problems happen in the space between the two.

Also, it would not necessarily be wrong to say that rand is a subset of the 2nd version of math vs. math/rand has a 2nd version. Either of those could be logically correct.

[Merovius]

@mikeschinkel By the same logic it would be correct to say that math/rand is a subset of the 2nd version of the standard library, thus making v2/math/rand the logical choice.

It also seems to me that math/v2/rand doesn't solve the problem you are describing. It just moves it from importing math/rand to importing math.

[mikeschinkel]

@Merovius — Yes to your first point, but that is a distinction without a difference.

And no, you are incorrect about the 2nd point.

But this is not a hill I want to die on. I have made my point, Russ can consider it and decide. To further bikeshed this would do a disservice to everyone on the thread.

[zephyrtronium]

I've studied random number generators and their applications as a hobby for quite a few years now, including making several packages in various languages related to the topic, not to mention being a bother plenty often in the previous math/rand proposals. I'm glad to see some movement here.

I think the thing I miss most from the suggested API is implementing encoding.BinaryMarshaler and its dual. This is generally straightforward and cheap for PRNGs and lends a lot of versatility, since it makes reproducing states trivial. Plus, x/exp/rand already implements it for PCG.

Another thing I would like to see is a second standard library generator optimized for quality rather than throughput. I've suggested this on #26263 and been turned down, but I don't think I clearly communicated at the time just how difficult it is to diagnose the problems that arise due to insufficient quality. The only people who will even think about it are the experts who will probably choose a different package in the first place. Mersenne Twister is the popular choice in this space, but others exist. (I don't think this is especially likely to happen, still.)

That said, if we don't intend to add more generators, I'm not sure it's terribly informative to put the name of the particular PRNG algorithm in the name of the constructor. Why NewPCG rather than just replacing NewSource and mentioning that the generator is PCG in documentation? I missed the explanation of this in the discussion post.

[zephyrtronium]

I wrote down the API I'd like if Rand were removed and it turned out to be nearly identical to @carlmjohnson's #60751 (comment). The differences would be that I would not add SafeSource; leave Default unspecified beyond "a Source which is safe to use from multiple goroutines;" keep Zipf but with its method changed to Uint64(Source); and include Int, Int32, Int64, Uint, Uint32, and Uint64 as the discussion post suggests.

Details

package rand

// Source is a source of a pseudo-random number stream.
type Source interface {
	Uint64() uint64
}

// PCG implements [Source] with the PCG-DXSM algorithm.
type PCG struct {
	// contains filtered or unexported fields
}

// NewPCG creates a new PCG with the given seed.
func NewPCG(s1, s2 uint64) *PCG
// Uint64 produces a new 64-bit uniformly distributed variate.
func (*PCG) Uint64() uint64
// Seed sets the internal state of the generator.
func (*PCG) Seed(s1, s2 uint64)

// Default returns a Source which is safe to use from multiple goroutines.
// (Internally, it may be a singleton wrapper around runtime.fastrand64.)
func Default() Source

func Exponential(Source) float64
func Float32(Source) float32
func Float64(Source) float64
func Int(Source) int
func Int32(Source) int32
func Int64(Source) int64
func Intn[T ~int | ~uint | ...](Source, T) T
func Normal(Source) float64
func Perm(Source, int) []int
func Shuffle(Source, int, func(int, int))
func Uint(Source) uint
func Uint32(Source) uint32
func Uint64(Source) uint64

type Zipf struct {
	// contains filtered or unexported fields
}

func NewZipf(float64, float64, uint64) *Zipf
func (*Zipf) Uint64(Source) uint64

[rsc]

@zephyrtronium I explicitly addressed additional sources in the comment at the top. Planning for addition of other sources is why the function is NewPCG and not NewSource.

[rsc]

Happy to add BinaryMarshaler / BinaryUnmarshaler to PCG (but not Source of course).

[rsc]

Regarding removing Rand and just having operations on Source, I mentioned that part of the goal for math/rand/v2 is to set the example for v2 generally. One example I think is critical to set is not to change APIs without a very good reason. Every detail we change just because we think it's a little bit nicer is a detail that anyone switching from the old package to the new package has to relearn, something that invalidates old blog posts, documentation, books, and so on.

We have experience with the existing API, and deviations should be well-motivated. This is why there are explanations on each of the changes above. Honestly, I am on the fence about renaming Int63n to Int64n. That might be a mistake. Designing a completely new API just because we think it might be better is definitely a mistake. It very well might be worse too, and we have no experience with it. The second-system effect is real, and we need to avoid it.

[kortschak]

The serialisation/deserialisation pair of methods was something that we really wanted in Gonum for being able to use PRNGs for simulations that can be reproduced, stopped and restarted.

[Merovius]

Remove Rand.Read and top-level Read. It is almost always a mistake to pretend that a pseudo-random generator is a good source of arbitrarily long byte sequences. The kinds of simulations and non-deterministic algorithms for which math/rand is a good choice almost never need byte sequences.

FWIW my one use case for this is that I tend to use Read for unique identifiers, by hex/base64/… encoding a 128 bit string or something like that. crypto/rand feels slow for that, but maybe that's misunderstood optimization. Just wanted to bring it up.

I fully agree that removing the top-level Read is beneficial to remove the ambiguity with crypto/rand, but I don't think that's a concern for Rand.Read, so we could decide to keep the latter?

[Merovius]

TBC I have zero passion to keep Read around - my use case works perfectly fine by doing fmt.Sprintf("%.16x%.16x", rand.Uint64(), rand.Uint64()) or something like that. I'm mainly curious to learn if it's that bad of an idea that it justifies removing Read.

[danderson]

FWIW, some order-of-magnitude numbers: modern cryptographically secure PRNGs produce secure random data at tens of gigabits/sec. If you don't need a reproducible (seedable) pseudorandom sequence, the system CSPRNG is almost always the correct source of random bits to use. And by using a CSPRNG, you will never have to ask yourself "does this UUID cause some obscure security issue if an attacker can steer the seeding of my RNG?".

IMO there is a lot of leftover impression in peoples's minds that CSPRNGs are slow/expensive due to how /dev/random behaved years ago on linux: it used algorithms that weren't very fast, and also deliberately slowed itself down even more due to the entropy estimation logic. At this point, the biggest bottleneck for RNG on modern systems is the syscall overhead, which is being worked on but can also be mitigated with larger reads + buffering in userspace if necessary (it's almost certainly not necessary).

This is obviously hyperbolic, but the Go standard library would do a great service to its users by making it as hard as possible to use any RNG other than crypto/rand. You should have to submit evidence that it's too slow or that you specifically need deterministic behavior before you're allowed to instantiate the non-cryptographically secure versions :P

[earthboundkid]

I think it would be convenient to have access to a goroutine safe Source other than the top level functions.

One way to do that would be to have a function that takes a Source and returns a safe Source that's wrapped in a mutex.

I also don't think Rand really pays for itself. Once you have all Source64s, there's not as much need to stash values in Rand.

So, putting it together, I would like to see something like this:

package rand

type Source {
  Uint64() uint64
}

func NewPCG(seed1, seed2 uint64) *PCG
type PCG struct { ... }
func (p *PCG) Uint64() uint64
func (p *PCG) Seed(seed1, seed2 uint64)

func SafeSource(Source) Source

func Default() Source { return default }

func init() { default = SafeSource(NewPCG(someseed(), someseed2()) }

func Float32(Source) float32
func Float64(Source) float64
func ExpFloat64(Source) float64
func NormFloat64(Source) float64

func Intn[Integer ...](Source, Integer) Integer

func Choice[T any](Source, []T) T

func Perm(Source, int) []int
func Shuffle[T any](Source, []T)
func ShuffleFunc(Source, int, func(int, int))

That's pretty tight and easy to understand.

[Merovius]

I think the only reason we currently need a magic locked source is Read. Without Read, a locked source can be provided without any issues by a third party library.

[zephyrtronium]

While I agree that Rand doesn't add much (I was intending to suggest a similar change of API), the correct way to use a PRNG in multiple goroutines is to use distinct PRNGs, not to share memory with a contended mutex. The current math/rand generator makes that hard because it has a very large state, but PCG's state is 16 bytes.

[earthboundkid]

Another use for the forbidden fruit of Go routine local values. :-)

[Merovius]

FWIW you can use new(maphash.Hash).Sum64() to get a seed from TLS (effectively this calls runtime.fastrand) and then use that to seed your artisanally crafted 3rd party local Source. My expectation is that we get a Default() Source that does effectively this. And if you want to use a shared custom Source you have to either implement locking or this scheme yourself.

I don't think it makes sense to provide a general mechanism for this. Though maybe, at some point, we get per-P sharding, which seems like a good way to address this.

[rsc]

I think it would be convenient to have access to a goroutine safe Source other than the top level functions.

type SourceFunc func() uint64
func (f SourceFunc) Uint64() uint64 { return f() }

var goroutineSafeSource rand.Source = SourceFunc(rand.Uint64)

We may also expose the source directly with a function named something like Default, although I would like to see the use cases a bit more clearly first.

[flyingmutant]

I've written an alternative rand package (pgregory.net/rand) that tries to improve upon speed and quality of math/rand and x/exp/rand. The most controversial change I've made is probably the removal of Source to make inlining possible, but that may be too radical for math/rand/v2. Most of the changes proposed here are great, but I would like to propose two additional ones:

• Use the v2 compatibility break to fix the bias in NormFloat64 and ExpFloat64 generators. I believe they were incorrectly ported from the original ziggurat algorithm that operates on float32 to float64: in x := float64(j) * float64(we[i]) the multipliers are not independent, which leads to detectable bias in generated numbers.

• Consider using 32.64 and 64.128 fixed-point multiplication instead of Lemire's algorithm in Intn, Uintn, Int32n, Uint32n, Int64n, Uint64n. These algorithms avoid both division and rejection sampling completely and guarantee that result is unbiased with either 1 - 2^-32 or 1 - 2^-64 probability (requiring at least 2^64 or 2^128 samples to detect the bias).

[rsc]

@flyingmutant, I looked at the 32.64 and 64.128 links and honestly it all seems a bit hand-wavy. Lemire's algorithm as implemented in my rand/v2 prototype looks like:

	hi, lo := bits.Mul64(r.Uint64(), n)
	if lo < n {
		thresh := (-n) % n
		for lo < thresh {
			hi, lo = bits.Mul64(r.Uint64(), n)
		}
	}
	return hi

The lo < n test is essentially always false, meaning it will predict very well, and so the overwhelmingly typical execution is a single Mul64 and a well-predicted branch. And it's less work than https://github.com/flyingmutant/rand/blob/d64f7d2a06fcd674bf70702ccefa12e9634094f0/rand.go#L274-L286. So it seems like a win-win: faster and unbiased. What am I missing?

Happy to look into breaking changes for the ExpFloat64 and NormFloat64 algorithms to use more efficient ones. Those did not get the same attention as the other methods in this prototype. Thanks for the links.

[flyingmutant]

Lemire's algorithm is definitely more of a tradeoff than a win-win. It guarantees that every result is unbiased, but leaves a lot of performance on the table because of unpredictability and easily triggered worst case behavior.

For Uint32n, there is no contest:

cpu: Intel(R) Core(TM) i7-10510U CPU @ 1.80GHz
BenchmarkUint32n/fp/small-8         	803663504	        1.428 ns/op
BenchmarkUint32n/fp/midrange-8      	825222320	        1.416 ns/op
BenchmarkUint32n/fp/big-8           	843441136	        1.360 ns/op
BenchmarkUint32n/fp/rand-8          	728327844	        1.485 ns/op
BenchmarkUint32n/lemire/small-8     	487817625	        3.202 ns/op
BenchmarkUint32n/lemire/midrange-8  	488942898	        3.059 ns/op
BenchmarkUint32n/lemire/big-8       	155839828	        7.826 ns/op
BenchmarkUint32n/lemire/rand-8      	100000000	       12.79 ns/op

You can't beat a single multiply

func Uint32n(n uint32) uint32 {
	res, _ := bits.Mul64(Uint64(), uint64(n))
	return uint32(res)
}

For Uint64n, there are cases where Lemire's will only perform a single multiplication while fixed-point version will perform two, but I'd say that in general fixed-point version is still much faster:

cpu: Intel(R) Core(TM) i7-10510U CPU @ 1.80GHz
BenchmarkUint64n/fp/small-8         	578032960	        2.059 ns/op
BenchmarkUint64n/fp/midrange-8      	342683163	        3.081 ns/op
BenchmarkUint64n/fp/big-8           	341557892	        3.292 ns/op
BenchmarkUint64n/fp/rand-8          	328838946	        3.631 ns/op
BenchmarkUint64n/lemire/small-8     	562186795	        2.173 ns/op
BenchmarkUint64n/lemire/midrange-8  	499154936	        2.226 ns/op
BenchmarkUint64n/lemire/big-8       	100000000	       10.50 ns/op
BenchmarkUint64n/lemire/rand-8      	74914030	       14.58 ns/op

The downside of the algorithms I am proposing is that on average you will get 1 biased result every 2^32 values for 32-bit n and 1 biased result every 2^64 values for 64-bit n. I consider this amount of bias to be non-detectable in practice, but that is of course debatable.

Benchmark code is here.

[randall77]

@flyingmutant Your benchmark code suffers from the fact that the results of the computation are unused, so some of the code you're trying to measure inlines and completely compiles away.

[rsc]

Indeed. See https://go-review.googlesource.com/c/go/+/502496/1/src/math/rand/v2/rand_test.go for what should be more accurate benchmark functions.

On top of that, if you're going to assume 64x64->128-bit multiplications are available and efficient on 32-bit systems, then a fair comparison of the Lemire Uint32n would use the 64x64, like (untested):

func Uint32nLemire(n uint32) uint32 {
	n64 := uint64(n)
	hi, lo := bits.Mul64(Uint64(), n64)
	if lo < n64 {
		thresh := (-n64) % n64
		for lo < thresh {
			hi, lo = bits.Mul64(Uint64(), n64)
		}
	}
	return uint32(hi)
}

That branch is almost never going to be true, even for large uint32 values of n. The code should run at basically the same speed as the biased one, just not biased.

[flyingmutant]

@randall77 I am well aware of the fact that result in theory can be compiled away, and try to work around it in my package benchmarks. I've left it out of the code above because of the reasoning that changing global wyrandState is enough to make sure that the code runs for real when benchmarking. Doing the sum + sink trick does not have any meaningful effect on results:

cpu: Intel(R) Core(TM) i7-10510U CPU @ 1.80GHz
BenchmarkUint64n/fp/small-8         	502715562	        2.113 ns/op
BenchmarkUint64n/fp/midrange-8      	311604848	        3.467 ns/op
BenchmarkUint64n/fp/big-8           	366345884	        3.217 ns/op
BenchmarkUint64n/fp/rand-8          	250678396	        4.527 ns/op
BenchmarkUint64n/lemire/small-8     	474199490	        2.573 ns/op
BenchmarkUint64n/lemire/midrange-8  	479039515	        2.245 ns/op
BenchmarkUint64n/lemire/big-8       	100000000	       11.04 ns/op
BenchmarkUint64n/lemire/rand-8      	73511844	       16.09 ns/op

@rsc I think that doing 64x64->128 multiplication (with or without rejection sampling) for Uint32n is definitely the way to go. For Uint64n, Lemire's algorithm does have an impossible to predict slow path that you can easily hit, I don't think there is a way around it.

Edit: Relying on wyrandState alone is definitely wrong. I've updated the gist to use sum + sink as well.

[leaxoy]

It is recommended to show the practice of this situation in other languages ​​for comparison.

[rsc]

For a new package, that would make sense. For a v2 of an existing package, the approach is (1) start with v1, and (2) make well-motivated, necessary changes only.

[earthboundkid]

https://fuglede.dk/en/blog/bias-in-net-rng/ is an interesting case. I don’t know what the follow up was though.

[zephyrtronium]

@carlmjohnson I think that article is a good contrast to illustrate why the proposed approach is correct.

The article finds a set of polynomials describing the problem scenario and says, "I wish I had a nice analytical proof of this relation, but it came about only through experimentation and searching for patterns." I have neither the experience nor the time to derive a full proof, but I am confident that the polynomials are derived directly from the multiplier that reduces the integer variate to floating point [0, 1). In particular, the degree of the polynomials is exactly the index of the only explicit nonzero bit in the mantissa, relative to the implicit 1. And, thinking about it in terms of "grade school" multiplication, it follows that that 22 bit has a visible effect on the multiplication result.

Moreover, in binary floating point arithmetic, multiplication by a (possibly negative) integer power of 2 is exact if and only if it doesn't overflow to infinity or underflow the normals. When the result is a normal number, it is equivalent to summing the exponent components of the numbers. So, generating integers uniform over [0, 2^53-1] (where 53 is the number of bits in a float64 mantissa) and dividing by 2^53, as suggested in the discussion post, cannot produce the problem described in the article.

[bagasme]

5. Use a more straightforward implementation in Float32 and Float64. Taking Float64 as an example, it originally used `float64(r.Int63()) / (1<<63)`, but this has the problem of occasionally rounding up to 1.0, which Float64 must not. We [tried changing it](https://codereview.appspot.com/22730043/diff/60001/src/pkg/math/rand/rand.go) to `float64(r.Int63n(1<<53) / (1<<53)`, which avoids the rounding problem, but we decided it was a backwards compatibile change to break the value streams that rand.Rand generators, and instead [added a retry loop](https://codereview.appspot.com/95460049/diff/50001/src/pkg/math/rand/rand.go) for the rare 1.0 case. Now we can make the breaking change, which is simpler and faster.
   Note that some people have observed that the simple division does not make use of all the possible float64 values in the range [0, 1). For example values like 1/(1<<54), 1/(1<<55), 3/(1<<55), and so on are not generated at all, both in today's math/rand and in this simpler algorithm. Only the values 0 and 1/(1<<53) are, not the ones in between. It is possible to introduce even more complex algorithms to spread out the low values more while preserving something that can be deemed a uniform distribution, but these algorithms seem like overkill. The simple division should continue to suffice.

For me, for completeness, there needs to be a way to generate these omitted values when simple division algorithm doesn't generate
them, but still preserving uniform distribution of the whole values.

[rsc]

That's easily provided by an importable third-party package if truly needed. I am not convinced it's wise in the standard package, for two reasons:

1. The simple division produces numbers that have uniform gaps between them. That is, the gaps near 1 are the same size as the gaps near zero. That specific uniformity is lost if you get clever about creating smaller gaps near zero.
2. Many times the result of Float32 or Float64 is scaled or translated or both (Float64()*X + Y). The +Y undoes all the work to create those gaps, and the *X scales the gaps so that in the new range it's still not true that every possible floating point value will be generated.

It does not seem worth the cost.

[zephyrtronium]

More strongly, generating the "omitted" values of float64 implies that the distribution is no longer numerically uniform when the generator is uniform. Rather, the density over same-size intervals increases exponentially toward zero.

The problem can be illustrated with a smaller domain. Say we have a hypothetical float8 with 1-bit sign, 4-bit exponent, and 3+1-bit mantissa. A positive finite float8 number is represented as a*2^e + b*2^(e-1) + c*2^(e-2) + d*2^(e-3) where ^ means exponentiation, e is an integer in -7...6, each of abcd are either 1 or 0, and a is 0 if and only if e is -7. Generating 4-bit integers and dividing by 2^4 can produce 0, 1/16, 1/8, 3/16, ..., 15/16, all with equal probability (assuming the integers are uniform). That means of float8's representable values, we're missing:

• all seven subnormals (1..7/512) and all eight normals for each exponent of -6 and -5, for a total of 23 numbers strictly between 0 and 1/16;
• 9/128, 5/64, 11/128, 3/32, 13/128, 7/64, 15/128 between 1/16 and 1/8;
• 9/64, 5/32, 11/64 between 1/8 and 3/16;
• 13/64, 7/32, 15/64 between 3/16 and 1/4;
• 9/32 between 1/4 and 5/16;
• 11/32 between 5/16 and 3/8;
• 13/32 between 3/8 and 7/16;
• 15/32 between 7/16 and 1/2. (All values between 1/2 and 1 can be produced by this procedure.)

Clearly these numbers are not distributed uniformly. In quartiles of all representable numbers, the third and fourth are equal sizes at 4 values each, but the second quartile has 8, and the first has 40.

A pattern I use commonly goes like

if rng.Float64() < p {
    // Do work with probability p.
}

If we generate every possible float64 value, that pattern performs much more work than probability p would imply, if p is not 0 or 1.

We could hypothetically use a non-uniform input to restore numerical uniformity to the output by selecting exponent and mantissa separately. To make the distribution uniform for our float8, we need 7 bits of entropy for the exponent, plus another 3 for the mantissa. For float64, we need 1023 bits for the exponent and 53 for the mantissa. PCG DXSM only provides 128 bits! So, your small values will be biased, non-uniform, and full of "gaps" anyway – but in a fashion that requires much fancier math specific to the PRNG algorithm to predict, rather than arithmetic – unless you use a different and far more expensive generator specifically for this purpose. (This is the quality issue I've mentioned elsewhere.)

[rsc]

@zephyrtronium For what it's worth, I don't think anyone is suggesting that Float64() return any of the values in the range [0,1) with equal probability. Doing so would, as you point out, skew too much toward smaller values. Instead the suggestion is to take some of the hits toward the low end and redistribute them to nearby floats, so that Float64() < p still has probability p but there are more possible results. I still don't believe that's worthwhile, but that's the suggestion.

[zephyrtronium]

I think that's a reasonable assumption, and I addressed why that is also not just "not worthwhile" but actually infeasible in my last paragraph. That's probably the more compelling argument, so maybe it should not have been a fifth the size of the other. 🙂

[szuecs]

Just for completeness as I think there's a lot of great input.

My normal use case for seeding and using math/rand is testing functions that are probability based.
If I don't get always the same numbers I will never know if I was unlucky with random values or something really broke.
As far as I understand the PCG would address it.

[rsc]

Thanks for the discussion everyone. I have filed a proposal: #61716.