From 15702648006ffacfaa8939d065472dad453b2b76 Mon Sep 17 00:00:00 2001 From: mschauer Date: Wed, 11 Nov 2015 21:26:20 +0100 Subject: [PATCH] Reintroduce random BigFloats removed in 4cc21ca --- base/gmp.jl | 39 +++++++++++++++++++++++++++++++++++++++ base/mpfr.jl | 11 ++++++++++- base/random.jl | 35 ++++++++++++++++++++++++++++++++++- 3 files changed, 83 insertions(+), 2 deletions(-) diff --git a/base/gmp.jl b/base/gmp.jl index 02d2016124d6c..b6e4e0b259e5f 100644 --- a/base/gmp.jl +++ b/base/gmp.jl @@ -49,6 +49,23 @@ type BigInt <: Integer end end +type GMPRandState + seed_alloc::Cint + seed_size::Cint + seed_d::Ptr{Void} + alg::Cint + alg_data::Ptr{Void} + + function GMPRandState() + randstate = new(zero(Cint), zero(Cint), C_NULL, zero(Cint), C_NULL) + # The default algorithm in GMP is currently MersenneTwister + ccall((:__gmp_randinit_default, :libgmp), Void, (Ptr{GMPRandState},), + &randstate) + finalizer(randstate, randclear) + randstate + end +end + _gmp_clear_func = C_NULL _mpfr_clear_func = C_NULL @@ -553,4 +570,26 @@ Base.checked_add(a::BigInt, b::BigInt) = a + b Base.checked_sub(a::BigInt, b::BigInt) = a - b Base.checked_mul(a::BigInt, b::BigInt) = a * b + +## GMPRandState + +function GMPRandState(seed::UInt64) + randstate = GMPRandState() + ccall((:__gmp_randseed_ui, :libgmp), Void, (Ptr{GMPRandState}, Culong), + &randstate, seed) + randstate +end + +function GMPRandState(seed::BigInt) + randstate = GMPRandState() + ccall((:__gmp_randseed, :libgmp), Void, (Ptr{GMPRandState}, Ptr{BigInt}), + &randstate, &seed) + randstate +end + + +function randclear(x::GMPRandState) + ccall((:__gmp_randclear, :libgmp), Void, (Ptr{GMPRandState},), &x) +end + end # module diff --git a/base/mpfr.jl b/base/mpfr.jl index 265d4cf990eef..6c3577853657b 100644 --- a/base/mpfr.jl +++ b/base/mpfr.jl @@ -22,7 +22,7 @@ import import Base.Rounding: rounding_raw, setrounding_raw -import Base.GMP: ClongMax, CulongMax, CdoubleMax +import Base.GMP: ClongMax, CulongMax, CdoubleMax, GMPRandState import Base.Math.lgamma_r @@ -880,4 +880,13 @@ get_emin_max() = ccall((:mpfr_get_emin_max, :libmpfr), Clong, ()) set_emax!(x) = ccall((:mpfr_set_emax, :libmpfr), Void, (Clong,), x) set_emin!(x) = ccall((:mpfr_set_emin, :libmpfr), Void, (Clong,), x) + +function urandomb(randstate::GMPRandState) + z = BigFloat() + ccall((:mpfr_urandomb,:libmpfr), Int32, + (Ptr{BigFloat}, Ptr{GMPRandState}), + &z, &randstate) + z +end + end #module diff --git a/base/random.jl b/base/random.jl index 9eaceb0c0a51b..aa01513249d58 100644 --- a/base/random.jl +++ b/base/random.jl @@ -3,7 +3,8 @@ module Random using Base.dSFMT -using Base.GMP: GMP_VERSION, Limb +using Base.GMP: GMP_VERSION, Limb, GMPRandState +using Base.MPFR: urandomb export srand, rand, rand!, @@ -126,6 +127,36 @@ function randjump(mt::MersenneTwister, jumps::Integer, jumppoly::AbstractString) end randjump(r::MersenneTwister, jumps::Integer) = randjump(r, jumps, dSFMT.JPOLY1e21) +## BigFloat using GMPRandState + +immutable BigFloatRNG <: AbstractRNG + randstate::GMPRandState + seed::Vector{UInt32} + function BigFloatRNG(seed_::Vector{UInt32}) + s = big(0) + for i in seed_ + s = s << 32 + i + end + randstate = GMPRandState(s) + gmprng = new(randstate, seed_) + end + BigFloatRNG() = BigFloatRNG([UInt32(0)]) +end + +function srand(r::BigFloatRNG, seed::Vector{UInt32}) + r.seed = seed + s = big(0) + for i in seed + s = s << 32 + i + end + r.randstate = GMPRandState(s) + return r +end + +function rand(r::BigFloatRNG, ::Type{BigFloat}) + urandomb(r.randstate) +end + ## initialization function __init__() @@ -1364,4 +1395,6 @@ function randcycle(r::AbstractRNG, n::Integer) end randcycle(n::Integer) = randcycle(GLOBAL_RNG, n) + + end # module