WIP: Reintroduce random BigFloats removed in 1a8885b

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

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

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