rand(BFloat16, n) appears shifted from [0, 1) into (0, 1]

[JeffreySarnoff]

julia> using Random
julia> Random.seed!(1618); # TaskLocalRNG()

julia> f16s = extrema(sort!(rand(Float16, 10_000_000)))
(Float16(0.0), Float16(0.9995))

julia> Random.seed!(1618); # TaskLocalRNG()

julia> bf16s = extrema(sort!(rand(BFloat16, 10_000_000)))
(BFloat16(3.0733645f-8), BFloat16(1.0))

[milankl]

Is that because for BFloat16 the actual rand is done in Float32?

julia> f32s = extrema(sort!(rand(Float32, 10_000_000)))
(0.0f0, 0.99999994f0)

julia> extrema(sort!(BFloat16.(rand(Float32, 10_000_000))))
(BFloat16(0.0), BFloat16(1.0))

julia> extrema(sort!(BFloat16.(rand(Float32, 10_000_000))))
(BFloat16(1.1920929f-7), BFloat16(1.0))

julia> extrema(sort!(BFloat16.(rand(Float32, 10_000_000))))
(BFloat16(5.9604645f-8), BFloat16(1.0))

But that's not the full story because 100m times would also not result in a 0

julia> bf16s = extrema(sort!(rand(BFloat16, 100_000_000)))
(BFloat16(1.0535587f-8), BFloat16(1.0))

[christiangnrd]

Seems like the rand is done in Float64. If you change it to Float32, the domain appears to be [0,1] due to the package rounding to the nearest even when converting from Float32 to BFloat16.

julia> extrema(sort!(BFloat16.(rand(Float32, 10_000_000))))
(BFloat16(0.0), BFloat16(1.0))

[milankl]

We should define it as

using BFloat16s, Random
import Random: AbstractRNG, Sampler

function Random.rand(rng::AbstractRNG, ::Sampler{BFloat16})
    u = reinterpret(UInt16, one(BFloat16))
    u |= rand(rng, UInt16) >> 9                     # u in [1,2)
    return reinterpret(BFloat16, u) - one(BFloat16) # -1 for [0,1)
end

which is also how rand is done for Float16, Float32, Float64 I believe. There's also RandomNumbers.randfloat which is statistically better, but a little slower, see https://juliarandom.github.io/RandomNumbers.jl/dev/man/basics/#conversion-to-float

Using the Sampler immediately also implements rand(BFloat16, dims...) now we'd have

julia> extrema(sort!(rand(BFloat16, 10_000_000)))
(BFloat16(0.0), BFloat16(0.9921875))

🎉