Use only DSFMT. Remove random stuff from random.

Allow seeding entropy pool from an array of seeds, and from a file (/dev/urandom).
JuliaLang · Jul 9, 2011 · b13d037 · b13d037
1 parent 37c543e
commit b13d037
Show file tree

Hide file tree

Showing 3 changed files with 33 additions and 113 deletions.
diff --git a/external/Makefile b/external/Makefile
@@ -127,7 +127,7 @@ clean-fdlibm:
 	rm -f $(FDLIBM_OBJ_SOURCE) $(FDLIBM_OBJ_TARGET)
 distclean-fdlibm: clean-fdlibm
 
-## MT ##
+## DSFMT ##
 
 DSFMT_OBJ_TARGET = $(EXTROOTLIB)/libMT.$(SHLIB_EXT)
 DSFMT_OBJ_SOURCE = random/libMT.$(SHLIB_EXT)
@@ -137,13 +137,12 @@ install-dsfmt: $(DSFMT_OBJ_TARGET)
 
 random/dsfmt-$(DSFMT_VER).tar.gz:
 	cd random && \
-	curl -o dsfmt-$(DSFMT_VER).tar.gz http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/dSFMT-src-$(DSFMT_VER).tar.gz -O http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/mt19937-64.tgz 
+	curl -o dsfmt-$(DSFMT_VER).tar.gz http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/dSFMT-src-$(DSFMT_VER).tar.gz 
 random/jl_random.c: random/dsfmt-$(DSFMT_VER).tar.gz
 	cd random && \
 	mkdir -p dsfmt-$(DSFMT_VER) && \
 	tar -C dsfmt-$(DSFMT_VER) --strip-components 1 -xf dsfmt-$(DSFMT_VER).tar.gz && \
-	cp dSFMT.h dsfmt-$(DSFMT_VER) && \
-	tar zxf mt19937-64.tgz
+	cp dSFMT.h dsfmt-$(DSFMT_VER)
 	touch $@
 $(DSFMT_OBJ_SOURCE): random/jl_random.c
 	cd random && \
@@ -156,7 +155,7 @@ $(DSFMT_OBJ_TARGET): $(DSFMT_OBJ_SOURCE)
 clean-dsfmt:
 	rm -f random/libMT.$(SHLIB_EXT)
 distclean-dsfmt: clean-dsfmt
-	cd random && rm -rf *.tgz *.tar.gz dsfmt-$(DSFMT_VER) mt19937-64
+	cd random && rm -rf *.tar.gz dsfmt-$(DSFMT_VER)
 
 ## OpenBLAS ##
 

diff --git a/external/random/jl_random.c b/external/random/jl_random.c
@@ -3,83 +3,24 @@
 */
 
 #include "jl_random.h"
-
-#include "mt19937ar.c"
 #include "dsfmt-2.1/dSFMT.c"
 
-double rand_double()
-{
-    union ieee754_double d;
-
-    d.ieee.mantissa0 = genrand_int32();
-    d.ieee.mantissa1 = genrand_int32();
-    d.ieee.negative = 0;
-    d.ieee.exponent = IEEE754_DOUBLE_BIAS + 0;    /* 2^0 */
-    return d.d - 1.0;
-}
-
-float rand_float()
-{
-    union ieee754_float f;
+#define RANDN_RESET -99999999
 
-    f.ieee.mantissa = genrand_int32();
-    f.ieee.negative = 0;
-    f.ieee.exponent = IEEE754_FLOAT_BIAS + 0;     /* 2^0 */
-    return f.f - 1.0;
-}
+static double randn_next = RANDN_RESET;
 
-static double randn_next = -42;
-double randn()
+void dsfmt_randn_reset()
 {
-    double s, vre, vim, ure, uim;
-
-    if (randn_next != -42) {
-        s = randn_next;
-        randn_next = -42;
-        return s;
-    }
-    do {
-        ure = rand_double();
-        uim = rand_double();
-        vre = 2*ure - 1;
-        vim = 2*uim - 1;
-        s = vre*vre + vim*vim;
-    } while (s >= 1);
-    s = sqrt(-2*log(s)/s);
-    randn_next = s * vre;
-    return s * vim;
-}
-
-void randomseed32(uint32_t s)
-{
-    init_genrand(s);
-    randn_next = -42;
-}
-
-void randomseed64(uint64_t s)
-{
-    init_by_array((uint32_t*)&s, 2);
-    randn_next = -42;
-}
-
-void randomize()
-{
-    struct timeval now;
-    uint64_t a;
-
-    gettimeofday(&now, NULL);
-    a = (((uint64_t)now.tv_sec)<<32) + (uint64_t)now.tv_usec;
-
-    randomseed64(a);
+    randn_next = RANDN_RESET;
 }
 
 double dsfmt_randn()
 {
     double s, vre, vim, ure, uim;
 
-    if (randn_next != -42) {
+    if (randn_next != RANDN_RESET) {
         s = randn_next;
-        randn_next = -42;
+        randn_next = RANDN_RESET;
         return s;
     }
     do {

diff --git a/j/random.j b/j/random.j
@@ -1,16 +1,24 @@
 libmt = dlopen("libMT")
 
 function mt_init()
-    #randomize()
-    srand(0)
+    srand(uint32(clock()))
+#    try
+#        srand("/dev/urandom", 4)
+#    catch
+#        srand(uint32(clock()))
+#    end
 end
 
-### DSFMT ###
-
 dsfmt_get_min_array_size() = ccall(dlsym(libmt, :dsfmt_get_min_array_size), Int32, ())
 
-srand(seed::Union(Int32,Uint32)) = ccall(dlsym(libmt, :dsfmt_gv_init_gen_rand), 
-                                         Void, (Uint32, ), uint32(seed))
+dsfmt_randn_reset() = ccall(dlsym(libmt, :dsfmt_randn_reset), Void, ())
+
+srand(seed::Uint32) = (ccall(dlsym(libmt, :dsfmt_gv_init_gen_rand), Void, (Uint32, ), seed);
+                       dsfmt_randn_reset())
+
+srand(seed::Vector{Uint32}) = (ccall(dlsym(libmt, :dsfmt_gv_init_by_array),
+                                     Void, (Vector{Uint32}, Int32), seed, length(seed));
+                               dsfmt_randn_reset())
 
 randf() = float32(rand())
 
@@ -37,43 +45,14 @@ function dsfmt_fill_array_open_open(A::Array{Float64})
     return A
 end
 
-# jl_randn_next = -42.0
-# function randn()
-#     global jl_randn_next
-#
-#     if (jl_randn_next != -42.0)
-#         s = jl_randn_next
-#         jl_randn_next = -42.0
-#         return s
-#     end
-#
-#     s = 1.0
-#     vre = 0.0
-#     vim = 0.0
-#     while (s >= 1.0)
-#         ure = rand()
-#         uim = rand()
-#         vre = 2.0*ure - 1.0
-#         vim = 2.0*uim - 1.0
-#         s = vre*vre + vim*vim
-#     end
-#
-#     s = sqrt(-2.0*log(s)/s)
-#     jl_randn_next = s * vre
-#     return s * vim
-# end
-
-
-### MT ###
-# This is the old code based on the original Mersenne Twister
-
-#randomize() = ccall(dlsym(libmt, :randomize), Void, ())
-#rand()     = ccall(dlsym(libmt, :rand_double),   Float64, ())
-#randf()    = ccall(dlsym(libmt, :rand_float),    Float32, ())
-#randui32() = ccall(dlsym(libmt, :genrand_int32), Uint32,  ())
-#randn()    = ccall(dlsym(libmt, :randn),         Float64, ())
-#srand(s::Union(Int32,Uint32)) = ccall(dlsym(libmt, :randomseed32), Void, (Uint32,), uint32(s))
-#srand(s::Union(Int64,Uint64)) = ccall(dlsym(libmt, :randomseed64), Void, (Uint64,), uint64(s))
+## Seed from a file
+
+function srand(fname::String, n::Int32)
+    fid = open(fname)
+    a = Array(Uint32, n)
+    read(fid, a)
+    srand(a)
+end
 
 ## Random integers
 
@@ -138,3 +117,4 @@ end # macro
 @rand_matrix_builder Float64 randn
 @rand_matrix_builder Float32 randf
 @rand_matrix_builder Uint32 randui32
+@rand_matrix_builder Uint64 randui64