Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Make rand work with AbstractArray instead of only with Range #8309

Merged
merged 9 commits into from
Oct 1, 2014
Merged

Make rand work with AbstractArray instead of only with Range #8309

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
ea2bd4b
restrict RandIntGen to generate random numbers in the interval [1,n]
rfourquet Sep 10, 2014
e2b14a1
remove redundant rand! method
rfourquet Sep 10, 2014
6d329ce
make rand work with AbstractArray instead of only with Range
rfourquet Sep 10, 2014
9d0282e
update the docs for rand/rand! and add some related tests
rfourquet Sep 11, 2014
57bd2b7
make RandIntGen stop handling "full range" random generation
rfourquet Sep 14, 2014
91aa94f
factor out common code in rand(::RandIntGen) methods
rfourquet Sep 14, 2014
d9814ff
remove rand!(::AbstractArray, ::AbstractArray) until #8246 is fixed
rfourquet Sep 16, 2014
38bef62
use only non-overflowing length in test/random.jl
rfourquet Sep 30, 2014
c754a60
RandIntGen: improve error message
rfourquet Sep 30, 2014
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
103 changes: 45 additions & 58 deletions base/random.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -146,91 +146,78 @@ 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
## Generate random integer within a range 1:n

# 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)
# 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<:Unsigned}
a::T # first element of the range
k::U # range length or zero for full range
k::U # range length
u::U # rejection threshold
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we still get any benefit from having both T and U as type parameters?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The random generation logic is implemented only with 3 types (parameter U): Uint32, Uint64, Uint128, and the generation for all types (param. T) is implemented in term of those. An instance r of RandIntGen needs to know what type of number to produce, in a call like rand(r), this is given by T. It's probably possible to have only one type parameter, but I felt that was beyond the scope of this PR.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That explains it. Thanks!

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)
@assert U in (Uint32, Uint64, Uint128) # respectively 32, mixed 32/64 and 128 bits of entropy
k < 1 && error("No integers in the range [1, $k]. Did you try to pick an element from a empty collection?")
new(k, U == 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
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3946,7 +3946,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 @@ -3960,9 +3960,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
View file
Open in desktop
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