Best way to get a vector of BigFloat values drawn from a uniform distribution

General Usage
,
1

Is there some better and simple way to implement this for BigFloat:

function rand_uniform(lo::F, hi::F, n::Signed) where {F <: AbstractFloat}
  ret = Vector{F}(undef, n)
  dist = Distributions.Uniform(lo, hi)
  for i ∈ eachindex(ret)
    ret[i] = rand(dist)
  end
  ret
end

This approach fails because of a Distributions bug, it returns Vector{Float64} instead of Vector{BigFloat} even though lo and hi are BigFloat:

rand(Distributions.Uniform(lo, hi), n)

This, on the other hand, fails with a stack overflow because of recursion, for some reason:

rand!(zeros(F, n), Distributions.Uniform(lo, hi))
1 Like
2

If you only need uniformly distributed numbers, than it might be as simple as

rand_uniform(lo::T, hi::T, n) where {T <: AbstractFloat} = (hi - lo) * rand(T, n) .+ lo

Example result:

julia> rand_uniform(big(2.0), big(3.0), 10)
10-element Vector{BigFloat}:
 2.9174284555376025201145330738918913272402954502191116535401251098289736926065
 2.19631750020517507697116801083470708819543746756621169816997226270517883434895
 2.436503101481603637612779871925157875792727427507907758741273987186075483228166
 2.252714788297513525073447778551048335730249363760619306742365774731777906788439
 2.30163769835545213520563708849814348263120048796123561822508355898440209243495
 2.11121906581318486235373868765400026140501833929437774286345600289310872262491
 2.608123684200581891039323989360851860777985059179252353217119384669158234991661
 2.148349095950641373356503905512936681549144318191796575473727217069682413323919
 2.085010092947450646897480351115790914412328410397893567808432852935416907217017
 2.909521362406743263806281684702704618299470210785351171484376026730330458982956

Edit: The code assumes that rand(T, n) returns numbers between 0 and 1, so it only works for floating-point inputs.

1 Like
3

Why not just use the built-in rand function?

rand_uniform(lo::F, hi::F, n::Integer) where {F <: AbstractFloat} =
    [rand(F) * (hi - lo) + lo for _=1:n]

(I used a comprehension, unlike @barucden’s solution, to avoid allocating extra temporary arrays.)

2 Likes
4

Out of curiosity, I wanted to see how big of a difference this makes. I am observing the opposite effect though!

julia> using BenchmarkTools

julia> rand_uniform1(lo::F, hi::F, n::Integer) where {F <: AbstractFloat} =
           [rand(F) * (hi - lo) + lo for _=1:n]
rand_uniform1 (generic function with 1 method)

julia> rand_uniform2(lo::T, hi::T, n) where {T <: AbstractFloat} = (hi - lo) * rand(T, n) .+ lo
rand_uniform2 (generic function with 1 method)

julia> @btime rand_uniform1(big(1.0), big(10.0), 1000);
  464.910 μs (9005 allocations: 508.14 KiB)

julia> @btime rand_uniform2(big(1.0), big(10.0), 1000);
  337.539 μs (6010 allocations: 328.90 KiB)

julia> @btime rand_uniform1(big(1.0), big(10.0), 10000);
  4.763 ms (90006 allocations: 4.96 MiB)

julia> @btime rand_uniform2(big(1.0), big(10.0), 10000);
  3.481 ms (60013 allocations: 3.20 MiB)

The version with a temporary array is faster and allocates less. Why is that?

1 Like
5

Yeah, this is getting interesting. Both your and Steven’s function’s presented some opportunities for further optimization:

function rand_uniform_nsajko(lo::F, hi::F, n::Signed) where {F <: AbstractFloat}
  ret = Vector{F}(undef, n)
  dist = Distributions.Uniform(lo, hi)
  for i ∈ eachindex(ret)
    ret[i] = rand(dist)
  end
  ret
end

rand_uniform_barucden(lo::T, hi::T, n) where {T <: AbstractFloat} =
  (hi - lo) * rand(T, n) .+ lo

rand_uniform_barucden_opt(lo::T, hi::T, n) where {T <: AbstractFloat} =
  (hi - lo) .* rand(T, n) .+ lo

rand_uniform_stevengj(lo::F, hi::F, n::Integer) where {F <: AbstractFloat} =
  [rand(F) * (hi - lo) + lo for _=1:n]

rand_uniform_stevengj_opt(lo::F, hi::F, n::Integer) where {F <: AbstractFloat} =
  let d = hi - lo
    [rand(F) * d + lo for _=1:n]
  end

# Almost as fast as rand_uniform_barucden_opt
rand_uniform_map(lo::F, hi::F, n::Signed) where {F <: AbstractFloat} =
  let v = rand(F, n)
    map!(
      let d = hi - lo, l = lo
        x -> x*d + l
      end,
      v,
      v,
    )
  end

The optimized version of your function, rand_uniform_barucden_opt, seems to be the fastest of all. I guess this could be because your approach uses rand(::Type, ::Int), which could be specially optimized in the Julia stdlib.