Skip to content

Commit

Permalink
re-add JuliaLang#8309 (which was reverted)
Browse files Browse the repository at this point in the history
This reverts commit 2abb59b.
This reverts commit 6117328.
  • Loading branch information
rfourquet committed Oct 10, 2014
1 parent 1d929d7 commit d52ebff
Show file tree
Hide file tree
Showing 3 changed files with 75 additions and 71 deletions.
108 changes: 47 additions & 61 deletions base/random.jl
Original file line number Diff line number Diff line change
Expand Up @@ -146,91 +146,77 @@ rand{T<:Number}(::Type{T}) = error("no random number generator for type $T; try
rand{T<:Number}(::Type{T}, dims::Int...) = rand(T, dims)


## Generate random integer within a range

# remainder function according to Knuth, where rem_knuth(a, 0) = a
rem_knuth(a::Uint, b::Uint) = a % (b + (b == 0)) + a * (b == 0)
rem_knuth{T<:Unsigned}(a::T, b::T) = b != 0 ? a % b : a

# maximum multiple of k <= 2^bits(T) decremented by one,
# that is 0xFFFFFFFF if k = typemax(T) - typemin(T) with intentional underflow
maxmultiple(k::Uint32) = itrunc(Uint32, div(0x0000000100000000,k + (k == 0))*k - 1)
maxmultiple(k::Uint64) = itrunc(Uint64, div(0x00000000000000010000000000000000, k + (k == 0))*k - 1)
# maximum multiple of k within 1:typemax(Uint128)
maxmultiple(k::Uint128) = div(typemax(Uint128), k + (k == 0))*k - 1
# maximum multiple of k within 1:2^32 or 1:2^64, depending on size
maxmultiplemix(k::Uint64) = itrunc(Uint64, div((k >> 32 != 0)*0x0000000000000000FFFFFFFF00000000 + 0x0000000100000000, k + (k == 0))*k - 1)

immutable RandIntGen{T<:Integer, U<:Unsigned}
a::T # first element of the range
k::U # range length or zero for full range
## Generate random integer within a range 1:n

# maxmultiple: precondition: k>0
# maximum multiple of k <= 2^numbits(k), decremented by one
maxmultiple(k::Uint32) = itrunc(Uint32, div(0x0000000100000000, k)*k - 1)
maxmultiple(k::Uint64) = itrunc(Uint64, div(0x00000000000000010000000000000000, k)*k - 1)
# maximum multiple of k <= typemax(Uint128), decremented by one
maxmultiple(k::Uint128) = div(typemax(Uint128), k)*k - 1
# 32 or 64 bits version depending on size
maxmultiplemix(k::Uint64) = k >> 32 == 0 ? uint64(maxmultiple(itrunc(Uint32, k))) : maxmultiple(k)

immutable RandIntGen{T<:Integer, U<:Union(Uint32, Uint64, Uint128)}
k::U # range length
u::U # rejection threshold
end
# generators with 32, 128 bits entropy
RandIntGen{T, U<:Union(Uint32, Uint128)}(a::T, k::U) = RandIntGen{T, U}(a, k, maxmultiple(k))
# mixed 32/64 bits entropy generator
RandIntGen{T}(a::T, k::Uint64) = RandIntGen{T,Uint64}(a, k, maxmultiplemix(k))


# generator for ranges
RandIntGen{T<:Unsigned}(r::UnitRange{T}) = isempty(r) ? error("range must be non-empty") : RandIntGen(first(r), last(r) - first(r) + one(T))
function RandIntGen(k::U)
k < 1 && error("No integers in the range [1, $k]. Did you try to pick an element from a empty collection?")
new(k, isa(k, Uint64) ? maxmultiplemix(k) : maxmultiple(k))
end
end

# generator API
# randintgen(k) returns a helper object for generating random integers in the range 1:k
randintgen{T<:Unsigned}(k::T) = RandIntGen{T, T}(k)
# specialized versions
for (T, U) in [(Uint8, Uint32), (Uint16, Uint32),
(Int8, Uint32), (Int16, Uint32), (Int32, Uint32), (Int64, Uint64), (Int128, Uint128),
(Bool, Uint32), (Char, Uint32)]

@eval RandIntGen(r::UnitRange{$T}) = isempty(r) ? error("range must be non-empty") : RandIntGen(first(r), convert($U, unsigned(last(r) - first(r)) + one($U))) # overflow ok
@eval randintgen(k::$T) = RandIntGen{$T,$U}(convert($U, k))
end

# this function uses 32 bit entropy for small ranges of length <= typemax(Uint32) + 1
# RandIntGen is responsible for providing the right value of k
function rand{T<:Union(Uint64, Int64)}(g::RandIntGen{T,Uint64})
local x::Uint64
if (g.k - 1) >> 32 == 0
x = rand(Uint32)
while x > g.u
x = rand(Uint32)
end
else
x = rand(Uint64)
while x > g.u
x = rand(Uint64)
end
@inline function rand_lessthan{U}(u::U)
while true
x = rand(U)
x <= u && return x
end
return reinterpret(T, reinterpret(Uint64, g.a) + rem_knuth(x, g.k))
end

function rand{T<:Integer, U<:Unsigned}(g::RandIntGen{T,U})
x = rand(U)
while x > g.u
x = rand(U)
end
itrunc(T, unsigned(g.a) + rem_knuth(x, g.k))
# this function uses 32 bit entropy for small ranges of length <= typemax(Uint32) + 1
# RandIntGen is responsible for providing the right value of k
function rand{T}(g::RandIntGen{T,Uint64})
x = (g.k - 1) >> 32 == 0 ?
uint64(rand_lessthan(itrunc(Uint32, g.u))) :
rand_lessthan(g.u)
return reinterpret(T, one(Uint64) + x % g.k)
end

rand{T<:Union(Signed,Unsigned,Bool,Char)}(r::UnitRange{T}) = rand(RandIntGen(r))
rand{T}(r::Range{T}) = r[rand(1:(length(r)))]
rand{T,U}(g::RandIntGen{T,U}) = itrunc(T, one(U) + rand_lessthan(g.u) % g.k)

function rand!(g::RandIntGen, A::AbstractArray)
for i = 1 : length(A)
@inbounds A[i] = rand(g)
end
return A
end

rand!{T<:Union(Signed,Unsigned,Bool,Char)}(r::UnitRange{T}, A::AbstractArray) = rand!(RandIntGen(r), A)
## Randomly draw a sample from an AbstractArray r
# (e.g. r is a range 0:2:8 or a vector [2, 3, 5, 7])

rand(r::AbstractArray) = r[rand(randintgen(length(r)))]

# Arrays of random integers

rand!(r::Range, A::AbstractArray) = rand!(r, A, ())

function rand!(r::Range, A::AbstractArray)
g = RandIntGen(1:(length(r)))
# TODO: this more general version is "disabled" until #8246 is resolved
function rand!(r::AbstractArray, A::AbstractArray, ::())
g = randintgen(length(r))
for i = 1 : length(A)
@inbounds A[i] = r[rand(g)]
end
return A
end

rand{T}(r::Range{T}, dims::Dims) = rand!(r, Array(T, dims))
rand(r::Range, dims::Int...) = rand(r, dims)
rand{T}(r::AbstractArray{T}, dims::Dims) = rand!(r, Array(T, dims), ())
rand(r::AbstractArray, dims::Int...) = rand(r, dims)


## random Bools
Expand Down
10 changes: 7 additions & 3 deletions doc/stdlib/base.rst
Original file line number Diff line number Diff line change
Expand Up @@ -3947,7 +3947,7 @@ Random number generation in Julia uses the `Mersenne Twister library <http://www

.. function:: rand!([rng], A)

Populate the array A with random number generated from the specified RNG.
Populate the array A with random numbers generated from the specified RNG.

.. function:: rand(rng::AbstractRNG, [dims...])

Expand All @@ -3961,9 +3961,13 @@ Random number generation in Julia uses the `Mersenne Twister library <http://www

Generate a random number or array of random numbes of the given type.

.. function:: rand(r, [dims...])
.. function:: rand(coll, [dims...])

Pick a random element or array of random elements from range ``r`` (for example, ``1:n`` or ``0:2:10``).
Pick a random element or array of random elements from the indexable collection ``coll`` (for example, ``1:n`` or ``['x','y','z']``).

.. function:: rand!(r, A)

Populate the array A with random values drawn uniformly from the range ``r``.

.. function:: randbool([dims...])

Expand Down
28 changes: 21 additions & 7 deletions test/random.jl
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,11 @@ rand!(MersenneTwister(0), A)
@test A == [0.8236475079774124 0.16456579813368521;
0.9103565379264364 0.17732884646626457]

@test rand(0:3:1000) in 0:3:1000
coll = {2, uint128(128), big(619), "string", 'c'}
@test rand(coll) in coll
@test issubset(rand(coll, 2, 3), coll)

# randn
@test randn(MersenneTwister(42)) == -0.5560268761463861
A = zeros(2, 2)
Expand All @@ -44,18 +49,27 @@ for T in (Int8, Uint8, Int16, Uint16, Int32, Uint32, Int64, Uint64, Int128, Uint
@test mod(r,2)==1

if T<:Integer && T!==Char
x = rand(typemin(T):typemax(T))
@test isa(x,T)
@test typemin(T) <= x <= typemax(T)
if T.size < Int.size
x = rand(typemin(T):typemax(T))
@test isa(x,T)
@test typemin(T) <= x <= typemax(T)
else # we bound the length of the range to typemax(T) so that there is no overflow
for (a, b) in ((typemin(T),typemin(T)+typemax(T)-one(T)),
(one(T),typemax(T)))
x = rand(a:b)
@test isa(x,T)
@test a <= x <= b
end
end
end
end

if sizeof(Int32) < sizeof(Int)
r = rand(int32(-1):typemax(Int32))
@test typeof(r) == Int32
@test -1 <= r <= typemax(Int32)
@test all([div(0x00010000000000000000,k)*k - 1 == Base.Random.RandIntGen(uint64(1:k)).u for k in 13 .+ int64(2).^(32:62)])
@test all([div(0x00010000000000000000,k)*k - 1 == Base.Random.RandIntGen(int64(1:k)).u for k in 13 .+ int64(2).^(32:61)])
@test all([div(0x00010000000000000000,k)*k - 1 == Base.Random.randintgen(uint64(k)).u for k in 13 .+ int64(2).^(32:62)])
@test all([div(0x00010000000000000000,k)*k - 1 == Base.Random.randintgen(int64(k)).u for k in 13 .+ int64(2).^(32:61)])

end

Expand Down Expand Up @@ -160,8 +174,8 @@ r = uint64(rand(uint32(97:122)))
srand(seed)
@test r == rand(uint64(97:122))

@test all([div(0x000100000000,k)*k - 1 == Base.Random.RandIntGen(uint64(1:k)).u for k in 13 .+ int64(2).^(1:30)])
@test all([div(0x000100000000,k)*k - 1 == Base.Random.RandIntGen(int64(1:k)).u for k in 13 .+ int64(2).^(1:30)])
@test all([div(0x000100000000,k)*k - 1 == Base.Random.randintgen(uint64(k)).u for k in 13 .+ int64(2).^(1:30)])
@test all([div(0x000100000000,k)*k - 1 == Base.Random.randintgen(int64(k)).u for k in 13 .+ int64(2).^(1:30)])

import Base.Random: uuid4, UUID

Expand Down

0 comments on commit d52ebff

Please sign in to comment.