Initial commit.

rebuild_random.py

@@ -0,0 +1,68 @@
+import gzip
+import random
+from itertools import imap
+
+
+print "Loading Magic"
+f = gzip.GzipFile("magic_data","r")
+magic = eval(f.read())
+f.close()
+print "Done."
+
+def rebuild_random(data):
+    assert len(data) >= 3360
+    vals = list(imap(ord, data))
+
+    def getbit(bit):
+        assert bit >= 0
+        return (vals[bit // 8] >> (7 - bit % 8)) & 1
+
+    state = []
+    for i in xrange(0, 624):
+        print "REBUILDING RANDOM-POOL ["+ ("#" * (i // 10)).ljust(62) +"]"
+        val = 0
+        data = magic[i % 2]
+        for bit in data:
+            val <<= 1
+            for b in bit:
+                val ^= getbit(b+(i//2)*8 - 8)
+        state.append(val)
+
+    state.append(0)
+    ran = random.Random()
+    ran.setstate((3, tuple(state),None))
+    for i in xrange(len(vals) - 3201 + 394):
+        _ = ran.randint(0,255)
+    return ran
+
+def random_string(length, random_module = random):
+    return "".join(chr(random_module.randint(0, 255)) for i in xrange(length))
+
+
+## OK... here's how it works
+
+# Shuffle the random-state a little bit
+random_string(random.randint(0,10000))
+
+# First we receive 3500 bytes from our random-function
+first_random_string = random_string(3360)
+
+# and put it into our magic function. It returns a Random-object
+# that is in the same state as the other random-object, reconstructed
+# out of the random-strings
+
+my_random = rebuild_random(first_random_string)
+
+# Now we expect this string
+expected_string = random_string(10000, my_random)
+
+# Let's see...
+second_random_string = random_string(10000)
+
+# if it matches.
+if expected_string == second_random_string:
+    print "RANDOM POOL SUCCESSFULLY REBUILT!"
+else:
+    print "Should not happen"
+
+

rebuild_random_multicore.py

@@ -0,0 +1,76 @@
+import gzip
+import random
+import multiprocessing
+from itertools import imap
+
+
+print "Loading Magic"
+f = gzip.GzipFile("magic_data","r")
+magic = eval(f.read())
+f.close()
+print "Done."
+
+def calc(i):
+    def getbit(bit):
+        assert bit >= 0
+        return (vals[bit // 8] >> (7 - bit % 8)) & 1
+    val = 0
+    data = magic[i % 2]
+    for bit in data:
+        val <<= 1
+        for b in bit:
+            val ^= getbit(b+(i//2)*8 - 8)
+    return val
+
+def rebuild_random(data):
+    global vals
+    assert len(data) >= 3360
+    vals = list(imap(ord, data))
+
+    state = []
+    #real_state = random._inst.getstate()[1]
+    pool = multiprocessing.Pool()
+    stepsize = multiprocessing.cpu_count() * 3
+
+    for i in xrange(0, 624, stepsize):
+        print "REBUILDING RANDOM-POOL ["+ ("#" * (i // 10)).ljust(62) +"]"
+        res = pool.map(calc, range(i, min(624, i + stepsize)))
+        state.extend(res)
+        #print val, real_state[i]
+
+    state.append(0)
+    ran = random.Random()
+    ran.setstate((3, tuple(state),None))
+    for i in xrange(len(vals) - 3201 + 394):
+        _ = ran.randint(0,255)
+    return ran
+
+def random_string(length, random_module = random):
+    return "".join(chr(random_module.randint(0, 255)) for i in xrange(length))
+
+
+## OK... here's how it works
+
+# Shuffle the random-state a little bit
+random_string(random.randint(0,10000))
+
+# First we receive 3500 bytes from our random-function
+first_random_string = random_string(3500)
+
+# and put it into our magic function. It returns a Random-object
+# that is in the same state as the other random-object, reconstructed
+# out of the random-strings
+
+my_random = rebuild_random(first_random_string)
+
+# Now we expect this string
+expected_string = random_string(10000, my_random)
+
+# Let's see...
+second_random_string = random_string(10000)
+
+# if it matches.
+if expected_string == second_random_string:
+    print "RANDOM POOL SUCCESSFULLY REBUILT!"
+else:
+    print "Should not happen"

twister.py

@@ -0,0 +1,106 @@
+"""
+Implementation of the Mersenne Twister that keeps track of every bit-manipulation
+"""
+
+from util import *
+from itertools import imap
+import random
+
+def _shiftr(v, anz):
+    vor = [frozenset()] * anz
+    return vor + v[:-anz]
+
+def _shiftl(v, anz):
+    vor = [frozenset()] * anz
+    return v[anz:] + vor
+
+def _xor(a, b):
+    return [c.symmetric_difference(d) for c, d in zip(a, b)]
+
+def _andint(v, mask):
+    neu = []
+    for i in xrange(32):
+        if mask & (1 << (31 - i)):
+            neu.append(v[i])
+        else:
+            neu.append(frozenset())
+    return neu
+
+class Twister(object):
+    def __init__(self):
+        self.tw = []
+        self.offset = 0
+        pos = 0
+        for i in xrange(N):
+            self.tw.append([frozenset((i,)) for i in xrange(pos, pos+32)])
+            pos += 32
+
+    def getint32(self):
+        mt = self.tw
+        if self.offset >= N:
+            # for (kk=0;kk<N-M;kk++)
+            for kk in xrange(0, N-M):
+                #y = (mt[kk]&UPPER_MASK)|(mt[kk+1]&LOWER_MASK);
+                #mt[kk] = mt[kk+M] ^ (y >> 1) ^ mag01[y & 0x1UL];
+                y = mt[kk][0:1]+mt[kk+1][1:]
+                t = _andint([y[-1]] * 32, MATRIX_A)
+                mt[kk] = _xor(_xor(mt[kk+M], _shiftr(y, 1)), t)
+
+            #for (;kk<N-1;kk++)
+            for kk in xrange(kk+1, N-1):
+                #y = (mt[kk]&UPPER_MASK)|(mt[kk+1]&LOWER_MASK);
+                #mt[kk] = mt[kk+(M-N)] ^ (y >> 1) ^ mag01[y & 0x1UL];
+                y = mt[kk][0:1]+mt[kk+1][1:]
+                t = _andint([y[-1]] * 32, MATRIX_A)
+                mt[kk] = _xor(_xor(mt[kk+(M-N)], _shiftr(y, 1)), t)
+
+            #y = (mt[N-1]&UPPER_MASK)|(mt[0]&LOWER_MASK);
+            #mt[N-1] = mt[M-1] ^ (y >> 1) ^ mag01[y & 0x1UL];
+            y = mt[N-1][0:1]+mt[0][1:]
+            t = _andint([y[-1]] * 32, MATRIX_A)
+            mt[N-1] = _xor(_xor(mt[M-1], _shiftr(y, 1)), t)
+
+            self.offset = 0
+
+        y = self.tw[self.offset]
+        y = _xor(y, _shiftr(y, 11))
+
+        t = _shiftl(y, 7)
+        t = _andint(t, 0x9d2c5680)
+        y = _xor(y, t)
+
+        t = _shiftl(y, 15)
+        t = _andint(t, 0xefc60000)
+        y = _xor(y, t)
+
+        y = _xor(y, _shiftr(y, 18))
+
+        self.offset += 1        
+        return y
+
+if __name__ == "__main__":
+    print "Verifying implementation..."
+    twister = Twister()
+    v = random.getrandbits(32)
+    state = random._inst.getstate()[1]
+    def getbit(i, state = state):
+        return (state[i//32]>>(31 - i % 32))&1
+    for byte in xrange(10000):
+        print "Byte %d" % (byte)
+        now_state = random._inst.getstate()[1]
+        if byte % 625 == 0:
+            falsch = False
+            for p, bits in enumerate(twister.tw):
+                for bn, bit in enumerate(bits):
+                    #assert sum(imap(getbit, bit)) & 1 == getbit(bn + p * 32, now_state), (p, bn)
+                    if sum(imap(getbit, bit)) & 1 != getbit(bn + p * 32, now_state):
+                        print (p, bn, bit)
+                        falsch = True
+
+            if falsch:
+                raise "Bad Bit"
+        w = twister.getint32()
+        for bit, d in enumerate(w):
+            assert sum(imap(getbit, d)) & 1 == (v >> (31 - bit)) & 1, byte
+        v = random.getrandbits(32)
+    print "OK"

util.py

@@ -0,0 +1,15 @@
+N = 624
+M = 397
+MATRIX_A = 0x9908b0df
+
+def y(y):
+    y ^= (y >> 11)
+    y ^= (y << 7) & 0x9d2c5680
+    y ^= (y << 15) & 0xefc60000
+    y ^= (y >> 18)
+    return y
+
+def random_random(a,b):
+    a=a>>5
+    b=b>>6;
+    return (a*67108864.0+b)*(1.0/9007199254740992.0)