From 94d90087fa40c57cf12d342dd10b966b49df6636 Mon Sep 17 00:00:00 2001
From: Rafael Fourquet <fourquet.rafael@gmail.com>
Date: Sun, 5 Oct 2014 16:49:12 +0530
Subject: [PATCH] re-add #8309 (which was reverted)

This reverts commit 2abb59bc443fa9ebf876ebc6a45f6bde18ac2c88.
This reverts commit 611732862e20adc8403adf0377981cdd6acf859d.
---
 base/random.jl      | 106 +++++++++++++++++++-------------------------
 doc/stdlib/base.rst |  10 +++--
 test/random.jl      |  28 +++++++++---
 3 files changed, 74 insertions(+), 70 deletions(-)

diff --git a/base/random.jl b/base/random.jl
index 3c71008283e97e..3d4e4d4bbcb7cc 100644
--- a/base/random.jl
+++ b/base/random.jl
@@ -24,7 +24,7 @@ type MersenneTwister <: AbstractRNG
     MersenneTwister(seed=0) = MersenneTwister(make_seed(seed))
 end
 
-function srand(r::MersenneTwister, seed)
+function srand(r::MersenneTwister, seed) 
     r.seed = seed
     dsfmt_init_gen_rand(r.state, seed)
     return r
@@ -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
+## 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
+# 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)
 
-# 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
+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
diff --git a/doc/stdlib/base.rst b/doc/stdlib/base.rst
index 5e1dfca01a0251..400109d639af76 100644
--- a/doc/stdlib/base.rst
+++ b/doc/stdlib/base.rst
@@ -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...])
 
@@ -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...])
 
diff --git a/test/random.jl b/test/random.jl
index 9e9970f2081724..6329470dc3dbeb 100644
--- a/test/random.jl
+++ b/test/random.jl
@@ -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)
@@ -44,9 +49,18 @@ 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
 
@@ -54,8 +68,8 @@ 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
 
@@ -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