Skip to content

Commit

Permalink
Xoshiro: allow any non-negative integer as a seed, via SHA2_256 (Juli…
Browse files Browse the repository at this point in the history
  • Loading branch information
rfourquet committed Sep 23, 2021
1 parent b55497d commit 2174ee1
Show file tree
Hide file tree
Showing 13 changed files with 53 additions and 94 deletions.
10 changes: 5 additions & 5 deletions doc/src/devdocs/subarrays.md
Original file line number Diff line number Diff line change
Expand Up @@ -19,14 +19,14 @@ julia> A = rand(2,3,4);
julia> S1 = view(A, :, 1, 2:3)
2×2 view(::Array{Float64, 3}, :, 1, 2:3) with eltype Float64:
0.342284 0.831961
0.237287 0.435938
0.839622 0.711389
0.967143 0.103929
julia> S2 = view(A, 1, :, 2:3)
3×2 view(::Array{Float64, 3}, 1, :, 2:3) with eltype Float64:
0.342284 0.831961
0.988944 0.927795
0.178426 0.404876
0.839622 0.711389
0.789764 0.806704
0.566704 0.962715
```
```@meta
DocTestSetup = nothing
Expand Down
20 changes: 10 additions & 10 deletions doc/src/manual/performance-tips.md
Original file line number Diff line number Diff line change
Expand Up @@ -77,12 +77,12 @@ julia> function sum_global()
end;
julia> @time sum_global()
0.010414 seconds (9.07 k allocations: 373.448 KiB, 98.40% compilation time)
493.6199223951192
0.011539 seconds (9.08 k allocations: 373.386 KiB, 98.69% compilation time)
523.0007221951678
julia> @time sum_global()
0.000108 seconds (3.49 k allocations: 70.156 KiB)
493.6199223951192
0.000091 seconds (3.49 k allocations: 70.156 KiB)
523.0007221951678
```

On the first call (`@time sum_global()`) the function gets compiled. (If you've not yet used [`@time`](@ref)
Expand Down Expand Up @@ -113,12 +113,12 @@ julia> function sum_arg(x)
end;
julia> @time sum_arg(x)
0.007971 seconds (3.96 k allocations: 200.171 KiB, 99.83% compilation time)
493.6199223951192
0.007551 seconds (3.98 k allocations: 200.548 KiB, 99.77% compilation time)
523.0007221951678
julia> @time sum_arg(x)
0.000003 seconds (1 allocation: 16 bytes)
493.6199223951192
0.000006 seconds (1 allocation: 16 bytes)
523.0007221951678
```

The 1 allocation seen is from running the `@time` macro itself in global scope. If we instead run
Expand All @@ -128,8 +128,8 @@ the timing in a function, we can see that indeed no allocations are performed:
julia> time_sum(x) = @time sum_arg(x);
julia> time_sum(x)
0.000001 seconds
493.6199223951192
0.000002 seconds
523.0007221951678
```

In some situations, your function may need to allocate memory as part of its operation, and this
Expand Down
8 changes: 4 additions & 4 deletions stdlib/LinearAlgebra/test/dense.jl
Original file line number Diff line number Diff line change
Expand Up @@ -22,10 +22,10 @@ Random.seed!(1234323)
@testset "for $elty" for elty in (Float32, Float64, ComplexF32, ComplexF64)
ainit = convert(Matrix{elty}, ainit)
for a in (copy(ainit), view(ainit, 1:n, 1:n))
@test cond(a,1) 122.15725126320953 atol=0.5
@test cond(a,2) 78.44837047684149 atol=0.5
@test cond(a,Inf) 174.10761543202744 atol=0.4
@test cond(a[:,1:5]) 6.7492840150789135 atol=0.01
@test cond(a,1) 198.3324294531168 atol=0.5
@test cond(a,2) 85.93920079319506 atol=0.5
@test cond(a,Inf) 149.7523084803039 atol=0.4
@test cond(a[:,1:5]) 8.319279144493297 atol=0.01
@test_throws ArgumentError cond(a,3)
end
end
Expand Down
2 changes: 1 addition & 1 deletion stdlib/LinearAlgebra/test/eigen.jl
Original file line number Diff line number Diff line change
Expand Up @@ -163,7 +163,7 @@ end
end

@testset "eigen of an Adjoint" begin
Random.seed!(1)
Random.seed!(4)
A = randn(3,3)
@test eigvals(A') == eigvals(copy(A'))
@test eigen(A') == eigen(copy(A'))
Expand Down
1 change: 1 addition & 0 deletions stdlib/Random/Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@ uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

[deps]
Serialization = "9e88b42a-f829-5b0c-bbe9-9e923198166b"
SHA = "ea8e919c-243c-51af-8825-aaa63cd721ce"

[extras]
Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
Expand Down
18 changes: 9 additions & 9 deletions stdlib/Random/docs/src/index.md
Original file line number Diff line number Diff line change
Expand Up @@ -151,22 +151,22 @@ Scalar and array methods for `Die` now work as expected:

```jldoctest Die; setup = :(Random.seed!(1))
julia> rand(Die)
Die(7)
Die(5)
julia> rand(MersenneTwister(0), Die)
Die(11)
julia> rand(Die, 3)
3-element Vector{Die}:
Die(13)
Die(8)
Die(20)
Die(9)
Die(15)
Die(14)
julia> a = Vector{Die}(undef, 3); rand!(a)
3-element Vector{Die}:
Die(4)
Die(14)
Die(10)
Die(19)
Die(7)
Die(17)
```

#### A simple sampler without pre-computed data
Expand All @@ -183,9 +183,9 @@ julia> rand(Die(4))
julia> rand(Die(4), 3)
3-element Vector{Any}:
3
2
4
3
3
```

Given a collection type `S`, it's currently assumed that if `rand(::S)` is defined, an object of type `eltype(S)` will be produced. In the last example, a `Vector{Any}` is produced; the reason is that `eltype(Die) == Any`. The remedy is to define `Base.eltype(::Type{Die}) = Int`.
Expand Down
2 changes: 1 addition & 1 deletion stdlib/Random/src/RNGs.jl
Original file line number Diff line number Diff line change
Expand Up @@ -384,7 +384,7 @@ end

seed!(rng::_GLOBAL_RNG, ::Nothing) = seed!(rng) # to resolve ambiguity

seed!(seed::Union{Nothing,Integer,Vector{UInt32},Vector{UInt64},NTuple{4,UInt64}}=nothing) =
seed!(seed::Union{Nothing,Integer,Vector{UInt32},Vector{UInt64}}=nothing) =
seed!(GLOBAL_RNG, seed)

rng_native_52(::_GLOBAL_RNG) = rng_native_52(default_rng())
Expand Down
1 change: 1 addition & 0 deletions stdlib/Random/src/Random.jl
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,7 @@ include("DSFMT.jl")
using .DSFMT
using Base.GMP.MPZ
using Base.GMP: Limb
import SHA

using Base: BitInteger, BitInteger_types, BitUnsigned, require_one_based_indexing

Expand Down
60 changes: 8 additions & 52 deletions stdlib/Random/src/Xoshiro.jl
Original file line number Diff line number Diff line change
Expand Up @@ -117,65 +117,21 @@ end

# Shared implementation between Xoshiro and TaskLocalRNG -- seeding

function seed!(x::Union{TaskLocalRNG,Xoshiro})
function seed!(rng::Union{TaskLocalRNG,Xoshiro})
# as we get good randomness from RandomDevice, we can skip hashing
parent = RandomDevice()
# Constants have nothing up their sleeve, see task.c
# 0x02011ce34bce797f == hash(UInt(1))|0x01
# 0x5a94851fb48a6e05 == hash(UInt(2))|0x01
# 0x3688cf5d48899fa7 == hash(UInt(3))|0x01
# 0x867b4bb4c42e5661 == hash(UInt(4))|0x01
setstate!(x,
0x02011ce34bce797f * rand(parent, UInt64),
0x5a94851fb48a6e05 * rand(parent, UInt64),
0x3688cf5d48899fa7 * rand(parent, UInt64),
0x867b4bb4c42e5661 * rand(parent, UInt64))
rd = RandomDevice()
setstate!(rng, rand(rd, UInt64), rand(rd, UInt64), rand(rd, UInt64), rand(rd, UInt64))
end

function seed!(rng::Union{TaskLocalRNG,Xoshiro}, seed::NTuple{4,UInt64})
# TODO: Consider a less ad-hoc construction
# We can afford burning a handful of cycles here, and we don't want any
# surprises with respect to bad seeds / bad interactions.

s0 = s = Base.hash_64_64(seed[1])
s1 = s += Base.hash_64_64(seed[2])
s2 = s += Base.hash_64_64(seed[3])
s3 = s += Base.hash_64_64(seed[4])

function seed!(rng::Union{TaskLocalRNG,Xoshiro}, seed::Union{Vector{UInt32}, Vector{UInt64}})
c = SHA.SHA2_256_CTX()
SHA.update!(c, reinterpret(UInt8, seed))
s0, s1, s2, s3 = reinterpret(UInt64, SHA.digest!(c))
setstate!(rng, s0, s1, s2, s3)

rand(rng, UInt64)
rand(rng, UInt64)
rand(rng, UInt64)
rand(rng, UInt64)
rng
end

function seed!(rng::Union{TaskLocalRNG, Xoshiro}, seed::UInt128)
seed0 = seed % UInt64
seed1 = (seed>>>64) % UInt64
seed!(rng, (seed0, seed1, zero(UInt64), zero(UInt64)))
end
seed!(rng::Union{TaskLocalRNG, Xoshiro}, seed::Integer) = seed!(rng, UInt128(seed))

function seed!(rng::Union{TaskLocalRNG, Xoshiro}, seed::AbstractVector{UInt64})
if length(seed) > 4
throw(ArgumentError("seed should have no more than 256 bits"))
end
seed0 = length(seed)>0 ? seed[1] : UInt64(0)
seed1 = length(seed)>1 ? seed[2] : UInt64(0)
seed2 = length(seed)>2 ? seed[3] : UInt64(0)
seed3 = length(seed)>3 ? seed[4] : UInt64(0)
seed!(rng, (seed0, seed1, seed2, seed3))
end
seed!(rng::Union{TaskLocalRNG, Xoshiro}, seed::Integer) = seed!(rng, make_seed(seed))

function seed!(rng::Union{TaskLocalRNG, Xoshiro}, seed::AbstractVector{UInt32})
if iseven(length(seed))
seed!(rng, reinterpret(UInt64, seed))
else
seed!(rng, UInt64[reinterpret(UInt64, @view(seed[begin:end-1])); seed[end] % UInt64])
end
end

@inline function rand(rng::Union{TaskLocalRNG, Xoshiro}, ::SamplerType{UInt128})
first = rand(rng, UInt64)
Expand Down
4 changes: 2 additions & 2 deletions stdlib/Random/src/misc.jl
Original file line number Diff line number Diff line change
Expand Up @@ -53,13 +53,13 @@ number generator, see [Random Numbers](@ref).
# Examples
```jldoctest
julia> Random.seed!(3); randstring()
"h8BzxSoS"
"Lxz5hUwn"
julia> randstring(MersenneTwister(3), 'a':'z', 6)
"ocucay"
julia> randstring("ACGT")
"CTTACTGC"
"TGCTCCTC"
```
!!! note
Expand Down
4 changes: 2 additions & 2 deletions stdlib/Random/test/runtests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -688,7 +688,7 @@ end
@testset "$RNG(seed) & Random.seed!(m::$RNG, seed) produce the same stream" for RNG=(MersenneTwister,Xoshiro)
seeds = Any[0, 1, 2, 10000, 10001, rand(UInt32, 8), rand(UInt128, 3)...]
if RNG == Xoshiro
push!(seeds, rand(UInt64, rand(1:4)), Tuple(rand(UInt64, 4)))
push!(seeds, rand(UInt64, rand(1:4)))
end
for seed=seeds
m = RNG(seed)
Expand All @@ -699,7 +699,7 @@ end
end

@testset "Random.seed!(seed) sets Random.GLOBAL_SEED" begin
seeds = Any[0, rand(UInt128), rand(UInt64, 4), Tuple(rand(UInt64, 4))]
seeds = Any[0, rand(UInt128), rand(UInt64, 4)]

for seed=seeds
Random.seed!(seed)
Expand Down
15 changes: 8 additions & 7 deletions stdlib/SparseArrays/src/sparsematrix.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1642,12 +1642,13 @@ argument specifies a random number generator, see [Random Numbers](@ref).
```jldoctest; setup = :(using Random; Random.seed!(1234))
julia> sprand(Bool, 2, 2, 0.5)
2×2 SparseMatrixCSC{Bool, Int64} with 2 stored entries:
1
1
1 1
julia> sprand(Float64, 3, 0.75)
3-element SparseVector{Float64, Int64} with 1 stored entry:
[3] = 0.787459
3-element SparseVector{Float64, Int64} with 2 stored entries:
[1] = 0.795547
[2] = 0.49425
```
"""
function sprand(r::AbstractRNG, m::Integer, n::Integer, density::AbstractFloat, rfn::Function, ::Type{T}=eltype(rfn(r, 1))) where T
Expand Down Expand Up @@ -1688,9 +1689,9 @@ argument specifies a random number generator, see [Random Numbers](@ref).
# Examples
```jldoctest; setup = :(using Random; Random.seed!(0))
julia> sprandn(2, 2, 0.75)
2×2 SparseMatrixCSC{Float64, Int64} with 1 stored entry:
⋅ 1.32078
2×2 SparseMatrixCSC{Float64, Int64} with 3 stored entries:
-1.20577
0.311817 -0.234641
```
"""
sprandn(r::AbstractRNG, m::Integer, n::Integer, density::AbstractFloat) =
Expand Down
2 changes: 1 addition & 1 deletion stdlib/SparseArrays/test/sparse.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1750,7 +1750,7 @@ end
A = guardseed(1234321) do
triu(sprand(10, 10, 0.2))
end
@test getcolptr(SparseArrays.droptol!(A, 0.01)) == [1, 1, 3, 4, 5, 6, 7, 11, 13, 15, 18]
@test getcolptr(SparseArrays.droptol!(A, 0.01)) == [1, 1, 1, 1, 2, 2, 2, 4, 4, 5, 5]
@test isequal(SparseArrays.droptol!(sparse([1], [1], [1]), 1), SparseMatrixCSC(1, 1, Int[1, 1], Int[], Int[]))
end

Expand Down

0 comments on commit 2174ee1

Please sign in to comment.