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.

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

RFC: Refactor of random.jl, and stick into modules. #1213

Closed
wants to merge 4 commits into from
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
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
89 changes: 44 additions & 45 deletions base/random.jl → base/RNG.jl
Original file line number Diff line number Diff line change
@@ -1,9 +1,14 @@
_jl_librandom = dlopen("librandom")
@windows_only _jl_advapi32 = dlopen("Advapi32")
module RNG

import Base.*
import Base.LibRandom.*

export librandom_init, rand, rand!, randn, randn!

## initialization

function _jl_librandom_init()
function librandom_init()

@unix_only begin
try
srand("/dev/urandom")
Expand All @@ -18,18 +23,20 @@ function _jl_librandom_init()
end
srand(seed)
end
_jl_randn_zig_init()
end

@windows_only begin
a=zeros(Uint32,2)
ccall(dlsym(_jl_advapi32,:SystemFunction036),stdcall,Uint8,(Ptr{Void},Uint64),convert(Ptr{Void},a),8)
a = zeros(Uint32, 2)
LibRandom.win32_SystemFunction036!(a)
srand(a)
end

LibRandom.randmtzig_create_ziggurat_tables()
end

# macros to generate random arrays

macro _jl_rand_matrix_builder(T, f)
macro rand_matrix_builder(T, f)
f! = esc(symbol("$(f)!"))
f = esc(f)
quote
Expand All @@ -44,7 +51,7 @@ macro _jl_rand_matrix_builder(T, f)
end
end

macro _jl_rand_matrix_builder_1arg(T, f)
macro rand_matrix_builder_1arg(T, f)
f! = esc(symbol("$(f)!"))
f = esc(f)
quote
Expand All @@ -63,14 +70,12 @@ end

function srand(seed::Uint32)
global RANDOM_SEED = seed
ccall(dlsym(_jl_librandom, :dsfmt_gv_init_gen_rand),
Void, (Uint32,), seed)
LibRandom.dsfmt_gv_init_gen_rand(seed)
end

function srand(seed::Vector{Uint32})
global RANDOM_SEED = seed
ccall(dlsym(_jl_librandom, :dsfmt_gv_init_by_array),
Void, (Ptr{Uint32}, Int32), seed, length(seed))
LibRandom.dsfmt_gv_init_by_array(seed)
end

srand(seed::Uint64) = srand([uint32(seed),uint32(seed>>32)])
Expand All @@ -88,20 +93,18 @@ srand(filename::String) = srand(filename, 4)

## rand()

rand() = ccall(dlsym(_jl_librandom, :dsfmt_gv_genrand_close_open), Float64, ())
rand() = LibRandom.dsfmt_gv_genrand_close_open()

const _jl_dsfmt_get_min_array_size =
ccall(dlsym(_jl_librandom, :dsfmt_get_min_array_size), Int32, ())
const dsfmt_min_array_size = LibRandom.dsfmt_get_min_array_size()

function rand!(A::Array{Float64})
n = numel(A)
if n <= _jl_dsfmt_get_min_array_size
if n <= dsfmt_min_array_size()
for i = 1:n
A[i] = rand()
end
else
ccall(dlsym(_jl_librandom, :dsfmt_gv_fill_array_close_open),
Void, (Ptr{Void}, Int32), A, n & 0xfffffffe)
A = LibRandom.dsfmt_gv_fill_array_close_open!(A, n & 0xfffffffe)
if isodd(n)
A[n] = rand()
end
Expand All @@ -114,15 +117,14 @@ rand(dims::Int...) = rand(dims)

## random integers

_jl_dsfmt_randui32() =
ccall(dlsym(_jl_librandom, :dsfmt_gv_genrand_uint32), Uint32, ())
dsfmt_randui32() = dsfmt_gv_genrand_uint32()

_jl_dsfmt_randui64() =
box(Uint64,or_int(zext64(unbox(Uint32,_jl_dsfmt_randui32())),
shl_int(zext64(unbox(Uint32,_jl_dsfmt_randui32())), 32)))
dsfmt_randui64() =
box(Uint64,or_int(zext64(unbox(Uint32,dsfmt_randui32())),
shl_int(zext64(unbox(Uint32,dsfmt_randui32())), 32)))

randi(::Type{Uint32}) = _jl_dsfmt_randui32()
randi(::Type{Uint64}) = _jl_dsfmt_randui64()
randi(::Type{Uint32}) = dsfmt_randui32()
randi(::Type{Uint64}) = dsfmt_randui64()
randi(::Type{Uint128}) = uint128(randi(Uint64))<<64 | randi(Uint64)

randi(::Type{Int32}) = int32(randi(Uint32)) & typemax(Int32)
Expand Down Expand Up @@ -175,39 +177,26 @@ randi(r::(Integer,Integer), dims::Int...) = randival(r[1], r[2], dims)

## random Bools

randbit() = int(_jl_dsfmt_randui32() & 1)
@_jl_rand_matrix_builder Int randbit
randbit() = int(dsfmt_randui32() & 1)
@rand_matrix_builder Int randbit

randbool() = randbit() == 1
@_jl_rand_matrix_builder Bool randbool
@rand_matrix_builder Bool randbool

## randn() - Normally distributed random numbers using Ziggurat algorithm

# The Ziggurat Method for generating random variables - Marsaglia and Tsang
# Paper and reference code: http://www.jstatsoft.org/v05/i08/

_jl_randn_zig_init() =
ccall(dlsym(_jl_librandom, :randmtzig_create_ziggurat_tables), Void, ())

function randn!(A::Array{Float64})
ccall(dlsym(_jl_librandom, :randmtzig_fill_randn),
Void, (Ptr{Float64}, Uint32), A, numel(A))
return A
end

randn() = ccall(dlsym(_jl_librandom, :randmtzig_randn), Float64, ())
randn() = LibRandom.randmtzig_randn()
randn!(A::Array{Float64}) = LibRandom.randmtzig_fill_randn!(A)
randn(dims::Dims) = randn!(Array(Float64, dims))
randn(dims::Int...) = randn(dims)

## randexp() - Exponentially distributed random numbers using Ziggurat algorithm

function randexp!(A::Array{Float64})
ccall(dlsym(_jl_librandom, :randmtzig_fill_exprnd),
Void, (Ptr{Float64}, Uint32), A, numel(A))
return A
end

randexp() = ccall(dlsym(_jl_librandom, :randmtzig_exprnd), Float64, ())
randexp() = LibRandom.randmtzig_exprnd()
randexp!(A::Array{Float64}) = LibRandom.randmtzig_fill_exprnd!(A)
randexp(dims::Dims) = randexp!(Array(Float64, dims))
randexp(dims::Int...) = randexp(dims)

Expand Down Expand Up @@ -247,16 +236,19 @@ function randg!(a::Real, A::Array{Float64})
end
A
end

function randg(a::Real)
if a <= 0. error("shape parameter a must be > 0") end
d = (a < 1. ? a + 1 : a) - 1.0/3.0
randg2(d, 1.0/sqrt(9.0d)) * (a > 1. ? 1. : rand()^(1./a))
end

randg(a::Real, dims::Dims) = randg!(a, Array(Float64, dims))
randg(a::Real, dims::Int...) = randg(a, dims)

## randchi2 - the distribution chi^2(df) is 2*gamma(df/2)
## for integer n, a chi^2(n) is the sum of n squared standard normals

function randchi2!(df::Real, A::Array{Float64})
if df == 1
for i in 1:numel(A)
Expand All @@ -269,9 +261,11 @@ function randchi2!(df::Real, A::Array{Float64})
for i in 1:numel(A) A[i] = 2.randg2(d,c) end
A
end

randchi2(df::Real) = df == 1 ? randn()^2 : 2.randg(df/2.)
randchi2(df::Real, dims::Dims) = randchi2!(df, Array(Float64, dims))
randchi2(df::Real, dims::Int...) = randchi2(df, dims)

const chi2rnd = randchi2 # alias chi2rnd

function randbeta!(alpha::Real, beta::Real, A::Array{Float64})
Expand All @@ -287,8 +281,13 @@ function randbeta!(alpha::Real, beta::Real, A::Array{Float64})
end

randbeta(alpha::Real, beta::Real) = (u=randg(alpha); u/(u + randg(beta)))

function randbeta(alpha::Real, beta::Real, dims::Dims)
randbeta!(alpha, beta, Array(Float64, dims))
end

randbeta(alpha::Real, beta::Real, dims::Int...) = randbeta(alpha, beta, dims)

const betarnd = randbeta

end # module
94 changes: 94 additions & 0 deletions base/RNG_librandom.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,94 @@
module LibRandom

import Base.*
export dsfmt_get_min_array_size, dsfmt_gv_init_gen_rand, dsfmt_gv_init_by_array,
dsfmt_gv_genrand_close_open, dsfmt_gv_genrand_uint32, dsfmt_gv_fill_array_close_open!,
randmtzig_create_ziggurat_tables, randmtzig_randn, randmtzig_fill_randn!,
randmtzig_exprnd, randmtzig_fill_exprnd!,
win32_SystemFunction036!

## DSFMT

function dsfmt_get_min_array_size()
min_array_size = ccall(dlsym(Base.librandom, :dsfmt_get_min_array_size),
Int32,
())
end

function dsfmt_gv_init_gen_rand(seed::Uint32)
ccall(dlsym(Base.librandom, :dsfmt_gv_init_gen_rand),
Void,
(Uint32,),
seed)
end

function dsfmt_gv_init_by_array(seed::Vector{Uint32})
ccall(dlsym(Base.librandom, :dsfmt_gv_init_by_array),
Void,
(Ptr{Uint32}, Int32),
seed, length(seed))
end

function dsfmt_gv_genrand_close_open()
r = ccall(dlsym(Base.librandom, :dsfmt_gv_genrand_close_open),
Float64,
())
end

function dsfmt_gv_genrand_uint32()
r = ccall(dlsym(Base.librandom, :dsfmt_gv_genrand_uint32),
Uint32,
())
end

function dsfmt_gv_fill_array_close_open!(A::Array{Float64})
ccall(dlsym(Base.librandom, :dsfmt_gv_fill_array_close_open),
Void,
(Ptr{Void}, Int32),
A, n)
return A
end

## randmtzig

function randmtzig_create_ziggurat_tables()
ccall(dlsym(Base.librandom, :randmtzig_create_ziggurat_tables),
Void,
())
end

function randmtzig_randn()
call(dlsym(Base.librandom, :randmtzig_randn),
Float64,
())
end

function randmtzig_fill_randn!(A)
ccall(dlsym(Base.librandom, :randmtzig_fill_randn),
Void,
(Ptr{Float64}, Uint32),
A, numel(A))
end

function randmtzig_exprnd()
ccall(dlsym(Base.librandom, :randmtzig_exprnd),
Float64,
())
end

function randmtzig_fill_exprnd!(A)
ccall(dlsym(Base.librandom, :randmtzig_fill_exprnd),
Void,
(Ptr{Float64}, Uint32),
A, numel(A))
end

## Windows entropy

@windows_only begin
function win32_SystemFunction036!(a::Array{Uint32})
ccall(dlsym(advapi32,:SystemFunction036),stdcall,Uint8,(Ptr{Void},Uint64),convert(Ptr{Void},a),8)
end
end

end # module
2 changes: 1 addition & 1 deletion base/client.jl
Original file line number Diff line number Diff line change
Expand Up @@ -198,7 +198,7 @@ function _start()
# set CPU core count
global const CPU_CORES = ccall(:jl_cpu_cores, Int32, ())

_jl_librandom_init()
RNG.librandom_init()

atexit(()->flush(stdout_stream))
try
Expand Down
2 changes: 2 additions & 0 deletions base/export.jl
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,8 @@ export
Grisu,
Printf,
PCRE,
LibRandom,
RNG,
FFTW,
DSP,

Expand Down
8 changes: 4 additions & 4 deletions base/start_image.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,12 +6,12 @@ _jl_lib = ccall(:jl_load_dynamic_library,Ptr{Void},(Ptr{None},),C_NULL)
@windows_only _jl_repl = ccall(:GetModuleHandleA,stdcall,Ptr{Void},(Ptr{Void},),C_NULL)

# Essential libraries
libpcre = dlopen("libpcre")
libgrisu = dlopen("libgrisu")
const libpcre = dlopen("libpcre")
const libgrisu = dlopen("libgrisu")
const librandom = dlopen("librandom");
@windows_only _jl_advapi32 = dlopen("Advapi32")
_jl_libm = dlopen("libm")
_jl_libfdm = dlopen("libfdm")
_jl_librandom = dlopen("librandom");
@windows_only _jl_advapi32 = dlopen("Advapi32")

# Optional libraries
const _jl_libblas = dlopen(_jl_libblas_name)
Expand Down
11 changes: 7 additions & 4 deletions base/sysimg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -111,7 +111,10 @@ include("combinatorics.jl")
include("statistics.jl")

# random number generation
include("random.jl")
include("RNG_librandom.jl")
include("RNG.jl")
import Base.RNG
import Base.RNG.*

# distributed arrays and memory-mapped arrays
include("darray.jl")
Expand Down Expand Up @@ -172,10 +175,10 @@ compile_hint(cmp, (Int32, Int32))
compile_hint(min, (Int32, Int32))
compile_hint(==, (ASCIIString, ASCIIString))
compile_hint(arg_gen, (ASCIIString,))
compile_hint(_jl_librandom_init, ())
compile_hint(srand, (ASCIIString, Int))
#compile_hint(librandom_init, ())
#compile_hint(srand, (ASCIIString, Int))
compile_hint(open, (ASCIIString, Bool, Bool, Bool, Bool))
compile_hint(srand, (Uint64,))
#compile_hint(srand, (Uint64,))
compile_hint(done, (IntSet, Int64))
compile_hint(next, (IntSet, Int64))
compile_hint(ht_keyindex, (Dict{Any,Any}, Int32))
Expand Down